QMCPACK
LatticeDeviationEstimator Class Reference

lattice deviation estimator More...

+ Inheritance diagram for LatticeDeviationEstimator:
+ Collaboration diagram for LatticeDeviationEstimator:

Public Member Functions

 LatticeDeviationEstimator (ParticleSet &P, ParticleSet &sP, const std::string &tgroup, const std::string &sgroup)
 
 ~LatticeDeviationEstimator () override
 
std::string getClassName () const override
 return class name More...
 
bool put (xmlNodePtr cur) override
 Read the input parameter. More...
 
bool get (std::ostream &os) const override
 write about the class More...
 
Return_t evaluate (ParticleSet &P) override
 Evaluate the local energy contribution of this component. More...
 
void addObservables (PropertySetType &plist, BufferType &collectables) override
 named values to the property list Default implementaton uses addValue(plist_) 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 registerCollectables (std::vector< ObservableHelper > &h5desc, hdf_archive &file) const override
 
void resetTargetParticleSet (ParticleSet &P) override
 Reset the data with the target ParticleSet. More...
 
std::unique_ptr< OperatorBasemakeClone (ParticleSet &qp, TrialWaveFunction &psi) final
 
- Public Member Functions inherited from OperatorBase
 OperatorBase ()
 Construct a new Operator Base object Default and unique empty constructor. More...
 
virtual ~OperatorBase ()=default
 
virtual bool dependsOnWaveFunction () const
 return true if this operator depends on a wavefunction More...
 
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 setParticlePropertyList (PropertySetType &plist, int offset)
 
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 setRandomGenerator (RandomBase< FullPrecRealType > *rng)
 Set the Random Generator object TODO: add docs. More...
 
virtual void add2Hamiltonian (ParticleSet &qp, TrialWaveFunction &psi, QMCHamiltonian &targetH)
 TODO: add docs. More...
 
virtual void getRequiredTraces (TraceManager &tm)
 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...
 

Private Attributes

SpeciesSettspecies
 
SpeciesSetsspecies
 
ParticleSettpset
 
ParticleSet spset
 
std::string tgroup
 
std::string sgroup
 
int num_sites
 
bool hdf5_out
 
int h5_index
 
bool per_xyz
 
std::vector< RealTypexyz2
 
xmlNodePtr input_xml
 
const int myTableID_
 

Additional Inherited Members

- 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 > >
 
- Protected Member Functions inherited from OperatorBase
virtual void contributeScalarQuantities ()
 
virtual void checkoutScalarQuantities (TraceManager &tm)
 
virtual void collectScalarQuantities ()
 
virtual void deleteScalarQuantities ()
 
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...
 
- 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_
 

Detailed Description

lattice deviation estimator

Compute deviation of species="tgroup" in target particle set from species="sgroup" in source particle set. The motivation is to observe the deviation of protons from their crystal sites in an all electron-proton simulation of hydrogen i.e. two-component system One can also create any reference point to measure the deviation of an up electron, for example <particleset name="wf_center"> <group name="origin" size="1"> <attrib name="position" datatype="posArray" condition="0"> 0.00000000 0.00000000 0.00000000 </attrib> </group> </particleset> <estimator type="latticedeviation" name="latdev" hdf5="yes" per_xyz="yes" source="wf_center" sgroup="origin" target="e" tgroup="u"> This estimator outputs to both scalar.dat and stat.h5. The scalar.dat entries are averaged over all particles, whereas the stat.h5 entries are particle-resolved. The two sets of outputs can be compared as a consistency check for the estimator.

Definition at line 39 of file LatticeDeviationEstimator.h.

Constructor & Destructor Documentation

◆ LatticeDeviationEstimator()

LatticeDeviationEstimator ( ParticleSet P,
ParticleSet sP,
const std::string &  tgroup,
const std::string &  sgroup 
)

Definition at line 17 of file LatticeDeviationEstimator.cpp.

References APP_ABORT, qmcplusplus::app_log(), SpeciesSet::findSpecies(), ParticleSet::first(), ParticleSet::last(), LatticeDeviationEstimator::num_sites, LatticeDeviationEstimator::sgroup, LatticeDeviationEstimator::spset, LatticeDeviationEstimator::sspecies, LatticeDeviationEstimator::tgroup, LatticeDeviationEstimator::tpset, and LatticeDeviationEstimator::tspecies.

