MADNESS  version 0.9
gfit.h
Go to the documentation of this file.
1 /*
2  This file is part of MADNESS.
3 
4  Copyright (C) 2007,2010 Oak Ridge National Laboratory
5 
6  This program is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2 of the License, or
9  (at your option) any later version.
10 
11  This program is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with this program; if not, write to the Free Software
18  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20  For more information please contact:
21 
22  Robert J. Harrison
23  Oak Ridge National Laboratory
24  One Bethel Valley Road
25  P.O. Box 2008, MS-6367
26 
27  email: harrisonrj@ornl.gov
28  tel: 865-241-3937
29  fax: 865-572-0680
30 
31 
32  $Id: key.h 2907 2012-06-14 10:15:05Z 3ru6ruWu $
33  */
34 
35 #ifndef MADNESS_MRA_GFIT_H__INCLUDED
36 #define MADNESS_MRA_GFIT_H__INCLUDED
37 
40 
41 //#include <iostream>
42 #include <madness/tensor/tensor.h>
43 #include <madness/constants.h>
44 namespace madness {
45 
46 template<typename T, std::size_t NDIM>
47 class GFit {
48 
49 public:
50 
52  GFit() {}
53 
55  static GFit CoulombFit(double lo, double hi, double eps, bool prnt=false) {
56  GFit fit=BSHFit(0.0,lo,hi,eps/(4.0*constants::pi),prnt);
57  fit.coeffs_.scale(4.0*constants::pi);
58  return fit;
59  }
60 
62 
70  static GFit BSHFit(double mu, double lo, double hi, double eps, bool prnt=false) {
71  GFit fit;
72  if (NDIM==3) bsh_fit(mu,lo,hi,eps,fit.coeffs_,fit.exponents_,prnt);
73  else bsh_fit_ndim(NDIM,mu,lo,hi,eps,fit.coeffs_,fit.exponents_,prnt);
74  return fit;
75  }
76 
78 
86  static GFit SlaterFit(double gamma, double lo, double hi, double eps, bool prnt=false) {
87  GFit fit;
88  slater_fit(gamma,lo,hi,eps,fit.coeffs_,fit.exponents_,prnt);
89  return fit;
90  }
91 
93 
96  static GFit GeneralFit() {
97  MADNESS_EXCEPTION("General GFit still to be implemented",1);
98  return GFit();
99  }
100 
102  Tensor<T> coeffs() const {return coeffs_;}
103 
105  Tensor<T> exponents() const {return exponents_;}
106 
107  void truncate_periodic_expansion(Tensor<double>& c, Tensor<double>& e,
108  double L, bool discardG0) const {
109  double tcut = 0.25/L/L;
110 
111  if (discardG0) {
112  // Relies on expnts being in decreasing order
113  for (int i=0; i<e.dim(0); ++i) {
114  if (e(i) < tcut) {
115  c = c(Slice(0,i));
116  e = e(Slice(0,i));
117  break;
118  }
119  }
120  } else {
121  // Relies on expnts being in decreasing order
122  int icut = -1;
123  for (int i=0; i<e.dim(0); ++i) {
124  if (e(i) < tcut) {
125  icut = i;
126  break;
127  }
128  }
129  if (icut > 0) {
130  for (int i=icut+1; i<e.dim(0); ++i) {
131  c(icut) += c(i);
132  }
133  c = c(Slice(0,icut));
134  e = e(Slice(0,icut));
135  }
136  }
137  }
138 
139 private:
140 
142 
145  template<typename funcT>
146  GFit(const funcT& f) {
147 
148  }
149 
151  Tensor<T> coeffs_;
152 
154  Tensor<T> exponents_;
155 
157 
165  static void bsh_fit(double mu, double lo, double hi, double eps,
166  Tensor<double>& pcoeff, Tensor<double>& pexpnt, bool prnt) {
167 
168  if (mu < 0.0) throw "cannot handle negative mu in bsh_fit";
169 
170  if (mu > 0) {
171  // Restrict hi according to the exponential decay
172  double r = -log(4*constants::pi*0.01*eps);
173  r = -log(r * 4*constants::pi*0.01*eps);
174  if (hi > r) hi = r;
175  }
176 
177  double TT;
178  double slo, shi;
179 
180  if (eps >= 1e-2) TT = 5;
181  else if (eps >= 1e-4) TT = 10;
182  else if (eps >= 1e-6) TT = 14;
183  else if (eps >= 1e-8) TT = 18;
184  else if (eps >= 1e-10) TT = 22;
185  else if (eps >= 1e-12) TT = 26;
186  else TT = 30;
187 
188  if (mu > 0) {
189  slo = -0.5*log(4.0*TT/(mu*mu));
190  }
191  else {
192  slo = log(eps/hi) - 1.0;
193  }
194  shi = 0.5*log(TT/(lo*lo));
195  if (shi <= slo) throw "bsh_fit: logic error in slo,shi";
196 
197  // Resolution required for quadrature over s
198  double h = 1.0/(0.2-.50*log10(eps)); // was 0.5 was 0.47
199 
200  // Truncate the number of binary digits in h's mantissa
201  // so that rounding does not occur when performing
202  // manipulations to determine the quadrature points and
203  // to limit the number of distinct values in case of
204  // multiple precisions being used at the same time.
205  h = floor(64.0*h)/64.0;
206 
207 
208  // Round shi/lo up/down to an integral multiple of quadrature points
209  shi = ceil(shi/h)*h;
210  slo = floor(slo/h)*h;
211 
212  long npt = long((shi-slo)/h+0.5);
213 
214  //if (prnt)
215  //std::cout << "mu " << mu << " slo " << slo << " shi " << shi << " npt " << npt << " h " << h << " " << eps << std::endl;
216 
217  Tensor<double> coeff(npt), expnt(npt);
218 
219  for (int i=0; i<npt; ++i) {
220  double s = slo + h*(npt-i); // i+1
221  coeff[i] = h*2.0/sqrt(constants::pi)*exp(-mu*mu*exp(-2.0*s)/4.0)*exp(s);
222  coeff[i] = coeff[i]/(4.0*constants::pi);
223  expnt[i] = exp(2.0*s);
224  }
225 
226 #if ONE_TERM
227  npt=1;
228  double s=1.0;
229  coeff[0]=1.0;
230  expnt[0] = exp(2.0*s);
231  coeff=coeff(Slice(0,0));
232  expnt=expnt(Slice(0,0));
233  print("only one term in gfit",s,coeff[0],expnt[0]);
234 
235 
236 #endif
237 
238  // Prune large exponents from the fit ... never necessary due to construction
239 
240  // Prune small exponents from Coulomb fit. Evaluate a gaussian at
241  // the range midpoint, and replace it there with the next most
242  // diffuse gaussian. Then examine the resulting error at the two
243  // end points ... if this error is less than the desired
244  // precision, can discard the diffuse gaussian.
245 
246  if (mu == 0.0) {
247  double mid = lo + (hi-lo)*0.5;
248  long i;
249  for (i=npt-1; i>0; --i) {
250  double cnew = coeff[i]*exp(-(expnt[i]-expnt[i-1])*mid*mid);
251  double errlo = coeff[i]*exp(-expnt[i]*lo*lo) -
252  cnew*exp(-expnt[i-1]*lo*lo);
253  double errhi = coeff[i]*exp(-expnt[i]*hi*hi) -
254  cnew*exp(-expnt[i-1]*hi*hi);
255  if (std::max(std::abs(errlo),std::abs(errhi)) > 0.03*eps) break;
256  npt--;
257  coeff[i-1] = coeff[i-1] + cnew;
258  }
259  coeff = coeff(Slice(0,npt-1));
260  expnt = expnt(Slice(0,npt-1));
261  }
262 
263  // Modify the coeffs of the largest exponents to satisfy the moment conditions
264  //
265  // SETTING NMOM>1 TURNS OUT TO BE A BAD IDEA (AS CURRENTLY IMPLEMENTED)
266  // [It is accurate and efficient for a one-shot application but it seems to
267  // introduce fine-scale noise that amplifies during iterative solution of
268  // the SCF and DFT equations ... the symptom is negative coeffs in the fit]
269  //
270  // SET NMOM=0 or 1 (1 recommended) unless you are doing a one-shot application
271  //
272  // Determine the effective range of the four largest exponents and compute
273  // moments of the exact and remainder of the fit. Then adjust the coeffs
274  // to reproduce the exact moments in that volume.
275  //
276  // If nmom!=4 we assume that we will eliminate n=-1 which is stored first
277  // in the moment list
278  //
279  // <r^i|gj> cj = <r^i|exact-remainder>
280  const long nmom = 1;
281  if (nmom > 0) {
282  Tensor<double> q(4), qg(4);
283  double range = sqrt(-log(1e-6)/expnt[nmom-1]);
284  if (prnt) print("exponent(nmom-1)",expnt[nmom-1],"has range", range);
285 
286  bsh_spherical_moments(mu, range, q);
287  Tensor<double> M(nmom,nmom);
288  for (int i=nmom; i<npt; ++i) {
289  Tensor<double> qt(4);
290  gaussian_spherical_moments(expnt[i], range, qt);
291  qg += qt*coeff[i];
292  }
293  if (nmom != 4) {
294  q = q(Slice(1,nmom));
295  qg = qg(Slice(1,nmom));
296  }
297  if (prnt) {
298  print("moments", q);
299  print("moments", qg);
300  }
301  q = q - qg;
302  for (int j=0; j<nmom; ++j) {
303  Tensor<double> qt(4);
304  gaussian_spherical_moments(expnt[j], range, qt);
305  if (nmom != 4) qt = qt(Slice(1,nmom));
306  for (int i=0; i<nmom; ++i) {
307  M(i,j) = qt[i];
308  }
309  }
310  Tensor<double> ncoeff;
311  gesv(M, q, ncoeff);
312  if (prnt) {
313  print("M\n",M);
314  print("old coeffs", coeff(Slice(0,nmom-1)));
315  print("new coeffs", ncoeff);
316  }
317 
318  coeff(Slice(0,nmom-1)) = ncoeff;
319  }
320 
321  if (prnt) {
322  for (int i=0; i<npt; ++i)
323  std::cout << i << " " << coeff[i] << " " << expnt[i] << std::endl;
324 
325  long npt = 300;
326  //double hi = 1.0;
327  //if (mu) hi = min(1.0,30.0/mu);
328  std::cout << " x value abserr relerr" << std::endl;
329  std::cout << " ------------ ------- -------- -------- " << std::endl;
330  double step = exp(log(hi/lo)/(npt+1));
331  for (int i=0; i<=npt; ++i) {
332  double r = lo*(pow(step,i+0.5));
333  double exact = exp(-mu*r)/r/4.0/constants::pi;
334  double test = 0.0;
335  for (int j=0; j<coeff.dim(0); ++j)
336  test += coeff[j]*exp(-r*r*expnt[j]);
337  double err = 0.0;
338  if (exact) err = (exact-test)/exact;
339  printf(" %.6e %8.1e %8.1e %8.1e\n",r, exact, exact-test, err);
340  }
341  }
342  pcoeff = coeff;
343  pexpnt = expnt;
344  }
345 
347 
352  static void slater_fit(double gamma, double lo, double hi, double eps,
353  Tensor<double>& pcoeff, Tensor<double>& pexpnt, bool prnt) {
354 
355  // empirical number TT for the upper integration limit
356  double TT;
357  if (eps >= 1e-2) TT = 5;
358  else if (eps >= 1e-4) TT = 10;
359  else if (eps >= 1e-6) TT = 14;
360  else if (eps >= 1e-8) TT = 18;
361  else if (eps >= 1e-10) TT = 22;
362  else if (eps >= 1e-12) TT = 26;
363  else TT = 30;
364 
365  // integration limits for quadrature over s: slo and shi
366  double slo=0.5 * log(eps) - 1.0;
367  double shi=log(TT/(lo*lo))*0.5;
368 
369  // Resolution required for quadrature over s
370  double h = 1.0/(0.2-.5*log10(eps)); // was 0.5 was 0.47
371 
372  // Truncate the number of binary digits in h's mantissa
373  // so that rounding does not occur when performing
374  // manipulations to determine the quadrature points and
375  // to limit the number of distinct values in case of
376  // multiple precisions being used at the same time.
377  h = floor(64.0*h)/64.0;
378 
379  // Round shi/lo up/down to an integral multiple of quadrature points
380  shi = ceil(shi/h)*h;
381  slo = floor(slo/h)*h;
382 
383  long npt = long((shi-slo)/h+0.5);
384 
385  Tensor<double> coeff(npt), expnt(npt);
386 
387  for (int i=0; i<npt; ++i) {
388  const double s = slo + h*(npt-i); // i+1
389  coeff[i] = h*exp(-gamma*gamma*exp(2.0*s) + s);
390  coeff[i]*=2.0*gamma/sqrt(constants::pi);
391  expnt[i] = 0.25*exp(-2.0*s);
392  }
393 
394  // Prune large exponents from the fit ... never necessary due to construction
395 
396  // Prune small exponents from Coulomb fit. Evaluate a gaussian at
397  // the range midpoint, and replace it there with the next most
398  // diffuse gaussian. Then examine the resulting error at the two
399  // end points ... if this error is less than the desired
400  // precision, can discard the diffuse gaussian.
401 
402  if (gamma == 0.0) {
403  double mid = lo + (hi-lo)*0.5;
404  long i;
405  for (i=npt-1; i>0; --i) {
406  double cnew = coeff[i]*exp(-(expnt[i]-expnt[i-1])*mid*mid);
407  double errlo = coeff[i]*exp(-expnt[i]*lo*lo) -
408  cnew*exp(-expnt[i-1]*lo*lo);
409  double errhi = coeff[i]*exp(-expnt[i]*hi*hi) -
410  cnew*exp(-expnt[i-1]*hi*hi);
411  if (std::max(std::abs(errlo),std::abs(errhi)) > 0.03*eps) break;
412  npt--;
413  coeff[i-1] = coeff[i-1] + cnew;
414  }
415  coeff = coeff(Slice(0,npt-1));
416  expnt = expnt(Slice(0,npt-1));
417  }
418 
419 
420  if (prnt) {
421  std::cout << "weights and roots for a Slater function with gamma=" << gamma << std::endl;
422  for (int i=0; i<npt; ++i)
423  std::cout << i << " " << coeff[i] << " " << expnt[i] << std::endl;
424 
425  long npt = 300;
426  //double hi = 1.0;
427  //if (mu) hi = min(1.0,30.0/mu);
428  std::cout << " x value abserr relerr" << std::endl;
429  std::cout << " ------------ ------- -------- -------- " << std::endl;
430  double step = exp(log(hi/lo)/(npt+1));
431  for (int i=0; i<=npt; ++i) {
432  double r = lo*(pow(step,i+0.5));
433  double exact = exp(-gamma*r);
434  double test = 0.0;
435  for (int j=0; j<coeff.dim(0); ++j)
436  test += coeff[j]*exp(-r*r*expnt[j]);
437  double err = 0.0;
438  if (exact) err = (exact-test)/exact;
439  printf(" %.6e %8.1e %8.1e %8.1e\n",r, exact, exact-test, err);
440  }
441  }
442  pcoeff = coeff;
443  pexpnt = expnt;
444  }
445 
446  void static bsh_fit_ndim(int ndim, double mu, double lo, double hi, double eps,
447  Tensor<double>& pcoeff, Tensor<double>& pexpnt, bool prnt) {
448 
449  if (mu > 0) {
450  // Restrict hi according to the exponential decay
451  double r = -log(4*constants::pi*0.01*eps);
452  r = -log(r * 4*constants::pi*0.01*eps);
453  if (hi > r) hi = r;
454  }
455 
456 
457  // Determine range of quadrature by estimating when
458  // kernel drops to exp(-100)
459 
460  double slo, shi;
461  if (mu > 0) {
462  slo = -0.5*log(4.0*100.0/(mu*mu));
463  slo = -0.5*log(4.0*(slo*ndim - 2.0*slo + 100.0)/(mu*mu));
464  }
465  else {
466  slo = log(eps/hi) - 1.0;
467  }
468  shi = 0.5*log(100.0/(lo*lo));
469 
470  // Resolution required for quadrature over s
471  double h = 1.0/(0.2-.50*log10(eps)); // was 0.5 was 0.47
472 
473  // Truncate the number of binary digits in h's mantissa
474  // so that rounding does not occur when performing
475  // manipulations to determine the quadrature points and
476  // to limit the number of distinct values in case of
477  // multiple precisions being used at the same time.
478  h = floor(64.0*h)/64.0;
479 
480 
481  // Round shi/lo up/down to an integral multiple of quadrature points
482  shi = ceil(shi/h)*h;
483  slo = floor(slo/h)*h;
484 
485  long npt = long((shi-slo)/h+0.5);
486 
487  if (prnt)
488  std::cout << "bsh: mu " << mu << " lo " << lo << " hi " << hi
489  << " eps " << eps << " slo " << slo << " shi " << shi
490  << " npt " << npt << " h " << h << std::endl;
491 
492 
493  // Compute expansion pruning small coeffs and large exponents
494  Tensor<double> coeff(npt), expnt(npt);
495  int nnpt=0;
496  for (int i=0; i<npt; ++i) {
497  double s = slo + h*(npt-i); // i+1
498  double c = exp(-0.25*mu*mu*exp(-2.0*s)+(ndim-2)*s)*0.5/pow(constants::pi,0.5*ndim);
499  double p = exp(2.0*s);
500  c = c*h;
501  if (c*exp(-p*lo*lo) > eps) {
502  coeff(nnpt) = c;
503  expnt(nnpt) = p;
504  ++nnpt;
505  }
506  }
507  npt = nnpt;
508 #if ONE_TERM
509  npt=1;
510  double s=1.0;
511  coeff[0]=1.0;
512  expnt[0] = exp(2.0*s);
513  coeff=coeff(Slice(0,0));
514  expnt=expnt(Slice(0,0));
515  print("only one term in gfit",s,coeff[0],expnt[0]);
516 
517 #endif
518 
519 
520  // Prune small exponents from Coulomb fit. Evaluate a gaussian at
521  // the range midpoint, and replace it there with the next most
522  // diffuse gaussian. Then examine the resulting error at the two
523  // end points ... if this error is less than the desired
524  // precision, can discard the diffuse gaussian.
525 
526  if (mu == 0.0) {
527  double mid = lo + (hi-lo)*0.5;
528  long i;
529  for (i=npt-1; i>0; --i) {
530  double cnew = coeff[i]*exp(-(expnt[i]-expnt[i-1])*mid*mid);
531  double errlo = coeff[i]*exp(-expnt[i]*lo*lo) -
532  cnew*exp(-expnt[i-1]*lo*lo);
533  double errhi = coeff[i]*exp(-expnt[i]*hi*hi) -
534  cnew*exp(-expnt[i-1]*hi*hi);
535  if (std::max(std::abs(errlo),std::abs(errhi)) > 0.03*eps) break;
536  npt--;
537  coeff[i-1] = coeff[i-1] + cnew;
538  }
539  }
540 
541  // Shrink array to correct size
542  coeff = coeff(Slice(0,npt-1));
543  expnt = expnt(Slice(0,npt-1));
544 
545 
546  if (prnt) {
547  for (int i=0; i<npt; ++i)
548  std::cout << i << " " << coeff[i] << " " << expnt[i] << std::endl;
549 
550  long npt = 300;
551  std::cout << " x value" << std::endl;
552  std::cout << " ------------ ---------------------" << std::endl;
553  double step = exp(log(hi/lo)/(npt+1));
554  for (int i=0; i<=npt; ++i) {
555  double r = lo*(pow(step,i+0.5));
556  double test = 0.0;
557  for (int j=0; j<coeff.dim(0); ++j)
558  test += coeff[j]*exp(-r*r*expnt[j]);
559  printf(" %.6e %20.10e\n",r, test);
560  }
561  }
562 
563  pcoeff = coeff;
564  pexpnt = expnt;
565  }
566 
567  // Returns in q[0..4] int(r^2(n+1)*exp(-alpha*r^2),r=0..R) n=-1,0,1,2
568  static void gaussian_spherical_moments(double alpha, double R, Tensor<double>& q) {
569  q[0] = -(-0.1e1 + exp(-alpha * R*R)) / alpha / 0.2e1;
570  q[1] = (-0.2e1 * R * pow(alpha, 0.3e1 / 0.2e1) + sqrt(constants::pi)
571  * erf(R * sqrt(alpha)) * alpha * exp(alpha * R*R))
572  * pow(alpha, -0.5e1 / 0.2e1) * exp(-alpha * R*R) / 0.4e1;
573  q[2] = -(-0.1e1 + exp(-alpha * R*R) + exp(-alpha * R*R) * alpha * R*R)
574  * pow(alpha, -0.2e1) / 0.2e1;
575  q[3] = -(-0.3e1 * sqrt(constants::pi) * erf(R * sqrt(alpha)) * pow(alpha, 0.2e1)
576  * exp(alpha * R*R) + 0.6e1 * R * pow(alpha, 0.5e1 / 0.2e1)
577  + 0.4e1 * pow(R, 0.3e1) * pow(alpha, 0.7e1 / 0.2e1))
578  * pow(alpha, -0.9e1 / 0.2e1) * exp(-alpha * R*R) / 0.8e1;
579  }
580 
581  // Returns in q[0..4] int(r^2(n+1)*exp(-mu*r)/(4*constants::pi*r),r=0..R) n=-1,0,1,2
582  static void bsh_spherical_moments(double mu, double R, Tensor<double>& q) {
583  if (mu == 0.0) {
584  q[0] = R / constants::pi / 0.4e1;
585  q[1] = pow(R, 0.2e1) / constants::pi / 0.8e1;
586  q[2] = pow(R, 0.3e1) / constants::pi / 0.12e2;
587  q[3] = pow(R, 0.4e1) / constants::pi / 0.16e2;
588  }
589  else {
590  q[0] = (exp(mu * R) - 0.1e1) / mu * exp(-mu * R) / constants::pi / 0.4e1;
591  q[1] = -(-exp(mu * R) + 0.1e1 + mu * R) * pow(mu, -0.2e1) / constants::pi
592  * exp(-mu * R) / 0.4e1;
593  q[2] = -(-0.2e1 * exp(mu * R) + 0.2e1 + 0.2e1 * mu * R + R*R *
594  pow(mu, 0.2e1))*pow(mu, -0.3e1) / constants::pi * exp(-mu * R) / 0.4e1;
595  q[3] = -(-0.6e1 * exp(mu * R) + 0.6e1 + 0.6e1 * mu * R + 0.3e1 * R*R
596  * pow(mu, 0.2e1) + pow(R, 0.3e1) * pow(mu, 0.3e1))
597  * pow(mu, -0.4e1) / constants::pi * exp(-mu * R) / 0.4e1;
598  }
599  }
600 
601 };
602 
603 
604 }
605 
606 #endif // MADNESS_MRA_GFIT_H__INCLUDED
607 
madness::gesv
void gesv(const Tensor< T > &a, const Tensor< T > &b, Tensor< T > &x)
Solve Ax = b for general A using the LAPACK *gesv routines.
Definition: lapack.cc:339
mpfr::exp
const mpreal exp(const mpreal &v, mp_rnd_t rnd_mode)
Definition: mpreal.h:2234
madness::constants::pi
const double pi
Mathematical constant pi.
Definition: constants.h:44
mpfr::erf
const mpreal erf(const mpreal &v, mp_rnd_t rnd_mode)
Definition: mpreal.h:2457
mpfr::log
const mpreal log(const mpreal &v, mp_rnd_t rnd_mode)
Definition: mpreal.h:2213
madness::GFit::SlaterFit
static GFit SlaterFit(double gamma, double lo, double hi, double eps, bool prnt=false)
return a fit for the Slater function
Definition: gfit.h:86
e1
FLOAT e1(FLOAT z)
Definition: y1.cc:144
mpfr::pow
const mpreal pow(const mpreal &a, const unsigned int b, mp_rnd_t rnd_mode=mpreal::default_rnd)
Definition: mpreal.h:2788
NDIM
const int NDIM
Definition: tdse1.cc:44
constants.h
Defines common mathematical and physical constants.
madness::f
NDIM & f
Definition: mra.h:2179
madness::GFit::coeffs
Tensor< T > coeffs() const
return the coefficients of the fit
Definition: gfit.h:102
fit
FLOAT fit(const FLOAT &x, const vector< FLOAT > &p)
Definition: y.cc:326
mpfr::ceil
const mpreal ceil(const mpreal &v)
Definition: mpreal.h:2596
tensor.h
Defines and implements most of Tensor.
madness::GFit::exponents
Tensor< T > exponents() const
return the exponents of the fit
Definition: gfit.h:105
max
#define max(a, b)
Definition: lda.h:53
madness::GFit::BSHFit
static GFit BSHFit(double mu, double lo, double hi, double eps, bool prnt=false)
return a fit for the bound-state Helmholtz function
Definition: gfit.h:70
madness::GFit::GeneralFit
static GFit GeneralFit()
return a fit for a general isotropic function
Definition: gfit.h:96
mpfr::log10
const mpreal log10(const mpreal &v, mp_rnd_t rnd_mode)
Definition: mpreal.h:2227
madness::GFit::truncate_periodic_expansion
void truncate_periodic_expansion(Tensor< double > &c, Tensor< double > &e, double L, bool discardG0) const
Definition: gfit.h:107
mu
Namespace for mathematical applications.
Definition: muParser.cpp:47
sqrt
tensorT sqrt(const tensorT &s, double tol=1e-8)
Computes matrix square root (not used any more?)
Definition: DFcode/moldft.cc:446
mpfr::gamma
const mpreal gamma(const mpreal &v, mp_rnd_t rnd_mode)
Definition: mpreal.h:2429
madness::GFit::CoulombFit
static GFit CoulombFit(double lo, double hi, double eps, bool prnt=false)
return a fit for the Coulomb function
Definition: gfit.h:55
test
void test(World &world, bool doloadbal=false)
Definition: dataloadbal.cc:225
madness::print
void print(const A &a)
Print a single item to std::cout terminating with new line.
Definition: print.h:122
madness::GFit
Definition: gfit.h:47
madness::GFit::GFit
GFit()
default ctor does nothing
Definition: gfit.h:52
MADNESS_EXCEPTION
#define MADNESS_EXCEPTION(msg, value)
Definition: worldexc.h:88
madness::Slice
A slice defines a sub-range or patch of a dimension.
Definition: slice.h:103
madness
Holds machinery to set up Functions/FuncImpls using various Factories and Interfaces.
Definition: chem/atomutil.cc:45
c
const double c
Definition: gfit.cc:200
mpfr::floor
const mpreal floor(const mpreal &v)
Definition: mpreal.h:2604
Generated on Fri Sep 26 2014 10:38:51 for MADNESS by   doxygen 1.8.8