Skip to content

Commit

Permalink
Fix failing complex build
Browse files Browse the repository at this point in the history
Signed-off-by: Steven Hahn <[email protected]>
  • Loading branch information
quantumsteve authored and williamfgc committed Sep 6, 2023
1 parent f78c3b0 commit c03543c
Show file tree
Hide file tree
Showing 5 changed files with 280 additions and 7 deletions.
2 changes: 1 addition & 1 deletion src/QMCWaveFunctions/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -87,7 +87,7 @@ if(OHMMS_DIM MATCHES 3)
LCAO/AOBasisBuilder.cpp
LCAO/SoaLocalizedBasisSet.cpp)
if(QMC_COMPLEX)
set(FERMION_SRCS ${FERMION_SRCS} LCAO/LCAOSpinorBuilder.cpp)
set(FERMION_SRCS ${FERMION_SRCS} LCAO/LCAOSpinorBuilder.cpp LCAO/LCAOSpinorBuilder.cpp)
else(QMC_COMPLEX)
#LCAO cusp correction is not ready for complex
set(FERMION_SRCS ${FERMION_SRCS}
Expand Down
207 changes: 207 additions & 0 deletions src/QMCWaveFunctions/LCAO/LCAOSpinorBuilderT.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,207 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2020 QMCPACK developers.
//
// File developed by: Cody A. Melton, [email protected], Sandia National Laboratories
//
// File created by: Cody A. Melton, [email protected], Sandia National Laboratories
//////////////////////////////////////////////////////////////////////////////////////

#include "LCAOSpinorBuilder.h"
#include "QMCWaveFunctions/SpinorSet.h"
#include "OhmmsData/AttributeSet.h"
#include "Utilities/ProgressReportEngine.h"
#include "hdf/hdf_archive.h"
#include "Message/CommOperators.h"

namespace qmcplusplus
{
template<class T>
LCAOSpinorBuilderT<T>::LCAOSpinorBuilderT(ParticleSet& els, ParticleSet& ions, Communicate* comm, xmlNodePtr cur)
: LCAOrbitalBuilder(els, ions, comm, cur)
{
ClassName = "LCAOSpinorBuilder";

if (h5_path == "")
myComm->barrier_and_abort("LCAOSpinorBuilder only works with href");
}

template<class T>
std::unique_ptr<SPOSetT<T>> LCAOSpinorBuilderT<T>::createSPOSetFromXML(xmlNodePtr cur)
{
ReportEngine PRE(ClassName, "createSPO(xmlNodePtr)");
std::string spo_name(""), optimize("no");
std::string basisset_name("LCAOBSet");
OhmmsAttributeSet spoAttrib;
spoAttrib.add(spo_name, "name");
spoAttrib.add(optimize, "optimize");
spoAttrib.add(basisset_name, "basisset");
spoAttrib.put(cur);

BasisSet_t* myBasisSet = nullptr;
if (basisset_map_.find(basisset_name) == basisset_map_.end())
myComm->barrier_and_abort("basisset \"" + basisset_name + "\" cannot be found\n");
else
myBasisSet = basisset_map_[basisset_name].get();

if (optimize == "yes")
app_log() << " SPOSet " << spo_name << " is optimizable\n";

std::unique_ptr<LCAOrbitalSet> upspo =
std::make_unique<LCAOrbitalSet>(spo_name + "_up", std::unique_ptr<BasisSet_t>(myBasisSet->makeClone()));
std::unique_ptr<LCAOrbitalSet> dnspo =
std::make_unique<LCAOrbitalSet>(spo_name + "_dn", std::unique_ptr<BasisSet_t>(myBasisSet->makeClone()));

loadMO(*upspo, *dnspo, cur);

//create spinor and register up/dn
auto spinor_set = std::make_unique<SpinorSet>(spo_name);
spinor_set->set_spos(std::move(upspo), std::move(dnspo));
return spinor_set;
}

template<class T>
bool LCAOSpinorBuilderT<T>::loadMO(LCAOrbitalSetT<T>& up, LCAOrbitalSetT<T>& dn, xmlNodePtr cur)
{
bool PBC = false;
int norb = up.getBasisSetSize();
std::string debugc("no");
OhmmsAttributeSet aAttrib;
aAttrib.add(norb, "size");
aAttrib.add(debugc, "debug");
aAttrib.put(cur);

up.setOrbitalSetSize(norb);
dn.setOrbitalSetSize(norb);

xmlNodePtr occ_ptr = nullptr;
cur = cur->xmlChildrenNode;
while (cur != nullptr)
{
std::string cname((const char*)(cur->name));
if (cname == "occupation")
{
occ_ptr = cur;
}
cur = cur->next;
}

hdf_archive hin(myComm);
if (myComm->rank() == 0)
{
if (!hin.open(h5_path, H5F_ACC_RDONLY))
myComm->barrier_and_abort("LCAOSpinorBuilder::loadMO missing or incorrect path to H5 file.");
hin.push("PBC");
PBC = false;
hin.read(PBC, "PBC");
hin.close();
}
myComm->bcast(PBC);
if (PBC)
myComm->barrier_and_abort("LCAOSpinorBuilder::loadMO lcao spinors not implemented in PBC");

bool success = putFromH5(up, dn, occ_ptr);


if (debugc == "yes")
{
app_log() << "UP: Single-particle orbital coefficients dims=" << up.C->rows() << " x " << up.C->cols()
<< std::endl;
app_log() << *up.C << std::endl;
app_log() << "DN: Single-particle orbital coefficients dims=" << dn.C->rows() << " x " << dn.C->cols()
<< std::endl;
app_log() << *dn.C << std::endl;
}
return success;
}

template<class T>
bool LCAOSpinorBuilderT<T>::putFromH5(LCAOrbitalSetT<T>& up, LCAOrbitalSetT<T>& dn, xmlNodePtr occ_ptr)
{
if (up.getBasisSetSize() == 0 || dn.getBasisSetSize() == 0)
{
myComm->barrier_and_abort("LCASpinorBuilder::loadMO detected ZERO BasisSetSize");
return false;
}

bool success = true;
hdf_archive hin(myComm);
if (myComm->rank() == 0)
{
istd::string setname = "/Super_Twist/eigenset_0";
readRealMatrixFromH5(hin, setname, upReal);
setname += "_imag";
readRealMatrixFromH5(hin, setname, upImag);

af(!hin.open(h5_path, H5F_ACC_RDONLY))
myComm->barrier_and_abort("LCAOSpinorBuilder::putFromH5 missing or incorrect path to H5 file");

Matrix<RealType> upReal;
Matrix<RealType> upImag;
ssert(upReal.rows() == upImag.rows());
assert(upReal.cols() == upImag.cols());

Matrix<ValueType> upTemp(upReal.rows(), upReal.cols());
for (int i = 0; i < upTemp.rows(); i++)
{
for (int j = 0; j < upTemp.cols(); j++)
{
upTemp[i][j] = ValueType(upReal[i][j], upImag[i][j]);
}
}

Matrix<RealType> dnReal;
Matrix<RealType> dnImag;
setname = "/Super_Twist/eigenset_1";
readRealMatrixFromH5(hin, setname, dnReal);
setname += "_imag";
readRealMatrixFromH5(hin, setname, dnImag);

assert(dnReal.rows() == dnImag.rows());
assert(dnReal.cols() == dnImag.cols());

Matrix<ValueType> dnTemp(dnReal.rows(), dnReal.cols());
for (int i = 0; i < dnTemp.rows(); i++)
{
for (int j = 0; j < dnTemp.cols(); j++)
{
dnTemp[i][j] = ValueType(dnReal[i][j], dnImag[i][j]);
}
}

assert(upReal.rows() == dnReal.rows());
assert(upReal.cols() == dnReal.cols());

Occ.resize(upReal.rows());
success = putOccupation(up, occ_ptr);

int norbs = up.getOrbitalSetSize();

int n = 0, i = 0;
while (i < norbs)
{
if (Occ[n] > 0.0)
{
std::copy(upTemp[n], upTemp[n + 1], (*up.C)[i]);
std::copy(dnTemp[n], dnTemp[n + 1], (*dn.C)[i]);
i++;
}
n++;
}

hin.close();
}

#ifdef HAVE_MPI
myComm->comm.broadcast_n(up.C->data(), up.C->size());
myComm->comm.broadcast_n(dn.C->data(), dn.C->size());
#endif

return success;
}

template class LCAOSpinorBuilderT<std::complex<double>>;
template class LCAOSpinorBuilderT<std::complex<float>>;
} // namespace qmcplusplus
64 changes: 64 additions & 0 deletions src/QMCWaveFunctions/LCAO/LCAOSpinorBuilderT.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
//////////////////////////////////////////////////////////////////////////////////////
// This file is distributed under the University of Illinois/NCSA Open Source License.
// See LICENSE file in top directory for details.
//
// Copyright (c) 2020 QMCPACK developers.
//
// File developed by: Cody A. Melton, [email protected], Sandia National Laboratories
//
// File created by: Cody A. Melton, [email protected], Sandia National Laboratories
//////////////////////////////////////////////////////////////////////////////////////


#ifndef QMCPLUSPLUS_SOA_LCAO_SPINOR_BUILDER_H
#define QMCPLUSPLUS_SOA_LCAO_SPINOR_BUILDER_H

#include "QMCWaveFunctions/LCAO/LCAOrbitalBuilderT.h"

namespace qmcplusplus
{
/** @file LCAOSpinorBuidler.h
*
* Derives from LCAOrbitalBuilder.h. Overrides createSPOSetFromXML method to read up and
* down channel from HDF5 and construct SpinorSet
*
*/
template<class T>
class LCAOSpinorBuilderT : public LCAOrbitalBuilderT<T>
{
public:
/** constructor
* \param els reference to the electrons
* \param ions reference to the ions
*
* Derives from LCAOrbitalBuilder, but will require an h5_path to be set
*/
LCAOSpinorBuilderT(ParticleSet& els, ParticleSet& ions, Communicate* comm, xmlNodePtr cur);

/** creates and returns SpinorSet
*
* Creates an up and down LCAOrbitalSet
* calls LCAOSpinorBuilder::loadMO to build up and down from the H5 file
* registers up and down into a SpinorSet and returns
*/
std::unique_ptr<SPOSetT<T>> createSPOSetFromXML(xmlNodePtr cur) override;

private:
/** load the up and down MO sets
*
* checks to make sure not PBC and initialize the Occ vector.
* call putFromH5 to parse the up and down MO coefficients
*/
bool loadMO(LCAOrbitalSetT<T>& up, LCAOrbitalSetT<T>& dn, xmlNodePtr cur);

/** parse h5 file for spinor info
*
* assumes the h5 file has KPTS_0/eigenset_0(_imag) for the real/imag part of up component of spinor
* assumes the h5 file as KPTS_0/eigenset_1(_imag) for the real/imag part of dn component of spinor
* reads the various coefficient matricies and broadcast
* after this, we have up/dn LCAOrbitalSet that can be registered to the SpinorSet
*/
bool putFromH5(LCAOrbitalSetT<T>& up, LCAOrbitalSetT<T>& dn, xmlNodePtr);
};
} // namespace qmcplusplus
#endif
2 changes: 1 addition & 1 deletion src/QMCWaveFunctions/LCAO/LCAOrbitalBuilderT.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -563,7 +563,7 @@ LCAOrbitalBuilderT<T>::createSPOSetFromXML(xmlNodePtr cur)
#ifdef HAVE_MPI
for (int orb_idx = 0; orb_idx < orbital_set_size; orb_idx++)
for (int center_idx = 0; center_idx < num_centers; center_idx++)
broadcastCuspInfo(
CuspCorrectionConstructionT<T>::broadcastCuspInfo(
info(center_idx, orb_idx), *this->myComm, 0);
#endif
}
Expand Down
12 changes: 7 additions & 5 deletions src/QMCWaveFunctions/SPOSetBuilderFactoryT.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@

#if defined(QMC_COMPLEX)
#include "QMCWaveFunctions/EinsplineSpinorSetBuilder.h"
#include "QMCWaveFunctions/LCAO/LCAOSpinorBuilder.h"
#include "QMCWaveFunctions/LCAO/LCAOSpinorBuilderT.h"
#endif

#if defined(HAVE_EINSPLINE)
Expand Down Expand Up @@ -154,7 +154,7 @@ SPOSetBuilderFactoryT<T>::createSPOSetBuilder(xmlNodePtr rootNode)
ions = pit->second.get();
if (targetPtcl.isSpinor())
#ifdef QMC_COMPLEX
bb = std::make_unique<LCAOSpinorBuilder>(
bb = std::make_unique<LCAOSpinorBuilderT<T>>(
targetPtcl, *ions, myComm, rootNode);
#else
PRE.error("Use of lcao spinors requires QMC_COMPLEX=1. Rebuild "
Expand Down Expand Up @@ -253,9 +253,11 @@ SPOSetBuilderFactoryT<T>::addSPOSet(std::unique_ptr<SPOSetT<T>> spo)
template <typename T>
std::string SPOSetBuilderFactoryT<T>::basisset_tag = "basisset";

template class SPOSetBuilderFactoryT<double>;
template class SPOSetBuilderFactoryT<float>;
#ifdef QMC_COMPLEX
template class SPOSetBuilderFactoryT<std::complex<double>>;
template class SPOSetBuilderFactoryT<std::complex<float>>;

#else
template class SPOSetBuilderFactoryT<double>;
template class SPOSetBuilderFactoryT<float>;
#endif
} // namespace qmcplusplus

0 comments on commit c03543c

Please sign in to comment.