scnc
/
qmcpack
mirror of https://github.com/QMCPACK/qmcpack.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

clean up nofk input block, add nexus support

This commit is contained in:
Jaron Krogel 2017-09-08 14:22:56 -04:00
parent 8c006b2841
commit fb65b4e84f
2 changed files with 96 additions and 102 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
nexus/library/qmcpack_input.py
View File

@ -2205,6 +2205,13 @@ class flux(QIxml):
attributes = ['type','name']
identifier = 'name'
#end class flux
class momentum(QIxml):
tag = 'estimator'
attributes = ['type','name','grid','samples','hdf5','wavefunction']
identifier = 'name'
write_types = obj(hdf5=yesno)
#end class momentum
estimator = QIxmlFactory(
name = 'estimator',
@ -2226,6 +2233,7 @@ estimator = QIxmlFactory(
skall = skall,
gofr = gofr,
flux = flux,
momentum = momentum,
),
typekey = 'type',
typekey2 = 'name'
@ -2457,7 +2465,7 @@ classes = [ #standard classes
header,local,force,forwardwalking,observable,record,rmc,pressure,dmccorrection,
nofk,mpc_est,flux,distancetable,cpp,element,spline,setparams,
backflow,transformation,cubicgrid,molecular_orbital_builder,cmc,sk,skall,gofr,
host,date,user,rpa_jastrow
host,date,user,rpa_jastrow,momentum
]
types = dict( #simple types and factories
#host = param,
@ -2641,6 +2649,9 @@ force.defaults.set(
pressure.defaults.set(
type='Pressure'
)
momentum.defaults.set(
type='momentum'
)
linear.defaults.set(

185
src/QMCHamiltonians/MomentumEstimator.cpp
View File

@ -166,124 +166,107 @@ bool MomentumEstimator::putSpecial(xmlNodePtr cur, ParticleSet& elns, bool rootN
OhmmsAttributeSet pAttrib;
std::string hdf5_flag="yes";
pAttrib.add(hdf5_flag,"hdf5");
pAttrib.add(kgrid,"grid");
pAttrib.add(M,"samples");
pAttrib.put(cur);
hdf5_out = (hdf5_flag=="yes");
// app_log()<<" MomentumEstimator::putSpecial "<< std::endl;
xmlNodePtr kids=cur->children;
while (kids!=NULL)
{
std::string cname((const char*)(kids->name));
// app_log()<<" MomentumEstimator::cname : "<<cname<< std::endl;
if (cname=="kpoints")
{
std::string ctype("manual");
OhmmsAttributeSet pAttrib;
pAttrib.add(ctype,"mode");
pAttrib.add(kgrid,"grid");
pAttrib.put(kids);
#if OHMMS_DIM==3
int numqtwists(6*kgrid+3);
std::vector<int> qk(0);
mappedQtonofK.resize(numqtwists,qk);
compQ.resize(numqtwists);
RealType qn(4.0*M_PI*M_PI*std::pow(Lattice.Volume,-2.0/3.0));
mappedQnorms.resize(numqtwists,qn*0.5/RealType(M));
if (twist[0]==0)
mappedQnorms[kgrid]=qn/RealType(M);
if (twist[1]==0)
mappedQnorms[3*kgrid+1]=qn/RealType(M);
if (twist[2]==0)
mappedQnorms[5*kgrid+2]=qn/RealType(M);
// app_log()<<" Jnorm="<<qn<< std::endl;
Q.resize(numqtwists);
for (int i=-kgrid; i<(kgrid+1); i++)
int numqtwists(6*kgrid+3);
std::vector<int> qk(0);
mappedQtonofK.resize(numqtwists,qk);
compQ.resize(numqtwists);
RealType qn(4.0*M_PI*M_PI*std::pow(Lattice.Volume,-2.0/3.0));
mappedQnorms.resize(numqtwists,qn*0.5/RealType(M));
if (twist[0]==0)
mappedQnorms[kgrid]=qn/RealType(M);
if (twist[1]==0)
mappedQnorms[3*kgrid+1]=qn/RealType(M);
if (twist[2]==0)
mappedQnorms[5*kgrid+2]=qn/RealType(M);
Q.resize(numqtwists);
for (int i=-kgrid; i<(kgrid+1); i++)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=i-twist[1];
kpt[2]=i-twist[2];
kpt=Lattice.k_cart(kpt);
Q[i+kgrid]=std::abs(kpt[0]);
Q[i+kgrid+(2*kgrid+1)]=std::abs(kpt[1]);
Q[i+kgrid+(4*kgrid+2)]=std::abs(kpt[2]);
}
app_log()<<" Using all k-space points with (nx^2+ny^2+nz^2)^0.5 < "<< kgrid <<" for Momentum Distribution."<< std::endl;
app_log()<<" My twist is:"<<twist[0]<<" "<<twist[1]<<" "<<twist[2]<< std::endl;
int indx(0);
int kgrid_squared=kgrid*kgrid;
for (int i=-kgrid; i<(kgrid+1); i++)
{
for (int j=-kgrid; j<(kgrid+1); j++)
{
for (int k=-kgrid; k<(kgrid+1); k++)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=i-twist[1];
kpt[2]=i-twist[2];
kpt=Lattice.k_cart(kpt);
Q[i+kgrid]=std::abs(kpt[0]);
Q[i+kgrid+(2*kgrid+1)]=std::abs(kpt[1]);
Q[i+kgrid+(4*kgrid+2)]=std::abs(kpt[2]);
}
app_log()<<" Using all k-space points with (nx^2+ny^2+nz^2)^0.5 < "<< kgrid <<" for Momentum Distribution."<< std::endl;
app_log()<<" My twist is:"<<twist[0]<<" "<<twist[1]<<" "<<twist[2]<< std::endl;
int indx(0);
int kgrid_squared=kgrid*kgrid;
for (int i=-kgrid; i<(kgrid+1); i++)
{
for (int j=-kgrid; j<(kgrid+1); j++)
if (i*i+j*j+k*k<=kgrid_squared) //if (std::sqrt(i*i+j*j+k*k)<=kgrid)
{
for (int k=-kgrid; k<(kgrid+1); k++)
{
if (i*i+j*j+k*k<=kgrid_squared) //if (std::sqrt(i*i+j*j+k*k)<=kgrid)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=j-twist[1];
kpt[2]=k-twist[2];
//convert to Cartesian: note that 2Pi is multiplied
kpt=Lattice.k_cart(kpt);
kPoints.push_back(kpt);
mappedQtonofK[i+kgrid].push_back(indx);
mappedQtonofK[j+kgrid+(2*kgrid+1)].push_back(indx);
mappedQtonofK[k+kgrid+(4*kgrid+2)].push_back(indx);
indx++;
}
}
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=j-twist[1];
kpt[2]=k-twist[2];
//convert to Cartesian: note that 2Pi is multiplied
kpt=Lattice.k_cart(kpt);
kPoints.push_back(kpt);
mappedQtonofK[i+kgrid].push_back(indx);
mappedQtonofK[j+kgrid+(2*kgrid+1)].push_back(indx);
mappedQtonofK[k+kgrid+(4*kgrid+2)].push_back(indx);
indx++;
}
}
}
}
#endif
#if OHMMS_DIM==2
int numqtwists(4*kgrid+2);
std::vector<int> qk(0);
mappedQtonofK.resize(numqtwists,qk);
compQ.resize(numqtwists);
RealType qn(2.0*M_PI/std::sqrt(Lattice.Volume));
mappedQnorms.resize(numqtwists,qn*0.5/RealType(M));
if (twist[0]==0)
mappedQnorms[kgrid]=qn/RealType(M);
if (twist[1]==0)
mappedQnorms[3*kgrid+1]=qn/RealType(M);
// app_log()<<" Jnorm="<<qn<< std::endl;
Q.resize(numqtwists);
for (int i=-kgrid; i<(kgrid+1); i++)
int numqtwists(4*kgrid+2);
std::vector<int> qk(0);
mappedQtonofK.resize(numqtwists,qk);
compQ.resize(numqtwists);
RealType qn(2.0*M_PI/std::sqrt(Lattice.Volume));
mappedQnorms.resize(numqtwists,qn*0.5/RealType(M));
if (twist[0]==0)
mappedQnorms[kgrid]=qn/RealType(M);
if (twist[1]==0)
mappedQnorms[3*kgrid+1]=qn/RealType(M);
Q.resize(numqtwists);
for (int i=-kgrid; i<(kgrid+1); i++)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=i-twist[1];
kpt=Lattice.k_cart(kpt);
Q[i+kgrid]=std::abs(kpt[0]);
Q[i+kgrid+(2*kgrid+1)]=std::abs(kpt[1]);
}
app_log()<<" Using all k-space points with (nx^2+ny^2)^0.5 < "<< kgrid <<" for Momentum Distribution."<< std::endl;
app_log()<<" My twist is:"<<twist[0]<<" "<<twist[1]<< std::endl;
int indx(0);
int kgrid_squared=kgrid*kgrid;
for (int i=-kgrid; i<(kgrid+1); i++)
{
for (int j=-kgrid; j<(kgrid+1); j++)
{
if (i*i+j*j<=kgrid_squared) //if (std::sqrt(i*i+j*j+k*k)<=kgrid)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=i-twist[1];
kpt[1]=j-twist[1];
//convert to Cartesian: note that 2Pi is multiplied
kpt=Lattice.k_cart(kpt);
Q[i+kgrid]=std::abs(kpt[0]);
Q[i+kgrid+(2*kgrid+1)]=std::abs(kpt[1]);
kPoints.push_back(kpt);
mappedQtonofK[i+kgrid].push_back(indx);
mappedQtonofK[j+kgrid+(2*kgrid+1)].push_back(indx);
indx++;
}
app_log()<<" Using all k-space points with (nx^2+ny^2)^0.5 < "<< kgrid <<" for Momentum Distribution."<< std::endl;
app_log()<<" My twist is:"<<twist[0]<<" "<<twist[1]<< std::endl;
int indx(0);
int kgrid_squared=kgrid*kgrid;
for (int i=-kgrid; i<(kgrid+1); i++)
{
for (int j=-kgrid; j<(kgrid+1); j++)
{
if (i*i+j*j<=kgrid_squared) //if (std::sqrt(i*i+j*j+k*k)<=kgrid)
{
PosType kpt;
kpt[0]=i-twist[0];
kpt[1]=j-twist[1];
//convert to Cartesian: note that 2Pi is multiplied
kpt=Lattice.k_cart(kpt);
kPoints.push_back(kpt);
mappedQtonofK[i+kgrid].push_back(indx);
mappedQtonofK[j+kgrid+(2*kgrid+1)].push_back(indx);
indx++;
}
}
}
#endif
}
kids=kids->next;
}
#endif
if (rootNode)
{
std::stringstream sstr;