QMCPACK
DensityMatrices1B Class Reference
+ Inheritance diagram for DensityMatrices1B:
+ Collaboration diagram for DensityMatrices1B:

Public Types

enum  { DIM = OHMMS_DIM }
 
enum  integrators { uniform_grid = 0, uniform, density, no_integrator }
 
enum  evaluators { loop = 0, matrix, no_evaluator }
 
enum  samplings { volume_based = 0, metropolis, no_sampling }
 
using Value_t = ValueType
 
using Grad_t = GradType
 
using ValueVector = SPOSet::ValueVector
 
using GradVector = SPOSet::GradVector
 
using Lattice_t = ParticleSet::ParticleLayout
 
using Vector_t = Vector< Value_t >
 
using Matrix_t = Matrix< Value_t >
 
using pts_t = std::vector< PosType >
 
using dens_t = std::vector< RealType >
 
- Public Types inherited from OperatorBase
enum  EnergyDomains { KINETIC = 0, POTENTIAL, NO_ENERGY_DOMAIN }
 enum to denote energy domain of operators More...
 
enum  QuantumDomains {
  NO_QUANTUM_DOMAIN = 0, CLASSICAL, QUANTUM, CLASSICAL_CLASSICAL,
  QUANTUM_CLASSICAL, QUANTUM_QUANTUM
}
 
enum  {
  PRIMARY = 0, OPTIMIZABLE = 1, RATIOUPDATE = 2, PHYSICAL = 3,
  COLLECTABLE = 4, NONLOCAL = 5
}
 enum for update_mode More...
 
using Return_t = FullPrecRealType
 type of return value of evaluate More...
 
using ValueMatrix = SPOSet::ValueMatrix
 For fast derivative evaluation. More...
 
using GradMatrix = SPOSet::GradMatrix
 
using BufferType = ParticleSet::Buffer_t
 typedef for the serialized buffer More...
 
using Walker_t = ParticleSet::Walker_t
 typedef for the walker More...
 
using ParticleScalar = ParticleSet::Scalar_t
 typedef for the ParticleScalar More...
 
using SPOMap = SPOSet::SPOMap
 typedef for SPOMap More...
 
- Public Types inherited from QMCTraits
enum  { DIM = OHMMS_DIM, DIM_VGL = OHMMS_DIM + 2 }
 
using QTBase = QMCTypes< OHMMS_PRECISION, DIM >
 
using QTFull = QMCTypes< OHMMS_PRECISION_FULL, DIM >
 
using RealType = QTBase::RealType
 
using ComplexType = QTBase::ComplexType
 
using ValueType = QTBase::ValueType
 
using PosType = QTBase::PosType
 
using GradType = QTBase::GradType
 
using TensorType = QTBase::TensorType
 
using IndexType = OHMMS_INDEXTYPE
 define other types More...
 
using FullPrecRealType = QTFull::RealType
 
using FullPrecValueType = QTFull::ValueType
 
using PropertySetType = RecordNamedProperty< FullPrecRealType >
 define PropertyList_t More...
 
using PtclGrpIndexes = std::vector< std::pair< int, int > >
 

Public Member Functions

 DensityMatrices1B (ParticleSet &P, TrialWaveFunction &psi, ParticleSet *Pcl)
 
 DensityMatrices1B (DensityMatrices1B &master, ParticleSet &P, TrialWaveFunction &psi)
 
 ~DensityMatrices1B () override
 
bool dependsOnWaveFunction () const override
 return true if this operator depends on a wavefunction More...
 
std::string getClassName () const override
 return class name More...
 
std::unique_ptr< OperatorBasemakeClone (ParticleSet &P, TrialWaveFunction &psi) final
 
bool put (xmlNodePtr cur) override
 Read the input parameter. More...
 
Return_t evaluate (ParticleSet &P) override
 Evaluate the local energy contribution of this component. More...
 
void getRequiredTraces (TraceManager &tm) override
 TODO: add docs. More...
 
void setRandomGenerator (RandomBase< FullPrecRealType > *rng) override
 Set the Random Generator object TODO: add docs. More...
 
void addObservables (PropertySetType &plist, BufferType &olist) override
 named values to the property list Default implementaton uses addValue(plist_) More...
 
void registerCollectables (std::vector< ObservableHelper > &h5desc, hdf_archive &file) const override
 
void resetTargetParticleSet (ParticleSet &P) override
 Reset the data with the target ParticleSet. More...
 
void setObservables (PropertySetType &plist) override
 Set the values evaluated by this object to plist Default implementation is to assign Value which is updated by evaluate function using my_index_. More...
 
void setParticlePropertyList (PropertySetType &plist, int offset) override
 
void contributeScalarQuantities () override
 
void checkoutScalarQuantities (TraceManager &tm) override
 
void collectScalarQuantities () override
 
void deleteScalarQuantities () override
 
bool get (std::ostream &os) const override
 write about the class More...
 
void reset ()
 
void set_state (xmlNodePtr cur)
 
void set_state (DensityMatrices1B &master)
 
void initialize ()
 
void finalize ()
 
void normalize ()
 
void report (const std::string &pad="")
 
void warmup_sampling ()
 
void generate_samples (RealType weight, int steps=0)
 
void generate_uniform_grid (RandomBase< FullPrecRealType > &rng)
 
void generate_uniform_samples (RandomBase< FullPrecRealType > &rng)
 
void generate_density_samples (bool save, int steps, RandomBase< FullPrecRealType > &rng)
 
void diffusion (RealType sqt, PosType &diff)
 
void density_only (const PosType &r, RealType &dens)
 
void density_drift (const PosType &r, RealType &dens, PosType &drift)
 
void get_energies (std::vector< Vector_t *> &E_n)
 
void generate_sample_basis (Matrix_t &Phi_mb)
 
void generate_sample_ratios (std::vector< Matrix_t *> Psi_nm)
 
void generate_particle_basis (ParticleSet &P, std::vector< Matrix_t *> &Phi_nb)
 
void update_basis (const PosType &r)
 
void update_basis_d012 (const PosType &r)
 
void test_overlap ()
 
void test_derivatives ()
 
void integrate (ParticleSet &P, int n)
 
Return_t evaluate_loop (ParticleSet &P)
 
Return_t evaluate_check (ParticleSet &P)
 
Return_t evaluate_matrix (ParticleSet &P)
 
bool match (Value_t e1, Value_t e2, RealType tol=1e-12)
 
bool same (Vector_t &v1, Vector_t &v2, RealType tol=1e-6)
 
bool same (Matrix_t &m1, Matrix_t &m2, RealType tol=1e-6)
 
void compare (const std::string &name, Vector_t &v1, Vector_t &v2, bool write=false, bool diff_only=true)
 
void compare (const std::string &name, Matrix_t &m1, Matrix_t &m2, bool write=false, bool diff_only=true)
 
- Public Member Functions inherited from OperatorBase
 OperatorBase ()
 Construct a new Operator Base object Default and unique empty constructor. More...
 
virtual ~OperatorBase ()=default
 
std::bitset< 8 > & getUpdateMode () noexcept
 get update_mode_ reference More...
 
Return_t getValue () const noexcept
 get a copy of value_ More...
 
std::string getName () const noexcept
 getter a copy of my_name_, rvalue small string optimization More...
 
void setName (const std::string name) noexcept
 Set my_name member, uses small string optimization (pass by value) More...
 
TraceRequestgetRequest () noexcept
 Get request_ member. More...
 
virtual void registerObservables (std::vector< ObservableHelper > &h5desc, hdf_archive &file) const
 add to observable descriptor for hdf5 The default implementation is to register a scalar for this->value_ More...
 
virtual void setHistories (Walker_t &ThisWalker)
 
virtual Return_t evaluateDeterministic (ParticleSet &P)
 Evaluate the local energy contribution of this component, deterministically based on current state. More...
 
virtual void mw_evaluate (const RefVectorWithLeader< OperatorBase > &o_list, const RefVectorWithLeader< TrialWaveFunction > &wf_list, const RefVectorWithLeader< ParticleSet > &p_list) const
 Evaluate the contribution of this component of multiple walkers. More...
 
virtual void mw_evaluatePerParticle (const RefVectorWithLeader< OperatorBase > &o_list, const RefVectorWithLeader< TrialWaveFunction > &wf_list, const RefVectorWithLeader< ParticleSet > &p_list, const std::vector< ListenerVector< RealType >> &listeners, const std::vector< ListenerVector< RealType >> &listeners_ions) const
 Evaluate the contribution of this component of multiple walkers per particle and report to registerd listeners from objects in Estimators. More...
 
virtual void mw_evaluateWithParameterDerivatives (const RefVectorWithLeader< OperatorBase > &o_list, const RefVectorWithLeader< ParticleSet > &p_list, const opt_variables_type &optvars, const RecordArray< ValueType > &dlogpsi, RecordArray< ValueType > &dhpsioverpsi) const
 TODO: add docs. More...
 
virtual Return_t rejectedMove (ParticleSet &P)
 TODO: add docs. More...
 
virtual Return_t evaluateWithToperator (ParticleSet &P)
 Evaluate the local energy contribution of this component with Toperators updated if requested. More...
 
virtual void mw_evaluateWithToperator (const RefVectorWithLeader< OperatorBase > &o_list, const RefVectorWithLeader< TrialWaveFunction > &wf_list, const RefVectorWithLeader< ParticleSet > &p_list) const
 Evaluate the contribution of this component of multiple walkers. More...
 
virtual void mw_evaluatePerParticleWithToperator (const RefVectorWithLeader< OperatorBase > &o_list, const RefVectorWithLeader< TrialWaveFunction > &wf_list, const RefVectorWithLeader< ParticleSet > &p_list, const std::vector< ListenerVector< RealType >> &listeners, const std::vector< ListenerVector< RealType >> &listeners_ions) const
 Evaluate the contribution of this component of multiple walkers per particle and report to registerd listeners from objects in Estimators. More...
 
virtual Return_t evaluateValueAndDerivatives (ParticleSet &P, const opt_variables_type &optvars, const Vector< ValueType > &dlogpsi, Vector< ValueType > &dhpsioverpsi)
 Evaluate value and derivatives wrt the optimizables. More...
 
virtual Return_t evaluateWithIonDerivs (ParticleSet &P, ParticleSet &ions, TrialWaveFunction &psi, ParticleSet::ParticlePos &hf_term, ParticleSet::ParticlePos &pulay_term)
 Evaluate contribution to local energy and derivatives w.r.t ionic coordinates from OperatorBase. More...
 
virtual Return_t evaluateWithIonDerivsDeterministic (ParticleSet &P, ParticleSet &ions, TrialWaveFunction &psi, ParticleSet::ParticlePos &hf_term, ParticleSet::ParticlePos &pulay_term)
 Evaluate contribution to local energy and derivatives w.r.t ionic coordinates from OperatorBase. More...
 
virtual void evaluateOneBodyOpMatrix (ParticleSet &P, const TWFFastDerivWrapper &psi, std::vector< ValueMatrix > &B)
 Evaluate "B" matrix for observable. More...
 
virtual void evaluateOneBodyOpMatrixForceDeriv (ParticleSet &P, ParticleSet &source, const TWFFastDerivWrapper &psi, const int iat, std::vector< std::vector< ValueMatrix >> &Bforce)
 Evaluate "dB/dR" matrices for observable. More...
 
virtual void updateSource (ParticleSet &s)
 Update data associated with a particleset. More...
 
virtual Return_t getEnsembleAverage ()
 Return an average value by collective operation. More...
 
virtual void createResource (ResourceCollection &collection) const
 Initialize a shared resource and hand it to a collection. More...
 
virtual void acquireResource (ResourceCollection &collection, const RefVectorWithLeader< OperatorBase > &o_list) const
 Acquire a shared resource from a collection. More...
 
virtual void releaseResource (ResourceCollection &collection, const RefVectorWithLeader< OperatorBase > &o_list) const
 Return a shared resource to a collection. More...
 
virtual void add2Hamiltonian (ParticleSet &qp, TrialWaveFunction &psi, QMCHamiltonian &targetH)
 TODO: add docs. More...
 
virtual void informOfPerParticleListener ()
 
bool isClassical () const noexcept
 
bool isQuantum () const noexcept
 
bool isClassicalClassical () const noexcept
 
bool isQuantumClassical () const noexcept
 
bool isQuantumQuantum () const noexcept
 
bool getMode (const int i) const noexcept
 Return the mode i. More...
 
bool isNonLocal () const noexcept
 TODO: add docs. More...
 
bool hasListener () const noexcept
 
void contributeTraceQuantities ()
 Make trace quantities available. More...
 
void checkoutTraceQuantities (TraceManager &tm)
 Checkout trace arrays Derived classes must guard individual checkouts using request info. More...
 
void collectScalarTraces ()
 Collect scalar trace data. More...
 
void deleteTraceQuantities ()
 delete trace arrays More...
 

Public Attributes

bool energy_mat
 
CompositeSPOSet basis_functions
 
ValueVector basis_values
 
ValueVector basis_norms
 
GradVector basis_gradients
 
ValueVector basis_laplacians
 
ValueVector integrated_values
 
bool warmed_up
 
std::vector< PosTypersamples
 
Vector< RealTypesample_weights
 
std::vector< ValueTypepsi_ratios
 
RealType dens
 
PosType drift
 
int nindex
 
int eindex
 
const Lattice_tlattice_
 
TrialWaveFunctionPsi
 
ParticleSetPq
 
const ParticleSetPc
 
TraceSample< TraceReal > * w_trace
 
TraceSample< TraceComp > * T_trace
 
