Skip to content

Commit

Permalink
Fix pw unit tests.
Browse files Browse the repository at this point in the history
  • Loading branch information
ye-luo committed Sep 28, 2023
1 parent 6ebe8f6 commit 5c1bc51
Show file tree
Hide file tree
Showing 4 changed files with 34 additions and 67 deletions.
6 changes: 2 additions & 4 deletions src/QMCWaveFunctions/PlaneWave/PWBasis.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -34,10 +34,8 @@ int PWBasis::readbasis(hdf_archive& h5basisgroup,
h5basisgroup.read(gvecs, "/electrons/kpoint_0/gvectors");
NumPlaneWaves = std::max(gvecs.size(), kplusgvecs_cart.size());
if (NumPlaneWaves == 0)
{
app_error() << " PWBasis::readbasis Basis is missing. Abort " << std::endl;
abort(); //FIX_ABORT
}
throw std::runtime_error(" PWBasis::readbasis Basis is missing.");

if (kplusgvecs_cart.empty())
{
kplusgvecs_cart.resize(NumPlaneWaves);
Expand Down
23 changes: 10 additions & 13 deletions src/QMCWaveFunctions/PlaneWave/PWOrbitalSetBuilder.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -34,13 +34,10 @@ PWOrbitalSetBuilder::PWOrbitalSetBuilder(const ParticleSet& p, Communicate* comm
//check the current href
bool success = getH5(cur, "href");
//Move through the XML tree and read basis information
cur = cur->children;
while (cur != nullptr)
{
std::string cname((const char*)(cur->name));
processChildren(cur, [&](const std::string& cname, const xmlNodePtr element) {
if (cname == "basisset")
{
const std::string a(getXMLAttributeValue(cur, "ecut"));
const std::string a(getXMLAttributeValue(element, "ecut"));
if (!a.empty())
myParam->Ecut = std::stod(a);
}
Expand All @@ -49,13 +46,14 @@ PWOrbitalSetBuilder::PWOrbitalSetBuilder(const ParticleSet& p, Communicate* comm
//close
if (success)
hfile.close();
success = getH5(cur, "hdata");
success = getH5(element, "hdata");
}
cur = cur->next;
}
});

if (!success)
throw std::runtime_error("h5 cannot be open for creating PW basis!");
//create PW Basis
createPWBasis(cur);
createPWBasis();
}

PWOrbitalSetBuilder::~PWOrbitalSetBuilder() = default;
Expand All @@ -82,10 +80,10 @@ std::unique_ptr<SPOSet> PWOrbitalSetBuilder::createSPOSetFromXML(xmlNodePtr cur)
* -- maximum_ecut
* - basis
*/
bool PWOrbitalSetBuilder::createPWBasis(xmlNodePtr cur)
bool PWOrbitalSetBuilder::createPWBasis()
{
//recycle int and double reader
int idata;
int idata = 0;
//start of parameters
hfile.read(idata, "electrons/number_of_kpoints");
int nkpts = idata;
Expand All @@ -110,9 +108,8 @@ bool PWOrbitalSetBuilder::createPWBasis(xmlNodePtr cur)
hfile.read(TwistAngle_DP, "/electrons/kpoint_0/reduced_k");
TwistAngle = TwistAngle_DP;
if (!myBasisSet)
{
myBasisSet = std::make_unique<PWBasis>(TwistAngle);
}

//Read the planewave basisset.
//Note that the same data is opened here for each twist angle-avoids duplication in the
//h5 file (which may become very large).
Expand Down
2 changes: 1 addition & 1 deletion src/QMCWaveFunctions/PlaneWave/PWOrbitalSetBuilder.h
Original file line number Diff line number Diff line change
Expand Up @@ -65,7 +65,7 @@ class PWOrbitalSetBuilder : public SPOSetBuilder

private:
bool getH5(xmlNodePtr cur, const char* aname);
bool createPWBasis(xmlNodePtr cur);
bool createPWBasis();
std::unique_ptr<SPOSet> createPW(xmlNodePtr cur, const std::string& objname, int spinIndex);
#if defined(QMC_COMPLEX)
void transform2GridData(PWBasis::GIndex_t& nG, int spinIndex, PWOrbitalSet& pwFunc);
Expand Down
70 changes: 21 additions & 49 deletions src/QMCWaveFunctions/tests/test_pw.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -16,9 +16,7 @@
#include "OhmmsPETE/OhmmsMatrix.h"
#include "Particle/ParticleSet.h"
#include "Particle/ParticleSetPool.h"
#include "QMCWaveFunctions/WaveFunctionComponent.h"
#include "QMCWaveFunctions/PlaneWave/PWOrbitalBuilder.h"
#include "QMCWaveFunctions/Fermion/SlaterDet.h"
#include "PlaneWave/PWOrbitalSetBuilder.h"


#include <stdio.h>
Expand Down Expand Up @@ -70,40 +68,27 @@ TEST_CASE("PlaneWave SPO from HDF for BCC H", "[wavefunction]")
elec.update();

//BCC H
const char* particles = R"(<tmp>
<determinantset type="PW" href="bccH.pwscf.h5" tilematrix="1 0 0 0 1 0 0 0 1" twistnum="0" source="ion">
<slaterdeterminant>
<determinant id="updet" size="1">
<occupation mode="ground" spindataset="0"/>
</determinant>
<determinant id="downdet" size="1">
<occupation mode="ground" spindataset="0"/>
</determinant>
</slaterdeterminant>
</determinantset>
</tmp>
const char* particles = R"(
<sposet_collection type="PW" href="bccH.pwscf.h5" tilematrix="1 0 0 0 1 0 0 0 1" twistnum="0" source="ion">
<sposet name="updet" size="1" spindataset="0">
<occupation mode="ground"/>
</sposet>
</sposet_collection>
)";

Libxml2Document doc;
bool okay = doc.parseFromString(particles);
REQUIRE(okay);

xmlNodePtr root = doc.getRoot();

xmlNodePtr pw1 = xmlFirstElementChild(root);


PWOrbitalBuilder pw_builder(c, elec, ptcl.getPool());
auto orb = pw_builder.buildComponent(pw1);
SlaterDet* sd = dynamic_cast<SlaterDet*>(orb.get());
REQUIRE(sd != nullptr);
REQUIRE(sd->Dets.size() == 2);
SPOSetPtr spo = sd->getPhi(0);
REQUIRE(spo != nullptr);
//SPOSet *spo = einSet.createSPOSetFromXML(ein1);
//REQUIRE(spo != nullptr);
PWOrbitalSetBuilder pw_builder(elec, c, root);
auto spo = pw_builder.createSPOSet(pw1);
REQUIRE(spo);

int orbSize = spo->getOrbitalSetSize();
const int orbSize = spo->getOrbitalSetSize();
elec.update();
SPOSet::ValueVector orbs(orbSize);
spo->evaluateValue(elec, 0, orbs);
Expand Down Expand Up @@ -189,40 +174,27 @@ TEST_CASE("PlaneWave SPO from HDF for LiH arb", "[wavefunction]")
elec.update();

//diamondC_1x1x1
const char* particles = R"(<tmp>
<determinantset type="PW" href="LiH-arb.pwscf.h5" tilematrix="1 0 0 0 1 0 0 0 1" twistnum="0" source="ion">
<slaterdeterminant>
<determinant id="updet" size="2">
<occupation mode="ground" spindataset="0"/>
</determinant>
<determinant id="downdet" size="2">
<occupation mode="ground" spindataset="0"/>
</determinant>
</slaterdeterminant>
</determinantset>
</tmp>
const char* particles = R"(
<sposet_collection type="PW" href="LiH-arb.pwscf.h5" tilematrix="1 0 0 0 1 0 0 0 1" twistnum="0" source="ion">
<sposet name="updet" size="2" spindataset="0">
<occupation mode="ground"/>
</sposet>
</sposet_collection>
)";

Libxml2Document doc;
bool okay = doc.parseFromString(particles);
REQUIRE(okay);

xmlNodePtr root = doc.getRoot();

xmlNodePtr pw1 = xmlFirstElementChild(root);


PWOrbitalBuilder pw_builder(c, elec, ptcl.getPool());
auto orb = pw_builder.buildComponent(pw1);
SlaterDet* sd = dynamic_cast<SlaterDet*>(orb.get());
REQUIRE(sd != nullptr);
REQUIRE(sd->Dets.size() == 2);
SPOSetPtr spo = sd->getPhi(0);
REQUIRE(spo != nullptr);
//SPOSet *spo = einSet.createSPOSetFromXML(ein1);
//REQUIRE(spo != nullptr);
PWOrbitalSetBuilder pw_builder(elec, c, root);
auto spo = pw_builder.createSPOSet(pw1);
REQUIRE(spo);

int orbSize = spo->getOrbitalSetSize();
const int orbSize = spo->getOrbitalSetSize();
elec.update();
SPOSet::ValueVector orbs(orbSize);
spo->evaluateValue(elec, 0, orbs);
Expand Down

0 comments on commit 5c1bc51

Please sign in to comment.