Path Integral Quantum Monte Carlo
worm.h
Go to the documentation of this file.
1 
8 #ifndef WORM_H
9 #define WORM_H
10 
11 #include "common.h"
12 #include "constants.h"
13 
14 class Path;
15 // ========================================================================
16 // Worm Class
17 // ========================================================================
25 class Worm {
26 
27  public:
28  Worm(int);
29  ~Worm();
30 
35  double maxWormCost;
37  int length;
38  int gap;
40 
41  /* Safe get methods for beads */
42  inline int beadOn(int,int) const;
43  inline int beadOn(const beadLocator&) const;
44 
45  /* Get the state of a bead */
46  beadState getState (const beadLocator &) const;
47 
48  /* The worm-trajectory factor for interactions and external potentials*/
49  double factor(const beadState, const beadLocator&) const;
50  double factor(const beadState state1) const { return 1.0 - 0.5*(state1 != NONE);};
51 
52  /* Safely add/delete beads */
53  inline void delBead(int,int);
54  inline void addBead(int,int);
55  inline void delBead(const beadLocator&);
56  inline void addBead(const beadLocator&);
57 
58  /* Reset all worm parameters to null */
59  void reset();
60 
61  /* Update all worm parameters */
62  void update(Path &, const beadLocator &, const beadLocator &);
63 
64  /* Determine if a given beadLocator is on a worm */
65  bool foundBead(const Path &, const beadLocator &);
66 
68  inline bool tooCostly() {
69  if (gap == 0)
70  return true;
71  else
72  return ((dot(sep,sep) * constants()->fourLambdaTauInv() / (1.0*gap)) > maxWormCost);
73  }
75  inline bool tooCostly(const dVec& _sep, int _gap) {
76  if (_gap == 0)
77  return true;
78  else {
79  return ((dot(_sep,_sep) * constants()->fourLambdaTauInv() / (1.0*_gap)) > maxWormCost);
80  }
81  }
82 
84  const Array <unsigned int, 2> & getBeads() const { return beads; }
85 
86  /* Test whether a given bead is on a worm */
88  int getNumBeadsOn() const {return numBeadsOn;}
90  void incNumBeadsOn() {++numBeadsOn;}
92  void decNumBeadsOn() {--numBeadsOn;}
94  void resetNumBeadsOn() {numBeadsOn = sum(beads);}
95 
96  friend class Path; // Path needs access to beads
97  friend class PathIntegralMonteCarlo; // Friends for I/O
98 
99  private:
100  Array <unsigned int,2> beads; // Is a bead present?
101  int numBeadsOn; // How many beads are present
102 };
103 
104 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
105 // INLINE FUNCTION DEFINITIONS
106 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
107 
108 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
110 inline beadState Worm::getState(const beadLocator &beadIndex) const {
111  if (all(beadIndex==head) || all(beadIndex==tail))
112  return HEADTAIL;
113  else if (all(beadIndex==special1) || all(beadIndex==special2))
114  return SPECIAL;
115  else
116  return NONE;
117 }
118 
119 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
120 
122 inline void Worm::delBead(int slice, int ptcl) {
123  PIMC_ASSERT((slice >= 0) && (slice < constants()->numTimeSlices()));
124  beads(slice,ptcl) = 0;
125 }
126 
128 inline void Worm::delBead(const beadLocator &beadIndex) {
129  PIMC_ASSERT((beadIndex[0] >= 0) && (beadIndex[0] < constants()->numTimeSlices()));
130  beads(beadIndex) = 0;
131 }
132 
133 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
135 inline void Worm::addBead(int slice, int ptcl) {
136  PIMC_ASSERT((slice >= 0) && (slice < constants()->numTimeSlices()));
137  beads(slice,ptcl) = 1;
138 }
139 
141 inline void Worm::addBead(const beadLocator &beadIndex) {
142  PIMC_ASSERT((beadIndex[0] >= 0) && (beadIndex[0] < constants()->numTimeSlices()));
143  beads(beadIndex) = 1;
144 }
145 
146 // -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=
148 inline int Worm::beadOn(int slice, int ptcl) const {
149  PIMC_ASSERT((slice >=0) && (slice < constants()->numTimeSlices()));
150  return beads(slice,ptcl);
151 }
152 
154 inline int Worm::beadOn(const beadLocator &beadIndex) const {
155  return beads(beadIndex);
156 }
157 
158 #endif
The main driver class for the entire path integral monte carlo program.
Definition: pimc.h:37
The space-time trajectories.
Definition: path.h:29
Contains information on the worm.
Definition: worm.h:25
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
Worm(int)
Constructor.
Definition: worm.cpp:22
dVec sep
The spatial separation between head and tail.
Definition: worm.h:36
beadLocator head
The coordinates of the worm head.
Definition: worm.h:31
beadLocator special2
Special bead, used in move updates.
Definition: worm.h:34
void delBead(int, int)
Safely delete a bead (int indexed)
Definition: worm.h:122
bool tooCostly(const dVec &_sep, int _gap)
Return true if the worm is too costly.
Definition: worm.h:75
double maxWormCost
The maximum 'cost' of inserting a worm.
Definition: worm.h:35
void update(Path &, const beadLocator &, const beadLocator &)
We update all worm properties for a new head and tail.
Definition: worm.cpp:69
int beadOn(int, int) const
Safely get a bead (int indexed)
Definition: worm.h:148
void incNumBeadsOn()
Increment the number of active beads.
Definition: worm.h:90
bool tooCostly()
Return true if the worm is too costly.
Definition: worm.h:68
beadState getState(const beadLocator &) const
Get the state of the supplied bead?
Definition: worm.h:110
bool foundBead(const Path &, const beadLocator &)
Test to see if a supplied bead is located on a worm.
Definition: worm.cpp:124
void resetNumBeadsOn()
Reset the number of active beads.
Definition: worm.h:94
void decNumBeadsOn()
Decrement the number of active beads.
Definition: worm.h:92
double factor(const beadState, const beadLocator &) const
Compute the value of the potential action trajectory factor.
Definition: worm.cpp:109
void reset()
Reset the worm to a null state.
Definition: worm.cpp:53
int getNumBeadsOn() const
Return the number of active beads.
Definition: worm.h:88
~Worm()
Destructor.
Definition: worm.cpp:46
const Array< unsigned int, 2 > & getBeads() const
Return the bead list.
Definition: worm.h:84
void addBead(int, int)
Safely add a bead (int indexed)
Definition: worm.h:135
beadLocator special1
Special bead, used in move updates.
Definition: worm.h:33
int length
The length of the worm.
Definition: worm.h:37
Global common header with shared dependencies and methods.
TinyVector< int, 2 > beadLocator
time-slice,bead-number world line index
Definition: common.h:117
beadState
Each bead can have three possible states.
Definition: common.h:129
TinyVector< double, NDIM > dVec
A NDIM-vector of type double.
Definition: common.h:111
#define PIMC_ASSERT(X)
Rename assert method.
Definition: common.h:64
ConstantParameters class definition.
ConstantParameters * constants()
Global public access to the constants.
Definition: constants.h:201