CombinedTraceSample< TraceReal > * Vq_trace
 
CombinedTraceSample< TraceReal > * Vc_trace
 
CombinedTraceSample< TraceReal > * Vqq_trace
 
CombinedTraceSample< TraceReal > * Vqc_trace
 
CombinedTraceSample< TraceReal > * Vcc_trace
 
std::vector< Value_tE_samp
 
CombinedTraceSample< TraceReal > * E_trace
 
bool initialized
 
bool normalized
 
bool volume_normed
 
int basis_size
 
int samples
 
int nparticles
 
int nspecies
 
std::vector< int > species_size
 
std::vector< std::string > species_name
 
std::vector< Vector_t * > E_N
 
std::vector< Matrix_t * > Phi_NB
 
std::vector< Matrix_t * > Psi_NM
 
std::vector< Matrix_t * > Phi_Psi_NB
 
std::vector< Matrix_t * > N_BB
 
std::vector< Matrix_t * > E_BB
 
Matrix_t Phi_MB
 
bool check_overlap
 
bool check_derivatives
 
integrators integrator
 
evaluators evaluator
 
samplings sampling
 
int points
 
RealType scale
 
PosType center
 
PosType rcorner
 
RealType volume
 
bool periodic
 
int warmup
 
RealType timestep
 
bool use_drift
 
int nmoves
 
int naccepted
 
RealType acceptance_ratio
 
bool write_acceptance_ratio
 
bool write_rstats
 
int ind_dims [DIM]
 
RealType metric
 
PosType rpcur
 
PosType dpcur
 
RealType rhocur
 
RandomBase< FullPrecRealType > * uniform_random
 

Protected Attributes

TimerList_t timers
 
- Protected Attributes inherited from OperatorBase
std::bitset< 8 > update_mode_
 set the current update mode More...
 
Return_t value_
 current value More...
 
std::string name_
 name of this object More...
 
TraceRequest request_
 whether traces are being collected More...
 
int my_index_
 starting index of this object More...
 
Return_t new_value_
 a new value for a proposed move More...
 
Walker_tt_walker_
 reference to the current walker More...
 
bool streaming_particles_
 
bool have_required_traces_
 

Additional Inherited Members

- Protected Member Functions inherited from OperatorBase
virtual void contributeParticleQuantities ()
 
virtual void checkoutParticleQuantities (TraceManager &tm)
 
virtual void deleteParticleQuantities ()
 
virtual void setComputeForces (bool compute)
 
void setEnergyDomain (EnergyDomains edomain)
 Set the Energy Domain. More...
 
void setQuantumDomain (QuantumDomains qdomain)
 set quantum domain More...
 
void oneBodyQuantumDomain (const ParticleSet &P)
 set quantum domain for one-body operator More...
 
void twoBodyQuantumDomain (const ParticleSet &P)
 set quantum domain for two-body operator More...
 
void twoBodyQuantumDomain (const ParticleSet &P1, const ParticleSet &P2)
 set quantum domain for two-body operator More...
 
void addValue (PropertySetType &plist)
 named values to the property list More...
 

Detailed Description

Definition at line 23 of file DensityMatrices1B.h.

Member Typedef Documentation

◆ dens_t

using dens_t = std::vector<RealType>

Definition at line 42 of file DensityMatrices1B.h.

◆ Grad_t

using Grad_t = GradType

Definition at line 35 of file DensityMatrices1B.h.

◆ GradVector

Definition at line 37 of file DensityMatrices1B.h.

◆ Lattice_t

Definition at line 38 of file DensityMatrices1B.h.

◆ Matrix_t

Definition at line 40 of file DensityMatrices1B.h.

◆ pts_t

using pts_t = std::vector<PosType>

Definition at line 41 of file DensityMatrices1B.h.

◆ Value_t

using Value_t = ValueType

Definition at line 34 of file DensityMatrices1B.h.

◆ ValueVector

Definition at line 36 of file DensityMatrices1B.h.

◆ Vector_t

Definition at line 39 of file DensityMatrices1B.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
Enumerator
DIM 

Definition at line 29 of file DensityMatrices1B.h.

30  {
31  DIM = OHMMS_DIM
32  };
#define OHMMS_DIM
Definition: config.h:64

◆ evaluators

◆ integrators

◆ samplings

Constructor & Destructor Documentation

◆ DensityMatrices1B() [1/2]

Definition at line 48 of file DensityMatrices1B.cpp.

References DensityMatrices1B::reset().

50  basis_functions("DensityMatrices1B::basis"),
51  lattice_(P.getLattice()),
52  Psi(psi),
53  Pq(P),
54  Pc(Pcl)
55 {
56  reset();
57 }
static const TimerNameList_t< DMTimers > DMTimerNames
TimerManager< NewTimer > & getGlobalTimerManager()

◆ DensityMatrices1B() [2/2]

Definition at line 59 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::initialize(), DensityMatrices1B::reset(), and DensityMatrices1B::set_state().

60  : OperatorBase(master),
62  basis_functions(master.basis_functions),
63  lattice_(P.getLattice()),
64  Psi(psi),
65  Pq(P),
66  Pc(master.Pc)
67 {
68  reset();
69  set_state(master);
70  initialize();
71  for (int i = 0; i < basis_size; ++i)
72  basis_norms[i] = master.basis_norms[i];
73 }
OperatorBase()
Construct a new Operator Base object Default and unique empty constructor.
static const TimerNameList_t< DMTimers > DMTimerNames
TimerManager< NewTimer > & getGlobalTimerManager()

◆ ~DensityMatrices1B()

Member Function Documentation

◆ addObservables()

void addObservables ( PropertySetType plist,
BufferType collectables 
)
overridevirtual

named values to the property list Default implementaton uses addValue(plist_)

Parameters
plistRecordNameProperty
collectablesObservables that are accumulated by evaluate

Reimplemented from OperatorBase.

Definition at line 522 of file DensityMatrices1B.cpp.

References PooledData< T >::add(), DensityMatrices1B::basis_size, PooledData< T >::current(), DensityMatrices1B::eindex, DensityMatrices1B::energy_mat, OperatorBase::my_index_, DensityMatrices1B::nindex, and DensityMatrices1B::nspecies.

523 {
524 #if defined(QMC_COMPLEX)
525  int nentries = 2 * basis_size * basis_size * nspecies;
526 #else
527  int nentries = basis_size * basis_size * nspecies;
528 #endif
529  my_index_ = collectables.current();
530  nindex = my_index_;
531  std::vector<RealType> ntmp(nentries);
532  collectables.add(ntmp.begin(), ntmp.end());
533  if (energy_mat)
534  {
535  eindex = nindex + nentries;
536  std::vector<RealType> etmp(nentries);
537  collectables.add(etmp.begin(), etmp.end());
538  }
539 }
int my_index_
starting index of this object
Definition: OperatorBase.h:535

◆ checkoutScalarQuantities()

void checkoutScalarQuantities ( TraceManager tm)
inlineoverridevirtual

Reimplemented from OperatorBase.

Definition at line 171 of file DensityMatrices1B.h.

171 {}

◆ collectScalarQuantities()

void collectScalarQuantities ( )
inlineoverridevirtual

Reimplemented from OperatorBase.

Definition at line 172 of file DensityMatrices1B.h.

172 {}

◆ compare() [1/2]

void compare ( const std::string &  name,
Vector_t v1,
Vector_t v2,
bool  write = false,
bool  diff_only = true 
)

Definition at line 1424 of file DensityMatrices1B.cpp.

References qmcplusplus::app_log(), qmcplusplus::real(), DensityMatrices1B::same(), and Vector< T, Alloc >::size().

Referenced by DensityMatrices1B::evaluate_matrix().

1425 {
1426  bool sm = same(v1, v2);
1427  std::string result = "differ";
1428  if (sm)
1429  result = "agree";
1430  app_log() << name << " " << result << std::endl;
1431  if (write && !sm)
1432  for (int i = 0; i < v1.size(); ++i)
1433  app_log() << " " << i << " " << real(v1[i]) << " " << real(v2[i]) << " " << real(v1[i] / v2[i]) << " "
1434  << real(v2[i] / v1[i]) << std::endl;
1435 }
std::ostream & app_log()
Definition: OutputManager.h:65
float real(const float &c)
real part of a scalar. Cannot be replaced by std::real due to AFQMC specific needs.
bool same(Vector_t &v1, Vector_t &v2, RealType tol=1e-6)

◆ compare() [2/2]

void compare ( const std::string &  name,
Matrix_t m1,
Matrix_t m2,
bool  write = false,
bool  diff_only = true 
)

Definition at line 1437 of file DensityMatrices1B.cpp.

References qmcplusplus::app_log(), Matrix< T, Alloc >::cols(), DensityMatrices1B::match(), qmcplusplus::real(), Matrix< T, Alloc >::rows(), and DensityMatrices1B::same().

1438 {
1439  bool sm = same(m1, m2);
1440  std::string result = "differ";
1441  if (sm)
1442  result = "agree";
1443  app_log() << name << " " << result << std::endl;
1444  if (write && !sm)
1445  for (int i = 0; i < m1.rows(); ++i)
1446  for (int j = 0; j < m1.cols(); ++j)
1447  if (!diff_only || !match(m1(i, j), m2(i, j)))
1448  app_log() << " " << i << " " << j << " " << real(m1(i, j)) << " " << real(m2(i, j)) << " "
1449  << real(m1(i, j) / m2(i, j)) << std::endl;
1450 }
std::ostream & app_log()
Definition: OutputManager.h:65
float real(const float &c)
real part of a scalar. Cannot be replaced by std::real due to AFQMC specific needs.
bool match(Value_t e1, Value_t e2, RealType tol=1e-12)
bool same(Vector_t &v1, Vector_t &v2, RealType tol=1e-6)

◆ contributeScalarQuantities()

void contributeScalarQuantities ( )
inlineoverridevirtual

Reimplemented from OperatorBase.

Definition at line 170 of file DensityMatrices1B.h.

170 {}

◆ deleteScalarQuantities()

void deleteScalarQuantities ( )
inlineoverridevirtual

Reimplemented from OperatorBase.

Definition at line 173 of file DensityMatrices1B.h.

173 {}

◆ density_drift()

void density_drift ( const PosType r,
RealType dens,
PosType drift 
)
inline

Definition at line 1068 of file DensityMatrices1B.cpp.

References qmcplusplus::abs(), DensityMatrices1B::basis_gradients, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, qmcplusplus::conj(), DensityMatrices1B::dens, DensityMatrices1B::DIM, DensityMatrices1B::drift, DensityMatrices1B::timestep, and DensityMatrices1B::update_basis_d012().

Referenced by DensityMatrices1B::generate_density_samples(), DensityMatrices1B::test_derivatives(), and DensityMatrices1B::warmup_sampling().

1069 {
1070  update_basis_d012(r);
1071  dens = 0.0;
1072  drift = 0.0;
1073  for (int i = 0; i < basis_size; ++i)
1074  {
1075  const Grad_t& bg = basis_gradients[i];
1076  Value_t b = basis_values[i];
1077  Value_t bc = qmcplusplus::conj(b);
1078  dens += std::abs(bc * b);
1079  for (int d = 0; d < DIM; ++d)
1080  drift[d] += std::real(bc * bg[d]);
1081  }
1082  drift *= timestep / dens;
1083  dens /= basis_size;
1084 }
QMCTraits::RealType real
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
void update_basis_d012(const PosType &r)

◆ density_only()

void density_only ( const PosType r,
RealType dens 
)
inline

Definition at line 1055 of file DensityMatrices1B.cpp.

References qmcplusplus::abs(), DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, qmcplusplus::conj(), DensityMatrices1B::dens, and DensityMatrices1B::update_basis().

Referenced by DensityMatrices1B::generate_density_samples().

1056 {
1057  update_basis(r);
1058  dens = 0.0;
1059  for (int i = 0; i < basis_size; ++i)
1060  {
1061  Value_t b = basis_values[i];
1062  dens += std::abs(qmcplusplus::conj(b) * b);
1063  }
1064  dens /= basis_size;
1065 }
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
void update_basis(const PosType &r)

◆ dependsOnWaveFunction()

bool dependsOnWaveFunction ( ) const
inlineoverridevirtual

return true if this operator depends on a wavefunction

Reimplemented from OperatorBase.

Definition at line 151 of file DensityMatrices1B.h.

151 { return true; }

◆ diffusion()

void diffusion ( RealType  sqt,
PosType diff 
)
inline

Definition at line 1048 of file DensityMatrices1B.cpp.

References qmcplusplus::assignGaussRand(), DensityMatrices1B::DIM, and DensityMatrices1B::uniform_random.

Referenced by DensityMatrices1B::generate_density_samples(), and DensityMatrices1B::warmup_sampling().

1049 {
1050  assignGaussRand(&diff[0], DIM, *uniform_random);
1051  diff *= sqt;
1052 }
void assignGaussRand(T *restrict a, unsigned n, RG &rng)
RandomBase< FullPrecRealType > * uniform_random

◆ evaluate()

DensityMatrices1B::Return_t evaluate ( ParticleSet P)
overridevirtual

Evaluate the local energy contribution of this component.

Parameters
Pinput configuration containing N particles
Returns
the value of the Hamiltonian component

Implements OperatorBase.

Definition at line 598 of file DensityMatrices1B.cpp.

References APP_ABORT, DensityMatrices1B::check_derivatives, qmcplusplus::DM_eval, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluate_loop(), DensityMatrices1B::evaluate_matrix(), DensityMatrices1B::evaluator, OperatorBase::have_required_traces_, DensityMatrices1B::loop, DensityMatrices1B::matrix, DensityMatrices1B::test_derivatives(), and DensityMatrices1B::timers.

