5 #ifndef MADNESS_CHEM_GTH_PSEUDOPOTENTIAL_H__INCLUDED
6 #define MADNESS_CHEM_GTH_PSEUDOPOTENTIAL_H__INCLUDED
9 #include <madness/external/tinyxml/tinyxml.h>
55 const double x = r[0]-
center[0];
const double y = r[1]-
center[1];
const double z = r[2]-
center[2];
56 double rr = std::sqrt(x*x + y*y + z*z);
58 double rs2 = rs*rs;
double rs4 = rs2*rs2;
double rs6 = rs2*rs4;
59 return -(
Zeff/rr)*erf(rs/std::sqrt(2.0))
101 double rsq = x*x + y*y + z*z;
103 if (rsq > 40.0)
return 0.0;
105 double rr = std::sqrt(rsq);
110 rval *= std::sqrt(2);
112 else if (
itmp2 == 1) {
113 rval *= rr*std::sqrt(2);
115 else if (
itmp2 == 2) {
116 rval *= rsq*std::sqrt(2);
118 else if (
itmp2 == 3) {
119 rval *= rr*rsq*std::sqrt(2);
121 else if (
itmp2 == 4) {
122 rval *= rsq*rsq*std::sqrt(2);
124 else if (
itmp2 == 5) {
125 rval *= rr*rsq*rsq*std::sqrt(2);
127 else if (
itmp2 == 6) {
128 rval *= rsq*rsq*rsq*std::sqrt(2);
130 else if (
itmp2 == 7) {
131 rval *= rr*rsq*rsq*rsq*std::sqrt(2);
135 rval *= (1./2.)*std::sqrt(1./
PI);
138 rval *= std::sqrt(3./4./
PI)*x;
141 rval *= std::sqrt(3./4./
PI)*y;
144 rval *= std::sqrt(3./4./
PI)*z;
151 rval *= (1./4.)*std::sqrt(5./
PI)*(-x*x - y*y + 2*z*z);
154 rval *= (1./2.)*std::sqrt(15./
PI)*(y*z);
157 rval *= (1./2.)*std::sqrt(15./
PI)*(x*z);
160 rval *= (1./2.)*std::sqrt(15./
PI)*(x*y);
163 rval *= (1./4.)*std::sqrt(15./
PI)*(x*x - y*y);
176 double x1 = c1[0];
double y1 = c1[1];
double z1 = c1[2];
177 double x2 = c2[0];
double y2 = c2[1];
double z2 = c2[2];
191 int ii = -1;
int jj = -1;
int kk = -1;
192 for (
int i = 0; i <= 1; i++) {
193 for (
int j = 0; j <= 1; j++) {
194 for (
int k = 0;
k <= 1;
k++) {
195 double x = (i == 0) ? c1[0] : c2[0];
196 double y = (j == 0) ? c1[1] : c2[1];
197 double z = (
k == 0) ? c1[2] : c2[2];
199 double rsq = rr[0]*rr[0]+rr[1]*rr[1]+rr[2]*rr[2];
203 ii = i; jj = j; kk =
k;
210 if ((ii < 0) || (jj < 0) || (kk < 0))
MADNESS_EXCEPTION(
"GTH_Pseudopotential: failed to find suitable minimum point\n", 0);
211 double x = (ii == 0) ? c1[0] : c2[0];
212 double y = (jj == 0) ? c1[1] : c2[1];
213 double z = (kk == 0) ? c1[2] : c2[2];
215 if (fabs(fval) < ftol)
return true;
222 double* x =
new double[npts];
223 double* y =
new double[npts];
224 double* z =
new double[npts];
225 double* rsq =
new double[npts];
226 double* rr =
new double[npts];
228 double* x1 = xvals[0];
double* x2 = xvals[1];
double* x3 = xvals[2];
229 for (
int i = 0; i < npts; i++) {
233 rsq[i] = x[i]*x[i] + y[i]*y[i] + z[i]*z[i];
234 rr[i] = std::sqrt(rsq[i]);
242 for (
int i = 0; i < npts; i++) {
243 fvals[i] *= std::sqrt(2);
246 else if (
itmp2 == 1) {
247 for (
int i = 0; i < npts; i++) {
248 fvals[i] *= rr[i]*std::sqrt(2);
251 else if (
itmp2 == 2) {
252 for (
int i = 0; i < npts; i++) {
253 fvals[i] *= rsq[i]*std::sqrt(2);
256 else if (
itmp2 == 3) {
257 for (
int i = 0; i < npts; i++) {
258 fvals[i] *= rr[i]*rsq[i]*std::sqrt(2);
261 else if (
itmp2 == 4) {
262 for (
int i = 0; i < npts; i++) {
263 fvals[i] *= rsq[i]*rsq[i]*std::sqrt(2);
266 else if (
itmp2 == 5) {
267 for (
int i = 0; i < npts; i++) {
268 fvals[i] *= rr[i]*rsq[i]*rsq[i]*std::sqrt(2);
271 else if (
itmp2 == 6) {
272 for (
int i = 0; i < npts; i++) {
273 fvals[i] *= rsq[i]*rsq[i]*rsq[i]*std::sqrt(2);
276 else if (
itmp2 == 7) {
277 for (
int i = 0; i < npts; i++) {
278 fvals[i] *= rr[i]*rsq[i]*rsq[i]*rsq[i];
283 for (
int i = 0; i < npts; i++) {
284 fvals[i] *= (1./2.)*std::sqrt(1./
PI)*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
288 for (
int i = 0; i < npts; i++) {
289 fvals[i] *= std::sqrt(3./4./
PI)*x[i]*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
293 for (
int i = 0; i < npts; i++) {
294 fvals[i] *= std::sqrt(3./4./
PI)*y[i]*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
298 for (
int i = 0; i < npts; i++) {
299 fvals[i] *= std::sqrt(3./4./
PI)*z[i]*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
307 for (
int i = 0; i < npts; i++) {
308 fvals[i] *= (1./4.)*std::sqrt(5./
PI)*(-x[i]*x[i] - y[i]*y[i] + 2*z[i]*z[i])*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
312 for (
int i = 0; i < npts; i++) {
313 fvals[i] *= (1./2.)*std::sqrt(15./
PI)*(y[i]*z[i])*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
317 for (
int i = 0; i < npts; i++) {
318 fvals[i] *= (1./2.)*std::sqrt(15./
PI)*(x[i]*z[i])*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
322 for (
int i = 0; i < npts; i++) {
323 fvals[i] *= (1./2.)*std::sqrt(15./
PI)*(x[i]*y[i])*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
327 for (
int i = 0; i < npts; i++) {
328 fvals[i] *= (1./4.)*std::sqrt(15./
PI)*(x[i]*x[i] - y[i]*y[i])*std::exp(-0.5*(rsq[i]/
alpha/
alpha));
371 f1 =
real_factory_3d(world).functor(functor).truncate_on_project().nofence().truncate_mode(0);
385 template <
typename Q>
392 std::array<real_tensor,118>
hlij;
393 std::array<real_tensor,118>
klij;
414 if (
radii[atype-1].dim(0) > 0)
433 VLocalFunctor(atom_localp[0], atom_localp[1], atom_localp[2], atom_localp[3], atom_localp[4], atom_localp[5], center))).
453 if (!doc.LoadFile()) {
460 if (
debug && world.
rank() == 0) {printf(
"atom atomic_number = %d\n", atype);}
462 bool success =
false;
463 for (TiXmlElement* node=doc.FirstChildElement(); node && !success; node=node->NextSiblingElement()) {
464 if (strcmp(node->Value(),
"name") == 0) {
465 std::string
name = node->GetText();
466 if (
debug && world.
rank() == 0) std::cout <<
"Loading pseudopotential file " <<
name << std::endl;
468 else if (strcmp(node->Value(),
"atom") == 0) {
469 const char* symbol = node->Attribute(
"symbol");
473 if (
debug && world.
rank() == 0) std::cout <<
" found atomic pseudopotential " << symbol << std::endl;
475 node->Attribute(
"lmax", &lmax);
476 if (
debug && world.
rank() == 0) std::cout <<
" maximum L is " << lmax << std::endl;
478 real_tensor t_hlij((
long)lmax+1, (
long)3, (
long)3);
479 real_tensor t_klij((
long)lmax+1, (
long)3, (
long)3);
481 TiXmlElement* xmlVLocal = node->FirstChildElement();
483 double zeff = 0.0; xmlVLocal->Attribute(
"Zeff", &zeff); t_localp[0] = zeff;
484 double lradius = 0.0; xmlVLocal->Attribute(
"radius", &lradius); t_localp(1) = lradius;
485 double C1 = 0.0; xmlVLocal->Attribute(
"C1", &C1); t_localp[2] = C1;
486 double C2 = 0.0; xmlVLocal->Attribute(
"C2", &C2); t_localp[3] = C2;
487 double C3 = 0.0; xmlVLocal->Attribute(
"C3", &C3); t_localp[4] = C3;
488 double C4 = 0.0; xmlVLocal->Attribute(
"C4", &C4); t_localp[5] = C4;
490 for (TiXmlElement* xmlLnlproj = xmlVLocal->NextSiblingElement();
491 xmlLnlproj; xmlLnlproj=xmlLnlproj->NextSiblingElement()) {
492 int lvalue = -1; xmlLnlproj->Attribute(
"l", &lvalue);
493 double radius = 0.0; xmlLnlproj->Attribute(
"radius", &radius); t_radii[lvalue] = radius;
494 double h00 = 0.0; xmlLnlproj->Attribute(
"h00", &h00); t_hlij(lvalue, 0, 0) = h00;
495 double h11 = 0.0; xmlLnlproj->Attribute(
"h11", &h11); t_hlij(lvalue, 1, 1) = h11;
496 double h22 = 0.0; xmlLnlproj->Attribute(
"h22", &h22); t_hlij(lvalue, 2, 2) = h22;
497 double k00 = 0.0; xmlLnlproj->Attribute(
"k00", &k00); t_klij(lvalue, 0, 0) = k00;
498 double k11 = 0.0; xmlLnlproj->Attribute(
"k11", &k11); t_klij(lvalue, 1, 1) = k11;
499 double k22 = 0.0; xmlLnlproj->Attribute(
"k22", &k22); t_klij(lvalue, 2, 2) = k22;
503 t_hlij(0, 0, 1) = -1./2.*std::sqrt(3./5.)*t_hlij(0, 1, 1);
504 t_hlij(0, 1, 0) = t_hlij(0, 0, 1);
505 t_hlij(0, 0, 2) = 1./2.*std::sqrt(5./21.)*t_hlij(0, 2, 2);
506 t_hlij(0, 2, 0) = t_hlij(0, 0, 2);
507 t_hlij(0, 1, 2) = -1./2.*std::sqrt(100./63.)*t_hlij(0, 2, 2);
508 t_hlij(0, 2, 1) = t_hlij(0, 1, 2);
510 t_hlij(1, 0, 1) = -1./2.*std::sqrt(5./7.)*t_hlij(1, 1, 1);
511 t_hlij(1, 1, 0) = t_hlij(1, 0, 1);
512 t_hlij(1, 0, 2) = 1./6.*std::sqrt(35./11.)*t_hlij(1, 2, 2);
513 t_hlij(1, 2, 0) = t_hlij(1, 0, 2);
514 t_hlij(1, 1, 2) = -1./6.*14./std::sqrt(11.)*t_hlij(1, 2, 2);
515 t_hlij(1, 2, 1) = t_hlij(1, 1, 2);
517 t_hlij(2, 0, 1) = -1./2.*std::sqrt(7./9.)*t_hlij(2, 1, 1);
518 t_hlij(2, 1, 0) = t_hlij(2, 0, 1);
519 t_hlij(2, 0, 2) = 1./2.*std::sqrt(63./143.)*t_hlij(2, 2, 2);
520 t_hlij(2, 2, 0) = t_hlij(2, 0, 2);
521 t_hlij(2, 1, 2) = -1./2.*18./std::sqrt(143.)*t_hlij(2, 2, 2);
522 t_hlij(2, 2, 1) = t_hlij(2, 1, 2);
526 localp[atype-1] = t_localp;
527 radii[atype-1] = t_radii;
528 hlij[atype-1] = t_hlij;
529 klij[atype-1] = t_klij;
542 unsigned int norbs =
psi.size();
548 unsigned int maxLL = 0;
553 for (
unsigned int iatom = 0; iatom < natoms; iatom++) {
558 if (atom_radii.
dim(0) > 0)
559 maxLL =
std::max(maxLL,(
unsigned int)atom_radii.
dim(0)-1);
563 Tensor<int> Pilm_lookup((
unsigned int) natoms, (
unsigned long) 3, (
unsigned long) maxLL+1, (
unsigned long) 2*maxLL+1);
564 for (
unsigned int iatom = 0; iatom < natoms; iatom++) {
574 for (
unsigned int j = 1; j <= 3; j++) {
575 for (
unsigned int l = 0; l <= maxLL; l++) {
576 for (
unsigned int m = 0;
m < 2*maxLL+1;
m++) {
577 Pilm_lookup(iatom, j-1, l,
m) = lidx++;
578 if (
m < 2*l+1) localproj.push_back(prlmstore.
nlmproj(world,l,
m,j));
594 Pilm = Pilm.
reshape(natoms, 3, maxLL+1, 2*maxLL+1, norbs);
596 Tensor<Q> Qilm((
unsigned int) natoms, (
unsigned long) 3, (
unsigned long) maxLL+1, (
unsigned long) 2*maxLL+1, (
unsigned int) norbs);
597 for (
unsigned int iorb=0; iorb<
psi.size(); iorb++) {
598 for (
unsigned int iatom = 0; iatom < natoms; iatom++) {
604 int maxL = atom_radii.
dim(0)-1;
605 for (
unsigned int i = 1; i <= 3; i++) {
606 for (
int l = 0; l <= maxL; l++) {
607 for (
int m = 0;
m < 2*l+1;
m++) {
609 for (
unsigned int j = 1; j <= 3; j++) {
610 s += atom_hlij(l,i-1,j-1)*Pilm(iatom, j-1,l,
m,iorb);
612 Qilm(iatom, i-1,l,
m,iorb) = s;
618 Qilm = Qilm.
reshape(natoms*3*(maxLL+1)*(2*maxLL+1),norbs);
621 double trantol = vtol2 / std::min(30.0,
double(localproj.size()));
626 int nocc = occ.
size();
628 for(
int i = 0;i < nocc;++i){
629 enl += occ[i] * nlmat(i, i);
644 gaxpy(world, 1.0, vpsi, 1.0, dpsi);
657 for (
unsigned int iatom = 0; iatom < natoms; iatom++) {
664 unsigned int maxLL = atom_radii.
dim(0)-1;
665 for (
unsigned int l = 0; l <= maxLL; l++) {
666 for (
unsigned int m = 0;
m < 2*l+1;
m++) {
667 for (
unsigned int i = 1; i <= 3; i++) {
669 for (
unsigned int j = 1; j <= 3; j++) {
671 for (
unsigned int iorb = 0; iorb <
psi.size(); iorb++) {
672 double val = atom_hlij(l,i-1,j-1)*(fprojj.inner(
psi[iorb]));
674 vpsi[iorb] += val*fproji;
double potential(const coord_3d &r)
Definition: 3dharmonic.cc:132
Definition: molecule.h:58
unsigned int atomic_number
Atomic number.
Definition: molecule.h:61
madness::Vector< double, 3 > get_coords() const
Definition: molecule.h:99
bool pseudo_atom
Indicates if this atom uses a pseudopotential.
Definition: molecule.h:63
long dim(int i) const
Returns the size of dimension i.
Definition: basetensor.h:147
long size() const
Returns the number of elements in the tensor.
Definition: basetensor.h:138
FunctionDefaults holds default paramaters as static class members.
Definition: funcdefaults.h:204
static const double & get_thresh()
Returns the default threshold.
Definition: funcdefaults.h:279
Abstract base class interface required for functors used as input to Functions.
Definition: function_interface.h:68
const Function< T, NDIM > & compress(bool fence=true) const
Compresses the function, transforming into wavelet basis. Possible non-blocking comm.
Definition: mra.h:721
Function< T, NDIM > & gaxpy(const T &alpha, const Function< Q, NDIM > &other, const R &beta, bool fence=true)
Inplace, general bi-linear operation in wavelet basis. No communication except for optional fence.
Definition: mra.h:981
Definition: gth_pseudopotential.h:386
std::vector< Function< Q, 3 > > apply_potential(World &world, const real_function_3d &potential, const std::vector< Function< Q, 3 > > &psi, const tensorT &occ, Q &enl)
Definition: gth_pseudopotential.h:537
void make_pseudo_potential(World &world)
Definition: gth_pseudopotential.h:401
std::array< real_tensor, 118 > klij
Definition: gth_pseudopotential.h:393
void reproject(int k, double thresh)
Definition: gth_pseudopotential.h:445
real_function_3d vlocalp
Definition: gth_pseudopotential.h:394
std::vector< unsigned int > atoms_with_projectors
Definition: gth_pseudopotential.h:395
GTHPseudopotential(World &world, Molecule molecule)
Definition: gth_pseudopotential.h:399
std::array< real_tensor, 118 > hlij
Definition: gth_pseudopotential.h:392
std::array< real_tensor, 118 > localp
Definition: gth_pseudopotential.h:390
void load_pseudo_from_file(World &world, const std::string filename)
Definition: gth_pseudopotential.h:449
std::array< real_tensor, 118 > radii
Definition: gth_pseudopotential.h:391
real_function_3d vlocalpot()
Definition: gth_pseudopotential.h:441
Molecule molecule
Definition: gth_pseudopotential.h:389
std::vector< Function< Q, 3 > > apply_potential_simple(World &world, const real_function_3d &potential, const std::vector< Function< Q, 3 > > &psi, const tensorT &occ, Q &enl)
Definition: gth_pseudopotential.h:649
Definition: molecule.h:124
const Atom & get_atom(unsigned int i) const
Definition: molecule.cc:447
size_t natom() const
Definition: molecule.h:387
Definition: gth_pseudopotential.h:71
int itmp
Definition: gth_pseudopotential.h:81
static const double gamma_data[17]
Definition: gth_pseudopotential.h:88
virtual bool supports_vectorized() const
Does the interface support a vectorized operator()?
Definition: gth_pseudopotential.h:86
double t1
Definition: gth_pseudopotential.h:82
virtual bool screened(const coord_3d &c1, const coord_3d &c2) const
Definition: gth_pseudopotential.h:173
double alpha
Definition: gth_pseudopotential.h:73
std::vector< coord_3d > specialpts
Definition: gth_pseudopotential.h:76
ProjRLMFunctor(double alpha, int l, int m, int i, const coord_3d ¢er)
Definition: gth_pseudopotential.h:90
int l
Definition: gth_pseudopotential.h:74
double operator()(const coord_3d &r) const
Definition: gth_pseudopotential.h:99
double sqrtPI
Definition: gth_pseudopotential.h:80
virtual void operator()(const Vector< double *, 3 > &xvals, double *MADNESS_RESTRICT fvals, int npts) const
Definition: gth_pseudopotential.h:220
std::vector< coord_3d > special_points() const
Override this to return list of special points to be refined more deeply.
Definition: gth_pseudopotential.h:346
Level special_level()
Override this change level refinement for special points (default is 6)
Definition: gth_pseudopotential.h:348
int m
Definition: gth_pseudopotential.h:74
int itmp2
Definition: gth_pseudopotential.h:81
coord_3d center
Definition: gth_pseudopotential.h:75
Definition: gth_pseudopotential.h:353
coord_3d center
Definition: gth_pseudopotential.h:357
real_function_3d nlmproj(World &world, int l, int m, int i)
Definition: gth_pseudopotential.h:363
ProjRLMStore(const real_tensor &radii, const coord_3d ¢er)
Definition: gth_pseudopotential.h:360
real_tensor radii
Definition: gth_pseudopotential.h:356
ProjRLMFunctor nlmproj_functor(World &world, int l, int m, int i)
Definition: gth_pseudopotential.h:380
int maxL
Definition: gth_pseudopotential.h:355
Tensor< T > reshape(int ndimnew, const long *d)
Returns new view/tensor reshaping size/number of dimensions to conforming tensor.
Definition: tensor.h:1384
Definition: gth_pseudopotential.h:40
double operator()(const coord_3d &r) const
Definition: gth_pseudopotential.h:54
double C3
Definition: gth_pseudopotential.h:42
double C1
Definition: gth_pseudopotential.h:42
Level special_level()
Override this change level refinement for special points (default is 6)
Definition: gth_pseudopotential.h:66
double C4
Definition: gth_pseudopotential.h:42
double Zeff
Definition: gth_pseudopotential.h:42
double zi
Definition: gth_pseudopotential.h:42
double C2
Definition: gth_pseudopotential.h:42
coord_3d center
Definition: gth_pseudopotential.h:43
std::vector< coord_3d > special_points() const
Override this to return list of special points to be refined more deeply.
Definition: gth_pseudopotential.h:64
std::vector< coord_3d > specialpts
Definition: gth_pseudopotential.h:44
VLocalFunctor(double Zeff, double zi, double C1, double C2, double C3, double C4, const coord_3d ¢er)
Definition: gth_pseudopotential.h:47
void fence(bool debug=false)
Synchronizes all processes in communicator AND globally ensures no pending AM or tasks.
Definition: worldgop.cc:161
A parallel world class.
Definition: world.h:132
ProcessID rank() const
Returns the process rank in this World (same as MPI_Comm_rank()).
Definition: world.h:318
WorldGopInterface & gop
Global operations.
Definition: world.h:205
double(* f1)(const coord_3d &)
Definition: derivatives.cc:55
static bool debug
Definition: dirac-hatom.cc:16
const double m
Definition: gfit.cc:199
double psi(const Vector< double, 3 > &r)
Definition: hatom_energy.cc:78
static double pow(const double *a, const double *b)
Definition: lda.h:74
#define max(a, b)
Definition: lda.h:51
#define MADNESS_EXCEPTION(msg, value)
Macro for throwing a MADNESS exception.
Definition: madness_exception.h:119
Main include file for MADNESS and defines Function interface.
const double pi
Mathematical constant .
Definition: constants.h:48
File holds all helper structures necessary for the CC_Operator and CC2 class.
Definition: DFParameters.h:10
static const char * filename
Definition: legendre.cc:96
void truncate(World &world, response_space &v, double tol, bool fence)
Definition: basic_operators.cc:30
Function< T, NDIM > project(const Function< T, NDIM > &other, int k=FunctionDefaults< NDIM >::get_k(), double thresh=FunctionDefaults< NDIM >::get_thresh(), bool fence=true)
Definition: mra.h:2399
Tensor< double > tensorT
Definition: distpm.cc:21
void compress(World &world, const std::vector< Function< T, NDIM > > &v, bool fence=true)
Compress a vector of functions.
Definition: vmra.h:133
std::shared_ptr< FunctionFunctorInterface< double, 3 > > real_functor_3d
Definition: functypedefs.h:107
std::vector< real_function_3d > vector_real_function_3d
Definition: functypedefs.h:79
unsigned int symbol_to_atomic_number(const std::string &symbol)
Definition: atomutil.cc:173
int Level
Definition: key.h:55
FunctionFactory< double, 3 > real_factory_3d
Definition: functypedefs.h:93
Function< TENSOR_RESULT_TYPE(L, R), NDIM > mul_sparse(const Function< L, NDIM > &left, const Function< R, NDIM > &right, double tol, bool fence=true)
Sparse multiplication — left and right must be reconstructed and if tol!=0 have tree of norms already...
Definition: mra.h:1742
double get_charge_from_file(const std::string filename, unsigned int atype)
Definition: gth_pseudopotential.cc:8
Vector< double, 3 > coord_3d
Definition: funcplot.h:1042
std::string name(const FuncType &type, const int ex=-1)
Definition: ccpairfunction.h:28
std::vector< Function< TENSOR_RESULT_TYPE(T, R), NDIM > > transform(World &world, const std::vector< Function< T, NDIM > > &v, const Tensor< R > &c, bool fence=true)
Transforms a vector of functions according to new[i] = sum[j] old[j]*c[j,i].
Definition: vmra.h:689
void matrix_inner(DistributedMatrix< T > &A, const std::vector< Function< T, NDIM > > &f, const std::vector< Function< T, NDIM > > &g, bool sym=false)
Definition: distpm.cc:46
void gaxpy(const double a, ScalarResult< T > &left, const double b, const T &right, const bool fence=true)
the result type of a macrotask must implement gaxpy
Definition: macrotaskq.h:140
static long abs(long a)
Definition: tensor.h:218
double Q(double a)
Definition: relops.cc:20
static const double PI
Definition: relops.cc:12
static const double thresh
Definition: rk.cc:45
static const long k
Definition: rk.cc:44
void e()
Definition: test_sig.cc:75
static const int truncate_mode
Definition: testcosine.cc:14