Path Integral Quantum Monte Carlo
Public Member Functions | Static Public Attributes
OneBodyDensityMatrixEstimator Class Reference

Compute the one body density matrix n(r) which can be used to find the momentum distribution function and structure factor. More...

#include <estimator.h>

+ Inheritance diagram for OneBodyDensityMatrixEstimator:
+ Collaboration diagram for OneBodyDensityMatrixEstimator:

Public Member Functions

 OneBodyDensityMatrixEstimator (Path &, ActionBase *, const MTRand &, double, int _frequency=20, string _label="obdm")
 Constructor. More...
 
 ~OneBodyDensityMatrixEstimator ()
 Destructor.
 
void sample ()
 Sample the OBDM. More...
 
void outputFooter ()
 For the one body density matrix estimator, we would like to output the acceptance information for the accumulate trial move.
 
string getName () const
 Get the name of the estimator.
 
- Public Member Functions inherited from EstimatorBase
 EstimatorBase (const Path &_path, ActionBase *_actionPtr, const MTRand &_random, double _maxR, int _frequency=1, string _label="")
 Constructor. More...
 
virtual ~EstimatorBase ()
 Destructor.
 
void reset ()
 Reset numAccumulated and the estimator to 0.
 
void restart (const uint32, const uint32)
 Restart the measurment process from a previous state.
 
virtual void output ()
 Output the estimator value to disk. More...
 
virtual void outputFlat ()
 Output a flat estimator value to disk. More...
 
bool baseSample ()
 Determine the basic sampling condition. More...
 
uint32 getTotNumAccumulated () const
 Get the total number of accumulated measurments.
 
uint32 getNumAccumulated () const
 Get the number of accumulated measurements since the last reset.
 
uint32 getNumSampled () const
 Get the number of samples since the last reset.
 
void prepare ()
 Prepare the estimator for i/o. More...
 
void addEndLine ()
 Add a carriage return to estimator files.
 
void appendLabel (string append)
 Append to default label. More...
 
string getLabel () const
 Get the estimator label.
 

Static Public Attributes

static const string name
 

Additional Inherited Members

- Protected Member Functions inherited from EstimatorBase
void initialize (int)
 Initialize estimator. More...
 
void initialize (vector< string >)
 Initialize estimator. More...
 
- Protected Attributes inherited from EstimatorBase
const Pathpath
 A constant reference to the paths.
 
ActionBaseactionPtr
 A pointer to the action.
 
MTRand random
 
double maxR
 
fstream * outFilePtr
 The output fie.
 
map< string, int > estIndex
 Map estimator labels to indices.
 
Array< double, 1 > estimator
 The estimator array.
 
Array< double, 1 > norm
 The normalization factor for each estimator.
 
int numEst
 The number of individual quantities measured.
 
int frequency
 The frequency at which we accumulate.
 
int startSlice
 Where imaginary time averages begin.
 
int endSlice
 Where imaginary time averages end.
 
int endDiagSlice
 Where imaginary time averages end for diagonal estimiators.
 
vector< double > sliceFactor
 Used to properly incorporate end affects.
 
string label
 The label used for the output file.
 
uint32 numSampled
 The number of times we have sampled.
 
uint32 numAccumulated
 The number of accumulated values.
 
uint32 totNumAccumulated
 The total number of accumulated values.
 
int numBeads0
 The target number of beads for the canonical ensemble.
 
bool diagonal
 Is this a diagonal estimator?
 
bool endLine
 Should we output a carriage return?
 
bool canonical
 Are we in the canonical ensemble?
 
string header
 The data file header.
 

Detailed Description

Compute the one body density matrix n(r) which can be used to find the momentum distribution function and structure factor.

We use a reference to the global random number generator, and thus this estimator can effect final results. We also need a non-constant local reference to the Path which is used for temporary updates.

Definition at line 644 of file estimator.h.

Constructor & Destructor Documentation

◆ OneBodyDensityMatrixEstimator()

OneBodyDensityMatrixEstimator::OneBodyDensityMatrixEstimator ( Path _path,
ActionBase _actionPtr,
const MTRand &  _random,
double  _maxR,
int  _frequency = 20,
string  _label = "obdm" 
)

Constructor.

