Skip to content

Commit

Permalink
fixed failing periodic test
Browse files Browse the repository at this point in the history
  • Loading branch information
fbischoff committed Sep 21, 2023
1 parent 0aea1ca commit 716d963
Show file tree
Hide file tree
Showing 2 changed files with 5 additions and 20 deletions.
3 changes: 2 additions & 1 deletion src/madness/mra/testper.cc
Original file line number Diff line number Diff line change
Expand Up @@ -65,7 +65,8 @@ int test_per(World& world) {
expnt[0] = 10000.0;
coeff[0] = sqrt(expnt[0]/constants::pi);
print(coeff,expnt);
SeparatedConvolution<double,1> op(world, coeff, expnt);
double lo=1.e-6; thresh=1.e-4; // dummy values
SeparatedConvolution<double,1> op(world, coeff, expnt, lo, thresh);

Function<double,1> f = FunctionFactory<double,1>(world).f(constant).initial_level(3).norefine();

Expand Down
22 changes: 3 additions & 19 deletions src/madness/mra/vmra.h
Original file line number Diff line number Diff line change
Expand Up @@ -518,12 +518,9 @@ namespace madness {
return v;
}

auto sv=copy(ovlp);
rr_cholesky(ovlp,tol,piv,rank); // destroys ovlp and gives back Upper ∆ Matrix from CCD
// ovlp zeroed such that input = inner(transpose(output),output).


// rearrange and truncate the functions according to the pivoting of the rr_cholesky
// rearrange and truncate the functions according to the pivoting of the rr_cholesky
std::vector<Function<T,NDIM> > pv(rank);
for(integer i=0;i<rank;++i){
pv[i]=v[piv[i]];
Expand All @@ -534,21 +531,8 @@ namespace madness {
Tensor<T> Linv = inverse(L);
Tensor<T> U = transpose(Linv);

// // L L^T = ovlp
// // Linv ovlp LT inv = 1
// Tensor<T> LTinv=inverse(ovlp);
// Tensor<T> test=inner(LTinv, inner(sv,Linv));
// print("LTinv, inner(sv,Linv)");
// print(test);
//
// print("ovlp - LLt");
// print(sv-inner(L,L,1,0));
// print("Linv ovlp Ltinv",inner(Linv, inner(sv,LTinv)).normf());
//
//
// throw;
World& world=v.front().world();
return transform(world, pv, U,tol);
return transform(world, pv, U);
}

/// convenience routine for orthonromalize_cholesky: orthonromalize_cholesky without information on pivoting and rank
Expand All @@ -572,7 +556,7 @@ namespace madness {
}
// compute overlap
World& world=v.front().world();
Tensor<T> ovlp = matrix_inner(world, v, v,true);
Tensor<T> ovlp = matrix_inner(world, v, v);
return orthonormalize_rrcd(v,ovlp,tol);
}

Expand Down

0 comments on commit 716d963

Please sign in to comment.