33#ifndef MADNESS_MRA_VMRA_H__INCLUDED
34#define MADNESS_MRA_VMRA_H__INCLUDED
132 template <
typename T, std::
size_t NDIM>
155 template <
typename T, std::
size_t NDIM>
164 template <
typename T, std::
size_t NDIM>
170 template <
typename T, std::
size_t NDIM>
188 template <
typename T, std::
size_t NDIM>
207 template <
typename T, std::
size_t NDIM>
210 for (
const auto&
f :
vf)
f.refine(
false);
217 template <
typename T, std::
size_t NDIM>
223 std::vector<FunctionImpl<T,NDIM>*>
v_ptr;
226 for (
unsigned int i=0; i<
vf.size(); ++i) {
232 typename std::vector<FunctionImpl<T, NDIM>*>::iterator it;
236 std::vector< Tensor<T> >
c(
v_ptr.size());
238 if (fence)
v_ptr[0]->world.gop.fence();
244 template <
typename T, std::
size_t NDIM>
259 template <
typename T, std::
size_t NDIM>
276 template <
typename T, std::
size_t NDIM>
279 const bool fence=
true) {
280 if (
v.size()==0)
return v;
283 for (
const auto&
f :
v)
284 if (
f.is_initialized()) {
288 if (
not dummy.is_initialized())
return v;
334 template <
typename T, std::
size_t NDIM>
345 for (
unsigned int i=0; i<
v.size(); ++i) {
346 v[i].truncate(tol,
false);
355 template <
typename T, std::
size_t NDIM>
357 double tol=0.0,
bool fence=
true) {
358 if (
v.size()>0)
truncate(
v[0].world(),
v,tol,fence);
365 template <
typename T, std::
size_t NDIM>
367 double thresh=0.0,
bool fence=
true) {
368 if (
v.size()==0)
return v;
369 for (
auto& vv :
v) vv.reduce_rank(
thresh,
false);
370 if (fence)
v[0].world().gop.fence();
376 template <
typename T, std::
size_t NDIM>
377 std::vector< Function<T,NDIM> >
384 std::vector< Function<T,NDIM> >
df(
v.size());
385 for (
unsigned int i=0; i<
v.size(); ++i) {
386 df[i] =
D(
v[i],
false);
393 template <
typename T, std::
size_t NDIM>
394 std::vector< Function<T,NDIM> >
397 std::vector< Function<T,NDIM> > r(n);
398 for (
int i=0; i<n; ++i)
407 template <
typename T, std::
size_t NDIM>
408 std::vector< Function<T,NDIM> >
411 std::vector< Function<T,NDIM> > r(n);
412 for (
int i=0; i<n; ++i)
420 template<
typename T, std::
size_t NDIM>
422 if (
vf_in.size()==0)
return std::vector<Function<T,NDIM>>();
431 for (
int i=0; i<s.dim(0); ++i)
Q(i,i) += 1.5;
438 for (
int i=0; i<
Q.dim(0); ++i)
439 for (
int j=0; j<i; ++j)
455 template <
typename T, std::
size_t NDIM>
459 if(
v.empty())
return v;
466 for(
size_t i=0;i<
v.size();++i) s(i)=1.0/(
sqrt(s(i)));
470 for(
size_t i=0;i<
v.size();++i){
471 for(
size_t j=0;j<
v.size();++j){
478 World& world=
v.front().world();
485 template <
typename T, std::
size_t NDIM>
487 if(
v.empty())
return v;
490 World& world=
v.front().world();
500 template <
typename T, std::
size_t NDIM>
506 if(
v.empty())
return v;
516 for(
size_t i=0;i<
v.size();++i) {
527 for(
size_t i=0;i<
v.size();++i){
528 for(
size_t j=0;j<
v.size();++j){
529 U(i,j)=U(i,j)*(
sqrts(j));
534 World& world=
v.front().world();
542 template <
typename T, std::
size_t NDIM>
545 if(
v.empty())
return v;
547 World& world=
v.front().world();
556 template <
typename T, std::
size_t NDIM>
561 if (
v.empty())
return v;
569 World& world=
v.front().world();
577 template <
typename T, std::
size_t NDIM>
579 if(
v.empty())
return v;
581 World& world=
v.front().world();
593 template <
typename T, std::
size_t NDIM>
608 std::vector<Function<T,NDIM> >
pv(rank);
618 World& world=
v.front().world();
626 template <
typename T, std::
size_t NDIM>
636 template <
typename T, std::
size_t NDIM>
642 World& world=
v.front().world();
648 template <
typename T, std::
size_t NDIM>
650 std::vector<Function<T,NDIM> >
v=
lhs;
651 for (std::size_t i = 0; i <
rhs.size(); ++i)
v.push_back(
rhs[i]);
655 template <
typename T, std::
size_t NDIM>
657 std::vector<Function<T,NDIM> >result;
658 for(
const auto& x:vv) result=
append(result,x);
662 template<
typename T, std::
size_t NDIM>
664 std::vector<std::shared_ptr<FunctionImpl<T,NDIM>>> result;
665 for (
auto&
f :
v) result.push_back(
f.get_impl());
669 template<
typename T, std::
size_t NDIM>
672 for (std::size_t i=0; i<vimpl.size(); ++i)
v[i].
set_impl(vimpl[i]);
675 template<
typename T, std::
size_t NDIM>
677 std::vector<Function<T,NDIM>>
v(vimpl.size());
678 for (std::size_t i=0; i<vimpl.size(); ++i)
v[i].
set_impl(vimpl[i]);
687 template <
typename T,
typename R, std::
size_t NDIM>
703 for (
int i=0; i<
m; ++i) {
704 for (
int j=0; j<n; ++j) {
705 if (
c(j,i) !=
R(0.0))
vc[i].gaxpy(resultT(1.0),
v[j],resultT(
c(j,i)),
false);
716 template <
typename T,
typename R, std::
size_t NDIM>
732 vv.get_impl()->get_tree_state()==
reconstructed,
"trees have to be reconstructed in transform_reconstructed");
736 for (
int i=0; i<
m; ++i) {
738 for (
int j=0; j<n; ++j) {
739 if (
c(j,i) !=
R(0.0))
v[j].get_impl()->accumulate_trees(*(result[i].
get_impl()),resultT(
c(j,i)),
true);
747 for (
auto& r : result) r.get_impl()->finalize_sum();
755 template <
typename L,
typename R, std::
size_t NDIM>
770 template <
typename T,
typename R, std::
size_t NDIM>
786 c.copy_to_replicated(
tmp);
792 for (
int i=0; i<
m; ++i) {
793 for (
int j=0; j<n; ++j) {
794 if (
tmp(j,i) !=
R(0.0))
vc[i].gaxpy(1.0,
v[j],
tmp(j,i),
false);
804 template <
typename T,
typename Q, std::
size_t NDIM>
810 for (
unsigned int i=0; i<
v.size(); ++i)
v[i].
scale(
factors[i],
false);
815 template <
typename T,
typename Q, std::
size_t NDIM>
821 for (
unsigned int i=0; i<
v.size(); ++i)
v[i].
scale(factor,
false);
826 template <
typename T, std::
size_t NDIM>
830 std::vector<double>
norms(
v.size());
831 for (
unsigned int i=0; i<
v.size(); ++i)
norms[i] =
v[i].norm2sq_local();
838 template <
typename T, std::
size_t NDIM>
842 for (
unsigned int i = 0; i <
v.size(); ++i)
norms[i] =
v[i].norm2sq_local();
850 template <
typename T, std::
size_t NDIM>
853 if (
v.size()==0)
return 0.0;
854 std::vector<double>
norms(
v.size());
855 for (
unsigned int i=0; i<
v.size(); ++i)
norms[i] =
v[i].norm2sq_local();
857 for (
unsigned int i=1; i<
v.size(); ++i)
norms[0] +=
norms[i];
905 template <
typename T, std::
size_t NDIM>
922 std::vector< Function<T,NDIM> >
ivec(
f.begin()+ilo,
f.begin()+ihi);
925 std::vector< Function<T,NDIM> >
jvec(
g.begin()+jlo,
g.begin()+jhi);
928 A.copy_from_replicated_patch(ilo, ihi - 1, jlo, jhi - 1,
P);
939 template <
typename T,
typename R, std::
size_t NDIM>
951 std::vector<const FunctionImpl<T,NDIM>*> left(
f.size());
952 std::vector<const FunctionImpl<R,NDIM>*> right(
g.size());
953 for (
unsigned int i=0; i<
f.size(); i++) left[i] =
f[i].
get_impl().get();
954 for (
unsigned int i=0; i<
g.size(); i++) right[i]=
g[i].
get_impl().get();
969 template <
typename T,
typename R, std::
size_t NDIM>
975 long n=
f.size(),
m=
g.size();
981 if ((
void*)(&
f) != (
void*)(&
g))
compress(world,
g);
983 for (
long i=0; i<n; ++i) {
986 for (
long j=0; j<jtop; ++j) {
987 r(i,j) =
f[i].inner_local(
g[j]);
988 if (
sym) r(j,i) =
conj(r(i,j));
1011 template <
typename T,
typename R, std::
size_t NDIM>
1016 long n=
f.size(),
m=
g.size();
1023 for (
long i=0; i<n; ++i) {
1024 r(i) =
f[i].inner_local(
g[i]);
1035 template <
typename T,
typename R, std::
size_t NDIM>
1046 for (
long i=0; i<n; ++i) {
1047 r(i) =
f.inner_local(
g[i]);
1058 template <
typename T,
typename R, std::
size_t NDIM>
1062 if(
f.empty())
return 0.0;
1063 else return inner(
f[0].world(),
f,
g).sum();
1068 template <
typename T,
typename R, std::
size_t NDIM>
1075 a.reconstruct(
false);
1082 template <
typename T,
typename R, std::
size_t NDIM>
1090 a.reconstruct(
false);
1093 for (
unsigned int i=0; i<
v.size(); ++i) {
1094 v[i].norm_tree(
false);
1113 template <
typename T,
typename R, std::
size_t NDIM>
1120 bool symm =
false) {
1126 for (
auto&
ff :
f)
ff.norm_tree(
false);
1130 std::vector<std::vector<Function<R,NDIM> > >result(
f.size());
1131 std::vector<Function<R,NDIM>>
g_i;
1132 for (
int64_t i=
f.size()-1; i>=0; --i) {
1134 result[i]=
vmulXX(
f[i],
g, tol,
false);
1146 template <
typename T, std::
size_t NDIM>
1152 for (
unsigned int i=0; i<
v.size(); ++i) {
1153 v[i].norm_tree(
false);
1159 template <
typename T,
typename R, std::
size_t NDIM>
1171 for (
unsigned int i=0; i<
a.size(); ++i) {
1172 q[i] =
mul(
a[i],
b[i],
false);
1185 template<
typename T, std::
size_t NDIM, std::
size_t LDIM>
1190 std::vector<Function<T,NDIM> > result(
g.size());
1191 for (
auto& r : result) r.set_impl(
f,
false);
1199 for (std::size_t i=0; i<result.size(); ++i) {
1206 for (
auto&
ig :
g)
ig.get_impl()->undo_redundant(
false);
1211 template<
typename T, std::
size_t NDIM, std::
size_t LDIM>
1213 const std::tuple<int,int,int>
v) {
1219 template <
typename T, std::
size_t NDIM>
1220 std::vector< Function<T,NDIM> >
1235 template <
typename T, std::
size_t NDIM>
1236 std::vector< Function<typename Tensor<T>::scalar_type,
NDIM> >
1242 std::vector<Function<scalartype,NDIM> > result(
v.size());
1243 for (
size_t i=0; i<
v.size(); ++i) result[i]=
abs_square(
v[i],
false);
1250 template <
typename T, std::
size_t NDIM>
1252 for (
unsigned int j=0; j<
v.size(); ++j) {
1253 v[j].set_thresh(
thresh,
false);
1259 template <
typename T, std::
size_t NDIM>
1260 std::vector< Function<T,NDIM> >
1265 std::vector< Function<T,NDIM> > r =
copy(world,
v);
1266 for (
unsigned int i=0; i<
v.size(); ++i) {
1274 template <
typename T,
typename R, std::
size_t NDIM>
1278 std::vector< Function<R,NDIM> > r(
v.size());
1279 for (
unsigned int i=0; i<
v.size(); ++i) {
1288 template <
typename T, std::
size_t NDIM>
1289 std::vector< Function<T,NDIM> >
1294 std::vector< Function<T,NDIM> > r(
v.size());
1295 for (
unsigned int i=0; i<
v.size(); ++i) {
1296 r[i] =
copy(
v[i],
false);
1304 template <
typename T, std::
size_t NDIM>
1305 std::vector< Function<T,NDIM> >
1308 std::vector< Function<T,NDIM> > r(
v.size());
1309 if (
v.size()>0) r=
copy(
v.front().world(),
v,fence);
1314 template <
typename T, std::
size_t NDIM>
1315 std::vector< Function<T,NDIM> >
1318 const unsigned int n,
1321 std::vector< Function<T,NDIM> > r(n);
1322 for (
unsigned int i=0; i<n; ++i) {
1323 r[i] =
copy(
v,
false);
1338 template <
typename T, std::
size_t NDIM>
1342 bool fence =
true) {
1344 std::vector<Function<T, NDIM>> r(
v.size());
1345 for (
unsigned int i = 0; i <
v.size(); ++i) {
1346 r[i] =
copy(
v[i], pmap,
false);
1353 template <
typename T,
typename R, std::
size_t NDIM>
1365 for (
unsigned int i=0; i<
a.size(); ++i) {
1366 r[i] =
add(
a[i],
b[i],
false);
1373 template <
typename T,
typename R, std::
size_t NDIM>
1384 for (
unsigned int i=0; i<
b.size(); ++i) {
1385 r[i] =
add(
a,
b[i],
false);
1390 template <
typename T,
typename R, std::
size_t NDIM>
1396 return add(world,
a,
b, fence);
1400 template <
typename T,
typename R, std::
size_t NDIM>
1412 for (
unsigned int i=0; i<
a.size(); ++i) {
1413 r[i] =
sub(
a[i],
b[i],
false);
1420 template <
typename T, std::
size_t NDIM>
1427 for (
unsigned int i=0; i<
f.size(); ++i) r.
gaxpy(1.0,
f[i],1.0,
false);
1432 template <
typename T, std::
size_t NDIM>
1449 std::vector<Function<T, NDIM>>
ivec(
f.begin() + ilo,
f.begin() + ihi);
1452 std::vector<Function<T, NDIM>>
jvec(
g.begin() + jlo,
g.begin() + jhi);
1455 A.copy_from_replicated_patch(ilo, ihi - 1, jlo, jhi - 1,
P);
1466 template <
typename T,
typename R, std::
size_t NDIM>
1477 std::vector<const FunctionImpl<T, NDIM>*> left(
f.size());
1478 std::vector<const FunctionImpl<R, NDIM>*> right(
g.size());
1479 for (
unsigned int i = 0; i <
f.size(); i++) left[i] =
f[i].
get_impl().get();
1480 for (
unsigned int i = 0; i <
g.size(); i++) right[i] =
g[i].
get_impl().get();
1491 template <
typename T,
typename R, std::
size_t NDIM>
1498 return sum(world,
mul(world,
a,
b,
true),fence);
1503 template <
typename T,
typename Q,
typename R, std::
size_t NDIM>
1513 if (
a.size()==0)
return std::vector<Function<resultT,NDIM> >();
1515 World& world=
a[0].world();
1516 std::vector<Function<resultT,NDIM> > result(
a.size());
1520 for (
unsigned int i=0; i<
a.size(); ++i) {
1526 for (
unsigned int i=0; i<
a.size(); ++i) {
1531 if (fence) world.gop.fence();
1537 template <
typename T,
typename Q,
typename R, std::
size_t NDIM>
1546 if (
a.size()==0)
return std::vector<Function<resultT,NDIM> >();
1548 World& world=
a[0].world();
1551 std::vector<Function<resultT,NDIM> > result(
a.size());
1552 for (
unsigned int i=0; i<
a.size(); ++i) {
1555 if (fence) world.gop.fence();
1561 template <
typename T,
typename Q,
typename R, std::
size_t NDIM>
1563 if (
a.size() == 0)
return;
1564 World& world=
a.front().world();
1569 template <
typename T,
typename Q,
typename R, std::
size_t NDIM>
1585 for (
unsigned int i=0; i<
a.size(); ++i) {
1593 template <
typename opT,
typename R, std::
size_t NDIM>
1596 const std::vector< std::shared_ptr<opT> >&
op,
1602 std::vector< Function<R,NDIM> >&
ncf = *
const_cast< std::vector< Function<R,NDIM>
>* >(&
f);
1608 for (
unsigned int i=0; i<
f.size(); ++i) {
1624 template <
typename T,
typename R, std::
size_t NDIM, std::
size_t KDIM>
1633 template <
typename T,
typename R, std::
size_t NDIM, std::
size_t KDIM>
1640 std::vector< Function<R,NDIM> >&
ncf = *
const_cast< std::vector< Function<R,NDIM>
>* >(&
f);
1650 for (
unsigned int i=0; i<
f.size(); ++i) {
1657 if (
op.destructive()) {
1658 for (
auto&
ff :
ncf)
ff.clear(
false);
1666 for (
auto& r : result) r.get_impl()->finalize_apply();
1670 for (
auto& r : result) r.get_impl()->print_timer();
1679 template <
typename T, std::
size_t NDIM>
1682 std::vector<double>
nn =
norm2s(world,
v);
1683 for (
unsigned int i=0; i<
v.size(); ++i)
v[i].
scale(1.0/
nn[i],
false);
1687 template <
typename T, std::
size_t NDIM>
1690 if(world.
rank()==0) std::cout <<
"print_size: " << msg <<
" is empty" << std::endl;
1691 }
else if(
v.size()==1){
1692 v.front().print_size(msg);
1701 template <
typename T, std::
size_t NDIM>
1705 if (x.is_initialized()) size+=x.size_local();
1707 const double d=
sizeof(
T);
1708 const double fac=1024*1024*1024;
1713 template <
typename T, std::
size_t NDIM>
1720 template <
typename T, std::
size_t NDIM>
1723 if (
v.empty())
return 0.0;
1725 const double d=
sizeof(
T);
1726 const double fac=1024*1024*1024;
1729 for(
unsigned int i=0;i<
v.size();i++){
1730 if (
v[i].is_initialized()) size+=
v[i].size();
1738 template <
typename T, std::
size_t NDIM>
1740 const double d=
sizeof(
T);
1741 const double fac=1024*1024*1024;
1742 double size=
f.size();
1751 template <
typename T,
typename opT, std::
size_t NDIM>
1754 const bool fence=
true) {
1761 for (
auto& out :
vout) out.set_impl(
vin[0],
false);
1772 template <
typename T, std::
size_t NDIM>
1780 template <
typename T, std::
size_t NDIM>
1788 template <
typename T, std::
size_t NDIM>
1796 template <
typename T, std::
size_t NDIM>
1804 template <
typename T, std::
size_t NDIM>
1812 template <
typename T, std::
size_t NDIM>
1820 template <
typename T,
typename R, std::
size_t NDIM>
1824 std::vector<Function<T,NDIM> >
tmp=
copy(
rhs[0].world(),
rhs);
1831 template <
typename T,
typename R, std::
size_t NDIM>
1843 template <
typename T,
typename R, std::
size_t NDIM>
1846 if (
v.size()>0)
return mul(
v[0].world(),
a,
v,
true);
1852 template <
typename T,
typename R, std::
size_t NDIM>
1855 if (
v.size()>0)
return mul(
v[0].world(),
a,
v,
true);
1860 template <
typename T, std::
size_t NDIM>
1867 template <
typename T, std::
size_t NDIM>
1876 template <
typename T, std::
size_t NDIM>
1877 std::vector<Function<typename Tensor<T>::scalar_type,
NDIM> >
1879 std::vector<Function<typename Tensor<T>::scalar_type,
NDIM> > result(
v.size());
1880 for (std::size_t i=0; i<
v.size(); ++i) result[i]=
real(
v[i],
false);
1881 if (fence
and result.size()>0) result[0].world().gop.fence();
1886 template <
typename T, std::
size_t NDIM>
1887 std::vector<Function<typename Tensor<T>::scalar_type,
NDIM> >
1889 std::vector<Function<typename Tensor<T>::scalar_type,
NDIM> > result(
v.size());
1890 for (std::size_t i=0; i<
v.size(); ++i) result[i]=
imag(
v[i],
false);
1891 if (fence
and result.size()>0) result[0].world().gop.fence();
1902 template <
typename T, std::
size_t NDIM>
1904 bool refine=
false,
bool fence=
true) {
1910 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
1913 std::vector<Function<T,NDIM> > result(
NDIM);
1914 for (
size_t i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
1920 template <
typename T, std::
size_t NDIM>
1922 bool refine=
false,
bool fence=
true) {
1928 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
1932 for (
unsigned int i=0; i<
NDIM; ++i) (*
grad[i]).set_ble1();
1934 std::vector<Function<T,NDIM> > result(
NDIM);
1935 for (
unsigned int i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
1941 template <
typename T, std::
size_t NDIM>
1943 bool refine=
false,
bool fence=
true) {
1949 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
1953 for (
unsigned int i=0; i<
NDIM; ++i) (*
grad[i]).set_ble2();
1955 std::vector<Function<T,NDIM> > result(
NDIM);
1956 for (
unsigned int i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
1962 template <
typename T, std::
size_t NDIM>
1964 bool refine=
false,
bool fence=
true) {
1970 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
1974 for (
unsigned int i=0; i<
NDIM; ++i) (*
grad[i]).set_bspline1();
1976 std::vector<Function<T,NDIM> > result(
NDIM);
1977 for (
unsigned int i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
1983 template <
typename T, std::
size_t NDIM>
1985 bool refine=
false,
bool fence=
true) {
1991 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
1995 for (
unsigned int i=0; i<
NDIM; ++i) (*
grad[i]).set_bspline2();
1997 std::vector<Function<T,NDIM> > result(
NDIM);
1998 for (
unsigned int i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
2004 template <
typename T, std::
size_t NDIM>
2006 bool refine=
false,
bool fence=
true) {
2012 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
2016 for (
unsigned int i=0; i<
NDIM; ++i) (*
grad[i]).set_bspline3();
2018 std::vector<Function<T,NDIM> > result(
NDIM);
2019 for (
unsigned int i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
f,
false);
2034 template <
typename T, std::
size_t NDIM>
2036 bool do_refine=
false,
bool fence=
true) {
2039 World& world=
v[0].world();
2043 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
2046 std::vector<Function<T,NDIM> > result(
NDIM);
2047 for (
size_t i=0; i<
NDIM; ++i) result[i]=
apply(*(
grad[i]),
v[i],
false);
2049 return sum(world,result,fence);
2060 template <
typename T, std::
size_t NDIM>
2062 bool do_refine=
false,
bool fence=
true) {
2065 World& world=
v[0].world();
2069 std::vector< std::shared_ptr< Derivative<T,NDIM> > >
grad=
2081 d[0].gaxpy(1.0,
dd[0],-1.0,
false);
2082 d[1].gaxpy(1.0,
dd[1],-1.0,
false);
2083 d[2].gaxpy(1.0,
dd[2],-1.0,
false);
2097 template <
typename T,
typename R, std::
size_t NDIM>
2100 bool do_refine=
false,
bool fence=
true) {
2104 World& world=
f[0].world();
2110 d[0]=
mul(
f[1],
g[2],
false);
2111 d[1]=
mul(
f[2],
g[0],
false);
2112 d[2]=
mul(
f[0],
g[1],
false);
2122 d[0].gaxpy(1.0,
dd[0],-1.0,
false);
2123 d[1].gaxpy(1.0,
dd[1],-1.0,
false);
2124 d[2].gaxpy(1.0,
dd[2],-1.0,
false);
2130 template<
typename T, std::
size_t NDIM>
2136 return node.coeff().size();
2147 template<
typename T,
size_t NDIM>
2149 const std::string
name) {
2150 if (world.
rank()==0)
print(
"loading vector of functions",
name);
2152 std::size_t
fsize=0;
2155 for (std::size_t i=0; i<
fsize; ++i) ar &
f[i];
2159 template<
typename T,
size_t NDIM>
2162 World& world=
f.front().world();
2163 if (world.
rank()==0)
print(
"saving vector of functions",
name);
2165 std::size_t
fsize=
f.size();
2167 for (std::size_t i=0; i<
fsize; ++i) ar &
f[i];
double q(double t)
Definition DKops.h:18
Definition test_ar.cc:118
long size() const
Returns the number of elements in the tensor.
Definition basetensor.h:138
Implements derivatives operators with variety of boundary conditions on simulation domain.
Definition derivative.h:266
Definition distributed_matrix.h:68
Manages data associated with a row/column/block distributed array.
Definition distributed_matrix.h:388
static void redistribute(World &world, const std::shared_ptr< WorldDCPmapInterface< Key< NDIM > > > &newpmap)
Sets the default process map and redistributes all functions using the old map.
Definition funcdefaults.h:414
FunctionFactory implements the named-parameter idiom for Function.
Definition function_factory.h:86
FunctionImpl holds all Function state to facilitate shallow copy semantics.
Definition funcimpl.h:945
static Tensor< TENSOR_RESULT_TYPE(T, R) > inner_local(const std::vector< const FunctionImpl< T, NDIM > * > &left, const std::vector< const FunctionImpl< R, NDIM > * > &right, bool sym)
Definition funcimpl.h:5924
void undo_redundant(const bool fence)
convert this from redundant to standard reconstructed form
Definition mraimpl.h:1560
void multiply(const implT *f, const FunctionImpl< T, LDIM > *g, const int particle)
multiply f (a pair function of NDIM) with an orbital g (LDIM=NDIM/2)
Definition funcimpl.h:3560
void change_tree_state(const TreeState finalstate, bool fence=true)
change the tree state of this function, might or might not respect fence!
Definition mraimpl.h:1386
static Tensor< TENSOR_RESULT_TYPE(T, R)> dot_local(const std::vector< const FunctionImpl< T, NDIM > * > &left, const std::vector< const FunctionImpl< R, NDIM > * > &right, bool sym)
Definition funcimpl.h:5976
FunctionNode holds the coefficients, etc., at each node of the 2^NDIM-tree.
Definition funcimpl.h:127
A multiresolution adaptive numerical function.
Definition mra.h:139
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:1006
void set_impl(const std::shared_ptr< FunctionImpl< T, NDIM > > &impl)
Replace current FunctionImpl with provided new one.
Definition mra.h:646
Key is the index for a node of the 2^NDIM-tree.
Definition key.h:68
std::shared_ptr< WorldDCPmapInterface< keyT > > load_balance(double fac=1.0, bool printstuff=false)
Actually does the partitioning of the tree.
Definition lbdeux.h:320
void add_tree(const Function< T, NDIM > &f, const costT &costfn, bool fence=false)
Accumulates cost from a function.
Definition lbdeux.h:289
Convolutions in separated form (including Gaussian)
Definition operator.h:139
A slice defines a sub-range or patch of a dimension.
Definition slice.h:103
A tensor is a multidimensional array.
Definition tensor.h:317
TensorTypeData< T >::scalar_type scalar_type
C++ typename of the real type associated with a complex type.
Definition tensor.h:409
T * ptr()
Returns a pointer to the internal data.
Definition tensor.h:1824
A simple, fixed dimension vector.
Definition vector.h:64
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
void sum(T *buf, size_t nelem)
Inplace global sum while still processing AM & tasks.
Definition worldgop.h:870
void fence()
Returns after all local tasks have completed.
Definition world_task_queue.h:1384
A parallel world class.
Definition world.h:132
WorldTaskQueue & taskq
Task queue.
Definition world.h:206
ProcessID rank() const
Returns the process rank in this World (same as MPI_Comm_rank()).
Definition world.h:320
WorldGopInterface & gop
Global operations.
Definition world.h:207
An archive for storing local or parallel data wrapping a BinaryFstreamOutputArchive.
Definition parallel_archive.h:321
int integer
Definition crayio.c:25
static const double R
Definition csqrt.cc:46
Declaration and initialization of tree traversal functions and generic derivative.
static double lo
Definition dirac-hatom.cc:23
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
#define rot(x, k)
Definition lookup3.c:72
#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:182
#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:207
Main include file for MADNESS and defines Function interface.
static const bool VERIFY_TREE
Definition mra.h:57
Namespace for all elements and tools of MADNESS.
Definition DFParameters.h:10
void save_function(const std::vector< Function< T, NDIM > > &f, const std::string name)
save a vector of functions
Definition vmra.h:2160
void rr_cholesky(Tensor< T > &A, typename Tensor< T >::scalar_type tol, Tensor< integer > &piv, int &rank)
Compute the rank-revealing Cholesky factorization.
Definition lapack.cc:1203
void make_redundant(World &world, const std::vector< Function< T, NDIM > > &v, bool fence=true)
change tree_state of a vector of functions to redundant
Definition vmra.h:189
std::vector< Function< T, NDIM > > orthonormalize_rrcd(const std::vector< Function< T, NDIM > > &v, Tensor< T > &ovlp, const double tol, Tensor< integer > &piv, int &rank)
Definition vmra.h:594
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:2745
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:1949
Function< typename TensorTypeData< Q >::scalar_type, NDIM > abs_square(const Function< Q, NDIM > &func)
Definition complexfun.h:121
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
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:2713
response_space scale(response_space a, double b)
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:2003
std::vector< Function< T, NDIM > > reduce_rank(std::vector< Function< T, NDIM > > v, double thresh=0.0, bool fence=true)
reduces the tensor rank of the coefficient tensor (if applicable)
Definition vmra.h:366
Tensor< double > norm2s_T(World &world, const std::vector< Function< T, NDIM > > &v)
Computes the 2-norms of a vector of functions.
Definition vmra.h:839
std::vector< double > norm2s(World &world, const std::vector< Function< T, NDIM > > &v)
Computes the 2-norms of a vector of functions.
Definition vmra.h:827
std::vector< Function< T, NDIM > > grad_bspline_one(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
Definition vmra.h:1963
void set_impl(std::vector< Function< T, NDIM > > &v, const std::vector< std::shared_ptr< FunctionImpl< T, NDIM > > > vimpl)
Definition vmra.h:670
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:2064
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:1733
void truncate(World &world, response_space &v, double tol, bool fence)
Definition basic_operators.cc:30
std::vector< std::shared_ptr< FunctionImpl< T, NDIM > > > get_impl(const std::vector< Function< T, NDIM > > &v)
Definition vmra.h:663
Function< T, NDIM > div(const std::vector< Function< T, NDIM > > &v, bool do_refine=false, bool fence=true)
shorthand div operator
Definition vmra.h:2035
std::vector< Function< T, NDIM > > orthonormalize_cd(const std::vector< Function< T, NDIM > > &v, Tensor< T > &ovlp)
Definition vmra.h:557
void norm_tree(World &world, const std::vector< Function< T, NDIM > > &v, bool fence=true)
Makes the norm tree for all functions in a vector.
Definition vmra.h:1147
tensorT Q2(const tensorT &s)
Given overlap matrix, return rotation with 2nd order error to orthonormalize the vectors.
Definition SCF.cc:135
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
std::vector< Function< TENSOR_RESULT_TYPE(T, R), NDIM > > cross(const std::vector< Function< T, NDIM > > &f, const std::vector< Function< R, NDIM > > &g, bool do_refine=false, bool fence=true)
shorthand cross operator
Definition vmra.h:2098
TreeState
Definition funcdefaults.h:59
@ nonstandard_after_apply
s and d coeffs, state after operator application
Definition funcdefaults.h:64
@ redundant_after_merge
s coeffs everywhere, must be summed up to yield the result
Definition funcdefaults.h:66
@ reconstructed
s coeffs at the leaves only
Definition funcdefaults.h:60
@ nonstandard
s and d coeffs in internal nodes
Definition funcdefaults.h:62
@ compressed
d coeffs in internal nodes, s and d coeffs at the root
Definition funcdefaults.h:61
@ redundant
s coeffs everywhere
Definition funcdefaults.h:65
@ nonstandard_with_leaves
like nonstandard, with s coeffs at the leaves
Definition funcdefaults.h:63
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:2078
void cholesky(Tensor< T > &A)
Compute the Cholesky factorization.
Definition lapack.cc:1174
std::vector< std::vector< Function< TENSOR_RESULT_TYPE(T, R), NDIM > > > matrix_mul_sparse(World &world, const std::vector< Function< R, NDIM > > &f, const std::vector< Function< R, NDIM > > &g, double tol, bool fence=true, bool symm=false)
Outer product of a vector of functions with a vector of functions using sparsity.
Definition vmra.h:1115
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
void compress(World &world, const std::vector< Function< T, NDIM > > &v, bool fence=true)
Compress a vector of functions.
Definition vmra.h:133
const std::vector< Function< T, NDIM > > & reconstruct(const std::vector< Function< T, NDIM > > &v)
reconstruct a vector of functions
Definition vmra.h:156
response_space transpose(response_space &f)
Definition basic_operators.cc:10
std::vector< Function< T, NDIM > > impl2function(const std::vector< std::shared_ptr< FunctionImpl< T, NDIM > > > vimpl)
Definition vmra.h:676
std::vector< Function< T, NDIM > > grad_bpsline_two(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
Definition vmra.h:1984
void set_thresh(World &world, std::vector< Function< T, NDIM > > &v, double thresh, bool fence=true)
Sets the threshold in a vector of functions.
Definition vmra.h:1251
double norm2(World &world, const std::vector< Function< T, NDIM > > &v)
Computes the 2-norm of a vector of functions.
Definition vmra.h:851
std::vector< Function< T, NDIM > > flatten(const std::vector< std::vector< Function< T, NDIM > > > &vv)
Definition vmra.h:656
std::vector< Function< T, NDIM > > orthonormalize_canonical(const std::vector< Function< T, NDIM > > &v, const Tensor< T > &ovlp, double lindep)
Definition vmra.h:501
std::vector< CCPairFunction< T, NDIM > > operator*(const double fac, const std::vector< CCPairFunction< T, NDIM > > &arg)
Definition ccpairfunction.h:1089
static void verify_tree(World &world, const std::vector< Function< T, NDIM > > &v)
Definition SCF.cc:72
std::vector< Function< T, NDIM > > multi_to_multi_op_values(const opT &op, const std::vector< Function< T, NDIM > > &vin, const bool fence=true)
apply op on the input vector yielding an output vector of functions
Definition vmra.h:1752
static const Slice _(0,-1, 1)
Tensor< T > inverse(const Tensor< T > &a_in)
invert general square matrix A
Definition lapack.cc:832
void load_balance(const real_function_6d &f, const bool leaf)
do some load-balancing
Definition madness/chem/mp2.cc:70
std::vector< CCPairFunction< T, NDIM > > operator-(const std::vector< CCPairFunction< T, NDIM > > c1, const std::vector< CCPairFunction< T, NDIM > > &c2)
Definition ccpairfunction.h:1060
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:1774
std::vector< Function< T, NDIM > > partial_mul(const Function< T, NDIM > f, const std::vector< Function< T, LDIM > > g, const int particle)
multiply a high-dimensional function with a low-dimensional function
Definition vmra.h:1186
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:1967
std::vector< Function< TENSOR_RESULT_TYPE(T, R), NDIM > > transform_reconstructed(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:718
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
response_space apply(World &world, std::vector< std::vector< std::shared_ptr< real_convolution_3d > > > &op, response_space &f)
Definition basic_operators.cc:39
std::vector< Function< T, NDIM > > append(const std::vector< Function< T, NDIM > > &lhs, const std::vector< Function< T, NDIM > > &rhs)
combine two vectors
Definition vmra.h:649
@ TT_2D
Definition gentensor.h:120
void refine(World &world, const std::vector< Function< T, NDIM > > &vf, bool fence=true)
refine the functions according to the autorefine criteria
Definition vmra.h:208
void print_size(World &world, const std::vector< Function< T, NDIM > > &v, const std::string &msg="vectorfunction")
Definition vmra.h:1688
NDIM & f
Definition mra.h:2448
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:1959
Function< TENSOR_RESULT_TYPE(T, R), NDIM > dot(World &world, const std::vector< Function< T, NDIM > > &a, const std::vector< Function< R, NDIM > > &b, bool fence=true)
Multiplies and sums two vectors of functions r = \sum_i a[i] * b[i].
Definition vmra.h:1493
std::vector< CCPairFunction< T, NDIM > > & operator-=(std::vector< CCPairFunction< T, NDIM > > &rhs, const std::vector< CCPairFunction< T, NDIM > > &lhs)
Definition ccpairfunction.h:1082
std::vector< Function< T, NDIM > > orthonormalize(const std::vector< Function< T, NDIM > > &vf_in)
orthonormalize the vectors
Definition vmra.h:421
NDIM const Function< R, NDIM > & g
Definition mra.h:2448
double wall_time()
Returns the wall time in seconds relative to an arbitrary origin.
Definition timers.cc:48
std::vector< Function< T, NDIM > > grad_ble_one(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
Definition vmra.h:1921
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:2152
std::vector< Function< T, NDIM > > grad_bspline_three(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
Definition vmra.h:2005
std::vector< Function< T, NDIM > > zero_functions_compressed(World &world, int n, bool fence=true)
Generates a vector of zero functions (compressed)
Definition vmra.h:409
double inner(response_space &a, response_space &b)
Definition response_functions.h:442
double imag(double x)
Definition complexfun.h:56
void load_function(World &world, std::vector< Function< T, NDIM > > &f, const std::string name)
load a vector of functions
Definition vmra.h:2148
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:1839
void refine_to_common_level(World &world, std::vector< Function< T, NDIM > > &vf, bool fence=true)
refine all functions to a common (finest) level
Definition vmra.h:218
std::vector< Function< T, NDIM > > grad_ble_two(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
Definition vmra.h:1942
static bool print_timings
Definition SCF.cc:104
void normalize(World &world, std::vector< Function< T, NDIM > > &v, bool fence=true)
Normalizes a vector of functions — v[i] = v[i].scale(1.0/v[i].norm2())
Definition vmra.h:1680
std::vector< Function< T, NDIM > > zero_functions(World &world, int n, bool fence=true)
Generates a vector of zero functions (reconstructed)
Definition vmra.h:395
std::vector< CCPairFunction< T, NDIM > > operator+(const std::vector< CCPairFunction< T, NDIM > > c1, const std::vector< CCPairFunction< T, NDIM > > &c2)
Definition ccpairfunction.h:1052
std::vector< Function< T, NDIM > > grad(const Function< T, NDIM > &f, bool refine=false, bool fence=true)
shorthand gradient operator
Definition vmra.h:1903
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:2404
DistributedMatrix< T > matrix_dot(const DistributedMatrixDistribution &d, const std::vector< Function< T, NDIM > > &f, const std::vector< Function< T, NDIM > > &g, bool sym=false)
Definition vmra.h:1433
std::vector< Function< T, NDIM > > orthonormalize_symmetric(const std::vector< Function< T, NDIM > > &v, const Tensor< T > &ovlp)
symmetric orthonormalization (see e.g. Szabo/Ostlund)
Definition vmra.h:456
double real(double x)
Definition complexfun.h:52
std::vector< CCPairFunction< T, NDIM > > & operator+=(std::vector< CCPairFunction< T, NDIM > > &lhs, const CCPairFunction< T, NDIM > &rhs)
Definition ccpairfunction.h:1068
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
std::string name(const FuncType &type, const int ex=-1)
Definition ccpairfunction.h:28
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
double get_size(World &world, const std::vector< Function< T, NDIM > > &v)
Definition vmra.h:1721
double get_size_local(World &world, const std::vector< Function< T, NDIM > > &v)
return the size of a vector of functions for each rank
Definition vmra.h:1702
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:2034
void syev(const Tensor< T > &A, Tensor< T > &V, Tensor< typename Tensor< T >::scalar_type > &e)
Real-symmetric or complex-Hermitian eigenproblem.
Definition lapack.cc:969
Tensor< TENSOR_RESULT_TYPE(T, R) > matrix_inner_old(World &world, const std::vector< Function< T, NDIM > > &f, const std::vector< Function< R, NDIM > > &g, bool sym=false)
Computes the matrix inner product of two function vectors - q(i,j) = inner(f[i],g[j])
Definition vmra.h:970
void make_nonstandard(World &world, std::vector< Function< T, NDIM > > &v, bool fence=true)
Generates non-standard form of a vector of functions.
Definition vmra.h:245
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:246
int distance(const madness::Hash_private::HashIterator< hashT > &it, const madness::Hash_private::HashIterator< hashT > &jt)
Definition worldhashmap.h:616
static long abs(long a)
Definition tensor.h:218
static const double b
Definition nonlinschro.cc:119
static const double d
Definition nonlinschro.cc:121
static const double a
Definition nonlinschro.cc:118
double Q(double a)
Definition relops.cc:20
static const double c
Definition relops.cc:10
static const double m
Definition relops.cc:9
static const double L
Definition rk.cc:46
static const double thresh
Definition rk.cc:45
Definition test_ar.cc:204
double operator()(const Key< 6 > &key, const FunctionNode< double, 6 > &node) const
Definition mp2.h:70
Definition lowrankfunction.h:332
Definition dirac-hatom.cc:108
AtomicInt sum
Definition test_atomicint.cc:46
int P
Definition test_binsorter.cc:9
double aa
Definition testbsh.cc:68
static const double alpha
Definition testcosine.cc:10
constexpr std::size_t NDIM
Definition testgconv.cc:54
#define TENSOR_RESULT_TYPE(L, R)
This macro simplifies access to TensorResultType.
Definition type_data.h:205
#define PROFILE_FUNC
Definition worldprofile.h:209
#define PROFILE_BLOCK(name)
Definition worldprofile.h:208