2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2008,2009,2010,2011,2012,2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
36 /* The source code in this file should be thread-safe.
37 Please keep it that way. */
40 #include "checkpoint.h"
49 #if GMX_NATIVE_WINDOWS
51 #include <sys/locking.h>
57 #include "buildinfo.h"
58 #include "gromacs/compat/make_unique.h"
59 #include "gromacs/fileio/filetypes.h"
60 #include "gromacs/fileio/gmxfio.h"
61 #include "gromacs/fileio/gmxfio-xdr.h"
62 #include "gromacs/fileio/xdr_datatype.h"
63 #include "gromacs/fileio/xdrf.h"
64 #include "gromacs/gmxlib/network.h"
65 #include "gromacs/math/vec.h"
66 #include "gromacs/math/vecdump.h"
67 #include "gromacs/math/vectypes.h"
68 #include "gromacs/mdtypes/awh-correlation-history.h"
69 #include "gromacs/mdtypes/awh-history.h"
70 #include "gromacs/mdtypes/commrec.h"
71 #include "gromacs/mdtypes/df_history.h"
72 #include "gromacs/mdtypes/edsamhistory.h"
73 #include "gromacs/mdtypes/energyhistory.h"
74 #include "gromacs/mdtypes/inputrec.h"
75 #include "gromacs/mdtypes/md_enums.h"
76 #include "gromacs/mdtypes/observableshistory.h"
77 #include "gromacs/mdtypes/pullhistory.h"
78 #include "gromacs/mdtypes/state.h"
79 #include "gromacs/mdtypes/swaphistory.h"
80 #include "gromacs/trajectory/trajectoryframe.h"
81 #include "gromacs/utility/arrayref.h"
82 #include "gromacs/utility/baseversion.h"
83 #include "gromacs/utility/cstringutil.h"
84 #include "gromacs/utility/fatalerror.h"
85 #include "gromacs/utility/futil.h"
86 #include "gromacs/utility/gmxassert.h"
87 #include "gromacs/utility/int64_to_int.h"
88 #include "gromacs/utility/programcontext.h"
89 #include "gromacs/utility/smalloc.h"
90 #include "gromacs/utility/sysinfo.h"
91 #include "gromacs/utility/txtdump.h"
97 #define CPT_MAGIC1 171817
98 #define CPT_MAGIC2 171819
99 #define CPTSTRLEN 1024
101 /*! \brief Enum of values that describe the contents of a cpt file
102 * whose format matches a version number
104 * The enum helps the code be more self-documenting and ensure merges
105 * do not silently resolve when two patches make the same bump. When
106 * adding new functionality, add a new element just above cptv_Count
107 * in this enumeration, and write code below that does the right thing
108 * according to the value of file_version.
111 cptv_Unknown = 17, /**< Version before numbering scheme */
112 cptv_RemoveBuildMachineInformation, /**< remove functionality that makes mdrun builds non-reproducible */
113 cptv_ComPrevStepAsPullGroupReference, /**< Allow using COM of previous step as pull group PBC reference */
114 cptv_PullAverage, /**< Added possibility to output average pull force and position */
115 cptv_Count /**< the total number of cptv versions */
118 /*! \brief Version number of the file format written to checkpoint
119 * files by this version of the code.
121 * cpt_version should normally only be changed, via adding a new field
122 * to cptv enumeration, when the header or footer format changes.
124 * The state data format itself is backward and forward compatible.
125 * But old code can not read a new entry that is present in the file
126 * (but can read a new format when new entries are not present).
128 * The cpt_version increases whenever the file format in the main
129 * development branch changes, due to an extension of the cptv enum above.
130 * Backward compatibility for reading old run input files is maintained
131 * by checking this version number against that of the file and then using
132 * the correct code path. */
133 static const int cpt_version = cptv_Count - 1;
136 const char *est_names[estNR] =
139 "box", "box-rel", "box-v", "pres_prev",
140 "nosehoover-xi", "thermostat-integral",
141 "x", "v", "sdx-unsupported", "CGp", "LD-rng-unsupported", "LD-rng-i-unsupported",
142 "disre_initf", "disre_rm3tav",
143 "orire_initf", "orire_Dtav",
144 "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng-unsupported", "MC-rng-i-unsupported",
149 eeksEKIN_N, eeksEKINH, eeksDEKINDL, eeksMVCOS, eeksEKINF, eeksEKINO, eeksEKINSCALEF, eeksEKINSCALEH, eeksVSCALE, eeksEKINTOTAL, eeksNR
152 static const char *eeks_names[eeksNR] =
154 "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
155 "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
159 eenhENERGY_N, eenhENERGY_AVER, eenhENERGY_SUM, eenhENERGY_NSUM,
160 eenhENERGY_SUM_SIM, eenhENERGY_NSUM_SIM,
161 eenhENERGY_NSTEPS, eenhENERGY_NSTEPS_SIM,
162 eenhENERGY_DELTA_H_NN,
163 eenhENERGY_DELTA_H_LIST,
164 eenhENERGY_DELTA_H_STARTTIME,
165 eenhENERGY_DELTA_H_STARTLAMBDA,
170 epullhPULL_NUMCOORDINATES, epullhPULL_NUMGROUPS, epullhPULL_NUMVALUESINXSUM,
171 epullhPULL_NUMVALUESINFSUM, epullhNR
175 epullcoordh_VALUE_REF_SUM, epullcoordh_VALUE_SUM, epullcoordh_DR01_SUM,
176 epullcoordh_DR23_SUM, epullcoordh_DR45_SUM, epullcoordh_FSCAL_SUM,
177 epullcoordh_DYNAX_SUM, epullcoordh_NR
181 epullgrouph_X_SUM, epullgrouph_NR
184 static const char *eenh_names[eenhNR] =
186 "energy_n", "energy_aver", "energy_sum", "energy_nsum",
187 "energy_sum_sim", "energy_nsum_sim",
188 "energy_nsteps", "energy_nsteps_sim",
190 "energy_delta_h_list",
191 "energy_delta_h_start_time",
192 "energy_delta_h_start_lambda"
195 static const char *ePullhNames[epullhNR] =
197 "pullhistory_numcoordinates", "pullhistory_numgroups", "pullhistory_numvaluesinxsum",
198 "pullhistory_numvaluesinfsum"
201 /* free energy history variables -- need to be preserved over checkpoint */
203 edfhBEQUIL, edfhNATLAMBDA, edfhWLHISTO, edfhWLDELTA, edfhSUMWEIGHTS, edfhSUMDG, edfhSUMMINVAR, edfhSUMVAR,
204 edfhACCUMP, edfhACCUMM, edfhACCUMP2, edfhACCUMM2, edfhTIJ, edfhTIJEMP, edfhNR
206 /* free energy history variable names */
207 static const char *edfh_names[edfhNR] =
209 "bEquilibrated", "N_at_state", "Wang-Landau Histogram", "Wang-Landau Delta", "Weights", "Free Energies", "minvar", "variance",
210 "accumulated_plus", "accumulated_minus", "accumulated_plus_2", "accumulated_minus_2", "Tij", "Tij_empirical"
213 /* AWH biasing history variables */
216 eawhhEQUILIBRATEHISTOGRAM,
219 eawhhCOORDPOINT, eawhhUMBRELLAGRIDPOINT,
221 eawhhLOGSCALEDSAMPLEWEIGHT,
223 eawhhFORCECORRELATIONGRID,
227 static const char *eawhh_names[eawhhNR] =
230 "awh_equilibrateHistogram",
233 "awh_coordpoint", "awh_umbrellaGridpoint",
235 "awh_logScaledSampleWeight",
237 "awh_forceCorrelationGrid"
245 static const char *epull_prev_step_com_names[epullsNR] =
247 "Pull groups prev step COM"
251 //! Higher level vector element type, only used for formatting checkpoint dumps
252 enum class CptElementType
254 integer, //!< integer
255 real, //!< float or double, not linked to precision of type real
256 real3, //!< float[3] or double[3], not linked to precision of type real
257 matrix3x3 //!< float[3][3] or double[3][3], not linked to precision of type real
260 //! \brief Parts of the checkpoint state, only used for reporting
263 microState, //!< The microstate of the simulated system
264 kineticEnergy, //!< Kinetic energy, needed for T/P-coupling state
265 energyHistory, //!< Energy observable statistics
266 freeEnergyHistory, //!< Free-energy state and observable statistics
267 accWeightHistogram, //!< Accelerated weight histogram method state
268 pullState, //!< COM of previous step.
269 pullHistory //!< Pull history statistics (sums since last written output)
272 //! \brief Return the name of a checkpoint entry based on part and part entry
273 static const char *entryName(StatePart part, int ecpt)
277 case StatePart::microState: return est_names [ecpt];
278 case StatePart::kineticEnergy: return eeks_names[ecpt];
279 case StatePart::energyHistory: return eenh_names[ecpt];
280 case StatePart::freeEnergyHistory: return edfh_names[ecpt];
281 case StatePart::accWeightHistogram: return eawhh_names[ecpt];
282 case StatePart::pullState: return epull_prev_step_com_names[ecpt];
283 case StatePart::pullHistory: return ePullhNames[ecpt];
289 static void cp_warning(FILE *fp)
291 fprintf(fp, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
294 [[ noreturn ]] static void cp_error()
296 gmx_fatal(FARGS, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
299 static void do_cpt_string_err(XDR *xd, const char *desc, gmx::ArrayRef<char> s, FILE *list)
301 char *data = s.data();
302 if (xdr_string(xd, &data, s.size()) == 0)
308 fprintf(list, "%s = %s\n", desc, data);
312 static int do_cpt_int(XDR *xd, const char *desc, int *i, FILE *list)
314 if (xdr_int(xd, i) == 0)
320 fprintf(list, "%s = %d\n", desc, *i);
325 static int do_cpt_u_chars(XDR *xd, const char *desc, int n, unsigned char *i, FILE *list)
329 fprintf(list, "%s = ", desc);
332 for (int j = 0; j < n && res; j++)
334 res &= xdr_u_char(xd, &i[j]);
337 fprintf(list, "%02x", i[j]);
352 static void do_cpt_int_err(XDR *xd, const char *desc, int *i, FILE *list)
354 if (do_cpt_int(xd, desc, i, list) < 0)
360 static void do_cpt_bool_err(XDR *xd, const char *desc, bool *b, FILE *list)
362 int i = static_cast<int>(*b);
364 if (do_cpt_int(xd, desc, &i, list) < 0)
372 static void do_cpt_step_err(XDR *xd, const char *desc, int64_t *i, FILE *list)
374 char buf[STEPSTRSIZE];
376 if (xdr_int64(xd, i) == 0)
382 fprintf(list, "%s = %s\n", desc, gmx_step_str(*i, buf));
386 static void do_cpt_double_err(XDR *xd, const char *desc, double *f, FILE *list)
388 if (xdr_double(xd, f) == 0)
394 fprintf(list, "%s = %f\n", desc, *f);
398 static void do_cpt_real_err(XDR *xd, real *f)
401 bool_t res = xdr_double(xd, f);
403 bool_t res = xdr_float(xd, f);
411 static void do_cpt_n_rvecs_err(XDR *xd, const char *desc, int n, rvec f[], FILE *list)
413 for (int i = 0; i < n; i++)
415 for (int j = 0; j < DIM; j++)
417 do_cpt_real_err(xd, &f[i][j]);
423 pr_rvecs(list, 0, desc, f, n);
427 template <typename T>
435 static const int value = xdr_datatype_int;
439 struct xdr_type<float>
441 static const int value = xdr_datatype_float;
445 struct xdr_type<double>
447 static const int value = xdr_datatype_double;
450 //! \brief Returns size in byte of an xdr_datatype
451 static inline unsigned int sizeOfXdrType(int xdrType)
455 case xdr_datatype_int:
457 case xdr_datatype_float:
458 return sizeof(float);
459 case xdr_datatype_double:
460 return sizeof(double);
461 default: GMX_RELEASE_ASSERT(false, "XDR data type not implemented");
467 //! \brief Returns the XDR process function for i/o of an XDR type
468 static inline xdrproc_t xdrProc(int xdrType)
472 case xdr_datatype_int:
473 return reinterpret_cast<xdrproc_t>(xdr_int);
474 case xdr_datatype_float:
475 return reinterpret_cast<xdrproc_t>(xdr_float);
476 case xdr_datatype_double:
477 return reinterpret_cast<xdrproc_t>(xdr_double);
478 default: GMX_RELEASE_ASSERT(false, "XDR data type not implemented");
484 /*! \brief Lists or only reads an xdr vector from checkpoint file
486 * When list!=NULL reads and lists the \p nf vector elements of type \p xdrType.
487 * The header for the print is set by \p part and \p ecpt.
488 * The formatting of the printing is set by \p cptElementType.
489 * When list==NULL only reads the elements.
491 static bool_t listXdrVector(XDR *xd, StatePart part, int ecpt, int nf, int xdrType,
492 FILE *list, CptElementType cptElementType)
496 const unsigned int elemSize = sizeOfXdrType(xdrType);
497 std::vector<char> data(nf*elemSize);
498 res = xdr_vector(xd, data.data(), nf, elemSize, xdrProc(xdrType));
504 case xdr_datatype_int:
505 pr_ivec(list, 0, entryName(part, ecpt), reinterpret_cast<const int *>(data.data()), nf, TRUE);
507 case xdr_datatype_float:
509 if (cptElementType == CptElementType::real3)
511 pr_rvecs(list, 0, entryName(part, ecpt), reinterpret_cast<const rvec *>(data.data()), nf/3);
516 /* Note: With double precision code dumping a single precision rvec will produce float iso rvec print, but that's a minor annoyance */
517 pr_fvec(list, 0, entryName(part, ecpt), reinterpret_cast<const float *>(data.data()), nf, TRUE);
520 case xdr_datatype_double:
522 if (cptElementType == CptElementType::real3)
524 pr_rvecs(list, 0, entryName(part, ecpt), reinterpret_cast<const rvec *>(data.data()), nf/3);
529 /* Note: With single precision code dumping a double precision rvec will produce float iso rvec print, but that's a minor annoyance */
530 pr_dvec(list, 0, entryName(part, ecpt), reinterpret_cast<const double *>(data.data()), nf, TRUE);
533 default: GMX_RELEASE_ASSERT(false, "Data type not implemented for listing");
540 //! \brief Convert a double array, typed char*, to float
541 gmx_unused static void convertArrayRealPrecision(const char *c, float *v, int n)
543 const double *d = reinterpret_cast<const double *>(c);
544 for (int i = 0; i < n; i++)
546 v[i] = static_cast<float>(d[i]);
550 //! \brief Convert a float array, typed char*, to double
551 static void convertArrayRealPrecision(const char *c, double *v, int n)
553 const float *f = reinterpret_cast<const float *>(c);
554 for (int i = 0; i < n; i++)
556 v[i] = static_cast<double>(f[i]);
560 //! \brief Generate an error for trying to convert to integer
561 static void convertArrayRealPrecision(const char gmx_unused *c, int gmx_unused *v, int gmx_unused n)
563 GMX_RELEASE_ASSERT(false, "We only expect type mismatches between float and double, not integer");
566 /*! \brief Low-level routine for reading/writing a vector of reals from/to file.
568 * This is the only routine that does the actually i/o of real vector,
569 * all other routines are intermediate level routines for specific real
570 * data types, calling this routine.
571 * Currently this routine is (too) complex, since it handles both real *
572 * and std::vector<real>. Using real * is deprecated and this routine
573 * will simplify a lot when only std::vector needs to be supported.
575 * The routine is generic to vectors with different allocators,
576 * because as part of reading a checkpoint there are vectors whose
577 * size is not known until reading has progressed far enough, so a
578 * resize method must be called.
580 * When not listing, we use either v or vector, depending on which is !=NULL.
581 * If nval >= 0, nval is used; on read this should match the passed value.
582 * If nval n<0, *nptr (with v) or vector->size() is used. On read using v,
583 * the value is stored in nptr
585 template<typename T, typename AllocatorType>
586 static int doVectorLow(XDR *xd, StatePart part, int ecpt, int sflags,
588 T **v, std::vector<T, AllocatorType> *vector,
589 FILE *list, CptElementType cptElementType)
591 GMX_RELEASE_ASSERT(list != nullptr || (v != nullptr && vector == nullptr) || (v == nullptr && vector != nullptr), "Without list, we should have exactly one of v and vector != NULL");
595 int numElemInTheFile;
600 GMX_RELEASE_ASSERT(nptr == nullptr, "With nval>=0 we should have nptr==NULL");
601 numElemInTheFile = nval;
607 GMX_RELEASE_ASSERT(nptr != nullptr, "With nval<0 we should have nptr!=NULL");
608 numElemInTheFile = *nptr;
612 numElemInTheFile = vector->size();
616 /* Read/write the vector element count */
617 res = xdr_int(xd, &numElemInTheFile);
622 /* Read/write the element data type */
623 constexpr int xdrTypeInTheCode = xdr_type<T>::value;
624 int xdrTypeInTheFile = xdrTypeInTheCode;
625 res = xdr_int(xd, &xdrTypeInTheFile);
631 if (list == nullptr && (sflags & (1 << ecpt)))
635 if (numElemInTheFile != nval)
637 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", entryName(part, ecpt), nval, numElemInTheFile);
640 else if (nptr != nullptr)
642 *nptr = numElemInTheFile;
645 bool typesMatch = (xdrTypeInTheFile == xdrTypeInTheCode);
649 sprintf(buf, "mismatch for state entry %s, code precision is %s, file precision is %s",
650 entryName(part, ecpt),
651 xdr_datatype_names[xdrTypeInTheCode],
652 xdr_datatype_names[xdrTypeInTheFile]);
654 /* Matching int and real should never occur, but check anyhow */
655 if (xdrTypeInTheFile == xdr_datatype_int ||
656 xdrTypeInTheCode == xdr_datatype_int)
658 gmx_fatal(FARGS, "Type %s: incompatible checkpoint formats or corrupted checkpoint file.", buf);
667 snew(*v, numElemInTheFile);
673 /* This conditional ensures that we don't resize on write.
674 * In particular in the state where this code was written
675 * vector has a size of numElemInThefile and we
676 * don't want to lose that padding here.
678 if (vector->size() < static_cast<unsigned int>(numElemInTheFile))
680 vector->resize(numElemInTheFile);
688 vChar = reinterpret_cast<char *>(vp);
692 snew(vChar, numElemInTheFile*sizeOfXdrType(xdrTypeInTheFile));
694 res = xdr_vector(xd, vChar,
695 numElemInTheFile, sizeOfXdrType(xdrTypeInTheFile),
696 xdrProc(xdrTypeInTheFile));
704 /* In the old code float-double conversion came for free.
705 * In the new code we still support it, mainly because
706 * the tip4p_continue regression test makes use of this.
707 * It's an open question if we do or don't want to allow this.
709 convertArrayRealPrecision(vChar, vp, numElemInTheFile);
715 res = listXdrVector(xd, part, ecpt, numElemInTheFile, xdrTypeInTheFile,
716 list, cptElementType);
722 //! \brief Read/Write a std::vector, on read checks the number of elements matches \p numElements, if specified.
723 template <typename T>
724 static int doVector(XDR *xd, StatePart part, int ecpt, int sflags,
725 std::vector<T> *vector, FILE *list, int numElements = -1)
727 return doVectorLow<T>(xd, part, ecpt, sflags, numElements, nullptr, nullptr, vector, list, CptElementType::real);
730 //! \brief Read/Write an ArrayRef<real>.
731 static int doRealArrayRef(XDR *xd, StatePart part, int ecpt, int sflags,
732 gmx::ArrayRef<real> vector, FILE *list)
734 real *v_real = vector.data();
735 return doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags, vector.size(), nullptr, &v_real, nullptr, list, CptElementType::real);
738 //! Convert from view of RVec to view of real.
739 static gmx::ArrayRef<real>
740 realArrayRefFromRVecArrayRef(gmx::ArrayRef<gmx::RVec> ofRVecs)
742 return gmx::arrayRefFromArray<real>(reinterpret_cast<real *>(ofRVecs.data()), ofRVecs.size() * DIM);
745 //! \brief Read/Write a PaddedVector whose value_type is RVec.
746 template <typename PaddedVectorOfRVecType>
747 static int doRvecVector(XDR *xd, StatePart part, int ecpt, int sflags,
748 PaddedVectorOfRVecType *v, int numAtoms, FILE *list)
750 const int numReals = numAtoms*DIM;
754 GMX_RELEASE_ASSERT(sflags & (1 << ecpt), "When not listing, the flag for the entry should be set when requesting i/o");
755 GMX_RELEASE_ASSERT(v->size() == numAtoms, "v should have sufficient size for numAtoms");
757 return doRealArrayRef(xd, part, ecpt, sflags, realArrayRefFromRVecArrayRef(makeArrayRef(*v)), list);
761 // Use the rebind facility to change the value_type of the
762 // allocator from RVec to real.
763 using realAllocator = typename std::allocator_traits<typename PaddedVectorOfRVecType::allocator_type>::template rebind_alloc<real>;
764 return doVectorLow<real, realAllocator>(xd, part, ecpt, sflags, numReals, nullptr, nullptr, nullptr, list, CptElementType::real);
768 /* This function stores n along with the reals for reading,
769 * but on reading it assumes that n matches the value in the checkpoint file,
770 * a fatal error is generated when this is not the case.
772 static int do_cpte_reals(XDR *xd, StatePart part, int ecpt, int sflags,
773 int n, real **v, FILE *list)
775 return doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags, n, nullptr, v, nullptr, list, CptElementType::real);
778 /* This function does the same as do_cpte_reals,
779 * except that on reading it ignores the passed value of *n
780 * and stores the value read from the checkpoint file in *n.
782 static int do_cpte_n_reals(XDR *xd, StatePart part, int ecpt, int sflags,
783 int *n, real **v, FILE *list)
785 return doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags, -1, n, v, nullptr, list, CptElementType::real);
788 static int do_cpte_real(XDR *xd, StatePart part, int ecpt, int sflags,
791 return doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags, 1, nullptr, &r, nullptr, list, CptElementType::real);
794 static int do_cpte_ints(XDR *xd, StatePart part, int ecpt, int sflags,
795 int n, int **v, FILE *list)
797 return doVectorLow < int, std::allocator < int>>(xd, part, ecpt, sflags, n, nullptr, v, nullptr, list, CptElementType::integer);
800 static int do_cpte_int(XDR *xd, StatePart part, int ecpt, int sflags,
803 return do_cpte_ints(xd, part, ecpt, sflags, 1, &i, list);
806 static int do_cpte_bool(XDR *xd, StatePart part, int ecpt, int sflags,
809 int i = static_cast<int>(*b);
810 int ret = do_cpte_int(xd, part, ecpt, sflags, &i, list);
815 static int do_cpte_doubles(XDR *xd, StatePart part, int ecpt, int sflags,
816 int n, double **v, FILE *list)
818 return doVectorLow < double, std::allocator < double>>(xd, part, ecpt, sflags, n, nullptr, v, nullptr, list, CptElementType::real);
821 static int do_cpte_double(XDR *xd, StatePart part, int ecpt, int sflags,
822 double *r, FILE *list)
824 return do_cpte_doubles(xd, part, ecpt, sflags, 1, &r, list);
827 static int do_cpte_matrix(XDR *xd, StatePart part, int ecpt, int sflags,
828 matrix v, FILE *list)
834 ret = doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags,
835 DIM*DIM, nullptr, &vr, nullptr, nullptr, CptElementType::matrix3x3);
837 if (list && ret == 0)
839 pr_rvecs(list, 0, entryName(part, ecpt), v, DIM);
846 static int do_cpte_nmatrix(XDR *xd, StatePart part, int ecpt, int sflags,
847 int n, real **v, FILE *list)
851 char name[CPTSTRLEN];
858 for (i = 0; i < n; i++)
860 reti = doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags, n, nullptr, &(v[i]), nullptr, nullptr, CptElementType::matrix3x3);
861 if (list && reti == 0)
863 sprintf(name, "%s[%d]", entryName(part, ecpt), i);
864 pr_reals(list, 0, name, v[i], n);
874 static int do_cpte_matrices(XDR *xd, StatePart part, int ecpt, int sflags,
875 int n, matrix **v, FILE *list)
878 matrix *vp, *va = nullptr;
884 res = xdr_int(xd, &nf);
889 if (list == nullptr && nf != n)
891 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", entryName(part, ecpt), n, nf);
893 if (list || !(sflags & (1<<ecpt)))
906 snew(vr, nf*DIM*DIM);
907 for (i = 0; i < nf; i++)
909 for (j = 0; j < DIM; j++)
911 for (k = 0; k < DIM; k++)
913 vr[(i*DIM+j)*DIM+k] = vp[i][j][k];
917 ret = doVectorLow < real, std::allocator < real>>(xd, part, ecpt, sflags,
918 nf*DIM*DIM, nullptr, &vr, nullptr, nullptr,
919 CptElementType::matrix3x3);
920 for (i = 0; i < nf; i++)
922 for (j = 0; j < DIM; j++)
924 for (k = 0; k < DIM; k++)
926 vp[i][j][k] = vr[(i*DIM+j)*DIM+k];
932 if (list && ret == 0)
934 for (i = 0; i < nf; i++)
936 pr_rvecs(list, 0, entryName(part, ecpt), vp[i], DIM);
947 // TODO Expand this into being a container of all data for
948 // serialization of a checkpoint, which can be stored by the caller
949 // (e.g. so that mdrun doesn't have to open the checkpoint twice).
950 // This will separate issues of allocation from those of
951 // serialization, help separate comparison from reading, and have
952 // better defined transformation functions to/from trajectory frame
955 // Several fields were once written to checkpoint file headers, but
956 // have been removed. So that old files can continue to be read,
957 // the names of such fields contain the string "_UNUSED" so that it
958 // is clear they should not be used.
959 struct CheckpointHeaderContents
962 char version[CPTSTRLEN];
963 char btime_UNUSED[CPTSTRLEN];
964 char buser_UNUSED[CPTSTRLEN];
965 char bhost_UNUSED[CPTSTRLEN];
967 char fprog[CPTSTRLEN];
968 char ftime[CPTSTRLEN];
984 int flagsPullHistory;
991 static void do_cpt_header(XDR *xd, gmx_bool bRead, FILE *list,
992 CheckpointHeaderContents *contents)
1005 res = xdr_int(xd, &magic);
1008 gmx_fatal(FARGS, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
1010 if (magic != CPT_MAGIC1)
1012 gmx_fatal(FARGS, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
1013 "The checkpoint file is corrupted or not a checkpoint file",
1019 gmx_gethostname(fhost, 255);
1021 do_cpt_string_err(xd, "GROMACS version", contents->version, list);
1022 // The following fields are no longer ever written with meaningful
1023 // content, but because they precede the file version, there is no
1024 // good way for new code to read the old and new formats, nor a
1025 // good way for old code to avoid giving an error while reading a
1026 // new format. So we read and write a field that no longer has a
1028 do_cpt_string_err(xd, "GROMACS build time UNUSED", contents->btime_UNUSED, list);
1029 do_cpt_string_err(xd, "GROMACS build user UNUSED", contents->buser_UNUSED, list);
1030 do_cpt_string_err(xd, "GROMACS build host UNUSED", contents->bhost_UNUSED, list);
1031 do_cpt_string_err(xd, "generating program", contents->fprog, list);
1032 do_cpt_string_err(xd, "generation time", contents->ftime, list);
1033 contents->file_version = cpt_version;
1034 do_cpt_int_err(xd, "checkpoint file version", &contents->file_version, list);
1035 if (contents->file_version > cpt_version)
1037 gmx_fatal(FARGS, "Attempting to read a checkpoint file of version %d with code of version %d\n", contents->file_version, cpt_version);
1039 if (contents->file_version >= 13)
1041 do_cpt_int_err(xd, "GROMACS double precision", &contents->double_prec, list);
1045 contents->double_prec = -1;
1047 if (contents->file_version >= 12)
1049 do_cpt_string_err(xd, "generating host", fhost, list);
1051 do_cpt_int_err(xd, "#atoms", &contents->natoms, list);
1052 do_cpt_int_err(xd, "#T-coupling groups", &contents->ngtc, list);
1053 if (contents->file_version >= 10)
1055 do_cpt_int_err(xd, "#Nose-Hoover T-chains", &contents->nhchainlength, list);
1059 contents->nhchainlength = 1;
1061 if (contents->file_version >= 11)
1063 do_cpt_int_err(xd, "#Nose-Hoover T-chains for barostat ", &contents->nnhpres, list);
1067 contents->nnhpres = 0;
1069 if (contents->file_version >= 14)
1071 do_cpt_int_err(xd, "# of total lambda states ", &contents->nlambda, list);
1075 contents->nlambda = 0;
1077 do_cpt_int_err(xd, "integrator", &contents->eIntegrator, list);
1078 if (contents->file_version >= 3)
1080 do_cpt_int_err(xd, "simulation part #", &contents->simulation_part, list);
1084 contents->simulation_part = 1;
1086 if (contents->file_version >= 5)
1088 do_cpt_step_err(xd, "step", &contents->step, list);
1093 do_cpt_int_err(xd, "step", &idum, list);
1094 contents->step = static_cast<int64_t>(idum);
1096 do_cpt_double_err(xd, "t", &contents->t, list);
1097 do_cpt_int_err(xd, "#PP-ranks", &contents->nnodes, list);
1098 do_cpt_int_err(xd, "dd_nc[x]", &contents->dd_nc[XX], list);
1099 do_cpt_int_err(xd, "dd_nc[y]", &contents->dd_nc[YY], list);
1100 do_cpt_int_err(xd, "dd_nc[z]", &contents->dd_nc[ZZ], list);
1101 do_cpt_int_err(xd, "#PME-only ranks", &contents->npme, list);
1102 do_cpt_int_err(xd, "state flags", &contents->flags_state, list);
1103 if (contents->file_version >= 4)
1105 do_cpt_int_err(xd, "ekin data flags", &contents->flags_eks, list);
1106 do_cpt_int_err(xd, "energy history flags", &contents->flags_enh, list);
1110 contents->flags_eks = 0;
1111 contents->flags_enh = (contents->flags_state >> (estORIRE_DTAV+1));
1112 contents->flags_state = (contents->flags_state & ~((1<<(estORIRE_DTAV+1)) |
1113 (1<<(estORIRE_DTAV+2)) |
1114 (1<<(estORIRE_DTAV+3))));
1116 if (contents->file_version >= 14)
1118 do_cpt_int_err(xd, "df history flags", &contents->flags_dfh, list);
1122 contents->flags_dfh = 0;
1125 if (contents->file_version >= 15)
1127 do_cpt_int_err(xd, "ED data sets", &contents->nED, list);
1134 if (contents->file_version >= 16)
1136 do_cpt_int_err(xd, "swap", &contents->eSwapCoords, list);
1140 contents->eSwapCoords = eswapNO;
1143 if (contents->file_version >= 17)
1145 do_cpt_int_err(xd, "AWH history flags", &contents->flags_awhh, list);
1149 contents->flags_awhh = 0;
1152 if (contents->file_version >= 18)
1154 do_cpt_int_err(xd, "pull history flags", &contents->flagsPullHistory, list);
1158 contents->flagsPullHistory = 0;
1162 static int do_cpt_footer(XDR *xd, int file_version)
1167 if (file_version >= 2)
1170 res = xdr_int(xd, &magic);
1175 if (magic != CPT_MAGIC2)
1184 static int do_cpt_state(XDR *xd,
1185 int fflags, t_state *state,
1189 const StatePart part = StatePart::microState;
1190 const int sflags = state->flags;
1191 // If reading, state->natoms was probably just read, so
1192 // allocations need to be managed. If writing, this won't change
1193 // anything that matters.
1194 state_change_natoms(state, state->natoms);
1195 for (int i = 0; (i < estNR && ret == 0); i++)
1197 if (fflags & (1<<i))
1201 case estLAMBDA: ret = doRealArrayRef(xd, part, i, sflags, gmx::arrayRefFromArray<real>(state->lambda.data(), state->lambda.size()), list); break;
1202 case estFEPSTATE: ret = do_cpte_int (xd, part, i, sflags, &state->fep_state, list); break;
1203 case estBOX: ret = do_cpte_matrix(xd, part, i, sflags, state->box, list); break;
1204 case estBOX_REL: ret = do_cpte_matrix(xd, part, i, sflags, state->box_rel, list); break;
1205 case estBOXV: ret = do_cpte_matrix(xd, part, i, sflags, state->boxv, list); break;
1206 case estPRES_PREV: ret = do_cpte_matrix(xd, part, i, sflags, state->pres_prev, list); break;
1207 case estSVIR_PREV: ret = do_cpte_matrix(xd, part, i, sflags, state->svir_prev, list); break;
1208 case estFVIR_PREV: ret = do_cpte_matrix(xd, part, i, sflags, state->fvir_prev, list); break;
1209 case estNH_XI: ret = doVector<double>(xd, part, i, sflags, &state->nosehoover_xi, list); break;
1210 case estNH_VXI: ret = doVector<double>(xd, part, i, sflags, &state->nosehoover_vxi, list); break;
1211 case estNHPRES_XI: ret = doVector<double>(xd, part, i, sflags, &state->nhpres_xi, list); break;
1212 case estNHPRES_VXI: ret = doVector<double>(xd, part, i, sflags, &state->nhpres_vxi, list); break;
1213 case estTHERM_INT: ret = doVector<double>(xd, part, i, sflags, &state->therm_integral, list); break;
1214 case estBAROS_INT: ret = do_cpte_double(xd, part, i, sflags, &state->baros_integral, list); break;
1215 case estVETA: ret = do_cpte_real(xd, part, i, sflags, &state->veta, list); break;
1216 case estVOL0: ret = do_cpte_real(xd, part, i, sflags, &state->vol0, list); break;
1217 case estX: ret = doRvecVector(xd, part, i, sflags, &state->x, state->natoms, list); break;
1218 case estV: ret = doRvecVector(xd, part, i, sflags, &state->v, state->natoms, list); break;
1219 /* The RNG entries are no longer written,
1220 * the next 4 lines are only for reading old files.
1222 case estLD_RNG_NOTSUPPORTED: ret = do_cpte_ints(xd, part, i, sflags, 0, nullptr, list); break;
1223 case estLD_RNGI_NOTSUPPORTED: ret = do_cpte_ints(xd, part, i, sflags, 0, nullptr, list); break;
1224 case estMC_RNG_NOTSUPPORTED: ret = do_cpte_ints(xd, part, i, sflags, 0, nullptr, list); break;
1225 case estMC_RNGI_NOTSUPPORTED: ret = do_cpte_ints(xd, part, i, sflags, 0, nullptr, list); break;
1226 case estDISRE_INITF: ret = do_cpte_real (xd, part, i, sflags, &state->hist.disre_initf, list); break;
1227 case estDISRE_RM3TAV: ret = do_cpte_n_reals(xd, part, i, sflags, &state->hist.ndisrepairs, &state->hist.disre_rm3tav, list); break;
1228 case estORIRE_INITF: ret = do_cpte_real (xd, part, i, sflags, &state->hist.orire_initf, list); break;
1229 case estORIRE_DTAV: ret = do_cpte_n_reals(xd, part, i, sflags, &state->hist.norire_Dtav, &state->hist.orire_Dtav, list); break;
1230 case estPULLCOMPREVSTEP: ret = doVector<double>(xd, part, i, sflags, &state->pull_com_prev_step, list); break;
1232 gmx_fatal(FARGS, "Unknown state entry %d\n"
1233 "You are reading a checkpoint file written by different code, which is not supported", i);
1240 static int do_cpt_ekinstate(XDR *xd, int fflags, ekinstate_t *ekins,
1245 const StatePart part = StatePart::kineticEnergy;
1246 for (int i = 0; (i < eeksNR && ret == 0); i++)
1248 if (fflags & (1<<i))
1253 case eeksEKIN_N: ret = do_cpte_int(xd, part, i, fflags, &ekins->ekin_n, list); break;
1254 case eeksEKINH: ret = do_cpte_matrices(xd, part, i, fflags, ekins->ekin_n, &ekins->ekinh, list); break;
1255 case eeksEKINF: ret = do_cpte_matrices(xd, part, i, fflags, ekins->ekin_n, &ekins->ekinf, list); break;
1256 case eeksEKINO: ret = do_cpte_matrices(xd, part, i, fflags, ekins->ekin_n, &ekins->ekinh_old, list); break;
1257 case eeksEKINTOTAL: ret = do_cpte_matrix(xd, part, i, fflags, ekins->ekin_total, list); break;
1258 case eeksEKINSCALEF: ret = doVector<double>(xd, part, i, fflags, &ekins->ekinscalef_nhc, list); break;
1259 case eeksVSCALE: ret = doVector<double>(xd, part, i, fflags, &ekins->vscale_nhc, list); break;
1260 case eeksEKINSCALEH: ret = doVector<double>(xd, part, i, fflags, &ekins->ekinscaleh_nhc, list); break;
1261 case eeksDEKINDL: ret = do_cpte_real(xd, part, i, fflags, &ekins->dekindl, list); break;
1262 case eeksMVCOS: ret = do_cpte_real(xd, part, i, fflags, &ekins->mvcos, list); break;
1264 gmx_fatal(FARGS, "Unknown ekin data state entry %d\n"
1265 "You are probably reading a new checkpoint file with old code", i);
1274 static int do_cpt_swapstate(XDR *xd, gmx_bool bRead,
1275 int eSwapCoords, swaphistory_t *swapstate, FILE *list)
1277 int swap_cpt_version = 2;
1279 if (eSwapCoords == eswapNO)
1284 swapstate->bFromCpt = bRead;
1285 swapstate->eSwapCoords = eSwapCoords;
1287 do_cpt_int_err(xd, "swap checkpoint version", &swap_cpt_version, list);
1288 if (bRead && swap_cpt_version < 2)
1290 gmx_fatal(FARGS, "Cannot read checkpoint files that were written with old versions"
1291 "of the ion/water position swapping protocol.\n");
1294 do_cpt_int_err(xd, "swap coupling steps", &swapstate->nAverage, list);
1296 /* When reading, init_swapcoords has not been called yet,
1297 * so we have to allocate memory first. */
1298 do_cpt_int_err(xd, "number of ion types", &swapstate->nIonTypes, list);
1301 snew(swapstate->ionType, swapstate->nIonTypes);
1304 for (int ic = 0; ic < eCompNR; ic++)
1306 for (int ii = 0; ii < swapstate->nIonTypes; ii++)
1308 swapstateIons_t *gs = &swapstate->ionType[ii];
1312 do_cpt_int_err(xd, "swap requested atoms", &gs->nMolReq[ic], list);
1316 do_cpt_int_err(xd, "swap requested atoms p", gs->nMolReq_p[ic], list);
1321 do_cpt_int_err(xd, "swap influx net", &gs->inflow_net[ic], list);
1325 do_cpt_int_err(xd, "swap influx net p", gs->inflow_net_p[ic], list);
1328 if (bRead && (nullptr == gs->nMolPast[ic]) )
1330 snew(gs->nMolPast[ic], swapstate->nAverage);
1333 for (int j = 0; j < swapstate->nAverage; j++)
1337 do_cpt_int_err(xd, "swap past atom counts", &gs->nMolPast[ic][j], list);
1341 do_cpt_int_err(xd, "swap past atom counts p", &gs->nMolPast_p[ic][j], list);
1347 /* Ion flux per channel */
1348 for (int ic = 0; ic < eChanNR; ic++)
1350 for (int ii = 0; ii < swapstate->nIonTypes; ii++)
1352 swapstateIons_t *gs = &swapstate->ionType[ii];
1356 do_cpt_int_err(xd, "channel flux A->B", &gs->fluxfromAtoB[ic], list);
1360 do_cpt_int_err(xd, "channel flux A->B p", gs->fluxfromAtoB_p[ic], list);
1365 /* Ion flux leakage */
1368 do_cpt_int_err(xd, "flux leakage", &swapstate->fluxleak, list);
1372 do_cpt_int_err(xd, "flux leakage", swapstate->fluxleak_p, list);
1376 for (int ii = 0; ii < swapstate->nIonTypes; ii++)
1378 swapstateIons_t *gs = &swapstate->ionType[ii];
1380 do_cpt_int_err(xd, "number of ions", &gs->nMol, list);
1384 snew(gs->channel_label, gs->nMol);
1385 snew(gs->comp_from, gs->nMol);
1388 do_cpt_u_chars(xd, "channel history", gs->nMol, gs->channel_label, list);
1389 do_cpt_u_chars(xd, "domain history", gs->nMol, gs->comp_from, list);
1392 /* Save the last known whole positions to checkpoint
1393 * file to be able to also make multimeric channels whole in PBC */
1394 do_cpt_int_err(xd, "Ch0 atoms", &swapstate->nat[eChan0], list);
1395 do_cpt_int_err(xd, "Ch1 atoms", &swapstate->nat[eChan1], list);
1398 snew(swapstate->xc_old_whole[eChan0], swapstate->nat[eChan0]);
1399 snew(swapstate->xc_old_whole[eChan1], swapstate->nat[eChan1]);
1400 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], swapstate->xc_old_whole[eChan0], list);
1401 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], swapstate->xc_old_whole[eChan1], list);
1405 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], *swapstate->xc_old_whole_p[eChan0], list);
1406 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], *swapstate->xc_old_whole_p[eChan1], list);
1413 static int do_cpt_enerhist(XDR *xd, gmx_bool bRead,
1414 int fflags, energyhistory_t *enerhist,
1424 GMX_RELEASE_ASSERT(enerhist != nullptr, "With energy history, we need a valid enerhist pointer");
1426 /* This is stored/read for backward compatibility */
1427 int energyHistoryNumEnergies = 0;
1430 enerhist->nsteps = 0;
1432 enerhist->nsteps_sim = 0;
1433 enerhist->nsum_sim = 0;
1435 else if (enerhist != nullptr)
1437 energyHistoryNumEnergies = enerhist->ener_sum_sim.size();
1440 delta_h_history_t *deltaH = enerhist->deltaHForeignLambdas.get();
1441 const StatePart part = StatePart::energyHistory;
1442 for (int i = 0; (i < eenhNR && ret == 0); i++)
1444 if (fflags & (1<<i))
1448 case eenhENERGY_N: ret = do_cpte_int(xd, part, i, fflags, &energyHistoryNumEnergies, list); break;
1449 case eenhENERGY_AVER: ret = doVector<double>(xd, part, i, fflags, &enerhist->ener_ave, list); break;
1450 case eenhENERGY_SUM: ret = doVector<double>(xd, part, i, fflags, &enerhist->ener_sum, list); break;
1451 case eenhENERGY_NSUM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum, list); break;
1452 case eenhENERGY_SUM_SIM: ret = doVector<double>(xd, part, i, fflags, &enerhist->ener_sum_sim, list); break;
1453 case eenhENERGY_NSUM_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum_sim, list); break;
1454 case eenhENERGY_NSTEPS: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps, list); break;
1455 case eenhENERGY_NSTEPS_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps_sim, list); break;
1456 case eenhENERGY_DELTA_H_NN:
1459 if (!bRead && deltaH != nullptr)
1461 numDeltaH = deltaH->dh.size();
1463 do_cpt_int_err(xd, eenh_names[i], &numDeltaH, list);
1466 if (deltaH == nullptr)
1468 enerhist->deltaHForeignLambdas = gmx::compat::make_unique<delta_h_history_t>();
1469 deltaH = enerhist->deltaHForeignLambdas.get();
1471 deltaH->dh.resize(numDeltaH);
1472 deltaH->start_lambda_set = FALSE;
1476 case eenhENERGY_DELTA_H_LIST:
1477 for (auto dh : deltaH->dh)
1479 ret = doVector<real>(xd, part, i, fflags, &dh, list);
1482 case eenhENERGY_DELTA_H_STARTTIME:
1483 ret = do_cpte_double(xd, part, i, fflags, &(deltaH->start_time), list); break;
1484 case eenhENERGY_DELTA_H_STARTLAMBDA:
1485 ret = do_cpte_double(xd, part, i, fflags, &(deltaH->start_lambda), list); break;
1487 gmx_fatal(FARGS, "Unknown energy history entry %d\n"
1488 "You are probably reading a new checkpoint file with old code", i);
1493 if ((fflags & (1<<eenhENERGY_SUM)) && !(fflags & (1<<eenhENERGY_SUM_SIM)))
1495 /* Assume we have an old file format and copy sum to sum_sim */
1496 enerhist->ener_sum_sim = enerhist->ener_sum;
1499 if ( (fflags & (1<<eenhENERGY_NSUM)) &&
1500 !(fflags & (1<<eenhENERGY_NSTEPS)))
1502 /* Assume we have an old file format and copy nsum to nsteps */
1503 enerhist->nsteps = enerhist->nsum;
1505 if ( (fflags & (1<<eenhENERGY_NSUM_SIM)) &&
1506 !(fflags & (1<<eenhENERGY_NSTEPS_SIM)))
1508 /* Assume we have an old file format and copy nsum to nsteps */
1509 enerhist->nsteps_sim = enerhist->nsum_sim;
1515 static int doCptPullCoordHist(XDR *xd, PullCoordinateHistory *pullCoordHist,
1516 const StatePart part, FILE *list)
1521 flags |= ((1<<epullcoordh_VALUE_REF_SUM) | (1<<epullcoordh_VALUE_SUM) | (1<<epullcoordh_DR01_SUM) |
1522 (1<<epullcoordh_DR23_SUM) | (1<<epullcoordh_DR45_SUM) | (1<<epullcoordh_FSCAL_SUM));
1524 for (int i = 0; i < epullcoordh_NR && ret == 0; i++)
1528 case epullcoordh_VALUE_REF_SUM: ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->valueRef), list); break;
1529 case epullcoordh_VALUE_SUM: ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->value), list); break;
1530 case epullcoordh_DR01_SUM:
1531 for (int j = 0; j < DIM && ret == 0; j++)
1533 ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->dr01[j]), list);
1536 case epullcoordh_DR23_SUM:
1537 for (int j = 0; j < DIM && ret == 0; j++)
1539 ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->dr23[j]), list);
1542 case epullcoordh_DR45_SUM:
1543 for (int j = 0; j < DIM && ret == 0; j++)
1545 ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->dr45[j]), list);
1548 case epullcoordh_FSCAL_SUM: ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->scalarForce), list); break;
1549 case epullcoordh_DYNAX_SUM:
1550 for (int j = 0; j < DIM && ret == 0; j++)
1552 ret = do_cpte_double(xd, part, i, flags, &(pullCoordHist->dynaX[j]), list);
1561 static int doCptPullGroupHist(XDR *xd, PullGroupHistory *pullGroupHist,
1562 const StatePart part, FILE *list)
1567 flags |= ((1<<epullgrouph_X_SUM));
1569 for (int i = 0; i < epullgrouph_NR; i++)
1573 case epullgrouph_X_SUM:
1574 for (int j = 0; j < DIM && ret == 0; j++)
1576 ret = do_cpte_double(xd, part, i, flags, &(pullGroupHist->x[j]), list);
1586 static int doCptPullHist(XDR *xd, gmx_bool bRead,
1587 int fflags, PullHistory *pullHist,
1588 const StatePart part,
1592 int pullHistoryNumCoordinates = 0;
1593 int pullHistoryNumGroups = 0;
1595 /* Retain the number of terms in the sum and the number of coordinates (used for writing
1596 * average pull forces and coordinates) in the pull history, in temporary variables,
1597 * in case they cannot be read from the checkpoint, in order to have backward compatibility */
1600 pullHist->numValuesInXSum = 0;
1601 pullHist->numValuesInFSum = 0;
1603 else if (pullHist != nullptr)
1605 pullHistoryNumCoordinates = pullHist->pullCoordinateSums.size();
1606 pullHistoryNumGroups = pullHist->pullGroupSums.size();
1610 GMX_RELEASE_ASSERT(fflags == 0, "Without pull history, all flags should be off");
1613 for (int i = 0; (i < epullhNR && ret == 0); i++)
1615 if (fflags & (1<<i))
1619 case epullhPULL_NUMCOORDINATES: ret = do_cpte_int(xd, part, i, fflags, &pullHistoryNumCoordinates, list); break;
1620 case epullhPULL_NUMGROUPS: do_cpt_int_err(xd, eenh_names[i], &pullHistoryNumGroups, list); break;
1621 case epullhPULL_NUMVALUESINXSUM: do_cpt_int_err(xd, eenh_names[i], &pullHist->numValuesInXSum, list); break;
1622 case epullhPULL_NUMVALUESINFSUM: do_cpt_int_err(xd, eenh_names[i], &pullHist->numValuesInFSum, list); break;
1624 gmx_fatal(FARGS, "Unknown pull history entry %d\n"
1625 "You are probably reading a new checkpoint file with old code", i);
1631 pullHist->pullCoordinateSums.resize(pullHistoryNumCoordinates);
1632 pullHist->pullGroupSums.resize(pullHistoryNumGroups);
1634 if (pullHist->numValuesInXSum > 0 || pullHist->numValuesInFSum > 0)
1636 for (size_t i = 0; i < pullHist->pullCoordinateSums.size() && ret == 0; i++)
1638 ret = doCptPullCoordHist(xd, &(pullHist->pullCoordinateSums[i]), part, list);
1640 for (size_t i = 0; i < pullHist->pullGroupSums.size() && ret == 0; i++)
1642 ret = doCptPullGroupHist(xd, &(pullHist->pullGroupSums[i]), part, list);
1649 static int do_cpt_df_hist(XDR *xd, int fflags, int nlambda, df_history_t **dfhistPtr, FILE *list)
1658 if (*dfhistPtr == nullptr)
1660 snew(*dfhistPtr, 1);
1661 (*dfhistPtr)->nlambda = nlambda;
1662 init_df_history(*dfhistPtr, nlambda);
1664 df_history_t *dfhist = *dfhistPtr;
1666 const StatePart part = StatePart::freeEnergyHistory;
1667 for (int i = 0; (i < edfhNR && ret == 0); i++)
1669 if (fflags & (1<<i))
1673 case edfhBEQUIL: ret = do_cpte_bool(xd, part, i, fflags, &dfhist->bEquil, list); break;
1674 case edfhNATLAMBDA: ret = do_cpte_ints(xd, part, i, fflags, nlambda, &dfhist->n_at_lam, list); break;
1675 case edfhWLHISTO: ret = do_cpte_reals(xd, part, i, fflags, nlambda, &dfhist->wl_histo, list); break;
1676 case edfhWLDELTA: ret = do_cpte_real(xd, part, i, fflags, &dfhist->wl_delta, list); break;
1677 case edfhSUMWEIGHTS: ret = do_cpte_reals(xd, part, i, fflags, nlambda, &dfhist->sum_weights, list); break;
1678 case edfhSUMDG: ret = do_cpte_reals(xd, part, i, fflags, nlambda, &dfhist->sum_dg, list); break;
1679 case edfhSUMMINVAR: ret = do_cpte_reals(xd, part, i, fflags, nlambda, &dfhist->sum_minvar, list); break;
1680 case edfhSUMVAR: ret = do_cpte_reals(xd, part, i, fflags, nlambda, &dfhist->sum_variance, list); break;
1681 case edfhACCUMP: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->accum_p, list); break;
1682 case edfhACCUMM: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->accum_m, list); break;
1683 case edfhACCUMP2: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->accum_p2, list); break;
1684 case edfhACCUMM2: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->accum_m2, list); break;
1685 case edfhTIJ: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->Tij, list); break;
1686 case edfhTIJEMP: ret = do_cpte_nmatrix(xd, part, i, fflags, nlambda, dfhist->Tij_empirical, list); break;
1689 gmx_fatal(FARGS, "Unknown df history entry %d\n"
1690 "You are probably reading a new checkpoint file with old code", i);
1699 /* This function stores the last whole configuration of the reference and
1700 * average structure in the .cpt file
1702 static int do_cpt_EDstate(XDR *xd, gmx_bool bRead,
1703 int nED, edsamhistory_t *EDstate, FILE *list)
1710 EDstate->bFromCpt = bRead;
1713 /* When reading, init_edsam has not been called yet,
1714 * so we have to allocate memory first. */
1717 snew(EDstate->nref, EDstate->nED);
1718 snew(EDstate->old_sref, EDstate->nED);
1719 snew(EDstate->nav, EDstate->nED);
1720 snew(EDstate->old_sav, EDstate->nED);
1723 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1724 for (int i = 0; i < EDstate->nED; i++)
1728 /* Reference structure SREF */
1729 sprintf(buf, "ED%d # of atoms in reference structure", i+1);
1730 do_cpt_int_err(xd, buf, &EDstate->nref[i], list);
1731 sprintf(buf, "ED%d x_ref", i+1);
1734 snew(EDstate->old_sref[i], EDstate->nref[i]);
1735 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref[i], list);
1739 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref_p[i], list);
1742 /* Average structure SAV */
1743 sprintf(buf, "ED%d # of atoms in average structure", i+1);
1744 do_cpt_int_err(xd, buf, &EDstate->nav[i], list);
1745 sprintf(buf, "ED%d x_av", i+1);
1748 snew(EDstate->old_sav[i], EDstate->nav[i]);
1749 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav[i], list);
1753 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav_p[i], list);
1760 static int do_cpt_correlation_grid(XDR *xd, gmx_bool bRead, gmx_unused int fflags,
1761 gmx::CorrelationGridHistory *corrGrid,
1762 FILE *list, int eawhh)
1766 do_cpt_int_err(xd, eawhh_names[eawhh], &(corrGrid->numCorrelationTensors), list);
1767 do_cpt_int_err(xd, eawhh_names[eawhh], &(corrGrid->tensorSize), list);
1768 do_cpt_int_err(xd, eawhh_names[eawhh], &(corrGrid->blockDataListSize), list);
1772 initCorrelationGridHistory(corrGrid, corrGrid->numCorrelationTensors, corrGrid->tensorSize, corrGrid->blockDataListSize);
1775 for (gmx::CorrelationBlockDataHistory &blockData : corrGrid->blockDataBuffer)
1777 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.blockSumWeight), list);
1778 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.blockSumSquareWeight), list);
1779 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.blockSumWeightX), list);
1780 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.blockSumWeightY), list);
1781 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.sumOverBlocksSquareBlockWeight), list);
1782 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.sumOverBlocksBlockSquareWeight), list);
1783 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.sumOverBlocksBlockWeightBlockWeightX), list);
1784 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.sumOverBlocksBlockWeightBlockWeightY), list);
1785 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.blockLength), list);
1786 do_cpt_int_err(xd, eawhh_names[eawhh], &(blockData.previousBlockIndex), list);
1787 do_cpt_double_err(xd, eawhh_names[eawhh], &(blockData.correlationIntegral), list);
1793 static int do_cpt_awh_bias(XDR *xd, gmx_bool bRead,
1794 int fflags, gmx::AwhBiasHistory *biasHistory,
1799 gmx::AwhBiasStateHistory *state = &biasHistory->state;
1800 for (int i = 0; (i < eawhhNR && ret == 0); i++)
1802 if (fflags & (1<<i))
1806 case eawhhIN_INITIAL:
1807 do_cpt_bool_err(xd, eawhh_names[i], &state->in_initial, list); break;
1808 case eawhhEQUILIBRATEHISTOGRAM:
1809 do_cpt_bool_err(xd, eawhh_names[i], &state->equilibrateHistogram, list); break;
1811 do_cpt_double_err(xd, eawhh_names[i], &state->histSize, list); break;
1817 numPoints = biasHistory->pointState.size();
1819 do_cpt_int_err(xd, eawhh_names[i], &numPoints, list);
1822 biasHistory->pointState.resize(numPoints);
1826 case eawhhCOORDPOINT:
1827 for (auto &psh : biasHistory->pointState)
1829 do_cpt_double_err(xd, eawhh_names[i], &psh.target, list);
1830 do_cpt_double_err(xd, eawhh_names[i], &psh.free_energy, list);
1831 do_cpt_double_err(xd, eawhh_names[i], &psh.bias, list);
1832 do_cpt_double_err(xd, eawhh_names[i], &psh.weightsum_iteration, list);
1833 do_cpt_double_err(xd, eawhh_names[i], &psh.weightsum_covering, list);
1834 do_cpt_double_err(xd, eawhh_names[i], &psh.weightsum_tot, list);
1835 do_cpt_double_err(xd, eawhh_names[i], &psh.weightsum_ref, list);
1836 do_cpt_step_err(xd, eawhh_names[i], &psh.last_update_index, list);
1837 do_cpt_double_err(xd, eawhh_names[i], &psh.log_pmfsum, list);
1838 do_cpt_double_err(xd, eawhh_names[i], &psh.visits_iteration, list);
1839 do_cpt_double_err(xd, eawhh_names[i], &psh.visits_tot, list);
1842 case eawhhUMBRELLAGRIDPOINT:
1843 do_cpt_int_err(xd, eawhh_names[i], &(state->umbrellaGridpoint), list); break;
1844 case eawhhUPDATELIST:
1845 do_cpt_int_err(xd, eawhh_names[i], &(state->origin_index_updatelist), list);
1846 do_cpt_int_err(xd, eawhh_names[i], &(state->end_index_updatelist), list);
1848 case eawhhLOGSCALEDSAMPLEWEIGHT:
1849 do_cpt_double_err(xd, eawhh_names[i], &(state->logScaledSampleWeight), list);
1850 do_cpt_double_err(xd, eawhh_names[i], &(state->maxLogScaledSampleWeight), list);
1852 case eawhhNUMUPDATES:
1853 do_cpt_step_err(xd, eawhh_names[i], &(state->numUpdates), list);
1855 case eawhhFORCECORRELATIONGRID:
1856 ret = do_cpt_correlation_grid(xd, bRead, fflags, &biasHistory->forceCorrelationGrid, list, i);
1859 gmx_fatal(FARGS, "Unknown awh history entry %d\n", i);
1867 static int do_cpt_awh(XDR *xd, gmx_bool bRead,
1868 int fflags, gmx::AwhHistory *awhHistory,
1875 std::shared_ptr<gmx::AwhHistory> awhHistoryLocal;
1877 if (awhHistory == nullptr)
1879 GMX_RELEASE_ASSERT(bRead, "do_cpt_awh should not be called for writing without an AwhHistory");
1881 awhHistoryLocal = std::make_shared<gmx::AwhHistory>();
1882 awhHistory = awhHistoryLocal.get();
1885 /* To be able to read and write the AWH data several parameters determining the layout of the AWH structs need to be known
1886 (nbias, npoints, etc.). The best thing (?) would be to have these passed to this function. When writing to a checkpoint
1887 these parameters are available in awh_history (after calling init_awh_history). When reading from a checkpoint though, there
1888 is no initialized awh_history (it is initialized and set in this function). The AWH parameters have not always been read
1889 at the time when this function is called for reading so I don't know how to pass them as input. Here, this is solved by
1890 when writing a checkpoint, also storing parameters needed for future reading of the checkpoint.
1892 Another issue is that some variables that AWH checkpoints don't have a registered enum and string (e.g. nbias below).
1893 One difficulty is the multilevel structure of the data which would need to be represented somehow. */
1898 numBias = awhHistory->bias.size();
1900 do_cpt_int_err(xd, "awh_nbias", &numBias, list);
1904 awhHistory->bias.resize(numBias);
1906 for (auto &bias : awhHistory->bias)
1908 ret = do_cpt_awh_bias(xd, bRead, fflags, &bias, list);
1914 do_cpt_double_err(xd, "awh_potential_offset", &awhHistory->potentialOffset, list);
1919 static int do_cpt_files(XDR *xd, gmx_bool bRead,
1920 std::vector<gmx_file_position_t> *outputfiles,
1921 FILE *list, int file_version)
1924 gmx_off_t mask = 0xFFFFFFFFL;
1925 int offset_high, offset_low;
1926 std::array<char, CPTSTRLEN> buf;
1927 GMX_RELEASE_ASSERT(outputfiles, "Must have valid outputfiles");
1929 // Ensure that reading pre-allocates outputfiles, while writing
1930 // writes what is already there.
1931 int nfiles = outputfiles->size();
1932 if (do_cpt_int(xd, "number of output files", &nfiles, list) != 0)
1938 outputfiles->resize(nfiles);
1941 for (auto &outputfile : *outputfiles)
1943 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1946 do_cpt_string_err(xd, "output filename", buf, list);
1947 std::copy(std::begin(buf), std::end(buf), std::begin(outputfile.filename));
1949 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1953 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1957 outputfile.offset = (static_cast<gmx_off_t>(offset_high) << 32 ) | ( static_cast<gmx_off_t>(offset_low) & mask );
1961 do_cpt_string_err(xd, "output filename", outputfile.filename, list);
1963 offset = outputfile.offset;
1971 offset_low = static_cast<int>(offset & mask);
1972 offset_high = static_cast<int>((offset >> 32) & mask);
1974 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1978 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1983 if (file_version >= 8)
1985 if (do_cpt_int(xd, "file_checksum_size", &outputfile.checksumSize,
1990 if (do_cpt_u_chars(xd, "file_checksum", outputfile.checksum.size(), outputfile.checksum.data(), list) != 0)
1997 outputfile.checksumSize = -1;
2004 void write_checkpoint(const char *fn, gmx_bool bNumberAndKeep,
2005 FILE *fplog, const t_commrec *cr,
2006 ivec domdecCells, int nppnodes,
2007 int eIntegrator, int simulation_part,
2008 gmx_bool bExpanded, int elamstats,
2009 int64_t step, double t,
2010 t_state *state, ObservablesHistory *observablesHistory)
2013 char *fntemp; /* the temporary checkpoint file name */
2015 char buf[1024], suffix[5+STEPSTRSIZE], sbuf[STEPSTRSIZE];
2018 if (DOMAINDECOMP(cr))
2020 npmenodes = cr->npmenodes;
2028 /* make the new temporary filename */
2029 snew(fntemp, std::strlen(fn)+5+STEPSTRSIZE);
2030 std::strcpy(fntemp, fn);
2031 fntemp[std::strlen(fn) - std::strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
2032 sprintf(suffix, "_%s%s", "step", gmx_step_str(step, sbuf));
2033 std::strcat(fntemp, suffix);
2034 std::strcat(fntemp, fn+std::strlen(fn) - std::strlen(ftp2ext(fn2ftp(fn))) - 1);
2036 /* if we can't rename, we just overwrite the cpt file.
2037 * dangerous if interrupted.
2039 snew(fntemp, std::strlen(fn));
2040 std::strcpy(fntemp, fn);
2042 std::string timebuf = gmx_format_current_time();
2046 fprintf(fplog, "Writing checkpoint, step %s at %s\n\n",
2047 gmx_step_str(step, buf), timebuf.c_str());
2050 /* Get offsets for open files */
2051 auto outputfiles = gmx_fio_get_output_file_positions();
2053 fp = gmx_fio_open(fntemp, "w");
2056 if (state->ekinstate.bUpToDate)
2059 ((1<<eeksEKIN_N) | (1<<eeksEKINH) | (1<<eeksEKINF) |
2060 (1<<eeksEKINO) | (1<<eeksEKINSCALEF) | (1<<eeksEKINSCALEH) |
2061 (1<<eeksVSCALE) | (1<<eeksDEKINDL) | (1<<eeksMVCOS));
2068 energyhistory_t *enerhist = observablesHistory->energyHistory.get();
2070 if (enerhist != nullptr && (enerhist->nsum > 0 || enerhist->nsum_sim > 0))
2072 flags_enh |= (1<<eenhENERGY_N) | (1<<eenhENERGY_NSTEPS) | (1<<eenhENERGY_NSTEPS_SIM);
2073 if (enerhist->nsum > 0)
2075 flags_enh |= ((1<<eenhENERGY_AVER) | (1<<eenhENERGY_SUM) |
2076 (1<<eenhENERGY_NSUM));
2078 if (enerhist->nsum_sim > 0)
2080 flags_enh |= ((1<<eenhENERGY_SUM_SIM) | (1<<eenhENERGY_NSUM_SIM));
2082 if (enerhist->deltaHForeignLambdas != nullptr)
2084 flags_enh |= ( (1<< eenhENERGY_DELTA_H_NN) |
2085 (1<< eenhENERGY_DELTA_H_LIST) |
2086 (1<< eenhENERGY_DELTA_H_STARTTIME) |
2087 (1<< eenhENERGY_DELTA_H_STARTLAMBDA) );
2091 PullHistory *pullHist = observablesHistory->pullHistory.get();
2092 int flagsPullHistory = 0;
2093 if (pullHist != nullptr && (pullHist->numValuesInXSum > 0 || pullHist->numValuesInFSum > 0))
2095 flagsPullHistory |= (1<<epullhPULL_NUMCOORDINATES);
2096 flagsPullHistory |= ((1<<epullhPULL_NUMGROUPS) | (1<<epullhPULL_NUMVALUESINXSUM) | (1<<epullhPULL_NUMVALUESINFSUM));
2102 flags_dfh = ((1<<edfhBEQUIL) | (1<<edfhNATLAMBDA) | (1<<edfhSUMWEIGHTS) | (1<<edfhSUMDG) |
2103 (1<<edfhTIJ) | (1<<edfhTIJEMP));
2106 flags_dfh |= ((1<<edfhWLDELTA) | (1<<edfhWLHISTO));
2108 if ((elamstats == elamstatsMINVAR) || (elamstats == elamstatsBARKER) || (elamstats == elamstatsMETROPOLIS))
2110 flags_dfh |= ((1<<edfhACCUMP) | (1<<edfhACCUMM) | (1<<edfhACCUMP2) | (1<<edfhACCUMM2)
2111 | (1<<edfhSUMMINVAR) | (1<<edfhSUMVAR));
2120 if (state->awhHistory != nullptr && !state->awhHistory->bias.empty())
2122 flags_awhh |= ( (1 << eawhhIN_INITIAL) |
2123 (1 << eawhhEQUILIBRATEHISTOGRAM) |
2124 (1 << eawhhHISTSIZE) |
2125 (1 << eawhhNPOINTS) |
2126 (1 << eawhhCOORDPOINT) |
2127 (1 << eawhhUMBRELLAGRIDPOINT) |
2128 (1 << eawhhUPDATELIST) |
2129 (1 << eawhhLOGSCALEDSAMPLEWEIGHT) |
2130 (1 << eawhhNUMUPDATES) |
2131 (1 << eawhhFORCECORRELATIONGRID));
2134 /* We can check many more things now (CPU, acceleration, etc), but
2135 * it is highly unlikely to have two separate builds with exactly
2136 * the same version, user, time, and build host!
2139 int nlambda = (state->dfhist ? state->dfhist->nlambda : 0);
2141 edsamhistory_t *edsamhist = observablesHistory->edsamHistory.get();
2142 int nED = (edsamhist ? edsamhist->nED : 0);
2144 swaphistory_t *swaphist = observablesHistory->swapHistory.get();
2145 int eSwapCoords = (swaphist ? swaphist->eSwapCoords : eswapNO);
2147 CheckpointHeaderContents headerContents =
2149 0, {0}, {0}, {0}, {0}, GMX_DOUBLE, {0}, {0},
2150 eIntegrator, simulation_part, step, t, nppnodes,
2152 state->natoms, state->ngtc, state->nnhpres,
2153 state->nhchainlength, nlambda, state->flags, flags_eks, flags_enh,
2154 flagsPullHistory, flags_dfh, flags_awhh,
2157 std::strcpy(headerContents.version, gmx_version());
2158 std::strcpy(headerContents.fprog, gmx::getProgramContext().fullBinaryPath());
2159 std::strcpy(headerContents.ftime, timebuf.c_str());
2160 if (DOMAINDECOMP(cr))
2162 copy_ivec(domdecCells, headerContents.dd_nc);
2165 do_cpt_header(gmx_fio_getxdr(fp), FALSE, nullptr, &headerContents);
2167 if ((do_cpt_state(gmx_fio_getxdr(fp), state->flags, state, nullptr) < 0) ||
2168 (do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, nullptr) < 0) ||
2169 (do_cpt_enerhist(gmx_fio_getxdr(fp), FALSE, flags_enh, enerhist, nullptr) < 0) ||
2170 (doCptPullHist(gmx_fio_getxdr(fp), FALSE, flagsPullHistory, pullHist, StatePart::pullHistory, nullptr) < 0) ||
2171 (do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, nlambda, &state->dfhist, nullptr) < 0) ||
2172 (do_cpt_EDstate(gmx_fio_getxdr(fp), FALSE, nED, edsamhist, nullptr) < 0) ||
2173 (do_cpt_awh(gmx_fio_getxdr(fp), FALSE, flags_awhh, state->awhHistory.get(), nullptr) < 0) ||
2174 (do_cpt_swapstate(gmx_fio_getxdr(fp), FALSE, eSwapCoords, swaphist, nullptr) < 0) ||
2175 (do_cpt_files(gmx_fio_getxdr(fp), FALSE, &outputfiles, nullptr,
2176 headerContents.file_version) < 0))
2178 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2181 do_cpt_footer(gmx_fio_getxdr(fp), headerContents.file_version);
2183 /* we really, REALLY, want to make sure to physically write the checkpoint,
2184 and all the files it depends on, out to disk. Because we've
2185 opened the checkpoint with gmx_fio_open(), it's in our list
2187 ret = gmx_fio_all_output_fsync();
2193 "Cannot fsync '%s'; maybe you are out of disk space?",
2194 gmx_fio_getname(ret));
2196 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV) == nullptr)
2202 gmx_warning("%s", buf);
2206 if (gmx_fio_close(fp) != 0)
2208 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2211 /* we don't move the checkpoint if the user specified they didn't want it,
2212 or if the fsyncs failed */
2214 if (!bNumberAndKeep && !ret)
2218 /* Rename the previous checkpoint file */
2219 std::strcpy(buf, fn);
2220 buf[std::strlen(fn) - std::strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
2221 std::strcat(buf, "_prev");
2222 std::strcat(buf, fn+std::strlen(fn) - std::strlen(ftp2ext(fn2ftp(fn))) - 1);
2225 /* we copy here so that if something goes wrong between now and
2226 * the rename below, there's always a state.cpt.
2227 * If renames are atomic (such as in POSIX systems),
2228 * this copying should be unneccesary.
2230 gmx_file_copy(fn, buf, FALSE);
2231 /* We don't really care if this fails:
2232 * there's already a new checkpoint.
2237 gmx_file_rename(fn, buf);
2240 if (gmx_file_rename(fntemp, fn) != 0)
2242 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
2245 #endif /* GMX_NO_RENAME */
2250 /*code for alternate checkpointing scheme. moved from top of loop over
2252 fcRequestCheckPoint();
2253 if (fcCheckPointParallel( cr->nodeid, NULL, 0) == 0)
2255 gmx_fatal( 3, __FILE__, __LINE__, "Checkpoint error on step %d\n", step );
2257 #endif /* end GMX_FAHCORE block */
2260 static void check_int(FILE *fplog, const char *type, int p, int f, gmx_bool *mm)
2262 bool foundMismatch = (p != f);
2270 fprintf(fplog, " %s mismatch,\n", type);
2271 fprintf(fplog, " current program: %d\n", p);
2272 fprintf(fplog, " checkpoint file: %d\n", f);
2273 fprintf(fplog, "\n");
2277 static void check_string(FILE *fplog, const char *type, const char *p,
2278 const char *f, gmx_bool *mm)
2280 bool foundMismatch = (std::strcmp(p, f) != 0);
2288 fprintf(fplog, " %s mismatch,\n", type);
2289 fprintf(fplog, " current program: %s\n", p);
2290 fprintf(fplog, " checkpoint file: %s\n", f);
2291 fprintf(fplog, "\n");
2295 static void check_match(FILE *fplog, const t_commrec *cr, const ivec dd_nc,
2296 const CheckpointHeaderContents &headerContents,
2297 gmx_bool reproducibilityRequested)
2299 /* Note that this check_string on the version will also print a message
2300 * when only the minor version differs. But we only print a warning
2301 * message further down with reproducibilityRequested=TRUE.
2303 gmx_bool versionDiffers = FALSE;
2304 check_string(fplog, "Version", gmx_version(), headerContents.version, &versionDiffers);
2306 gmx_bool precisionDiffers = FALSE;
2307 check_int (fplog, "Double prec.", GMX_DOUBLE, headerContents.double_prec, &precisionDiffers);
2308 if (precisionDiffers)
2310 const char msg_precision_difference[] =
2311 "You are continuing a simulation with a different precision. Not matching\n"
2312 "single/double precision will lead to precision or performance loss.\n";
2315 fprintf(fplog, "%s\n", msg_precision_difference);
2319 gmx_bool mm = (versionDiffers || precisionDiffers);
2321 if (reproducibilityRequested)
2323 check_string(fplog, "Program name", gmx::getProgramContext().fullBinaryPath(), headerContents.fprog, &mm);
2325 check_int (fplog, "#ranks", cr->nnodes, headerContents.nnodes, &mm);
2328 if (cr->nnodes > 1 && reproducibilityRequested)
2330 check_int (fplog, "#PME-ranks", cr->npmenodes, headerContents.npme, &mm);
2332 int npp = cr->nnodes;
2333 if (cr->npmenodes >= 0)
2335 npp -= cr->npmenodes;
2337 int npp_f = headerContents.nnodes;
2338 if (headerContents.npme >= 0)
2340 npp_f -= headerContents.npme;
2344 check_int (fplog, "#DD-cells[x]", dd_nc[XX], headerContents.dd_nc[XX], &mm);
2345 check_int (fplog, "#DD-cells[y]", dd_nc[YY], headerContents.dd_nc[YY], &mm);
2346 check_int (fplog, "#DD-cells[z]", dd_nc[ZZ], headerContents.dd_nc[ZZ], &mm);
2352 /* Gromacs should be able to continue from checkpoints between
2353 * different patch level versions, but we do not guarantee
2354 * compatibility between different major/minor versions - check this.
2358 sscanf(gmx_version(), "%5d", &gmx_major);
2359 int ret = sscanf(headerContents.version, "%5d", &cpt_major);
2360 gmx_bool majorVersionDiffers = (ret < 1 || gmx_major != cpt_major);
2362 const char msg_major_version_difference[] =
2363 "The current GROMACS major version is not identical to the one that\n"
2364 "generated the checkpoint file. In principle GROMACS does not support\n"
2365 "continuation from checkpoints between different versions, so we advise\n"
2366 "against this. If you still want to try your luck we recommend that you use\n"
2367 "the -noappend flag to keep your output files from the two versions separate.\n"
2368 "This might also work around errors where the output fields in the energy\n"
2369 "file have changed between the different versions.\n";
2371 const char msg_mismatch_notice[] =
2372 "GROMACS patchlevel, binary or parallel settings differ from previous run.\n"
2373 "Continuation is exact, but not guaranteed to be binary identical.\n";
2375 if (majorVersionDiffers)
2379 fprintf(fplog, "%s\n", msg_major_version_difference);
2382 else if (reproducibilityRequested)
2384 /* Major & minor versions match at least, but something is different. */
2387 fprintf(fplog, "%s\n", msg_mismatch_notice);
2393 static void lockLogFile(FILE *fplog,
2394 t_fileio *chksum_file,
2395 const char *filename,
2398 /* Note that there are systems where the lock operation
2399 * will succeed, but a second process can also lock the file.
2400 * We should probably try to detect this.
2402 #if defined __native_client__
2405 #elif GMX_NATIVE_WINDOWS
2406 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX) == -1)
2408 struct flock fl; /* don't initialize here: the struct order is OS
2410 fl.l_type = F_WRLCK;
2411 fl.l_whence = SEEK_SET;
2416 if (fcntl(fileno(gmx_fio_getfp(chksum_file)), F_SETLK, &fl) == -1)
2419 if (errno == ENOSYS)
2423 gmx_fatal(FARGS, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2429 fprintf(fplog, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", filename);
2433 else if (errno == EACCES || errno == EAGAIN)
2435 gmx_fatal(FARGS, "Failed to lock: %s. Already running "
2436 "simulation?", filename);
2440 gmx_fatal(FARGS, "Failed to lock: %s. %s.",
2441 filename, std::strerror(errno));
2446 //! Check whether chksum_file output file has a checksum that matches \c outputfile from the checkpoint.
2447 static void checkOutputFile(t_fileio *chksum_file,
2448 const gmx_file_position_t &outputfile)
2450 /* compute md5 chksum */
2451 std::array<unsigned char, 16> digest;
2452 if (outputfile.checksumSize != -1)
2454 if (gmx_fio_get_file_md5(chksum_file, outputfile.offset,
2455 &digest) != outputfile.checksumSize) /*at the end of the call the file position is at the end of the file*/
2457 gmx_fatal(FARGS, "Can't read %d bytes of '%s' to compute checksum. The file has been replaced or its contents have been modified. Cannot do appending because of this condition.",
2458 outputfile.checksumSize,
2459 outputfile.filename);
2463 /* compare md5 chksum */
2464 if (outputfile.checksumSize != -1 &&
2465 digest != outputfile.checksum)
2469 fprintf(debug, "chksum for %s: ", outputfile.filename);
2470 for (int j = 0; j < 16; j++)
2472 fprintf(debug, "%02x", digest[j]);
2474 fprintf(debug, "\n");
2476 gmx_fatal(FARGS, "Checksum wrong for '%s'. The file has been replaced or its contents have been modified. Cannot do appending because of this condition.",
2477 outputfile.filename);
2481 static void read_checkpoint(const char *fn, t_fileio *logfio,
2482 const t_commrec *cr,
2485 int *init_fep_state,
2486 CheckpointHeaderContents *headerContents,
2487 t_state *state, gmx_bool *bReadEkin,
2488 ObservablesHistory *observablesHistory,
2489 gmx_bool bAppendOutputFiles, gmx_bool bForceAppend,
2490 gmx_bool reproducibilityRequested)
2494 char buf[STEPSTRSIZE];
2497 fp = gmx_fio_open(fn, "r");
2498 do_cpt_header(gmx_fio_getxdr(fp), TRUE, nullptr, headerContents);
2500 if (bAppendOutputFiles &&
2501 headerContents->file_version >= 13 && headerContents->double_prec != GMX_DOUBLE)
2503 gmx_fatal(FARGS, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
2506 // If we are appending, then we don't write to the open log file
2507 // because we still need to compute a checksum for it. Otherwise
2508 // we report to the new log file about the checkpoint file that we
2509 // are reading from.
2510 FILE *fplog = bAppendOutputFiles ? nullptr : gmx_fio_getfp(logfio);
2513 fprintf(fplog, "\n");
2514 fprintf(fplog, "Reading checkpoint file %s\n", fn);
2515 fprintf(fplog, " file generated by: %s\n", headerContents->fprog);
2516 fprintf(fplog, " file generated at: %s\n", headerContents->ftime);
2517 fprintf(fplog, " GROMACS double prec.: %d\n", headerContents->double_prec);
2518 fprintf(fplog, " simulation part #: %d\n", headerContents->simulation_part);
2519 fprintf(fplog, " step: %s\n", gmx_step_str(headerContents->step, buf));
2520 fprintf(fplog, " time: %f\n", headerContents->t);
2521 fprintf(fplog, "\n");
2524 if (headerContents->natoms != state->natoms)
2526 gmx_fatal(FARGS, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", headerContents->natoms, state->natoms);
2528 if (headerContents->ngtc != state->ngtc)
2530 gmx_fatal(FARGS, "Checkpoint file is for a system of %d T-coupling groups, while the current system consists of %d T-coupling groups", headerContents->ngtc, state->ngtc);
2532 if (headerContents->nnhpres != state->nnhpres)
2534 gmx_fatal(FARGS, "Checkpoint file is for a system of %d NH-pressure-coupling variables, while the current system consists of %d NH-pressure-coupling variables", headerContents->nnhpres, state->nnhpres);
2537 int nlambdaHistory = (state->dfhist ? state->dfhist->nlambda : 0);
2538 if (headerContents->nlambda != nlambdaHistory)
2540 gmx_fatal(FARGS, "Checkpoint file is for a system with %d lambda states, while the current system consists of %d lambda states", headerContents->nlambda, nlambdaHistory);
2543 init_gtc_state(state, state->ngtc, state->nnhpres, headerContents->nhchainlength); /* need to keep this here to keep the tpr format working */
2544 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
2546 if (headerContents->eIntegrator != eIntegrator)
2548 gmx_fatal(FARGS, "Cannot change integrator during a checkpoint restart. Perhaps you should make a new .tpr with grompp -f new.mdp -t %s", fn);
2551 if (headerContents->flags_state != state->flags)
2553 gmx_fatal(FARGS, "Cannot change a simulation algorithm during a checkpoint restart. Perhaps you should make a new .tpr with grompp -f new.mdp -t %s", fn);
2558 check_match(fplog, cr, dd_nc, *headerContents,
2559 reproducibilityRequested);
2562 ret = do_cpt_state(gmx_fio_getxdr(fp), headerContents->flags_state, state, nullptr);
2563 *init_fep_state = state->fep_state; /* there should be a better way to do this than setting it here.
2564 Investigate for 5.0. */
2569 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), headerContents->flags_eks, &state->ekinstate, nullptr);
2574 *bReadEkin = (((headerContents->flags_eks & (1<<eeksEKINH)) != 0) ||
2575 ((headerContents->flags_eks & (1<<eeksEKINF)) != 0) ||
2576 ((headerContents->flags_eks & (1<<eeksEKINO)) != 0) ||
2577 (((headerContents->flags_eks & (1<<eeksEKINSCALEF)) |
2578 (headerContents->flags_eks & (1<<eeksEKINSCALEH)) |
2579 (headerContents->flags_eks & (1<<eeksVSCALE))) != 0));
2581 if (headerContents->flags_enh && observablesHistory->energyHistory == nullptr)
2583 observablesHistory->energyHistory = gmx::compat::make_unique<energyhistory_t>();
2585 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2586 headerContents->flags_enh, observablesHistory->energyHistory.get(), nullptr);
2592 if (headerContents->flagsPullHistory)
2594 if (observablesHistory->pullHistory == nullptr)
2596 observablesHistory->pullHistory = gmx::compat::make_unique<PullHistory>();
2598 ret = doCptPullHist(gmx_fio_getxdr(fp), TRUE,
2599 headerContents->flagsPullHistory, observablesHistory->pullHistory.get(), StatePart::pullHistory, nullptr);
2606 if (headerContents->file_version < 6)
2608 gmx_fatal(FARGS, "Continuing from checkpoint files written before GROMACS 4.5 is not supported");
2611 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), headerContents->flags_dfh, headerContents->nlambda, &state->dfhist, nullptr);
2617 if (headerContents->nED > 0 && observablesHistory->edsamHistory == nullptr)
2619 observablesHistory->edsamHistory = gmx::compat::make_unique<edsamhistory_t>(edsamhistory_t {});
2621 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, headerContents->nED, observablesHistory->edsamHistory.get(), nullptr);
2627 if (headerContents->flags_awhh != 0 && state->awhHistory == nullptr)
2629 state->awhHistory = std::make_shared<gmx::AwhHistory>();
2631 ret = do_cpt_awh(gmx_fio_getxdr(fp), TRUE,
2632 headerContents->flags_awhh, state->awhHistory.get(), nullptr);
2638 if (headerContents->eSwapCoords != eswapNO && observablesHistory->swapHistory == nullptr)
2640 observablesHistory->swapHistory = gmx::compat::make_unique<swaphistory_t>(swaphistory_t {});
2642 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, headerContents->eSwapCoords, observablesHistory->swapHistory.get(), nullptr);
2648 std::vector<gmx_file_position_t> outputfiles;
2649 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, nullptr, headerContents->file_version);
2655 ret = do_cpt_footer(gmx_fio_getxdr(fp), headerContents->file_version);
2660 if (gmx_fio_close(fp) != 0)
2662 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2665 /* If the user wants to append to output files,
2666 * we use the file pointer positions of the output files stored
2667 * in the checkpoint file and truncate the files such that any frames
2668 * written after the checkpoint time are removed.
2669 * All files are md5sum checked such that we can be sure that
2670 * we do not truncate other (maybe imprortant) files.
2672 if (bAppendOutputFiles)
2674 if (outputfiles.empty())
2676 gmx_fatal(FARGS, "No names of output files were recorded in the checkpoint");
2678 if (fn2ftp(outputfiles[0].filename) != efLOG)
2680 /* make sure first file is log file so that it is OK to use it for
2683 gmx_fatal(FARGS, "The first output file should always be the log "
2684 "file but instead is: %s. Cannot do appending because of this condition.", outputfiles[0].filename);
2686 for (const auto &outputfile : outputfiles)
2688 if (outputfile.offset < 0)
2690 gmx_fatal(FARGS, "The original run wrote a file called '%s' which "
2691 "is larger than 2 GB, but mdrun did not support large file"
2692 " offsets. Can not append. Run mdrun with -noappend",
2693 outputfile.filename);
2697 // Handle the log file separately - it comes first in the
2698 // list, so we make sure that we retain a lock on the open
2699 // file that is never lifted after the checksum is calculated.
2701 const gmx_file_position_t &logOutputFile = outputfiles[0];
2704 lockLogFile(fplog, logfio, logOutputFile.filename, bForceAppend);
2705 checkOutputFile(logfio, logOutputFile);
2707 if (gmx_fio_seek(logfio, logOutputFile.offset))
2709 gmx_fatal(FARGS, "Seek error! Failed to truncate log-file: %s.", std::strerror(errno));
2715 // Now handle the remaining outputfiles
2716 gmx::ArrayRef<gmx_file_position_t> outputfilesView(outputfiles);
2717 for (const auto &outputfile : outputfilesView.subArray(1, outputfiles.size()-1))
2719 t_fileio *chksum_file = gmx_fio_open(outputfile.filename, "r+");
2720 checkOutputFile(chksum_file, outputfile);
2721 gmx_fio_close(chksum_file);
2723 if (!GMX_NATIVE_WINDOWS)
2725 /* For FAHCORE, we do this elsewhere*/
2726 rc = gmx_truncate(outputfile.filename, outputfile.offset);
2729 gmx_fatal(FARGS, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfile.filename);
2738 void load_checkpoint(const char *fn, t_fileio *logfio,
2739 const t_commrec *cr, const ivec dd_nc,
2740 t_inputrec *ir, t_state *state,
2741 gmx_bool *bReadEkin,
2742 ObservablesHistory *observablesHistory,
2743 gmx_bool bAppend, gmx_bool bForceAppend,
2744 gmx_bool reproducibilityRequested)
2746 CheckpointHeaderContents headerContents;
2749 /* Read the state from the checkpoint file */
2750 read_checkpoint(fn, logfio,
2752 ir->eI, &(ir->fepvals->init_fep_state),
2754 state, bReadEkin, observablesHistory,
2755 bAppend, bForceAppend,
2756 reproducibilityRequested);
2760 gmx_bcast(sizeof(headerContents.step), &headerContents.step, cr);
2761 gmx_bcast(sizeof(*bReadEkin), bReadEkin, cr);
2763 ir->bContinuation = TRUE;
2764 // TODO Should the following condition be <=? Currently if you
2765 // pass a checkpoint written by an normal completion to a restart,
2766 // mdrun will read all input, does some work but no steps, and
2767 // write successful output. But perhaps that is not desirable.
2768 if ((ir->nsteps >= 0) && (ir->nsteps < headerContents.step))
2770 // Note that we do not intend to support the use of mdrun
2771 // -nsteps to circumvent this condition.
2772 char nstepsString[STEPSTRSIZE], stepString[STEPSTRSIZE];
2773 gmx_step_str(ir->nsteps, nstepsString);
2774 gmx_step_str(headerContents.step, stepString);
2775 gmx_fatal(FARGS, "The input requested %s steps, however the checkpoint "
2776 "file has already reached step %s. The simulation will not "
2777 "proceed, because either your simulation is already complete, "
2778 "or your combination of input files don't match.",
2779 nstepsString, stepString);
2781 if (ir->nsteps >= 0)
2783 ir->nsteps += ir->init_step - headerContents.step;
2785 ir->init_step = headerContents.step;
2786 ir->simulation_part = headerContents.simulation_part + 1;
2789 void read_checkpoint_part_and_step(const char *filename,
2790 int *simulation_part,
2795 if (filename == nullptr ||
2796 !gmx_fexist(filename) ||
2797 ((fp = gmx_fio_open(filename, "r")) == nullptr))
2799 *simulation_part = 0;
2804 CheckpointHeaderContents headerContents;
2805 do_cpt_header(gmx_fio_getxdr(fp), TRUE, nullptr, &headerContents);
2807 *simulation_part = headerContents.simulation_part;
2808 *step = headerContents.step;
2811 static void read_checkpoint_data(t_fileio *fp, int *simulation_part,
2812 int64_t *step, double *t, t_state *state,
2813 std::vector<gmx_file_position_t> *outputfiles)
2817 CheckpointHeaderContents headerContents;
2818 do_cpt_header(gmx_fio_getxdr(fp), TRUE, nullptr, &headerContents);
2819 *simulation_part = headerContents.simulation_part;
2820 *step = headerContents.step;
2821 *t = headerContents.t;
2822 state->natoms = headerContents.natoms;
2823 state->ngtc = headerContents.ngtc;
2824 state->nnhpres = headerContents.nnhpres;
2825 state->nhchainlength = headerContents.nhchainlength;
2826 state->flags = headerContents.flags_state;
2828 do_cpt_state(gmx_fio_getxdr(fp), state->flags, state, nullptr);
2833 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), headerContents.flags_eks, &state->ekinstate, nullptr);
2839 energyhistory_t enerhist;
2840 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2841 headerContents.flags_enh, &enerhist, nullptr);
2846 PullHistory pullHist = {};
2847 ret = doCptPullHist(gmx_fio_getxdr(fp), TRUE,
2848 headerContents.flagsPullHistory, &pullHist, StatePart::pullHistory, nullptr);
2854 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), headerContents.flags_dfh, headerContents.nlambda, &state->dfhist, nullptr);
2860 edsamhistory_t edsamhist = {};
2861 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, headerContents.nED, &edsamhist, nullptr);
2867 ret = do_cpt_awh(gmx_fio_getxdr(fp), TRUE,
2868 headerContents.flags_awhh, state->awhHistory.get(), nullptr);
2874 swaphistory_t swaphist = {};
2875 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, headerContents.eSwapCoords, &swaphist, nullptr);
2881 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE,
2883 nullptr, headerContents.file_version);
2890 ret = do_cpt_footer(gmx_fio_getxdr(fp), headerContents.file_version);
2897 void read_checkpoint_trxframe(t_fileio *fp, t_trxframe *fr)
2900 int simulation_part;
2904 std::vector<gmx_file_position_t> outputfiles;
2905 read_checkpoint_data(fp, &simulation_part, &step, &t, &state, &outputfiles);
2907 fr->natoms = state.natoms;
2909 fr->step = int64_to_int(step,
2910 "conversion of checkpoint to trajectory");
2914 fr->lambda = state.lambda[efptFEP];
2915 fr->fep_state = state.fep_state;
2917 fr->bX = ((state.flags & (1<<estX)) != 0);
2920 fr->x = makeRvecArray(state.x, state.natoms);
2922 fr->bV = ((state.flags & (1<<estV)) != 0);
2925 fr->v = makeRvecArray(state.v, state.natoms);
2928 fr->bBox = ((state.flags & (1<<estBOX)) != 0);
2931 copy_mat(state.box, fr->box);
2935 void list_checkpoint(const char *fn, FILE *out)
2942 fp = gmx_fio_open(fn, "r");
2943 CheckpointHeaderContents headerContents;
2944 do_cpt_header(gmx_fio_getxdr(fp), TRUE, out, &headerContents);
2945 state.natoms = headerContents.natoms;
2946 state.ngtc = headerContents.ngtc;
2947 state.nnhpres = headerContents.nnhpres;
2948 state.nhchainlength = headerContents.nhchainlength;
2949 state.flags = headerContents.flags_state;
2950 ret = do_cpt_state(gmx_fio_getxdr(fp), state.flags, &state, out);
2955 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), headerContents.flags_eks, &state.ekinstate, out);
2961 energyhistory_t enerhist;
2962 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2963 headerContents.flags_enh, &enerhist, out);
2967 PullHistory pullHist = {};
2968 ret = doCptPullHist(gmx_fio_getxdr(fp), TRUE,
2969 headerContents.flagsPullHistory, &pullHist, StatePart::pullHistory, out);
2974 ret = do_cpt_df_hist(gmx_fio_getxdr(fp),
2975 headerContents.flags_dfh, headerContents.nlambda, &state.dfhist, out);
2980 edsamhistory_t edsamhist = {};
2981 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, headerContents.nED, &edsamhist, out);
2986 ret = do_cpt_awh(gmx_fio_getxdr(fp), TRUE,
2987 headerContents.flags_awhh, state.awhHistory.get(), out);
2992 swaphistory_t swaphist = {};
2993 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, headerContents.eSwapCoords, &swaphist, out);
2998 std::vector<gmx_file_position_t> outputfiles;
2999 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, out, headerContents.file_version);
3004 ret = do_cpt_footer(gmx_fio_getxdr(fp), headerContents.file_version);
3011 if (gmx_fio_close(fp) != 0)
3013 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
3017 /* This routine cannot print tons of data, since it is called before the log file is opened. */
3019 read_checkpoint_simulation_part_and_filenames(t_fileio *fp,
3020 int *simulation_part,
3021 std::vector<gmx_file_position_t> *outputfiles)
3027 read_checkpoint_data(fp, simulation_part, &step, &t, &state, outputfiles);
3028 if (gmx_fio_close(fp) != 0)
3030 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");