32 #ifndef MADNESS_MRA_MRA_H__INCLUDED
33 #define MADNESS_MRA_MRA_H__INCLUDED
48 #define FUNCTION_INSTANTIATE_1
49 #define FUNCTION_INSTANTIATE_2
50 #define FUNCTION_INSTANTIATE_3
51 #if !defined(HAVE_IBMBGP) || !defined(HAVE_IBMBGQ)
52 #define FUNCTION_INSTANTIATE_4
53 #define FUNCTION_INSTANTIATE_5
54 #define FUNCTION_INSTANTIATE_6
61 void startup(World& world,
int argc,
char** argv,
bool doprint=
false,
bool make_stdcout_nice_to_reals =
true);
79 template<
typename T, std::
size_t NDIM>
82 template<
typename T, std::
size_t NDIM>
85 template<
typename T, std::
size_t NDIM>
88 template<
typename T, std::
size_t NDIM>
89 class FunctionFactory;
91 template<
typename T, std::
size_t NDIM>
92 class FunctionFunctorInterface;
94 template<
typename T, std::
size_t NDIM>
97 template<
typename T, std::
size_t NDIM>
100 template<
typename T, std::
size_t NDIM>
103 template<
typename T, std::
size_t NDIM, std::
size_t LDIM,
typename opT>
106 template<
typename T, std::
size_t NDIM,
typename opT>
109 template<
typename T, std::
size_t NDIM>
121 template <
typename T, std::
size_t NDIM>
129 std::shared_ptr< FunctionImpl<T,NDIM> >
impl;
133 if(
impl==NULL)
return false;
173 if (
this != &
f)
impl =
f.impl;
188 const double eps=1
e-15;
195 for (std::size_t
d=0;
d<
NDIM; ++
d) {
196 if (xsim[
d] < -eps) {
199 else if (xsim[
d] < eps) {
203 if (xsim[
d] > 1.0+eps) {
206 else if (xsim[
d] > 1.0-eps) {
221 const double eps=1
e-15;
228 for (std::size_t
d=0;
d<
NDIM; ++
d) {
229 if (xsim[
d] < -eps) {
232 else if (xsim[
d] < eps) {
236 if (xsim[
d] > 1.0+eps) {
239 else if (xsim[
d] > 1.0-eps) {
243 return impl->eval_local_only(xsim,maxlevel);
254 const double eps=1
e-15;
261 for (std::size_t
d=0;
d<
NDIM; ++
d) {
262 if (xsim[
d] < -eps) {
265 else if (xsim[
d] < eps) {
269 if (xsim[
d] > 1.0+eps) {
272 else if (xsim[
d] > 1.0-eps) {
291 const double eps=1
e-15;
298 for (std::size_t
d=0;
d<
NDIM; ++
d) {
299 if (xsim[
d] < -eps) {
302 else if (xsim[
d] < eps) {
306 if (xsim[
d] > 1.0+eps) {
309 else if (xsim[
d] > 1.0-eps) {
331 const std::vector<long>& npt,
332 bool eval_refine =
false)
const {
335 const double eps=1
e-14;
339 for (std::size_t
d=0;
d<
NDIM; ++
d) {
340 simlo[
d] = cell(
d,0);
341 simhi[
d] = cell(
d,1);
348 for (std::size_t
d=0;
d<
NDIM; ++
d) {
353 double delta = eps*(simhi[
d]-simlo[
d]);
357 return impl->eval_plot_cube(simlo, simhi, npt, eval_refine);
376 if (
impl->world.rank() == 0) result =
eval(xuser).get();
377 impl->world.gop.broadcast(result);
385 T operator()(
double x,
double y=0,
double z=0,
double xx=0,
double yy=0,
double zz=0)
const {
388 if (
NDIM>=2) r[1] = y;
389 if (
NDIM>=3) r[2] = z;
390 if (
NDIM>=4) r[3] = xx;
391 if (
NDIM>=5) r[4] = yy;
392 if (
NDIM>=6) r[5] = zz;
412 impl->world.gop.broadcast(result);
424 template <
typename funcT>
439 template <
typename funcT>
448 impl->world.gop.fence();
465 return impl->is_compressed();
476 return impl->is_reconstructed();
486 return impl ?
impl->is_nonstandard() :
false;
494 return impl->tree_size();
500 print(
"function",
name,
"not assigned yet");
510 return impl->max_depth();
518 return impl->max_local_depth();
526 return impl->max_nodes();
533 return impl->min_nodes();
550 if (!
impl)
return true;
551 return impl->get_autorefine();
561 impl->set_autorefine(value);
569 if (!
impl)
return 0.0;
570 return impl->get_thresh();
580 impl->set_thresh(value);
589 return impl->get_k();
604 if (!
impl)
return *
this;
614 const std::shared_ptr< FunctionImpl<T,NDIM> >&
get_impl()
const {
632 this->impl->set_functor(functor);
633 print(
"set functor in mra.h");
641 template <
typename R>
643 impl = std::shared_ptr<implT>(
new implT(*
f.get_impl(),
f.get_pmap(), zero));
656 const std::shared_ptr< WorldDCPmapInterface< Key<NDIM> > >&
get_pmap()
const {
659 return impl->get_pmap();
671 impl->distribute(newmap);
681 return impl->norm2sq_local();
695 impl->world.gop.fence();
787 if (not
impl)
return *
this;
789 if (finalstate == current_state)
return *
this;
792 impl->change_tree_state(finalstate,
fence);
810 template <
typename opT>
821 return impl->autorefine_square_test(key, t);
824 template <
typename Archive>
void serialize (Archive& ar) {}
837 bool fence =
true)
const {
866 if (
impl)
impl->print_tree_json(os);
872 os <<
"digraph G {" << std::endl;
873 if (
impl)
impl->print_tree_graphviz(os);
874 os <<
"}" << std::endl;
891 template <
typename Archive>
void serialize(Archive& ar) {}
903 template <
typename opT>
913 template <
typename opT>
923 template <
typename opT>
952 template <
typename Q>
980 template <
typename Q,
typename R>
987 "gaxpy requires both functions to be in the same tree state");
991 if (not same_world) {
1007 template <
typename Q>
1020 return gaxpy(
T(1.0), other,
Q(1.0),
true);
1027 template <
typename Q>
1040 return gaxpy(
T(1.0), other,
Q(-1.0),
true);
1047 template <
typename Q>
1092 if (!
impl)
return 0.0;
1094 return impl->trace_local();
1101 if (!
impl)
return 0.0;
1104 impl->world.gop.fence();
1110 template <
typename R>
1117 return impl->inner_local(*(
g.get_impl()));
1124 template<
typename R>
1130 impl->get_coeffs().clear();
1142 template<
typename opT>
1146 impl->get_coeffs().clear();
1157 impl->get_coeffs().clear();
1165 template<
typename opT>
1169 impl->get_coeffs().clear();
1182 impl->get_coeffs().clear();
1193 template<
typename opT>
1197 impl->get_coeffs().clear();
1210 impl->get_coeffs().clear();
1220 template<
size_t LDIM,
size_t KDIM,
typename opT>
1228 impl->finalize_sum();
1234 template<
size_t LDIM,
size_t KDIM>
1241 impl->finalize_sum();
1250 template <
typename R>
1256 if (not
g.is_initialized())
return 0.0;
1268 if (
g.is_on_demand())
return this->inner_on_demand(
g);
1279 impl->world.gop.fence();
1331 impl->world.gop.fence();
1344 const bool leaf_refine=
true)
const {
1347 T local =
impl->inner_adaptive_local(
f, leaf_refine);
1349 impl->world.gop.fence();
1357 template <
typename L>
1372 template<
typename R>
1380 func->replicate_low_dim_functions(
true);
1397 template <
typename R,
size_t LDIM>
1403 static const size_t KDIM=
NDIM-
LDIM;
1412 this->
get_impl()->project_out(result.
get_impl().get(),
g.get_impl().get(),dim,
true);
1416 result.
get_impl()->trickle_down(
false);
1438 template <
typename Archive>
1442 long magic = 0l,
id = 0l, ndim = 0l,
k = 0l;
1443 ar & magic &
id & ndim &
k;
1459 template <
typename Archive>
1473 if (not
impl)
return;
1479 template <
typename Q,
typename opT>
1492 template <
typename Q, std::
size_t D>
1495 std::vector< std::shared_ptr< FunctionImpl<Q,D> > > r(
v.size());
1496 for (
unsigned int i=0; i<
v.size(); ++i) r[i] =
v[i].
get_impl();
1501 template <
typename opT>
1503 std::vector<implT*>
v(vf.size(),NULL);
1504 for (
unsigned int i=0; i<
v.size(); ++i) {
1520 template <
typename opT>
1524 const bool fence=
true) {
1525 std::vector<implT*> vimplin(vin.size(),NULL);
1526 for (
unsigned int i=0; i<vin.size(); ++i) {
1527 if (vin[i].
is_initialized()) vimplin[i] = vin[i].get_impl().get();
1529 std::vector<implT*> vimplout(vout.size(),NULL);
1530 for (
unsigned int i=0; i<vout.size(); ++i) {
1531 if (vout[i].
is_initialized()) vimplout[i] = vout[i].get_impl().get();
1534 impl->multi_to_multi_op_values(
op, vimplin, vimplout,
fence);
1541 template <
typename L,
typename R>
1549 std::vector<FunctionImpl<T,NDIM>*> vresult(right.size());
1550 std::vector<const FunctionImpl<R,NDIM>*> vright(right.size());
1551 for (
unsigned int i=0; i<right.size(); ++i) {
1552 result[i].set_impl(left,
false);
1553 vresult[i] = result[i].impl.get();
1554 vright[i] = right[i].get_impl().get();
1558 vresult[0]->mulXXvec(left.
get_impl().get(), vright, vresult, tol,
fence);
1564 template<
typename L,
typename R>
1574 impl->multiply(leaf_op1,gimpl,fimpl,
fence);
1577 impl->multiply(leaf_op1,fimpl,gimpl,
fence);
1582 template <
typename R,
typename Q>
1593 template <
typename L,
typename R>
1613 impl.reset(
new implT(*
f.impl,
f.get_pmap(),
false));
1627 for (std::size_t i=0; i<
NDIM; ++i)
MADNESS_ASSERT((mirrormap[i]==1) or (mirrormap[i]==-1));
1628 impl.reset(
new implT(*
f.impl,
f.get_pmap(),
false));
1643 const std::vector<long>& map,
const std::vector<long>&
mirror,
1649 for (std::size_t i=0; i<map.size(); ++i)
MADNESS_ASSERT(map[i]>=0 &&
static_cast<std::size_t
>(map[i])<
NDIM);
1651 impl.reset(
new implT(*
f.impl,
f.get_pmap(),
false));
1662 double local =
impl->check_symmetry_local();
1664 impl->world.gop.fence();
1665 double asy=sqrt(
local);
1666 if (this->
world().rank()==0)
print(
"asymmetry wrt particle",asy);
1683 impl->chop_at_level(n,
true);
1690 template <
typename T,
typename opT, std::
size_t NDIM>
1699 template <
typename Q,
typename T, std::
size_t NDIM>
1706 result.set_impl(
f,
false);
1707 result.get_impl()->scale_oop(
alpha,*
f.get_impl(),fence);
1713 template <
typename Q,
typename T, std::
size_t NDIM>
1724 template <
typename Q,
typename T, std::
size_t NDIM>
1733 template <
typename Q,
typename T, std::
size_t NDIM>
1740 template <
typename L,
typename R,std::
size_t NDIM>
1751 result.set_impl(left,
false);
1752 result.get_impl()->mulXX(left.
get_impl().get(), right.
get_impl().get(), tol, fence);
1757 template <
typename L,
typename R,std::
size_t NDIM>
1764 template <
typename L,
typename R,
typename opT, std::
size_t NDIM>
1772 result.set_impl(left,
false);
1773 result.get_impl()->binaryXX(left.
get_impl().get(), right.
get_impl().get(),
op, fence);
1778 template <
typename Q,
typename opT, std::
size_t NDIM>
1779 Function<typename opT::resultT, NDIM>
1781 if (!
func.is_reconstructed())
func.reconstruct();
1785 result.
get_impl()->unaryXXvalues(
func.get_impl().get(),
op, fence);
1791 template <
typename Q,
typename opT, std::
size_t NDIM>
1792 Function<typename opT::resultT, NDIM>
1794 if (!
func.is_reconstructed())
func.reconstruct();
1805 template <
typename L,
typename R, std::
size_t D>
1810 vresult[0].vmulXX(left, vright, vresult, tol, fence);
1818 template <
typename L,
typename R, std::
size_t NDIM>
1824 return mul(left,right,
true);
1830 template<
typename T, std::
size_t KDIM, std::
size_t LDIM>
1831 Function<T,KDIM+LDIM>
1846 std::vector<std::shared_ptr<FunctionImpl<T,KDIM>>> vleft=
get_impl(left);
1847 std::vector<std::shared_ptr<FunctionImpl<T,LDIM>>> vright=
get_impl(right);
1856 template<
typename T, std::
size_t KDIM, std::
size_t LDIM>
1857 Function<T,KDIM+LDIM>
1859 typedef std::vector<Function<T,KDIM>> vector;
1864 template<
typename T, std::
size_t KDIM, std::
size_t LDIM,
typename opT>
1865 Function<T,KDIM+LDIM>
1880 print(
"incomplete FunctionFactory in Function::hartree_product");
1889 std::vector<std::shared_ptr<FunctionImpl<T,KDIM>>> vleft;
1890 std::vector<std::shared_ptr<FunctionImpl<T,LDIM>>> vright;
1892 vright.push_back(right.
get_impl());
1896 if (not same) right.
standard(
false);
1905 template <
typename L,
typename R,std::
size_t NDIM>
1915 template <
typename L,
typename R,std::
size_t NDIM>
1921 return result.gaxpy_oop(
alpha, left,
beta, right, fence);
1925 template <
typename L,
typename R,std::
size_t NDIM>
1934 template<
typename T, std::
size_t NDIM>
1950 template <
typename L,
typename R, std::
size_t NDIM>
1964 return add(left,right,
true);
1969 template <
typename L,
typename R,std::
size_t NDIM>
1980 template <
typename L,
typename R, std::
size_t NDIM>
1992 return sub(left,right,
true);
2001 template <
typename T, std::
size_t NDIM>
2004 bool fence =
true) {
2009 result.
set_impl(std::shared_ptr<implT>(
new implT(*
f.get_impl(), pmap,
false)));
2010 result.
get_impl()->copy_coeffs(*
f.get_impl(), fence);
2016 template <
typename T, std::
size_t NDIM>
2019 return copy(
f,
f.get_pmap(), fence);
2031 template <
typename T,
typename Q, std::
size_t NDIM>
2037 result.
get_impl()->copy_coeffs(*
f.get_impl(), fence);
2045 template <
typename T, std::
size_t NDIM>
2049 return result.
conj(fence);
2062 template <
typename opT,
typename T, std::
size_t LDIM>
2066 World& world=
f1.front().world();
2069 typedef std::vector<Function<T,LDIM>> vecfuncL;
2071 vecfuncL& ff1 =
const_cast< vecfuncL&
>(
f1);
2072 vecfuncL& ff2 =
const_cast< vecfuncL&
>(
f2);
2074 bool same=(ff1[0].get_impl()==ff2[0].get_impl());
2081 for (
auto&
f :
f1)
f.make_nonstandard(
true,
false);
2082 for (
auto&
f :
f2)
f.make_nonstandard(
true,
false);
2094 for (
size_t i=0; i<
f1.size(); ++i)
2098 if (
op.print_timings) {
2103 result.
get_impl()->finalize_apply();
2105 if (
op.modified()) {
2106 result.
get_impl()->trickle_down(
true);
2108 result.
get_impl()->reconstruct(
true);
2111 if (not same)
standard(world,ff2,
false);
2118 template <
typename opT,
typename R, std::
size_t NDIM>
2123 constexpr std::size_t OPDIM=opT::opdim;
2124 constexpr
bool low_dim=(OPDIM*2==
NDIM);
2127 if (
NDIM <= 3 and (not low_dim)) {
2128 result.set_impl(
f,
false);
2129 result.get_impl()->apply(
op, *
f.get_impl(), fence);
2135 result.set_impl(
f,
false);
2136 r1.set_impl(
f,
false);
2138 result.get_impl()->reset_timer();
2141 result.get_impl()->apply_source_driven(
op, *
f.get_impl(), fence);
2164 template <
typename opT,
typename R, std::
size_t NDIM>
2179 if (
op.modified()) {
2184 ff.
get_impl()->undo_redundant(
false);
2185 result.
get_impl()->trickle_down(
true);
2195 fff.
get_impl()->timer_filter.print(
"filter");
2196 fff.
get_impl()->timer_compress_svd.print(
"compress_svd");
2205 double elapsed=result.
get_impl()->finalize_apply();
2206 if (
print_timings) printf(
"time in finalize_apply %8.2f\n",elapsed);
2213 result.
get_impl()->reconstruct(
true);
2216 if (
op.destructive()) {
2229 template <
typename opT,
typename R, std::
size_t NDIM>
2239 result.set_impl(ff,
false);
2240 result.get_impl()->apply_1d_realspace_push(
op, ff.
get_impl().get(),
axis, fence);
2257 template <
typename T, std::
size_t NDIM>
2262 return result.
mapdim(
f,map,fence);
2269 template <
typename T, std::
size_t NDIM>
2274 return result.
mirror(
f,mirrormap,fence);
2286 template <
typename T, std::
size_t NDIM>
2289 const std::vector<long>&
mirror,
bool fence=
true) {
2300 template <
typename T, std::
size_t NDIM>
2301 typename std::enable_if_t<NDIM%2==0, Function<T,NDIM>>
2304 std::vector<long> map(
NDIM);
2305 constexpr std::size_t LDIM=
NDIM/2;
2306 static_assert(LDIM*2==
NDIM);
2307 for (
auto d=0;
d<LDIM; ++
d) {
2326 template <
typename T, std::
size_t NDIM>
2332 std::vector<long> map(
NDIM);
2335 if (symmetry==
"sy_particle") {
2336 map[0]=3; map[1]=4; map[2]=5;
2337 map[3]=0; map[4]=1; map[5]=2;
2338 }
else if (symmetry==
"cx") {
2339 map[0]=0; map[1]=2; map[2]=1;
2340 map[3]=3; map[4]=5; map[5]=4;
2342 }
else if (symmetry==
"cy") {
2343 map[0]=2; map[1]=1; map[2]=0;
2344 map[3]=5; map[4]=4; map[5]=3;
2346 }
else if (symmetry==
"cz") {
2347 map[0]=1; map[1]=0; map[2]=2;
2348 map[3]=4; map[4]=3; map[5]=5;
2351 if (
f.world().rank()==0) {
2352 print(
"unknown parameter in symmetrize:",symmetry);
2358 result.
get_impl()->average(*
f.get_impl());
2371 template<
typename T, std::
size_t NDIM, std::
size_t LDIM>
2374 static_assert(LDIM+LDIM==
NDIM);
2397 template <
typename T, std::
size_t NDIM>
2415 template <
typename T,
typename R, std::
size_t NDIM>
2427 template<std::size_t
NDIM,
typename T, std::size_t LDIM,
typename R, std::size_t KDIM,
2428 std::size_t CDIM = (KDIM + LDIM -
NDIM) / 2>
2431 const std::array<int, CDIM> v2,
int task=0) {
2432 bool prepare = ((
task==0) or (
task==1));
2433 bool work = ((
task==0) or (
task==2));
2434 bool finish = ((
task==0) or (
task==3));
2436 static_assert((KDIM + LDIM -
NDIM) % 2 == 0,
"faulty dimensions in inner (partial version)");
2437 static_assert(KDIM + LDIM - 2 * CDIM ==
NDIM,
"faulty dimensions in inner (partial version)");
2440 for (
int i=0; i<CDIM-1; ++i)
MADNESS_CHECK((v1[i]+1)==v1[i+1]);
2443 for (
int i=0; i<CDIM-1; ++i)
MADNESS_CHECK((v2[i]+1)==v2[i+1]);
2458 f.get_impl()->compute_snorm_and_dnorm(
false);
2459 for (
auto&
g : vg)
g.get_impl()->compute_snorm_and_dnorm(
false);
2464 std::vector<Function<resultT,NDIM>> result(vg.size());
2467 for (
int i=0; i<vg.size(); ++i) {
2470 result[i].get_impl()->partial_inner(*
f.get_impl(),*(vg[i]).get_impl(),v1,v2);
2484 auto erase_list = [] (
const auto& funcimpl) {
2485 typedef typename std::decay_t<decltype(funcimpl)>::keyT keyTT;
2486 std::list<keyTT> to_be_erased;
2487 for (
auto it=funcimpl.get_coeffs().begin(); it!=funcimpl.get_coeffs().end(); ++it) {
2488 const auto& key=it->first;
2489 const auto& node=it->second;
2490 if (not node.has_children()) to_be_erased.push_back(key);
2492 return to_be_erased;
2497 for (
auto&
g : vg) {
2517 template<std::size_t
NDIM,
typename T, std::size_t LDIM,
typename R, std::size_t KDIM,
2518 std::size_t CDIM = (KDIM + LDIM -
NDIM) / 2>
2521 const std::array<int, CDIM> v2,
int task=0) {
2529 template <
typename T, std::
size_t LDIM,
typename R, std::
size_t KDIM>
2533 std::array<int,1>({std::get<0>(v1)}),
2534 std::array<int,1>({std::get<0>(v2)}));
2541 template <
typename T, std::
size_t LDIM,
typename R, std::
size_t KDIM>
2545 std::array<int,2>({std::get<0>(v1),std::get<1>(v1)}),
2546 std::array<int,2>({std::get<0>(v2),std::get<1>(v2)}));
2553 template <
typename T, std::
size_t LDIM,
typename R, std::
size_t KDIM>
2557 std::array<int,3>({std::get<0>(v1),std::get<1>(v1),std::get<2>(v1)}),
2558 std::array<int,3>({std::get<0>(v2),std::get<1>(v2),std::get<2>(v2)}));
2571 template <
typename T,
typename opT, std::
size_t NDIM>
2574 std::shared_ptr< FunctionFunctorInterface<double,3> >
func(
new opT(
g));
2586 template <
typename T,
typename opT, std::
size_t NDIM>
2591 template <
typename T,
typename R, std::
size_t NDIM>
2594 return (
f*
R(1.0)).add_scalar(r);
2597 template <
typename T,
typename R, std::
size_t NDIM>
2600 return (
f*
R(1.0)).add_scalar(r);
2603 template <
typename T,
typename R, std::
size_t NDIM>
2606 return (
f*
R(1.0)).add_scalar(-r);
2609 template <
typename T,
typename R, std::
size_t NDIM>
2612 return (
f*
R(-1.0)).add_scalar(r);
2616 template <std::
size_t NDIM>
2626 template <std::
size_t NDIM>
2636 template <std::
size_t NDIM>
2647 template <std::
size_t NDIM>
2661 template <std::
size_t NDIM>
2667 template <std::
size_t NDIM>
2673 template <std::
size_t NDIM>
2680 template <
typename T, std::
size_t NDIM>
2684 return result.
square(
true);
2688 template <
typename T, std::
size_t NDIM>
2692 return result.
abs(
true);
2696 template <
typename T, std::
size_t NDIM>
2697 typename std::enable_if<!TensorTypeData<T>::iscomplex, Function<T,NDIM> >
::type
2705 template <
typename T, std::
size_t NDIM>
2706 typename std::enable_if<TensorTypeData<T>::iscomplex, Function<typename Tensor<T>::scalar_type,
NDIM> >
::type
2712 template <std::
size_t NDIM>
2718 template <std::
size_t NDIM>
2729 template <
class archiveT,
class T, std::
size_t NDIM>
2736 template <
class archiveT,
class T, std::
size_t NDIM>
2744 template <
class T, std::
size_t NDIM>
2750 template <
class T, std::
size_t NDIM>
2763 template<
typename T, std::
size_t NDIM>
double q(double t)
Definition: DKops.h:18
This header should include pretty much everything needed for the parallel runtime.
long dim(int i) const
Returns the size of dimension i.
Definition: basetensor.h:147
This class is used to specify boundary conditions for all operators.
Definition: funcdefaults.h:101
CompositeFunctorInterface implements a wrapper of holding several functions and functors.
Definition: function_interface.h:165
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
FunctionFactory implements the named-parameter idiom for Function.
Definition: function_factory.h:86
FunctionFactory & empty()
Definition: function_factory.h:239
FunctionFactory & fence(bool fence=true)
Definition: function_factory.h:269
FunctionFactory & nofence()
Definition: function_factory.h:275
virtual FunctionFactory & thresh(double thresh)
Definition: function_factory.h:191
virtual FunctionFactory & k(int k)
Definition: function_factory.h:186
Abstract base class interface required for functors used as input to Functions.
Definition: function_interface.h:68
FunctionImpl holds all Function state to facilitate shallow copy semantics.
Definition: funcimpl.h:941
const dcT & get_coeffs() const
Definition: mraimpl.h:322
void reconstruct(bool fence)
reconstruct this tree – respects fence
Definition: mraimpl.h:1490
bool is_on_demand() const
Definition: mraimpl.h:262
FunctionNode holds the coefficients, etc., at each node of the 2^NDIM-tree.
Definition: funcimpl.h:124
A multiresolution adaptive numerical function.
Definition: mra.h:122
void print_tree_json(std::ostream &os=std::cout) const
Definition: mra.h:864
TENSOR_RESULT_TYPE(T, R) inner(const Function< R
Returns the inner product.
Function< T, NDIM > & gaxpy_oop(T alpha, const Function< L, NDIM > &left, T beta, const Function< R, NDIM > &right, bool fence)
This is replaced with alpha*left + beta*right ... private.
Definition: mra.h:1594
Function< T, NDIM > & fill_cuspy_tree(const bool fence=true)
Special refinement on 6D boxes where the electrons come close (meet)
Definition: mra.h:1179
void gaxpy_ext(const Function< L, NDIM > &left, T(*f)(const coordT &), T alpha, T beta, double tol, bool fence=true) const
Definition: mra.h:1358
const Function< T, NDIM > & refine(bool fence=true) const
Inplace autorefines the function using same test as for squaring.
Definition: mra.h:830
T inner_adaptive(const std::shared_ptr< FunctionFunctorInterface< T, NDIM > > f, const bool leaf_refine=true) const
Definition: mra.h:1343
void unaryop_coeff(const opT &op, bool fence=true)
Unary operation applied inplace to the coefficients.
Definition: mra.h:914
bool is_compressed() const
Returns true if compressed, false otherwise. No communication.
Definition: mra.h:462
Function< T, NDIM > & fill_tree(const opT &op, bool fence=true)
With this being an on-demand function, fill the MRA tree according to different criteria.
Definition: mra.h:1143
Function< T, NDIM > & fill_cuspy_tree(const opT &op, const bool fence=true)
Definition: mra.h:1166
void set_impl(const Function< R, NDIM > &f, bool zero=true)
Replace current FunctionImpl with a new one using the same parameters & map as f.
Definition: mra.h:642
bool autorefine() const
Retunrs.
Definition: mra.h:548
void print_size(const std::string name) const
print some info about this
Definition: mra.h:498
void broaden(const BoundaryConditions< NDIM > &bc=FunctionDefaults< NDIM >::get_bc(), bool fence=true) const
Inplace broadens support in scaling function basis.
Definition: mra.h:836
void print_info() const
Print a summary of the load balancing info.
Definition: mra.h:880
Future< Level > evaldepthpt(const coordT &xuser) const
Definition: mra.h:252
void norm_tree(bool fence=true) const
Initializes information about the function norm at all length scales.
Definition: mra.h:701
void load(World &world, Archive &ar)
Replaces this function with one loaded from an archive using the default processor map.
Definition: mra.h:1439
double norm2sq_local() const
Returns the square of the norm of the local function ... no communication.
Definition: mra.h:678
void sum_down(bool fence=true) const
Sums scaling coeffs down tree restoring state with coeffs only at leaves. Optional fence....
Definition: mra.h:798
Level depthpt(const coordT &xuser) const
Definition: mra.h:406
Function< T, NDIM > & abs(bool fence=true)
Returns *this for chaining.
Definition: mra.h:1068
const std::shared_ptr< WorldDCPmapInterface< Key< NDIM > > > & get_pmap() const
Returns a shared pointer to the process map.
Definition: mra.h:656
void set_autorefine(bool value, bool fence=true)
Sets the value of the autorefine flag. Optional global fence.
Definition: mra.h:558
Function< T, NDIM > & fill_tree(const Function< R, NDIM > &g, bool fence=true)
With this being an on-demand function, fill the MRA tree according to different criteria.
Definition: mra.h:1125
T trace() const
Returns global value of int(f(x),x) ... global comm required.
Definition: mra.h:1099
void print_tree_graphviz(std::ostream &os=std::cout) const
Process 0 prints a graphviz-formatted output of all nodes in the tree (collective)
Definition: mra.h:870
double norm2() const
Returns the 2-norm of the function ... global sum ... works in either basis.
Definition: mra.h:688
void change_tensor_type(const TensorArgs &targs, bool fence=true)
change the tensor type of the coefficients in the FunctionNode
Definition: mra.h:1472
World & world() const
Returns the world.
Definition: mra.h:648
TreeState gstate
Definition: mra.h:1275
TENSOR_RESULT_TYPE(T, R) inner_on_demand(const Function< R
Returns the inner product for one on-demand function.
T operator()(const coordT &xuser) const
Evaluates the function at a point in user coordinates. Collective operation.
Definition: mra.h:371
Function< T, NDIM > & map_and_mirror(const Function< T, NDIM > &f, const std::vector< long > &map, const std::vector< long > &mirror, bool fence)
This is replaced with mirror(map(f)) ... private.
Definition: mra.h:1642
Function< T, NDIM > & operator=(const Function< T, NDIM > &f)
Assignment is shallow. No communication, works in either basis.
Definition: mra.h:171
Function< T, NDIM > & multiop_values(const opT &op, const std::vector< Function< T, NDIM > > &vf)
This is replaced with op(vector of functions) ... private.
Definition: mra.h:1502
MADNESS_ASSERT(g.is_compressed())
static std::vector< std::shared_ptr< FunctionImpl< Q, D > > > vimpl(const std::vector< Function< Q, D > > &v)
Returns vector of FunctionImpl pointers corresponding to vector of functions.
Definition: mra.h:1493
const Function< T, NDIM > & change_tree_state(const TreeState finalstate, bool fence=true) const
changes tree state to given state
Definition: mra.h:785
Function< T, NDIM > & chop_at_level(const int n, const bool fence=true)
remove all nodes with level higher than n
Definition: mra.h:1680
void verify_tree() const
Verifies the tree data structure ... global sync implied.
Definition: mra.h:453
IsSupported< TensorTypeData< Q >, Function< T, NDIM > >::type & operator*=(const Q q)
Inplace scaling by a constant.
Definition: mra.h:1049
std::pair< bool, T > eval_local_only(const Vector< double, NDIM > &xuser, Level maxlevel) const
Evaluate function only if point is local returning (true,value); otherwise return (false,...
Definition: mra.h:220
Function< TENSOR_RESULT_TYPE(T, R), NDIM-LDIM > project_out(const Function< R, LDIM > &g, const int dim) const
project this on the low-dim function g: h(x) = <f(x,y) | g(y)>
Definition: mra.h:1398
int k() const
Returns the number of multiwavelets (k). No communication.
Definition: mra.h:586
double thresh() const
Returns value of truncation threshold. No communication.
Definition: mra.h:567
Function< T, NDIM > & add_scalar(T t, bool fence=true)
Inplace add scalar. No communication except for optional fence.
Definition: mra.h:963
void set_thresh(double value, bool fence=true)
Sets the value of the truncation threshold. Optional global fence.
Definition: mra.h:577
void distribute(std::shared_ptr< WorldDCPmapInterface< Key< NDIM > > > newmap) const
distribute this function according to newmap
Definition: mra.h:669
Function< T, NDIM > & mirror(const Function< T, NDIM > &f, const std::vector< long > &mirrormap, bool fence)
This is replaced with mirror(f) ... private.
Definition: mra.h:1623
Future< long > evalR(const coordT &xuser) const
Evaluates the function rank at a point in user coordinates. Possible non-blocking comm.
Definition: mra.h:289
void unaryop(const opT &op, bool fence=true)
Inplace unary operation on function values.
Definition: mra.h:904
void standard(bool fence=true)
Converts the function standard compressed form. Possible non-blocking comm.
Definition: mra.h:748
void unaryop_node(const opT &op, bool fence=true)
Unary operation applied inplace to the nodes.
Definition: mra.h:924
std::size_t size() const
Returns the number of coefficients in the function ... collective global sum.
Definition: mra.h:538
bool is_on_demand() const
Definition: mra.h:636
std::size_t max_nodes() const
Returns the max number of nodes on a processor.
Definition: mra.h:523
std::shared_ptr< FunctionImpl< T, NDIM > > impl
Definition: mra.h:129
Function< T, NDIM > conj(bool fence=true)
Inplace complex conjugate. No communication except for optional fence.
Definition: mra.h:942
Function< T, NDIM > & scale(const Q q, bool fence=true)
Inplace, scale the function by a constant. No communication except for optional fence.
Definition: mra.h:953
void replicate(bool fence=true) const
replicate this function, generating a unique pmap
Definition: mra.h:663
Function< typename opT::resultT, NDIM > & unary_op_coeffs(const Function< Q, NDIM > &func, const opT &op, bool fence)
This is replaced with left*right ... private.
Definition: mra.h:1480
T trace_local() const
Returns local contribution to int(f(x),x) ... no communication.
Definition: mra.h:1090
void mul_on_demand(const Function< L, NDIM > &f, const Function< R, NDIM > &g, bool fence=true)
Same as operator* but with optional fence and no automatic reconstruction.
Definition: mra.h:1565
Vector< double, NDIM > coordT
Type of vector holding coordinates.
Definition: mra.h:139
void store(Archive &ar) const
Stores the function to an archive.
Definition: mra.h:1460
Function< T, NDIM > & square(bool fence=true)
Inplace squaring of function ... global comm only if not reconstructed.
Definition: mra.h:1059
std::size_t max_local_depth() const
Returns the maximum local depth of the function tree ... no communications.
Definition: mra.h:515
Function< T, NDIM > & operator-=(const Function< Q, NDIM > &other)
Inplace subtraction of functions in the wavelet basis.
Definition: mra.h:1028
void vtransform(const std::vector< Function< R, NDIM > > &v, const Tensor< Q > &c, std::vector< Function< T, NDIM > > &vresult, double tol, bool fence=true)
sparse transformation of a vector of functions ... private
Definition: mra.h:1583
Function< T, NDIM > & fill_tree(bool fence=true)
With this being an on-demand function, fill the MRA tree according to different criteria.
Definition: mra.h:1154
Function< T, NDIM > & truncate(double tol=0.0, bool fence=true)
Truncate the function with optional fence. Compresses with fence if not compressed.
Definition: mra.h:602
impl world gop sum(local)
TENSOR_RESULT_TYPE(T, R) local
std::size_t max_depth() const
Returns the maximum depth of the function tree ... collective global sum.
Definition: mra.h:507
Function()
Default constructor makes uninitialized function. No communication.
Definition: mra.h:154
std::size_t tree_size() const
Returns the number of nodes in the function tree ... collective global sum.
Definition: mra.h:491
FunctionImpl< T, NDIM > implT
Definition: mra.h:136
void clear(bool fence=true)
Clears the function as if constructed uninitialized. Optional fence.
Definition: mra.h:847
void refine_general(const opT &op, bool fence=true) const
Inplace autorefines the function. Optional fence. Possible non-blocking comm.
Definition: mra.h:811
static void doconj(const Key< NDIM >, Tensor< T > &t)
Definition: mra.h:934
const std::shared_ptr< FunctionImpl< T, NDIM > > & get_impl() const
Returns a shared-pointer to the implementation.
Definition: mra.h:614
TreeState state
Definition: mra.h:1274
TENSOR_RESULT_TYPE(T, R) inner_local(const Function< R
Returns local part of inner product ... throws if both not compressed.
void set_functor(const std::shared_ptr< FunctionFunctorInterface< T, NDIM > > functor)
Replace the current functor with the provided new one.
Definition: mra.h:631
bool impl_initialized() const
Definition: mra.h:132
Future< T > eval(const coordT &xuser) const
Evaluates the function at a point in user coordinates. Possible non-blocking comm.
Definition: mra.h:186
const Function< T, NDIM > & reconstruct(bool fence=true) const
Reconstructs the function, transforming into scaling function basis. Possible non-blocking comm.
Definition: mra.h:775
return local
Definition: mra.h:1299
Function< T, NDIM > & fill_nuclear_cuspy_tree(const opT &op, const size_t particle, const bool fence=true)
Definition: mra.h:1194
auto func
Definition: mra.h:1377
~Function()
Destruction of any underlying implementation is deferred to next global fence.
Definition: mra.h:178
void do_hartree_product(const std::vector< std::shared_ptr< FunctionImpl< T, LDIM >>> left, const std::vector< std::shared_ptr< FunctionImpl< T, KDIM >>> right, const opT *op)
perform the hartree product of f*g, invoked by result
Definition: mra.h:1221
Function(const Function< T, NDIM > &f)
Copy constructor is shallow. No communication, works in either basis.
Definition: mra.h:165
void set_impl(const std::shared_ptr< FunctionImpl< T, NDIM > > &impl)
Replace current FunctionImpl with provided new one.
Definition: mra.h:621
T operator()(double x, double y=0, double z=0, double xx=0, double yy=0, double zz=0) const
Evaluates the function at a point in user coordinates. Collective operation.
Definition: mra.h:385
T inner_ext_local(const std::shared_ptr< FunctionFunctorInterface< T, NDIM > > f, const bool leaf_refine=true, const bool keep_redundant=false) const
Definition: mra.h:1310
std::size_t min_nodes() const
Returns the min number of nodes on a processor.
Definition: mra.h:530
constexpr std::size_t LDIM
Definition: mra.h:1376
return impl inner_local * g())
void do_hartree_product(const std::vector< std::shared_ptr< FunctionImpl< T, LDIM >>> left, const std::vector< std::shared_ptr< FunctionImpl< T, KDIM >>> right)
perform the hartree product of f*g, invoked by result
Definition: mra.h:1235
Function< T, NDIM > & operator+=(const Function< Q, NDIM > &other)
Inplace addition of functions in the wavelet basis.
Definition: mra.h:1008
bool is_nonstandard() const
Returns true if nonstandard-compressed, false otherwise. No communication.
Definition: mra.h:484
void verify() const
Asserts that the function is initialized.
Definition: mra.h:142
Function< T, NDIM/2 > dirac_convolution(const bool fence=true) const
Definition: mra.h:1422
double err(const funcT &func) const
Returns an estimate of the difference ||this-func|| ... global sum performed.
Definition: mra.h:440
Function< T, NDIM > & reduce_rank(const double thresh=0.0, const bool fence=true)
reduce the rank of the coefficient tensors
Definition: mra.h:1672
T inner_ext(const std::shared_ptr< FunctionFunctorInterface< T, NDIM > > f, const bool leaf_refine=true, const bool keep_redundant=false) const
Definition: mra.h:1326
void make_redundant(bool fence=true)
Converts the function to redundant form, i.e. sum coefficients on all levels.
Definition: mra.h:760
double check_symmetry() const
check symmetry of a function by computing the 2nd derivative
Definition: mra.h:1658
void multi_to_multi_op_values(const opT &op, const std::vector< Function< T, NDIM > > &vin, std::vector< Function< T, NDIM > > &vout, const bool fence=true)
apply op on the input vector yielding an output vector of functions
Definition: mra.h:1521
Function< T, NDIM > & fill_nuclear_cuspy_tree(const size_t particle, const bool fence=true)
Special refinement on 6D boxes for the nuclear potentials (regularized with cusp, non-regularized wit...
Definition: mra.h:1207
FunctionFactory< T, NDIM > factoryT
Definition: mra.h:138
bool is_initialized() const
Returns true if the function is initialized.
Definition: mra.h:147
change_tree_state(state, false)
bool is_reconstructed() const
Returns true if reconstructed, false otherwise. No communication.
Definition: mra.h:473
MADNESS_ASSERT(is_compressed())
void vmulXX(const Function< L, NDIM > &left, const std::vector< Function< R, NDIM > > &right, std::vector< Function< T, NDIM > > &result, double tol, bool fence)
Multiplication of function * vector of functions using recursive algorithm of mulxx.
Definition: mra.h:1542
Function< T, NDIM > & mapdim(const Function< T, NDIM > &f, const std::vector< long > &map, bool fence)
This is replaced with mapdim(f) ... private.
Definition: mra.h:1608
double errsq_local(const funcT &func) const
Returns an estimate of the difference ||this-func||^2 from local data.
Definition: mra.h:425
Function< T, NDIM > & abs_square(bool fence=true)
Returns *this for chaining.
Definition: mra.h:1077
void make_nonstandard(bool keepleaves, bool fence=true) const
Compresses the function retaining scaling function coeffs. Possible non-blocking comm.
Definition: mra.h:734
void print_tree(std::ostream &os=std::cout) const
Process 0 prints a summary of all nodes in the tree (collective)
Definition: mra.h:857
FunctionNode< T, NDIM > nodeT
Definition: mra.h:137
const Function< T, NDIM > & compress(bool fence=true) const
Compresses the function, transforming into wavelet basis. Possible non-blocking comm.
Definition: mra.h:721
Tensor< T > eval_cube(const Tensor< double > &cell, const std::vector< long > &npt, bool eval_refine=false) const
Evaluates a cube/slice of points (probably for plotting) ... collective but no fence necessary.
Definition: mra.h:330
Function(const factoryT &factory)
Constructor from FunctionFactory provides named parameter idiom. Possible non-blocking communication.
Definition: mra.h:158
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
void unaryop(T(*f)(T))
Inplace unary operation on function values.
Definition: mra.h:895
A future is a possibly yet unevaluated value.
Definition: future.h:373
remote_refT remote_ref(World &world) const
Returns a structure used to pass references to another process.
Definition: future.h:675
Key is the index for a node of the 2^NDIM-tree.
Definition: key.h:66
Traits class to specify support of numeric types.
Definition: type_data.h:56
A tensor is a multidimension array.
Definition: tensor.h:317
Tensor< T > & conj()
Inplace complex conjugate.
Definition: tensor.h:716
Tensor< T > & emul(const Tensor< T > &t)
Inplace multiply by corresponding elements of argument Tensor.
Definition: tensor.h:1798
void erase(const keyT &key)
Erases entry from container (non-blocking comm if remote)
Definition: worlddc.h:1105
Interface to be provided by any process map.
Definition: worlddc.h:82
void fence(bool debug=false)
Synchronizes all processes in communicator AND globally ensures no pending AM or tasks.
Definition: worldgop.cc:161
bool set_forbid_fence(bool value)
Set forbid_fence flag to new value and return old value.
Definition: worldgop.h:674
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
ProcessID size() const
Returns the number of processes in this World (same as MPI_Comm_size()).
Definition: world.h:328
unsigned long id() const
Definition: world.h:313
WorldGopInterface & gop
Global operations.
Definition: world.h:205
World * get_world() const
Returns a pointer to the world.
Definition: parallel_archive.h:130
An archive for storing local or parallel data wrapping a BinaryFstreamOutputArchive.
Definition: parallel_archive.h:321
Objects that implement their own parallel archive interface should derive from this class.
Definition: parallel_archive.h:58
static const double R
Definition: csqrt.cc:46
Declaration and initialization of tree traversal functions and generic derivative.
double(* f1)(const coord_3d &)
Definition: derivatives.cc:55
double(* f2)(const coord_3d &)
Definition: derivatives.cc:56
const double delta
Definition: dielectric_external_field.cc:119
Provides FunctionDefaults and utilities for coordinate transformation.
Provides FunctionCommonData, FunctionImpl and FunctionFactory.
Defines/implements plotting interface for functions.
Provides typedefs to hide use of templates and to increase interoperability.
auto T(World &world, response_space &f) -> response_space
Definition: global_functions.cc:34
const double beta
Definition: gygi_soltion.cc:62
static const double v
Definition: hatom_sf_dirac.cc:20
Tensor< double > op(const Tensor< double > &x)
Definition: kain.cc:508
Multidimension Key for MRA tree and associated iterators.
Implements (2nd generation) static load/data balancing for functions.
#define max(a, b)
Definition: lda.h:51
#define MADNESS_CHECK(condition)
Check a condition — even in a release build the condition is always evaluated so it can have side eff...
Definition: madness_exception.h:190
#define MADNESS_EXCEPTION(msg, value)
Macro for throwing a MADNESS exception.
Definition: madness_exception.h:119
#define MADNESS_ASSERT(condition)
Assert a condition that should be free of side-effects since in release builds this might be a no-op.
Definition: madness_exception.h:134
#define MADNESS_CHECK_THROW(condition, msg)
Check a condition — even in a release build the condition is always evaluated so it can have side eff...
Definition: madness_exception.h:210
Header to declare stuff which has not yet found a home.
static const bool VERIFY_TREE
Definition: mra.h:57
double norm(const T &t)
Definition: adquad.h:42
File holds all helper structures necessary for the CC_Operator and CC2 class.
Definition: DFParameters.h:10
Function< TENSOR_RESULT_TYPE(L, R), NDIM > sub(const Function< L, NDIM > &left, const Function< R, NDIM > &right, bool fence=true)
Same as operator- but with optional fence and no automatic compression.
Definition: mra.h:1971
Function< T, KDIM+LDIM > hartree_product(const std::vector< Function< T, KDIM >> &left, const std::vector< Function< T, LDIM >> &right)
Performs a Hartree/outer product on the two given low-dimensional function vectors.
Definition: mra.h:1832
std::vector< CCPairFunction< T, NDIM > > operator+(const std::vector< CCPairFunction< T, NDIM >> c1, const std::vector< CCPairFunction< T, NDIM > > &c2)
Definition: ccpairfunction.h:1047
double abs(double x)
Definition: complexfun.h:48
Function< TENSOR_RESULT_TYPE(L, R), NDIM > add(const Function< L, NDIM > &left, const Function< R, NDIM > &right, bool fence=true)
Same as operator+ but with optional fence and no automatic compression.
Definition: mra.h:1927
std::enable_if_t< NDIM%2==0, Function< T, NDIM > > swap_particles(const Function< T, NDIM > &f)
swap particles 1 and 2
Definition: mra.h:2302
Function< T, NDIM > mapdim(const Function< T, NDIM > &f, const std::vector< long > &map, bool fence=true)
Generate a new function by reordering dimensions ... optional fence.
Definition: mra.h:2259
Function< TENSOR_RESULT_TYPE(Q, T), NDIM > mul(const Q alpha, const Function< T, NDIM > &f, bool fence=true)
Returns new function equal to alpha*f(x) with optional fence.
Definition: mra.h:1701
std::vector< CCPairFunction< T, NDIM > > operator-(const std::vector< CCPairFunction< T, NDIM >> c1, const std::vector< CCPairFunction< T, NDIM > > &c2)
Definition: ccpairfunction.h:1055
Function< T, NDIM > conj(const Function< T, NDIM > &f, bool fence=true)
Return the complex conjugate of the input function with the same distribution and optional fence.
Definition: mra.h:2046
response_space apply(World &world, std::vector< std::vector< std::shared_ptr< real_convolution_3d >>> &op, response_space &f)
Definition: basic_operators.cc:39
Function< double, NDIM > abssq(const Function< double_complex, NDIM > &z, bool fence=true)
Returns a new function that is the square of the absolute value of the input.
Definition: mra.h:2713
Function< typename opT::resultT, NDIM > unary_op_coeffs(const Function< Q, NDIM > &func, const opT &op, bool fence=true)
Out of place application of unary operation to scaling function coefficients with optional fence.
Definition: mra.h:1793
Function< typename TensorTypeData< Q >::scalar_type, NDIM > abs_square(const Function< Q, NDIM > &func)
Definition: complexfun.h:121
Function< T, NDIM > multiply(const Function< T, NDIM > f, const Function< T, LDIM > g, const int particle, const bool fence=true)
multiply a high-dimensional function with a low-dimensional function
Definition: mra.h:2372
Function< T, NDIM > multiop_values(const opT &op, const std::vector< Function< T, NDIM > > &vf)
Definition: mra.h:1691
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
const std::vector< Function< T, NDIM > > & change_tree_state(const std::vector< Function< T, NDIM >> &v, const TreeState finalstate, const bool fence=true)
change tree state of the functions
Definition: vmra.h:277
TreeState
Definition: funcdefaults.h:58
@ nonstandard_after_apply
s and d coeffs, state after operator application
Definition: funcdefaults.h:63
@ redundant_after_merge
s coeffs everywhere, must be summed up to yield the result
Definition: funcdefaults.h:65
@ reconstructed
s coeffs at the leaves only
Definition: funcdefaults.h:59
@ nonstandard
s and d coeffs in internal nodes
Definition: funcdefaults.h:61
@ unknown
Definition: funcdefaults.h:67
@ compressed
d coeffs in internal nodes, s and d coeffs at the root
Definition: funcdefaults.h:60
@ redundant
s coeffs everywhere
Definition: funcdefaults.h:64
@ nonstandard_with_leaves
like nonstandard, with s coeffs at the leaves
Definition: funcdefaults.h:62
static void user_to_sim(const Vector< double, NDIM > &xuser, Vector< double, NDIM > &xsim)
Convert user coords (cell[][]) to simulation coords ([0,1]^ndim)
Definition: funcdefaults.h:524
void standard(World &world, std::vector< Function< T, NDIM > > &v, bool fence=true)
Generates standard form of a vector of functions.
Definition: vmra.h:260
Function< T, NDIM > copy(const Function< T, NDIM > &f, const std::shared_ptr< WorldDCPmapInterface< Key< NDIM > > > &pmap, bool fence=true)
Create a new copy of the function with different distribution and optional fence.
Definition: mra.h:2002
Function< Q, NDIM > convert(const Function< T, NDIM > &f, bool fence=true)
Type conversion implies a deep copy. No communication except for optional fence.
Definition: mra.h:2032
Function< T, NDIM > square(const Function< T, NDIM > &f, bool fence=true)
Create a new function that is the square of f - global comm only if not reconstructed.
Definition: mra.h:2681
std::vector< std::shared_ptr< FunctionImpl< T, NDIM > > > get_impl(const std::vector< Function< T, NDIM >> &v)
Definition: vmra.h:663
Function< T, NDIM > map_and_mirror(const Function< T, NDIM > &f, const std::vector< long > &map, const std::vector< long > &mirror, bool fence=true)
This is replaced with mirror(map(f)), optional fence.
Definition: mra.h:2288
int Level
Definition: key.h:55
std::vector< Function< TENSOR_RESULT_TYPE(L, R), D > > vmulXX(const Function< L, D > &left, const std::vector< Function< R, D > > &vright, double tol, bool fence=true)
Use the vmra/mul(...) interface instead.
Definition: mra.h:1807
std::string get_mra_data_dir()
Definition: startup.cc:208
Function< TENSOR_RESULT_TYPE(L, R), NDIM > gaxpy_oop(TENSOR_RESULT_TYPE(L, R) alpha, const Function< L, NDIM > &left, TENSOR_RESULT_TYPE(L, R) beta, const Function< R, NDIM > &right, bool fence=true)
Returns new function alpha*left + beta*right optional fence and no automatic compression.
Definition: mra.h:1917
std::shared_ptr< FunctionFunctorInterface< double, 3 > > func(new opT(g))
Function< TENSOR_RESULT_TYPE(L, R), NDIM > binary_op(const Function< L, NDIM > &left, const Function< R, NDIM > &right, const opT &op, bool fence=true)
Generate new function = op(left,right) where op acts on the function values.
Definition: mra.h:1766
TENSOR_RESULT_TYPE(T, R) inner(const Function< T
Computes the scalar/inner product between two functions.
Definition: vmra.h:1059
std::vector< Function< TENSOR_RESULT_TYPE(T, R), NDIM > > innerXX(const Function< T, LDIM > &f, const std::vector< Function< R, KDIM >> &vg, const std::array< int, CDIM > v1, const std::array< int, CDIM > v2, int task=0)
Computes the partial scalar/inner product between two functions, returns a low-dim function.
Definition: mra.h:2430
void print(const T &t, const Ts &... ts)
Print items to std::cout (items separated by spaces) and terminate with a new line.
Definition: print.h:225
@ TT_2D
Definition: gentensor.h:120
NDIM & f
Definition: mra.h:2416
Function< T, NDIM > symmetrize(const Function< T, NDIM > &f, const std::string symmetry, bool fence=true)
symmetrize a function
Definition: mra.h:2328
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
NDIM const Function< R, NDIM > & g
Definition: mra.h:2416
std::vector< CCPairFunction< T, NDIM > > operator*(const double fac, const std::vector< CCPairFunction< T, NDIM > > &arg)
Definition: ccpairfunction.h:1084
double inner(response_space &a, response_space &b)
Definition: response_functions.h:442
double imag(double x)
Definition: complexfun.h:56
void startup(World &world, int argc, char **argv, bool doprint=false, bool make_stdcout_nice_to_reals=true)
Definition: startup.cc:64
std::string type(const PairType &n)
Definition: PNOParameters.h:18
static bool print_timings
Definition: SCF.cc:104
Function< typename opT::resultT, NDIM > unary_op(const Function< Q, NDIM > &func, const opT &op, bool fence=true)
Out of place application of unary operation to function values with optional fence.
Definition: mra.h:1780
Function< T, NDIM > gaxpy_oop_reconstructed(const double alpha, const Function< T, NDIM > &left, const double beta, const Function< T, NDIM > &right, const bool fence=true)
Returns new function alpha*left + beta*right optional fence, having both addends reconstructed.
Definition: mra.h:1935
Function< T, NDIM > mirror(const Function< T, NDIM > &f, const std::vector< long > &mirrormap, bool fence=true)
Generate a new function by mirroring within the dimensions .. optional fence.
Definition: mra.h:2271
void load(Function< T, NDIM > &f, const std::string name)
Definition: mra.h:2751
double real(double x)
Definition: complexfun.h:52
Function< TENSOR_RESULT_TYPE(typename opT::opT, R), NDIM > apply_only(const opT &op, const Function< R, NDIM > &f, bool fence=true)
Apply operator ONLY in non-standard form - required other steps missing !!
Definition: mra.h:2120
std::string name(const FuncType &type, const int ex=-1)
Definition: ccpairfunction.h:28
void save(const Function< T, NDIM > &f, const std::string name)
Definition: mra.h:2745
Function< TENSOR_RESULT_TYPE(typename opT::opT, R), NDIM > apply_1d_realspace_push(const opT &op, const Function< R, NDIM > &f, int axis, bool fence=true)
Definition: mra.h:2231
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
const std::vector< Function< T, NDIM > > & reconstruct(const std::vector< Function< T, NDIM > > &v)
reconstruct a vector of functions
Definition: vmra.h:156
Implements most functionality of separated operators.
Implements ParallelInputArchive and ParallelOutputArchive for parallel serialization of data.
double Q(double a)
Definition: relops.cc:20
static const double c
Definition: relops.cc:10
static const double L
Definition: rk.cc:46
static const double thresh
Definition: rk.cc:45
static const long k
Definition: rk.cc:44
Definition: test_ar.cc:204
void serialize(Archive &ar)
Definition: mra.h:891
T(* f)(T)
Definition: mra.h:886
void operator()(const Key< NDIM > &key, Tensor< T > &t) const
Definition: mra.h:888
SimpleUnaryOpWrapper(T(*f)(T))
Definition: mra.h:887
void serialize(Archive &ar)
Definition: mra.h:824
bool operator()(implT *impl, const Key< NDIM > &key, const nodeT &t) const
Definition: mra.h:820
Definition: type_data.h:146
TensorArgs holds the arguments for creating a LowRankTensor.
Definition: gentensor.h:134
Default load of an object via serialize(ar, t).
Definition: archive.h:666
static void store(const ParallelOutputArchive< archiveT > &ar, const Function< T, NDIM > &f)
Definition: mra.h:2738
Default store of an object via serialize(ar, t).
Definition: archive.h:611
double resultT
Definition: mra.h:2649
void serialize(Archive &ar)
Definition: mra.h:2655
Tensor< double > operator()(const Key< NDIM > &key, const Tensor< double_complex > &t) const
Definition: mra.h:2650
void serialize(Archive &ar)
Definition: mra.h:2644
Tensor< double > operator()(const Key< NDIM > &key, const Tensor< double_complex > &t) const
Definition: mra.h:2639
double resultT
Definition: mra.h:2638
void serialize(Archive &ar)
Definition: mra.h:2633
Tensor< double > operator()(const Key< NDIM > &key, const Tensor< double_complex > &t) const
Definition: mra.h:2629
double resultT
Definition: mra.h:2628
Tensor< double > operator()(const Key< NDIM > &key, const Tensor< double_complex > &t) const
Definition: mra.h:2619
double resultT
Definition: mra.h:2618
void serialize(Archive &ar)
Definition: mra.h:2623
Definition: funcimpl.h:606
returns true if the result of a hartree_product is a leaf node (compute norm & error)
Definition: funcimpl.h:496
Definition: funcimpl.h:560
Definition: lowrankfunction.h:332
Defines and implements most of Tensor.
#define UNARY_OPTIMIZED_ITERATOR(X, x, exp)
Definition: tensor_macros.h:658
const double alpha
Definition: test_coulomb.cc:54
int task(int i)
Definition: test_runtime.cpp:4
void d()
Definition: test_sig.cc:79
void e()
Definition: test_sig.cc:75
static const std::size_t NDIM
Definition: testpdiff.cc:42
std::size_t axis
Definition: testpdiff.cc:59
Defines operations on vectors of Functions.
Implements WorldContainer.
#define PROFILE_FUNC
Definition: worldprofile.h:209
#define PROFILE_MEMBER_FUNC(classname)
Definition: worldprofile.h:210