mirror of https://github.com/QMCPACK/qmcpack.git
Merge branch 'develop' into extract_min
This commit is contained in:
Ye Luo
GitHub
commit
3dd9ba50c9
Binary file not shown.
|
|
@ -0,0 +1,9 @@
|
|||
# NEXUS User Guide
|
||||
|
||||
This directory contains the LaTeX source for the NEXUS user guide. The
|
||||
PDF guide is not currently created by the build process and must be
|
||||
created manually.
|
||||
|
||||
- A script, build_nexus_user_guide.sh provides the current "best" version of the user guide.
|
||||
- TeX Live 2017 or equivalent, or more recent versions are recommended
|
||||
- Python and the python module pygmentize are required for the source listings.
|
||||
|
|
@ -0,0 +1,24 @@
|
|||
#!/bin/sh
|
||||
if [ -z "$(command -v pdflatex)" ]; then
|
||||
echo "Error: pdflatex is required for build_nexus_user_guide.sh"
|
||||
exit 1
|
||||
fi
|
||||
|
||||
# From QMCPACK build_manual.sh, for eventual versioning
|
||||
#if [ $# -eq 1 ]; then
|
||||
# QMCPACK_VER="$1"
|
||||
# cp -p version.tex version.save.tex
|
||||
# sed -i "s/Development Version/$QMCPACK_VER/" version.tex
|
||||
#fi
|
||||
|
||||
|
||||
pdflatex -shell-escape nexus_user_guide.tex
|
||||
pdflatex -shell-escape nexus_user_guide.tex
|
||||
pdflatex -shell-escape nexus_user_guide.tex
|
||||
|
||||
# From QMCPACK build_manual.sh, for eventual versioning
|
||||
#if [ ! -z "$QMCPACK_VER" ]; then
|
||||
# mv version.save.tex version.tex
|
||||
# echo Added QMCPACK version $QMCPACK_VER
|
||||
#fi
|
||||
exit 0
|
||||
|
|
@ -1,6 +0,0 @@
|
|||
This is makeindex, version 2.15 [TeX Live 2009] (kpathsea + Thai support).
|
||||
Scanning input file project_suite_user_guide.idx....done (7 entries accepted, 0 rejected).
|
||||
Sorting entries....done (22 comparisons).
|
||||
Generating output file project_suite_user_guide.ind....done (23 lines written, 0 warnings).
|
||||
Output written in project_suite_user_guide.ind.
|
||||
Transcript written in project_suite_user_guide.ilg.
|
||||
|
|
@ -1,23 +0,0 @@
|
|||
\begin{theindex}
|
||||
|
||||
\item basic\_qmc, \hyperpage{8}
|
||||
|
||||
\indexspace
|
||||
|
||||
\item generate\_physical\_system, \hyperpage{8}
|
||||
|
||||
\indexspace
|
||||
|
||||
\item load\_analyzer\_image, \hyperpage{8}
|
||||
|
||||
\indexspace
|
||||
|
||||
\item run\_project, \hyperpage{8}
|
||||
|
||||
\indexspace
|
||||
|
||||
\item settings, \hyperpage{8}
|
||||
\item standard\_qmc, \hyperpage{8}
|
||||
\item Structure, \hyperpage{8}
|
||||
|
||||
\end{theindex}
|
||||
Binary file not shown.
|
|
@ -1,54 +0,0 @@
|
|||
\contentsline {chapter}{Contents}{ii}{section*.1}
|
||||
\contentsline {chapter}{\chapternumberline {1}Using this document}{1}{chapter.1}
|
||||
\contentsline {chapter}{\chapternumberline {2}Overview of Nexus}{2}{chapter.2}
|
||||
\contentsline {section}{\numberline {2.1}What Nexus is}{2}{section.2.1}
|
||||
\contentsline {section}{\numberline {2.2}What Nexus can do}{2}{section.2.2}
|
||||
\contentsline {section}{\numberline {2.3}How Nexus is used}{2}{section.2.3}
|
||||
\contentsline {chapter}{\chapternumberline {3}Nexus Installation}{3}{chapter.3}
|
||||
\contentsline {section}{\numberline {3.1}Setting environment variables}{3}{section.3.1}
|
||||
\contentsline {section}{\numberline {3.2}Installing Python dependencies}{3}{section.3.2}
|
||||
\contentsline {chapter}{\chapternumberline {4}Nexus User Scripts}{4}{chapter.4}
|
||||
\contentsline {section}{\numberline {4.1}Nexus imports}{4}{section.4.1}
|
||||
\contentsline {section}{\numberline {4.2}Nexus settings: global state and user-specific information}{5}{section.4.2}
|
||||
\contentsline {section}{\numberline {4.3}Physical system specificaton}{6}{section.4.3}
|
||||
\contentsline {section}{\numberline {4.4}Workflow specification}{7}{section.4.4}
|
||||
\contentsline {subsection}{Generating simulation objects}{7}{section*.2}
|
||||
\contentsline {subsection}{Composing workflows from simulation objects}{11}{section*.3}
|
||||
\contentsline {section}{\numberline {4.5}Workflow execution}{13}{section.4.5}
|
||||
\contentsline {section}{\numberline {4.6}Data analysis}{14}{section.4.6}
|
||||
\contentsline {chapter}{\chapternumberline {5}Complete Examples}{15}{chapter.5}
|
||||
\contentsline {section}{\numberline {5.1}Example 1: Bulk Diamond VMC}{17}{section.5.1}
|
||||
\contentsline {section}{\numberline {5.2}Example 2: Graphene Sheet DMC}{23}{section.5.2}
|
||||
\contentsline {section}{\numberline {5.3}Example 3: C 20 Molecule DMC}{33}{section.5.3}
|
||||
\contentsline {section}{\numberline {5.4}Example 4: Automated oxygen dimer binding curve}{41}{section.5.4}
|
||||
\contentsline {chapter}{\chapternumberline {6}Recommended Reading}{48}{chapter.6}
|
||||
\contentsline {section}{\numberline {6.1}Helpful Links for Installing Python Modules}{48}{section.6.1}
|
||||
\contentsline {section}{\numberline {6.2}Helpful Links for Installing Electronic Structure Codes}{48}{section.6.2}
|
||||
\contentsline {subsection}{PWSCF: pw.x, pw2qmcpack.x, pw2casino.x}{48}{section*.4}
|
||||
\contentsline {subsection}{QMCPACK: qmcapp, qmcapp\_complex, convert4qmc, ppconvert, sqd}{48}{section*.5}
|
||||
\contentsline {subsection}{Wfconvert: wfconvert}{48}{section*.6}
|
||||
\contentsline {subsection}{VASP}{49}{section*.7}
|
||||
\contentsline {subsection}{GAMESS}{49}{section*.8}
|
||||
\contentsline {section}{\numberline {6.3}Brushing Up On Python}{49}{section.6.3}
|
||||
\contentsline {subsection}{Python}{49}{section*.9}
|
||||
\contentsline {subsection}{NumPy}{49}{section*.10}
|
||||
\contentsline {subsection}{Matplotlib}{49}{section*.11}
|
||||
\contentsline {subsection}{Scipy and H5Py}{50}{section*.12}
|
||||
\contentsline {section}{\numberline {6.4}Quantum Monte Carlo: Theory and Practice}{50}{section.6.4}
|
||||
\contentsline {appendix}{\chapternumberline {A}Basic Python constructs}{51}{appendix.A}
|
||||
\contentsline {section}{\numberline {A.1}Intrinsic types: \texttt {int, float, str}}{51}{section.A.1}
|
||||
\contentsline {section}{\numberline {A.2}Container types: \texttt {tuple, list, array, dict, obj}}{51}{section.A.2}
|
||||
\contentsline {section}{\numberline {A.3}Conditional Statements: \texttt {if/elif/else}}{53}{section.A.3}
|
||||
\contentsline {section}{\numberline {A.4}Iteration: \texttt {for}}{54}{section.A.4}
|
||||
\contentsline {section}{\numberline {A.5}Functions: \texttt {def}, argument syntax}{55}{section.A.5}
|
||||
\contentsline {appendix}{\chapternumberline {B}QMC Practice in a Nutshell}{57}{appendix.B}
|
||||
\contentsline {section}{\numberline {B.1}VMC and DMC in the abstract}{57}{section.B.1}
|
||||
\contentsline {section}{\numberline {B.2}From expectation values to random walks}{58}{section.B.2}
|
||||
\contentsline {section}{\numberline {B.3}Quality orbitals: planewaves, cutoffs, splines, and meshes}{58}{section.B.3}
|
||||
\contentsline {section}{\numberline {B.4}Quality Jastrows: less variance = more efficient}{59}{section.B.4}
|
||||
\contentsline {section}{\numberline {B.5}Finite size effects: k-points, supercells, and corrections}{60}{section.B.5}
|
||||
\contentsline {section}{\numberline {B.6}Imaginary time discretization: the DMC timestep}{60}{section.B.6}
|
||||
\contentsline {section}{\numberline {B.7}Population control bias: safety in numbers}{61}{section.B.7}
|
||||
\contentsline {section}{\numberline {B.8}The fixed node/phase approximation: varying the nodes/phase}{62}{section.B.8}
|
||||
\contentsline {section}{\numberline {B.9}Pseudopotentials: theoretical dissonance, the locality approximation, and T-moves}{62}{section.B.9}
|
||||
\contentsline {section}{\numberline {B.10}Other approximations: what else is missing?}{63}{section.B.10}
|
||||
|
|
@ -66,6 +66,36 @@ public:
|
|||
nElms=nIntPnts;
|
||||
}
|
||||
|
||||
void setPhiAndEta(SPOSetBasePtr Phi, SPOSetBasePtr Eta)
|
||||
{
|
||||
Psi1 = Phi;
|
||||
Psi2 = Eta;
|
||||
int norb = Psi1->OrbitalSetSize;
|
||||
|
||||
val1.resize(norb);
|
||||
grad1.resize(norb);
|
||||
lapl1.resize(norb);
|
||||
|
||||
val2.resize(norb);
|
||||
grad2.resize(norb);
|
||||
lapl2.resize(norb);
|
||||
}
|
||||
|
||||
void allocateELspace()
|
||||
{
|
||||
ELideal.resize(nElms);
|
||||
ELorig.resize(nElms);
|
||||
ELcurr.resize(nElms);
|
||||
pos.resize(nElms);
|
||||
}
|
||||
|
||||
void computeValAtZero()
|
||||
{
|
||||
TinyVector<RealType,3> ddr=0;
|
||||
evaluate(Psi1,ddr,val1,grad1,lapl1);
|
||||
valAtZero = val1[0];
|
||||
}
|
||||
|
||||
~CuspCorr() {}
|
||||
|
||||
/*
|
||||
|
|
@ -285,6 +315,7 @@ private:
|
|||
|
||||
SPOSetBasePtr Psi1,Psi2;
|
||||
|
||||
public:
|
||||
// cutoff
|
||||
RealType beta0,DX,eta0, ELorigAtRc;
|
||||
RealType Rc_init,Rc,C,sg,Z,valAtZero,valAtRc,Rc_max;
|
||||
|
|
|
|||
|
|
@ -36,8 +36,8 @@ FOREACH(fname ${FILES_TO_COPY})
|
|||
EXECUTE_PROCESS(COMMAND ${CMAKE_COMMAND} -E copy "${CMAKE_CURRENT_SOURCE_DIR}/${fname}" ${UTEST_DIR})
|
||||
ENDFOREACH()
|
||||
|
||||
IF ((NOT QMC_COMPLEX) AND (NOT QMC_MIXED_PRECISION))
|
||||
SET(MO_SRCS test_MO.cpp)
|
||||
IF ((NOT QMC_COMPLEX) AND (NOT QMC_MIXED_PRECISION) AND (NOT QMC_CUDA))
|
||||
SET(MO_SRCS test_MO.cpp test_cusp_corr.cpp)
|
||||
ENDIF()
|
||||
|
||||
ADD_EXECUTABLE(${UTEST_EXE} test_wf.cpp test_bspline_jastrow.cpp test_einset.cpp test_pw.cpp
|
||||
|
|
|
|||
|
|
@ -0,0 +1,231 @@
|
|||
from __future__ import print_function
|
||||
|
||||
# Cusp corrections for gaussian orbitals
|
||||
|
||||
# From "Scheme for adding electron-nucleus cusps to Gaussian orbitals" A. Ma, D. Towler, N. D. Drummond, R. J. Needs, Journal of Chemical Physics 122, 224322(2005) https://doi.org/10.1063/1.1940588
|
||||
|
||||
# Also qmc_algorithms/Wavefunctions/CuspCorrection.ipynb
|
||||
|
||||
|
||||
from sympy import *
|
||||
import gaussian_orbitals
|
||||
import read_qmcpack
|
||||
|
||||
|
||||
|
||||
alpha = IndexedBase('alpha')
|
||||
rc = Symbol('r_c')
|
||||
X1,X2,X3,X4,X5 = symbols('X_1 X_2 X_3 X_4 X_5')
|
||||
r = Symbol('r')
|
||||
Zeff = Symbol('Z_eff')
|
||||
|
||||
|
||||
|
||||
|
||||
def solve_for_alpha():
|
||||
p = alpha[0] + alpha[1]*r + alpha[2]*r**2 + alpha[3]*r**3 + alpha[4]*r**4
|
||||
|
||||
|
||||
# Constraint equations
|
||||
# Value matches at r_c
|
||||
eq1 = Eq(p.subs(r,rc), X1)
|
||||
# 1st derivative matches at r_c
|
||||
eq2 = Eq(diff(p,r).subs(r,rc), X2)
|
||||
# 2nd derivative matches at r_c
|
||||
eq3 = Eq((diff(p,r,2)+diff(p,r)**2).subs(r,rc),X3)
|
||||
# Cusp condition
|
||||
eq4 = Eq(diff(p,r).subs(r,0),X4)
|
||||
# Value of phi tilde at r=0
|
||||
eq5 = Eq(p.subs(r,0),X5)
|
||||
|
||||
sln = solve([eq1, eq2, eq3, eq4, eq5],[alpha[0], alpha[1], alpha[2], alpha[3], alpha[4]])
|
||||
|
||||
print('Symbolic solution for alpha:')
|
||||
for i,s in enumerate(sln[0]):
|
||||
print(alpha[i],' = ',s)
|
||||
print()
|
||||
|
||||
return sln[0]
|
||||
|
||||
alpha_sln = solve_for_alpha()
|
||||
|
||||
|
||||
def solve_for_alpha(Xvals, rc_val):
|
||||
svals = {rc: rc_val, X1:Xvals[1], X2:Xvals[2], X3:Xvals[3], X4:Xvals[4], X5:Xvals[5]}
|
||||
alpha_vals = []
|
||||
for i,s in enumerate(alpha_sln):
|
||||
alpha_vals.append(s.subs(svals))
|
||||
return alpha_vals
|
||||
|
||||
def simple_X_vals():
|
||||
rc_val = 0.1
|
||||
Xvals = [0.0, 1.0, 2.0, 3.0, 4.0, 5.0] # first entry is unused
|
||||
|
||||
print('For X: ',Xvals[1:])
|
||||
alpha_vals = solve_for_alpha(Xvals, rc_val)
|
||||
for i,a in enumerate(alpha_vals):
|
||||
print(alpha[i], ' = ', a)
|
||||
#svals = {rc: rc_val, X1:Xvals[1], X2:Xvals[2], X3:Xvals[3], X4:Xvals[4], X5:Xvals[5]}
|
||||
#for i,s in enumerate(alpha_sln):
|
||||
# print(alpha[i], ' = ', s.subs(svals))
|
||||
|
||||
def evalX(gto, Z_val, rc_val):
|
||||
Xvals = [0.0]*6
|
||||
val_rc, grad_rc, lap_rc = [v[0] for v in gto.eval_vgl(rc_val, 0.0, 0.0)]
|
||||
val_zero, grad_zero, lap_zero = [v[0] for v in gto.eval_vgl(0.0, 0.0, 0.0)]
|
||||
Xvals[1] = log(val_rc)
|
||||
Xvals[2] = grad_rc[0]/val_rc
|
||||
Xvals[3] = (lap_rc - 2.0*grad_rc[0]/rc_val)/val_rc
|
||||
Xvals[4] = -Z_val
|
||||
Xvals[5] = log(abs(val_zero)) # initially use phi at 0
|
||||
return Xvals
|
||||
|
||||
|
||||
|
||||
return he_gto, Xvals, alpha_vals
|
||||
|
||||
def output_required_xvals_and_alphas(Xvals, alpha_vals):
|
||||
print('Xvals = ',Xvals)
|
||||
print(' // From gen_cusp_corr.py')
|
||||
for i in range(5):
|
||||
print(' REQUIRE(X[%d] == Approx(%.15f));'%(i,Xvals[i+1]))
|
||||
|
||||
print()
|
||||
print(' // From gen_cusp_corr.py')
|
||||
for i in range(5):
|
||||
print(' REQUIRE(cusp.alpha[%d] == Approx(%.15f));'%(i,alpha_vals[i]))
|
||||
|
||||
def output_array(v, name):
|
||||
print(' // From gen_cusp_corr.py')
|
||||
for i,a in enumerate(v):
|
||||
print(' REQUIRE(%s[%d] == Approx(%.15f));'%(name, i, a))
|
||||
print()
|
||||
|
||||
|
||||
def del_spherical(e, r):
|
||||
"""Compute Laplacian for expression e with respect to symbol r.
|
||||
Currently works only with radial dependence"""
|
||||
t1 = r*r*diff(e, r)
|
||||
t2 = diff(t1, r)/(r*r)
|
||||
return simplify(t2)
|
||||
|
||||
def get_symbolic_effective_local_energy():
|
||||
p = alpha[0] + alpha[1]*r + alpha[2]*r**2 + alpha[3]*r**3 + alpha[4]*r**4
|
||||
R_sym = exp(p)
|
||||
phi_tilde = R_sym
|
||||
effEl_sym = -S.Half * del_spherical(phi_tilde, r)/phi_tilde - Zeff/r
|
||||
return effEl_sym
|
||||
|
||||
|
||||
def eval_El(El_sym, r_val, Zeff_val, alpha_vals):
|
||||
slist = {alpha[0]:alpha_vals[0], alpha[1]:alpha_vals[1], alpha[2]: alpha_vals[2], alpha[3]:alpha_vals[3],
|
||||
alpha[4]:alpha_vals[4], Zeff:Zeff_val, r:r_val}
|
||||
val = El_sym.subs(slist).evalf()
|
||||
return val
|
||||
|
||||
|
||||
def get_grid():
|
||||
pos = []
|
||||
for i in range(10):
|
||||
rval = .012*(i+1)
|
||||
pos.append(rval)
|
||||
return pos
|
||||
|
||||
def get_original_local_energy(pos, gto, El_sym, alpha_vals, rc_val, Zeff_val):
|
||||
vals = []
|
||||
for rval in pos:
|
||||
val, grad, lap = gto.eval_vgl(rval, 0.0, 0.0)
|
||||
real_el = -.5*lap[0]/val[0] - Zeff_val/rval
|
||||
vals.append(real_el)
|
||||
|
||||
return vals
|
||||
|
||||
|
||||
def get_current_local_energy(pos, gto, El_sym, alpha_vals, rc_val, dE, Zeff_val):
|
||||
vals = []
|
||||
for rval in pos:
|
||||
el = eval_El(El_sym, rval, 2.0, alpha_vals)
|
||||
if rval < rc_val:
|
||||
vals.append(el + dE)
|
||||
else:
|
||||
val, grad, lap = gto.eval_vgl(rval, 0.0, 0.0)
|
||||
real_el = -.5*lap[0]/val[0] - Zeff_val/rval
|
||||
vals.append(real_el + dE)
|
||||
|
||||
return vals
|
||||
|
||||
|
||||
|
||||
def get_ideal_local_energy(pos, rc_val, beta0_val=0.0):
|
||||
ideal_EL = []
|
||||
beta0 = Symbol('beta_0')
|
||||
beta_vals = [beta0, 3.25819, -15.0126, 33.7308, -42.8705, 31.2276, -12.1316, 1.94692]
|
||||
El_terms = [beta_vals[n]*r**(n+1) for n in range(1,8)]
|
||||
Z = Symbol('Z')
|
||||
Z_val = 2.0
|
||||
EL_ideal_sym = Z*Z*(beta0 + sum(El_terms))
|
||||
|
||||
slist = {beta0: 0.0, Z:Z_val, r: rc_val}
|
||||
|
||||
idealEL_at_rc = EL_ideal_sym.subs(slist).evalf()
|
||||
#print('idealEL at rc = ',idealEL_at_rc)
|
||||
beta0_val = (-idealEL_at_rc)/Z_val/Z_val
|
||||
#print('beta0_val = ',beta0_val)
|
||||
|
||||
for rval in pos:
|
||||
slist = {beta0: beta0_val, Z:Z_val, r: rval}
|
||||
v = EL_ideal_sym.subs(slist).evalf()
|
||||
ideal_EL.append(v)
|
||||
|
||||
return ideal_EL
|
||||
|
||||
|
||||
def values_for_He():
|
||||
rc_val = 0.1
|
||||
Z_val = 2
|
||||
|
||||
basis_set,MO = read_qmcpack.parse_qmc_wf('he_sto3g.wfj.xml')
|
||||
he_gto = gaussian_orbitals.GTO(basis_set['He'])
|
||||
|
||||
Xvals = evalX(he_gto, Z_val, rc_val)
|
||||
|
||||
alpha_vals = solve_for_alpha(Xvals, rc_val)
|
||||
output_required_xvals_and_alphas(Xvals, alpha_vals)
|
||||
|
||||
El_sym = get_symbolic_effective_local_energy()
|
||||
print("El_sym = ",El_sym)
|
||||
|
||||
Zeff_val = 2.0
|
||||
el_at_rc = -eval_El(El_sym, rc_val, Zeff_val, alpha_vals)
|
||||
dE = el_at_rc
|
||||
print('el at rc_val = ',el_at_rc)
|
||||
|
||||
|
||||
pos = get_grid()
|
||||
|
||||
current_EL = get_current_local_energy(pos, he_gto, El_sym, alpha_vals, rc_val, dE, Zeff_val)
|
||||
original_EL = get_original_local_energy(pos, he_gto, El_sym, alpha_vals, rc_val, Zeff_val)
|
||||
#print('Current effective local energy')
|
||||
#for (p,v) in zip(pos,current_EL):
|
||||
# print(p,v)
|
||||
print(" // Grid for local energy evaluations")
|
||||
output_array(pos, "cusp.pos");
|
||||
|
||||
print(" // Original local energy")
|
||||
output_array(original_EL, "cusp.ELorig")
|
||||
|
||||
print(" // Current local energy")
|
||||
output_array(current_EL, "cusp.ELcurr")
|
||||
|
||||
print(" // Ideal local energy")
|
||||
ideal_EL = get_ideal_local_energy(pos,rc_val)
|
||||
output_array(ideal_EL, "cusp.ELideal")
|
||||
|
||||
chi2 = 0.0
|
||||
for rval, ideal, curr in zip(pos, ideal_EL, current_EL):
|
||||
chi2 += (ideal - curr)**2
|
||||
print(' REQUIRE(chi2 == Approx(%.10f)); '%chi2)
|
||||
|
||||
if __name__ == '__main__':
|
||||
#simple_X_vals()
|
||||
values_for_He()
|
||||
|
|
@ -0,0 +1,235 @@
|
|||
|
||||
//////////////////////////////////////////////////////////////////////////////////////
|
||||
// This file is distributed under the University of Illinois/NCSA Open Source License.
|
||||
// See LICENSE file in top directory for details.
|
||||
//
|
||||
// Copyright (c) 2018 Jeongnim Kim and QMCPACK developers.
|
||||
//
|
||||
// File developed by: Mark Dewing, mdewing@anl.gov, Argonne National Laboratory
|
||||
//
|
||||
// File created by: Mark Dewing, mdewing@anl.gov, Argonne National Laboratory
|
||||
//////////////////////////////////////////////////////////////////////////////////////
|
||||
|
||||
|
||||
#include "catch.hpp"
|
||||
|
||||
#include "Configuration.h"
|
||||
|
||||
#include "Message/Communicate.h"
|
||||
#include "Numerics/OneDimGridBase.h"
|
||||
#include "Particle/DistanceTableData.h"
|
||||
#include "ParticleIO/XMLParticleIO.h"
|
||||
#include "Numerics/GaussianBasisSet.h"
|
||||
|
||||
#include "QMCWaveFunctions/MolecularOrbitals/LocalizedBasisSet.h"
|
||||
#include "QMCWaveFunctions/MolecularOrbitals/SphericalBasisSet.h"
|
||||
#include "QMCWaveFunctions/MolecularOrbitals/NGOBuilder.h"
|
||||
#include "QMCWaveFunctions/MolecularOrbitals/CuspCorr.h"
|
||||
|
||||
#include "QMCWaveFunctions/SPOSetBuilderFactory.h"
|
||||
|
||||
namespace qmcplusplus
|
||||
{
|
||||
|
||||
TEST_CASE("CuspCorrection He", "[wavefunction]")
|
||||
{
|
||||
OHMMS::Controller->initialize(0, NULL);
|
||||
Communicate *c = OHMMS::Controller;
|
||||
|
||||
ParticleSet elec;
|
||||
std::vector<int> agroup(2);
|
||||
agroup[0] = 1;
|
||||
agroup[1] = 1;
|
||||
elec.setName("e");
|
||||
elec.create(agroup);
|
||||
elec.R[0] = 0.0;
|
||||
|
||||
SpeciesSet &tspecies = elec.getSpeciesSet();
|
||||
int upIdx = tspecies.addSpecies("u");
|
||||
int downIdx = tspecies.addSpecies("d");
|
||||
int massIdx = tspecies.addAttribute("mass");
|
||||
tspecies(massIdx, upIdx) = 1.0;
|
||||
tspecies(massIdx, downIdx) = 1.0;
|
||||
|
||||
ParticleSet ions;
|
||||
ions.setName("ion0");
|
||||
ions.create(1);
|
||||
ions.R[0] = 0.0;
|
||||
SpeciesSet &ispecies = ions.getSpeciesSet();
|
||||
int heIdx = ispecies.addSpecies("He");
|
||||
ions.update();
|
||||
|
||||
|
||||
#ifdef ENABLE_SOA
|
||||
elec.addTable(ions,DT_SOA);
|
||||
#else
|
||||
elec.addTable(ions,DT_AOS);
|
||||
#endif
|
||||
elec.update();
|
||||
|
||||
|
||||
Libxml2Document doc;
|
||||
bool okay = doc.parse("he_sto3g.wfj.xml");
|
||||
REQUIRE(okay);
|
||||
xmlNodePtr root = doc.getRoot();
|
||||
|
||||
TrialWaveFunction psi(c);
|
||||
|
||||
OrbitalBuilderBase::PtclPoolType particle_set_map;
|
||||
particle_set_map["e"] = &elec;
|
||||
particle_set_map["ion0"] = &ions;
|
||||
|
||||
|
||||
SPOSetBuilderFactory bf(elec, psi, particle_set_map);
|
||||
|
||||
OhmmsXPathObject MO_base("//determinantset", doc.getXPathContext());
|
||||
REQUIRE(MO_base.size() == 1);
|
||||
|
||||
SPOSetBuilder *bb = bf.createSPOSetBuilder(MO_base[0]);
|
||||
REQUIRE(bb != NULL);
|
||||
|
||||
OhmmsXPathObject slater_base("//determinant", doc.getXPathContext());
|
||||
bb->loadBasisSetFromXML(MO_base[0]);
|
||||
SPOSetBase *sposet = bb->createSPOSet(slater_base[0]);
|
||||
|
||||
|
||||
typedef OneDimGridBase<double> GridType;
|
||||
typedef LCOrbitalSet<LocalizedBasisSet<SphericalBasisSet<NGOrbital, GridType>>, false> OrbType;
|
||||
OrbType *lcob = dynamic_cast<OrbType *>(sposet);
|
||||
REQUIRE(lcob != NULL);
|
||||
|
||||
|
||||
typedef CuspCorr<LocalizedBasisSet<SphericalBasisSet<NGOrbital, GridType>>> CuspCorrType;
|
||||
typedef OrbitalSetTraits<double>::ValueVector_t ValueVector_t;
|
||||
|
||||
double rc = 0.1;
|
||||
int npts = 10;
|
||||
CuspCorrType cusp(rc, npts, &elec, &ions);
|
||||
|
||||
typedef NGOBuilder::CenteredOrbitalType COT;
|
||||
OrbType bs_phi(lcob->myBasisSet);
|
||||
bs_phi.setOrbitalSetSize(lcob->OrbitalSetSize);
|
||||
bs_phi.BasisSetSize = lcob->BasisSetSize;
|
||||
bs_phi.setIdentity(false);
|
||||
|
||||
*(bs_phi.C) = *(lcob->C);
|
||||
|
||||
OrbType bs_eta(lcob->myBasisSet);
|
||||
bs_eta.setOrbitalSetSize(lcob->OrbitalSetSize);
|
||||
bs_eta.BasisSetSize = lcob->BasisSetSize;
|
||||
bs_eta.setIdentity(false);
|
||||
*(bs_eta.C) = *(lcob->C);
|
||||
// For He in a minimal basis, there are only s-type orbitals
|
||||
(*bs_eta.C)(0,0) = 0.0;
|
||||
|
||||
cusp.setPhiAndEta(&bs_phi, &bs_eta);
|
||||
|
||||
cusp.curCenter = 0;
|
||||
cusp.curOrb = 0;
|
||||
cusp.Z = 2.0;
|
||||
cusp.computeValAtZero();
|
||||
|
||||
ValueVector_t X(5);
|
||||
cusp.evalX(X);
|
||||
|
||||
//std::cout << "X = " << X << std::endl;
|
||||
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(X[0] == Approx(-0.033436891110336));
|
||||
REQUIRE(X[1] == Approx(-0.653568722769692));
|
||||
REQUIRE(X[2] == Approx(-5.819488164002633));
|
||||
REQUIRE(X[3] == Approx(-2.000000000000000));
|
||||
REQUIRE(X[4] == Approx(-0.000396345019839));
|
||||
|
||||
cusp.X2alpha(X);
|
||||
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(cusp.alpha[0] == Approx(-0.000396345019839));
|
||||
REQUIRE(cusp.alpha[1] == Approx(-2.000000000000000));
|
||||
REQUIRE(cusp.alpha[2] == Approx(56.659413909100188));
|
||||
REQUIRE(cusp.alpha[3] == Approx(-599.993590267020409));
|
||||
REQUIRE(cusp.alpha[4] == Approx(2003.589050855219512));
|
||||
|
||||
|
||||
cusp.allocateELspace();
|
||||
|
||||
// compute original EL
|
||||
cusp.getELorig();
|
||||
|
||||
// compute current EL
|
||||
cusp.getELcurr();
|
||||
|
||||
// compute ideal EL
|
||||
cusp.getELideal(cusp.ELorigAtRc);
|
||||
|
||||
// Grid for local energy evaluations
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(cusp.pos[0] == Approx(0.012000000000000));
|
||||
REQUIRE(cusp.pos[1] == Approx(0.024000000000000));
|
||||
REQUIRE(cusp.pos[2] == Approx(0.036000000000000));
|
||||
REQUIRE(cusp.pos[3] == Approx(0.048000000000000));
|
||||
REQUIRE(cusp.pos[4] == Approx(0.060000000000000));
|
||||
REQUIRE(cusp.pos[5] == Approx(0.072000000000000));
|
||||
REQUIRE(cusp.pos[6] == Approx(0.084000000000000));
|
||||
REQUIRE(cusp.pos[7] == Approx(0.096000000000000));
|
||||
REQUIRE(cusp.pos[8] == Approx(0.108000000000000));
|
||||
REQUIRE(cusp.pos[9] == Approx(0.120000000000000));
|
||||
|
||||
// Original local energy
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(cusp.ELorig[0] == Approx(-156.654088753559449));
|
||||
REQUIRE(cusp.ELorig[1] == Approx(-73.346068180623860));
|
||||
REQUIRE(cusp.ELorig[2] == Approx(-45.610385939854496));
|
||||
REQUIRE(cusp.ELorig[3] == Approx(-31.780236703094037));
|
||||
REQUIRE(cusp.ELorig[4] == Approx(-23.522092887496903));
|
||||
REQUIRE(cusp.ELorig[5] == Approx(-18.057926774366479));
|
||||
REQUIRE(cusp.ELorig[6] == Approx(-14.196956436578184));
|
||||
REQUIRE(cusp.ELorig[7] == Approx(-11.343582162638119));
|
||||
REQUIRE(cusp.ELorig[8] == Approx(-9.166698588588746));
|
||||
REQUIRE(cusp.ELorig[9] == Approx(-7.467419605354912));
|
||||
|
||||
|
||||
// Current local energy
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(cusp.ELcurr[0] == Approx(-119.501377810849448));
|
||||
REQUIRE(cusp.ELcurr[1] == Approx(-84.586558430353278));
|
||||
REQUIRE(cusp.ELcurr[2] == Approx(-55.798846293994899));
|
||||
REQUIRE(cusp.ELcurr[3] == Approx(-32.804136849327627));
|
||||
REQUIRE(cusp.ELcurr[4] == Approx(-15.556451408317010));
|
||||
REQUIRE(cusp.ELcurr[5] == Approx(-4.108843171781601));
|
||||
REQUIRE(cusp.ELcurr[6] == Approx(1.506652537070609));
|
||||
REQUIRE(cusp.ELcurr[7] == Approx(1.330023830008599));
|
||||
REQUIRE(cusp.ELcurr[8] == Approx(1.387870101712975));
|
||||
REQUIRE(cusp.ELcurr[9] == Approx(3.087149084946809));
|
||||
|
||||
// Ideal local energy
|
||||
// From gen_cusp_corr.py
|
||||
REQUIRE(cusp.ELideal[0] == Approx(-0.080399129819489));
|
||||
REQUIRE(cusp.ELideal[1] == Approx(-0.075454680124444));
|
||||
REQUIRE(cusp.ELideal[2] == Approx(-0.067869571712510));
|
||||
REQUIRE(cusp.ELideal[3] == Approx(-0.058114416794904));
|
||||
REQUIRE(cusp.ELideal[4] == Approx(-0.046606832483210));
|
||||
REQUIRE(cusp.ELideal[5] == Approx(-0.033715675742608));
|
||||
REQUIRE(cusp.ELideal[6] == Approx(-0.019765043739301));
|
||||
REQUIRE(cusp.ELideal[7] == Approx(-0.005038047738043));
|
||||
REQUIRE(cusp.ELideal[8] == Approx(0.010219631429280));
|
||||
REQUIRE(cusp.ELideal[9] == Approx(0.025796398438142));
|
||||
|
||||
double chi2 = cusp.getchi2();
|
||||
REQUIRE(chi2 == Approx(25854.2846426019));
|
||||
|
||||
|
||||
double data[9];
|
||||
Vector<double> xgrid;
|
||||
Vector<double> rad_orb;
|
||||
xgrid.resize(1);
|
||||
xgrid[0] = 0.012;
|
||||
rad_orb.resize(1);
|
||||
|
||||
cusp.execute(0, 0, 2.0, &bs_phi, &bs_eta, xgrid, rad_orb, "none", rc, data);
|
||||
|
||||
SPOSetBuilderFactory::clear();
|
||||
|
||||
}
|
||||
|
||||
}
|
||||
Loading…
Reference in New Issue