599 {
602  {
603  if (check_derivatives)
605  if (evaluator == loop)
606  evaluate_loop(P);
607  else if (evaluator == matrix)
608  evaluate_matrix(P);
609  else
610  APP_ABORT("DensityMatrices1B::evaluate invalid evaluator");
611  }
612  return 0.0;
613 }
ScopeGuard< NewTimer > ScopedTimer
Definition: NewTimer.h:257
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
Return_t evaluate_loop(ParticleSet &P)
Return_t evaluate_matrix(ParticleSet &P)

◆ evaluate_check()

DensityMatrices1B::Return_t evaluate_check ( ParticleSet P)

Definition at line 769 of file DensityMatrices1B.cpp.

References APP_ABORT, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, TrialWaveFunction::calcRatio(), qmcplusplus::conj(), DensityMatrices1B::E_trace, DensityMatrices1B::energy_mat, DensityMatrices1B::integrated_values, qmcplusplus::Units::distance::m, ParticleSet::makeMove(), qmcplusplus::n, qmcplusplus::Units::time::ns, DensityMatrices1B::nspecies, DensityMatrices1B::Psi, ParticleSet::R, ParticleSet::rejectMove(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, TraceSample< T >::sample, DensityMatrices1B::sample_weights, DensityMatrices1B::samples, DensityMatrices1B::species_size, and DensityMatrices1B::update_basis().

Referenced by DensityMatrices1B::evaluate_matrix().

770 {
771 #ifdef DMCHECK
772  APP_ABORT("DensityMatrices1B::evaluate_check use of E_trace in this function needs to be replaces with "
773  "get_energies() and E_samp");
774  int n = 0;
775  for (int s = 0; s < nspecies; ++s)
776  {
777  Matrix_t& Phi_mb = Phi_MBtmp;
778  Matrix_t& Psi_nm = *Psi_NMtmp[s];
779  Matrix_t& Phi_Psi_nb = *Phi_Psi_NBtmp[s];
780  Matrix_t& Phi_nb = *Phi_NBtmp[s];
781  Vector_t& E_n = *E_Ntmp[s];
782  Matrix_t& N_bb = *N_BBtmp[s];
783  Matrix_t& E_bb = *E_BBtmp[s];
784 
785  for (int ij = 0; ij < basis_size * basis_size; ++ij)
786  N_bb(ij) = 0.0;
787  for (int ij = 0; ij < basis_size * basis_size; ++ij)
788  E_bb(ij) = 0.0;
789 
790  for (int ns = 0; ns < species_size[s]; ++ns, ++n)
791  {
792  std::fill(integrated_values.begin(), integrated_values.end(), 0.0);
793  for (int m = 0; m < samples; ++m)
794  {
795  PosType& rsamp = rsamples[m];
796  update_basis(rsamp);
797  PosType dr = rsamp - P.R[n];
798  P.makeMove(n, dr);
800  P.rejectMove(n);
801  for (int i = 0; i < basis_size; ++i)
802  {
803  integrated_values[i] += ratio * basis_values[i];
804  Phi_mb(m, i) = basis_values[i];
805  }
806  Psi_nm(ns, m) = ratio;
807  }
808  update_basis(P.R[n]);
809  for (int i = 0; i < basis_size; ++i)
810  Phi_Psi_nb(ns, i) = integrated_values[i];
811  for (int i = 0; i < basis_size; ++i)
812  Phi_nb(ns, i) = qmcplusplus::conj(basis_values[i]);
813  for (int i = 0; i < basis_size; ++i)
814  {
816  for (int j = 0; j < basis_size; ++j)
817  {
818  Value_t val = phi_i * integrated_values[j];
819  N_bb(i, j) += val;
820  }
821  }
822  if (energy_mat)
823  {
824  RealType e_n = E_trace->sample[n]; //replace this with traces access later
825  E_n[ns] = e_n;
826  for (int i = 0; i < basis_size; ++i)
827  {
828  Value_t ephi_i = e_n * qmcplusplus::conj(basis_values[i]);
829  Phi_nb(ns, i) = ephi_i;
830  for (int j = 0; j < basis_size; ++j)
831  {
832  Value_t val = ephi_i * integrated_values[j];
833  E_bb(i, j) += val;
834  }
835  }
836  }
837  }
838  }
839 #endif
840  return 0.0;
841 }
ValueType calcRatio(ParticleSet &P, int iat, ComputeType ct=ComputeType::ALL)
compute psi(R_new) / psi(R_current) ratio It returns a complex value if the wavefunction is complex...
QMCTraits::PosType PosType
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
QMCTraits::RealType RealType
CombinedTraceSample< TraceReal > * E_trace
void update_basis(const PosType &r)
std::vector< PosType > rsamples

◆ evaluate_loop()

DensityMatrices1B::Return_t evaluate_loop ( ParticleSet P)

Definition at line 844 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, ParticleSet::Collectables, qmcplusplus::conj(), DensityMatrices1B::E_trace, DensityMatrices1B::eindex, DensityMatrices1B::energy_mat, DensityMatrices1B::generate_samples(), ParticleSet::getTotalNum(), qmcplusplus::imag(), DensityMatrices1B::integrate(), DensityMatrices1B::integrated_values, DensityMatrices1B::metric, qmcplusplus::n, DensityMatrices1B::nindex, DensityMatrices1B::nparticles, qmcplusplus::Units::time::ns, DensityMatrices1B::nspecies, ParticleSet::R, qmcplusplus::real(), qmcplusplus::Units::time::s, TraceSample< T >::sample, DensityMatrices1B::species_size, OperatorBase::t_walker_, DensityMatrices1B::update_basis(), DensityMatrices1B::w_trace, DensityMatrices1B::warmed_up, DensityMatrices1B::warmup_sampling(), and Walker< t_traits, p_traits >::Weight.

Referenced by DensityMatrices1B::evaluate().

845 {
846  const int basis_size2 = basis_size * basis_size;
847  if (!warmed_up)
848  warmup_sampling();
849  RealType weight;
850  if (energy_mat)
851  weight = w_trace->sample[0] * metric;
852  else
853  weight = t_walker_->Weight * metric;
854  int nparticles = P.getTotalNum();
855  generate_samples(weight);
856  int n = 0;
857  for (int s = 0; s < nspecies; ++s)
858  {
859  for (int ns = 0; ns < species_size[s]; ++ns, ++n)
860  {
861  integrate(P, n);
862  update_basis(P.R[n]);
863  int ij = nindex + s * basis_size2;
864  for (int i = 0; i < basis_size; ++i)
865  {
867  for (int j = 0; j < basis_size; ++j)
868  {
869  Value_t val = phi_i * integrated_values[j];
870  P.Collectables[ij] += real(val);
871  ij++;
872 #if defined(QMC_COMPLEX)
873  P.Collectables[ij] += imag(val);
874  ij++;
875 #endif
876  }
877  }
878  if (energy_mat)
879  {
880  RealType e_n = E_trace->sample[n]; //replace this with traces access later
881  int ij = eindex + s * basis_size2;
882  for (int i = 0; i < basis_size; ++i)
883  {
884  Value_t ephi_i = e_n * qmcplusplus::conj(basis_values[i]);
885  for (int j = 0; j < basis_size; ++j)
886  {
887  Value_t val = ephi_i * integrated_values[j];
888  P.Collectables[ij] += real(val);
889  ij++;
890 #if defined(QMC_COMPLEX)
891  P.Collectables[ij] += imag(val);
892  ij++;
893 #endif
894  }
895  }
896  }
897  }
898  }
899  return 0.0;
900 }
float real(const float &c)
real part of a scalar. Cannot be replaced by std::real due to AFQMC specific needs.
void integrate(ParticleSet &P, int n)
TraceSample< TraceReal > * w_trace
float imag(const float &c)
imaginary part of a scalar. Cannot be replaced by std::imag due to AFQMC specific needs...
Walker_t * t_walker_
reference to the current walker
Definition: OperatorBase.h:541
void generate_samples(RealType weight, int steps=0)
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
QMCTraits::RealType RealType
CombinedTraceSample< TraceReal > * E_trace
FullPrecRealType Weight
Weight of the walker.
Definition: Walker.h:102
void update_basis(const PosType &r)

◆ evaluate_matrix()

DensityMatrices1B::Return_t evaluate_matrix ( ParticleSet P)

Definition at line 616 of file DensityMatrices1B.cpp.

References APP_ABORT, qmcplusplus::app_log(), DensityMatrices1B::basis_size, ParticleSet::Collectables, DensityMatrices1B::compare(), qmcplusplus::MatrixOperators::diag_product(), qmcplusplus::DM_accumulate, qmcplusplus::DM_matrix_products, DensityMatrices1B::E_BB, DensityMatrices1B::E_N, DensityMatrices1B::eindex, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluate_check(), DensityMatrices1B::generate_particle_basis(), DensityMatrices1B::generate_sample_basis(), DensityMatrices1B::generate_sample_ratios(), DensityMatrices1B::generate_samples(), DensityMatrices1B::get_energies(), qmcplusplus::imag(), DensityMatrices1B::metric, qmcplusplus::n, DensityMatrices1B::N_BB, DensityMatrices1B::nindex, DensityMatrices1B::nspecies, DensityMatrices1B::Phi_MB, DensityMatrices1B::Phi_NB, DensityMatrices1B::Phi_Psi_NB, qmcplusplus::MatrixOperators::product(), qmcplusplus::MatrixOperators::product_AtB(), DensityMatrices1B::Psi_NM, qmcplusplus::real(), DensityMatrices1B::report(), qmcplusplus::Units::time::s, TraceSample< T >::sample, DensityMatrices1B::sample_weights, OperatorBase::t_walker_, DensityMatrices1B::timers, DensityMatrices1B::w_trace, DensityMatrices1B::warmed_up, DensityMatrices1B::warmup_sampling(), and Walker< t_traits, p_traits >::Weight.

Referenced by DensityMatrices1B::evaluate().

617 {
618  //perform warmup sampling the first time
619  if (!warmed_up)
620  warmup_sampling();
621  // get weight and single particle energy trace data
622  RealType weight;
623  if (energy_mat)
624  weight = w_trace->sample[0] * metric;
625  else
626  weight = t_walker_->Weight * metric;
627 
628  if (energy_mat)
629  get_energies(E_N); // energies : particles x 1
630  // compute sample positions (monte carlo or deterministic)
631  generate_samples(weight);
632  // compute basis and wavefunction ratio values in matrix form
633  generate_sample_basis(Phi_MB); // basis : samples x basis_size
634  generate_sample_ratios(Psi_NM); // conj(Psi ratio) : particles x samples
635  generate_particle_basis(P, Phi_NB); // conj(basis) : particles x basis_size
636  // perform integration via matrix products
637  {
638  ScopedTimer local_timer(timers[DM_matrix_products]);
639  for (int s = 0; s < nspecies; ++s)
640  {
641  Matrix_t& Psi_nm = *Psi_NM[s];
642  Matrix_t& Phi_Psi_nb = *Phi_Psi_NB[s];
643  Matrix_t& Phi_nb = *Phi_NB[s];
644  diag_product(Psi_nm, sample_weights, Psi_nm);
645  product(Psi_nm, Phi_MB, Phi_Psi_nb); // ratio*basis : particles x basis_size
646  product_AtB(Phi_nb, Phi_Psi_nb, *N_BB[s]); // conj(basis)^T*ratio*basis : basis_size^2
647  if (energy_mat)
648  {
649  Vector_t& E = *E_N[s];
650  diag_product(E, Phi_nb, Phi_nb); // diag(energies)*qmcplusplus::conj(basis)
651  product_AtB(Phi_nb, Phi_Psi_nb, *E_BB[s]); // (energies*conj(basis))^T*ratio*basis
652  }
653  }
654  }
655  // accumulate data into collectables
656  {
657  ScopedTimer local_timer(timers[DM_accumulate]);
658  const int basis_size2 = basis_size * basis_size;
659  int ij = nindex;
660  for (int s = 0; s < nspecies; ++s)
661  {
662  //int ij=nindex; // for testing
663  const Matrix_t& NDM = *N_BB[s];
664  for (int n = 0; n < basis_size2; ++n)
665  {
666  Value_t val = NDM(n);
667  P.Collectables[ij] += real(val);
668  ij++;
669 #if defined(QMC_COMPLEX)
670  P.Collectables[ij] += imag(val);
671  ij++;
672 #endif
673  }
674  }
675  if (energy_mat)
676  {
677  int ij = eindex;
678  for (int s = 0; s < nspecies; ++s)
679  {
680  //int ij=eindex; // for testing
681  const Matrix_t& EDM = *E_BB[s];
682  for (int n = 0; n < basis_size2; ++n)
683  {
684  Value_t val = EDM(n);
685  P.Collectables[ij] += real(val);
686  ij++;
687 #if defined(QMC_COMPLEX)
688  P.Collectables[ij] += imag(val);
689  ij++;
690 #endif
691  }
692  }
693  }
694  }
695 
696 
697  // jtk come back to this
698  // why are matrices so similar across species?
699  // O atom, einspline, points=20, blocks=1, steps=1
700  //
701  //app_log()<<" ntraces = ";
702  //for(int s=0;s<nspecies;++s)
703  //{
704  // const Matrix_t& NDM = *N_BB[s];
705  // Value_t ntrace = 0.0;
706  // for(int i=0;i<basis_size;++i)
707  // ntrace+=NDM(i,i);
708  // app_log()<<ntrace<<" ";
709  //}
710  //app_log()<< std::endl;
711  //
712  //app_log()<<"nmatrices"<< std::endl;
713  //for(int s=0;s<nspecies;++s)
714  //{
715  // app_log()<<" species "<<s<< std::endl;
716  // const Matrix_t& NDM = *N_BB[s];
717  // for(int i=0;i<basis_size;++i)
718  // {
719  // for(int j=0;j<basis_size;++j)
720  // app_log()<<" "<<real(NDM(i,j));
721  // app_log()<< std::endl;
722  // }
723  //}
724  //
725  //app_log()<<"positions"<< std::endl;
726  //int p=0;
727  //for(int s=0;s<nspecies;++s)
728  //{
729  // app_log()<<" species "<<s<< std::endl;
730  // for(int ps=0;ps<species_size[s];++ps,++p)
731  // app_log()<<" "<<p<<" "<<P.R[p]-P.getLattice().Center<< std::endl;
732  //}
733  //
734  ////app_log()<<"basis_values"<< std::endl;
735  ////int p=0;
736  ////for(int s=0;s<nspecies;++s)
737  ////{
738  //// app_log()<<" species "<<s<< std::endl;
739  //// for(int ps=0;ps<species_size[s];++ps,++p)
740  //// app_log()<<" "<<p<<" "<<P.R[p]<< std::endl;
741  ////}
742 
743 
744 #ifdef DMCHECK
745  report();
746  app_log() << "DM Check" << std::endl;
747  evaluate_check(P);
748  compare(" Phi_MB", Phi_MB, Phi_MBtmp);
749  for (int s = 0; s < nspecies; ++s)
750  {
751  app_log() << " species " << s << std::endl;
752  compare(" E_N ", *E_N[s], *E_Ntmp[s]);
753  compare(" Phi_NB ", *Phi_NB[s], *Phi_NBtmp[s]);
754  compare(" Psi_NM ", *Psi_NM[s], *Psi_NMtmp[s]);
755  compare(" Phi_Psi_NB", *Phi_Psi_NB[s], *Phi_Psi_NBtmp[s], true);
756  compare(" N_BB ", *N_BB[s], *N_BBtmp[s], true);
757  if (energy_mat)
758  compare(" E_BB ", *E_BB[s], *E_BBtmp[s], true);
759  }
760  app_log() << "end DM Check" << std::endl;
761  APP_ABORT("DM Check");
762 #endif
763 
764 
765  return 0.0;
766 }
std::vector< Matrix_t * > Phi_NB
std::ostream & app_log()
Definition: OutputManager.h:65
std::vector< Vector_t * > E_N
std::vector< Matrix_t * > Psi_NM
Return_t evaluate_check(ParticleSet &P)
float real(const float &c)
real part of a scalar. Cannot be replaced by std::real due to AFQMC specific needs.
ScopeGuard< NewTimer > ScopedTimer
Definition: NewTimer.h:257
TraceSample< TraceReal > * w_trace
float imag(const float &c)
imaginary part of a scalar. Cannot be replaced by std::imag due to AFQMC specific needs...
std::vector< Matrix_t * > E_BB
void product(const Matrix< T > &A, const Matrix< T > &B, Matrix< T > &C)
static function to perform C=AB for real matrices
void get_energies(std::vector< Vector_t *> &E_n)
void generate_sample_basis(Matrix_t &Phi_mb)
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
Walker_t * t_walker_
reference to the current walker
Definition: OperatorBase.h:541
void generate_samples(RealType weight, int steps=0)
QMCTraits::RealType RealType
void report(const std::string &pad="")
FullPrecRealType Weight
Weight of the walker.
Definition: Walker.h:102
std::vector< Matrix_t * > Phi_Psi_NB
void compare(const std::string &name, Vector_t &v1, Vector_t &v2, bool write=false, bool diff_only=true)
void generate_particle_basis(ParticleSet &P, std::vector< Matrix_t *> &Phi_nb)
void generate_sample_ratios(std::vector< Matrix_t *> Psi_nm)
void product_AtB(const Matrix< T > &A, const Matrix< T > &B, Matrix< T > &C)
std::vector< Matrix_t * > N_BB
void diag_product(const Matrix< T1 > &A, const Vector< T2 > &B, Matrix< T3 > &C)
C = A*diag(B)

◆ finalize()

void finalize ( )

Definition at line 395 of file DensityMatrices1B.cpp.

References qmcplusplus::delete_iter(), DensityMatrices1B::E_BB, DensityMatrices1B::E_N, DensityMatrices1B::energy_mat, DensityMatrices1B::N_BB, DensityMatrices1B::Phi_NB, DensityMatrices1B::Phi_Psi_NB, and DensityMatrices1B::Psi_NM.

Referenced by DensityMatrices1B::~DensityMatrices1B().

396 {
397  delete_iter(Phi_NB.begin(), Phi_NB.end());
398  delete_iter(Psi_NM.begin(), Psi_NM.end());
399  delete_iter(Phi_Psi_NB.begin(), Phi_Psi_NB.end());
400  delete_iter(N_BB.begin(), N_BB.end());
401  if (energy_mat)
402  {
403  delete_iter(E_N.begin(), E_N.end());
404  delete_iter(E_BB.begin(), E_BB.end());
405  }
406 
407 #ifdef DMCHECK
408  delete_iter(Phi_NBtmp.begin(), Phi_NBtmp.end());
409  delete_iter(Psi_NMtmp.begin(), Psi_NMtmp.end());
410  delete_iter(Phi_Psi_NBtmp.begin(), Phi_Psi_NBtmp.end());
411  delete_iter(N_BBtmp.begin(), N_BBtmp.end());
412  if (energy_mat)
413  {
414  delete_iter(E_Ntmp.begin(), E_Ntmp.end());
415  delete_iter(E_BBtmp.begin(), E_BBtmp.end());
416  }
417 #endif
418 }
std::vector< Matrix_t * > Phi_NB
void delete_iter(IT first, IT last)
delete the pointers in [first,last)
std::vector< Vector_t * > E_N
std::vector< Matrix_t * > Psi_NM
std::vector< Matrix_t * > E_BB
std::vector< Matrix_t * > Phi_Psi_NB
std::vector< Matrix_t * > N_BB

◆ generate_density_samples()

void generate_density_samples ( bool  save,
int  steps,
RandomBase< FullPrecRealType > &  rng 
)
inline

Definition at line 996 of file DensityMatrices1B.cpp.

References DensityMatrices1B::acceptance_ratio, qmcplusplus::app_log(), DensityMatrices1B::density_drift(), DensityMatrices1B::density_only(), DensityMatrices1B::diffusion(), qmcplusplus::dot(), DensityMatrices1B::dpcur, qmcplusplus::exp(), qmcplusplus::n, DensityMatrices1B::naccepted, DensityMatrices1B::nmoves, omp_get_thread_num(), DensityMatrices1B::rhocur, DensityMatrices1B::rpcur, DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::sample_weights, qmcplusplus::sqrt(), DensityMatrices1B::timestep, DensityMatrices1B::use_drift, and DensityMatrices1B::write_acceptance_ratio.

Referenced by DensityMatrices1B::generate_samples().

997 {
998  RealType sqt = std::sqrt(timestep);
999  RealType ot = 1.0 / timestep;
1000  PosType r = rpcur; //current position
1001  PosType d = dpcur; //current drift
1002  RealType rho = rhocur; //current density
1003  for (int s = 0; s < steps; ++s)
1004  {
1005  nmoves++;
1006  PosType n, rp, dp, ds; //diffusion, trial pos/drift, drift sum
1007  RealType rhop, ratio, Pacc; //trial density, dens ratio, acc prob
1008  diffusion(sqt, n); //get diffusion
1009  if (use_drift)
1010  {
1011  rp = r + n + d; //update trial position
1012  density_drift(rp, rhop, dp); //get trial drift and density
1013  ratio = rhop / rho; //density ratio
1014  ds = dp + d; //drift sum
1015  Pacc = ratio * std::exp(-ot * (dot(n, ds) + .5 * dot(ds, ds))); //acceptance probability
1016  }
1017  else
1018  {
1019  rp = r + n; //update trial position
1020  density_only(rp, rhop); //get trial density
1021  ratio = rhop / rho; //density ratio
1022  Pacc = ratio; //acceptance probability
1023  }
1024  if (rng() < Pacc)
1025  { //accept move
1026  r = rp;
1027  d = dp;
1028  rho = rhop;
1029  naccepted++;
1030  }
1031  if (save)
1032  {
1033  rsamples[s] = r;
1034  sample_weights[s] = 1.0 / rho;
1035  }
1036  }
1038 
1040  app_log() << "dm1b acceptance_ratio = " << acceptance_ratio << std::endl;
1041 
1042  rpcur = r;
1043  dpcur = d;
1044  rhocur = rho;
1045 }
void density_drift(const PosType &r, RealType &dens, PosType &drift)
QTBase::RealType RealType
Definition: Configuration.h:58
std::ostream & app_log()
Definition: OutputManager.h:65
QMCTraits::PosType PosType
omp_int_t omp_get_thread_num()
Definition: OpenMP.h:25
void diffusion(RealType sqt, PosType &diff)
MakeReturn< UnaryNode< FnExp, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t exp(const Vector< T1, C1 > &l)
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)
QMCTraits::RealType RealType
Tensor< typename BinaryReturn< T1, T2, OpMultiply >::Type_t, D > dot(const AntiSymTensor< T1, D > &lhs, const AntiSymTensor< T2, D > &rhs)
void density_only(const PosType &r, RealType &dens)
std::vector< PosType > rsamples

