Skip to content

Commit

Permalink
Implement SPOSetBuilderFactoryT and most required builders
Browse files Browse the repository at this point in the history
  • Loading branch information
@PhilipFackler @williamfgc
PhilipFackler authored and williamfgc committed Sep 6, 2023
1 parent 0eef653 commit 47bc03b
Show file tree
Hide file tree
Showing 20 changed files with 3,755 additions and 480 deletions.
26 changes: 21 additions & 5 deletions src/QMCWaveFunctions/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -67,21 +67,36 @@ set(JASTROW_SRCS
set(JASTROW_OMPTARGET_SRCS
Jastrow/TwoBodyJastrow.cpp
Jastrow/BsplineFunctor.cpp)
set(FERMION_SRCS ${FERMION_SRCS} ElectronGas/FreeOrbital.cpp ElectronGas/FreeOrbitalT.cpp ElectronGas/FreeOrbitalBuilder.cpp)
set(FERMION_SRCS ${FERMION_SRCS}
ElectronGas/FreeOrbital.cpp
ElectronGas/FreeOrbitalT.cpp
ElectronGas/FreeOrbitalBuilder.cpp
ElectronGas/FreeOrbitalBuilderT.cpp)

# wavefunctions only availbale to 3-dim problems
if(OHMMS_DIM MATCHES 3)

set(JASTROW_SRCS ${JASTROW_SRCS} Jastrow/eeI_JastrowBuilder.cpp Jastrow/CountingJastrowBuilder.cpp)

set(FERMION_SRCS ${FERMION_SRCS} LCAO/LCAOrbitalSet.cpp LCAO/LCAOrbitalSetT.cpp LCAO/LCAOrbitalBuilder.cpp LCAO/MultiQuinticSpline1D.cpp
LCAO/AOBasisBuilder.cpp LCAO/SoaLocalizedBasisSet.cpp)
set(FERMION_SRCS ${FERMION_SRCS}
LCAO/LCAOrbitalSet.cpp
LCAO/LCAOrbitalSetT.cpp
LCAO/LCAOrbitalBuilder.cpp
LCAO/LCAOrbitalBuilderT.cpp
LCAO/MultiQuinticSpline1D.cpp
LCAO/AOBasisBuilder.cpp
LCAO/SoaLocalizedBasisSet.cpp)
if(QMC_COMPLEX)
set(FERMION_SRCS ${FERMION_SRCS} LCAO/LCAOSpinorBuilder.cpp)
else(QMC_COMPLEX)
#LCAO cusp correction is not ready for complex
set(FERMION_SRCS ${FERMION_SRCS} LCAO/LCAOrbitalSetWithCorrection.cpp LCAO/LCAOrbitalSetWithCorrectionT.cpp
LCAO/CuspCorrectionConstruction.cpp LCAO/SoaCuspCorrection.cpp LCAO/SoaCuspCorrectionT.cpp)
set(FERMION_SRCS ${FERMION_SRCS}
LCAO/LCAOrbitalSetWithCorrection.cpp
LCAO/LCAOrbitalSetWithCorrectionT.cpp
LCAO/CuspCorrectionConstruction.cpp
LCAO/CuspCorrectionConstructionT.cpp
LCAO/SoaCuspCorrection.cpp
LCAO/SoaCuspCorrectionT.cpp)
endif(QMC_COMPLEX)

if(HAVE_EINSPLINE)
Expand Down Expand Up @@ -134,6 +149,7 @@ set(FERMION_SRCS
Fermion/DiracDeterminantWithBackflow.cpp
Fermion/SlaterDetWithBackflow.cpp
SPOSetBuilderFactory.cpp
SPOSetBuilderFactoryT.cpp
TrialWaveFunction.cpp
TWFdispatcher.cpp
TWFFastDerivWrapper.cpp
Expand Down
196 changes: 114 additions & 82 deletions src/QMCWaveFunctions/CompositeSPOSetT.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,10 +14,9 @@
// at Urbana-Champaign
//////////////////////////////////////////////////////////////////////////////////////

#include "CompositeSPOSetT.h"
#include "QMCWaveFunctions/CompositeSPOSetT.h"

#include "OhmmsData/AttributeSet.h"
#include "QMCWaveFunctions/SPOSetBuilderFactory.h"
#include "Utilities/IteratorUtility.h"

#include <algorithm>
Expand All @@ -38,27 +37,27 @@ template <typename MAT1, typename MAT2>
inline void
insert_columns(const MAT1& small, MAT2& big, int offset_c)
{
const int c = small.cols();
for (int i = 0; i < small.rows(); ++i)
std::copy(small[i], small[i] + c, big[i] + offset_c);
const int c = small.cols();
for (int i = 0; i < small.rows(); ++i)
std::copy(small[i], small[i] + c, big[i] + offset_c);
}
} // namespace MatrixOperators

template <typename T>
CompositeSPOSetT<T>::CompositeSPOSetT(const std::string& my_name) :
SPOSetT<T>(my_name)
SPOSetT<T>(my_name)
{
this->OrbitalSetSize = 0;
component_offsets.reserve(4);
this->OrbitalSetSize = 0;
component_offsets.reserve(4);
}

template <typename T>
CompositeSPOSetT<T>::CompositeSPOSetT(const CompositeSPOSetT<T>& other) :
SPOSetT<T>(other)
SPOSetT<T>(other)
{
for (auto& element : other.components) {
this->add(element->makeClone());
}
for (auto& element : other.components) {
this->add(element->makeClone());
}
}

template <typename T>
Expand All @@ -68,120 +67,120 @@ template <typename T>
void
CompositeSPOSetT<T>::add(std::unique_ptr<SPOSetT<T>> component)
{
if (components.empty())
component_offsets.push_back(0); // add 0
if (components.empty())
component_offsets.push_back(0); // add 0

int norbs = component->size();
components.push_back(std::move(component));
component_values.emplace_back(norbs);
component_gradients.emplace_back(norbs);
component_laplacians.emplace_back(norbs);
int norbs = component->size();
components.push_back(std::move(component));
component_values.emplace_back(norbs);
component_gradients.emplace_back(norbs);
component_laplacians.emplace_back(norbs);

this->OrbitalSetSize += norbs;
component_offsets.push_back(this->OrbitalSetSize);
this->OrbitalSetSize += norbs;
component_offsets.push_back(this->OrbitalSetSize);
}

template <typename T>
void
CompositeSPOSetT<T>::report()
{
app_log() << "CompositeSPOSetT" << std::endl;
app_log() << " ncomponents = " << components.size() << std::endl;
app_log() << " components" << std::endl;
for (int i = 0; i < components.size(); ++i) {
app_log() << " " << i << std::endl;
components[i]->basic_report(" ");
}
app_log() << "CompositeSPOSetT" << std::endl;
app_log() << " ncomponents = " << components.size() << std::endl;
app_log() << " components" << std::endl;
for (int i = 0; i < components.size(); ++i) {
app_log() << " " << i << std::endl;
components[i]->basic_report(" ");
}
}

template <typename T>
std::unique_ptr<SPOSetT<T>>
CompositeSPOSetT<T>::makeClone() const
{
return std::make_unique<CompositeSPOSetT<T>>(*this);
return std::make_unique<CompositeSPOSetT<T>>(*this);
}

template <typename T>
void
CompositeSPOSetT<T>::evaluateValue(
const ParticleSet& P, int iat, ValueVector& psi)
const ParticleSet& P, int iat, ValueVector& psi)
{
int n = 0;
for (int c = 0; c < components.size(); ++c) {
SPOSetT<T>& component = *components[c];
ValueVector& values = component_values[c];
component.evaluateValue(P, iat, values);
std::copy(values.begin(), values.end(), psi.begin() + n);
n += component.size();
}
int n = 0;
for (int c = 0; c < components.size(); ++c) {
SPOSetT<T>& component = *components[c];
ValueVector& values = component_values[c];
component.evaluateValue(P, iat, values);
std::copy(values.begin(), values.end(), psi.begin() + n);
n += component.size();
}
}

template <typename T>
void
CompositeSPOSetT<T>::evaluateVGL(const ParticleSet& P, int iat,
ValueVector& psi, GradVector& dpsi, ValueVector& d2psi)
ValueVector& psi, GradVector& dpsi, ValueVector& d2psi)
{
int n = 0;
for (int c = 0; c < components.size(); ++c) {
SPOSetT<T>& component = *components[c];
ValueVector& values = component_values[c];
GradVector& gradients = component_gradients[c];
ValueVector& laplacians = component_laplacians[c];
component.evaluateVGL(P, iat, values, gradients, laplacians);
std::copy(values.begin(), values.end(), psi.begin() + n);
std::copy(gradients.begin(), gradients.end(), dpsi.begin() + n);
std::copy(laplacians.begin(), laplacians.end(), d2psi.begin() + n);
n += component.size();
}
int n = 0;
for (int c = 0; c < components.size(); ++c) {
SPOSetT<T>& component = *components[c];
ValueVector& values = component_values[c];
GradVector& gradients = component_gradients[c];
ValueVector& laplacians = component_laplacians[c];
component.evaluateVGL(P, iat, values, gradients, laplacians);
std::copy(values.begin(), values.end(), psi.begin() + n);
std::copy(gradients.begin(), gradients.end(), dpsi.begin() + n);
std::copy(laplacians.begin(), laplacians.end(), d2psi.begin() + n);
n += component.size();
}
}

