Path Integral Quantum Monte Carlo
estimator.h
Go to the documentation of this file.
1 
9 #include "common.h"
10 
11 #ifndef ESTIMATOR_H
12 #define ESTIMATOR_H
13 
14 class Path;
15 class ActionBase;
16 class Potential;
17 
18 // ========================================================================
19 // EstimatorBase Class
20 // ========================================================================
29 
30  public:
31  EstimatorBase(const Path & _path, ActionBase *_actionPtr, const MTRand &_random,
32  double _maxR, int _frequency=1, string _label="");
33  virtual ~EstimatorBase();
34 
35  /* Sample the estimator */
36  virtual void sample();
37 
38  /* Reset the estimator */
39  void reset();
40 
41  /* Restart the estimator */
42  void restart(const uint32, const uint32);
43 
44  /* Output the estimator */
45  virtual void output();
46 
47  /* Output a flat average estimator */
48  virtual void outputFlat();
49 
51  virtual void outputFooter() {};
52 
53  /* Evaluate the basic conditions for sampling */
54  bool baseSample();
55 
58 
61 
63  uint32 getNumSampled() const { return numSampled; }
64 
66  virtual string getName() const { return "base"; }
67 
69  void prepare();
70 
72  void addEndLine() {endLine = true;};
73 
75  void appendLabel(string append);
76 
78  string getLabel() const {return label;};
79 
80  protected:
81  const Path &path;
83  MTRand random; // A local copy of the random number generator
84  double maxR; // An estimator cutoff radius
85 
86  fstream *outFilePtr;
87 
88  map<string,int> estIndex;
89 
90  Array<double,1> estimator;
91  Array<double,1> norm;
92 
93  int numEst;
94  int frequency;
95  int startSlice;
96  int endSlice;
98  vector<double> sliceFactor;
99  string label;
100 
104  int numBeads0;
105 
106  bool diagonal;
107  bool endLine;
108  bool canonical;
109 
110  string header;
111 
113  virtual void accumulate() {}
114 
115  /* Initialize the estimator */
116  void initialize(int);
117  void initialize(vector<string>);
118 };
119 
120 // ========================================================================
121 // Time Estimator Class
122 // ========================================================================
128 class TimeEstimator : public EstimatorBase {
129 
130  public:
131  TimeEstimator(const Path &, ActionBase *, const MTRand &, double,
132  int _frequency=1, string _label="estimator");
133  ~TimeEstimator() {};
134 
135  void sample(); // Overload to always-on sampling
136 
137  static const string name;
138  string getName() const {return name;}
139 
140  void output(); // overload the output
141 
142  private:
143  void accumulate(); // Accumulate values
144  std::chrono::high_resolution_clock::time_point time_begin;
145  std::chrono::high_resolution_clock::time_point time_end;
146 
147 };
148 
149 
150 // ========================================================================
151 // Energy Estimator Class
152 // ========================================================================
163 
164  public:
165  EnergyEstimator (const Path &, ActionBase *, const MTRand &, double,
166  int _frequency=1, string _label="estimator");
168 
169  static const string name;
170  string getName() const {return name;}
171 
172  private:
173  void accumulate(); // Accumulate values
174 
175 };
176 
177 
178 // ========================================================================
179 // Thermodynamic and Virial Energy Estimator Class
180 // ========================================================================
200 
201  public:
202  VirialEnergyEstimator (const Path &, ActionBase *, const MTRand &, double,
203  int _frequency=1, string _label="virial");
205 
206  static const string name;
207  string getName() const {return name;}
208 
209  private:
210  void accumulate(); // Accumulate values
211 
212 };
213 
214 // ========================================================================
215 // Number Particles Estimator Class
216 // ========================================================================
221 
222  public:
223  NumberParticlesEstimator(const Path &, ActionBase *, const MTRand &, double,
224  int _frequency=1, string _label="estimator");
226 
227  static const string name;
228  string getName() const {return name;}
229 
230  private:
231  void accumulate(); // Accumulate values
232 };
233 
234 // ========================================================================
235 // Particle Density Estimator Class
236 // ========================================================================
242 
243  public:
244  ParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double,
245  int _frequency=1, string _label="position");
247 
248  static const string name;
249  string getName() const {return name;}
250 
251  void output() {outputFlat();} // overload the output
252 
253  private:
254  void accumulate(); // Accumulate values
255  vector<string> diffLabels;
256 };
257 
258 // ========================================================================
259 // BIPARTITION DENSITY ESTIMATOR CLASS
260 // ========================================================================
266 
267  public:
268  BipartitionDensityEstimator(const Path &, ActionBase *, const MTRand &, double,
269  int _frequency=1, string _label="bipart_dens");
271 
272  static const string name;
273  string getName() const {return name;}
274 
275  private:
276  void accumulate(); // Accumulate values
277 };
278 
279 // ========================================================================
280 // Linear Particle Density Estimator Class
281 // ========================================================================
287 
288  public:
289  LinearParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double,
290  int _frequency=1, string _label="lineardensity");
292 
293  static const string name;
294  string getName() const {return name;}
295 
296  private:
297  int numGrid; // The number of grid points
298  double dz; // The size of each spatial bin in the z-direction
299  void accumulate(); // Accumulate values
300  dVec side; // Local copy of container geometry
301 };
302 
303 // ========================================================================
304 // Plane Particle Density Estimator Class
305 // ========================================================================
311 
312  public:
313  PlaneParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double,
314  int _frequency=1, string _label="planedensity");
316 
317  static const string name;
318  string getName() const {return name;}
319 
320  private:
321  int numGrid; // The total number of grid boxes
322  int numLinearGrid; // The linear number of grid boxes
323  dVec dl; // The linear size of each spatial bin
324  void accumulate(); // Accumulate values
325  dVec side; // Local copy of container geometry
326 };
327 
328 // ========================================================================
329 // Averaged Plane Particle Density Estimator Class
330 // ========================================================================
336 
337  public:
338  PlaneParticleAveragePositionEstimator(const Path &, ActionBase *, const MTRand &, double,
339  int _frequency=1, string _label="planeavedensity");
341 
342  static const string name;
343  string getName() const {return name;}
344 
345  void output() {outputFlat();} // overload the output
346 
347  private:
348  int numGrid; // The total number of grid boxes
349  int numLinearGrid; // The linear number of grid boxes
350  dVec dl; // The linear size of each spatial bin
351  void accumulate(); // Accumulate values
352  dVec side; // Local copy of container geometry
353 };
354 
355 // ========================================================================
356 // Averaged Plane External Potential Estimator Class
357 // ========================================================================
363 
364  public:
366  const MTRand &, double, int _frequency=1, string _label="planeaveVext");
368 
369  static const string name;
370  string getName() const {return name;}
371 
372  void output(); // overload the output
373 
374 
375  private:
376  int numGrid; // The total number of grid boxes
377  int numLinearGrid; // The linear number of grid boxes
378  dVec dl; // The linear size of each spatial bin
379  void accumulate(); // Accumulate values
380  dVec side; // Local copy of container geometry
381 };
382 
383 
384 // ========================================================================
385 // Number Distribution Estimator Class
386 // ========================================================================
391 
392  public:
393  NumberDistributionEstimator(const Path &, ActionBase *, const MTRand &, double,
394  int _frequency=1, string _label="number");
396 
397  static const string name;
398  string getName() const {return name;}
399 
400  private:
401  int startParticleNumber; // The smallest number of particles
402  int endParticleNumber; // The largest number of particles
403  int maxNumParticles; // The particle Range
404  int particleShift; // The shift in the particle array
405 
406  void accumulate(); // Accumulate values
407 };
408 
409 // ========================================================================
410 // Superfluid Fraction Estimator Class
411 // ========================================================================
417 
418  public:
419  SuperfluidFractionEstimator(const Path &, ActionBase *, const MTRand &, double,
420  int _frequency=1, string _label="super");
422 
423  static const string name;
424  string getName() const {return name;}
425 
426  private:
427  int windMax; // The maximum winding number considered
428 
429  void accumulate(); // Accumulate values
430 
431 };
432 
433 // ========================================================================
434 // Local Superfluid Density Estimator Class
435 // ========================================================================
440 
441  public:
442  LocalSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double,
443  int _frequency=1, string _label="locsuper");
445 
446  void output();
447 
448  static const string name;
449  string getName() const {return name;}
450 
451  private:
452  int numGrid; // The number of grid points
453  double dR; // The size of the radial bin
454  Array <double,1> locAz; // The local path area estimator
455  Array <double,1> locA2; // The local area squared
456  Array <double,1> locWz; // The local winding number estimator
457 
458  void accumulate(); // Accumulate values
459 
460 };
461 
462 // ========================================================================
463 // Plane Winding Superfluid Density Estimator Class
464 // ========================================================================
469 
470  public:
471  PlaneWindingSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &,
472  double, int _frequency=1, string _label="planewind");
474 
475  static const string name;
476  string getName() const {return name;}
477 
478  private:
479  dVec side; // A local copy of the dimensions
480  double dx; // The linear x-size of the spatial bin
481  double dy; // The linear y-size of the spatial bin
482  int numGrid; // The number of grid points
483  Array <double,1> locWz; // The local winding number estimator
484  void accumulate(); // Accumulate values
485 
486 };
487 
488 // ========================================================================
489 // Plane Area Superfluid Density Estimator Class
490 // ========================================================================
495 
496  public:
497  PlaneAreaSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &,
498  double, int _frequency=1, string _label="planearea");
500 
501  static const string name;
502  string getName() const {return name;}
503 
504  private:
505  dVec side; // A local copy of the dimensions
506  double dx; // The linear x-size of the spatial bin
507  double dy; // The linear y-size of the spatial bin
508  int numGrid; // The number of grid points
509  Array <double,1> locAz; // The local area estimator
510  void accumulate(); // Accumulate values
511 
512 };
513 
514 // ========================================================================
515 // Radial Winding Superfluid Density Estimator Class
516 // ========================================================================
521 
522  public:
523  RadialWindingSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &,
524  double, int _frequency=1, string _label="radwind");
526 
527  static const string name;
528  string getName() const {return name;}
529 
530  private:
531  double dR; // The size of the radial bin
532  int numGrid; // The number of grid points
533  Array <double,1> locWz; // The local winding number estimator
534  void accumulate(); // Accumulate values
535 
536 };
537 
538 // ========================================================================
539 // Radial Area Superfluid Density Estimator Class
540 // ========================================================================
545 
546  public:
547  RadialAreaSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &,
548  double, int _frequency=1, string _label="radarea");
550 
551  static const string name;
552  string getName() const {return name;}
553 
554  private:
555  double dR; // The size of the radial bin
556  int numGrid; // The number of grid points
557  Array <double,1> locAz; // The local winding number estimator
558  void accumulate(); // Accumulate values
559 
560 };
561 
562 // ========================================================================
563 // Diagonal Fraction Estimator Class
564 // ========================================================================
569 
570  public:
571  DiagonalFractionEstimator(const Path &, ActionBase *, const MTRand &,
572  double, int _frequency=1, string _label="estimator");
574 
575  void sample(); // Overload to always-on sampling
576 
577  static const string name;
578  string getName() const {return name;}
579 
580  private:
581  void accumulate(); // Accumulate values
582 };
583 
584 // ========================================================================
585 // Permutation Cycle Estimator Class
586 // ========================================================================
591 
592  public:
593  PermutationCycleEstimator(const Path &, ActionBase *, const MTRand &,
594  double, int _frequency=1, string _label="pcycle");
596 
597  static const string name;
598  string getName() const {return name;}
599 
600  private:
601  Array <bool,1> doBead; // Used for ensuring we don't double count beads
602  int maxNumCycles; // The maximum number of cycles to consider
603  void accumulate(); // Accumulate values
604 };
605 
606 // ========================================================================
607 // Local Permutation Cycle Estimator Class
608 // ========================================================================
614 
615  public:
616  LocalPermutationEstimator(const Path &, ActionBase *, const MTRand &,
617  double, int _frequency=1, string _label="locperm");
619 
620  void output() {outputFlat();}
621 
622  static const string name;
623  string getName() const {return name;}
624 
625  private:
626  Array <int, 1> numBeadInGrid;
627  Array <bool,1> doBead; // Used for ensuring we don't double count beads
628  int maxNumCycles; // The maximum number of cycles to consider
629  void accumulate(); // Accumulate values
630 };
631 
632 
633 // ========================================================================
634 // One Body Density Matrix Estimator Class
635 // ========================================================================
645 
646  public:
647  OneBodyDensityMatrixEstimator(Path &, ActionBase *, const MTRand &,
648  double, int _frequency=20, string _label="obdm");
650 
651  void sample(); // Sample the estimator
652  void outputFooter(); // Output the acceptance footer to disk
653 
654  static const string name;
655  string getName() const {return name;}
656 
657  private:
658  Path &lpath; // A non-constant local reference to the path
659 
660  double dR; // The discretization
661  int numReps; // The number of measurments reps
662  uint32 numAccepted; // The number of moves accepted
663  uint32 numAttempted; // The number of moves attempted
664 
665 
666  dVec newTailPos,oldTailPos; // The new and old tail position
667  dVec newHeadPos; // The new head position
668  dVec newRanPos,neighborPos; // The random shift
669 
670  double sqrt2LambdaTau; // sqrt(2 * lambda * tau)
671  double rho0Norm; // Free density matrix
672  double oldAction,newAction; // The old and new action
673 
674  /* Get a random vector */
675  dVec getRandomVector(const double);
676 
677  /* Returns a position used by the staging algorithm */
678  dVec newStagingPosition(const beadLocator &, const beadLocator &, const int, const int);
679 
680  /* Accumulate values */
681  void accumulate();
682 };
683 
684 // ========================================================================
685 // Pair Correlation Estimator Class
686 // ========================================================================
691 
692  public:
693  PairCorrelationEstimator(const Path &, ActionBase *, const MTRand &,
694  double, int _frequency=1, string _label="pair");
696 
697  static const string name;
698  string getName() const {return name;}
699 
700  private:
701  void accumulate(); // Accumulate values
702  double dR; // The discretization
703 };
704 
705 // ========================================================================
706 // Static Structure Factor Estimator Class
707 // ========================================================================
712 
713  public:
715  const MTRand &, double, int _frequency=1, string _label="ssf");
717 
718  static const string name;
719  string getName() const {return name;}
720 
721  private:
722  void accumulate(); // Accumulate values
723  Array <double,1> sf; // structure factor
724  Array <double,1> q; // the q values
725  int numq; // number of q values
726 };
727 
728 // ========================================================================
729 // Intermediate Scattering Function Estimator Class
730 // ========================================================================
735 
736  public:
738  const MTRand &, double, int _frequency=1, string _label="isf");
740 
741  static const string name;
742  string getName() const {return name;}
743 
744  private:
745  void accumulate(); // Accumulate values
746  Array <double,1> isf; // local intermediate scattering function
747 
748  int numq; // the number of q-magnitudes
749  Array <int,1> numqVecs; // the number of q-vectors with a given magnitude
750  vector <vector<dVec> > q; // the q-vectors
751 };
752 
753 // ========================================================================
754 // Radial Density Estimator Class
755 // ========================================================================
760 
761  public:
762  RadialDensityEstimator(const Path &, ActionBase *, const MTRand &,
763  double, int _frequency=1, string _label="radial");
765 
766  static const string name;
767  string getName() const {return name;}
768 
769  private:
770  void accumulate(); // Accumulate values
771  double dR; // The discretization
772 };
773 
774 // ========================================================================
775 // Worm Properties Estimator Class
776 // ========================================================================
781 
782  public:
783  WormPropertiesEstimator(const Path &, ActionBase *, const MTRand &,
784  double, int _frequency=1, string _label="worm");
786 
787  static const string name;
788  string getName() const {return name;}
789 
790  private:
791  dVec sep; // head-tail separation
792  void accumulate(); // Accumulate values
793 };
794 
797 // CYLINDER ESTIMATORS ///////////////////////////////////////////////////
800 
801 // ========================================================================
802 // Cylinder Energy Estimator Class
803 // ========================================================================
814 
815  public:
816  CylinderEnergyEstimator(const Path &, ActionBase *, const MTRand &,
817  double, int _frequency=1, string _label="cyl_estimator");
819 
820  static const string name;
821  string getName() const {return name;}
822 
823  private:
824 
825  void accumulate(); // Accumulate values
826 
827 };
828 
829 // ========================================================================
830 // Cylinder Number Particles Estimator Class
831 // ========================================================================
836 
837  public:
838  CylinderNumberParticlesEstimator(const Path &, ActionBase *, const MTRand &,
839  double, int _frequency=1, string _label="cyl_estimator");
841 
842  static const string name;
843  string getName() const {return name;}
844 
845  private:
846  void accumulate(); // Accumulate values
847 };
848 
849 // ========================================================================
850 // Cylinder Number Distribution Estimator Class
851 // ========================================================================
856 
857  public:
858  CylinderNumberDistributionEstimator(const Path &, ActionBase *, const MTRand &,
859  double, int _frequency=1, string _label="cyl_number");
861 
862  static const string name;
863  string getName() const {return name;}
864 
865  private:
866  int maxNumParticles; // The maximum number considered
867 
868  void accumulate(); // Accumulate values
869 };
870 
871 // ========================================================================
872 // Cylinder Linear Density Estimator Class
873 // ========================================================================
878 
879  public:
880  CylinderLinearDensityEstimator(const Path &, ActionBase *, const MTRand &,
881  double, int _frequency=1, string _label="cyl_linedensity");
883 
884  static const string name;
885  string getName() const {return name;}
886 
887  private:
888  double dz; // The bin-size in the z-direction
889  double Lz; // The length of the cylinder
890 
891  void accumulate(); // Accumulate values
892 };
893 
894 // ========================================================================
895 // Cylinder Superfluid Fraction Estimator Class
896 // ========================================================================
902 
903  public:
904  CylinderSuperfluidFractionEstimator(const Path &, ActionBase *, const MTRand &,
905  double, int _frequency=1, string _label="cyl_super");
907 
908  static const string name;
909  string getName() const {return name;}
910 
911  private:
912  Array <bool,1> doBead; // Used for ensuring we don't double count beads
913  int windMax; // The maximum winding number considered
914 
915  void accumulate(); // Accumulate values
916 
917 };
918 
919 // ========================================================================
920 // Cylinder One Body Density Matrix Estimator Class
921 // ========================================================================
931 
932  public:
934  double, int _frequency=20, string _label="cyl_obdm");
936 
937  void sample(); // Sample the estimator
938 
939  static const string name;
940  string getName() const {return name;}
941 
942  private:
943  Path &lpath; // A non-constant local reference to the path
944 
945  double dR; // The discretization
946  int numReps; // The number of measurments reps
947  uint32 numAccepted; // The number of moves accepted
948  uint32 numAttempted; // The number of moves attempted
949 
950 
951  dVec newTailPos,oldTailPos; // The new and old tail position
952  dVec newHeadPos; // The new head position
953  dVec newRanPos,neighborPos; // The random shift
954 
955  double sqrt2LambdaTau; // sqrt(2 * lambda * tau)
956  double rho0Norm; // Free density matrix
957  double oldAction,newAction; // The old and new action
958 
959  /* Get a random vector */
960  dVec getRandomVector(const double);
961 
962  /* Returns a position used by the staging algorithm */
963  dVec newStagingPosition(const beadLocator &, const beadLocator &, const int, const int);
964 
965  /* Accumulate values */
966  void accumulate();
967 };
968 
969 // ========================================================================
970 // Pair Correlation Estimator Class
971 // ========================================================================
976 
977  public:
978  CylinderPairCorrelationEstimator(const Path &, ActionBase *, const MTRand &,
979  double, int _frequency=1, string _label="cyl_pair");
981 
982  void sample(); // Sample the estimator
983 
984  static const string name;
985  string getName() const {return name;}
986 
987  private:
988  void accumulate(); // Accumulate values
989  double dR; // The discretization
990 };
991 
992 // ========================================================================
993 // Cylinder Linear Potential Estimator Class
994 // ========================================================================
999 
1000  public:
1001  CylinderLinearPotentialEstimator(const Path &, ActionBase *, const MTRand &,
1002  double, int _frequency=1, string _label="cyl_linepotential");
1004 
1005  static const string name;
1006  string getName() const {return name;}
1007 
1008  private:
1009 
1010  double dz; // The discretization
1011  double Lz; // The length of the cylinder
1012 
1013  void accumulate(); // Accumulate values
1014 };
1015 
1016 // ========================================================================
1017 // Cylinder Radial Potential Estimator Class
1018 // ========================================================================
1023 
1024  public:
1025  CylinderRadialPotentialEstimator(const Path &, ActionBase *, const MTRand &,
1026  double, int _frequency=1, string _label="cyl_potential");
1028 
1029  static const string name;
1030  string getName() const {return name;}
1031 
1032  private:
1033 
1034  double dR; // The discretization
1035  Array <double,1> radPot; // Used for normalization
1036 
1037  void accumulate(); // Accumulate values
1038  void accumulate1(); // Accumulate values
1039 };
1040 
1041 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1042 // BEGIN PIGS ESTIMATORS ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1043 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1044 
1045 // ========================================================================
1046 // Potential Energy Estimator Class
1047 // ========================================================================
1053 
1054  public:
1055  PotentialEnergyEstimator(const Path &, ActionBase *, const MTRand &,
1056  double, int _frequency=1, string _label="potential");
1058 
1059  static const string name;
1060  string getName() const {return name;}
1061 
1062  private:
1063  void accumulate(); // Accumulate values
1064 
1065 };
1066 
1067 // ========================================================================
1068 // Kinetic Energy Estimator Class
1069 // ========================================================================
1074 
1075  public:
1076  KineticEnergyEstimator(const Path &, ActionBase *, const MTRand &,
1077  double, int _frequency=1, string _label="kinetic");
1079 
1080  static const string name;
1081  string getName() const {return name;}
1082 
1083  private:
1084  void accumulate(); // Accumulate values
1085 
1086 };
1087 
1088 // ========================================================================
1089 // Total Energy Estimator Class
1090 // ========================================================================
1095 
1096  public:
1097  TotalEnergyEstimator(const Path &, ActionBase *, const MTRand &,
1098  double, int _frequency=1, string _label="energy");
1100 
1101  static const string name;
1102  string getName() const {return name;}
1103 
1104  private:
1105  void accumulate(); // Accumulate values
1106 
1107 };
1108 
1109 // ========================================================================
1110 // Themodynamic Potential Energy Estimator Class
1111 // ========================================================================
1116 
1117  public:
1118  ThermoPotentialEnergyEstimator(const Path &, ActionBase *, const MTRand &,
1119  double, int _frequency=1, string _label="thpotential");
1121 
1122  static const string name;
1123  string getName() const {return name;}
1124 
1125  private:
1126  void accumulate(); // Accumulate values
1127 
1128 };
1129 
1130 // ========================================================================
1131 // Position Estimator Class
1132 // ========================================================================
1138 
1139  public:
1140  PositionEstimator(const Path &, ActionBase *, const MTRand &,
1141  double, int _frequency=1, string _label="position");
1143 
1144  static const string name;
1145  string getName() const {return name;}
1146 
1147  private:
1148  void accumulate(); // Accumulate values
1149 
1150 };
1151 
1152 // ========================================================================
1153 // Particle Resolved Position Estimator Class
1154 // ========================================================================
1160 
1161  public:
1162  ParticleResolvedPositionEstimator(const Path &, ActionBase *, const MTRand &,
1163  double, int _frequency=1, string _label="prposition");
1165 
1166  static const string name;
1167  string getName() const {return name;}
1168 
1169  private:
1170  void accumulate(); // Accumulate values
1171 
1172 };
1173 
1174 
1175 // ========================================================================
1176 // Particle Correlation Estimator Class
1177 // ========================================================================
1183 
1184  public:
1185  ParticleCorrelationEstimator(const Path &, ActionBase *, const MTRand &,
1186  double, int _frequency=1, string _label="prcorrelation");
1188 
1189  static const string name;
1190  string getName() const {return name;}
1191 
1192 
1193  private:
1194  void accumulate(); // Accumulate values
1195 
1196 };
1197 
1198 
1199 // ========================================================================
1200 // Velocity Estimator Class
1201 // ========================================================================
1207 
1208  public:
1209  VelocityEstimator(const Path &, ActionBase *, const MTRand &,
1210  double, int _frequency=1, string _label="velocity");
1212 
1213  static const string name;
1214  string getName() const {return name;}
1215 
1216  private:
1217  void accumulate(); // Accumulate values
1218 
1219 };
1220 
1221 // ========================================================================
1222 // SubregionOccupation Estimator Class
1223 // ========================================================================
1229 
1230  public:
1231  SubregionOccupationEstimator(const Path &, ActionBase *, const MTRand &,
1232  double, int _frequency=1, string _label="subregionocc");
1234 
1235  static const string name;
1236  string getName() const {return name;}
1237 
1238  private:
1239  void accumulate(); // Accumulate values
1240 
1241 };
1242 // ========================================================================
1243 // One Body Density Matrix Estimator Class
1244 // ========================================================================
1254 
1255  public:
1256  PIGSOneBodyDensityMatrixEstimator(Path &, ActionBase *, const MTRand &,
1257  double, int _frequency=1, string _label="obdm");
1259 
1260  void sample(); // Sample the estimator
1261  void outputFooter(); // Output the acceptance footer to disk
1262 
1263  static const string name;
1264  string getName() const {return name;}
1265 
1266  private:
1267  Path &lpath; // A non-constant local reference to the path
1268 
1269  double dR; // The discretization
1270  int numReps; // The number of measurments reps
1271  uint32 numAccepted; // The number of moves accepted
1272  uint32 numAttempted; // The number of moves attempted
1273 
1274 
1275  dVec newTailPos,oldTailPos; // The new and old tail position
1276  dVec newRanPos,neighborPos; // The random shift
1277 
1278  double sqrt2LambdaTau; // sqrt(2 * lambda * tau)
1279  double oldAction,newAction; // The old and new action
1280 
1281  /* Get a random vector */
1282  dVec getRandomVector(const double);
1283 
1284  /* Accumulate values */
1285  void accumulate();
1286 };
1287 
1288 
1289 // ========================================================================
1290 // Doubled Estimator Base Class
1291 // ========================================================================
1296 
1297  public:
1298  DoubledEstimator(const Path &, const Path &, ActionBase*, ActionBase*,
1299  const MTRand &, double, int _frequency=1, string _label="");
1301 
1302  string getName() const {return "doubled base";}
1303 
1304  protected:
1305  const Path &path2;
1306 };
1307 
1308 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1309 // END DOUBLED ESTIMATOR BASE CLASS ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1310 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1311 
1312 
1313 // ========================================================================
1314 // Swap Estimator Class
1315 // ========================================================================
1320 
1321  public:
1323  const MTRand &, double, int _frequency=1, string _label="swap");
1324  ~SwapEstimator();
1325 
1326  static const string name;
1327  string getName() const {return name;}
1328 
1329  private:
1330  Path &lpath; // A non-constant local reference to path 1
1331  Path &lpath2; // A non-constant local reference to path 2
1332  ActionBase *actionPtr; // The action pointer for path 1
1333  ActionBase *actionPtr2; // The action pointer for path 2
1334 
1335  void accumulate(); // Accumulate values
1336  void accumulateOpen(); // Accumulate values for open paths
1337  void accumulateClosed(); // Accumulate values for open paths
1338 };
1339 
1340 
1341 // ========================================================================
1342 // Swap Estimator Class
1343 // ========================================================================
1348 
1349  public:
1351  const MTRand &, double, int _frequency=1, string _label="entpart");
1353 
1354  static const string name;
1355  string getName() const {return name;}
1356 
1357  private:
1358  Path &lpath; // A non-constant local reference to path 1
1359  Path &lpath2; // A non-constant local reference to path 2
1360  ActionBase *actionPtr; // The action pointer for path 1
1361  ActionBase *actionPtr2; // The action pointer for path 2
1362 
1363  void accumulate(); // Accumulate values
1364 };
1365 
1366 
1367 
1368 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1369 // END PIGS ESTIMATORS ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1370 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
1371 #endif
Holds a base class that all action classes will be derived from.
Definition: action.h:29
Compute density inside film and in bulk separately for Excluded volume potentials.
Definition: estimator.h:265
~BipartitionDensityEstimator()
Destructor.
Definition: estimator.cpp:894
string getName() const
Get the name of the estimator.
Definition: estimator.h:273
BipartitionDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="bipart_dens")
Constructor.
Definition: estimator.cpp:878
Computes the total energy via the thermodynamic estimator.
Definition: estimator.h:813
CylinderEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_estimator")
Constructor.
Definition: estimator.cpp:2919
~CylinderEnergyEstimator()
Destructor.
Definition: estimator.cpp:2939
string getName() const
Get the name of the estimator.
Definition: estimator.h:821
Computes the density as a function of distance along the cylinder axis.
Definition: estimator.h:877
CylinderLinearDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_linedensity")
Constructor.
Definition: estimator.cpp:3120
~CylinderLinearDensityEstimator()
Destructor.
Definition: estimator.cpp:3145
string getName() const
Get the name of the estimator.
Definition: estimator.h:885
Compute the effective linear potential along the axis of the cylinder.
Definition: estimator.h:998
~CylinderLinearPotentialEstimator()
Destructor.
Definition: estimator.cpp:3628
string getName() const
Get the name of the estimator.
Definition: estimator.h:1006
CylinderLinearPotentialEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_linepotential")
Constructor.
Definition: estimator.cpp:3605
Computes the probability distribution function for the number of particles.
Definition: estimator.h:855
~CylinderNumberDistributionEstimator()
Destructor.
Definition: estimator.cpp:3093
string getName() const
Get the name of the estimator.
Definition: estimator.h:863
CylinderNumberDistributionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_number")
Constructor.
Definition: estimator.cpp:3075
Computes the average number of particles, as well as density.
Definition: estimator.h:835
~CylinderNumberParticlesEstimator()
Destructor.
Definition: estimator.cpp:3048
string getName() const
Get the name of the estimator.
Definition: estimator.h:843
CylinderNumberParticlesEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_estimator")
Constructor.
Definition: estimator.cpp:3032
Compute the one body density matrix n(r) which can be used to find the momentum distribution function...
Definition: estimator.h:930
CylinderOneBodyDensityMatrixEstimator(Path &, ActionBase *, const MTRand &, double, int _frequency=20, string _label="cyl_obdm")
Constructor.
Definition: estimator.cpp:3330
string getName() const
Get the name of the estimator.
Definition: estimator.h:940
Compute the two-body pair-correlation function, g(r) ~ <rho(r)rho(0)>.
Definition: estimator.h:975
string getName() const
Get the name of the estimator.
Definition: estimator.h:985
~CylinderPairCorrelationEstimator()
Destructor.
Definition: estimator.cpp:3559
void sample()
Sample the estimator.
Definition: estimator.cpp:3583
CylinderPairCorrelationEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_pair")
Constructor.
Definition: estimator.cpp:3536
Compute the effective radial potential in a cylinder.
Definition: estimator.h:1022
CylinderRadialPotentialEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_potential")
Constructor.
Definition: estimator.cpp:3692
~CylinderRadialPotentialEstimator()
Destructor.
Definition: estimator.cpp:3712
string getName() const
Get the name of the estimator.
Definition: estimator.h:1030
Compute the superfluid fraction, as well as the winding number probability distribution.
Definition: estimator.h:901
~CylinderSuperfluidFractionEstimator()
Destructor.
Definition: estimator.cpp:3211
string getName() const
Get the name of the estimator.
Definition: estimator.h:909
CylinderSuperfluidFractionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="cyl_super")
Constructor.
Definition: estimator.cpp:3186
Compute the fraction of time we spend in the diagonal ensemble.
Definition: estimator.h:568
void sample()
Overload sampling to make sure it is always done, regardless of ensemble.
Definition: estimator.cpp:1944
string getName() const
Get the name of the estimator.
Definition: estimator.h:578
~DiagonalFractionEstimator()
Destructor.
Definition: estimator.cpp:1929
DiagonalFractionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="estimator")
Constructor.
Definition: estimator.cpp:1915
Base class for estimators that use two paths.
Definition: estimator.h:1295
const Path & path2
A constant reference to the paths.
Definition: estimator.h:1305
~DoubledEstimator()
Destructor.
Definition: estimator.cpp:4576
string getName() const
Get the name of the estimator.
Definition: estimator.h:1302
DoubledEstimator(const Path &, const Path &, ActionBase *, ActionBase *, const MTRand &, double, int _frequency=1, string _label="")
Constructor.
Definition: estimator.cpp:4566
Computes the total energy via the thermodynamic estimator.
Definition: estimator.h:162
string getName() const
Get the name of the estimator.
Definition: estimator.h:170
EnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="estimator")
Constructor.
Definition: estimator.cpp:409
~EnergyEstimator()
Destructor.
Definition: estimator.cpp:420
Computes the Swap Estimator between two paths.
Definition: estimator.h:1347
EntPartEstimator(Path &, Path &, ActionBase *, ActionBase *, const MTRand &, double, int _frequency=1, string _label="entpart")
Constructor.
Definition: estimator.cpp:4889
string getName() const
Get the name of the estimator.
Definition: estimator.h:1355
~EntPartEstimator()
Destructor.
Definition: estimator.cpp:4916
The base class that all estimator classes will be derived from.
Definition: estimator.h:28
virtual void outputFlat()
Output a flat estimator value to disk.
Definition: estimator.cpp:305
int endSlice
Where imaginary time averages end.
Definition: estimator.h:96
bool diagonal
Is this a diagonal estimator?
Definition: estimator.h:106
void reset()
Reset numAccumulated and the estimator to 0.
Definition: estimator.cpp:261
Array< double, 1 > estimator
The estimator array.
Definition: estimator.h:90
virtual void outputFooter()
Ouptut the fooder to disk.
Definition: estimator.h:51
ActionBase * actionPtr
A pointer to the action.
Definition: estimator.h:82
uint32 getNumAccumulated() const
Get the number of accumulated measurements since the last reset.
Definition: estimator.h:60
virtual ~EstimatorBase()
Destructor.
Definition: estimator.cpp:146
Array< double, 1 > norm
The normalization factor for each estimator.
Definition: estimator.h:91
bool canonical
Are we in the canonical ensemble?
Definition: estimator.h:108
int endDiagSlice
Where imaginary time averages end for diagonal estimiators.
Definition: estimator.h:97
virtual void output()
Output the estimator value to disk.
Definition: estimator.cpp:283
vector< double > sliceFactor
Used to properly incorporate end affects.
Definition: estimator.h:98
void restart(const uint32, const uint32)
Restart the measurment process from a previous state.
Definition: estimator.cpp:269
void prepare()
Prepare the estimator for i/o.
Definition: estimator.cpp:240
void addEndLine()
Add a carriage return to estimator files.
Definition: estimator.h:72
string header
The data file header.
Definition: estimator.h:110
bool endLine
Should we output a carriage return?
Definition: estimator.h:107
int numEst
The number of individual quantities measured.
Definition: estimator.h:93
int startSlice
Where imaginary time averages begin.
Definition: estimator.h:95
int frequency
The frequency at which we accumulate.
Definition: estimator.h:94
void appendLabel(string append)
Append to default label.
Definition: estimator.cpp:325
map< string, int > estIndex
Map estimator labels to indices.
Definition: estimator.h:88
uint32 totNumAccumulated
The total number of accumulated values.
Definition: estimator.h:103
virtual string getName() const
Get the name of the estimator.
Definition: estimator.h:66
uint32 getTotNumAccumulated() const
Get the total number of accumulated measurments.
Definition: estimator.h:57
virtual void sample()
Sample the estimator.
Definition: estimator.cpp:187
EstimatorBase(const Path &_path, ActionBase *_actionPtr, const MTRand &_random, double _maxR, int _frequency=1, string _label="")
Constructor.
Definition: estimator.cpp:103
string getLabel() const
Get the estimator label.
Definition: estimator.h:78
uint32 numSampled
The number of times we have sampled.
Definition: estimator.h:101
virtual void accumulate()
Accumulate the estimator.
Definition: estimator.h:113
uint32 numAccumulated
The number of accumulated values.
Definition: estimator.h:102
uint32 getNumSampled() const
Get the number of samples since the last reset.
Definition: estimator.h:63
bool baseSample()
Determine the basic sampling condition.
Definition: estimator.cpp:161
const Path & path
A constant reference to the paths.
Definition: estimator.h:78
int numBeads0
The target number of beads for the canonical ensemble.
Definition: estimator.h:104
fstream * outFilePtr
The output fie.
Definition: estimator.h:86
void initialize(int)
Initialize estimator.
Definition: estimator.cpp:201
string label
The label used for the output file.
Definition: estimator.h:99
Compute the intermediate scattering function F(q,\tau)
Definition: estimator.h:734
string getName() const
Get the name of the estimator.
Definition: estimator.h:742
IntermediateScatteringFunctionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="isf")
Constructor.
Definition: estimator.cpp:2668
Computes the total energy using a mixed estimator.
Definition: estimator.h:1073
string getName() const
Get the name of the estimator.
Definition: estimator.h:1081
KineticEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="kinetic")
Constructor.
Definition: estimator.cpp:3890
~KineticEnergyEstimator()
Destructor.
Definition: estimator.cpp:3908
Create a 1d histogram of particle positions in the z-direction.
Definition: estimator.h:286
~LinearParticlePositionEstimator()
Destructor.
Definition: estimator.cpp:986
LinearParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="lineardensity")
Constructor.
Definition: estimator.cpp:955
string getName() const
Get the name of the estimator.
Definition: estimator.h:294
Particle permutation number density histogram.
Definition: estimator.h:613
string getName() const
Get the name of the estimator.
Definition: estimator.h:623
void output()
Output the estimator value to disk.
Definition: estimator.h:620
LocalPermutationEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="locperm")
Constructor.
Definition: estimator.cpp:2117
~LocalPermutationEstimator()
Destructor.
Definition: estimator.cpp:2140
Compute the local superfluid density.
Definition: estimator.h:439
LocalSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="locsuper")
Constructor.
Definition: estimator.cpp:1769
~LocalSuperfluidDensityEstimator()
Destructor.
Definition: estimator.cpp:1808
string getName() const
Get the name of the estimator.
Definition: estimator.h:449
void output()
overload the output
Definition: estimator.cpp:1818
Computes the probability distribution function for the number of particles.
Definition: estimator.h:390
~NumberDistributionEstimator()
Destructor.
Definition: estimator.cpp:788
string getName() const
Get the name of the estimator.
Definition: estimator.h:398
NumberDistributionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="number")
Constructor.
Definition: estimator.cpp:759
Computes the average number of particles, as well as density.
Definition: estimator.h:220
NumberParticlesEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="estimator")
Constructor.
Definition: estimator.cpp:719
~NumberParticlesEstimator()
Destructor.
Definition: estimator.cpp:733
string getName() const
Get the name of the estimator.
Definition: estimator.h:228
Compute the one body density matrix n(r) which can be used to find the momentum distribution function...
Definition: estimator.h:644
string getName() const
Get the name of the estimator.
Definition: estimator.h:655
void outputFooter()
For the one body density matrix estimator, we would like to output the acceptance information for the...
Definition: estimator.cpp:2478
OneBodyDensityMatrixEstimator(Path &, ActionBase *, const MTRand &, double, int _frequency=20, string _label="obdm")
Constructor.
Definition: estimator.cpp:2256
~OneBodyDensityMatrixEstimator()
Destructor.
Definition: estimator.cpp:2284
void sample()
Sample the OBDM.
Definition: estimator.cpp:2296
Compute the one body density matrix n(r) which can be used to find the momentum distribution function...
Definition: estimator.h:1253
void outputFooter()
For the one body density matrix estimator, we would like to output the acceptance information for the...
Definition: estimator.cpp:4550
~PIGSOneBodyDensityMatrixEstimator()
Destructor.
Definition: estimator.cpp:4404
string getName() const
Get the name of the estimator.
Definition: estimator.h:1264
PIGSOneBodyDensityMatrixEstimator(Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="obdm")
Constructor.
Definition: estimator.cpp:4376
void sample()
Sample the OBDM.
Definition: estimator.cpp:4414
Compute the two-body pair-correlation function, g(r) ~ <rho(r)rho(0)>.
Definition: estimator.h:690
PairCorrelationEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="pair")
Constructor.
Definition: estimator.cpp:2497
~PairCorrelationEstimator()
Destructor.
Definition: estimator.cpp:2545
string getName() const
Get the name of the estimator.
Definition: estimator.h:698
Computes the average position of each particle in 1D at the center time slice.
Definition: estimator.h:1182
~ParticleCorrelationEstimator()
Destructor.
Definition: estimator.cpp:4209
ParticleCorrelationEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="prcorrelation")
Constructor.
Definition: estimator.cpp:4191
string getName() const
Get the name of the estimator.
Definition: estimator.h:1190
Create histogram of particle positions.
Definition: estimator.h:241
ParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="position")
Constructor.
Definition: estimator.cpp:816
void output()
Output the estimator value to disk.
Definition: estimator.h:251
string getName() const
Get the name of the estimator.
Definition: estimator.h:249
~ParticlePositionEstimator()
Destructor.
Definition: estimator.cpp:843
Computes the average position of each particle in 1D at the center time slice.
Definition: estimator.h:1159
ParticleResolvedPositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="prposition")
Constructor.
Definition: estimator.cpp:4134
~ParticleResolvedPositionEstimator()
Destructor.
Definition: estimator.cpp:4155
string getName() const
Get the name of the estimator.
Definition: estimator.h:1167
The space-time trajectories.
Definition: path.h:29
Computes the particle permutation cycle probability distribution.
Definition: estimator.h:590
~PermutationCycleEstimator()
Destructor.
Definition: estimator.cpp:2037
PermutationCycleEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="pcycle")
Constructor.
Definition: estimator.cpp:2015
string getName() const
Get the name of the estimator.
Definition: estimator.h:598
Compute the radially averaged local superfluid density.
Definition: estimator.h:494
~PlaneAreaSuperfluidDensityEstimator()
Destructor.
Definition: estimator.cpp:1530
string getName() const
Get the name of the estimator.
Definition: estimator.h:502
PlaneAreaSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="planearea")
Constructor.
Definition: estimator.cpp:1500
Create a 2d histogram of particle positions but only store the average.
Definition: estimator.h:362
PlaneAverageExternalPotentialEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="planeaveVext")
Constructor.
Definition: estimator.cpp:1191
string getName() const
Get the name of the estimator.
Definition: estimator.h:370
void output()
Output a flat estimator value to disk.
Definition: estimator.cpp:1261
Create a 2d histogram of particle positions but only store the average.
Definition: estimator.h:335
void output()
Output the estimator value to disk.
Definition: estimator.h:345
PlaneParticleAveragePositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="planeavedensity")
Constructor.
Definition: estimator.cpp:1112
string getName() const
Get the name of the estimator.
Definition: estimator.h:343
Create a 2d histogram of particle positions.
Definition: estimator.h:310
string getName() const
Get the name of the estimator.
Definition: estimator.h:318
PlaneParticlePositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="planedensity")
Constructor.
Definition: estimator.cpp:1029
~PlaneParticlePositionEstimator()
Destructor.
Definition: estimator.cpp:1066
Compute the radially averaged local superfluid density.
Definition: estimator.h:468
PlaneWindingSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="planewind")
Constructor.
Definition: estimator.cpp:1411
string getName() const
Get the name of the estimator.
Definition: estimator.h:476
Computes the average value of the position in 1D.
Definition: estimator.h:1137
string getName() const
Get the name of the estimator.
Definition: estimator.h:1145
~PositionEstimator()
Destructor.
Definition: estimator.cpp:4102
PositionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="position")
Constructor.
Definition: estimator.cpp:4084
Computes the potential energy along the worldline.
Definition: estimator.h:1052
PotentialEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="potential")
Constructor.
Definition: estimator.cpp:3844
string getName() const
Get the name of the estimator.
Definition: estimator.h:1060
~PotentialEnergyEstimator()
Destructor.
Definition: estimator.cpp:3861
Compute the radially averaged local superfluid density.
Definition: estimator.h:544
string getName() const
Get the name of the estimator.
Definition: estimator.h:552
RadialAreaSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="radarea")
Constructor.
Definition: estimator.cpp:1681
~RadialAreaSuperfluidDensityEstimator()
Destructor.
Definition: estimator.cpp:1710
Compute the density as a function of position in the radial direction.
Definition: estimator.h:759
string getName() const
Get the name of the estimator.
Definition: estimator.h:767
~RadialDensityEstimator()
Destructor.
Definition: estimator.cpp:2853
RadialDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="radial")
Constructor.
Definition: estimator.cpp:2829
Compute the radially averaged local superfluid density.
Definition: estimator.h:520
string getName() const
Get the name of the estimator.
Definition: estimator.h:528
RadialWindingSuperfluidDensityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="radwind")
Constructor.
Definition: estimator.cpp:1595
Compute the static structure factor S(q)
Definition: estimator.h:711
string getName() const
Get the name of the estimator.
Definition: estimator.h:719
StaticStructureFactorEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="ssf")
Constructor.
Definition: estimator.cpp:2577
~StaticStructureFactorEstimator()
Destructor.
Definition: estimator.cpp:2607
Computes the imaginary time resolved "velocity" for the first particle .
Definition: estimator.h:1228
~SubregionOccupationEstimator()
Destructor.
Definition: estimator.cpp:4322
SubregionOccupationEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="subregionocc")
Constructor.
Definition: estimator.cpp:4304
string getName() const
Get the name of the estimator.
Definition: estimator.h:1236
Compute the superfluid fraction, as well as the winding number probability distribution.
Definition: estimator.h:416
SuperfluidFractionEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="super")
Constructor.
Definition: estimator.cpp:1297
~SuperfluidFractionEstimator()
Destructor.
Definition: estimator.cpp:1324
string getName() const
Get the name of the estimator.
Definition: estimator.h:424
Computes the Swap Estimator between two paths.
Definition: estimator.h:1319
SwapEstimator(Path &, Path &, ActionBase *, ActionBase *, const MTRand &, double, int _frequency=1, string _label="swap")
Constructor.
Definition: estimator.cpp:4589
~SwapEstimator()
Destructor.
Definition: estimator.cpp:4606
string getName() const
Get the name of the estimator.
Definition: estimator.h:1327
Computes the total energy using a mixed estimator.
Definition: estimator.h:1115
~ThermoPotentialEnergyEstimator()
Destructor.
Definition: estimator.cpp:4055
string getName() const
Get the name of the estimator.
Definition: estimator.h:1123
ThermoPotentialEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="thpotential")
Constructor.
Definition: estimator.cpp:4037
An estimator which tracks the ammount of time between bins, summing them into a total at the end.
Definition: estimator.h:128
void sample()
Overload sampling to make sure it is always done, regardless of ensemble.
Definition: estimator.cpp:354
TimeEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="estimator")
Constructor.
Definition: estimator.cpp:340
string getName() const
Get the name of the estimator.
Definition: estimator.h:138
void output()
Grab the final time and write to disk.
Definition: estimator.cpp:376
Computes the total energy using a mixed estimator.
Definition: estimator.h:1094
~TotalEnergyEstimator()
Destructor.
Definition: estimator.cpp:3985
TotalEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="energy")
Constructor.
Definition: estimator.cpp:3967
string getName() const
Get the name of the estimator.
Definition: estimator.h:1102
Computes the imaginary time resolved "velocity" for the first particle .
Definition: estimator.h:1206
~VelocityEstimator()
Destructor.
Definition: estimator.cpp:4267
VelocityEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="velocity")
Constructor.
Definition: estimator.cpp:4249
string getName() const
Get the name of the estimator.
Definition: estimator.h:1214
Computes the total energy via the thermodynamic estimator.
Definition: estimator.h:199
VirialEnergyEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="virial")
Constructor.
Definition: estimator.cpp:523
string getName() const
Get the name of the estimator.
Definition: estimator.h:207
~VirialEnergyEstimator()
Destructor.
Definition: estimator.cpp:539
Compute various properties related to the worm in the simulation.
Definition: estimator.h:780
string getName() const
Get the name of the estimator.
Definition: estimator.h:788
WormPropertiesEstimator(const Path &, ActionBase *, const MTRand &, double, int _frequency=1, string _label="worm")
Constructor.
Definition: estimator.cpp:1968
~WormPropertiesEstimator()
Destructor.
Definition: estimator.cpp:1985
Global common header with shared dependencies and methods.
unsigned long uint32
Unsigned integer type, at least 32 bits.
Definition: common.h:105
TinyVector< int, 2 > beadLocator
time-slice,bead-number world line index
Definition: common.h:117
TinyVector< double, NDIM > dVec
A NDIM-vector of type double.
Definition: common.h:111