◆ generate_particle_basis()

void generate_particle_basis ( ParticleSet P,
std::vector< Matrix_t *> &  Phi_nb 
)

Definition at line 1194 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, qmcplusplus::conj(), qmcplusplus::DM_gen_particle_basis, qmcplusplus::n, DensityMatrices1B::nspecies, ParticleSet::R, qmcplusplus::Units::time::s, DensityMatrices1B::species_size, DensityMatrices1B::timers, and DensityMatrices1B::update_basis().

Referenced by DensityMatrices1B::evaluate_matrix().

1195 {
1197  int p = 0;
1198  for (int s = 0; s < nspecies; ++s)
1199  {
1200  int nb = 0;
1201  Matrix_t& P_nb = *Phi_nb[s];
1202  for (int n = 0; n < species_size[s]; ++n, ++p)
1203  {
1204  update_basis(P.R[p]);
1205  for (int b = 0; b < basis_size; ++b, ++nb)
1206  P_nb(nb) = qmcplusplus::conj(basis_values[b]);
1207  }
1208  }
1209 }
ScopeGuard< NewTimer > ScopedTimer
Definition: NewTimer.h:257
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
void update_basis(const PosType &r)

◆ generate_sample_basis()

void generate_sample_basis ( Matrix_t Phi_mb)

◆ generate_sample_ratios()

void generate_sample_ratios ( std::vector< Matrix_t *>  Psi_nm)

Definition at line 1171 of file DensityMatrices1B.cpp.

References qmcplusplus::conj(), qmcplusplus::DM_gen_sample_ratios, TrialWaveFunction::evaluateRatiosAlltoOne(), qmcplusplus::Units::distance::m, ParticleSet::makeVirtualMoves(), qmcplusplus::n, DensityMatrices1B::nspecies, DensityMatrices1B::Pq, DensityMatrices1B::Psi, DensityMatrices1B::psi_ratios, DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::samples, DensityMatrices1B::species_size, and DensityMatrices1B::timers.

Referenced by DensityMatrices1B::evaluate_matrix().

1172 {
1174  for (int m = 0; m < samples; ++m)
1175  {
1176  // get N ratios for the current sample point
1179 
1180  // collect ratios into per-species matrices
1181  int p = 0;
1182  for (int s = 0; s < nspecies; ++s)
1183  {
1184  Matrix_t& P_nm = *Psi_nm[s];
1185  for (int n = 0; n < species_size[s]; ++n, ++p)
1186  {
1187  P_nm(n, m) = qmcplusplus::conj(psi_ratios[p]);
1188  }
1189  }
1190  }
1191 }
void evaluateRatiosAlltoOne(ParticleSet &P, std::vector< ValueType > &ratios)
ScopeGuard< NewTimer > ScopedTimer
Definition: NewTimer.h:257
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
std::vector< ValueType > psi_ratios
std::vector< PosType > rsamples
void makeVirtualMoves(const SingleParticlePos &newpos)
Handles virtual moves for all the particles to a single newpos.

◆ generate_samples()

void generate_samples ( RealType  weight,
int  steps = 0 
)
inline

Definition at line 903 of file DensityMatrices1B.cpp.

