1 ////////////////////////////////////////////////////////////////////////////////////// 2 // This file is distributed under the University of Illinois/NCSA Open Source License. 3 // See LICENSE file in top directory for details. 4 // 5 // Copyright (c) 2016 Jeongnim Kim and QMCPACK developers. 6 // 7 // File developed by: Jeremy McMinnis, jmcminis@gmail.com, University of Illinois at Urbana-Champaign 8 // Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign 9 // Cynthia Gu, zg1@ornl.gov, Oak Ridge National Laboratory 10 // Jaron T. Krogel, krogeljt@ornl.gov, Oak Ridge National Laboratory 11 // Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory 12 // 13 // File created by: Jeongnim Kim, jeongnim.kim@gmail.com, University of Illinois at Urbana-Champaign 14 ////////////////////////////////////////////////////////////////////////////////////// 15 16 17 /** 18 * @file QMCDriver.h 19 * @brief Declaration of QMCDriver 20 */ 21 #ifndef QMCPLUSPLUS_QMCDRIVER_H 22 #define QMCPLUSPLUS_QMCDRIVER_H 23 24 #include "Configuration.h" 25 #include "OhmmsData/ParameterSet.h" 26 #include "Utilities/PooledData.h" 27 #include "Utilities/TimerManager.h" 28 #include "Utilities/ScopedProfiler.h" 29 #include "QMCWaveFunctions/TrialWaveFunction.h" 30 #include "QMCWaveFunctions/WaveFunctionPool.h" 31 #include "QMCHamiltonians/QMCHamiltonian.h" 32 #include "Estimators/EstimatorManagerBase.h" 33 #include "QMCDrivers/DriverTraits.h" 34 #include "QMCDrivers/QMCDriverInterface.h" 35 #include "QMCDrivers/GreenFunctionModifiers/DriftModifierBase.h" 36 #include "QMCDrivers/SimpleFixedNodeBranch.h" 37 #include "QMCDrivers/BranchIO.h" 38 class Communicate; 39 40 namespace qmcplusplus 41 { 42 /** @defgroup QMCDrivers QMC Driver group 43 * QMC drivers that implement QMC algorithms 44 */ 45 46 /** @defgroup WalkerByWalker QMC Drivers using walker-by-walker update 47 * @brief Move all the particles for each walker 48 */ 49 50 /** @defgroup ParticleByParticle QMC Drivers using particle-by-particle update 51 * @brief Move particle by particle 52 */ 53 54 /** @defgroup MultiplePsi QMC Drivers for energy differences 55 * @brief Umbrella sampling over multiple H/Psi 56 * 57 * This class of QMC drivers are suitable to evaluate 58 * the energy differences of multiple H-Psi pairs. 59 */ 60 61 //forward declarations: Do not include headers if not needed 62 class MCWalkerConfiguration; 63 class HDFWalkerOutput; 64 class TraceManager; 65 66 /** @ingroup QMCDrivers 67 * @{ 68 * @brief abstract base class for QMC engines 69 */ 70 class QMCDriver : public QMCDriverInterface, public QMCTraits, public MPIObjectBase 71 { 72 public: 73 /** enumeration coupled with QMCMode */ 74 enum 75 { 76 QMC_UPDATE_MODE, 77 QMC_MULTIPLE, 78 QMC_OPTIMIZE, 79 QMC_WARMUP 80 }; 81 82 typedef MCWalkerConfiguration::Walker_t Walker_t; 83 typedef Walker_t::Buffer_t Buffer_t; 84 /** bits to classify QMCDriver 85 * 86 * - qmc_driver_mode[QMC_UPDATE_MODE]? particle-by-particle: walker-by-walker 87 * - qmc_driver_mode[QMC_MULTIPLE]? multiple H/Psi : single H/Psi 88 * - qmc_driver_mode[QMC_OPTIMIZE]? optimization : vmc/dmc/rmc 89 */ 90 std::bitset<QMC_MODE_MAX> qmc_driver_mode; 91 92 /// whether to allow traces 93 bool allow_traces; 94 /// traces xml 95 xmlNodePtr traces_xml; 96 97 /// Constructor. 98 QMCDriver(MCWalkerConfiguration& w, 99 TrialWaveFunction& psi, 100 QMCHamiltonian& h, 101 Communicate* comm, 102 const std::string& QMC_driver_type, 103 bool enable_profiling = false); 104 105 ///Copy Constructor (disabled). 106 QMCDriver(const QMCDriver&) = delete; 107 ///Copy operator (disabled). 108 QMCDriver& operator=(const QMCDriver&) = delete; 109 110 virtual ~QMCDriver() override; 111 112 ///return current step current()113 inline int current() const { return CurrentStep; } 114 115 /** set the update mode 116 * @param pbyp if true, use particle-by-particle update 117 */ setUpdateMode(bool pbyp)118 inline void setUpdateMode(bool pbyp) override { qmc_driver_mode[QMC_UPDATE_MODE] = pbyp; } 119 120 /** Set the status of the QMCDriver 121 * @param aname the root file name 122 * @param h5name root name of the master hdf5 file containing previous qmcrun 123 * @param append if true, the run is a continuation of the previous qmc 124 * 125 * All output files will be of 126 * the form "aname.s00X.suffix", where "X" is number 127 * of previous QMC runs for the simulation and "suffix" 128 * is the suffix for the output file. 129 */ 130 void setStatus(const std::string& aname, const std::string& h5name, bool append) override; 131 132 /** add QMCHamiltonian/TrialWaveFunction pair for multiple 133 * @param h QMCHamiltonian 134 * @param psi TrialWaveFunction 135 * 136 * *Multiple* drivers use multiple H/Psi pairs to perform correlated sampling 137 * for energy difference evaluations. 138 */ 139 void add_H_and_Psi(QMCHamiltonian* h, TrialWaveFunction* psi) override; 140 141 /** initialize with xmlNode 142 */ 143 void process(xmlNodePtr cur) override; 144 145 /** return a xmlnode with update **/ 146 xmlNodePtr getQMCNode(); 147 148 void putWalkers(std::vector<xmlNodePtr>& wset) override; 149 putTraces(xmlNodePtr txml)150 inline void putTraces(xmlNodePtr txml) override { traces_xml = txml; } 151 requestTraces(bool traces)152 inline void requestTraces(bool traces) override { allow_traces = traces; } 153 getEngineName()154 std::string getEngineName() override { return QMCType; } 155 156 template<class PDT> setValue(const std::string & aname,PDT x)157 void setValue(const std::string& aname, PDT x) 158 { 159 m_param.setValue(aname, x); 160 } 161 162 ///set the BranchEngineType setBranchEngine(std::unique_ptr<BranchEngineType> && be)163 void setBranchEngine(std::unique_ptr<BranchEngineType>&& be) override { branchEngine = std::move(be); } 164 165 ///return BranchEngineType* getBranchEngine()166 std::unique_ptr<BranchEngineType> getBranchEngine() override { return std::move(branchEngine); } 167 addObservable(const std::string & aname)168 int addObservable(const std::string& aname) 169 { 170 if (Estimators) 171 return Estimators->addObservable(aname.c_str()); 172 else 173 return -1; 174 } 175 getObservable(int i)176 RealType getObservable(int i) { return Estimators->getObservable(i); } 177 setTau(RealType i)178 void setTau(RealType i) { Tau = i; } 179 180 ///set global offsets of the walkers 181 void setWalkerOffsets(); 182 183 ///Observables manager 184 EstimatorManagerBase* Estimators; 185 186 ///Traces manager 187 std::unique_ptr<TraceManager> Traces; 188 189 ///return the random generators getRng()190 inline std::vector<RandomGenerator_t*>& getRng() { return Rng; } 191 192 ///return the i-th random generator getRng(int i)193 inline RandomGenerator_t& getRng(int i) override { return (*Rng[i]); } 194 getDriverMode()195 unsigned long getDriverMode() override { return qmc_driver_mode.to_ulong(); } 196 197 protected: 198 ///branch engine 199 std::unique_ptr<BranchEngineType> branchEngine; 200 ///drift modifer 201 DriftModifierBase* DriftModifier; 202 ///randomize it 203 bool ResetRandom; 204 ///flag to append or restart the run 205 bool AppendRun; 206 ///flag to turn off dumping configurations 207 bool DumpConfig; 208 ///true, if it is a real QMC engine 209 bool IsQMCDriver; 210 /** the number of times this QMCDriver is executed 211 * 212 * MyCounter is initialized to zero by the constructor and is incremented 213 * whenever a run is completed by calling finalize(int block) or 214 * using MyCounter++ as in RQMC. 215 */ 216 int MyCounter; 217 ///the number to delay updates by 218 int kDelay; 219 /** period of dumping walker configurations and everything else for restart 220 * 221 * The unit is a block. 222 */ 223 int Period4CheckPoint; 224 /** period of dumping walker positions and IDs for Forward Walking 225 * 226 * The unit is in steps. 227 */ 228 int storeConfigs; 229 230 ///Period to recalculate the walker properties from scratch. 231 int Period4CheckProperties; 232 233 /** period of recording walker configurations 234 * 235 * Default is 0 indicating that only the last configuration will be saved. 236 */ 237 int Period4WalkerDump; 238 239 /** period of recording walker positions and IDs for forward walking afterwards 240 * 241 */ 242 int Period4ConfigDump; 243 244 ///current step 245 IndexType CurrentStep; 246 247 ///maximum number of blocks 248 IndexType nBlocks; 249 250 ///maximum number of steps 251 IndexType nSteps; 252 253 ///number of steps between a step: VMCs do not evaluate energies 254 IndexType nSubSteps; 255 256 ///number of warmup steps 257 IndexType nWarmupSteps; 258 259 ///counter for number of moves accepted 260 IndexType nAccept; 261 262 ///counter for number of moves /rejected 263 IndexType nReject; 264 265 /// the number of blocks between recomptePsi 266 IndexType nBlocksBetweenRecompute; 267 268 ///the number of walkers 269 IndexType nTargetWalkers; 270 ///the number of saved samples 271 IndexType nTargetSamples; 272 ///alternate method of setting QMC run parameters 273 IndexType nStepsBetweenSamples; 274 ///samples per thread 275 RealType nSamplesPerThread; 276 ///target population 277 RealType nTargetPopulation; 278 279 280 ///timestep 281 RealType Tau; 282 283 ///maximum cpu in secs 284 int MaxCPUSecs; 285 286 ///Time-step factor \f$ 1/(2\tau)\f$ 287 RealType m_oneover2tau; 288 ///Time-step factor \f$ \sqrt{\tau}\f$ 289 RealType m_sqrttau; 290 291 ///pointer to qmc node in xml file 292 xmlNodePtr qmcNode; 293 294 ///type of QMC driver 295 const std::string QMCType; 296 ///the root of h5File 297 std::string h5FileRoot; 298 ///root of all the output files 299 std::string RootName; 300 301 ///store any parameter that has to be read from a file 302 ParameterSet m_param; 303 304 ///walker ensemble 305 MCWalkerConfiguration& W; 306 307 ///trial function 308 TrialWaveFunction& Psi; 309 310 ///Hamiltonian 311 QMCHamiltonian& H; 312 313 ///record engine for walkers 314 HDFWalkerOutput* wOut; 315 316 ///a list of TrialWaveFunctions for multiple method 317 std::vector<TrialWaveFunction*> Psi1; 318 319 ///a list of QMCHamiltonians for multiple method 320 std::vector<QMCHamiltonian*> H1; 321 322 ///Random number generators 323 std::vector<RandomGenerator_t*> Rng; 324 325 ///a list of mcwalkerset element 326 std::vector<xmlNodePtr> mcwalkerNodePtr; 327 328 ///temporary storage for drift 329 ParticleSet::ParticlePos_t drift; 330 331 ///temporary storage for random displacement 332 ParticleSet::ParticlePos_t deltaR; 333 334 ///turn on spin moves 335 std::string SpinMoves; 336 RealType SpinMass; 337 338 bool putQMCInfo(xmlNodePtr cur); 339 340 void addWalkers(int nwalkers); 341 342 /** record the state of the block 343 * @param block current block 344 * 345 * virtual function with a default implementation 346 */ 347 virtual void recordBlock(int block) override; 348 349 /** finalize a qmc section 350 * @param block current block 351 * @param dumpwalkers if true, dump walkers 352 * 353 * Accumulate energy and weight is written to a hdf5 file. 354 * Finialize the estimators 355 */ 356 bool finalize(int block, bool dumpwalkers = true); 357 358 int rotation; 359 std::string getRotationName(std::string RootName); 360 std::string getLastRotationName(std::string RootName); get_root_name()361 const std::string& get_root_name() const override { return RootName; } 362 363 private: 364 NewTimer* checkpointTimer; 365 ///time the driver lifetime 366 ScopedTimer driver_scope_timer_; 367 ///profile the driver lifetime 368 ScopedProfiler driver_scope_profiler_; 369 }; 370 /**@}*/ 371 } // namespace qmcplusplus 372 373 #endif 374