template <typename T>
void
CompositeSPOSetT<T>::evaluate_notranspose(const ParticleSet& P, int first,
int last, ValueMatrix& logdet, GradMatrix& dlogdet, ValueMatrix& d2logdet)
int last, ValueMatrix& logdet, GradMatrix& dlogdet, ValueMatrix& d2logdet)
{
const int nat = last - first;
for (int c = 0; c < components.size(); ++c) {
int norb = components[c]->size();
ValueMatrix v(nat, norb);
GradMatrix g(nat, norb);
ValueMatrix l(nat, norb);
components[c]->evaluate_notranspose(P, first, last, v, g, l);
int n = component_offsets[c];
MatrixOperators::insert_columns(v, logdet, n);
MatrixOperators::insert_columns(g, dlogdet, n);
MatrixOperators::insert_columns(l, d2logdet, n);
}
const int nat = last - first;
for (int c = 0; c < components.size(); ++c) {
int norb = components[c]->size();
ValueMatrix v(nat, norb);
GradMatrix g(nat, norb);
ValueMatrix l(nat, norb);
components[c]->evaluate_notranspose(P, first, last, v, g, l);
int n = component_offsets[c];
MatrixOperators::insert_columns(v, logdet, n);
MatrixOperators::insert_columns(g, dlogdet, n);
MatrixOperators::insert_columns(l, d2logdet, n);
}
}

template <typename T>
void
CompositeSPOSetT<T>::evaluate_notranspose(const ParticleSet& P, int first,
int last, ValueMatrix& logdet, GradMatrix& dlogdet,
HessMatrix& grad_grad_logdet)
int last, ValueMatrix& logdet, GradMatrix& dlogdet,
HessMatrix& grad_grad_logdet)
{
const int nat = last - first;
for (int c = 0; c < components.size(); ++c) {
int norb = components[c]->size();
ValueMatrix v(nat, norb);
GradMatrix g(nat, norb);
HessMatrix h(nat, norb);
components[c]->evaluate_notranspose(P, first, last, v, g, h);
int n = component_offsets[c];
MatrixOperators::insert_columns(v, logdet, n);
MatrixOperators::insert_columns(g, dlogdet, n);
MatrixOperators::insert_columns(h, grad_grad_logdet, n);
}
const int nat = last - first;
for (int c = 0; c < components.size(); ++c) {
int norb = components[c]->size();
ValueMatrix v(nat, norb);
GradMatrix g(nat, norb);
HessMatrix h(nat, norb);
components[c]->evaluate_notranspose(P, first, last, v, g, h);
int n = component_offsets[c];
MatrixOperators::insert_columns(v, logdet, n);
MatrixOperators::insert_columns(g, dlogdet, n);
MatrixOperators::insert_columns(h, grad_grad_logdet, n);
}
}

template <typename T>
void
CompositeSPOSetT<T>::evaluate_notranspose(const ParticleSet& P, int first,
int last, ValueMatrix& logdet, GradMatrix& dlogdet,
HessMatrix& grad_grad_logdet, GGGMatrix& grad_grad_grad_logdet)
int last, ValueMatrix& logdet, GradMatrix& dlogdet,
HessMatrix& grad_grad_logdet, GGGMatrix& grad_grad_grad_logdet)
{
not_implemented(
"evaluate_notranspose(P,first,last,logdet,dlogdet,ddlogdet,dddlogdet)");
not_implemented(
"evaluate_notranspose(P,first,last,logdet,dlogdet,ddlogdet,dddlogdet)");
}

// Class concrete types from ValueType
Expand All @@ -190,4 +189,37 @@ template class CompositeSPOSetT<float>;
template class CompositeSPOSetT<std::complex<double>>;
template class CompositeSPOSetT<std::complex<float>>;

template <typename T>
std::unique_ptr<SPOSetT<T>>
CompositeSPOSetBuilderT<T>::createSPOSetFromXML(xmlNodePtr cur)
{
std::vector<std::string> spolist;
putContent(spolist, cur);
if (spolist.empty()) {
return nullptr;
}

auto spo_now = std::make_unique<CompositeSPOSetT<T>>(
getXMLAttributeValue(cur, "name"));
for (int i = 0; i < spolist.size(); ++i) {
const SPOSetT<T>* spo = sposet_builder_factory_.getSPOSet(spolist[i]);
if (spo)
spo_now->add(spo->makeClone());
}
return (spo_now->size()) ? std::unique_ptr<SPOSetT<T>>{std::move(spo_now)} :
nullptr;
}

template <typename T>
std::unique_ptr<SPOSetT<T>>
CompositeSPOSetBuilderT<T>::createSPOSet(xmlNodePtr cur, SPOSetInputInfo& input)
{
return createSPOSetFromXML(cur);
}

template class CompositeSPOSetBuilderT<double>;
template class CompositeSPOSetBuilderT<float>;
template class CompositeSPOSetBuilderT<std::complex<double>>;
template class CompositeSPOSetBuilderT<std::complex<float>>;

} // namespace qmcplusplus
Loading

0 comments on commit 47bc03b

Please sign in to comment.