References qmcplusplus::app_log(), Vector< T, Alloc >::begin(), DensityMatrices1B::density, DensityMatrices1B::DIM, qmcplusplus::DM_gen_samples, Vector< T, Alloc >::end(), DensityMatrices1B::generate_density_samples(), DensityMatrices1B::generate_uniform_grid(), DensityMatrices1B::generate_uniform_samples(), DensityMatrices1B::integrator, DensityMatrices1B::metropolis, omptarget::min(), omp_get_thread_num(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::sample_weights, DensityMatrices1B::samples, DensityMatrices1B::sampling, qmcplusplus::sqrt(), DensityMatrices1B::timers, DensityMatrices1B::uniform, DensityMatrices1B::uniform_grid, DensityMatrices1B::uniform_random, and DensityMatrices1B::write_rstats.

Referenced by DensityMatrices1B::evaluate_loop(), DensityMatrices1B::evaluate_matrix(), DensityMatrices1B::test_derivatives(), and DensityMatrices1B::warmup_sampling().

904 {
906  auto& rng = *uniform_random;
907  bool save = false;
908  if (steps == 0)
909  {
910  save = true;
911  steps = samples;
912  }
913  if (integrator == uniform_grid)
915  else if (integrator == uniform)
917  else if (integrator == density)
918  generate_density_samples(save, steps, rng);
919 
920  if (save)
921  {
922  if (sampling == metropolis)
923  {
924  sample_weights *= weight;
925  }
926  else
927  {
928  std::fill(sample_weights.begin(), sample_weights.end(), weight);
929  }
930  }
931 
932 
933  // temporary check
934  if (write_rstats && omp_get_thread_num() == 0)
935  {
936  PosType rmin = std::numeric_limits<RealType>::max();
937  PosType rmax = -std::numeric_limits<RealType>::max();
938  PosType rmean = 0.0;
939  PosType rstd = 0.0;
940  for (int s = 0; s < rsamples.size(); ++s)
941  for (int d = 0; d < DIM; ++d)
942  {
943  RealType rd = rsamples[s][d];
944  rmin[d] = std::min(rmin[d], rd);
945  rmax[d] = std::max(rmax[d], rd);
946  rmean[d] += rd;
947  rstd[d] += rd * rd;
948  }
949  rmean /= rsamples.size();
950  rstd /= rsamples.size();
951  for (int d = 0; d < DIM; ++d)
952  rstd[d] = std::sqrt(rstd[d] - rmean[d] * rmean[d]);
953  app_log() << "\nrsamples properties:" << std::endl;
954  app_log() << " rmin = " << rmin << std::endl;
955  app_log() << " rmax = " << rmax << std::endl;
956  app_log() << " rmean = " << rmean << std::endl;
957  app_log() << " rstd = " << rstd << std::endl;
958  }
959 }
void generate_density_samples(bool save, int steps, RandomBase< FullPrecRealType > &rng)
std::ostream & app_log()
Definition: OutputManager.h:65
ScopeGuard< NewTimer > ScopedTimer
Definition: NewTimer.h:257
T min(T a, T b)
QMCTraits::PosType PosType
omp_int_t omp_get_thread_num()
Definition: OpenMP.h:25
void generate_uniform_grid(RandomBase< FullPrecRealType > &rng)
RandomBase< FullPrecRealType > * uniform_random
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)
QMCTraits::RealType RealType
void generate_uniform_samples(RandomBase< FullPrecRealType > &rng)
std::vector< PosType > rsamples

◆ generate_uniform_grid()

void generate_uniform_grid ( RandomBase< FullPrecRealType > &  rng)
inline

Definition at line 962 of file DensityMatrices1B.cpp.

References DensityMatrices1B::DIM, DensityMatrices1B::ind_dims, DensityMatrices1B::lattice_, DensityMatrices1B::points, DensityMatrices1B::rcorner, DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::samples, and DensityMatrices1B::scale.

Referenced by DensityMatrices1B::generate_samples().

963 {
964  PosType rp;
965  PosType ushift = 0.0;
966  RealType du = scale / points;
967  for (int d = 0; d < DIM; ++d)
968  ushift[d] += rng() * du;
969  for (int s = 0; s < samples; ++s)
970  {
971  int nrem = s;
972  for (int d = 0; d < DIM - 1; ++d)
973  {
974  int ind = nrem / ind_dims[d];
975  rp[d] = ind * du + ushift[d];
976  nrem -= ind * ind_dims[d];
977  }
978  rp[DIM - 1] = nrem * du + ushift[DIM - 1];
979  rsamples[s] = lattice_.toCart(rp) + rcorner;
980  }
981 }
QMCTraits::PosType PosType
QMCTraits::RealType RealType
std::vector< PosType > rsamples

◆ generate_uniform_samples()

void generate_uniform_samples ( RandomBase< FullPrecRealType > &  rng)
inline

◆ get()

bool get ( std::ostream &  os) const
inlineoverridevirtual

write about the class

Implements OperatorBase.

Definition at line 176 of file DensityMatrices1B.h.

176 { return false; }

◆ get_energies()

void get_energies ( std::vector< Vector_t *> &  E_n)

Definition at line 1129 of file DensityMatrices1B.cpp.

References qmcplusplus::accum_constant(), qmcplusplus::accum_sample(), DensityMatrices1B::E_samp, DensityMatrices1B::nparticles, DensityMatrices1B::nspecies, qmcplusplus::Units::time::ps, qmcplusplus::Units::time::s, Vector< T, Alloc >::size(), DensityMatrices1B::T_trace, DensityMatrices1B::Vc_trace, DensityMatrices1B::Vcc_trace, DensityMatrices1B::Vq_trace, DensityMatrices1B::Vqc_trace, and DensityMatrices1B::Vqq_trace.

Referenced by DensityMatrices1B::evaluate_matrix().

1130 {
1131  Value_t Vc = 0;
1132  Vc += accum_constant(Vc_trace);
1133  Vc += accum_constant(Vcc_trace);
1134  Vc /= nparticles;
1135  std::fill(E_samp.begin(), E_samp.end(), Vc);
1139  accum_sample(E_samp, Vqc_trace, 2.0);
1140 
1141  int p = 0;
1142  for (int s = 0; s < nspecies; ++s)
1143  {
1144  Vector_t& E = *E_n[s];
1145  for (int ps = 0; ps < E.size(); ++ps, ++p)
1146  E[ps] = E_samp[p];
1147  }
1148 
1149  //E_trace->combine();
1150  //RealType E = 0.0;
1151  //for(int p=0;p<nparticles;++p)
1152  // E += E_samp[p];
1153  //app_log()<<" E = "<<E<<" "<<E_trace->sample[0]<< std::endl;
1154  //APP_ABORT("dm1b::get_energies check sp traces");
1155 }
std::vector< Value_t > E_samp
RealType accum_constant(CombinedTraceSample< TraceReal > *etrace, RealType weight=1.0)
TraceSample< TraceComp > * T_trace
void accum_sample(std::vector< Value_t > &E_samp, CombinedTraceSample< T > *etrace, RealType weight=1.0)
CombinedTraceSample< TraceReal > * Vc_trace
CombinedTraceSample< TraceReal > * Vq_trace
CombinedTraceSample< TraceReal > * Vcc_trace
CombinedTraceSample< TraceReal > * Vqq_trace
CombinedTraceSample< TraceReal > * Vqc_trace

◆ getClassName()

std::string getClassName ( ) const
inlineoverridevirtual

return class name

Implements OperatorBase.

Definition at line 152 of file DensityMatrices1B.h.

152 { return "DensityMatrices1B"; }

◆ getRequiredTraces()

void getRequiredTraces ( TraceManager tm)
overridevirtual

TODO: add docs.

Parameters
tm

Reimplemented from OperatorBase.

Definition at line 498 of file DensityMatrices1B.cpp.

References DensityMatrices1B::E_samp, DensityMatrices1B::energy_mat, TraceManager::get_complex_trace(), TraceManager::get_real_combined_trace(), TraceManager::get_real_trace(), OperatorBase::have_required_traces_, DensityMatrices1B::nparticles, DensityMatrices1B::Pc, DensityMatrices1B::Pq, DensityMatrices1B::T_trace, DensityMatrices1B::Vc_trace, DensityMatrices1B::Vcc_trace, DensityMatrices1B::Vq_trace, DensityMatrices1B::Vqc_trace, DensityMatrices1B::Vqq_trace, and DensityMatrices1B::w_trace.

499 {
500  w_trace = tm.get_real_trace("weight");
501  if (energy_mat)
502  {
503  T_trace = tm.get_complex_trace(Pq, "Kinetic_complex");
504  Vq_trace = tm.get_real_combined_trace(Pq, "Vq");
505  Vqq_trace = tm.get_real_combined_trace(Pq, "Vqq");
506  Vqc_trace = tm.get_real_combined_trace(Pq, "Vqc");
507  if (Pc)
508  {
509  Vc_trace = tm.get_real_combined_trace(*Pc, "Vc");
510  Vcc_trace = tm.get_real_combined_trace(*Pc, "Vcc");
511  }
512 
513  E_samp.resize(nparticles);
514  }
515  have_required_traces_ = true;
516 }
std::vector< Value_t > E_samp
TraceSample< TraceReal > * w_trace
TraceSample< TraceComp > * T_trace
CombinedTraceSample< TraceReal > * Vc_trace
CombinedTraceSample< TraceReal > * Vq_trace
CombinedTraceSample< TraceReal > * Vcc_trace
CombinedTraceSample< TraceReal > * Vqq_trace
CombinedTraceSample< TraceReal > * Vqc_trace

◆ initialize()

void initialize ( )

Definition at line 319 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_gradients, DensityMatrices1B::basis_laplacians, DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, DensityMatrices1B::E_BB, DensityMatrices1B::E_N, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluator, ParticleSet::getSpeciesSet(), ParticleSet::getTotalNum(), ParticleSet::groupsize(), DensityMatrices1B::initialized, DensityMatrices1B::integrated_values, DensityMatrices1B::matrix, DensityMatrices1B::metropolis, DensityMatrices1B::N_BB, DensityMatrices1B::normalize(), DensityMatrices1B::normalized, DensityMatrices1B::nparticles, DensityMatrices1B::nspecies, SpeciesSet::numAttributes(), DensityMatrices1B::Phi_MB, DensityMatrices1B::Phi_NB, DensityMatrices1B::Phi_Psi_NB, DensityMatrices1B::Pq, DensityMatrices1B::Psi_NM, DensityMatrices1B::psi_ratios, Matrix< T, Alloc >::resize(), Vector< T, Alloc >::resize(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::sample_weights, DensityMatrices1B::samples, DensityMatrices1B::sampling, SpeciesSet::size(), DensityMatrices1B::species_name, DensityMatrices1B::species_size, SpeciesSet::speciesName, qmcplusplus::sqrt(), DensityMatrices1B::volume, and DensityMatrices1B::volume_normed.

Referenced by DensityMatrices1B::DensityMatrices1B(), and DensityMatrices1B::put().

320 {
321  // get particle information
322  SpeciesSet& species = Pq.getSpeciesSet();
324  nspecies = species.size();
325  int natt = species.numAttributes();
326  for (int s = 0; s < nspecies; ++s)
327  species_size.push_back(Pq.groupsize(s));
328  for (int s = 0; s < nspecies; ++s)
329  species_name.push_back(species.speciesName[s]);
330 
331  // allocate space
332  basis_values.resize(basis_size);
334  basis_norms.resize(basis_size);
335  RealType bn_standard = 1.0;
336  if (volume_normed)
337  bn_standard = 1.0 / std::sqrt(volume);
338  for (int i = 0; i < basis_size; ++i)
339  basis_norms[i] = bn_standard;
340 
341  rsamples.resize(samples);
343  psi_ratios.resize(nparticles);
344 
345  if (evaluator == matrix)
346  {
348  for (int s = 0; s < nspecies; ++s)
349  {
350  int specs_size = species_size[s];
351  Phi_NB.push_back(new Matrix_t(specs_size, basis_size));
352  Psi_NM.push_back(new Matrix_t(specs_size, samples));
353  Phi_Psi_NB.push_back(new Matrix_t(specs_size, basis_size));
354  N_BB.push_back(new Matrix_t(basis_size, basis_size));
355  if (energy_mat)
356  {
357  E_N.push_back(new Vector_t(specs_size));
358  E_BB.push_back(new Matrix_t(basis_size, basis_size));
359  }
360  }
361 
362 #ifdef DMCHECK
363  Phi_MBtmp.resize(samples, basis_size);
364  for (int s = 0; s < nspecies; ++s)
365  {
366  int specs_size = species_size[s];
367  Phi_NBtmp.push_back(new Matrix_t(specs_size, basis_size));
368  Psi_NMtmp.push_back(new Matrix_t(specs_size, samples));
369  Phi_Psi_NBtmp.push_back(new Matrix_t(specs_size, basis_size));
370  N_BBtmp.push_back(new Matrix_t(basis_size, basis_size));
371  if (energy_mat)
372  {
373  E_Ntmp.push_back(new Vector_t(specs_size));
374  E_BBtmp.push_back(new Matrix_t(basis_size, basis_size));
375  }
376  }
377 #endif
378  }
379 
380  if (sampling == metropolis)
381  {
382  basis_gradients.resize(basis_size);
384  }
385 
386  if (!normalized)
387  {
388  normalize();
389  }
390 
391  initialized = true;
392 }
void resize(size_type n, Type_t val=Type_t())
Resize the container.
Definition: OhmmsVector.h:166
std::vector< Matrix_t * > Phi_NB
size_t getTotalNum() const
Definition: ParticleSet.h:493
std::vector< Vector_t * > E_N
std::vector< Matrix_t * > Psi_NM
void resize(size_type n, size_type m)
Resize the container.
Definition: OhmmsMatrix.h:99
std::vector< Matrix_t * > E_BB
int groupsize(int igroup) const
return the size of a group
Definition: ParticleSet.h:527
std::vector< std::string > species_name
SpeciesSet & getSpeciesSet()
retrun the SpeciesSet of this particle set
Definition: ParticleSet.h:231
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)
QMCTraits::RealType RealType
std::vector< Matrix_t * > Phi_Psi_NB
std::vector< ValueType > psi_ratios
std::vector< Matrix_t * > N_BB
std::vector< PosType > rsamples

◆ integrate()

void integrate ( ParticleSet P,
int  n 
)
inline