The one body density matrix estimator is initialized. We measure NOBDMSEP positions, out to the maximum separation in the sample (which may depend on the type of simulation cell).

Definition at line 2256 of file estimator.cpp.

2258  :
2259  EstimatorBase(_path,_actionPtr,_random,_maxR,_frequency,_label),
2260  lpath(_path)
2261 {
2262 
2263  sqrt2LambdaTau = sqrt(2.0 * constants()->lambda() * constants()->tau());
2264 
2265  /* We chooose the maximum separation to be sqrt(NDIM)*min(L)/2 */
2266  dR = 0.5*sqrt(sum(path.boxPtr->periodic))*(blitz::min(path.boxPtr->side)) / (1.0*NOBDMSEP);
2267 
2268  /* This is an off-diagonal estimator*/
2270  diagonal = false;
2271 
2272  /* The header is the first line which contains the spatial separations */
2273  header = str(format("#%15.3E") % 0.0);
2274  for (int n = 1; n < NOBDMSEP; n++)
2275  header.append(str(format("%16.3E") % (n*dR)));
2276 
2277  numReps = 5;
2278  norm = 1.0 / (1.0*numReps);
2279 }
TinyVector< unsigned int, NDIM > periodic
Determines which dimensions have periodic bc.
Definition: container.h:29
dVec side
The linear dimensions of the box.
Definition: container.h:31
bool diagonal
Is this a diagonal estimator?
Definition: estimator.h:106
Array< double, 1 > norm
The normalization factor for each estimator.
Definition: estimator.h:91
string header
The data file header.
Definition: estimator.h:110
EstimatorBase(const Path &_path, ActionBase *_actionPtr, const MTRand &_random, double _maxR, int _frequency=1, string _label="")
Constructor.
Definition: estimator.cpp:103
const Path & path
A constant reference to the paths.
Definition: estimator.h:78
void initialize(int)
Initialize estimator.
Definition: estimator.cpp:201
const Container * boxPtr
A constant reference to the container class.
Definition: path.h:43
#define NOBDMSEP
Spatial separations to be used in the one body density matrix.
Definition: common.h:91
Ttype & min(Ttype &x, Ttype &y)
Minimum of two inputs.
Definition: common.h:142
ConstantParameters * constants()
Global public access to the constants.
Definition: constants.h:201
+ Here is the call graph for this function:

Member Function Documentation

◆ sample()

void OneBodyDensityMatrixEstimator::sample ( )
virtual

Sample the OBDM.

We overload the sample method for the one body density matrix as we only want to measure when the gap is not too large, otherwise we will be dominated by tiny close probabilities.

!!NB!! We only measure the OBDM when the tail is on an even time slice

Reimplemented from EstimatorBase.

Definition at line 2296 of file estimator.cpp.

2296  {
2297  numSampled++;
2298  if ( frequency &&
2299  ((numSampled % frequency) == 0) &&
2301  (path.worm.gap > 0) && (path.worm.gap <= constants()->Mbar()) &&
2302  ( (lpath.worm.tail[0] % 2) == 0) ) {
2303 
2304  /* If we are canonical, we want the closed configuration to be close
2305  * to our ideal one */
2306  if ( (!canonical) ||
2307  (abs(path.worm.getNumBeadsOn()+path.worm.gap-numBeads0) <= 2) ) {
2309  numAccumulated++;
2310  accumulate();
2311  }
2312  }
2313 }
bool canonical
Are we in the canonical ensemble?
Definition: estimator.h:108
int frequency
The frequency at which we accumulate.
Definition: estimator.h:94
uint32 totNumAccumulated
The total number of accumulated values.
Definition: estimator.h:103
uint32 numSampled
The number of times we have sampled.
Definition: estimator.h:101
uint32 numAccumulated
The number of accumulated values.
Definition: estimator.h:102
int numBeads0
The target number of beads for the canonical ensemble.
Definition: estimator.h:104
Worm worm
Details on the worm.
Definition: path.h:44
bool isConfigDiagonal
Stores the diagonality of the configuration.
Definition: worm.h:39
int gap
numTimeSlices - length
Definition: worm.h:38
beadLocator tail
The coordinates of the worm tail.
Definition: worm.h:32
int getNumBeadsOn() const
Return the number of active beads.
Definition: worm.h:88
+ Here is the call graph for this function:

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