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

A derived class which performs a staging move, which exactly samples the kinetic action. More...

#include <move.h>

+ Inheritance diagram for StagingMove:
+ Collaboration diagram for StagingMove:

Public Member Functions

 StagingMove (Path &, ActionBase *, MTRand &, ensemble _operateOnConfig=ANY)
 Constructor.
 
 ~StagingMove ()
 Destructor.
 
bool attemptMove ()
 Perform the staging move. More...
 
string getName ()
 return the move name
 
- Public Member Functions inherited from MoveBase
 MoveBase (Path &, ActionBase *, MTRand &, ensemble _operateOnConfig=ANY, bool _varLength=false)
 Move naming conventions: More...
 
virtual ~MoveBase ()
 Destructor.
 
double getAcceptanceRatio ()
 Get the acceptance ratio.
 
double getTotAcceptanceRatio ()
 Get the total acceptance ratio.
 
double getAcceptanceRatioLevel (int n)
 Get the acceptance ratio by level.
 
int getNumAttempted ()
 Get the number of moves attempted.
 
int getNumAccepted ()
 Get the number of moves accepted.
 
int getNumAttemptedLevel (int n)
 Get the number of moves attempted by level.
 
int getNumAcceptedLevel (int n)
 Get the number of moves accepted by level.
 
void resetTotAccept ()
 Reset the total accepted counter.
 
void resetAccept ()
 Reset the number accepted counter.
 

Static Public Attributes

static const string name
 

Additional Inherited Members

- Data Fields inherited from MoveBase
ensemble operateOnConfig
 What configurations do we operate on?
 
bool variableLength
 Does the move have a variable length?
 
string name1
 
- Protected Member Functions inherited from MoveBase
virtual void keepMove ()
 Keep the move. More...
 
dVec newStagingPosition (const beadLocator &, const beadLocator &, const int, const int)
 Returns a new staging position which will exactly sample the kinetic action. More...
 
dVec newStagingPosition (const beadLocator &, const beadLocator &, const int, const int, iVec &)
 Returns a new staging position which will exactly sample the kinetic action in different winding sectors. More...
 
iVec sampleWindingSector (const beadLocator &, const beadLocator &, const int, double &)
 Obtain a winding sector for a stage-like move. More...
 
iVec getWindingNumber (const beadLocator &, const beadLocator &)
 Find the winding number for a path between two beads. More...
 
dVec newFreeParticlePosition (const beadLocator &)
 Generates a new position, which exactly samples the free particle density matrix. More...
 
dVec newBisectionPosition (const beadLocator &, const int)
 Returns a new bisection position which will exactly sample the kinetic action. More...
 
void printMoveState (string)
 
void checkMove (int, double)
 
- Protected Attributes inherited from MoveBase
Pathpath
 A reference to the paths.
 
ActionBaseactionPtr
 A base pointer to the action.
 
MTRand & random
 A reference to the RNG.
 
bool success
 Did we sucessfully perform a move?
 
uint32 numAccepted
 The number of accepted moves.
 
uint32 numAttempted
 The number of attempted moves.
 
int numToMove
 The number of particles moved.
 
int numLevels
 
Array< uint32, 1 > numAcceptedLevel
 The number of moves accepted at each level.
 
Array< uint32, 1 > numAttemptedLevel
 The number of moves attempted at each level.
 
Array< dVec, 1 > originalPos
 The original particle positions.
 
Array< dVec, 1 > newPos
 New particle positions.
 
vector< iVecwinding
 The winding vectors

 
vector< int > windingSector
 Used to index different winding sectors.
 
vector< double > cumrho0
 Used for tower-sampling winding sectors.
 
int maxWind
 The largest winding number.
 
int numWind
 The total number of winding vectors.
 
double oldAction
 The original potential action.
 
double newAction
 The new potential action.
 
double deltaAction
 The action difference.
 
double sqrt2LambdaTau
 sqrt(2 * Lambda * tau)
 
double sqrtLambdaTau
 sqrt(Lambda * tau)
 
beadLocator nBeadIndex
 Neighbor bead index.
 
dVec neighborPos
 Staging neighbor position.
 
dVec newRanPos
 Staing random position.
 
double newK
 
double oldK
 The old and new kinetic action.
 
double newV
 
double oldV
 The old and new potential action.
 
- Static Protected Attributes inherited from MoveBase
static uint32 totAccepted = 0
 The total number of moves accepted.
 
static uint32 totAttempted = 0
 The total number of moves attempted.
 

Detailed Description

A derived class which performs a staging move, which exactly samples the kinetic action.

Definition at line 271 of file move.h.

Member Function Documentation

◆ attemptMove()

bool StagingMove::attemptMove ( )
virtual

Perform the staging move.

The staging move is non-local in imaginary time and exactly samples the kinetic density matrix.

Implements MoveBase.

Definition at line 1316 of file move.cpp.

1316  {
1317 
1318  success = false;
1319 
1320  /* Only perform a move if we have beads */
1321  if (path.worm.getNumBeadsOn() == 0)
1322  return success;
1323 
1324  checkMove(0,0.0);
1325 
1326  /* We don't bother staging if we only have a worm, as this will be taken care
1327  * of with other moves. It is also possible that we could get stuck in an
1328  * infinite loop here.*/
1329  if (path.getTrueNumParticles()==0)
1330  return false;
1331 
1332  /* Randomly select the start bead of the stage */
1333 #if PIGS
1334  startBead[0] = random.randInt((path.numTimeSlices-1)-constants()->Mbar());
1335 #else
1336  startBead[0] = random.randInt(path.numTimeSlices-1);
1337 #endif
1338 
1339  /* We need to worry about the possibility of an empty slice for small
1340  * numbers of particles. */
1341  if (path.numBeadsAtSlice(startBead[0]) == 0)
1342  return false;
1343  startBead[1] = random.randInt(path.numBeadsAtSlice(startBead[0])-1);
1344 
1345  /* Do we perform varialble length staging updates? */
1346  if (constants()->varUpdates())
1347  stageLength = 2 + random.randInt(constants()->Mbar()-2);
1348 
1349  /* Now we have to make sure that we are moving an active trajectory,
1350  * otherwise we exit immediatly */
1351  beadLocator beadIndex;
1352  beadIndex = startBead;
1353  for (int k = 0; k < (stageLength); k++) {
1354  if (!path.worm.beadOn(beadIndex) || all(path.next(beadIndex)==XXX))
1355  return false;
1356  beadIndex = path.next(beadIndex);
1357  }
1358  endBead = beadIndex;
1359 
1360  /* If we haven't found the worm head, and all beads are on, try to perform the move */
1361 
1362  /* Increment the attempted counters */
1363  numAttempted++;
1364  totAttempted++;
1365 
1366  double totalrho0;
1367  iVec wind;
1368  wind = sampleWindingSector(startBead,endBead,stageLength,totalrho0);
1369 
1370  /* Get the current action for the path segment to be updated */
1371  oldAction = actionPtr->potentialAction(startBead,path.prev(endBead));
1372 
1373  /* Perform the staging update, generating the new path and updating bead
1374  * positions, while storing the old one */
1375  beadIndex = startBead;
1376  int k = 0;
1377  dVec pos;
1378  bool movedIntoSubRegionA = false;
1379  do {
1380  beadIndex = path.next(beadIndex);
1381  originalPos(k) = path(beadIndex);
1382  path.updateBead(beadIndex,
1383  newStagingPosition(path.prev(beadIndex),endBead,stageLength,k,wind));
1384  if (!movedIntoSubRegionA){
1385  movedIntoSubRegionA = path.inSubregionA(beadIndex);
1386  }
1387  ++k;
1388  } while (!all(beadIndex==path.prev(endBead)));
1389 
1390  if ( !movedIntoSubRegionA ) {
1391  /* Get the new action for the updated path segment */
1392  newAction = actionPtr->potentialAction(startBead,path.prev(endBead));
1393 
1394  /* The actual Metropolis test */
1395  if ( random.rand() < exp(-(newAction-oldAction)) ) {
1396  keepMove();
1397  checkMove(1,newAction-oldAction);
1398  }
1399  else {
1400  undoMove();
1401  checkMove(2,0.0);
1402  }
1403  }
1404  else {
1405  undoMove();
1406  }
1407 
1408  return success;
1409 }
virtual double potentialAction()
The effective potential inter-ACTION for various pass conditions.
Definition: action.h:48
uint32 numAttempted
The number of attempted moves.
Definition: move.h:87
Array< dVec, 1 > originalPos
The original particle positions.
Definition: move.h:97
Path & path
A reference to the paths.
Definition: move.h:80
dVec newStagingPosition(const beadLocator &, const beadLocator &, const int, const int)
Returns a new staging position which will exactly sample the kinetic action.
Definition: move.cpp:268
iVec sampleWindingSector(const beadLocator &, const beadLocator &, const int, double &)
Obtain a winding sector for a stage-like move.
Definition: move.cpp:357
MTRand & random
A reference to the RNG.
Definition: move.h:82
ActionBase * actionPtr
A base pointer to the action.
Definition: move.h:81
virtual void keepMove()
Keep the move.
Definition: move.cpp:248
double newAction
The new potential action.
Definition: move.h:108
bool success
Did we sucessfully perform a move?
Definition: move.h:84
double oldAction
The original potential action.
Definition: move.h:107
static uint32 totAttempted
The total number of moves attempted.
Definition: move.h:92
bool inSubregionA(const beadLocator &) const
Checks to see if bead is in subregion A/B at break slice + 1.
Definition: path.cpp:518
const int numTimeSlices
A local constant copy of the number of time slices.
Definition: path.h:37
Array< int, 1 > numBeadsAtSlice
The number of active beads at a given time slice.
Definition: path.h:48
void updateBead(const beadLocator &, const dVec &)
Update the position of a bead in the worldine configuration.
Definition: path.cpp:182
beadLocator & next(int slice, int ptcl)
Move one link forward in imaginary time.
Definition: path.h:86
Worm worm
Details on the worm.
Definition: path.h:44
beadLocator & prev(int slice, int ptcl)
Move one link backward in imaginary time.
Definition: path.h:95
int getTrueNumParticles() const
The number of active particles.
Definition: path.h:54
int beadOn(int, int) const
Safely get a bead (int indexed)
Definition: worm.h:148
int getNumBeadsOn() const
Return the number of active beads.
Definition: worm.h:88
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
TinyVector< int, NDIM > iVec
A NDIM-vector of type integer.
Definition: common.h:114
#define XXX
Used to refer to a nonsense beadIndex.
Definition: common.h:98
ConstantParameters * constants()
Global public access to the constants.
Definition: constants.h:201
+ Here is the call graph for this function:

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