Definition at line 1212 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, TrialWaveFunction::calcRatio(), qmcplusplus::conj(), DensityMatrices1B::integrated_values, ParticleSet::makeMove(), qmcplusplus::n, DensityMatrices1B::Psi, ParticleSet::R, ParticleSet::rejectMove(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::sample_weights, DensityMatrices1B::samples, and DensityMatrices1B::update_basis().

Referenced by DensityMatrices1B::evaluate_loop().

1213 {
1214  std::fill(integrated_values.begin(), integrated_values.end(), 0.0);
1215  for (int s = 0; s < samples; ++s)
1216  {
1217  PosType& rsamp = rsamples[s];
1218  update_basis(rsamp);
1219  P.makeMove(n, rsamp - P.R[n]);
1221  P.rejectMove(n);
1222  for (int i = 0; i < basis_size; ++i)
1223  integrated_values[i] += ratio * basis_values[i];
1224  }
1225 }
ValueType calcRatio(ParticleSet &P, int iat, ComputeType ct=ComputeType::ALL)
compute psi(R_new) / psi(R_current) ratio It returns a complex value if the wavefunction is complex...
QMCTraits::PosType PosType
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
void update_basis(const PosType &r)
std::vector< PosType > rsamples

◆ makeClone()

std::unique_ptr< OperatorBase > makeClone ( ParticleSet P,
TrialWaveFunction psi 
)
finalvirtual

Implements OperatorBase.

Definition at line 83 of file DensityMatrices1B.cpp.

84 {
85  return std::make_unique<DensityMatrices1B>(*this, P, psi);
86 }

◆ match()

bool match ( Value_t  e1,
Value_t  e2,
RealType  tol = 1e-12 
)

Definition at line 1396 of file DensityMatrices1B.cpp.

References qmcplusplus::abs().

Referenced by DensityMatrices1B::compare(), and DensityMatrices1B::same().

1397 {
1398  return std::abs(e1 - e2) < tol;
1399  //return std::abs(2*(e1-e2)/(e1+e2)) < tol;
1400 }
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)

◆ normalize()

void normalize ( )
inline

Definition at line 1252 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, qmcplusplus::conj(), DensityMatrices1B::DIM, DensityMatrices1B::lattice_, DensityMatrices1B::normalized, DensityMatrices1B::points, qmcplusplus::pow(), DensityMatrices1B::rcorner, qmcplusplus::real(), DensityMatrices1B::scale, qmcplusplus::sqrt(), DensityMatrices1B::update_basis(), and DensityMatrices1B::volume.

Referenced by DensityMatrices1B::initialize().

1253 {
1254  int ngrid = std::max(200, points);
1255  int ngtot = pow(ngrid, DIM);
1256  RealType du = scale / ngrid;
1257  RealType dV = volume / ngtot;
1258  PosType rp;
1259  ValueVector bnorms;
1260  int gdims[DIM];
1261  gdims[0] = pow(ngrid, DIM - 1);
1262  for (int d = 1; d < DIM; ++d)
1263  gdims[d] = gdims[d - 1] / ngrid;
1264  bnorms.resize(basis_size);
1265  for (int i = 0; i < basis_size; ++i)
1266  bnorms[i] = 0.0;
1267  std::fill(basis_norms.begin(), basis_norms.end(), 1.0);
1268  for (int p = 0; p < ngtot; ++p)
1269  {
1270  int nrem = p;
1271  for (int d = 0; d < DIM - 1; ++d)
1272  {
1273  int ind = nrem / gdims[d];
1274  rp[d] = ind * du + du / 2;
1275  nrem -= ind * gdims[d];
1276  }
1277  rp[DIM - 1] = nrem * du + du / 2;
1278  rp = lattice_.toCart(rp) + rcorner;
1279  update_basis(rp);
1280  for (int i = 0; i < basis_size; ++i)
1281  bnorms[i] += qmcplusplus::conj(basis_values[i]) * basis_values[i] * dV;
1282  }
1283  for (int i = 0; i < basis_size; ++i)
1284  basis_norms[i] = 1.0 / std::sqrt(real(bnorms[i]));
1285  normalized = true;
1286 }
float real(const float &c)
real part of a scalar. Cannot be replaced by std::real due to AFQMC specific needs.
QMCTraits::PosType PosType
MakeReturn< BinaryNode< FnPow, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t, typename CreateLeaf< Vector< T2, C2 > >::Leaf_t > >::Expression_t pow(const Vector< T1, C1 > &l, const Vector< T2, C2 > &r)
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)
QMCTraits::RealType RealType
void update_basis(const PosType &r)

◆ put()

bool put ( xmlNodePtr  cur)
overridevirtual

Read the input parameter.

Parameters
curxml node for a OperatorBase object

Implements OperatorBase.

Definition at line 138 of file DensityMatrices1B.cpp.

References DensityMatrices1B::check_overlap, DensityMatrices1B::initialize(), DensityMatrices1B::report(), DensityMatrices1B::set_state(), and DensityMatrices1B::test_overlap().

139 {
140  // read in parameters from the input xml
141  set_state(cur);
142 
143  // resize local data and perform warmup sampling, if necessary
144  initialize();
145 
146  // write a report to the log
147  report();
148 
149  if (check_overlap)
150  test_overlap();
151 
152  return true;
153 }
void report(const std::string &pad="")

◆ registerCollectables()

void registerCollectables ( std::vector< ObservableHelper > &  h5desc,
hdf_archive file 
) const
overridevirtual

Reimplemented from OperatorBase.

Definition at line 542 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::eindex, DensityMatrices1B::energy_mat, OperatorBase::name_, DensityMatrices1B::nindex, DensityMatrices1B::nspecies, qmcplusplus::oh, qmcplusplus::Units::time::s, ObservableHelper::set_dimensions(), and DensityMatrices1B::species_name.

543 {
544 #if defined(QMC_COMPLEX)
545  std::vector<int> ng(3);
546  ng[0] = basis_size;
547  ng[1] = basis_size;
548  ng[2] = 2;
549  int nentries = ng[0] * ng[1] * ng[2];
550 #else
551  std::vector<int> ng(2);
552  ng[0] = basis_size;
553  ng[1] = basis_size;
554  int nentries = ng[0] * ng[1];
555 #endif
556 
557  hdf_path hdf_name{name_};
558  hdf_name /= "number_matrix";
559  for (int s = 0; s < nspecies; ++s)
560  {
561  h5desc.emplace_back(hdf_name / species_name[s]);
562  auto& oh = h5desc.back();
563  oh.set_dimensions(ng, nindex + s * nentries);
564  }
565 
566  if (energy_mat)
567  {
568  hdf_name.replace_subgroup("energy_matrix");
569  for (int s = 0; s < nspecies; ++s)
570  {
571  h5desc.emplace_back(hdf_name / species_name[s]);
572  auto& oh = h5desc.back();
573  oh.set_dimensions(ng, eindex + s * nentries);
574  }
575  }
576 }
std::string name_
name of this object
Definition: OperatorBase.h:527
std::vector< std::string > species_name
void set_dimensions(const std::vector< int > &dims, int first)
set the shape of this observable
ObservableHelper oh

◆ report()

void report ( const std::string &  pad = "")

Definition at line 421 of file DensityMatrices1B.cpp.

