Simpatico  v1.10
MdEnergyAnalyzer.cpp
1 /*
2 * Simpatico - Simulation Package for Polymeric and Molecular Liquids
3 *
4 * Copyright 2010 - 2017, The Regents of the University of Minnesota
5 * Distributed under the terms of the GNU General Public License.
6 */
7 
8 #include "MdEnergyAnalyzer.h"
9 #include <mcMd/mdSimulation/MdSimulation.h>
10 #include <mcMd/mdSimulation/MdSystem.h>
11 
12 #ifndef SIMP_NOPAIR
13 #include <mcMd/potentials/pair/MdPairPotential.h>
14 #endif
15 #ifdef SIMP_BOND
16 #include <mcMd/potentials/bond/BondPotential.h>
17 #endif
18 #ifdef SIMP_ANGLE
19 #include <mcMd/potentials/angle/AnglePotential.h>
20 #endif
21 #ifdef SIMP_DIHEDRAL
22 #include <mcMd/potentials/dihedral/DihedralPotential.h>
23 #endif
24 #ifdef SIMP_COULOMB
25 #include <mcMd/potentials/coulomb/MdCoulombPotential.h>
26 #endif
27 #ifdef SIMP_EXTERNAL
28 #include <mcMd/potentials/external/ExternalPotential.h>
29 #endif
30 
31 #include <util/accumulators/Average.h>
32 #include <util/format/Int.h>
33 #include <util/format/Dbl.h>
34 #include <util/mpi/MpiLoader.h>
35 #include <util/misc/ioUtil.h>
36 
37 #include <sstream>
38 
39 namespace McMd
40 {
41 
42  using namespace Util;
43 
44  /*
45  * Constructor.
46  */
48  : SystemAnalyzer<MdSystem>(system),
49  totalAveragePtr_(0),
50  kineticAveragePtr_(0),
51  potentialAveragePtr_(0),
52  #ifndef SIMP_NOPAIR
53  pairAveragePtr_(0),
54  #endif
55  #ifdef SIMP_BOND
56  bondAveragePtr_(0),
57  #endif
58  #ifdef SIMP_ANGLE
59  angleAveragePtr_(0),
60  #endif
61  #ifdef SIMP_DIHEDRAL
62  dihedralAveragePtr_(0),
63  #endif
64  #ifdef SIMP_COULOMB
65  coulombRSpaceAveragePtr_(0),
66  coulombKSpaceAveragePtr_(0),
67  coulombAveragePtr_(0),
68  coulombComponents_(false),
69  #endif
70  #ifdef SIMP_EXTERNAL
71  externalAveragePtr_(0),
72  #endif
73  nSamplePerBlock_(0),
74  isInitialized_(false)
75  { setClassName("MdEnergyAnalyzer"); }
76 
77  /*
78  * Read interval and outputFileName.
79  */
80  void MdEnergyAnalyzer::readParameters(std::istream& in)
81  {
82  MdSystem& sys = system();
83 
84  readInterval(in);
86  nSamplePerBlock_ = 0;
87  readOptional<int>(in, "nSamplePerBlock", nSamplePerBlock_);
88  #ifdef SIMP_COULOMB
89  if (sys.hasCoulombPotential()) {
90  coulombComponents_ = false;
91  readOptional<bool>(in, "coulombComponents", coulombComponents_);
92  }
93  #endif
94 
95  totalAveragePtr_ = new Average;
96  totalAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
97  kineticAveragePtr_ = new Average;
98  kineticAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
99  potentialAveragePtr_ = new Average;
100  potentialAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
101  #ifndef SIMP_NOPAIR
102  pairAveragePtr_ = new Average;
103  pairAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
104  #endif
105  #ifdef SIMP_BOND
106  if (sys.hasBondPotential()) {
107  bondAveragePtr_ = new Average;
108  bondAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
109  }
110  #endif
111  #ifdef SIMP_ANGLE
112  if (sys.hasAnglePotential()) {
113  angleAveragePtr_ = new Average;
114  angleAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
115  }
116  #endif
117  #ifdef SIMP_DIHEDRAL
118  if (sys.hasDihedralPotential()) {
119  dihedralAveragePtr_ = new Average;
120  dihedralAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
121  }
122  #endif
123  #ifdef SIMP_COULOMB
124  if (sys.hasCoulombPotential()) {
125  if (coulombComponents_) {
126  coulombRSpaceAveragePtr_ = new Average;
127  coulombRSpaceAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
128  coulombKSpaceAveragePtr_ = new Average;
129  coulombKSpaceAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
130  }
131  coulombAveragePtr_ = new Average;
132  coulombAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
133  }
134  #endif
135  #ifdef SIMP_EXTERNAL
136  if (sys.hasExternalPotential()) {
137  externalAveragePtr_ = new Average;
138  externalAveragePtr_->setNSamplePerBlock(nSamplePerBlock_);
139  }
140  #endif
141 
142  isInitialized_ = true;
143  }
144 
145  /*
146  * Load internal state from an archive.
147  */
149  {
150  MdSystem& sys = system();
151 
152  // Load interval and outputFileName
154 
155  nSamplePerBlock_ = 0; // default value
156  bool isRequired = false;
157  loadParameter<int>(ar, "nSamplePerBlock", nSamplePerBlock_, isRequired);
158  #ifdef SIMP_COULOMB
159  if (sys.hasCoulombPotential()) {
160  loadParameter<bool>(ar, "coulombComponents", coulombComponents_, isRequired);
161  }
162  #endif
163 
164  // Load Average accumulators
165  totalAveragePtr_ = new Average;
166  ar >> *totalAveragePtr_;
167  kineticAveragePtr_ = new Average;
168  ar >> *kineticAveragePtr_;
169  potentialAveragePtr_ = new Average;
170  ar >> *potentialAveragePtr_;
171  #ifndef SIMP_NOPAIR
172  pairAveragePtr_ = new Average;
173  ar >> *pairAveragePtr_;
174  #endif
175  #ifdef SIMP_BOND
176  if (sys.hasBondPotential()) {
177  bondAveragePtr_ = new Average;
178  ar >> *bondAveragePtr_;
179  }
180  #endif
181  #ifdef SIMP_ANGLE
182  if (sys.hasAnglePotential()) {
183  angleAveragePtr_ = new Average;
184  ar >> *angleAveragePtr_;
185  }
186  #endif
187  #ifdef SIMP_DIHEDRAL
188  if (sys.hasDihedralPotential()) {
189  dihedralAveragePtr_ = new Average;
190  ar >> *dihedralAveragePtr_;
191  }
192  #endif
193  #ifdef SIMP_COULOMB
194  if (sys.hasCoulombPotential()) {
195  if (coulombComponents_) {
196  coulombRSpaceAveragePtr_ = new Average;
197  ar >> *coulombRSpaceAveragePtr_;
198  coulombKSpaceAveragePtr_ = new Average;
199  ar >> *coulombKSpaceAveragePtr_;
200  }
201  coulombAveragePtr_ = new Average;
202  ar >> *coulombAveragePtr_;
203  }
204  #endif
205  #ifdef SIMP_EXTERNAL
206  if (sys.hasExternalPotential()) {
207  externalAveragePtr_ = new Average;
208  ar >> *externalAveragePtr_;
209  }
210  #endif
211 
212  // If nSamplePerBlock != 0, open an output file for block averages.
213  if (nSamplePerBlock_ != 0) {
214  fileMaster().openOutputFile(Analyzer::outputFileName(".dat"), outputFile_);
215  }
216 
217  isInitialized_ = true;
218  }
219 
220  /*
221  * Save internal state to an archive.
222  */
224  {
225  MdSystem& sys = system();
226 
227  Analyzer::save(ar);
228  bool isActive = bool(nSamplePerBlock_);
229  Parameter::saveOptional(ar, nSamplePerBlock_, isActive);
230  #ifdef SIMP_COULOMB
231  if (sys.hasCoulombPotential()) {
232  isActive = coulombComponents_;
233  Parameter::saveOptional(ar, coulombComponents_, isActive);
234  }
235  #endif
236 
237  // Save average accumulators
238  ar << *totalAveragePtr_;
239  ar << *kineticAveragePtr_;
240  ar << *potentialAveragePtr_;
241  #ifndef SIMP_NOPAIR
242  ar << *pairAveragePtr_;
243  #endif
244  #ifdef SIMP_BOND
245  if (sys.hasBondPotential()) {
246  ar << *bondAveragePtr_;
247  }
248  #endif
249  #ifdef SIMP_ANGLE
250  if (sys.hasAnglePotential()) {
251  assert(angleAveragePtr_);
252  ar << *angleAveragePtr_;
253  }
254  #endif
255  #ifdef SIMP_DIHEDRAL
256  if (sys.hasDihedralPotential()) {
257  assert(dihedralAveragePtr_);
258  ar << *dihedralAveragePtr_;
259  }
260  #endif
261  #ifdef SIMP_COULOMB
262  if (sys.hasCoulombPotential()) {
263  if (coulombComponents_) {
264  assert(coulombRSpaceAveragePtr_);
265  ar << *coulombRSpaceAveragePtr_;
266  assert(coulombKSpaceAveragePtr_);
267  ar << *coulombKSpaceAveragePtr_;
268  }
269  assert(coulombAveragePtr_);
270  ar << *coulombAveragePtr_;
271  }
272  #endif
273  #ifdef SIMP_EXTERNAL
274  if (sys.hasExternalPotential()) {
275  assert(externalAveragePtr_);
276  ar << *externalAveragePtr_;
277  }
278  #endif
279  }
280 
281  /*
282  * Reset nSample.
283  */
285  {
286  MdSystem& sys = system();
287 
288  totalAveragePtr_->clear();
289  kineticAveragePtr_->clear();
290  potentialAveragePtr_->clear();
291  #ifndef SIMP_NOPAIR
292  pairAveragePtr_->clear();
293  #endif
294  #ifdef SIMP_BOND
295  if (sys.hasBondPotential()) {
296  bondAveragePtr_->clear();
297  }
298  #endif
299  #ifdef SIMP_ANGLE
300  if (sys.hasAnglePotential()) {
301  assert(angleAveragePtr_);
302  angleAveragePtr_->clear();
303  }
304  #endif
305  #ifdef SIMP_DIHEDRAL
306  if (sys.hasDihedralPotential()) {
307  assert(dihedralAveragePtr_);
308  dihedralAveragePtr_->clear();
309  }
310  #endif
311  #ifdef SIMP_COULOMB
312  if (sys.hasCoulombPotential()) {
313  if (coulombComponents_) {
314  assert(coulombRSpaceAveragePtr_);
315  coulombRSpaceAveragePtr_->clear();
316  assert(coulombKSpaceAveragePtr_);
317  coulombKSpaceAveragePtr_->clear();
318  }
319  assert(coulombAveragePtr_);
320  coulombAveragePtr_->clear();
321  }
322  #endif
323  #ifdef SIMP_EXTERNAL
324  if (sys.hasExternalPotential()) {
325  assert(externalAveragePtr_);
326  externalAveragePtr_->clear();
327  }
328  #endif
329  }
330 
331  /*
332  * Open outputfile
333  */
335  {
336  MdSystem& sys = system();
337  Simulation& sim = sys.simulation();
338 
339  if (outputFile_.is_open()) {
340  outputFile_.close();
341  }
342 
343  if (nSamplePerBlock_ > 0) {
344 
345  // Open *.fmt file with format of *.dat data file
346  std::string filename;
347  filename = outputFileName(".fmt");
348  sim.fileMaster().openOutputFile(filename, outputFile_);
349 
350  outputFile_ << " iStep";
351  outputFile_ << " Kinetic";
352  #ifndef SIMP_NOPAIR
353  outputFile_ << " Pair";
354  #endif
355  #ifdef SIMP_BOND
356  if (sys.hasBondPotential()) {
357  outputFile_ << " Bond";
358  }
359  #endif
360  #ifdef SIMP_ANGLE
361  if (sys.hasAnglePotential()) {
362  outputFile_ << " Angle";
363  }
364  #endif
365  #ifdef SIMP_DIHEDRAL
366  if (sys.hasDihedralPotential()) {
367  outputFile_ << " Dihedral";
368  }
369  #endif
370  #ifdef SIMP_COULOMB
371  if (sys.hasCoulombPotential()) {
372  outputFile_ << " Coulomb";
373  }
374  #endif
375  #ifdef SIMP_EXTERNAL
376  if (sys.hasExternalPotential()) {
377  outputFile_ << " External";
378  }
379  #endif
380  outputFile_ << " Total";
381  outputFile_.close();
382 
383  // Open *.dat data file, leave open for writing
384  filename = outputFileName(".dat");
385  sim.fileMaster().openOutputFile(filename, outputFile_);
386 
387  }
388 
389  }
390 
391  /*
392  * Output energy to file
393  */
394  void MdEnergyAnalyzer::sample(long iStep)
395  {
396  if (isAtInterval(iStep)) {
397  MdSystem& sys = system();
398 
399  //outputFile_ << Int(iStep, 10);
400  double kinetic = sys.kineticEnergy();
401  kineticAveragePtr_->sample(kinetic);
402  //outputFile_ << Dbl(kinetic, 15);
403 
404  double potential = 0.0;
405  #ifndef SIMP_NOPAIR
406  double pair = sys.pairPotential().energy();
407  potential += pair;
408  pairAveragePtr_->sample(pair);
409  //outputFile_ << Dbl(pair, 15);
410  #endif
411  #ifdef SIMP_BOND
412  if (sys.hasBondPotential()) {
413  double bond = sys.bondPotential().energy();
414  potential += bond;
415  bondAveragePtr_->sample(bond);
416  //outputFile_ << Dbl(bond, 15);
417  }
418  #endif
419  #ifdef SIMP_ANGLE
420  if (sys.hasAnglePotential()) {
421  double angle = sys.anglePotential().energy();
422  potential += angle;
423  assert(angleAveragePtr_);
424  angleAveragePtr_->sample(angle);
425  //outputFile_ << Dbl(angle, 15);
426  }
427  #endif
428  #ifdef SIMP_DIHEDRAL
429  if (sys.hasDihedralPotential()) {
430  double dihedral = sys.dihedralPotential().energy();
431  potential += dihedral;
432  assert(dihedralAveragePtr_);
433  dihedralAveragePtr_->sample(dihedral);
434  // outputFile_ << Dbl(dihedral, 15);
435  }
436  #endif
437  #ifdef SIMP_COULOMB
438  if (sys.hasCoulombPotential()) {
439  if (coulombComponents_) {
440  // R-space contribution
441  double coulombRSpace = sys.coulombPotential().rSpaceEnergy();
442  assert(coulombRSpaceAveragePtr_);
443  coulombRSpaceAveragePtr_->sample(coulombRSpace);
444  // K-space contribution
445  double coulombKSpace = sys.coulombPotential().kSpaceEnergy();
446  assert(coulombKSpaceAveragePtr_);
447  coulombKSpaceAveragePtr_->sample(coulombKSpace);
448  }
449  // Total
450  double coulomb = sys.coulombPotential().energy();
451  assert(coulombAveragePtr_);
452  coulombAveragePtr_->sample(coulomb);
453  potential += coulomb;
454  // outputFile_ << Dbl(dihedral, 15);
455  }
456  #endif
457  #ifdef SIMP_EXTERNAL
458  if (sys.hasExternalPotential()) {
459  double external = sys.externalPotential().energy();
460  potential += external;
461  assert(externalAveragePtr_);
462  externalAveragePtr_->sample(external);
463  // outputFile_ << Dbl(external, 15);
464  }
465  #endif
466 
467  assert(potentialAveragePtr_);
468  potentialAveragePtr_->sample(potential);
469  double total = kinetic + potential;
470  totalAveragePtr_->sample(total);
471  // outputFile_ << Dbl(total, 20)
472  // << std::endl;
473 
474  // Output block averages, if needed
475  if (nSamplePerBlock_ > 0 && totalAveragePtr_->isBlockComplete()) {
476  int beginStep = iStep - (nSamplePerBlock_ - 1)*interval();
477  outputFile_ << Int(beginStep, 12);
478  outputFile_ << Dbl(kineticAveragePtr_->blockAverage());
479  #ifndef SIMP_NOPAIR
480  outputFile_ << Dbl(pairAveragePtr_->blockAverage());
481  #endif
482  #ifdef SIMP_BOND
483  if (sys.hasBondPotential()) {
484  outputFile_ << Dbl(bondAveragePtr_->blockAverage());
485  }
486  #endif
487  #ifdef SIMP_ANGLE
488  if (sys.hasAnglePotential()) {
489  assert(angleAveragePtr_);
490  outputFile_ << Dbl(angleAveragePtr_->blockAverage());
491  }
492  #endif
493  #ifdef SIMP_DIHEDRAL
494  if (sys.hasDihedralPotential()) {
495  assert(dihedralAveragePtr_);
496  outputFile_ << Dbl(dihedralAveragePtr_->blockAverage());
497  }
498  #endif
499  #ifdef SIMP_COULOMB
500  if (sys.hasCoulombPotential()) {
501  assert(coulombAveragePtr_);
502  outputFile_ << Dbl(coulombAveragePtr_->blockAverage());
503  }
504  #endif
505  #ifdef SIMP_EXTERNAL
506  if (sys.hasExternalPotential()) {
507  assert(externalAveragePtr_);
508  outputFile_ << Dbl(externalAveragePtr_->blockAverage());
509  }
510  #endif
511  outputFile_ << Dbl(totalAveragePtr_->blockAverage());
512  outputFile_ << "\n";
513  }
514 
515  }
516  }
517 
518  /*
519  * Output results to file after simulation is completed.
520  */
522  {
523  MdSystem& sys = system();
524  Simulation& sim = sys.simulation();
525 
526  // Close data (*.dat) file, if any
527  if (outputFile_.is_open()) {
528  outputFile_.close();
529  }
530 
531  // Write parameter (*.prm) file
532  sim.fileMaster().openOutputFile(outputFileName(".prm"), outputFile_);
533  ParamComposite::writeParam(outputFile_);
534  outputFile_.close();
535 
536  // Write average (*.ave) file
537  sim.fileMaster().openOutputFile(outputFileName(".ave"), outputFile_);
538  double ave, err;
539  ave = kineticAveragePtr_->average();
540  err = kineticAveragePtr_->blockingError();
541  outputFile_ << "Kinetic " << Dbl(ave)
542  << " +- " << Dbl(err, 9, 2) << "\n";
543  #ifndef SIMP_NOPAIR
544  ave = pairAveragePtr_->average();
545  err = pairAveragePtr_->blockingError();
546  outputFile_ << "Pair " << Dbl(ave)
547  << " +- " << Dbl(err, 9, 2) << "\n";
548  #endif
549  #ifdef SIMP_BOND
550  if (sys.hasBondPotential()) {
551  ave = bondAveragePtr_->average();
552  err = bondAveragePtr_->blockingError();
553  outputFile_ << "Bond " << Dbl(ave)
554  << " +- " << Dbl(err, 9, 2) << "\n";
555  }
556  #endif
557  #ifdef SIMP_ANGLE
558  if (sys.hasAnglePotential()) {
559  assert(angleAveragePtr_);
560  ave = angleAveragePtr_->average();
561  err = angleAveragePtr_->blockingError();
562  outputFile_ << "Angle " << Dbl(ave)
563  << " +- " << Dbl(err, 9, 2) << "\n";
564  }
565  #endif
566  #ifdef SIMP_DIHEDRAL
567  if (sys.hasDihedralPotential()) {
568  assert(dihedralAveragePtr_);
569  ave = dihedralAveragePtr_->average();
570  err = dihedralAveragePtr_->blockingError();
571  outputFile_ << "Dihedral " << Dbl(ave)
572  << " +- " << Dbl(err, 9, 2) << "\n";
573  }
574  #endif
575  #ifdef SIMP_COULOMB
576  if (sys.hasCoulombPotential()) {
577 
578  assert(coulombAveragePtr_);
579  ave = coulombAveragePtr_->average();
580  err = coulombAveragePtr_->blockingError();
581  outputFile_ << "Coulomb " << Dbl(ave)
582  << " +- " << Dbl(err, 9, 2) << "\n";
583 
584  if (coulombComponents_) {
585  assert(coulombRSpaceAveragePtr_);
586  ave = coulombRSpaceAveragePtr_->average();
587  err = coulombRSpaceAveragePtr_->blockingError();
588  outputFile_ << "Coulomb(R) " << Dbl(ave)
589  << " +- " << Dbl(err, 9, 2) << "\n";
590  assert(coulombKSpaceAveragePtr_);
591  ave = coulombKSpaceAveragePtr_->average();
592  err = coulombKSpaceAveragePtr_->blockingError();
593  outputFile_ << "Coulomb(K) " << Dbl(ave) << " +- "
594  << Dbl(err, 9, 2) << "\n";
595  }
596  }
597  #endif
598  #ifdef SIMP_EXTERNAL
599  if (sys.hasExternalPotential()) {
600  assert(externalAveragePtr_);
601  ave = externalAveragePtr_->average();
602  err = externalAveragePtr_->blockingError();
603  outputFile_ << "External " << Dbl(ave) << " +- " << Dbl(err, 9, 2) << "\n";
604  }
605  #endif
606 
607  assert(potentialAveragePtr_);
608  ave = potentialAveragePtr_->average();
609  err = potentialAveragePtr_->blockingError();
610  outputFile_ << "Potential " << Dbl(ave) << " +- " << Dbl(err, 9, 2) << "\n";
611 
612  ave = totalAveragePtr_->average();
613  err = totalAveragePtr_->blockingError();
614  outputFile_ << "Total " << Dbl(ave) << " +- " << Dbl(err, 9, 2) << "\n";
615 
616  outputFile_.close();
617 
618  // Write error analysis (*.aer) file
619  sim.fileMaster().openOutputFile(outputFileName(".aer"), outputFile_);
620  outputFile_ <<
621  "---------------------------------------------------------------------------------\n";
622  outputFile_ << "Kinetic:\n\n";
623  kineticAveragePtr_->output(outputFile_);
624  outputFile_ <<
625  "---------------------------------------------------------------------------------\n";
626  outputFile_ << "Pair:\n\n";
627  #ifndef SIMP_NOPAIR
628  pairAveragePtr_->output(outputFile_);
629  #endif
630  #ifdef SIMP_BOND
631  if (sys.hasBondPotential()) {
632  outputFile_ <<
633  "---------------------------------------------------------------------------------\n";
634  outputFile_ << "Bond:\n\n";
635  bondAveragePtr_->output(outputFile_);
636  }
637  #endif
638  #ifdef SIMP_ANGLE
639  if (sys.hasAnglePotential()) {
640  outputFile_ <<
641  "---------------------------------------------------------------------------------\n";
642  outputFile_ << "Angle:\n\n";
643  assert(angleAveragePtr_);
644  angleAveragePtr_->output(outputFile_);
645  }
646  #endif
647  #ifdef SIMP_DIHEDRAL
648  if (sys.hasDihedralPotential()) {
649  outputFile_ <<
650  "---------------------------------------------------------------------------------\n";
651  outputFile_ << "Dihedral:\n\n";
652  assert(dihedralAveragePtr_);
653  dihedralAveragePtr_->output(outputFile_);
654  }
655  #endif
656  #ifdef SIMP_COULOMB
657  if (sys.hasCoulombPotential()) {
658  outputFile_ <<
659  "---------------------------------------------------------------------------------\n";
660  outputFile_ << "Coulomb:\n\n";
661  assert(coulombAveragePtr_);
662  coulombAveragePtr_->output(outputFile_);
663  }
664  #endif
665  #ifdef SIMP_EXTERNAL
666  if (sys.hasExternalPotential()) {
667  outputFile_ <<
668  "---------------------------------------------------------------------------------\n";
669  outputFile_ << "External:\n\n";
670  assert(externalAveragePtr_);
671  externalAveragePtr_->output(outputFile_);
672  }
673  #endif
674 
675  outputFile_ <<
676  "---------------------------------------------------------------------------------\n";
677  outputFile_ << "Potential:\n\n";
678  assert(potentialAveragePtr_);
679  potentialAveragePtr_->output(outputFile_);
680 
681  outputFile_ <<
682  "---------------------------------------------------------------------------------\n";
683  outputFile_ << "Total:\n\n";
684  assert(totalAveragePtr_);
685  totalAveragePtr_->output(outputFile_);
686 
687  outputFile_.close();
688  }
689 
690 
691 }
virtual void save(Serializable::OArchive &ar)
Load parameters to archive.
virtual void readParameters(std::istream &in)
Read interval and output file name.
void clear()
Clear all accumulators, set to empty initial state.
Definition: Average.cpp:42
DihedralPotential & dihedralPotential() const
Return DihedralPotential by reference.
Definition: MdSystem.h:571
bool isBlockComplete() const
Is the current block average complete?
Definition: Average.h:232
virtual double energy(double rsq, int iAtomType, int jAtomType) const =0
Return pair energy for a single pair.
double rSpaceEnergy()
Return short-range r-space part of Coulomb energy.
bool isRequired() const
Is this ParamComposite required in the input file?
virtual void save(Serializable::OArchive &ar)
Save internal state to an output archive.
Calculates the average and variance of a sampled property.
Definition: Average.h:43
double average() const
Return the average of all sampled values.
void openOutputFile(const std::string &filename, std::ofstream &out, std::ios_base::openmode mode=std::ios_base::out) const
Open an output file.
Definition: FileMaster.cpp:290
virtual void loadParameters(Serializable::IArchive &ar)
Load parameters from archive.
virtual void output()
Write file averages and error analysis to file.
MdSystem & system()
Return reference to parent system.
Wrapper for a double precision number, for formatted ostream output.
Definition: Dbl.h:39
double blockingError() const
Return estimated error on average from blocking analysis.
Definition: Average.cpp:133
void readOutputFileName(std::istream &in)
Read outputFileName from file.
bool hasBondPotential() const
Does a bond potential exist?.
Definition: MdSystem.h:531
ExternalPotential & externalPotential() const
Return ExternalPotential by reference.
Definition: MdSystem.h:605
The main object in a simulation, which coordinates others.
bool isActive() const
Is this parameter active?
bool hasExternalPotential() const
Does an external potential exist?.
Definition: MdSystem.h:599
Saving / output archive for binary ostream.
virtual double energy(const Vector &R1, const Vector &R2, const Vector &R3, int type) const =0
Returns potential energy for one dihedral.
MdPairPotential & pairPotential() const
Return MdPairPotential by reference.
Definition: MdSystem.h:520
virtual double energy(const Vector &position, int i) const =0
Returns external potential energy of a single particle.
double kSpaceEnergy()
Get long-range k-space part of Coulomb energy.
Simulation & simulation() const
Get the parent Simulation by reference.
Definition: System.h:1055
virtual void writeParam(std::ostream &out)
Write all parameters to an output stream.
void readInterval(std::istream &in)
Read interval from file, with error checking.
virtual void clear()
Clear nSample counter and all accumulators.
Utility classes for scientific computation.
Definition: accumulators.mod:1
virtual void loadParameters(Serializable::IArchive &ar)
Load internal state from an input archive.
virtual void sample(long iStep)
Write energy to file.
bool hasAnglePotential() const
Does angle potential exist?.
Definition: MdSystem.h:548
AnglePotential & anglePotential() const
Return AnglePotential by reference.
Definition: MdSystem.h:554
Wrapper for an int, for formatted ostream output.
Definition: Int.h:36
bool hasCoulombPotential() const
Does a Coulomb potential exist?.
Definition: MdSystem.h:582
double kineticEnergy() const
Compute and return total kinetic energy.
Definition: MdSystem.cpp:1003
Template for Analyzer associated with one System.
void output(std::ostream &out) const
Output final statistical properties to file.
Definition: Average.cpp:178
virtual double energy(double cosTheta, int type) const =0
Returns potential energy for one angle.
virtual double energy(double rSq, int type) const =0
Returns potential energy for one bond.
void setNSamplePerBlock(int nSamplePerBlock)
Set nSamplePerBlock.
Definition: Average.cpp:63
bool hasDihedralPotential() const
Does a dihedral potential exist?.
Definition: MdSystem.h:565
Saving archive for binary istream.
void sample(double value)
Add a sampled value to the ensemble.
Definition: Average.cpp:94
Single-processor Monte Carlo (MC) and molecular dynamics (MD).
MdEnergyAnalyzer(MdSystem &system)
Constructor.
void setClassName(const char *className)
Set class name string.
virtual void setup()
Setup before main loop.
FileMaster & fileMaster()
Get the FileMaster by reference.
A System for Molecular Dynamics simulation.
Definition: MdSystem.h:68
bool isAtInterval(long counter) const
Return true iff counter is a multiple of the interval.
const std::string & outputFileName() const
Return outputFileName string.
int interval() const
Get interval value.
static void saveOptional(Serializable::OArchive &ar, Type &value, bool isActive)
Save an optional parameter value to an output archive.
Definition: Parameter.h:224
double energy()
Get total Coulomb energy.
BondPotential & bondPotential() const
Return BondPotential by reference.
Definition: MdSystem.h:537
MdCoulombPotential & coulombPotential() const
Return CoulombPotential by reference.
Definition: MdSystem.h:588
FileMaster & fileMaster()
Get the FileMaster object.