mirror of https://github.com/QMCPACK/qmcpack.git
RMC now supports resizing reptiles.
git-svn-id: https://subversion.assembla.com/svn/qmcdev/trunk@6374 e5b18d87-469d-4833-9cc0-8cdfa06e9491
This commit is contained in:
Raymond Clay
parent
c452e0fb6b
commit
22ff3bc463
|
|
@ -17,21 +17,21 @@
|
|||
|
||||
namespace qmcplusplus
|
||||
{
|
||||
class MCWalkerConfiguration;
|
||||
class MCWalkerConfiguration;
|
||||
|
||||
class Reptile: public QMCTraits
|
||||
{
|
||||
class Reptile:public QMCTraits
|
||||
{
|
||||
|
||||
|
||||
public:
|
||||
public:
|
||||
typedef MCWalkerConfiguration::Walker_t Walker_t;
|
||||
//typedef Walker_t::Buffer_t Buffer_t;
|
||||
// typedef MCWalkerConfiguration::Walker_t Walker_t;
|
||||
typedef MCWalkerConfiguration::iterator WalkerIter_t;
|
||||
typedef vector<Walker_t::ParticlePos_t> ReptileConfig_t;
|
||||
typedef vector < Walker_t::ParticlePos_t > ReptileConfig_t;
|
||||
|
||||
std::vector<IndexType> Action;
|
||||
std::vector<IndexType> TransProb;
|
||||
std::vector < IndexType > Action;
|
||||
std::vector < IndexType > TransProb;
|
||||
|
||||
RealType forwardprob;
|
||||
RealType backwardprob;
|
||||
|
|
@ -41,235 +41,274 @@ public:
|
|||
RealType tau;
|
||||
|
||||
IndexType direction, headindex, nbeads;
|
||||
MCWalkerConfiguration& w;
|
||||
Walker_t* prophead;
|
||||
MCWalkerConfiguration & w;
|
||||
Walker_t *prophead;
|
||||
WalkerIter_t repstart, repend;
|
||||
|
||||
inline Reptile(MCWalkerConfiguration& W, WalkerIter_t start, WalkerIter_t end):
|
||||
w(W),repstart(start),repend(end),direction(1),headindex(0),prophead(0) //, r2prop(0.0), r2accept(0.0),tau(0.0)
|
||||
inline Reptile (MCWalkerConfiguration & W, WalkerIter_t start, WalkerIter_t end):w (W), repstart (start), repend (end), direction (1), headindex (0), prophead (0) //, r2prop(0.0), r2accept(0.0),tau(0.0)
|
||||
{
|
||||
w=W;
|
||||
Action.resize(3);
|
||||
Action[0]=w.addProperty("ActionBackward");
|
||||
Action[1]=w.addProperty("ActionForward");
|
||||
Action[2]=w.addProperty("ActionLocal");
|
||||
TransProb.resize(2);
|
||||
TransProb[0]=w.addProperty("TransProbBackward");
|
||||
TransProb[1]=w.addProperty("TransProbForward");
|
||||
w = W;
|
||||
Action.resize (3);
|
||||
Action[0] = w.addProperty ("ActionBackward");
|
||||
Action[1] = w.addProperty ("ActionForward");
|
||||
Action[2] = w.addProperty ("ActionLocal");
|
||||
TransProb.resize (2);
|
||||
TransProb[0] = w.addProperty ("TransProbBackward");
|
||||
TransProb[1] = w.addProperty ("TransProbForward");
|
||||
|
||||
nbeads=repend-repstart;
|
||||
nbeads = repend - repstart;
|
||||
}
|
||||
|
||||
~Reptile() {}
|
||||
|
||||
inline IndexType size(){return nbeads;}
|
||||
|
||||
inline Walker_t& operator[](IndexType i)
|
||||
~Reptile ()
|
||||
{
|
||||
return getWalker(getBeadIndex(i)) ;
|
||||
}
|
||||
|
||||
inline IndexType wrapIndex(IndexType repindex)
|
||||
inline IndexType size ()
|
||||
{
|
||||
return (repindex%nbeads + nbeads)%nbeads;
|
||||
return nbeads;
|
||||
}
|
||||
|
||||
inline Walker_t& getWalker(IndexType i)
|
||||
inline Walker_t & operator[] (IndexType i)
|
||||
{
|
||||
WalkerIter_t bead = repstart + wrapIndex(i);
|
||||
return getWalker (getBeadIndex (i));
|
||||
}
|
||||
|
||||
inline IndexType wrapIndex (IndexType repindex)
|
||||
{
|
||||
return (repindex % nbeads + nbeads) % nbeads;
|
||||
}
|
||||
|
||||
inline Walker_t & getWalker (IndexType i)
|
||||
{
|
||||
WalkerIter_t bead = repstart + wrapIndex (i);
|
||||
return **bead;
|
||||
}
|
||||
|
||||
inline IndexType getBeadIndex(IndexType i)
|
||||
inline IndexType getBeadIndex (IndexType i)
|
||||
{
|
||||
return wrapIndex(headindex + direction*i) ;
|
||||
return wrapIndex (headindex + direction * i);
|
||||
}
|
||||
inline Walker_t& getBead(IndexType i)
|
||||
inline Walker_t & getBead (IndexType i)
|
||||
{
|
||||
return getWalker(getBeadIndex(i)) ;
|
||||
return getWalker (getBeadIndex (i));
|
||||
}
|
||||
inline Walker_t& getHead()
|
||||
inline Walker_t & getHead ()
|
||||
{
|
||||
return getWalker(getBeadIndex(0)) ;
|
||||
return getWalker (getBeadIndex (0));
|
||||
}
|
||||
inline Walker_t& getTail()
|
||||
inline Walker_t & getTail ()
|
||||
{
|
||||
return getWalker(getBeadIndex(nbeads-1)) ;
|
||||
return getWalker (getBeadIndex (nbeads - 1));
|
||||
}
|
||||
inline Walker_t& getNext()
|
||||
inline Walker_t & getNext ()
|
||||
{
|
||||
return getWalker(getBeadIndex(nbeads-2)) ;
|
||||
return getWalker (getBeadIndex (nbeads - 2));
|
||||
}
|
||||
inline Walker_t& getCenter()
|
||||
inline Walker_t & getCenter ()
|
||||
{
|
||||
return getWalker(getBeadIndex( (nbeads-1)/2 ));
|
||||
return getWalker (getBeadIndex ((nbeads - 1) / 2));
|
||||
}
|
||||
//inline void setProposedHead(){
|
||||
|
||||
inline void flip()
|
||||
inline void flip ()
|
||||
{
|
||||
// direction*=-1;
|
||||
// headindex = getBeadIndex(nbeads-1);
|
||||
headindex = wrapIndex(headindex - direction);
|
||||
direction*=-1;
|
||||
headindex = wrapIndex (headindex - direction);
|
||||
direction *= -1;
|
||||
}
|
||||
|
||||
inline void setDirection(IndexType dir)
|
||||
inline void setDirection (IndexType dir)
|
||||
{
|
||||
direction=dir;
|
||||
direction = dir;
|
||||
}
|
||||
|
||||
inline void setBead(Walker_t& walker, IndexType i)
|
||||
inline void setBead (Walker_t & walker, IndexType i)
|
||||
{
|
||||
IndexType index = getBeadIndex(i);
|
||||
Walker_t& newbead(getWalker(index));
|
||||
newbead=walker; //This should be a hard copy
|
||||
IndexType index = getBeadIndex (i);
|
||||
Walker_t & newbead (getWalker (index));
|
||||
newbead = walker; //This should be a hard copy
|
||||
}
|
||||
|
||||
inline void setHead(Walker_t& overwrite)
|
||||
inline void setHead (Walker_t & overwrite)
|
||||
{
|
||||
//overwrite last element.
|
||||
headindex = getBeadIndex(nbeads-1); //sets to position of tail.
|
||||
Walker_t& newhead(getBead(0));
|
||||
newhead=overwrite;
|
||||
headindex = getBeadIndex (nbeads - 1); //sets to position of tail.
|
||||
Walker_t & newhead (getBead (0));
|
||||
newhead = overwrite;
|
||||
}
|
||||
//This function does two things: 1.) Moves the reptile forward 1 step. 2.) Returns the new head.
|
||||
inline Walker_t& getNewHead()
|
||||
inline Walker_t & getNewHead ()
|
||||
{
|
||||
//overwrite last element.
|
||||
headindex = getBeadIndex(nbeads-1); //sets to position of tail.
|
||||
return getWalker(headindex);
|
||||
headindex = getBeadIndex (nbeads - 1); //sets to position of tail.
|
||||
return getWalker (headindex);
|
||||
}
|
||||
|
||||
void saveAction(Walker_t& walker, IndexType d, RealType val, IndexType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType actionindex= 2;
|
||||
if (direction!=0) actionindex=(1-d*direction)/2;
|
||||
walker.Properties(nPsi,Action[actionindex])=val;
|
||||
}
|
||||
|
||||
RealType getDirectionalAction(Walker_t& walker, IndexType d, IndexType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType actionindex= 2;
|
||||
if (d!=0) actionindex=(1-direction*d)/2;
|
||||
|
||||
return walker.Properties(nPsi,Action[actionindex]);
|
||||
}
|
||||
|
||||
RealType getLinkAction(Walker_t& new_walker, Walker_t& old_walker, IndexType d, IndexType nPsi=0)
|
||||
void saveAction (Walker_t & walker, IndexType d, RealType val,
|
||||
IndexType nPsi = 0)
|
||||
{
|
||||
RealType af=getDirectionalAction(old_walker, +1, nPsi);
|
||||
RealType ab=getDirectionalAction(new_walker, -1, nPsi);
|
||||
RealType a0=getDirectionalAction(old_walker, 0, nPsi) + getDirectionalAction(new_walker,0,nPsi);
|
||||
return af+ab+a0;
|
||||
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType actionindex = 2;
|
||||
if (direction != 0)
|
||||
actionindex = (1 - d * direction) / 2;
|
||||
walker.Properties (nPsi, Action[actionindex]) = val;
|
||||
}
|
||||
|
||||
void saveTransProb(Walker_t& walker, IndexType d, RealType val,IndexType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex =(1-d*direction)/2;
|
||||
walker.Properties(nPsi,TransProb[transindex])=val;
|
||||
}
|
||||
|
||||
void saveTransProb(ParticleSet& W, IndexType d, RealType val,IndexType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex =(1-d*direction)/2;
|
||||
W.Properties(nPsi,TransProb[transindex])=val;
|
||||
}
|
||||
RealType getTransProb(Walker_t& walker, IndexType d, RealType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex =(1-d*direction)/2;
|
||||
return walker.Properties(nPsi,TransProb[transindex]);
|
||||
}
|
||||
RealType getTransProb(ParticleSet& W, IndexType d, RealType nPsi=0){
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex =(1-d*direction)/2;
|
||||
return W.Properties(nPsi,TransProb[transindex]);
|
||||
}
|
||||
|
||||
inline void printState()
|
||||
RealType getDirectionalAction (Walker_t & walker, IndexType d,
|
||||
IndexType nPsi = 0)
|
||||
{
|
||||
app_log()<<"********PRINT REPTILE STATE*********\n";
|
||||
app_log()<<"Direction="<<direction<<" Headindex="<<headindex<<" tail="<<getBeadIndex(nbeads-1)<<"\n next="<<getBeadIndex(nbeads-2)<<" nbeads="<<nbeads<<endl;
|
||||
app_log()<<"BeadIndex\tWrapIndex\tEnergy\tAction[0]\tAction[1]\tAction[2]\t\n";
|
||||
for( int i=0; i<nbeads; i++)
|
||||
{
|
||||
app_log()<<i<<"\t"<<getBeadIndex(i)<<"\t"<<getBead(i).Properties(LOCALENERGY)<<"\t"<<getBead(i).Properties(Action[0])<<"\t"<<getBead(i).Properties(Action[1])<<"\t"<<getBead(i).Properties(Action[2])<<"\n";
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType actionindex = 2;
|
||||
if (d != 0)
|
||||
actionindex = (1 - direction * d) / 2;
|
||||
|
||||
return walker.Properties (nPsi, Action[actionindex]);
|
||||
}
|
||||
app_log()<<"POSITIONS===============:\n";
|
||||
for( int i=0; i<nbeads; i++)
|
||||
|
||||
RealType getLinkAction (Walker_t & new_walker, Walker_t & old_walker,
|
||||
IndexType d, IndexType nPsi = 0)
|
||||
{
|
||||
RealType af = getDirectionalAction (old_walker, +1, nPsi);
|
||||
RealType ab = getDirectionalAction (new_walker, -1, nPsi);
|
||||
RealType a0 =
|
||||
getDirectionalAction (old_walker, 0,
|
||||
nPsi) + getDirectionalAction (new_walker, 0,
|
||||
nPsi);
|
||||
return af + ab + a0;
|
||||
|
||||
}
|
||||
|
||||
void saveTransProb (Walker_t & walker, IndexType d, RealType val,
|
||||
IndexType nPsi = 0)
|
||||
{
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex = (1 - d * direction) / 2;
|
||||
walker.Properties (nPsi, TransProb[transindex]) = val;
|
||||
}
|
||||
|
||||
void saveTransProb (ParticleSet & W, IndexType d, RealType val,
|
||||
IndexType nPsi = 0)
|
||||
{
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex = (1 - d * direction) / 2;
|
||||
W.Properties (nPsi, TransProb[transindex]) = val;
|
||||
}
|
||||
RealType getTransProb (Walker_t & walker, IndexType d, RealType nPsi = 0)
|
||||
{
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex = (1 - d * direction) / 2;
|
||||
return walker.Properties (nPsi, TransProb[transindex]);
|
||||
}
|
||||
RealType getTransProb (ParticleSet & W, IndexType d, RealType nPsi = 0)
|
||||
{
|
||||
//IndexType repdirection=circbuffer.get_direction();
|
||||
IndexType transindex = (1 - d * direction) / 2;
|
||||
return W.Properties (nPsi, TransProb[transindex]);
|
||||
}
|
||||
|
||||
inline void printState ()
|
||||
{
|
||||
app_log () << "********PRINT REPTILE STATE*********\n";
|
||||
app_log () << "Direction=" << direction << " Headindex=" << headindex
|
||||
<< " tail=" << getBeadIndex (nbeads -
|
||||
1) << "\n next=" <<
|
||||
getBeadIndex (nbeads - 2) << " nbeads=" << nbeads << endl;
|
||||
app_log () <<
|
||||
"BeadIndex\tWrapIndex\tEnergy\tAction[0]\tAction[1]\tAction[2]\t\n";
|
||||
for (int i = 0; i < nbeads; i++)
|
||||
{
|
||||
app_log () << i << "\t" << getBeadIndex (i) << "\t" << getBead (i).
|
||||
Properties (LOCALENERGY) << "\t" << getBead (i).
|
||||
Properties (Action[0]) << "\t" << getBead (i).
|
||||
Properties (Action[1]) << "\t" << getBead (i).
|
||||
Properties (Action[2]) << "\n";
|
||||
}
|
||||
app_log () << "POSITIONS===============:\n";
|
||||
for (int i = 0; i < nbeads; i++)
|
||||
{
|
||||
// app_log()<<i<<"\t1"<<1<<"\t"<<getBead(i).R[0]<<"\n";
|
||||
// app_log()<<i<<"\t2"<<2<<"\t"<<getBead(i).R[1]<<"\n";
|
||||
app_log()<<"BEAD #"<<i<<" tau = "<<tau*i<<endl;
|
||||
app_log()<<getBead(i).R<<endl;
|
||||
app_log () << "BEAD #" << i << " tau = " << tau * i << endl;
|
||||
app_log () << getBead (i).R << endl;
|
||||
}
|
||||
app_log()<<"GVECS===============:\n";
|
||||
for( int i=0; i<nbeads; i++)
|
||||
app_log () << "GVECS===============:\n";
|
||||
for (int i = 0; i < nbeads; i++)
|
||||
{
|
||||
// app_log()<<i<<"\t1"<<1<<"\t"<<getBead(i).G[0]<<"\n";
|
||||
// app_log()<<i<<"\t2"<<2<<"\t"<<getBead(i).G[1]<<"\n";
|
||||
app_log()<<"BEAD #"<<i<<" tau = "<<tau*i<<endl;
|
||||
app_log()<<getBead(i).G<<endl;
|
||||
app_log () << "BEAD #" << i << " tau = " << tau * i << endl;
|
||||
app_log () << getBead (i).G << endl;
|
||||
}
|
||||
app_log()<<"************************************\n";
|
||||
app_log () << "************************************\n";
|
||||
}
|
||||
inline RealType getTau ()
|
||||
{
|
||||
return tau;
|
||||
}
|
||||
inline void setTau (RealType t)
|
||||
{
|
||||
tau = t;
|
||||
}
|
||||
inline RealType getTau(){return tau;}
|
||||
inline void setTau(RealType t){tau=t;}
|
||||
|
||||
|
||||
//This takes a value of imaginary time "t" and returns a 3N particle position vector, corresponding to a time slice extrapolated
|
||||
// from the current reptile. If t>length of reptile, then return the last bead. if t<0; return the first bead.
|
||||
inline Walker_t::ParticlePos_t linearInterp(RealType t)
|
||||
inline Walker_t::ParticlePos_t linearInterp (RealType t)
|
||||
{
|
||||
IndexType nbead=IndexType(t/tau); //Calculate the lower bound on the timeslice. t is between binnum*Tau and (binnum+1)Tau
|
||||
RealType beadfrac=t/tau - nbead; //the fractional coordinate between n and n+1 bead
|
||||
if(nbead<=0)
|
||||
IndexType nbead = IndexType (t / tau); //Calculate the lower bound on the timeslice. t is between binnum*Tau and (binnum+1)Tau
|
||||
RealType beadfrac = t / tau - nbead; //the fractional coordinate between n and n+1 bead
|
||||
if (nbead <= 0)
|
||||
{
|
||||
ParticleSet::ParticlePos_t result=getHead().R;
|
||||
ParticleSet::ParticlePos_t result = getHead ().R;
|
||||
return result;
|
||||
}
|
||||
else if(nbead>=nbeads-1)
|
||||
else if (nbead >= nbeads - 1)
|
||||
{
|
||||
ParticleSet::ParticlePos_t result=getTail().R;
|
||||
ParticleSet::ParticlePos_t result = getTail ().R;
|
||||
return result;
|
||||
}
|
||||
|
||||
else
|
||||
{
|
||||
Walker_t::ParticlePos_t dR(getBead(nbead+1).R), interpR(getBead(nbead).R);
|
||||
dR=dR-getBead(nbead).R;
|
||||
Walker_t::ParticlePos_t dR (getBead (nbead + 1).R),
|
||||
interpR (getBead (nbead).R);
|
||||
dR = dR - getBead (nbead).R;
|
||||
|
||||
interpR = getBead(nbead).R + beadfrac*dR;
|
||||
interpR = getBead (nbead).R + beadfrac * dR;
|
||||
return interpR;
|
||||
}
|
||||
|
||||
}
|
||||
inline ReptileConfig_t getReptileSlicePositions(RealType tau, RealType beta)
|
||||
inline ReptileConfig_t getReptileSlicePositions (RealType tau,
|
||||
RealType beta)
|
||||
{
|
||||
IndexType nbeads_new=IndexType(beta/tau);
|
||||
ReptileConfig_t new_reptile_coords(0);
|
||||
IndexType nbeads_new = IndexType (beta / tau);
|
||||
ReptileConfig_t new_reptile_coords (0);
|
||||
|
||||
for (IndexType i=0; i<nbeads_new; i++) new_reptile_coords.push_back(linearInterp(tau*i));
|
||||
for (IndexType i = 0; i < nbeads_new; i++)
|
||||
new_reptile_coords.push_back (linearInterp (tau * i));
|
||||
|
||||
return new_reptile_coords;
|
||||
}
|
||||
|
||||
inline void setReptileSlicePositions(ReptileConfig_t& rept)
|
||||
inline void setReptileSlicePositions (ReptileConfig_t & rept)
|
||||
{
|
||||
if (rept.size() == nbeads)
|
||||
if (rept.size () == nbeads)
|
||||
{
|
||||
for (int i=0; i<nbeads; i++) getBead(i).R=rept[i];
|
||||
for (int i = 0; i < nbeads; i++)
|
||||
getBead (i).R = rept[i];
|
||||
}
|
||||
else;
|
||||
}
|
||||
|
||||
inline void setReptileSlicePositions(Walker_t::ParticlePos_t R)
|
||||
inline void setReptileSlicePositions (Walker_t::ParticlePos_t R)
|
||||
{
|
||||
for (int i=0; i<nbeads; i++) getBead(i).R=R;
|
||||
for (int i = 0; i < nbeads; i++)
|
||||
getBead (i).R = R;
|
||||
}
|
||||
|
||||
};
|
||||
};
|
||||
|
||||
|
||||
}
|
||||
|
|
|
|||
|
|
@ -14,89 +14,96 @@ namespace qmcplusplus
|
|||
{
|
||||
|
||||
/// Constructor.
|
||||
RMCSingleOMP::RMCSingleOMP(MCWalkerConfiguration& w, TrialWaveFunction& psi, QMCHamiltonian& h,
|
||||
HamiltonianPool& hpool, WaveFunctionPool& ppool):
|
||||
QMCDriver(w,psi,h,ppool), CloneManager(hpool),prestepsVMC(-1), rescaleDrift("no"), beta(-1), beads(-1),fromScratch(true)
|
||||
{
|
||||
RMCSingleOMP::RMCSingleOMP (MCWalkerConfiguration & w,
|
||||
TrialWaveFunction & psi, QMCHamiltonian & h,
|
||||
HamiltonianPool & hpool,
|
||||
WaveFunctionPool & ppool):QMCDriver (w, psi, h,
|
||||
ppool),
|
||||
CloneManager (hpool), prestepsVMC (-1), rescaleDrift ("no"), beta (-1),
|
||||
beads (-1), fromScratch (true)
|
||||
{
|
||||
RootName = "rmc";
|
||||
QMCType ="RMCSingleOMP";
|
||||
QMCDriverMode.set(QMC_UPDATE_MODE,1);
|
||||
QMCDriverMode.set(QMC_WARMUP,0);
|
||||
m_param.add(rescaleDrift,"drift","string");
|
||||
m_param.add(beta,"beta","double");
|
||||
m_param.add(beads,"beads","int");
|
||||
m_param.add(resizeReptile,"resize","int");
|
||||
m_param.add(prestepsVMC,"vmcpresteps","int");
|
||||
QMCType = "RMCSingleOMP";
|
||||
QMCDriverMode.set (QMC_UPDATE_MODE, 1);
|
||||
QMCDriverMode.set (QMC_WARMUP, 0);
|
||||
m_param.add (rescaleDrift, "drift", "string");
|
||||
m_param.add (beta, "beta", "double");
|
||||
m_param.add (beads, "beads", "int");
|
||||
m_param.add (resizeReptile, "resize", "int");
|
||||
m_param.add (prestepsVMC, "vmcpresteps", "int");
|
||||
// if (w.reptile_new==0){w.reptile_new = new Reptile_new(&w);};
|
||||
|
||||
// w.reptile_new->loadFromWalkerList(w.WalkerList);
|
||||
// app_log()<<" reptile_new.size()="<<w.reptile_new->size()<<endl;
|
||||
// app_log()<<" WalkerList.size()="<<w.WalkerList.size()<<endl;
|
||||
|
||||
Action.resize(3);
|
||||
Action[0]=w.addProperty("ActionBackward");
|
||||
Action[1]=w.addProperty("ActionForward");
|
||||
Action[2]=w.addProperty("ActionLocal");
|
||||
Action.resize (3);
|
||||
Action[0] = w.addProperty ("ActionBackward");
|
||||
Action[1] = w.addProperty ("ActionForward");
|
||||
Action[2] = w.addProperty ("ActionLocal");
|
||||
// app_log()<<Action[0]<<" "<<Action[1]<<" "<<Action[2]<<endl;
|
||||
TransProb.resize(2);
|
||||
TransProb[0]=w.addProperty("TransProbBackward");
|
||||
TransProb[1]=w.addProperty("TransProbForward");
|
||||
TransProb.resize (2);
|
||||
TransProb[0] = w.addProperty ("TransProbBackward");
|
||||
TransProb[1] = w.addProperty ("TransProbForward");
|
||||
// app_log()<<TransProb[0]<<" "<<TransProb[1]<<endl;
|
||||
}
|
||||
}
|
||||
|
||||
bool RMCSingleOMP::run()
|
||||
{
|
||||
resetRun();
|
||||
bool RMCSingleOMP::run ()
|
||||
{
|
||||
resetRun ();
|
||||
//start the main estimator
|
||||
// app_log()<<"Starting Estimators\n";
|
||||
Estimators->start(nBlocks);
|
||||
Estimators->start (nBlocks);
|
||||
// app_log()<<"Starting Movers\n";
|
||||
for (int ip=0; ip<NumThreads; ++ip)
|
||||
Movers[ip]->startRun(nBlocks,false);
|
||||
Traces->startRun(nBlocks,traceClones);
|
||||
const bool has_collectables=W.Collectables.size();
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
Movers[ip]->startRun (nBlocks, false);
|
||||
Traces->startRun (nBlocks, traceClones);
|
||||
const bool has_collectables = W.Collectables.size ();
|
||||
|
||||
for (int block=0; block<nBlocks; ++block)
|
||||
for (int block = 0; block < nBlocks; ++block)
|
||||
{
|
||||
// app_log()<<"Block "<<block<<endl;
|
||||
#pragma omp parallel
|
||||
#pragma omp parallel
|
||||
{
|
||||
int ip=omp_get_thread_num();
|
||||
int ip = omp_get_thread_num ();
|
||||
//IndexType updatePeriod=(QMCDriverMode[QMC_UPDATE_MODE])?Period4CheckProperties:(nBlocks+1)*nSteps;
|
||||
IndexType updatePeriod=(QMCDriverMode[QMC_UPDATE_MODE])?Period4CheckProperties:0;
|
||||
IndexType updatePeriod =
|
||||
(QMCDriverMode[QMC_UPDATE_MODE]) ? Period4CheckProperties : 0;
|
||||
//assign the iterators and resuse them
|
||||
MCWalkerConfiguration::iterator wit(W.begin()+wPerNode[ip]), wit_end(W.begin()+wPerNode[ip+1]);
|
||||
MCWalkerConfiguration::iterator wit (W.begin () + wPerNode[ip]),
|
||||
wit_end (W.begin () + wPerNode[ip + 1]);
|
||||
// MCWalkerConfiguration::iterator wit(first);
|
||||
Movers[ip]->startBlock(nSteps);
|
||||
int now_loc=CurrentStep;
|
||||
Movers[ip]->startBlock (nSteps);
|
||||
int now_loc = CurrentStep;
|
||||
|
||||
RealType cnorm=1.0/static_cast<RealType>(wPerNode[ip+1]-wPerNode[ip]);
|
||||
RealType cnorm =
|
||||
1.0 / static_cast < RealType > (wPerNode[ip + 1] - wPerNode[ip]);
|
||||
|
||||
|
||||
|
||||
for (int step=0; step<nSteps; ++step)
|
||||
for (int step = 0; step < nSteps; ++step)
|
||||
{
|
||||
//collectables are reset, it is accumulated while advancing walkers
|
||||
wClones[ip]->resetCollectables();
|
||||
wClones[ip]->resetCollectables ();
|
||||
// wit=first;
|
||||
Movers[ip]->advanceWalkers(wit,wit_end,false);
|
||||
if(has_collectables)
|
||||
Movers[ip]->advanceWalkers (wit, wit_end, false);
|
||||
if (has_collectables)
|
||||
wClones[ip]->Collectables *= cnorm;
|
||||
// wit=first;
|
||||
Movers[ip]->accumulate(wit,wit_end);
|
||||
Movers[ip]->accumulate (wit, wit_end);
|
||||
|
||||
++now_loc;
|
||||
//if (updatePeriod&& now_loc%updatePeriod==0) Movers[ip]->updateWalkers(wit,wit_end);
|
||||
if (Period4WalkerDump&& now_loc%myPeriod4WalkerDump==0)
|
||||
wClones[ip]->saveEnsemble(wit,wit_end);
|
||||
if (Period4WalkerDump && now_loc % myPeriod4WalkerDump == 0)
|
||||
wClones[ip]->saveEnsemble (wit, wit_end);
|
||||
|
||||
branchEngine->collect(CurrentStep, W, branchClones); //Ray Clay: For now, collects and syncs based on first reptile. Need a better way to do this.
|
||||
branchEngine->collect (CurrentStep, W, branchClones); //Ray Clay: For now, collects and syncs based on first reptile. Need a better way to do this.
|
||||
}
|
||||
// wClones[ip]->reptile->calcTauScaling();
|
||||
// wClones[ip]->reptile->calcERun();
|
||||
//branchEngine->collect(CurrentStep, W, branchClones);
|
||||
// wClones[ip]->reptile->resetR2Avg();
|
||||
Movers[ip]->stopBlock(false);
|
||||
Movers[ip]->stopBlock (false);
|
||||
// if (block > 2*wClones[ip]->reptile->nbeads){
|
||||
// wClones[ip]->reptile->printState();
|
||||
// ((RMCUpdateAllWithDrift*)Movers[ip])->checkReptile(wit,wit_end);
|
||||
|
|
@ -104,65 +111,74 @@ bool RMCSingleOMP::run()
|
|||
// if (block>1) branchEngine->updateReptileStats(*wClones[ip]);
|
||||
|
||||
//}
|
||||
}//end-of-parallel for
|
||||
} //end-of-parallel for
|
||||
//Estimators->accumulateCollectables(wClones,nSteps);
|
||||
CurrentStep+=nSteps;
|
||||
CurrentStep += nSteps;
|
||||
// branchEngine->collect(CurrentStep, W, branchClones);
|
||||
Estimators->stopBlock(estimatorClones);
|
||||
Estimators->stopBlock (estimatorClones);
|
||||
//why was this commented out? Are checkpoints stored some other way?
|
||||
if(storeConfigs)
|
||||
recordBlock(block);
|
||||
}//block
|
||||
Estimators->stop(estimatorClones);
|
||||
if (storeConfigs)
|
||||
recordBlock (block);
|
||||
} //block
|
||||
Estimators->stop (estimatorClones);
|
||||
//copy back the random states
|
||||
for (int ip=0; ip<NumThreads; ++ip)
|
||||
*(RandomNumberControl::Children[ip])=*(Rng[ip]);
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
*(RandomNumberControl::Children[ip]) = *(Rng[ip]);
|
||||
//return nbeads and stuff to its orginal unset state;
|
||||
resetVars();
|
||||
app_log()<<"Reptile State of Thread #0\n";
|
||||
W.ReptileList[0]->printState();
|
||||
resetVars ();
|
||||
app_log () << "Reptile State of Thread #0\n";
|
||||
W.ReptileList[0]->printState ();
|
||||
//finalize a qmc section
|
||||
return finalize(nBlocks);
|
||||
}
|
||||
return finalize (nBlocks);
|
||||
}
|
||||
|
||||
void RMCSingleOMP::resetRun()
|
||||
{
|
||||
m_param.put(qmcNode);
|
||||
//For now, assume that nReptiles=NumThreads;
|
||||
nReptiles=NumThreads;
|
||||
|
||||
if(beads<1)
|
||||
beads=beta/Tau;
|
||||
else
|
||||
beta=beads*Tau;
|
||||
|
||||
app_log()<<"Projection time: "<<beta<<" Ha^-1"<<endl;
|
||||
//Calculate the number of VMC presteps if not given:
|
||||
if (prestepsVMC==-1 && fromScratch==true) prestepsVMC=beads+2;
|
||||
//Check to see if the MCWalkerConfiguration is in a state suitable for reptation
|
||||
if (!W.ReptileList.empty())
|
||||
void RMCSingleOMP::resetRun ()
|
||||
{
|
||||
fromScratch=false;
|
||||
m_param.put (qmcNode);
|
||||
//For now, assume that nReptiles=NumThreads;
|
||||
nReptiles = NumThreads;
|
||||
|
||||
app_log()<<"Previous RMC reptiles detected...\n";
|
||||
if(Tau==W.ReptileList[0]->getTau() && beads==W.ReptileList[0]->size()) app_log()<<" Using current reptiles\n"; //do nothing
|
||||
if (beads < 1)
|
||||
beads = beta / Tau;
|
||||
else
|
||||
beta = beads * Tau;
|
||||
|
||||
app_log () << "Projection time: " << beta << " Ha^-1" << endl;
|
||||
//Calculate the number of VMC presteps if not given:
|
||||
if (prestepsVMC == -1 && fromScratch == true)
|
||||
prestepsVMC = beads + 2;
|
||||
//Check to see if the MCWalkerConfiguration is in a state suitable for reptation
|
||||
if (!W.ReptileList.empty ())
|
||||
{
|
||||
fromScratch = false;
|
||||
|
||||
app_log () << "Previous RMC reptiles detected...\n";
|
||||
if (Tau == W.ReptileList[0]->getTau ()
|
||||
&& beads == W.ReptileList[0]->size ())
|
||||
app_log () << " Using current reptiles\n"; //do nothing
|
||||
else //we need to extrapolate off of the current reptile set.
|
||||
{
|
||||
//pull the reptile configurations out
|
||||
app_log()<<" Previous Tau/Beta: "<<W.ReptileList[0]->getTau()<<"/"<<W.ReptileList[0]->getTau()*W.ReptileList[0]->size()<<endl;
|
||||
app_log()<<" New Tau/Beta: "<<Tau<<"/"<<beta<<endl;
|
||||
app_log()<<" Linear interpolation to get new reptile.\n";
|
||||
vector<ReptileConfig_t> repSamps(0);
|
||||
for (IndexType sampid=0; sampid<W.ReptileList.size() && sampid<nReptiles; sampid++)
|
||||
repSamps.push_back(W.ReptileList[sampid]->getReptileSlicePositions(Tau, beta));
|
||||
app_log () << " Previous Tau/Beta: " << W.ReptileList[0]->
|
||||
getTau () << "/" << W.ReptileList[0]->getTau () *
|
||||
W.ReptileList[0]->size () << endl;
|
||||
app_log () << " New Tau/Beta: " << Tau << "/" << beta <<
|
||||
endl;
|
||||
app_log () << " Linear interpolation to get new reptile.\n";
|
||||
vector < ReptileConfig_t > repSamps (0);
|
||||
for (IndexType sampid = 0;
|
||||
sampid < W.ReptileList.size () && sampid < nReptiles;
|
||||
sampid++)
|
||||
repSamps.push_back (W.ReptileList[sampid]->
|
||||
getReptileSlicePositions (Tau, beta));
|
||||
|
||||
//In the event of a disparity in the number of requested reptiles and the ones received.... just copy
|
||||
//Copies cyclically. First iteration copies the first entry, second the second, and so on. So we don't replicate just one config.
|
||||
for (IndexType copyid=0; repSamps.size()<nReptiles; copyid++)
|
||||
repSamps.push_back(repSamps[copyid]);
|
||||
for (IndexType copyid = 0; repSamps.size () < nReptiles; copyid++)
|
||||
repSamps.push_back (repSamps[copyid]);
|
||||
|
||||
|
||||
resetReptiles(repSamps, Tau);
|
||||
resetReptiles (repSamps, Tau);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
@ -170,19 +186,20 @@ void RMCSingleOMP::resetRun()
|
|||
else
|
||||
{
|
||||
//Initialize on whatever walkers are in MCWalkerConfiguration.
|
||||
app_log()<<"Using walkers from previous non-RMC run.\n";
|
||||
vector<ParticlePos_t> wSamps(0);
|
||||
MCWalkerConfiguration::iterator wit(W.begin()), wend(W.end());
|
||||
for (IndexType sampid=0; wit!=wend && sampid<nReptiles; wit++)
|
||||
wSamps.push_back((**wit).R);
|
||||
app_log () << "Using walkers from previous non-RMC run.\n";
|
||||
vector < ParticlePos_t > wSamps (0);
|
||||
MCWalkerConfiguration::iterator wit (W.begin ()), wend (W.end ());
|
||||
for (IndexType sampid = 0; wit != wend && sampid < nReptiles; wit++)
|
||||
wSamps.push_back ((**wit).R);
|
||||
|
||||
for (IndexType copyid=0; wSamps.size()<nReptiles; copyid++)
|
||||
wSamps.push_back(wSamps[copyid]);
|
||||
resetReptiles(wSamps,beads, Tau);
|
||||
for (IndexType copyid = 0; wSamps.size () < nReptiles; copyid++)
|
||||
wSamps.push_back (wSamps[copyid]);
|
||||
resetReptiles (wSamps, beads, Tau);
|
||||
}
|
||||
|
||||
//Now that we know if we're starting from scratch... decide whether to force VMC warmup.
|
||||
if (prestepsVMC==-1 && fromScratch==true) prestepsVMC=beads+2;
|
||||
if (prestepsVMC == -1 && fromScratch == true)
|
||||
prestepsVMC = beads + 2;
|
||||
/* if(resizeReptile<2)
|
||||
{
|
||||
// compute number of beads you need. either by looking at the beads parameter, or by looking at beta and tau.
|
||||
|
|
@ -206,46 +223,58 @@ void RMCSingleOMP::resetRun()
|
|||
///Norm
|
||||
app_log()<<"resizing the reptile not yet implemented."<<endl;
|
||||
}*/
|
||||
makeClones(W,Psi,H);
|
||||
myPeriod4WalkerDump=(Period4WalkerDump>0)?Period4WalkerDump:(nBlocks+1)*nSteps;
|
||||
makeClones (W, Psi, H);
|
||||
myPeriod4WalkerDump =
|
||||
(Period4WalkerDump > 0) ? Period4WalkerDump : (nBlocks + 1) * nSteps;
|
||||
|
||||
if (Movers.empty())
|
||||
if (Movers.empty ())
|
||||
{
|
||||
Movers.resize(NumThreads,0);
|
||||
branchClones.resize(NumThreads,0);
|
||||
estimatorClones.resize(NumThreads,0);
|
||||
traceClones.resize(NumThreads,0);
|
||||
Rng.resize(NumThreads,0);
|
||||
Movers.resize (NumThreads, 0);
|
||||
branchClones.resize (NumThreads, 0);
|
||||
estimatorClones.resize (NumThreads, 0);
|
||||
traceClones.resize (NumThreads, 0);
|
||||
Rng.resize (NumThreads, 0);
|
||||
// W.loadEnsemble(wClones);
|
||||
branchEngine->initReptile(W);
|
||||
branchEngine->initReptile (W);
|
||||
#if !defined(BGP_BUG)
|
||||
#pragma omp parallel for
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for(int ip=0; ip<NumThreads; ++ip)
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
{
|
||||
ostringstream os;
|
||||
estimatorClones[ip]= new EstimatorManager(*Estimators);//,*hClones[ip]);
|
||||
estimatorClones[ip]->resetTargetParticleSet(*wClones[ip]);
|
||||
estimatorClones[ip]->setCollectionMode(false);
|
||||
Rng[ip]=new RandomGenerator_t(*(RandomNumberControl::Children[ip]));
|
||||
traceClones[ip] = Traces->makeClone();
|
||||
hClones[ip]->setRandomGenerator(Rng[ip]);
|
||||
branchClones[ip] = new BranchEngineType(*branchEngine);
|
||||
estimatorClones[ip] = new EstimatorManager (*Estimators); //,*hClones[ip]);
|
||||
estimatorClones[ip]->resetTargetParticleSet (*wClones[ip]);
|
||||
estimatorClones[ip]->setCollectionMode (false);
|
||||
Rng[ip] =
|
||||
new RandomGenerator_t (*(RandomNumberControl::Children[ip]));
|
||||
traceClones[ip] = Traces->makeClone ();
|
||||
hClones[ip]->setRandomGenerator (Rng[ip]);
|
||||
branchClones[ip] = new BranchEngineType (*branchEngine);
|
||||
// branchClones[ip]->initReptile(W);
|
||||
if (QMCDriverMode[QMC_UPDATE_MODE])
|
||||
{
|
||||
os <<" PbyP moves with drift, using RMCUpdatePbyPWithDriftFast"<<endl;
|
||||
Movers[ip]=new RMCUpdatePbyPWithDrift(*wClones[ip],*psiClones[ip],*hClones[ip],*Rng[ip],Action,TransProb);
|
||||
os <<
|
||||
" PbyP moves with drift, using RMCUpdatePbyPWithDriftFast"
|
||||
<< endl;
|
||||
Movers[ip] =
|
||||
new RMCUpdatePbyPWithDrift (*wClones[ip], *psiClones[ip],
|
||||
*hClones[ip], *Rng[ip], Action,
|
||||
TransProb);
|
||||
|
||||
}
|
||||
else
|
||||
{
|
||||
os <<" walker moves with drift, using RMCUpdateAllWithDriftFast"<<endl;
|
||||
Movers[ip]=new RMCUpdateAllWithDrift(*wClones[ip],*psiClones[ip],*hClones[ip],*Rng[ip],Action,TransProb);
|
||||
os <<
|
||||
" walker moves with drift, using RMCUpdateAllWithDriftFast"
|
||||
<< endl;
|
||||
Movers[ip] =
|
||||
new RMCUpdateAllWithDrift (*wClones[ip], *psiClones[ip],
|
||||
*hClones[ip], *Rng[ip], Action,
|
||||
TransProb);
|
||||
}
|
||||
Movers[ip]->nSubSteps=nSubSteps;
|
||||
if(ip==0)
|
||||
app_log() << os.str() << endl;
|
||||
Movers[ip]->nSubSteps = nSubSteps;
|
||||
if (ip == 0)
|
||||
app_log () << os.str () << endl;
|
||||
}
|
||||
}
|
||||
else
|
||||
|
|
@ -253,50 +282,55 @@ void RMCSingleOMP::resetRun()
|
|||
#if !defined(BGP_BUG)
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for(int ip=0; ip<NumThreads; ++ip)
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
{
|
||||
traceClones[ip]->transfer_state_from(*Traces);
|
||||
traceClones[ip]->transfer_state_from (*Traces);
|
||||
}
|
||||
}
|
||||
app_log().flush();
|
||||
app_log ().flush ();
|
||||
#if !defined(BGP_BUG)
|
||||
#pragma omp parallel for
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for(int ip=0; ip<NumThreads; ++ip)
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
{
|
||||
Movers[ip]->put(qmcNode);
|
||||
Movers[ip]->resetRun(branchClones[ip],estimatorClones[ip], traceClones[ip]);
|
||||
Movers[ip]->put (qmcNode);
|
||||
Movers[ip]->resetRun (branchClones[ip], estimatorClones[ip],
|
||||
traceClones[ip]);
|
||||
// wClones[ip]->reptile = new Reptile(*wClones[ip], W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1]);
|
||||
wClones[ip]->reptile=W.ReptileList[ip];
|
||||
wClones[ip]->reptile = W.ReptileList[ip];
|
||||
//app_log()<<"Thread # "<<ip<<endl;
|
||||
// printf(" Thread# %d WalkerList.size()=%d \n",ip,wClones[ip]->WalkerList.size());
|
||||
|
||||
// wClones[ip]->reptile->printState();
|
||||
wClones[ip]->activeBead= 0;
|
||||
wClones[ip]->activeBead = 0;
|
||||
wClones[ip]->direction = +1;
|
||||
|
||||
if (QMCDriverMode[QMC_UPDATE_MODE])
|
||||
{
|
||||
app_log()<<ip<<" initWalkers for pbyp...\n";
|
||||
Movers[ip]->initWalkersForPbyP(W.ReptileList[ip]->repstart, W.ReptileList[ip]->repend);
|
||||
app_log () << ip << " initWalkers for pbyp...\n";
|
||||
Movers[ip]->initWalkersForPbyP (W.ReptileList[ip]->repstart,
|
||||
W.ReptileList[ip]->repend);
|
||||
}
|
||||
else
|
||||
{
|
||||
Movers[ip]->initWalkers(W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1]);
|
||||
Movers[ip]->initWalkers (W.begin () + wPerNode[ip],
|
||||
W.begin () + wPerNode[ip + 1]);
|
||||
// wClones[ip]->reptile = new Reptile(*wClones[ip], W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1]);
|
||||
}
|
||||
|
||||
app_log()<<"ResetRun()::Reptile state after resizing and init for thread #0\n";
|
||||
W.ReptileList[0]->printState();
|
||||
app_log () <<
|
||||
"ResetRun()::Reptile state after resizing and init for thread #0\n";
|
||||
W.ReptileList[0]->printState ();
|
||||
//this will "unroll" the reptile according to forced VMC steps (no bounce). See beginning of function for logic of setting prestepVMC.
|
||||
for (IndexType prestep=0; prestep<prestepsVMC; prestep++)
|
||||
for (IndexType prestep = 0; prestep < prestepsVMC; prestep++)
|
||||
{
|
||||
app_log()<<"prestep# "<<prestep<<endl;
|
||||
Movers[ip]->advanceWalkers(W.begin(),W.begin(),true);
|
||||
app_log () << "prestep# " << prestep << endl;
|
||||
Movers[ip]->advanceWalkers (W.begin (), W.begin (), true);
|
||||
}
|
||||
|
||||
//set up initial action and transprob.
|
||||
MCWalkerConfiguration::iterator wit(W.begin()+wPerNode[ip]), wit_end(W.begin()+wPerNode[ip+1]);
|
||||
MCWalkerConfiguration::iterator wit (W.begin () + wPerNode[ip]),
|
||||
wit_end (W.begin () + wPerNode[ip + 1]);
|
||||
|
||||
//IndexType initsteps = wClones[ip]->reptile->nbeads * 2;
|
||||
|
||||
|
|
@ -311,91 +345,102 @@ void RMCSingleOMP::resetRun()
|
|||
// Movers[ip]->updateWalkers(W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1]);
|
||||
}
|
||||
//app_log()<<"Check\n";
|
||||
branchEngine->checkParameters(W);
|
||||
branchEngine->checkParameters (W);
|
||||
|
||||
#if !defined(BGP_BUG)
|
||||
#pragma omp parallel for
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for(int ip=0; ip<NumThreads; ++ip)
|
||||
for (int ip = 0; ip < NumThreads; ++ip)
|
||||
{
|
||||
for (int prestep=0; prestep<nWarmupSteps; ++prestep)
|
||||
for (int prestep = 0; prestep < nWarmupSteps; ++prestep)
|
||||
{
|
||||
// app_log()<<"Advance Walkers\n";
|
||||
Movers[ip]->advanceWalkers(W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1],false);
|
||||
Movers[ip]->advanceWalkers (W.begin () + wPerNode[ip],
|
||||
W.begin () + wPerNode[ip + 1], false);
|
||||
// app_log()<<"Collect\n";
|
||||
branchEngine->collect(CurrentStep, W, branchClones);
|
||||
branchEngine->collect (CurrentStep, W, branchClones);
|
||||
}
|
||||
Movers[ip]->updateWalkers(W.begin()+wPerNode[ip],W.begin()+wPerNode[ip+1]);
|
||||
Movers[ip]->updateWalkers (W.begin () + wPerNode[ip],
|
||||
W.begin () + wPerNode[ip + 1]);
|
||||
|
||||
}
|
||||
|
||||
fromScratch=false;
|
||||
}
|
||||
fromScratch = false;
|
||||
}
|
||||
|
||||
bool
|
||||
RMCSingleOMP::put(xmlNodePtr q)
|
||||
{
|
||||
m_param.put(q);
|
||||
bool RMCSingleOMP::put (xmlNodePtr q)
|
||||
{
|
||||
m_param.put (q);
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
//This will resize the MCWalkerConfiguration and initialize the ReptileList. Does not care for previous runs.
|
||||
void RMCSingleOMP::resetReptiles(int nReptiles_in, int nbeads_in, RealType tau)
|
||||
{
|
||||
for( MCWalkerConfiguration::ReptileList_t::iterator it=W.ReptileList.begin(); it!=W.ReptileList.end(); it++) delete *it;
|
||||
W.ReptileList.clear();
|
||||
void RMCSingleOMP::resetReptiles (int nReptiles_in, int nbeads_in,
|
||||
RealType tau)
|
||||
{
|
||||
for (MCWalkerConfiguration::ReptileList_t::iterator it =
|
||||
W.ReptileList.begin (); it != W.ReptileList.end (); it++)
|
||||
delete *it;
|
||||
W.ReptileList.clear ();
|
||||
// Maybe we should be more vigorous in cleaning the MCWC WalkerList?
|
||||
vector<int> repWalkerSlice;
|
||||
int nwtot= nbeads_in*nReptiles_in;
|
||||
FairDivideLow(nwtot,nReptiles_in,repWalkerSlice);
|
||||
if(W.getActiveWalkers()-nwtot !=0)
|
||||
addWalkers(nwtot-W.getActiveWalkers());
|
||||
vector < int >repWalkerSlice;
|
||||
int nwtot = nbeads_in * nReptiles_in;
|
||||
FairDivideLow (nwtot, nReptiles_in, repWalkerSlice);
|
||||
if (W.getActiveWalkers () - nwtot != 0)
|
||||
addWalkers (nwtot - W.getActiveWalkers ());
|
||||
|
||||
for(int i=0; i<nReptiles_in; i++)
|
||||
for (int i = 0; i < nReptiles_in; i++)
|
||||
{
|
||||
W.ReptileList.push_back(new Reptile(W, W.begin()+repWalkerSlice[i],W.begin()+repWalkerSlice[i+1]));
|
||||
W.ReptileList[i]->setTau(tau);
|
||||
W.ReptileList.
|
||||
push_back (new
|
||||
Reptile (W, W.begin () + repWalkerSlice[i],
|
||||
W.begin () + repWalkerSlice[i + 1]));
|
||||
W.ReptileList[i]->setTau (tau);
|
||||
}
|
||||
|
||||
|
||||
}
|
||||
}
|
||||
//This will resize the MCWalkerConfiguration and initialize Reptile list. It will then reinitialize the MCWC with a list of Reptile coordinates
|
||||
void RMCSingleOMP::resetReptiles(vector<ReptileConfig_t>& reptile_samps, RealType tau)
|
||||
{
|
||||
if(reptile_samps.empty())
|
||||
void RMCSingleOMP::resetReptiles (vector < ReptileConfig_t > &reptile_samps,
|
||||
RealType tau)
|
||||
{
|
||||
APP_ABORT("RMCSingleOMP::resetReptiles(vector< ReptileConfig_t > reptile_samps): No samples!\n");
|
||||
if (reptile_samps.empty ())
|
||||
{
|
||||
APP_ABORT
|
||||
("RMCSingleOMP::resetReptiles(vector< ReptileConfig_t > reptile_samps): No samples!\n");
|
||||
}
|
||||
else
|
||||
{
|
||||
IndexType nReptiles_in=reptile_samps.size();
|
||||
IndexType nBeads_in = reptile_samps[0].size();
|
||||
resetReptiles(nReptiles_in,nBeads_in,tau);
|
||||
IndexType nReptiles_in = reptile_samps.size ();
|
||||
IndexType nBeads_in = reptile_samps[0].size ();
|
||||
resetReptiles (nReptiles_in, nBeads_in, tau);
|
||||
|
||||
for ( IndexType i=0; i<W.ReptileList.size(); i++)
|
||||
for (IndexType i = 0; i < W.ReptileList.size (); i++)
|
||||
{
|
||||
W.ReptileList[i]->setReptileSlicePositions(reptile_samps[i]);
|
||||
W.ReptileList[i]->setReptileSlicePositions (reptile_samps[i]);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
//For # of walker samples, create that many reptiles with nbeads each. Initialize each reptile to have the value of the walker "seed".
|
||||
void RMCSingleOMP::resetReptiles(vector<ParticlePos_t>& walker_samps, int nBeads_in, RealType tau)
|
||||
{
|
||||
if(walker_samps.empty())
|
||||
void RMCSingleOMP::resetReptiles (vector < ParticlePos_t > &walker_samps,
|
||||
int nBeads_in, RealType tau)
|
||||
{
|
||||
APP_ABORT("RMCSingleOMP::resetReptiles(vector< ParticlePos_t > walker_samps): No samples!\n");
|
||||
if (walker_samps.empty ())
|
||||
{
|
||||
APP_ABORT
|
||||
("RMCSingleOMP::resetReptiles(vector< ParticlePos_t > walker_samps): No samples!\n");
|
||||
}
|
||||
else
|
||||
{
|
||||
IndexType nReptiles_in=walker_samps.size();
|
||||
resetReptiles(nReptiles_in,nBeads_in, tau);
|
||||
IndexType nReptiles_in = walker_samps.size ();
|
||||
resetReptiles (nReptiles_in, nBeads_in, tau);
|
||||
|
||||
for ( IndexType i=0; i<W.ReptileList.size(); i++)
|
||||
for (IndexType i = 0; i < W.ReptileList.size (); i++)
|
||||
{
|
||||
W.ReptileList[i]->setReptileSlicePositions(walker_samps[i]);
|
||||
W.ReptileList[i]->setReptileSlicePositions (walker_samps[i]);
|
||||
}
|
||||
}
|
||||
|
||||
}
|
||||
}
|
||||
|
||||
};
|
||||
|
|
|
|||
Loading…
Reference in New Issue