References qmcplusplus::app_log(), DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::center, Matrix< T, Alloc >::cols(), DensityMatrices1B::E_BB, DensityMatrices1B::E_N, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluator, DensityMatrices1B::initialized, DensityMatrices1B::integrator, DensityMatrices1B::lattice_, DensityMatrices1B::matrix, DensityMatrices1B::metric, DensityMatrices1B::metropolis, DensityMatrices1B::N_BB, DensityMatrices1B::normalized, DensityMatrices1B::nparticles, DensityMatrices1B::nspecies, DensityMatrices1B::periodic, DensityMatrices1B::Phi_MB, DensityMatrices1B::Phi_NB, DensityMatrices1B::Phi_Psi_NB, DensityMatrices1B::points, DensityMatrices1B::Psi_NM, DensityMatrices1B::rcorner, Matrix< T, Alloc >::rows(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::samples, DensityMatrices1B::sampling, DensityMatrices1B::scale, DensityMatrices1B::species_size, DensityMatrices1B::timestep, DensityMatrices1B::volume, DensityMatrices1B::volume_based, DensityMatrices1B::volume_normed, and DensityMatrices1B::warmup.

Referenced by DensityMatrices1B::evaluate_matrix(), and DensityMatrices1B::put().

422 {
423  std::vector<std::string> integrator_list;
424  integrator_list.push_back("uniform_grid");
425  integrator_list.push_back("uniform");
426  integrator_list.push_back("density");
427  integrator_list.push_back("no_integrator");
428  std::vector<std::string> evaluator_list;
429  evaluator_list.push_back("loop");
430  evaluator_list.push_back("matrix");
431  evaluator_list.push_back("no_evaluator");
432  std::vector<std::string> sampling_list;
433  sampling_list.push_back("volume_based");
434  sampling_list.push_back("metropolis");
435  sampling_list.push_back("no_sampling");
436 
437  std::ostream& out = app_log();
438  //std::ostream& out = std::cout;
439 
440  out << pad << "DensityMatrices1B" << std::endl;
441 
442  out << pad << " integrator = " << integrator_list[(int)integrator] << std::endl;
443  out << pad << " sampling = " << sampling_list[(int)sampling] << std::endl;
444  out << pad << " evaluator = " << evaluator_list[(int)evaluator] << std::endl;
445  out << pad << " periodic = " << periodic << std::endl;
446  if (sampling == volume_based)
447  {
448  PosType rmax = rcorner + 2 * scale * lattice_.Center;
449  out << pad << " points = " << points << std::endl;
450  out << pad << " scale = " << scale << std::endl;
451  out << pad << " center = " << center << std::endl;
452  out << pad << " rmin = " << rcorner << std::endl;
453  out << pad << " rmax = " << rmax << std::endl;
454  out << pad << " volume = " << volume << std::endl;
455  }
456  else if (sampling == metropolis)
457  {
458  out << pad << " warmup = " << warmup << std::endl;
459  out << pad << " timestep = " << timestep << std::endl;
460  }
461  out << pad << " metric = " << metric << std::endl;
462  out << pad << " nparticles = " << nparticles << std::endl;
463  out << pad << " nspecies = " << nspecies << std::endl;
464  for (int s = 0; s < nspecies; ++s)
465  out << pad << " species " << s << " = " << species_size[s] << std::endl;
466  out << pad << " basis_size = " << basis_size << std::endl;
467  out << pad << " samples = " << samples << std::endl;
468  out << pad << " energy_mat = " << energy_mat << std::endl;
469  out << pad << " initialized = " << initialized << std::endl;
470  out << pad << " normalized = " << normalized << std::endl;
471  out << pad << " volume_normed = " << volume_normed << std::endl;
472  out << pad << " rsamples : " << rsamples.size() << std::endl;
473  if (evaluator == matrix)
474  {
475  out << pad << " Phi_MB : " << Phi_MB.rows() << " " << Phi_MB.cols() << std::endl;
476  for (int s = 0; s < nspecies; ++s)
477  {
478  out << pad << " matrices/vectors for species " << s << std::endl;
479  if (energy_mat)
480  if (energy_mat)
481  out << pad << " E_N : " << E_N[s]->size() << std::endl;
482  out << pad << " Phi_NB : " << Phi_NB[s]->rows() << " " << Phi_NB[s]->cols() << std::endl;
483  out << pad << " Psi_NM : " << Psi_NM[s]->rows() << " " << Psi_NM[s]->cols() << std::endl;
484  out << pad << " Phi_Psi_NB : " << Phi_Psi_NB[s]->rows() << " " << Phi_Psi_NB[s]->cols() << std::endl;
485  out << pad << " N_BB : " << N_BB[s]->rows() << " " << N_BB[s]->cols() << std::endl;
486  if (energy_mat)
487  if (energy_mat)
488  out << pad << " E_BB : " << E_BB[s]->rows() << " " << E_BB[s]->cols() << std::endl;
489  }
490  }
491  out << pad << " basis_norms" << std::endl;
492  for (int i = 0; i < basis_size; ++i)
493  out << pad << " " << i << " " << basis_norms[i] << std::endl;
494  out << pad << "end DensityMatrices1B" << std::endl;
495 }
std::vector< Matrix_t * > Phi_NB
std::ostream & app_log()
Definition: OutputManager.h:65
std::vector< Vector_t * > E_N
std::vector< Matrix_t * > Psi_NM
size_type cols() const
Definition: OhmmsMatrix.h:78
QMCTraits::PosType PosType
std::vector< Matrix_t * > E_BB
size_type rows() const
Definition: OhmmsMatrix.h:77
std::vector< Matrix_t * > Phi_Psi_NB
std::vector< Matrix_t * > N_BB
std::vector< PosType > rsamples

◆ reset()

void reset ( )

Definition at line 89 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::center, DensityMatrices1B::check_derivatives, DensityMatrices1B::check_overlap, OperatorBase::COLLECTABLE, DensityMatrices1B::eindex, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluator, DensityMatrices1B::initialized, DensityMatrices1B::integrator, DensityMatrices1B::loop, DensityMatrices1B::metric, DensityMatrices1B::nindex, DensityMatrices1B::normalized, DensityMatrices1B::points, OperatorBase::request_, TraceRequest::request_array(), TraceRequest::request_scalar(), DensityMatrices1B::samples, DensityMatrices1B::sampling, DensityMatrices1B::scale, DensityMatrices1B::T_trace, DensityMatrices1B::timestep, DensityMatrices1B::uniform_grid, DensityMatrices1B::uniform_random, OperatorBase::update_mode_, DensityMatrices1B::use_drift, DensityMatrices1B::Vc_trace, DensityMatrices1B::Vcc_trace, DensityMatrices1B::volume_based, DensityMatrices1B::volume_normed, DensityMatrices1B::Vq_trace, DensityMatrices1B::Vqc_trace, DensityMatrices1B::Vqq_trace, DensityMatrices1B::w_trace, DensityMatrices1B::warmed_up, DensityMatrices1B::warmup, DensityMatrices1B::write_acceptance_ratio, and DensityMatrices1B::write_rstats.

Referenced by DensityMatrices1B::DensityMatrices1B().

90 {
91  // uninitialized data
92  w_trace = NULL;
93  T_trace = NULL;
94  Vq_trace = NULL;
95  Vc_trace = NULL;
96  Vqq_trace = NULL;
97  Vqc_trace = NULL;
98  Vcc_trace = NULL;
99  basis_size = -1;
100  nindex = -1;
101  eindex = -1;
102  uniform_random = NULL;
103  // basic HamiltonianBase info
104  update_mode_.set(COLLECTABLE, 1);
105  // default values
106  energy_mat = false;
108  evaluator = loop;
110  scale = 1.0;
111  center = 0.0;
112  points = 10;
113  samples = 10;
114  warmup = 30;
115  timestep = .5;
116  use_drift = false;
117  warmed_up = false;
118  metric = 1.0;
119  write_acceptance_ratio = false;
120  write_rstats = false;
121  normalized = true;
122  volume_normed = true;
123  check_overlap = false;
124  check_derivatives = false;
125  // trace data is required
126  request_.request_scalar("weight");
127  request_.request_array("Kinetic_complex");
128  request_.request_array("Vq");
129  request_.request_array("Vc");
130  request_.request_array("Vqq");
131  request_.request_array("Vqc");
132  request_.request_array("Vcc");
133  // has not been initialized
134  initialized = false;
135 }
TraceSample< TraceReal > * w_trace
TraceRequest request_
whether traces are being collected
Definition: OperatorBase.h:531
TraceSample< TraceComp > * T_trace
RandomBase< FullPrecRealType > * uniform_random
void request_array(const std::string &name, bool write=false)
Definition: TraceManager.h:233
CombinedTraceSample< TraceReal > * Vc_trace
CombinedTraceSample< TraceReal > * Vq_trace
CombinedTraceSample< TraceReal > * Vcc_trace
CombinedTraceSample< TraceReal > * Vqq_trace
void request_scalar(const std::string &name, bool write=false)
Definition: TraceManager.h:224
CombinedTraceSample< TraceReal > * Vqc_trace
std::bitset< 8 > update_mode_
set the current update mode
Definition: OperatorBase.h:521

◆ resetTargetParticleSet()

void resetTargetParticleSet ( ParticleSet P)
inlineoverridevirtual

Reset the data with the target ParticleSet.

Parameters
Pnew target ParticleSet

Implements OperatorBase.

Definition at line 167 of file DensityMatrices1B.h.

167 {}

◆ same() [1/2]

bool same ( Vector_t v1,
Vector_t v2,
RealType  tol = 1e-6 
)

Definition at line 1403 of file DensityMatrices1B.cpp.

References APP_ABORT, DensityMatrices1B::match(), and Vector< T, Alloc >::size().

Referenced by DensityMatrices1B::compare().

1404 {
1405  if (v1.size() != v2.size())
1406  APP_ABORT("DensityMatrices1B::same(vector) vectors differ in size");
1407  bool sm = true;
1408  for (int i = 0; i < v1.size(); ++i)
1409  sm &= match(v1[i], v2[i]);
1410  return sm;
1411 }
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
bool match(Value_t e1, Value_t e2, RealType tol=1e-12)

◆ same() [2/2]

bool same ( Matrix_t m1,
Matrix_t m2,
RealType  tol = 1e-6 
)

Definition at line 1413 of file DensityMatrices1B.cpp.

References APP_ABORT, Matrix< T, Alloc >::cols(), DensityMatrices1B::match(), qmcplusplus::n, and Matrix< T, Alloc >::rows().

1414 {
1415  if (m1.rows() != m2.rows() || m1.cols() != m2.cols())
1416  APP_ABORT("DensityMatrices1B::same(matrix) matrices differ in size");
1417  bool sm = true;
1418  int n = m1.rows() * m1.cols();
1419  for (int i = 0; i < n; ++i)
1420  sm &= match(m1(i), m2(i));
1421  return sm;
1422 }
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
bool match(Value_t e1, Value_t e2, RealType tol=1e-12)

◆ set_state() [1/2]

void set_state ( xmlNodePtr  cur)

Definition at line 156 of file DensityMatrices1B.cpp.

References CompositeSPOSet::add(), APP_ABORT, DensityMatrices1B::basis_functions, DensityMatrices1B::basis_size, DensityMatrices1B::center, DensityMatrices1B::check_derivatives, DensityMatrices1B::check_overlap, DensityMatrices1B::density, DensityMatrices1B::DIM, qmcplusplus::Units::charge::e, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluator, qmcplusplus::exp(), TrialWaveFunction::getSPOMap(), getXMLAttributeValue(), DensityMatrices1B::ind_dims, DensityMatrices1B::integrator, DensityMatrices1B::lattice_, qmcplusplus::log(), DensityMatrices1B::loop, DensityMatrices1B::matrix, DensityMatrices1B::metric, DensityMatrices1B::metropolis, DensityMatrices1B::normalized, DensityMatrices1B::periodic, DensityMatrices1B::points, qmcplusplus::pow(), DensityMatrices1B::Psi, putContent(), DensityMatrices1B::rcorner, DensityMatrices1B::samples, DensityMatrices1B::sampling, DensityMatrices1B::scale, SPOSet::size(), qmcplusplus::spomap, qmcplusplus::SUPERCELL_OPEN, DensityMatrices1B::timestep, DensityMatrices1B::uniform, DensityMatrices1B::uniform_grid, DensityMatrices1B::use_drift, DensityMatrices1B::volume, DensityMatrices1B::volume_based, DensityMatrices1B::volume_normed, DensityMatrices1B::warmup, DensityMatrices1B::write_acceptance_ratio, and DensityMatrices1B::write_rstats.

Referenced by DensityMatrices1B::DensityMatrices1B(), and DensityMatrices1B::put().

157 {
158  bool center_defined = false;
159  std::string emstr = "no";
160  std::string igstr = "uniform_grid";
161  std::string evstr = "loop";
162  std::string costr = "no";
163  std::string cdstr = "no";
164  std::string arstr = "no";
165  std::string udstr = "no";
166  std::string wrstr = "no";
167  std::string nmstr = "yes";
168  std::string vnstr = "yes";
169  std::vector<std::string> sposets;
170 
171  xmlNodePtr element = cur->xmlChildrenNode;
172  while (element != NULL)
173  {
174  std::string ename((const char*)element->name);
175  if (ename == "parameter")
176  {
177  const std::string name(getXMLAttributeValue(element, "name"));
178  if (name == "basis")
179  putContent(sposets, element);
180  else if (name == "energy_matrix")
181  putContent(emstr, element);
182  else if (name == "integrator")
183  putContent(igstr, element);
184  else if (name == "evaluator")
185  putContent(evstr, element);
186  else if (name == "scale")
187  putContent(scale, element);
188  else if (name == "center")
189  {
190  putContent(center, element);
191  center_defined = true;
192  }
193  else if (name == "points")
194  putContent(points, element);
195  else if (name == "samples")
196  putContent(samples, element);
197  else if (name == "warmup")
198  putContent(warmup, element);
199  else if (name == "timestep")
200  putContent(timestep, element);
201  else if (name == "use_drift")
202  putContent(udstr, element);
203  else if (name == "check_overlap")
204  putContent(costr, element);
205  else if (name == "check_derivatives")
206  putContent(cdstr, element);
207  else if (name == "acceptance_ratio")
208  putContent(arstr, element);
209  else if (name == "rstats")
210  putContent(wrstr, element);
211  else if (name == "normalized")
212  putContent(nmstr, element);
213  else if (name == "volume_normed")
214  putContent(vnstr, element);
215  }
216  element = element->next;
217  }
218 
219  if (scale > 1.0 + 1e-10)
220  {
221  APP_ABORT("DensityMatrices1B::put scale must be less than one");
222  }
223  else if (scale < 0.0 - 1e-10)
224  APP_ABORT("DensityMatrices1B::put scale must be greater than zero");
225 
226  // get volume and cell information
227  if (!center_defined)
228  center = lattice_.Center;
229  volume = lattice_.Volume * std::exp(DIM * std::log(scale));
230  periodic = lattice_.SuperCellEnum != SUPERCELL_OPEN;
231  rcorner = center - scale * lattice_.Center;
232 
233  energy_mat = emstr == "yes";
234  if (igstr == "uniform_grid")
235  {
238  samples = pow(points, DIM);
239  metric = volume / samples;
240  ind_dims[0] = pow(points, DIM - 1);
241  for (int d = 1; d < DIM; ++d)
242  ind_dims[d] = ind_dims[d - 1] / points;
243  }
244  else if (igstr == "uniform")
245  {
248  metric = volume / samples;
249  }
250  else if (igstr == "density")
251  {
254  metric = 1.0 / samples;
255  }
256  else
257  throw std::runtime_error(
258  "DensityMatrices1B::set_state invalid integrator\n valid options are: uniform_grid, uniform, density");
259 
260  if (evstr == "loop")
261  evaluator = loop;
262  else if (evstr == "matrix")
263  evaluator = matrix;
264  else
265  throw std::runtime_error("DensityMatrices1B::set_state invalid evaluator\n valid options are: loop, matrix");
266 
267  normalized = nmstr == "yes";
268  volume_normed = vnstr == "yes";
269  use_drift = udstr == "yes";
270  check_overlap = costr == "yes";
271  check_derivatives = cdstr == "yes";
272  write_rstats = wrstr == "yes";
273  write_acceptance_ratio = arstr == "yes";
274 
275 
276  // get the sposets that form the basis
277  if (sposets.size() == 0)
278  throw std::runtime_error("DensityMatrices1B::put basis must have at least one sposet");
279 
280  for (int i = 0; i < sposets.size(); ++i)
281  {
282  auto& spomap = Psi.getSPOMap();
283  auto spo_it = spomap.find(sposets[i]);
284  if (spo_it == spomap.end())
285  throw std::runtime_error("DensityMatrices1B::put sposet " + sposets[i] + " does not exist.");
286  basis_functions.add(spo_it->second->makeClone());
287  }
289 
290  if (basis_size < 1)
291  throw std::runtime_error("DensityMatrices1B::put basis_size must be greater than one");
292 }
int size() const
return the size of the orbital set Ye: this needs to be replaced by getOrbitalSetSize(); ...
Definition: SPOSet.h:75
MakeReturn< BinaryNode< FnPow, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t, typename CreateLeaf< Vector< T2, C2 > >::Leaf_t > >::Expression_t pow(const Vector< T1, C1 > &l, const Vector< T2, C2 > &r)
const SPOMap & getSPOMap() const
spomap_ reference accessor
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
MakeReturn< UnaryNode< FnExp, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t exp(const Vector< T1, C1 > &l)
void add(std::unique_ptr< SPOSet > component)
add a sposet component to this composite sposet
std::string getXMLAttributeValue(const xmlNodePtr cur, const std::string_view name)
get the value string for attribute name if name is unfound in cur you get an empty string back this i...
bool putContent(T &a, xmlNodePtr cur)
replaces a&#39;s value with the first "element" in the "string" returned by XMLNodeString{cur}.
Definition: libxmldefs.h:88
MakeReturn< UnaryNode< FnLog, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t log(const Vector< T1, C1 > &l)

◆ set_state() [2/2]

void set_state ( DensityMatrices1B master)

Definition at line 295 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_size, DensityMatrices1B::DIM, DensityMatrices1B::energy_mat, DensityMatrices1B::evaluator, DensityMatrices1B::ind_dims, DensityMatrices1B::integrator, DensityMatrices1B::metric, DensityMatrices1B::normalized, DensityMatrices1B::periodic, DensityMatrices1B::points, DensityMatrices1B::rcorner, DensityMatrices1B::samples, DensityMatrices1B::sampling, DensityMatrices1B::scale, DensityMatrices1B::timestep, DensityMatrices1B::use_drift, DensityMatrices1B::volume, DensityMatrices1B::volume_normed, and DensityMatrices1B::warmup.

296 {
297  basis_size = master.basis_size;
298  energy_mat = master.energy_mat;
299  integrator = master.integrator;
300  evaluator = master.evaluator;
301  sampling = master.sampling;
302  scale = master.scale;
303  points = master.points;
304  samples = master.samples;
305  warmup = master.warmup;
306  timestep = master.timestep;
307  use_drift = master.use_drift;
308  volume = master.volume;
309  periodic = master.periodic;
310  metric = master.metric;
311  rcorner = master.rcorner;
312  normalized = master.normalized;
313  volume_normed = master.volume_normed;
314  for (int d = 0; d < DIM; ++d)
315  ind_dims[d] = master.ind_dims[d];
316 }

◆ setObservables()

void setObservables ( PropertySetType plist)
inlineoverridevirtual

Set the values evaluated by this object to plist Default implementation is to assign Value which is updated by evaluate function using my_index_.

Parameters
plistRecordNameProperty

Reimplemented from OperatorBase.

Definition at line 168 of file DensityMatrices1B.h.

168 {}

◆ setParticlePropertyList()

void setParticlePropertyList ( PropertySetType plist,
int  offset 
)
inlineoverridevirtual

Reimplemented from OperatorBase.

Definition at line 169 of file DensityMatrices1B.h.

169 {}

◆ setRandomGenerator()

void setRandomGenerator ( RandomBase< FullPrecRealType > *  rng)
overridevirtual

Set the Random Generator object TODO: add docs.

Parameters
rng

Reimplemented from OperatorBase.

Definition at line 519 of file DensityMatrices1B.cpp.

References DensityMatrices1B::uniform_random.

519 { uniform_random = rng; }
RandomBase< FullPrecRealType > * uniform_random

◆ test_derivatives()

void test_derivatives ( )

Definition at line 1350 of file DensityMatrices1B.cpp.

References APP_ABORT, qmcplusplus::app_log(), DensityMatrices1B::dens, DensityMatrices1B::density_drift(), DensityMatrices1B::DIM, DensityMatrices1B::drift, qmcplusplus::Units::charge::e, DensityMatrices1B::generate_samples(), DensityMatrices1B::rsamples, qmcplusplus::Units::time::s, DensityMatrices1B::timestep, and DensityMatrices1B::warmup_sampling().

Referenced by DensityMatrices1B::evaluate().

1351 {
1352  app_log() << "DensityMatrices1B::test_derivatives checking drift" << std::endl;
1353 
1354  PosType r, rtmp;
1355 
1356  RealType delta = 1e-5;
1357 
1358  RealType dens, densp, densm;
1359  PosType drift, driftn, drifttmp;
1360 
1361  app_log() << " warming up" << std::endl;
1362  warmup_sampling();
1363  app_log() << " generating samples" << std::endl;
1364  generate_samples(1.0);
1365 
1366  app_log() << " testing derivatives at sample points" << std::endl;
1367  for (int s = 0; s < rsamples.size(); ++s)
1368  {
1369  r = rsamples[s];
1370 
1371  density_drift(r, dens, drift);
1372 
1373  for (int d = 0; d < DIM; ++d)
1374  {
1375  rtmp = r;
1376 
1377  rtmp[d] = r[d] + delta;
1378  density_drift(rtmp, densp, drifttmp);
1379 
1380  rtmp[d] = r[d] - delta;
1381  density_drift(rtmp, densm, drifttmp);
1382 
1383  driftn[d] = (densp - densm) / (2 * delta);
1384  }
1385  driftn *= .5 * timestep / dens;
1386 
1387  app_log() << s << std::endl;
1388  app_log() << " " << driftn << std::endl;
1389  app_log() << " " << drift << std::endl;
1390  }
1391 
1392  APP_ABORT("DensityMatrices1B::test_derivatives");
1393 }
void density_drift(const PosType &r, RealType &dens, PosType &drift)
std::ostream & app_log()
Definition: OutputManager.h:65
QMCTraits::PosType PosType
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
void generate_samples(RealType weight, int steps=0)
QMCTraits::RealType RealType
std::vector< PosType > rsamples

◆ test_overlap()

void test_overlap ( )
inline

Definition at line 1289 of file DensityMatrices1B.cpp.

References qmcplusplus::abs(), APP_ABORT, qmcplusplus::app_log(), DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, qmcplusplus::conj(), DensityMatrices1B::DIM, DensityMatrices1B::lattice_, omptarget::min(), DensityMatrices1B::points, qmcplusplus::pow(), DensityMatrices1B::rcorner, Array< T, D, ALLOC >::resize(), DensityMatrices1B::scale, DensityMatrices1B::update_basis(), and DensityMatrices1B::volume.

Referenced by DensityMatrices1B::put().

1290 {
1291  //int ngrid = std::max(200,points);
1292  int ngrid = std::max(50, points);
1293  int ngtot = pow(ngrid, DIM);
1294 
1295  PosType rp;
1296  RealType du = scale / ngrid;
1297  RealType dV = volume / ngtot;
1298 
1299  PosType rmin = std::numeric_limits<RealType>::max();
1300  PosType rmax = -std::numeric_limits<RealType>::max();
1301  int gdims[DIM];
1302  gdims[0] = pow(ngrid, DIM - 1);
1303  for (int d = 1; d < DIM; ++d)
1304  gdims[d] = gdims[d - 1] / ngrid;
1305 
1306  Array<Value_t, 2> omat;
1307  omat.resize(basis_size, basis_size);
1308  for (int i = 0; i < basis_size; ++i)
1309  for (int j = 0; j < basis_size; ++j)
1310  omat(i, j) = 0.0;
1311 
1312  for (int p = 0; p < ngtot; ++p)
1313  {
1314  int nrem = p;
1315  for (int d = 0; d < DIM - 1; ++d)
1316  {
1317  int ind = nrem / gdims[d];
1318  rp[d] = ind * du + du / 2;
1319  nrem -= ind * gdims[d];
1320  }
1321  rp[DIM - 1] = nrem * du + du / 2;
1322  rp = lattice_.toCart(rp) + rcorner;
1323  update_basis(rp);
1324  for (int i = 0; i < basis_size; ++i)
1325  for (int j = 0; j < basis_size; ++j)
1326  omat(i, j) += qmcplusplus::conj(basis_values[i]) * basis_values[j] * dV;
1327  for (int d = 0; d < DIM; ++d)
1328  {
1329  rmin[d] = std::min(rmin[d], rp[d]);
1330  rmax[d] = std::max(rmax[d], rp[d]);
1331  }
1332  }
1333 
1334  app_log() << "DensityMatrices1B::test_overlap checking overlap matrix" << std::endl;
1335  app_log() << " rmin = " << rmin << std::endl;
1336  app_log() << " rmax = " << rmax << std::endl;
1337  app_log() << " overlap scale " << std::abs(omat(0, 0)) << std::endl;
1338  app_log() << " overlap matrix:" << std::endl;
1339  for (int i = 0; i < basis_size; ++i)
1340  {
1341  app_log() << std::endl;
1342  for (int j = 0; j < basis_size; ++j)
1343  app_log() << std::abs(omat(i, j)) / std::abs(omat(0, 0)) << " ";
1344  }
1345  app_log() << std::endl;
1346  APP_ABORT("DensityMatrices1B::test_overlap");
1347 }
MakeReturn< UnaryNode< FnFabs, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t abs(const Vector< T1, C1 > &l)
std::ostream & app_log()
Definition: OutputManager.h:65
void resize(const std::array< SIZET, D > &dims)
Resize the container.
Definition: OhmmsArray.h:65
T min(T a, T b)
QMCTraits::PosType PosType
MakeReturn< BinaryNode< FnPow, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t, typename CreateLeaf< Vector< T2, C2 > >::Leaf_t > >::Expression_t pow(const Vector< T1, C1 > &l, const Vector< T2, C2 > &r)
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
float conj(const float &c)
Workaround to allow conj on scalar to return real instead of complex.
QMCTraits::RealType RealType
A D-dimensional Array class based on PETE.
Definition: OhmmsArray.h:25
void update_basis(const PosType &r)

◆ update_basis()

void update_basis ( const PosType r)
inline

Definition at line 1228 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_functions, DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, CompositeSPOSet::evaluateValue(), ParticleSet::makeMove(), DensityMatrices1B::Pq, ParticleSet::R, and ParticleSet::rejectMove().

Referenced by DensityMatrices1B::density_only(), DensityMatrices1B::evaluate_check(), DensityMatrices1B::evaluate_loop(), DensityMatrices1B::generate_particle_basis(), DensityMatrices1B::generate_sample_basis(), DensityMatrices1B::integrate(), DensityMatrices1B::normalize(), and DensityMatrices1B::test_overlap().

1229 {
1230  Pq.makeMove(0, r - Pq.R[0]);
1232  Pq.rejectMove(0);
1233  for (int i = 0; i < basis_size; ++i)
1234  basis_values[i] *= basis_norms[i];
1235 }
void rejectMove(Index_t iat)
reject a proposed move in regular mode
ParticlePos R
Position.
Definition: ParticleSet.h:79
void evaluateValue(const ParticleSet &P, int iat, ValueVector &psi) override
evaluate the values of this single-particle orbital set
void makeMove(Index_t iat, const SingleParticlePos &displ, bool maybe_accept=true)
move the iat-th particle to active_pos_

◆ update_basis_d012()

void update_basis_d012 ( const PosType r)
inline

Definition at line 1238 of file DensityMatrices1B.cpp.

References DensityMatrices1B::basis_functions, DensityMatrices1B::basis_gradients, DensityMatrices1B::basis_laplacians, DensityMatrices1B::basis_norms, DensityMatrices1B::basis_size, DensityMatrices1B::basis_values, CompositeSPOSet::evaluateVGL(), ParticleSet::makeMove(), DensityMatrices1B::Pq, ParticleSet::R, and ParticleSet::rejectMove().

Referenced by DensityMatrices1B::density_drift().

1239 {
1240  Pq.makeMove(0, r - Pq.R[0]);
1242  Pq.rejectMove(0);
1243  for (int i = 0; i < basis_size; ++i)
1244  basis_values[i] *= basis_norms[i];
1245  for (int i = 0; i < basis_size; ++i)
1246  basis_gradients[i] *= basis_norms[i];
1247  for (int i = 0; i < basis_size; ++i)
1248  basis_laplacians[i] *= basis_norms[i];
1249 }
void rejectMove(Index_t iat)
reject a proposed move in regular mode
ParticlePos R
Position.
Definition: ParticleSet.h:79
void evaluateVGL(const ParticleSet &P, int iat, ValueVector &psi, GradVector &dpsi, ValueVector &d2psi) override
evaluate the values, gradients and laplacians of this single-particle orbital set ...
void makeMove(Index_t iat, const SingleParticlePos &displ, bool maybe_accept=true)
move the iat-th particle to active_pos_

◆ warmup_sampling()

void warmup_sampling ( )

Definition at line 579 of file DensityMatrices1B.cpp.

References APP_ABORT, DensityMatrices1B::center, DensityMatrices1B::density, DensityMatrices1B::density_drift(), DensityMatrices1B::diffusion(), DensityMatrices1B::dpcur, DensityMatrices1B::generate_samples(), DensityMatrices1B::integrator, DensityMatrices1B::metropolis, DensityMatrices1B::rhocur, DensityMatrices1B::rpcur, DensityMatrices1B::sampling, qmcplusplus::sqrt(), DensityMatrices1B::timestep, DensityMatrices1B::warmed_up, and DensityMatrices1B::warmup.

Referenced by DensityMatrices1B::evaluate_loop(), DensityMatrices1B::evaluate_matrix(), and DensityMatrices1B::test_derivatives().

580 {
581  if (sampling == metropolis)
582  {
583  if (!warmed_up)
584  {
586  rpcur += center;
587  if (integrator == density)
589  else
590  APP_ABORT("DensityMatrices1B::warmup_sampling invalid integrator");
591  }
592  generate_samples(1.0, warmup);
593  warmed_up = true;
594  }
595 }
void density_drift(const PosType &r, RealType &dens, PosType &drift)
void diffusion(RealType sqt, PosType &diff)
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
void generate_samples(RealType weight, int steps=0)
MakeReturn< UnaryNode< FnSqrt, typename CreateLeaf< Vector< T1, C1 > >::Leaf_t > >::Expression_t sqrt(const Vector< T1, C1 > &l)

Member Data Documentation

◆ acceptance_ratio

RealType acceptance_ratio

Definition at line 132 of file DensityMatrices1B.h.

Referenced by DensityMatrices1B::generate_density_samples().

◆ basis_functions

◆ basis_gradients

◆ basis_laplacians

ValueVector basis_laplacians

◆ basis_norms

◆ basis_size

◆ basis_values

◆ center

◆ check_derivatives

bool check_derivatives

◆ check_overlap

bool check_overlap

◆ dens

◆ dpcur

◆ drift

◆ E_BB

◆ E_N

◆ E_samp

std::vector<Value_t> E_samp

◆ E_trace

◆ eindex

◆ energy_mat

◆ evaluator

◆ ind_dims

int ind_dims[DIM]

◆ initialized

◆ integrated_values

◆ integrator

◆ lattice_

◆ metric

◆ N_BB

◆ naccepted

int naccepted

Definition at line 131 of file DensityMatrices1B.h.

Referenced by DensityMatrices1B::generate_density_samples().

◆ nindex

◆ nmoves

int nmoves

Definition at line 130 of file DensityMatrices1B.h.

Referenced by DensityMatrices1B::generate_density_samples().

◆ normalized

◆ nparticles

◆ nspecies

◆ Pc

const ParticleSet* Pc

Definition at line 85 of file DensityMatrices1B.h.

Referenced by DensityMatrices1B::getRequiredTraces().

◆ periodic

bool periodic

Definition at line 126 of file DensityMatrices1B.h.

Referenced by DensityMatrices1B::report(), and DensityMatrices1B::set_state().

◆ Phi_MB

◆ Phi_NB

◆ Phi_Psi_NB

◆ points

◆ Pq

◆ Psi

◆ Psi_NM

◆ psi_ratios

std::vector<ValueType> psi_ratios

◆ rcorner

◆ rhocur

◆ rpcur

◆ rsamples

◆ sample_weights

◆ samples

◆ sampling

◆ scale

◆ species_name

std::vector<std::string> species_name

◆ species_size

◆ T_trace

◆ timers

◆ timestep

◆ uniform_random

◆ use_drift

◆ Vc_trace

◆ Vcc_trace

◆ volume

◆ volume_normed

◆ Vq_trace

◆ Vqc_trace

◆ Vqq_trace

◆ w_trace

◆ warmed_up

◆ warmup

◆ write_acceptance_ratio

bool write_acceptance_ratio

◆ write_rstats


The documentation for this class was generated from the following files: