Merge pull request #4414 from rcclay/constant_sposet_for_debug

Create ConstantSPOSet for Unit Testing

src/QMCWaveFunctions/tests/CMakeLists.txt

@@ -100,6 +100,7 @@ set(TRIALWF_SRC
     test_TWFGrads.cpp)
 
 set(SPOSET_SRC
+    ConstantSPOSet.cpp
     test_spo_collection_input_spline.cpp
     test_spo_collection_input_LCAO_xml.cpp
     test_spo_collection_input_MSD_LCAO_h5.cpp
@@ -108,6 +109,7 @@ set(SPOSET_SRC
     test_CompositeSPOSet.cpp
     test_hybridrep.cpp
     test_pw.cpp
+    test_ConstantSPOSet.cpp
     ${MO_SRCS})
 if(NiO_a16_H5_FOUND)
   set(SPOSET_SRC ${SPOSET_SRC} test_einset_NiO_a16.cpp)

src/QMCWaveFunctions/tests/ConstantSPOSet.cpp

@@ -0,0 +1,98 @@
+//////////////////////////////////////////////////////////////////////////////////////
+// This file is distributed under the University of Illinois/NCSA Open Source License.
+// See LICENSE file in top directory for details.
+//
+// Copyright (c) 2023 Raymond Clay and QMCPACK developers.
+//
+// File developed by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//
+// File created by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//////////////////////////////////////////////////////////////////////////////////////
+
+#include "QMCWaveFunctions/tests/ConstantSPOSet.h"
+
+namespace qmcplusplus
+{
+ConstantSPOSet::ConstantSPOSet(const std::string& my_name, const int nparticles, const int norbitals)
+    : SPOSet(my_name), numparticles_(nparticles)
+{
+  OrbitalSetSize = norbitals;
+  ref_psi_.resize(numparticles_, OrbitalSetSize);
+  ref_egrad_.resize(numparticles_, OrbitalSetSize);
+  ref_elapl_.resize(numparticles_, OrbitalSetSize);
+
+  ref_psi_   = 0.0;
+  ref_egrad_ = 0.0;
+  ref_elapl_ = 0.0;
+};
+
+std::unique_ptr<SPOSet> ConstantSPOSet::makeClone() const
+{
+  auto myclone = std::make_unique<ConstantSPOSet>(my_name_, numparticles_, OrbitalSetSize);
+  myclone->setRefVals(ref_psi_);
+  myclone->setRefEGrads(ref_egrad_);
+  myclone->setRefELapls(ref_elapl_);
+  return myclone;
+};
+
+std::string ConstantSPOSet::getClassName() const { return "ConstantSPOSet"; };
+
+void ConstantSPOSet::checkOutVariables(const opt_variables_type& active)
+{
+  APP_ABORT("ConstantSPOSet should not call checkOutVariables");
+};
+
+void ConstantSPOSet::setOrbitalSetSize(int norbs) { APP_ABORT("ConstantSPOSet should not call setOrbitalSetSize()"); }
+
+void ConstantSPOSet::setRefVals(const ValueMatrix& vals)
+{
+  assert(vals.cols() == OrbitalSetSize);
+  assert(vals.rows() == numparticles_);
+  ref_psi_ = vals;
+};
+void ConstantSPOSet::setRefEGrads(const GradMatrix& grads)
+{
+  assert(grads.cols() == OrbitalSetSize);
+  assert(grads.rows() == numparticles_);
+  ref_egrad_ = grads;
+};
+void ConstantSPOSet::setRefELapls(const ValueMatrix& lapls)
+{
+  assert(lapls.cols() == OrbitalSetSize);
+  assert(lapls.rows() == numparticles_);
+  ref_elapl_ = lapls;
+};
+
+void ConstantSPOSet::evaluateValue(const ParticleSet& P, int iat, ValueVector& psi)
+{
+  assert(psi.size() == OrbitalSetSize);
+  for (int iorb = 0; iorb < OrbitalSetSize; iorb++)
+    psi[iorb] = ref_psi_(iat, iorb);
+};
+
+void ConstantSPOSet::evaluateVGL(const ParticleSet& P, int iat, ValueVector& psi, GradVector& dpsi, ValueVector& d2psi)
+{
+  for (int iorb = 0; iorb < OrbitalSetSize; iorb++)
+  {
+    psi[iorb]   = ref_psi_(iat, iorb);
+    dpsi[iorb]  = ref_egrad_(iat, iorb);
+    d2psi[iorb] = ref_elapl_(iat, iorb);
+  }
+};
+
+void ConstantSPOSet::evaluate_notranspose(const ParticleSet& P,
+                                          int first,
+                                          int last,
+                                          ValueMatrix& logdet,
+                                          GradMatrix& dlogdet,
+                                          ValueMatrix& d2logdet)
+{
+  for (int iat = first, i = 0; iat < last; ++iat, ++i)
+  {
+    ValueVector v(logdet[i], logdet.cols());
+    GradVector g(dlogdet[i], dlogdet.cols());
+    ValueVector l(d2logdet[i], d2logdet.cols());
+    evaluateVGL(P, iat, v, g, l);
+  }
+}
+} //namespace qmcplusplus

src/QMCWaveFunctions/tests/ConstantSPOSet.h

@@ -0,0 +1,87 @@
+//////////////////////////////////////////////////////////////////////////////////////
+// This file is distributed under the University of Illinois/NCSA Open Source License.
+// See LICENSE file in top directory for details.
+//
+// Copyright (c) 2023 Raymond Clay and QMCPACK developers.
+//
+// File developed by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//
+// File created by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//////////////////////////////////////////////////////////////////////////////////////
+
+
+#ifndef QMCPLUSPLUS_CONSTANTSPOSET_H
+#define QMCPLUSPLUS_CONSTANTSPOSET_H
+
+#include "QMCWaveFunctions/SPOSet.h"
+
+namespace qmcplusplus
+{
+/** Constant SPOSet for testing purposes.  Fixed N_elec x N_orb matrices storing value, gradients, and laplacians, etc.,
+   *  These values are accessed through standard SPOSet calls like evaluateValue, evaluateVGL, etc.
+   *  Exists to provide deterministic and known output to objects requiring SPOSet evaluations.      
+   *
+   */
+class ConstantSPOSet : public SPOSet
+{
+public:
+  ConstantSPOSet(const std::string& my_name) = delete;
+
+  //Constructor needs number of particles and number of orbitals.  This is the minimum
+  //amount of information needed to sanely construct all data members and perform size
+  //checks later.
+  ConstantSPOSet(const std::string& my_name, const int nparticles, const int norbitals);
+
+  std::unique_ptr<SPOSet> makeClone() const override;
+
+  std::string getClassName() const override;
+
+  void checkOutVariables(const opt_variables_type& active) override;
+
+  void setOrbitalSetSize(int norbs) override;
+
+  /**
+  * @brief Setter method to set \phi_j(r_i). Stores input matrix in ref_psi_.
+  * @param Nelec x Nion ValueType matrix of \phi_j(r_i)
+  * @return void
+  */ 
+  void setRefVals(const ValueMatrix& vals);
+  /**
+  * @brief Setter method to set \nabla_i \phi_j(r_i). Stores input matrix in ref_egrad_.
+  * @param Nelec x Nion GradType matrix of \grad_i \phi_j(r_i)
+  * @return void
+  */ 
+  void setRefEGrads(const GradMatrix& grads);
+  /**
+  * @brief Setter method to set \nabla^2_i \phi_j(r_i). Stores input matrix in ref_elapl_.
+  * @param Nelec x Nion GradType matrix of \grad^2_i \phi_j(r_i)
+  * @return void
+  */ 
+  void setRefELapls(const ValueMatrix& lapls);
+
+  void evaluateValue(const ParticleSet& P, int iat, ValueVector& psi) override;
+
+  void evaluateVGL(const ParticleSet& P, int iat, ValueVector& psi, GradVector& dpsi, ValueVector& d2psi) override;
+
+  void evaluate_notranspose(const ParticleSet& P,
+                            int first,
+                            int last,
+                            ValueMatrix& logdet,
+                            GradMatrix& dlogdet,
+                            ValueMatrix& d2logdet) override;
+
+private:
+  const int numparticles_; /// evaluate_notranspose arrays are nparticle x norb matrices.
+                           /// To ensure consistent array sizing and enforcement,
+                           /// we agree at construction how large these matrices will be.
+                           /// norb is stored in SPOSet::OrbitalSetSize.  
+
+  //Value, electron gradient, and electron laplacian at "reference configuration".
+  //i.e. before any attempted moves.
+
+  ValueMatrix ref_psi_;
+  GradMatrix ref_egrad_;
+  ValueMatrix ref_elapl_;
+};
+} // namespace qmcplusplus
+#endif

src/QMCWaveFunctions/tests/test_ConstantSPOSet.cpp

@@ -0,0 +1,136 @@
+//////////////////////////////////////////////////////////////////////////////////////
+// This file is distributed under the University of Illinois/NCSA Open Source License.
+// See LICENSE file in top directory for details.
+//
+// Copyright (c) 2023 QMCPACK developers.
+//
+// File developed by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//
+// File created by: Raymond Clay, rclay@sandia.gov, Sandia National Laboratories
+//////////////////////////////////////////////////////////////////////////////////////
+
+
+#include "catch.hpp"
+#include "Configuration.h"
+#include "QMCWaveFunctions/WaveFunctionTypes.hpp"
+#include "QMCWaveFunctions/tests/ConstantSPOSet.h"
+#include "Utilities/for_testing/checkMatrix.hpp"
+namespace qmcplusplus
+{
+//Ray:  Figure out how to template me on value type.
+TEST_CASE("ConstantSPOSet", "[wavefunction]")
+{
+  //For now, do a small square case.
+  const int nelec   = 2;
+  const int norb    = 2;
+  using WF          = WaveFunctionTypes<QMCTraits::ValueType, QMCTraits::FullPrecValueType>;
+  using Real        = WF::Real;
+  using Value       = WF::Value;
+  using Grad        = WF::Grad;
+  using ValueVector = Vector<Value>;
+  using GradVector  = Vector<Grad>;
+  using ValueMatrix = Matrix<Value>;
+  using GradMatrix  = Matrix<Grad>;
+
+  ValueVector row0{Value(0.92387953), Value(0.92387953)};
+  ValueVector row1{Value(0.29131988), Value(0.81078057)};
+
+  GradVector grow0{Grad({-2.22222, -1.11111, 0.33333}), Grad({8.795388, -0.816057, -0.9238793})};
+  GradVector grow1{Grad({2.22222, 1.11111, -0.33333}), Grad({-8.795388, 0.816057, 0.9238793})};
+
+  ValueVector lrow0{Value(-0.2234545), Value(0.72340234)};
+  ValueVector lrow1{Value(-12.291810), Value(6.879057)};
+
+
+  ValueMatrix spomat;
+  GradMatrix gradspomat;
+  ValueMatrix laplspomat;
+
+  spomat.resize(nelec, norb);
+  gradspomat.resize(nelec, norb);
+  laplspomat.resize(nelec, norb);
+
+  for (int iorb = 0; iorb < norb; iorb++)
+  {
+    spomat(0, iorb) = row0[iorb];
+    spomat(1, iorb) = row1[iorb];
+
+    gradspomat(0, iorb) = grow0[iorb];
+    gradspomat(1, iorb) = grow1[iorb];
+
+    laplspomat(0, iorb) = lrow0[iorb];
+    laplspomat(1, iorb) = lrow1[iorb];
+  }
+
+
+  const SimulationCell simulation_cell;
+  ParticleSet elec(simulation_cell);
+
+  elec.create({nelec});
+
+  ValueVector psiV = {0.0, 0.0};
+  ValueVector psiL = {0.0, 0.0};
+  GradVector psiG;
+  psiG.resize(norb);
+
+  //Test of value only constructor.
+  auto sposet = std::make_unique<ConstantSPOSet>("constant_spo", nelec, norb);
+  sposet->setRefVals(spomat);
+  sposet->setRefEGrads(gradspomat);
+  sposet->setRefELapls(laplspomat);
+
+  sposet->evaluateValue(elec, 0, psiV);
+
+  CHECK(psiV[0] == row0[0]);
+  CHECK(psiV[1] == row0[1]);
+
+
+  psiV = 0.0;
+
+  sposet->evaluateValue(elec, 1, psiV);
+  CHECK(psiV[0] == row1[0]);
+  CHECK(psiV[1] == row1[1]);
+
+  psiV = 0.0;
+
+  sposet->evaluateVGL(elec, 1, psiV, psiG, psiL);
+
+  for (int iorb = 0; iorb < norb; iorb++)
+  {
+    CHECK(psiV[iorb] == row1[iorb]);
+    CHECK(psiL[iorb] == lrow1[iorb]);
+
+    for (int idim = 0; idim < OHMMS_DIM; idim++)
+      CHECK(psiG[iorb][idim] == grow1[iorb][idim]);
+  }
+  //Test of evaluate_notranspose.
+  ValueMatrix phimat, lphimat;
+  GradMatrix gphimat;
+  phimat.resize(nelec, norb);
+  gphimat.resize(nelec, norb);
+  lphimat.resize(nelec, norb);
+
+  const int first_index = 0; //Only 2 electrons in this case.
+  const int last_index  = 2;
+  sposet->evaluate_notranspose(elec, first_index, last_index, phimat, gphimat, lphimat);
+
+  checkMatrix(phimat, spomat);
+  checkMatrix(lphimat, laplspomat);
+
+  //Test of makeClone()
+  auto sposet_vgl2 = sposet->makeClone();
+  phimat           = 0.0;
+  gphimat          = 0.0;
+  lphimat          = 0.0;
+
+  sposet_vgl2->evaluate_notranspose(elec, first_index, last_index, phimat, gphimat, lphimat);
+
+  checkMatrix(phimat, spomat);
+  checkMatrix(lphimat, laplspomat);
+
+  //Lastly, check if name is correct.
+  std::string myname = sposet_vgl2->getClassName();
+  std::string targetstring("ConstantSPOSet");
+  CHECK(myname == targetstring);
+}
+} // namespace qmcplusplus