21  : tspecies(P.getSpeciesSet()),
22  sspecies(sP.getSpeciesSet()),
23  tpset(P),
24  spset(sP),
25  tgroup(tgroup_in),
26  sgroup(sgroup_in),
27  hdf5_out(false),
28  per_xyz(false),
29  myTableID_(P.addTable(sP))
30 {
31  // calculate number of source particles to use as lattice sites
32  int src_species_id = sspecies.findSpecies(sgroup);
33  num_sites = spset.last(src_species_id) - spset.first(src_species_id);
34  int tar_species_id = tspecies.findSpecies(tgroup);
35  int num_tars = tpset.last(tar_species_id) - tpset.first(tar_species_id);
36  if (num_tars != num_sites)
37  {
38  app_log() << "number of target particles = " << num_tars << std::endl;
39  app_log() << "number of source particles = " << num_sites << std::endl;
40  APP_ABORT("nsource != ntargets in LatticeDeviationEstimator");
41  }
42 }
std::ostream & app_log()
Definition: OutputManager.h:65
int first(int igroup) const
return the first index of a group i
Definition: ParticleSet.h:514
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
int last(int igroup) const
return the last index of a group i
Definition: ParticleSet.h:517
int findSpecies(const std::string &name) const
if the input species is not found, add a new species
Definition: SpeciesSet.h:114

◆ ~LatticeDeviationEstimator()

~LatticeDeviationEstimator ( )
inlineoverride

Definition at line 43 of file LatticeDeviationEstimator.h.

43 {}

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 169 of file LatticeDeviationEstimator.cpp.

References PooledData< T >::add(), RecordNamedProperty< T >::add(), LatticeDeviationEstimator::h5_index, LatticeDeviationEstimator::hdf5_out, OperatorBase::my_index_, OperatorBase::name_, LatticeDeviationEstimator::num_sites, OHMMS_DIM, LatticeDeviationEstimator::per_xyz, PooledData< T >::size(), and RecordNamedProperty< T >::size().

170 {
171  // get myIndex for scalar.dat
172  if (per_xyz)
173  {
174  my_index_ = plist.size();
175  for (int idir = 0; idir < OHMMS_DIM; idir++)
176  {
177  std::stringstream ss;
178  ss << idir;
179  plist.add(name_ + "_dir" + ss.str());
180  }
181  }
182  else
183  {
184  my_index_ = plist.add(name_); // same as OperatorBase::addObservables
185  }
186 
187  // get h5_index for stat.h5
188  if (hdf5_out)
189  {
190  h5_index = collectables.size();
191  std::vector<RealType> tmp;
192  if (per_xyz)
193  {
194  tmp.resize(num_sites * OHMMS_DIM);
195  }
196  else
197  {
198  tmp.resize(num_sites);
199  }
200  collectables.add(tmp.begin(), tmp.end());
201  }
202 }
int my_index_
starting index of this object
Definition: OperatorBase.h:535
#define OHMMS_DIM
Definition: config.h:64
std::string name_
name of this object
Definition: OperatorBase.h:527

◆ evaluate()

LatticeDeviationEstimator::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 96 of file LatticeDeviationEstimator.cpp.

References APP_ABORT, qmcplusplus::app_log(), ParticleSet::Collectables, ParticleSet::getDistTableAB(), ParticleSet::getTotalNum(), ParticleSet::GroupID, LatticeDeviationEstimator::h5_index, LatticeDeviationEstimator::hdf5_out, LatticeDeviationEstimator::myTableID_, LatticeDeviationEstimator::num_sites, OHMMS_DIM, LatticeDeviationEstimator::per_xyz, LatticeDeviationEstimator::sgroup, SpeciesSet::speciesName, LatticeDeviationEstimator::spset, LatticeDeviationEstimator::sspecies, OperatorBase::t_walker_, LatticeDeviationEstimator::tgroup, LatticeDeviationEstimator::tpset, LatticeDeviationEstimator::tspecies, OperatorBase::value_, Walker< t_traits, p_traits >::Weight, and LatticeDeviationEstimator::xyz2.

97 { // calculate <r^2> averaged over lattice sites
98  value_ = 0.0;
99  std::fill(xyz2.begin(), xyz2.end(), 0.0);
100 
101  RealType wgt = t_walker_->Weight;
102  const auto& d_table = P.getDistTableAB(myTableID_);
103 
104  // temp variables
105  RealType r, r2;
106  PosType dr;
107 
108  int nsite(0); // site index
109  int cur_jat(-1); // target particle index
110  for (int iat = 0; iat < spset.getTotalNum(); iat++)
111  { // for each desired source particle
112  if (sspecies.speciesName[spset.GroupID[iat]] == sgroup)
113  { // find desired species
114  for (int jat = cur_jat + 1; jat < tpset.getTotalNum(); jat++)
115  { // find corresponding (!!!! assume next) target particle
116  if (tspecies.speciesName[tpset.GroupID[jat]] == tgroup)
117  {
118  // distance between particle iat in source pset, and jat in target pset
119  r = d_table.getDistRow(jat)[iat];
120  r2 = r * r;
121  value_ += r2;
122 
123  if (hdf5_out & !per_xyz)
124  { // store deviration for each lattice site if h5 file is available
125  P.Collectables[h5_index + nsite] = wgt * r2;
126  }
127 
128  if (per_xyz)
129  {
130  dr = d_table.getDisplRow(jat)[iat];
131  for (int idir = 0; idir < OHMMS_DIM; idir++)
132  {
133  RealType dir2 = dr[idir] * dr[idir];
134  xyz2[idir] += dir2;
135  if (hdf5_out)
136  {
137  P.Collectables[h5_index + nsite * OHMMS_DIM + idir] = wgt * dir2;
138  }
139  }
140  }
141 
142  cur_jat = jat;
143  break;
144  }
145  }
146  nsite += 1; // count the number of sites, for checking only
147  } // found desired species (source particle)
148  }
149 
150  if (nsite != num_sites)
151  {
152  app_log() << "num_sites = " << num_sites << " nsites = " << nsite << std::endl;
153  APP_ABORT("Mismatch in LatticeDeivationEstimator.");
154  }
155 
156  // average per site
157  value_ /= num_sites;
158  if (per_xyz)
159  {
160  for (int idir = 0; idir < OHMMS_DIM; idir++)
161  {
162  xyz2[idir] /= num_sites;
163  }
164  }
165 
166  return value_;
167 }
size_t getTotalNum() const
Definition: ParticleSet.h:493
std::ostream & app_log()
Definition: OutputManager.h:65
#define OHMMS_DIM
Definition: config.h:64
ParticleIndex GroupID
Species ID.
Definition: ParticleSet.h:77
QMCTraits::PosType PosType
#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
std::vector< std::string > speciesName
Species name list.
Definition: SpeciesSet.h:44
QMCTraits::RealType RealType
Return_t value_
current value
Definition: OperatorBase.h:524
FullPrecRealType Weight
Weight of the walker.
Definition: Walker.h:102

◆ get()

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

write about the class

Implements OperatorBase.

Definition at line 90 of file LatticeDeviationEstimator.cpp.

References OhmmsElementBase::getName(), OperatorBase::name_, and LatticeDeviationEstimator::spset.

91 { // class description
92  os << "LatticeDeviationEstimator: " << name_ << "lattice = " << spset.getName();
93  return true;
94 }
const std::string & getName() const
return the name
std::string name_
name of this object
Definition: OperatorBase.h:527

◆ getClassName()

std::string getClassName ( ) const
inlineoverridevirtual

return class name

Implements OperatorBase.

Definition at line 45 of file LatticeDeviationEstimator.h.

45 { return "LatticeDeviationEstimator"; }

◆ makeClone()

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

Implements OperatorBase.

Definition at line 218 of file LatticeDeviationEstimator.cpp.

References LatticeDeviationEstimator::input_xml, LatticeDeviationEstimator::sgroup, LatticeDeviationEstimator::spset, and LatticeDeviationEstimator::tgroup.

219 {
220  // default constructor does not work with threads
221  //LatticeDeviationEstimator* myclone = new LatticeDeviationEstimator(*this);
222  std::unique_ptr<LatticeDeviationEstimator> myclone =
223  std::make_unique<LatticeDeviationEstimator>(qp, spset, tgroup, sgroup);
224  myclone->put(input_xml);
225  return myclone;
226 }

◆ put()

bool put ( xmlNodePtr  cur)
overridevirtual

Read the input parameter.

Parameters
curxml node for a OperatorBase object

Implements OperatorBase.

Definition at line 44 of file LatticeDeviationEstimator.cpp.

References OhmmsAttributeSet::add(), APP_ABORT, OperatorBase::COLLECTABLE, LatticeDeviationEstimator::hdf5_out, LatticeDeviationEstimator::input_xml, OHMMS_DIM, LatticeDeviationEstimator::per_xyz, OhmmsAttributeSet::put(), OperatorBase::update_mode_, and LatticeDeviationEstimator::xyz2.

45 {
46  input_xml = cur;
47  std::string hdf5_flag = "no";
48  std::string xyz_flag = "no";
49  OhmmsAttributeSet attrib;
50  attrib.add(hdf5_flag, "hdf5");
51  attrib.add(xyz_flag, "per_xyz");
52  attrib.put(cur);
53 
54  if (hdf5_flag == "yes")
55  {
56  hdf5_out = true;
57  }
58  else if (hdf5_flag == "no")
59  {
60  hdf5_out = false;
61  }
62  else
63  {
64  APP_ABORT("LatticeDeviationEstimator::put() - Please choose \"yes\" or \"no\" for hdf5 flag");
65  } // end if hdf5_flag
66  if (hdf5_out)
67  {
68  // change the default behavior of addValue() called by addObservables()
69  // YY: does this still matter if I have overwritten addObservables()?
70  update_mode_.set(COLLECTABLE, 1);
71  }
72 
73  if (xyz_flag == "yes")
74  {
75  per_xyz = true;
76  xyz2.resize(OHMMS_DIM);
77  }
78  else if (xyz_flag == "no")
79  {
80  per_xyz = false;
81  }
82  else
83  {
84  APP_ABORT("LatticeDeviationEstimator::put() - Please choose \"yes\" or \"no\" for per_xyz flag");
85  } // end if xyz_flag
86 
87  return true;
88 }
bool put(xmlNodePtr cur)
assign attributes to the set
Definition: AttributeSet.h:55
#define OHMMS_DIM
Definition: config.h:64
class to handle a set of attributes of an xmlNode
Definition: AttributeSet.h:24
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
void add(PDT &aparam, const std::string &aname, std::vector< PDT > candidate_values={}, TagStatus status=TagStatus::OPTIONAL)
add a new attribute
Definition: AttributeSet.h:42
std::bitset< 8 > update_mode_
set the current update mode
Definition: OperatorBase.h:521

◆ registerCollectables()

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

Reimplemented from OperatorBase.

Definition at line 228 of file LatticeDeviationEstimator.cpp.

References LatticeDeviationEstimator::h5_index, LatticeDeviationEstimator::hdf5_out, OperatorBase::name_, LatticeDeviationEstimator::num_sites, OHMMS_DIM, and LatticeDeviationEstimator::per_xyz.

229 {
230  if (!hdf5_out)
231  return;
232 
233  // one scalar per lattice site (i.e. source particle)
234  std::vector<int> ndim(1, num_sites);
235  if (per_xyz)
236  {
237  // one scalar per lattice site per dimension
238  ndim[0] = num_sites * OHMMS_DIM;
239  }
240 
241  // open hdf5 entry and resize
242  h5desc.emplace_back(hdf_path{name_});
243  auto& h5o = h5desc.back();
244  h5o.set_dimensions(ndim, h5_index);
245 }
#define OHMMS_DIM
Definition: config.h:64
std::string name_
name of this object
Definition: OperatorBase.h:527

◆ resetTargetParticleSet()

void resetTargetParticleSet ( ParticleSet P)
overridevirtual

Reset the data with the target ParticleSet.

Parameters
Pnew target ParticleSet

Implements OperatorBase.

Definition at line 216 of file LatticeDeviationEstimator.cpp.

216 {}

◆ setObservables()

void setObservables ( PropertySetType plist)
overridevirtual

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 204 of file LatticeDeviationEstimator.cpp.

References RecordNamedProperty< T >::begin(), copy(), OperatorBase::my_index_, LatticeDeviationEstimator::per_xyz, OperatorBase::value_, and LatticeDeviationEstimator::xyz2.

205 {
206  if (per_xyz)
207  {
208  copy(xyz2.begin(), xyz2.end(), plist.begin() + my_index_);
209  }
210  else
211  {
212  plist[my_index_] = value_; // default behavior
213  }
214 }
int my_index_
starting index of this object
Definition: OperatorBase.h:535
void copy(const Array< T1, 3 > &src, Array< T2, 3 > &dest)
Definition: Blitz.h:639
Return_t value_
current value
Definition: OperatorBase.h:524

Member Data Documentation

◆ h5_index

◆ hdf5_out

◆ input_xml

xmlNodePtr input_xml
private

◆ myTableID_

const int myTableID_
private

Definition at line 76 of file LatticeDeviationEstimator.h.

Referenced by LatticeDeviationEstimator::evaluate().

◆ num_sites

◆ per_xyz

◆ sgroup

◆ spset

◆ sspecies

◆ tgroup

◆ tpset

◆ tspecies

◆ xyz2


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