1#ifndef MADNESS_APPS_MOLDFT_DKOPS_H_INCLUDED
2#define MADNESS_APPS_MOLDFT_DKOPS_H_INCLUDED
22double w(
double t,
double eps){
23 double c = 137.0359895;
25 return 1.0/(
c*sqrt(
PI))*exp(-1.0/2.0*t-
c*
c*exp(-t))-(1.0+eps/(
c*
c))*exp((2.0*eps+(eps*eps)/(
c*
c))*exp(-t)-t)*(erfc((
c+eps/
c)*(exp(-t/2.0)))-2.0);
28double Ebark(
double t,
double eps){
30 return pow(1.0/(2.0*
q(t)), 1.5)*
w(t, eps)*exp(-1.0/(4.0*
q(t))*
R*
R);
37 std::vector<double>
cvec;
38 std::vector<double>
tvec;
59 for(
int i = 0; i <
cvec.size(); i++){
113 std::vector<double>
c;
114 std::vector<double> t;
117 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/Pbar_t.csv");
119 if(world.
rank() == 0) std::cerr <<
"Unable to open Pbar_t.csv" << std::endl;
131 inf.open(
"/gpfs/home/jscanderson/DKproject/Pbar_c.csv");
133 if(world.
rank() == 0) std::cerr <<
"Unable to open Pbar_c.csv" << std::endl;
160 for(
int i = 0; i < n; i++){
172 std::vector<double>
c;
173 std::vector<double> t;
176 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/Tbar_t.csv");
178 if(world.
rank() == 0) std::cerr <<
"Unable to open Tbar_t.csv" << std::endl;
190 inf.open(
"/gpfs/home/jscanderson/DKproject/Tbar_c.csv");
192 if(world.
rank() == 0) std::cerr <<
"Unable to open Tbar_c.csv" << std::endl;
219 for(
int i = 0; i < n; i++){
232 std::vector<double>
c, t;
235 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/A_t.csv");
237 if(world.
rank() == 0) std::cerr <<
"Unable to open A_t.csv" << std::endl;
249 inf.open(
"/gpfs/home/jscanderson/DKproject/A_c.csv");
251 if(world.
rank() == 0) std::cerr <<
"Unable to open A_c.csv" << std::endl;
275 for(
int i = 0; i < n; i++){
288 std::vector<double>
c, t;
291 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/A_t.csv");
293 if(world.
rank() == 0) std::cerr <<
"Unable to open A_t.csv" << std::endl;
305 inf.open(
"/gpfs/home/jscanderson/DKproject/A_c.csv");
307 if(world.
rank() == 0) std::cerr <<
"Unable to open A_c.csv" << std::endl;
331 for(
int i = 0; i < n; i++){
344 std::vector<double>
c, t;
347 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/PbarA_t.csv");
349 if(world.
rank() == 0) std::cerr <<
"Unable to open PbarA_t.csv" << std::endl;
361 inf.open(
"/gpfs/home/jscanderson/DKproject/PbarA_c.csv");
363 if(world.
rank() == 0) std::cerr <<
"Unable to open PbarA_c.csv" << std::endl;
387 for(
int i = 0; i < n; i++){
399 std::vector<double>
c, t;
401 std::ifstream
inf(
"/gpfs/home/jscanderson/DKproject/PbarA_t.csv");
403 if(world.
rank() == 0) std::cerr <<
"Unable to open PbarA_t.csv" << std::endl;
415 inf.open(
"/gpfs/home/jscanderson/DKproject/PbarA_c.csv");
417 if(world.
rank() == 0) std::cerr <<
"Unable to open PbarA_c.csv" << std::endl;
439 double fac =
pow(2.0*
pi,1.5);
440 for(
int i = 0; i < n; i++){
446 double hi = width.
normf();
448 std::vector<real_convolution_3d_ptr>
gradG(3);
451 std::vector< ConvolutionND<double,3> > ops(n);
452 for (
int mu=0;
mu<n;
mu++) {
455 double cc = std::pow(t[
mu]/
pi,1.5);
458 for (
int d=0;
d<3;
d++) {
real_convolution_3d Tbar(World &world)
Definition DKops.h:170
real_convolution_3d A2(World &world)
Definition DKops.h:286
real_convolution_3d Ebar(World &world, double eps)
Definition DKops.h:33
double w(double t, double eps)
Definition DKops.h:22
real_convolution_3d Pbar(World &world)
Definition DKops.h:111
real_convolution_3d PbarA(World &world)
Definition DKops.h:342
static const double opthresh
Definition DKops.h:16
double q(double t)
Definition DKops.h:18
std::vector< real_convolution_3d_ptr > gradPbarA(World &world)
Definition DKops.h:397
double Ebark(double t, double eps)
Definition DKops.h:28
Definition test_ar.cc:118
static int get_k()
Returns the default wavelet order.
Definition funcdefaults.h:163
static const Tensor< double > & get_cell_width()
Returns the width of each user cell dimension.
Definition funcdefaults.h:369
A tensor is a multidimensional array.
Definition tensor.h:317
float_scalar_type normf() const
Returns the Frobenius norm of the tensor.
Definition tensor.h:1726
A simple, fixed dimension vector.
Definition vector.h:64
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:320
Defines common mathematical and physical constants.
static const double R
Definition csqrt.cc:46
static double pow(const double *a, const double *b)
Definition lda.h:74
Main include file for MADNESS and defines Function interface.
constexpr double pi
Mathematical constant .
Definition constants.h:48
Namespace for all elements and tools of MADNESS.
Definition DFParameters.h:10
std::shared_ptr< real_convolution_3d > real_convolution_3d_ptr
Definition functypedefs.h:150
SeparatedConvolution< double, 3 > real_convolution_3d
Definition functypedefs.h:136
constexpr Vector< T, sizeof...(Ts)+1 > vec(T t, Ts... ts)
Factory function for creating a madness::Vector.
Definition vector.h:750
static XNonlinearSolver< std::vector< Function< T, NDIM > >, T, vector_function_allocator< T, NDIM > > nonlinear_vector_solver(World &world, const long nvec)
Definition nonlinsol.h:284
static const string dir
Definition corepotential.cc:249
const double cc
Definition navstokes_cosines.cc:107
const double mu
Definition navstokes_cosines.cc:95
static const double d
Definition nonlinschro.cc:121
Implements most functionality of separated operators.
static const double c
Definition relops.cc:10
static const double PI
Definition relops.cc:12
static const long k
Definition rk.cc:44
Definition convolution1d.h:990
void e()
Definition test_sig.cc:75
static const double pi
Definition testcosine.cc:6