2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2008,2009,2010,2011,2012 by the GROMACS development team.
5 * Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team.
6 * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
38 /* The source code in this file should be thread-safe.
39 Please keep it that way. */
42 #include "checkpoint.h"
53 #include "buildinfo.h"
54 #include "gromacs/fileio/filetypes.h"
55 #include "gromacs/fileio/gmxfio.h"
56 #include "gromacs/fileio/gmxfio_xdr.h"
57 #include "gromacs/fileio/xdr_datatype.h"
58 #include "gromacs/fileio/xdrf.h"
59 #include "gromacs/gmxlib/network.h"
60 #include "gromacs/math/vec.h"
61 #include "gromacs/math/vecdump.h"
62 #include "gromacs/math/vectypes.h"
63 #include "gromacs/mdtypes/awh_correlation_history.h"
64 #include "gromacs/mdtypes/awh_history.h"
65 #include "gromacs/mdtypes/commrec.h"
66 #include "gromacs/mdtypes/df_history.h"
67 #include "gromacs/mdtypes/edsamhistory.h"
68 #include "gromacs/mdtypes/energyhistory.h"
69 #include "gromacs/mdtypes/inputrec.h"
70 #include "gromacs/mdtypes/md_enums.h"
71 #include "gromacs/mdtypes/mdrunoptions.h"
72 #include "gromacs/mdtypes/observableshistory.h"
73 #include "gromacs/mdtypes/pullhistory.h"
74 #include "gromacs/mdtypes/state.h"
75 #include "gromacs/mdtypes/swaphistory.h"
76 #include "gromacs/trajectory/trajectoryframe.h"
77 #include "gromacs/utility/arrayref.h"
78 #include "gromacs/utility/baseversion.h"
79 #include "gromacs/utility/cstringutil.h"
80 #include "gromacs/utility/fatalerror.h"
81 #include "gromacs/utility/futil.h"
82 #include "gromacs/utility/gmxassert.h"
83 #include "gromacs/utility/int64_to_int.h"
84 #include "gromacs/utility/keyvaluetree.h"
85 #include "gromacs/utility/keyvaluetreebuilder.h"
86 #include "gromacs/utility/keyvaluetreeserializer.h"
87 #include "gromacs/utility/mdmodulenotification.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"
93 #define CPT_MAGIC1 171817
94 #define CPT_MAGIC2 171819
96 /*! \brief Enum of values that describe the contents of a cpt file
97 * whose format matches a version number
99 * The enum helps the code be more self-documenting and ensure merges
100 * do not silently resolve when two patches make the same bump. When
101 * adding new functionality, add a new element just above cptv_Count
102 * in this enumeration, and write code below that does the right thing
103 * according to the value of file_version.
107 cptv_Unknown
= 17, /**< Version before numbering scheme */
108 cptv_RemoveBuildMachineInformation
, /**< remove functionality that makes mdrun builds non-reproducible */
109 cptv_ComPrevStepAsPullGroupReference
, /**< Allow using COM of previous step as pull group PBC reference */
110 cptv_PullAverage
, /**< Added possibility to output average pull force and position */
111 cptv_MdModules
, /**< Added checkpointing for MdModules */
112 cptv_Count
/**< the total number of cptv versions */
115 /*! \brief Version number of the file format written to checkpoint
116 * files by this version of the code.
118 * cpt_version should normally only be changed, via adding a new field
119 * to cptv enumeration, when the header or footer format changes.
121 * The state data format itself is backward and forward compatible.
122 * But old code can not read a new entry that is present in the file
123 * (but can read a new format when new entries are not present).
125 * The cpt_version increases whenever the file format in the main
126 * development branch changes, due to an extension of the cptv enum above.
127 * Backward compatibility for reading old run input files is maintained
128 * by checking this version number against that of the file and then using
129 * the correct code path. */
130 static const int cpt_version
= cptv_Count
- 1;
133 const char* est_names
[estNR
] = { "FE-lambda",
139 "thermostat-integral",
144 "LD-rng-unsupported",
145 "LD-rng-i-unsupported",
158 "MC-rng-unsupported",
159 "MC-rng-i-unsupported",
160 "barostat-integral" };
177 static const char* eeks_names
[eeksNR
] = { "Ekin_n", "Ekinh", "dEkindlambda",
178 "mv_cos", "Ekinf", "Ekinh_old",
179 "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC",
191 eenhENERGY_NSTEPS_SIM
,
192 eenhENERGY_DELTA_H_NN
,
193 eenhENERGY_DELTA_H_LIST
,
194 eenhENERGY_DELTA_H_STARTTIME
,
195 eenhENERGY_DELTA_H_STARTLAMBDA
,
201 epullhPULL_NUMCOORDINATES
,
202 epullhPULL_NUMGROUPS
,
203 epullhPULL_NUMVALUESINXSUM
,
204 epullhPULL_NUMVALUESINFSUM
,
210 epullcoordh_VALUE_REF_SUM
,
211 epullcoordh_VALUE_SUM
,
212 epullcoordh_DR01_SUM
,
213 epullcoordh_DR23_SUM
,
214 epullcoordh_DR45_SUM
,
215 epullcoordh_FSCAL_SUM
,
216 epullcoordh_DYNAX_SUM
,
226 static const char* eenh_names
[eenhNR
] = { "energy_n",
235 "energy_delta_h_list",
236 "energy_delta_h_start_time",
237 "energy_delta_h_start_lambda" };
239 static const char* ePullhNames
[epullhNR
] = { "pullhistory_numcoordinates", "pullhistory_numgroups",
240 "pullhistory_numvaluesinxsum",
241 "pullhistory_numvaluesinfsum" };
243 /* free energy history variables -- need to be preserved over checkpoint */
262 /* free energy history variable names */
263 static const char* edfh_names
[edfhNR
] = { "bEquilibrated",
265 "Wang-Landau Histogram",
273 "accumulated_plus_2",
274 "accumulated_minus_2",
278 /* AWH biasing history variables */
282 eawhhEQUILIBRATEHISTOGRAM
,
286 eawhhUMBRELLAGRIDPOINT
,
288 eawhhLOGSCALEDSAMPLEWEIGHT
,
290 eawhhFORCECORRELATIONGRID
,
294 static const char* eawhh_names
[eawhhNR
] = { "awh_in_initial", "awh_equilibrateHistogram",
295 "awh_histsize", "awh_npoints",
296 "awh_coordpoint", "awh_umbrellaGridpoint",
297 "awh_updatelist", "awh_logScaledSampleWeight",
298 "awh_numupdates", "awh_forceCorrelationGrid" };
306 static const char* epull_prev_step_com_names
[epullsNR
] = { "Pull groups prev step COM" };
309 //! Higher level vector element type, only used for formatting checkpoint dumps
310 enum class CptElementType
312 integer
, //!< integer
313 real
, //!< float or double, not linked to precision of type real
314 real3
, //!< float[3] or double[3], not linked to precision of type real
315 matrix3x3
//!< float[3][3] or double[3][3], not linked to precision of type real
318 //! \brief Parts of the checkpoint state, only used for reporting
321 microState
, //!< The microstate of the simulated system
322 kineticEnergy
, //!< Kinetic energy, needed for T/P-coupling state
323 energyHistory
, //!< Energy observable statistics
324 freeEnergyHistory
, //!< Free-energy state and observable statistics
325 accWeightHistogram
, //!< Accelerated weight histogram method state
326 pullState
, //!< COM of previous step.
327 pullHistory
//!< Pull history statistics (sums since last written output)
330 //! \brief Return the name of a checkpoint entry based on part and part entry
331 static const char* entryName(StatePart part
, int ecpt
)
335 case StatePart::microState
: return est_names
[ecpt
];
336 case StatePart::kineticEnergy
: return eeks_names
[ecpt
];
337 case StatePart::energyHistory
: return eenh_names
[ecpt
];
338 case StatePart::freeEnergyHistory
: return edfh_names
[ecpt
];
339 case StatePart::accWeightHistogram
: return eawhh_names
[ecpt
];
340 case StatePart::pullState
: return epull_prev_step_com_names
[ecpt
];
341 case StatePart::pullHistory
: return ePullhNames
[ecpt
];
347 static void cp_warning(FILE* fp
)
349 fprintf(fp
, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
352 [[noreturn
]] static void cp_error()
354 gmx_fatal(FARGS
, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
357 static void do_cpt_string_err(XDR
* xd
, const char* desc
, gmx::ArrayRef
<char> s
, FILE* list
)
359 char* data
= s
.data();
360 if (xdr_string(xd
, &data
, s
.size()) == 0)
366 fprintf(list
, "%s = %s\n", desc
, data
);
370 static int do_cpt_int(XDR
* xd
, const char* desc
, int* i
, FILE* list
)
372 if (xdr_int(xd
, i
) == 0)
378 fprintf(list
, "%s = %d\n", desc
, *i
);
383 static int do_cpt_u_chars(XDR
* xd
, const char* desc
, int n
, unsigned char* i
, FILE* list
)
387 fprintf(list
, "%s = ", desc
);
390 for (int j
= 0; j
< n
&& res
; j
++)
392 res
&= xdr_u_char(xd
, &i
[j
]);
395 fprintf(list
, "%02x", i
[j
]);
410 static void do_cpt_int_err(XDR
* xd
, const char* desc
, int* i
, FILE* list
)
412 if (do_cpt_int(xd
, desc
, i
, list
) < 0)
418 static void do_cpt_bool_err(XDR
* xd
, const char* desc
, bool* b
, FILE* list
)
420 int i
= static_cast<int>(*b
);
422 if (do_cpt_int(xd
, desc
, &i
, list
) < 0)
430 static void do_cpt_step_err(XDR
* xd
, const char* desc
, int64_t* i
, FILE* list
)
432 char buf
[STEPSTRSIZE
];
434 if (xdr_int64(xd
, i
) == 0)
440 fprintf(list
, "%s = %s\n", desc
, gmx_step_str(*i
, buf
));
444 static void do_cpt_double_err(XDR
* xd
, const char* desc
, double* f
, FILE* list
)
446 if (xdr_double(xd
, f
) == 0)
452 fprintf(list
, "%s = %f\n", desc
, *f
);
456 static void do_cpt_real_err(XDR
* xd
, real
* f
)
459 bool_t res
= xdr_double(xd
, f
);
461 bool_t res
= xdr_float(xd
, f
);
469 static void do_cpt_n_rvecs_err(XDR
* xd
, const char* desc
, int n
, rvec f
[], FILE* list
)
471 for (int i
= 0; i
< n
; i
++)
473 for (int j
= 0; j
< DIM
; j
++)
475 do_cpt_real_err(xd
, &f
[i
][j
]);
481 pr_rvecs(list
, 0, desc
, f
, n
);
493 static const int value
= xdr_datatype_int
;
497 struct xdr_type
<float>
499 static const int value
= xdr_datatype_float
;
503 struct xdr_type
<double>
505 static const int value
= xdr_datatype_double
;
508 //! \brief Returns size in byte of an xdr_datatype
509 static inline unsigned int sizeOfXdrType(int xdrType
)
513 case xdr_datatype_int
: return sizeof(int);
514 case xdr_datatype_float
: return sizeof(float);
515 case xdr_datatype_double
: return sizeof(double);
516 default: GMX_RELEASE_ASSERT(false, "XDR data type not implemented");
522 //! \brief Returns the XDR process function for i/o of an XDR type
523 static inline xdrproc_t
xdrProc(int xdrType
)
527 case xdr_datatype_int
: return reinterpret_cast<xdrproc_t
>(xdr_int
);
528 case xdr_datatype_float
: return reinterpret_cast<xdrproc_t
>(xdr_float
);
529 case xdr_datatype_double
: return reinterpret_cast<xdrproc_t
>(xdr_double
);
530 default: GMX_RELEASE_ASSERT(false, "XDR data type not implemented");
536 /*! \brief Lists or only reads an xdr vector from checkpoint file
538 * When list!=NULL reads and lists the \p nf vector elements of type \p xdrType.
539 * The header for the print is set by \p part and \p ecpt.
540 * The formatting of the printing is set by \p cptElementType.
541 * When list==NULL only reads the elements.
543 static bool_t
listXdrVector(XDR
* xd
, StatePart part
, int ecpt
, int nf
, int xdrType
, FILE* list
, CptElementType cptElementType
)
547 const unsigned int elemSize
= sizeOfXdrType(xdrType
);
548 std::vector
<char> data(nf
* elemSize
);
549 res
= xdr_vector(xd
, data
.data(), nf
, elemSize
, xdrProc(xdrType
));
555 case xdr_datatype_int
:
556 pr_ivec(list
, 0, entryName(part
, ecpt
), reinterpret_cast<const int*>(data
.data()), nf
, TRUE
);
558 case xdr_datatype_float
:
560 if (cptElementType
== CptElementType::real3
)
562 pr_rvecs(list
, 0, entryName(part
, ecpt
),
563 reinterpret_cast<const rvec
*>(data
.data()), nf
/ 3);
568 /* Note: With double precision code dumping a single precision rvec will produce float iso rvec print, but that's a minor annoyance */
569 pr_fvec(list
, 0, entryName(part
, ecpt
),
570 reinterpret_cast<const float*>(data
.data()), nf
, TRUE
);
573 case xdr_datatype_double
:
575 if (cptElementType
== CptElementType::real3
)
577 pr_rvecs(list
, 0, entryName(part
, ecpt
),
578 reinterpret_cast<const rvec
*>(data
.data()), nf
/ 3);
583 /* Note: With single precision code dumping a double precision rvec will produce float iso rvec print, but that's a minor annoyance */
584 pr_dvec(list
, 0, entryName(part
, ecpt
),
585 reinterpret_cast<const double*>(data
.data()), nf
, TRUE
);
588 default: GMX_RELEASE_ASSERT(false, "Data type not implemented for listing");
595 //! \brief Convert a double array, typed char*, to float
596 gmx_unused
static void convertArrayRealPrecision(const char* c
, float* v
, int n
)
598 const double* d
= reinterpret_cast<const double*>(c
);
599 for (int i
= 0; i
< n
; i
++)
601 v
[i
] = static_cast<float>(d
[i
]);
605 //! \brief Convert a float array, typed char*, to double
606 static void convertArrayRealPrecision(const char* c
, double* v
, int n
)
608 const float* f
= reinterpret_cast<const float*>(c
);
609 for (int i
= 0; i
< n
; i
++)
611 v
[i
] = static_cast<double>(f
[i
]);
615 //! \brief Generate an error for trying to convert to integer
616 static void convertArrayRealPrecision(const char gmx_unused
* c
, int gmx_unused
* v
, int gmx_unused n
)
618 GMX_RELEASE_ASSERT(false,
619 "We only expect type mismatches between float and double, not integer");
622 /*! \brief Low-level routine for reading/writing a vector of reals from/to file.
624 * This is the only routine that does the actually i/o of real vector,
625 * all other routines are intermediate level routines for specific real
626 * data types, calling this routine.
627 * Currently this routine is (too) complex, since it handles both real *
628 * and std::vector<real>. Using real * is deprecated and this routine
629 * will simplify a lot when only std::vector needs to be supported.
631 * The routine is generic to vectors with different allocators,
632 * because as part of reading a checkpoint there are vectors whose
633 * size is not known until reading has progressed far enough, so a
634 * resize method must be called.
636 * When not listing, we use either v or vector, depending on which is !=NULL.
637 * If nval >= 0, nval is used; on read this should match the passed value.
638 * If nval n<0, *nptr (with v) or vector->size() is used. On read using v,
639 * the value is stored in nptr
641 template<typename T
, typename AllocatorType
>
642 static int doVectorLow(XDR
* xd
,
649 std::vector
<T
, AllocatorType
>* vector
,
651 CptElementType cptElementType
)
653 GMX_RELEASE_ASSERT(list
!= nullptr || (v
!= nullptr && vector
== nullptr)
654 || (v
== nullptr && vector
!= nullptr),
655 "Without list, we should have exactly one of v and vector != NULL");
659 int numElemInTheFile
;
664 GMX_RELEASE_ASSERT(nptr
== nullptr, "With nval>=0 we should have nptr==NULL");
665 numElemInTheFile
= nval
;
671 GMX_RELEASE_ASSERT(nptr
!= nullptr, "With nval<0 we should have nptr!=NULL");
672 numElemInTheFile
= *nptr
;
676 numElemInTheFile
= vector
->size();
680 /* Read/write the vector element count */
681 res
= xdr_int(xd
, &numElemInTheFile
);
686 /* Read/write the element data type */
687 constexpr int xdrTypeInTheCode
= xdr_type
<T
>::value
;
688 int xdrTypeInTheFile
= xdrTypeInTheCode
;
689 res
= xdr_int(xd
, &xdrTypeInTheFile
);
695 if (list
== nullptr && (sflags
& (1 << ecpt
)))
699 if (numElemInTheFile
!= nval
)
702 "Count mismatch for state entry %s, code count is %d, file count is %d\n",
703 entryName(part
, ecpt
), nval
, numElemInTheFile
);
706 else if (nptr
!= nullptr)
708 *nptr
= numElemInTheFile
;
711 bool typesMatch
= (xdrTypeInTheFile
== xdrTypeInTheCode
);
715 sprintf(buf
, "mismatch for state entry %s, code precision is %s, file precision is %s",
716 entryName(part
, ecpt
), xdr_datatype_names
[xdrTypeInTheCode
],
717 xdr_datatype_names
[xdrTypeInTheFile
]);
719 /* Matching int and real should never occur, but check anyhow */
720 if (xdrTypeInTheFile
== xdr_datatype_int
|| xdrTypeInTheCode
== xdr_datatype_int
)
723 "Type %s: incompatible checkpoint formats or corrupted checkpoint file.", buf
);
732 snew(*v
, numElemInTheFile
);
738 /* This conditional ensures that we don't resize on write.
739 * In particular in the state where this code was written
740 * vector has a size of numElemInThefile and we
741 * don't want to lose that padding here.
743 if (vector
->size() < static_cast<unsigned int>(numElemInTheFile
))
745 vector
->resize(numElemInTheFile
);
753 vChar
= reinterpret_cast<char*>(vp
);
757 snew(vChar
, numElemInTheFile
* sizeOfXdrType(xdrTypeInTheFile
));
759 res
= xdr_vector(xd
, vChar
, numElemInTheFile
, sizeOfXdrType(xdrTypeInTheFile
),
760 xdrProc(xdrTypeInTheFile
));
768 /* In the old code float-double conversion came for free.
769 * In the new code we still support it, mainly because
770 * the tip4p_continue regression test makes use of this.
771 * It's an open question if we do or don't want to allow this.
773 convertArrayRealPrecision(vChar
, vp
, numElemInTheFile
);
779 res
= listXdrVector(xd
, part
, ecpt
, numElemInTheFile
, xdrTypeInTheFile
, list
, cptElementType
);
785 //! \brief Read/Write a std::vector, on read checks the number of elements matches \p numElements, if specified.
788 doVector(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, std::vector
<T
>* vector
, FILE* list
, int numElements
= -1)
790 return doVectorLow
<T
>(xd
, part
, ecpt
, sflags
, numElements
, nullptr, nullptr, vector
, list
,
791 CptElementType::real
);
794 //! \brief Read/Write an ArrayRef<real>.
795 static int doRealArrayRef(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, gmx::ArrayRef
<real
> vector
, FILE* list
)
797 real
* v_real
= vector
.data();
798 return doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, vector
.size(), nullptr,
799 &v_real
, nullptr, list
, CptElementType::real
);
802 //! Convert from view of RVec to view of real.
803 static gmx::ArrayRef
<real
> realArrayRefFromRVecArrayRef(gmx::ArrayRef
<gmx::RVec
> ofRVecs
)
805 return gmx::arrayRefFromArray
<real
>(reinterpret_cast<real
*>(ofRVecs
.data()), ofRVecs
.size() * DIM
);
808 //! \brief Read/Write a PaddedVector whose value_type is RVec.
809 template<typename PaddedVectorOfRVecType
>
811 doRvecVector(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, PaddedVectorOfRVecType
* v
, int numAtoms
, FILE* list
)
813 const int numReals
= numAtoms
* DIM
;
818 sflags
& (1 << ecpt
),
819 "When not listing, the flag for the entry should be set when requesting i/o");
820 GMX_RELEASE_ASSERT(v
->size() == numAtoms
, "v should have sufficient size for numAtoms");
822 return doRealArrayRef(xd
, part
, ecpt
, sflags
, realArrayRefFromRVecArrayRef(makeArrayRef(*v
)), list
);
826 // Use the rebind facility to change the value_type of the
827 // allocator from RVec to real.
828 using realAllocator
=
829 typename
std::allocator_traits
<typename
PaddedVectorOfRVecType::allocator_type
>::template rebind_alloc
<real
>;
830 return doVectorLow
<real
, realAllocator
>(xd
, part
, ecpt
, sflags
, numReals
, nullptr, nullptr,
831 nullptr, list
, CptElementType::real
);
835 /* This function stores n along with the reals for reading,
836 * but on reading it assumes that n matches the value in the checkpoint file,
837 * a fatal error is generated when this is not the case.
839 static int do_cpte_reals(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int n
, real
** v
, FILE* list
)
841 return doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, n
, nullptr, v
, nullptr,
842 list
, CptElementType::real
);
845 /* This function does the same as do_cpte_reals,
846 * except that on reading it ignores the passed value of *n
847 * and stores the value read from the checkpoint file in *n.
849 static int do_cpte_n_reals(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int* n
, real
** v
, FILE* list
)
851 return doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, -1, n
, v
, nullptr, list
,
852 CptElementType::real
);
855 static int do_cpte_real(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, real
* r
, FILE* list
)
857 return doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, 1, nullptr, &r
, nullptr,
858 list
, CptElementType::real
);
861 static int do_cpte_ints(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int n
, int** v
, FILE* list
)
863 return doVectorLow
<int, std::allocator
<int>>(xd
, part
, ecpt
, sflags
, n
, nullptr, v
, nullptr,
864 list
, CptElementType::integer
);
867 static int do_cpte_int(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int* i
, FILE* list
)
869 return do_cpte_ints(xd
, part
, ecpt
, sflags
, 1, &i
, list
);
872 static int do_cpte_bool(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, bool* b
, FILE* list
)
874 int i
= static_cast<int>(*b
);
875 int ret
= do_cpte_int(xd
, part
, ecpt
, sflags
, &i
, list
);
880 static int do_cpte_doubles(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int n
, double** v
, FILE* list
)
882 return doVectorLow
<double, std::allocator
<double>>(xd
, part
, ecpt
, sflags
, n
, nullptr, v
,
883 nullptr, list
, CptElementType::real
);
886 static int do_cpte_double(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, double* r
, FILE* list
)
888 return do_cpte_doubles(xd
, part
, ecpt
, sflags
, 1, &r
, list
);
891 static int do_cpte_matrix(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, matrix v
, FILE* list
)
897 ret
= doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, DIM
* DIM
, nullptr, &vr
,
898 nullptr, nullptr, CptElementType::matrix3x3
);
900 if (list
&& ret
== 0)
902 pr_rvecs(list
, 0, entryName(part
, ecpt
), v
, DIM
);
909 static int do_cpte_nmatrix(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int n
, real
** v
, FILE* list
)
913 char name
[CPTSTRLEN
];
920 for (i
= 0; i
< n
; i
++)
922 reti
= doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, n
, nullptr, &(v
[i
]),
923 nullptr, nullptr, CptElementType::matrix3x3
);
924 if (list
&& reti
== 0)
926 sprintf(name
, "%s[%d]", entryName(part
, ecpt
), i
);
927 pr_reals(list
, 0, name
, v
[i
], n
);
937 static int do_cpte_matrices(XDR
* xd
, StatePart part
, int ecpt
, int sflags
, int n
, matrix
** v
, FILE* list
)
940 matrix
*vp
, *va
= nullptr;
946 res
= xdr_int(xd
, &nf
);
951 if (list
== nullptr && nf
!= n
)
953 gmx_fatal(FARGS
, "Count mismatch for state entry %s, code count is %d, file count is %d\n",
954 entryName(part
, ecpt
), n
, nf
);
956 if (list
|| !(sflags
& (1 << ecpt
)))
969 snew(vr
, nf
* DIM
* DIM
);
970 for (i
= 0; i
< nf
; i
++)
972 for (j
= 0; j
< DIM
; j
++)
974 for (k
= 0; k
< DIM
; k
++)
976 vr
[(i
* DIM
+ j
) * DIM
+ k
] = vp
[i
][j
][k
];
980 ret
= doVectorLow
<real
, std::allocator
<real
>>(xd
, part
, ecpt
, sflags
, nf
* DIM
* DIM
, nullptr,
981 &vr
, nullptr, nullptr, CptElementType::matrix3x3
);
982 for (i
= 0; i
< nf
; i
++)
984 for (j
= 0; j
< DIM
; j
++)
986 for (k
= 0; k
< DIM
; k
++)
988 vp
[i
][j
][k
] = vr
[(i
* DIM
+ j
) * DIM
+ k
];
994 if (list
&& ret
== 0)
996 for (i
= 0; i
< nf
; i
++)
998 pr_rvecs(list
, 0, entryName(part
, ecpt
), vp
[i
], DIM
);
1009 static void do_cpt_header(XDR
* xd
, gmx_bool bRead
, FILE* list
, CheckpointHeaderContents
* contents
)
1022 res
= xdr_int(xd
, &magic
);
1026 "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
1028 if (magic
!= CPT_MAGIC1
)
1031 "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
1032 "The checkpoint file is corrupted or not a checkpoint file",
1038 gmx_gethostname(fhost
, 255);
1040 do_cpt_string_err(xd
, "GROMACS version", contents
->version
, list
);
1041 // The following fields are no longer ever written with meaningful
1042 // content, but because they precede the file version, there is no
1043 // good way for new code to read the old and new formats, nor a
1044 // good way for old code to avoid giving an error while reading a
1045 // new format. So we read and write a field that no longer has a
1047 do_cpt_string_err(xd
, "GROMACS build time UNUSED", contents
->btime_UNUSED
, list
);
1048 do_cpt_string_err(xd
, "GROMACS build user UNUSED", contents
->buser_UNUSED
, list
);
1049 do_cpt_string_err(xd
, "GROMACS build host UNUSED", contents
->bhost_UNUSED
, list
);
1050 do_cpt_string_err(xd
, "generating program", contents
->fprog
, list
);
1051 do_cpt_string_err(xd
, "generation time", contents
->ftime
, list
);
1052 contents
->file_version
= cpt_version
;
1053 do_cpt_int_err(xd
, "checkpoint file version", &contents
->file_version
, list
);
1054 if (contents
->file_version
> cpt_version
)
1057 "Attempting to read a checkpoint file of version %d with code of version %d\n",
1058 contents
->file_version
, cpt_version
);
1060 if (contents
->file_version
>= 13)
1062 do_cpt_int_err(xd
, "GROMACS double precision", &contents
->double_prec
, list
);
1066 contents
->double_prec
= -1;
1068 if (contents
->file_version
>= 12)
1070 do_cpt_string_err(xd
, "generating host", fhost
, list
);
1072 do_cpt_int_err(xd
, "#atoms", &contents
->natoms
, list
);
1073 do_cpt_int_err(xd
, "#T-coupling groups", &contents
->ngtc
, list
);
1074 if (contents
->file_version
>= 10)
1076 do_cpt_int_err(xd
, "#Nose-Hoover T-chains", &contents
->nhchainlength
, list
);
1080 contents
->nhchainlength
= 1;
1082 if (contents
->file_version
>= 11)
1084 do_cpt_int_err(xd
, "#Nose-Hoover T-chains for barostat ", &contents
->nnhpres
, list
);
1088 contents
->nnhpres
= 0;
1090 if (contents
->file_version
>= 14)
1092 do_cpt_int_err(xd
, "# of total lambda states ", &contents
->nlambda
, list
);
1096 contents
->nlambda
= 0;
1098 do_cpt_int_err(xd
, "integrator", &contents
->eIntegrator
, list
);
1099 if (contents
->file_version
>= 3)
1101 do_cpt_int_err(xd
, "simulation part #", &contents
->simulation_part
, list
);
1105 contents
->simulation_part
= 1;
1107 if (contents
->file_version
>= 5)
1109 do_cpt_step_err(xd
, "step", &contents
->step
, list
);
1114 do_cpt_int_err(xd
, "step", &idum
, list
);
1115 contents
->step
= static_cast<int64_t>(idum
);
1117 do_cpt_double_err(xd
, "t", &contents
->t
, list
);
1118 do_cpt_int_err(xd
, "#PP-ranks", &contents
->nnodes
, list
);
1119 do_cpt_int_err(xd
, "dd_nc[x]", &contents
->dd_nc
[XX
], list
);
1120 do_cpt_int_err(xd
, "dd_nc[y]", &contents
->dd_nc
[YY
], list
);
1121 do_cpt_int_err(xd
, "dd_nc[z]", &contents
->dd_nc
[ZZ
], list
);
1122 do_cpt_int_err(xd
, "#PME-only ranks", &contents
->npme
, list
);
1123 do_cpt_int_err(xd
, "state flags", &contents
->flags_state
, list
);
1124 if (contents
->file_version
>= 4)
1126 do_cpt_int_err(xd
, "ekin data flags", &contents
->flags_eks
, list
);
1127 do_cpt_int_err(xd
, "energy history flags", &contents
->flags_enh
, list
);
1131 contents
->flags_eks
= 0;
1132 contents
->flags_enh
= (contents
->flags_state
>> (estORIRE_DTAV
+ 1));
1133 contents
->flags_state
= (contents
->flags_state
1134 & ~((1 << (estORIRE_DTAV
+ 1)) | (1 << (estORIRE_DTAV
+ 2))
1135 | (1 << (estORIRE_DTAV
+ 3))));
1137 if (contents
->file_version
>= 14)
1139 do_cpt_int_err(xd
, "df history flags", &contents
->flags_dfh
, list
);
1143 contents
->flags_dfh
= 0;
1146 if (contents
->file_version
>= 15)
1148 do_cpt_int_err(xd
, "ED data sets", &contents
->nED
, list
);
1155 if (contents
->file_version
>= 16)
1157 do_cpt_int_err(xd
, "swap", &contents
->eSwapCoords
, list
);
1161 contents
->eSwapCoords
= eswapNO
;
1164 if (contents
->file_version
>= 17)
1166 do_cpt_int_err(xd
, "AWH history flags", &contents
->flags_awhh
, list
);
1170 contents
->flags_awhh
= 0;
1173 if (contents
->file_version
>= 18)
1175 do_cpt_int_err(xd
, "pull history flags", &contents
->flagsPullHistory
, list
);
1179 contents
->flagsPullHistory
= 0;
1183 static int do_cpt_footer(XDR
* xd
, int file_version
)
1188 if (file_version
>= 2)
1191 res
= xdr_int(xd
, &magic
);
1196 if (magic
!= CPT_MAGIC2
)
1205 static int do_cpt_state(XDR
* xd
, int fflags
, t_state
* state
, FILE* list
)
1208 const StatePart part
= StatePart::microState
;
1209 const int sflags
= state
->flags
;
1210 // If reading, state->natoms was probably just read, so
1211 // allocations need to be managed. If writing, this won't change
1212 // anything that matters.
1213 state_change_natoms(state
, state
->natoms
);
1214 for (int i
= 0; (i
< estNR
&& ret
== 0); i
++)
1216 if (fflags
& (1 << i
))
1221 ret
= doRealArrayRef(
1222 xd
, part
, i
, sflags
,
1223 gmx::arrayRefFromArray
<real
>(state
->lambda
.data(), state
->lambda
.size()),
1227 ret
= do_cpte_int(xd
, part
, i
, sflags
, &state
->fep_state
, list
);
1229 case estBOX
: ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->box
, list
); break;
1231 ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->box_rel
, list
);
1233 case estBOXV
: ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->boxv
, list
); break;
1235 ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->pres_prev
, list
);
1238 ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->svir_prev
, list
);
1241 ret
= do_cpte_matrix(xd
, part
, i
, sflags
, state
->fvir_prev
, list
);
1244 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->nosehoover_xi
, list
);
1247 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->nosehoover_vxi
, list
);
1250 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->nhpres_xi
, list
);
1253 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->nhpres_vxi
, list
);
1256 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->therm_integral
, list
);
1259 ret
= do_cpte_double(xd
, part
, i
, sflags
, &state
->baros_integral
, list
);
1261 case estVETA
: ret
= do_cpte_real(xd
, part
, i
, sflags
, &state
->veta
, list
); break;
1262 case estVOL0
: ret
= do_cpte_real(xd
, part
, i
, sflags
, &state
->vol0
, list
); break;
1264 ret
= doRvecVector(xd
, part
, i
, sflags
, &state
->x
, state
->natoms
, list
);
1267 ret
= doRvecVector(xd
, part
, i
, sflags
, &state
->v
, state
->natoms
, list
);
1269 /* The RNG entries are no longer written,
1270 * the next 4 lines are only for reading old files.
1272 case estLD_RNG_NOTSUPPORTED
:
1273 ret
= do_cpte_ints(xd
, part
, i
, sflags
, 0, nullptr, list
);
1275 case estLD_RNGI_NOTSUPPORTED
:
1276 ret
= do_cpte_ints(xd
, part
, i
, sflags
, 0, nullptr, list
);
1278 case estMC_RNG_NOTSUPPORTED
:
1279 ret
= do_cpte_ints(xd
, part
, i
, sflags
, 0, nullptr, list
);
1281 case estMC_RNGI_NOTSUPPORTED
:
1282 ret
= do_cpte_ints(xd
, part
, i
, sflags
, 0, nullptr, list
);
1284 case estDISRE_INITF
:
1285 ret
= do_cpte_real(xd
, part
, i
, sflags
, &state
->hist
.disre_initf
, list
);
1287 case estDISRE_RM3TAV
:
1288 ret
= do_cpte_n_reals(xd
, part
, i
, sflags
, &state
->hist
.ndisrepairs
,
1289 &state
->hist
.disre_rm3tav
, list
);
1291 case estORIRE_INITF
:
1292 ret
= do_cpte_real(xd
, part
, i
, sflags
, &state
->hist
.orire_initf
, list
);
1295 ret
= do_cpte_n_reals(xd
, part
, i
, sflags
, &state
->hist
.norire_Dtav
,
1296 &state
->hist
.orire_Dtav
, list
);
1298 case estPULLCOMPREVSTEP
:
1299 ret
= doVector
<double>(xd
, part
, i
, sflags
, &state
->pull_com_prev_step
, list
);
1303 "Unknown state entry %d\n"
1304 "You are reading a checkpoint file written by different code, which "
1313 static int do_cpt_ekinstate(XDR
* xd
, int fflags
, ekinstate_t
* ekins
, FILE* list
)
1317 const StatePart part
= StatePart::kineticEnergy
;
1318 for (int i
= 0; (i
< eeksNR
&& ret
== 0); i
++)
1320 if (fflags
& (1 << i
))
1326 ret
= do_cpte_int(xd
, part
, i
, fflags
, &ekins
->ekin_n
, list
);
1329 ret
= do_cpte_matrices(xd
, part
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinh
, list
);
1332 ret
= do_cpte_matrices(xd
, part
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinf
, list
);
1335 ret
= do_cpte_matrices(xd
, part
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinh_old
, list
);
1338 ret
= do_cpte_matrix(xd
, part
, i
, fflags
, ekins
->ekin_total
, list
);
1340 case eeksEKINSCALEF
:
1341 ret
= doVector
<double>(xd
, part
, i
, fflags
, &ekins
->ekinscalef_nhc
, list
);
1344 ret
= doVector
<double>(xd
, part
, i
, fflags
, &ekins
->vscale_nhc
, list
);
1346 case eeksEKINSCALEH
:
1347 ret
= doVector
<double>(xd
, part
, i
, fflags
, &ekins
->ekinscaleh_nhc
, list
);
1350 ret
= do_cpte_real(xd
, part
, i
, fflags
, &ekins
->dekindl
, list
);
1352 case eeksMVCOS
: ret
= do_cpte_real(xd
, part
, i
, fflags
, &ekins
->mvcos
, list
); break;
1355 "Unknown ekin data state entry %d\n"
1356 "You are probably reading a new checkpoint file with old code",
1366 static int do_cpt_swapstate(XDR
* xd
, gmx_bool bRead
, int eSwapCoords
, swaphistory_t
* swapstate
, FILE* list
)
1368 int swap_cpt_version
= 2;
1370 if (eSwapCoords
== eswapNO
)
1375 swapstate
->bFromCpt
= bRead
;
1376 swapstate
->eSwapCoords
= eSwapCoords
;
1378 do_cpt_int_err(xd
, "swap checkpoint version", &swap_cpt_version
, list
);
1379 if (bRead
&& swap_cpt_version
< 2)
1382 "Cannot read checkpoint files that were written with old versions"
1383 "of the ion/water position swapping protocol.\n");
1386 do_cpt_int_err(xd
, "swap coupling steps", &swapstate
->nAverage
, list
);
1388 /* When reading, init_swapcoords has not been called yet,
1389 * so we have to allocate memory first. */
1390 do_cpt_int_err(xd
, "number of ion types", &swapstate
->nIonTypes
, list
);
1393 snew(swapstate
->ionType
, swapstate
->nIonTypes
);
1396 for (int ic
= 0; ic
< eCompNR
; ic
++)
1398 for (int ii
= 0; ii
< swapstate
->nIonTypes
; ii
++)
1400 swapstateIons_t
* gs
= &swapstate
->ionType
[ii
];
1404 do_cpt_int_err(xd
, "swap requested atoms", &gs
->nMolReq
[ic
], list
);
1408 do_cpt_int_err(xd
, "swap requested atoms p", gs
->nMolReq_p
[ic
], list
);
1413 do_cpt_int_err(xd
, "swap influx net", &gs
->inflow_net
[ic
], list
);
1417 do_cpt_int_err(xd
, "swap influx net p", gs
->inflow_net_p
[ic
], list
);
1420 if (bRead
&& (nullptr == gs
->nMolPast
[ic
]))
1422 snew(gs
->nMolPast
[ic
], swapstate
->nAverage
);
1425 for (int j
= 0; j
< swapstate
->nAverage
; j
++)
1429 do_cpt_int_err(xd
, "swap past atom counts", &gs
->nMolPast
[ic
][j
], list
);
1433 do_cpt_int_err(xd
, "swap past atom counts p", &gs
->nMolPast_p
[ic
][j
], list
);
1439 /* Ion flux per channel */
1440 for (int ic
= 0; ic
< eChanNR
; ic
++)
1442 for (int ii
= 0; ii
< swapstate
->nIonTypes
; ii
++)
1444 swapstateIons_t
* gs
= &swapstate
->ionType
[ii
];
1448 do_cpt_int_err(xd
, "channel flux A->B", &gs
->fluxfromAtoB
[ic
], list
);
1452 do_cpt_int_err(xd
, "channel flux A->B p", gs
->fluxfromAtoB_p
[ic
], list
);
1457 /* Ion flux leakage */
1460 do_cpt_int_err(xd
, "flux leakage", &swapstate
->fluxleak
, list
);
1464 do_cpt_int_err(xd
, "flux leakage", swapstate
->fluxleak_p
, list
);
1468 for (int ii
= 0; ii
< swapstate
->nIonTypes
; ii
++)
1470 swapstateIons_t
* gs
= &swapstate
->ionType
[ii
];
1472 do_cpt_int_err(xd
, "number of ions", &gs
->nMol
, list
);
1476 snew(gs
->channel_label
, gs
->nMol
);
1477 snew(gs
->comp_from
, gs
->nMol
);
1480 do_cpt_u_chars(xd
, "channel history", gs
->nMol
, gs
->channel_label
, list
);
1481 do_cpt_u_chars(xd
, "domain history", gs
->nMol
, gs
->comp_from
, list
);
1484 /* Save the last known whole positions to checkpoint
1485 * file to be able to also make multimeric channels whole in PBC */
1486 do_cpt_int_err(xd
, "Ch0 atoms", &swapstate
->nat
[eChan0
], list
);
1487 do_cpt_int_err(xd
, "Ch1 atoms", &swapstate
->nat
[eChan1
], list
);
1490 snew(swapstate
->xc_old_whole
[eChan0
], swapstate
->nat
[eChan0
]);
1491 snew(swapstate
->xc_old_whole
[eChan1
], swapstate
->nat
[eChan1
]);
1492 do_cpt_n_rvecs_err(xd
, "Ch0 whole x", swapstate
->nat
[eChan0
], swapstate
->xc_old_whole
[eChan0
], list
);
1493 do_cpt_n_rvecs_err(xd
, "Ch1 whole x", swapstate
->nat
[eChan1
], swapstate
->xc_old_whole
[eChan1
], list
);
1497 do_cpt_n_rvecs_err(xd
, "Ch0 whole x", swapstate
->nat
[eChan0
],
1498 *swapstate
->xc_old_whole_p
[eChan0
], list
);
1499 do_cpt_n_rvecs_err(xd
, "Ch1 whole x", swapstate
->nat
[eChan1
],
1500 *swapstate
->xc_old_whole_p
[eChan1
], list
);
1507 static int do_cpt_enerhist(XDR
* xd
, gmx_bool bRead
, int fflags
, energyhistory_t
* enerhist
, FILE* list
)
1516 GMX_RELEASE_ASSERT(enerhist
!= nullptr,
1517 "With energy history, we need a valid enerhist pointer");
1519 /* This is stored/read for backward compatibility */
1520 int energyHistoryNumEnergies
= 0;
1523 enerhist
->nsteps
= 0;
1525 enerhist
->nsteps_sim
= 0;
1526 enerhist
->nsum_sim
= 0;
1528 else if (enerhist
!= nullptr)
1530 energyHistoryNumEnergies
= enerhist
->ener_sum_sim
.size();
1533 delta_h_history_t
* deltaH
= enerhist
->deltaHForeignLambdas
.get();
1534 const StatePart part
= StatePart::energyHistory
;
1535 for (int i
= 0; (i
< eenhNR
&& ret
== 0); i
++)
1537 if (fflags
& (1 << i
))
1542 ret
= do_cpte_int(xd
, part
, i
, fflags
, &energyHistoryNumEnergies
, list
);
1544 case eenhENERGY_AVER
:
1545 ret
= doVector
<double>(xd
, part
, i
, fflags
, &enerhist
->ener_ave
, list
);
1547 case eenhENERGY_SUM
:
1548 ret
= doVector
<double>(xd
, part
, i
, fflags
, &enerhist
->ener_sum
, list
);
1550 case eenhENERGY_NSUM
:
1551 do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsum
, list
);
1553 case eenhENERGY_SUM_SIM
:
1554 ret
= doVector
<double>(xd
, part
, i
, fflags
, &enerhist
->ener_sum_sim
, list
);
1556 case eenhENERGY_NSUM_SIM
:
1557 do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsum_sim
, list
);
1559 case eenhENERGY_NSTEPS
:
1560 do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsteps
, list
);
1562 case eenhENERGY_NSTEPS_SIM
:
1563 do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsteps_sim
, list
);
1565 case eenhENERGY_DELTA_H_NN
:
1568 if (!bRead
&& deltaH
!= nullptr)
1570 numDeltaH
= deltaH
->dh
.size();
1572 do_cpt_int_err(xd
, eenh_names
[i
], &numDeltaH
, list
);
1575 if (deltaH
== nullptr)
1577 enerhist
->deltaHForeignLambdas
= std::make_unique
<delta_h_history_t
>();
1578 deltaH
= enerhist
->deltaHForeignLambdas
.get();
1580 deltaH
->dh
.resize(numDeltaH
);
1581 deltaH
->start_lambda_set
= FALSE
;
1585 case eenhENERGY_DELTA_H_LIST
:
1586 for (auto dh
: deltaH
->dh
)
1588 ret
= doVector
<real
>(xd
, part
, i
, fflags
, &dh
, list
);
1591 case eenhENERGY_DELTA_H_STARTTIME
:
1592 ret
= do_cpte_double(xd
, part
, i
, fflags
, &(deltaH
->start_time
), list
);
1594 case eenhENERGY_DELTA_H_STARTLAMBDA
:
1595 ret
= do_cpte_double(xd
, part
, i
, fflags
, &(deltaH
->start_lambda
), list
);
1599 "Unknown energy history entry %d\n"
1600 "You are probably reading a new checkpoint file with old code",
1606 if ((fflags
& (1 << eenhENERGY_SUM
)) && !(fflags
& (1 << eenhENERGY_SUM_SIM
)))
1608 /* Assume we have an old file format and copy sum to sum_sim */
1609 enerhist
->ener_sum_sim
= enerhist
->ener_sum
;
1612 if ((fflags
& (1 << eenhENERGY_NSUM
)) && !(fflags
& (1 << eenhENERGY_NSTEPS
)))
1614 /* Assume we have an old file format and copy nsum to nsteps */
1615 enerhist
->nsteps
= enerhist
->nsum
;
1617 if ((fflags
& (1 << eenhENERGY_NSUM_SIM
)) && !(fflags
& (1 << eenhENERGY_NSTEPS_SIM
)))
1619 /* Assume we have an old file format and copy nsum to nsteps */
1620 enerhist
->nsteps_sim
= enerhist
->nsum_sim
;
1626 static int doCptPullCoordHist(XDR
* xd
, PullCoordinateHistory
* pullCoordHist
, const StatePart part
, FILE* list
)
1631 flags
|= ((1 << epullcoordh_VALUE_REF_SUM
) | (1 << epullcoordh_VALUE_SUM
)
1632 | (1 << epullcoordh_DR01_SUM
) | (1 << epullcoordh_DR23_SUM
)
1633 | (1 << epullcoordh_DR45_SUM
) | (1 << epullcoordh_FSCAL_SUM
));
1635 for (int i
= 0; i
< epullcoordh_NR
&& ret
== 0; i
++)
1639 case epullcoordh_VALUE_REF_SUM
:
1640 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->valueRef
), list
);
1642 case epullcoordh_VALUE_SUM
:
1643 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->value
), list
);
1645 case epullcoordh_DR01_SUM
:
1646 for (int j
= 0; j
< DIM
&& ret
== 0; j
++)
1648 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->dr01
[j
]), list
);
1651 case epullcoordh_DR23_SUM
:
1652 for (int j
= 0; j
< DIM
&& ret
== 0; j
++)
1654 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->dr23
[j
]), list
);
1657 case epullcoordh_DR45_SUM
:
1658 for (int j
= 0; j
< DIM
&& ret
== 0; j
++)
1660 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->dr45
[j
]), list
);
1663 case epullcoordh_FSCAL_SUM
:
1664 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->scalarForce
), list
);
1666 case epullcoordh_DYNAX_SUM
:
1667 for (int j
= 0; j
< DIM
&& ret
== 0; j
++)
1669 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullCoordHist
->dynaX
[j
]), list
);
1678 static int doCptPullGroupHist(XDR
* xd
, PullGroupHistory
* pullGroupHist
, const StatePart part
, FILE* list
)
1683 flags
|= ((1 << epullgrouph_X_SUM
));
1685 for (int i
= 0; i
< epullgrouph_NR
; i
++)
1689 case epullgrouph_X_SUM
:
1690 for (int j
= 0; j
< DIM
&& ret
== 0; j
++)
1692 ret
= do_cpte_double(xd
, part
, i
, flags
, &(pullGroupHist
->x
[j
]), list
);
1702 static int doCptPullHist(XDR
* xd
, gmx_bool bRead
, int fflags
, PullHistory
* pullHist
, const StatePart part
, FILE* list
)
1705 int pullHistoryNumCoordinates
= 0;
1706 int pullHistoryNumGroups
= 0;
1708 /* Retain the number of terms in the sum and the number of coordinates (used for writing
1709 * average pull forces and coordinates) in the pull history, in temporary variables,
1710 * in case they cannot be read from the checkpoint, in order to have backward compatibility */
1713 pullHist
->numValuesInXSum
= 0;
1714 pullHist
->numValuesInFSum
= 0;
1716 else if (pullHist
!= nullptr)
1718 pullHistoryNumCoordinates
= pullHist
->pullCoordinateSums
.size();
1719 pullHistoryNumGroups
= pullHist
->pullGroupSums
.size();
1723 GMX_RELEASE_ASSERT(fflags
== 0, "Without pull history, all flags should be off");
1726 for (int i
= 0; (i
< epullhNR
&& ret
== 0); i
++)
1728 if (fflags
& (1 << i
))
1732 case epullhPULL_NUMCOORDINATES
:
1733 ret
= do_cpte_int(xd
, part
, i
, fflags
, &pullHistoryNumCoordinates
, list
);
1735 case epullhPULL_NUMGROUPS
:
1736 do_cpt_int_err(xd
, eenh_names
[i
], &pullHistoryNumGroups
, list
);
1738 case epullhPULL_NUMVALUESINXSUM
:
1739 do_cpt_int_err(xd
, eenh_names
[i
], &pullHist
->numValuesInXSum
, list
);
1741 case epullhPULL_NUMVALUESINFSUM
:
1742 do_cpt_int_err(xd
, eenh_names
[i
], &pullHist
->numValuesInFSum
, list
);
1746 "Unknown pull history entry %d\n"
1747 "You are probably reading a new checkpoint file with old code",
1754 pullHist
->pullCoordinateSums
.resize(pullHistoryNumCoordinates
);
1755 pullHist
->pullGroupSums
.resize(pullHistoryNumGroups
);
1757 if (pullHist
->numValuesInXSum
> 0 || pullHist
->numValuesInFSum
> 0)
1759 for (size_t i
= 0; i
< pullHist
->pullCoordinateSums
.size() && ret
== 0; i
++)
1761 ret
= doCptPullCoordHist(xd
, &(pullHist
->pullCoordinateSums
[i
]), part
, list
);
1763 for (size_t i
= 0; i
< pullHist
->pullGroupSums
.size() && ret
== 0; i
++)
1765 ret
= doCptPullGroupHist(xd
, &(pullHist
->pullGroupSums
[i
]), part
, list
);
1772 static int do_cpt_df_hist(XDR
* xd
, int fflags
, int nlambda
, df_history_t
** dfhistPtr
, FILE* list
)
1781 if (*dfhistPtr
== nullptr)
1783 snew(*dfhistPtr
, 1);
1784 (*dfhistPtr
)->nlambda
= nlambda
;
1785 init_df_history(*dfhistPtr
, nlambda
);
1787 df_history_t
* dfhist
= *dfhistPtr
;
1789 const StatePart part
= StatePart::freeEnergyHistory
;
1790 for (int i
= 0; (i
< edfhNR
&& ret
== 0); i
++)
1792 if (fflags
& (1 << i
))
1797 ret
= do_cpte_bool(xd
, part
, i
, fflags
, &dfhist
->bEquil
, list
);
1800 ret
= do_cpte_ints(xd
, part
, i
, fflags
, nlambda
, &dfhist
->n_at_lam
, list
);
1803 ret
= do_cpte_reals(xd
, part
, i
, fflags
, nlambda
, &dfhist
->wl_histo
, list
);
1806 ret
= do_cpte_real(xd
, part
, i
, fflags
, &dfhist
->wl_delta
, list
);
1808 case edfhSUMWEIGHTS
:
1809 ret
= do_cpte_reals(xd
, part
, i
, fflags
, nlambda
, &dfhist
->sum_weights
, list
);
1812 ret
= do_cpte_reals(xd
, part
, i
, fflags
, nlambda
, &dfhist
->sum_dg
, list
);
1815 ret
= do_cpte_reals(xd
, part
, i
, fflags
, nlambda
, &dfhist
->sum_minvar
, list
);
1818 ret
= do_cpte_reals(xd
, part
, i
, fflags
, nlambda
, &dfhist
->sum_variance
, list
);
1821 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->accum_p
, list
);
1824 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->accum_m
, list
);
1827 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->accum_p2
, list
);
1830 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->accum_m2
, list
);
1833 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->Tij
, list
);
1836 ret
= do_cpte_nmatrix(xd
, part
, i
, fflags
, nlambda
, dfhist
->Tij_empirical
, list
);
1841 "Unknown df history entry %d\n"
1842 "You are probably reading a new checkpoint file with old code",
1852 /* This function stores the last whole configuration of the reference and
1853 * average structure in the .cpt file
1855 static int do_cpt_EDstate(XDR
* xd
, gmx_bool bRead
, int nED
, edsamhistory_t
* EDstate
, FILE* list
)
1862 EDstate
->bFromCpt
= bRead
;
1865 /* When reading, init_edsam has not been called yet,
1866 * so we have to allocate memory first. */
1869 snew(EDstate
->nref
, EDstate
->nED
);
1870 snew(EDstate
->old_sref
, EDstate
->nED
);
1871 snew(EDstate
->nav
, EDstate
->nED
);
1872 snew(EDstate
->old_sav
, EDstate
->nED
);
1875 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1876 for (int i
= 0; i
< EDstate
->nED
; i
++)
1880 /* Reference structure SREF */
1881 sprintf(buf
, "ED%d # of atoms in reference structure", i
+ 1);
1882 do_cpt_int_err(xd
, buf
, &EDstate
->nref
[i
], list
);
1883 sprintf(buf
, "ED%d x_ref", i
+ 1);
1886 snew(EDstate
->old_sref
[i
], EDstate
->nref
[i
]);
1887 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nref
[i
], EDstate
->old_sref
[i
], list
);
1891 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nref
[i
], EDstate
->old_sref_p
[i
], list
);
1894 /* Average structure SAV */
1895 sprintf(buf
, "ED%d # of atoms in average structure", i
+ 1);
1896 do_cpt_int_err(xd
, buf
, &EDstate
->nav
[i
], list
);
1897 sprintf(buf
, "ED%d x_av", i
+ 1);
1900 snew(EDstate
->old_sav
[i
], EDstate
->nav
[i
]);
1901 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nav
[i
], EDstate
->old_sav
[i
], list
);
1905 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nav
[i
], EDstate
->old_sav_p
[i
], list
);
1912 static int do_cpt_correlation_grid(XDR
* xd
,
1914 gmx_unused
int fflags
,
1915 gmx::CorrelationGridHistory
* corrGrid
,
1921 do_cpt_int_err(xd
, eawhh_names
[eawhh
], &(corrGrid
->numCorrelationTensors
), list
);
1922 do_cpt_int_err(xd
, eawhh_names
[eawhh
], &(corrGrid
->tensorSize
), list
);
1923 do_cpt_int_err(xd
, eawhh_names
[eawhh
], &(corrGrid
->blockDataListSize
), list
);
1927 initCorrelationGridHistory(corrGrid
, corrGrid
->numCorrelationTensors
, corrGrid
->tensorSize
,
1928 corrGrid
->blockDataListSize
);
1931 for (gmx::CorrelationBlockDataHistory
& blockData
: corrGrid
->blockDataBuffer
)
1933 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.blockSumWeight
), list
);
1934 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.blockSumSquareWeight
), list
);
1935 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.blockSumWeightX
), list
);
1936 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.blockSumWeightY
), list
);
1937 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.sumOverBlocksSquareBlockWeight
), list
);
1938 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.sumOverBlocksBlockSquareWeight
), list
);
1939 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.sumOverBlocksBlockWeightBlockWeightX
), list
);
1940 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.sumOverBlocksBlockWeightBlockWeightY
), list
);
1941 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.blockLength
), list
);
1942 do_cpt_int_err(xd
, eawhh_names
[eawhh
], &(blockData
.previousBlockIndex
), list
);
1943 do_cpt_double_err(xd
, eawhh_names
[eawhh
], &(blockData
.correlationIntegral
), list
);
1949 static int do_cpt_awh_bias(XDR
* xd
, gmx_bool bRead
, int fflags
, gmx::AwhBiasHistory
* biasHistory
, FILE* list
)
1953 gmx::AwhBiasStateHistory
* state
= &biasHistory
->state
;
1954 for (int i
= 0; (i
< eawhhNR
&& ret
== 0); i
++)
1956 if (fflags
& (1 << i
))
1960 case eawhhIN_INITIAL
:
1961 do_cpt_bool_err(xd
, eawhh_names
[i
], &state
->in_initial
, list
);
1963 case eawhhEQUILIBRATEHISTOGRAM
:
1964 do_cpt_bool_err(xd
, eawhh_names
[i
], &state
->equilibrateHistogram
, list
);
1967 do_cpt_double_err(xd
, eawhh_names
[i
], &state
->histSize
, list
);
1974 numPoints
= biasHistory
->pointState
.size();
1976 do_cpt_int_err(xd
, eawhh_names
[i
], &numPoints
, list
);
1979 biasHistory
->pointState
.resize(numPoints
);
1983 case eawhhCOORDPOINT
:
1984 for (auto& psh
: biasHistory
->pointState
)
1986 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.target
, list
);
1987 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.free_energy
, list
);
1988 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.bias
, list
);
1989 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.weightsum_iteration
, list
);
1990 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.weightsum_covering
, list
);
1991 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.weightsum_tot
, list
);
1992 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.weightsum_ref
, list
);
1993 do_cpt_step_err(xd
, eawhh_names
[i
], &psh
.last_update_index
, list
);
1994 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.log_pmfsum
, list
);
1995 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.visits_iteration
, list
);
1996 do_cpt_double_err(xd
, eawhh_names
[i
], &psh
.visits_tot
, list
);
1999 case eawhhUMBRELLAGRIDPOINT
:
2000 do_cpt_int_err(xd
, eawhh_names
[i
], &(state
->umbrellaGridpoint
), list
);
2002 case eawhhUPDATELIST
:
2003 do_cpt_int_err(xd
, eawhh_names
[i
], &(state
->origin_index_updatelist
), list
);
2004 do_cpt_int_err(xd
, eawhh_names
[i
], &(state
->end_index_updatelist
), list
);
2006 case eawhhLOGSCALEDSAMPLEWEIGHT
:
2007 do_cpt_double_err(xd
, eawhh_names
[i
], &(state
->logScaledSampleWeight
), list
);
2008 do_cpt_double_err(xd
, eawhh_names
[i
], &(state
->maxLogScaledSampleWeight
), list
);
2010 case eawhhNUMUPDATES
:
2011 do_cpt_step_err(xd
, eawhh_names
[i
], &(state
->numUpdates
), list
);
2013 case eawhhFORCECORRELATIONGRID
:
2014 ret
= do_cpt_correlation_grid(xd
, bRead
, fflags
,
2015 &biasHistory
->forceCorrelationGrid
, list
, i
);
2017 default: gmx_fatal(FARGS
, "Unknown awh history entry %d\n", i
);
2025 static int do_cpt_awh(XDR
* xd
, gmx_bool bRead
, int fflags
, gmx::AwhHistory
* awhHistory
, FILE* list
)
2031 std::shared_ptr
<gmx::AwhHistory
> awhHistoryLocal
;
2033 if (awhHistory
== nullptr)
2035 GMX_RELEASE_ASSERT(bRead
,
2036 "do_cpt_awh should not be called for writing without an AwhHistory");
2038 awhHistoryLocal
= std::make_shared
<gmx::AwhHistory
>();
2039 awhHistory
= awhHistoryLocal
.get();
2042 /* To be able to read and write the AWH data several parameters determining the layout of the AWH structs need to be known
2043 (nbias, npoints, etc.). The best thing (?) would be to have these passed to this function. When writing to a checkpoint
2044 these parameters are available in awh_history (after calling init_awh_history). When reading from a checkpoint though, there
2045 is no initialized awh_history (it is initialized and set in this function). The AWH parameters have not always been read
2046 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
2047 when writing a checkpoint, also storing parameters needed for future reading of the checkpoint.
2049 Another issue is that some variables that AWH checkpoints don't have a registered enum and string (e.g. nbias below).
2050 One difficulty is the multilevel structure of the data which would need to be represented somehow. */
2055 numBias
= awhHistory
->bias
.size();
2057 do_cpt_int_err(xd
, "awh_nbias", &numBias
, list
);
2061 awhHistory
->bias
.resize(numBias
);
2063 for (auto& bias
: awhHistory
->bias
)
2065 ret
= do_cpt_awh_bias(xd
, bRead
, fflags
, &bias
, list
);
2071 do_cpt_double_err(xd
, "awh_potential_offset", &awhHistory
->potentialOffset
, list
);
2076 static void do_cpt_mdmodules(int fileVersion
,
2077 t_fileio
* checkpointFileHandle
,
2078 const gmx::MdModulesNotifier
& mdModulesNotifier
)
2080 if (fileVersion
>= cptv_MdModules
)
2082 gmx::FileIOXdrSerializer
serializer(checkpointFileHandle
);
2083 gmx::KeyValueTreeObject mdModuleCheckpointParameterTree
=
2084 gmx::deserializeKeyValueTree(&serializer
);
2085 gmx::MdModulesCheckpointReadingDataOnMaster mdModuleCheckpointReadingDataOnMaster
= {
2086 mdModuleCheckpointParameterTree
, fileVersion
2088 mdModulesNotifier
.notifier_
.notify(mdModuleCheckpointReadingDataOnMaster
);
2092 static int do_cpt_files(XDR
* xd
, gmx_bool bRead
, std::vector
<gmx_file_position_t
>* outputfiles
, FILE* list
, int file_version
)
2095 gmx_off_t mask
= 0xFFFFFFFFL
;
2096 int offset_high
, offset_low
;
2097 std::array
<char, CPTSTRLEN
> buf
;
2098 GMX_RELEASE_ASSERT(outputfiles
, "Must have valid outputfiles");
2100 // Ensure that reading pre-allocates outputfiles, while writing
2101 // writes what is already there.
2102 int nfiles
= outputfiles
->size();
2103 if (do_cpt_int(xd
, "number of output files", &nfiles
, list
) != 0)
2109 outputfiles
->resize(nfiles
);
2112 for (auto& outputfile
: *outputfiles
)
2114 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
2117 do_cpt_string_err(xd
, "output filename", buf
, list
);
2118 std::copy(std::begin(buf
), std::end(buf
), std::begin(outputfile
.filename
));
2120 if (do_cpt_int(xd
, "file_offset_high", &offset_high
, list
) != 0)
2124 if (do_cpt_int(xd
, "file_offset_low", &offset_low
, list
) != 0)
2128 outputfile
.offset
= (static_cast<gmx_off_t
>(offset_high
) << 32)
2129 | (static_cast<gmx_off_t
>(offset_low
) & mask
);
2133 do_cpt_string_err(xd
, "output filename", outputfile
.filename
, list
);
2135 offset
= outputfile
.offset
;
2143 offset_low
= static_cast<int>(offset
& mask
);
2144 offset_high
= static_cast<int>((offset
>> 32) & mask
);
2146 if (do_cpt_int(xd
, "file_offset_high", &offset_high
, list
) != 0)
2150 if (do_cpt_int(xd
, "file_offset_low", &offset_low
, list
) != 0)
2155 if (file_version
>= 8)
2157 if (do_cpt_int(xd
, "file_checksum_size", &outputfile
.checksumSize
, list
) != 0)
2161 if (do_cpt_u_chars(xd
, "file_checksum", outputfile
.checksum
.size(),
2162 outputfile
.checksum
.data(), list
)
2170 outputfile
.checksumSize
= -1;
2176 static void mpiBarrierBeforeRename(const bool applyMpiBarrierBeforeRename
, MPI_Comm mpiBarrierCommunicator
)
2178 if (applyMpiBarrierBeforeRename
)
2181 MPI_Barrier(mpiBarrierCommunicator
);
2183 GMX_RELEASE_ASSERT(false, "Should not request a barrier without MPI");
2184 GMX_UNUSED_VALUE(mpiBarrierCommunicator
);
2189 void write_checkpoint(const char* fn
,
2190 gmx_bool bNumberAndKeep
,
2192 const t_commrec
* cr
,
2196 int simulation_part
,
2202 ObservablesHistory
* observablesHistory
,
2203 const gmx::MdModulesNotifier
& mdModulesNotifier
,
2204 bool applyMpiBarrierBeforeRename
,
2205 MPI_Comm mpiBarrierCommunicator
)
2208 char* fntemp
; /* the temporary checkpoint file name */
2210 char buf
[1024], suffix
[5 + STEPSTRSIZE
], sbuf
[STEPSTRSIZE
];
2213 if (DOMAINDECOMP(cr
))
2215 npmenodes
= cr
->npmenodes
;
2223 /* make the new temporary filename */
2224 snew(fntemp
, std::strlen(fn
) + 5 + STEPSTRSIZE
);
2225 std::strcpy(fntemp
, fn
);
2226 fntemp
[std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1] = '\0';
2227 sprintf(suffix
, "_%s%s", "step", gmx_step_str(step
, sbuf
));
2228 std::strcat(fntemp
, suffix
);
2229 std::strcat(fntemp
, fn
+ std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1);
2231 /* if we can't rename, we just overwrite the cpt file.
2232 * dangerous if interrupted.
2234 snew(fntemp
, std::strlen(fn
));
2235 std::strcpy(fntemp
, fn
);
2237 std::string timebuf
= gmx_format_current_time();
2241 fprintf(fplog
, "Writing checkpoint, step %s at %s\n\n", gmx_step_str(step
, buf
), timebuf
.c_str());
2244 /* Get offsets for open files */
2245 auto outputfiles
= gmx_fio_get_output_file_positions();
2247 fp
= gmx_fio_open(fntemp
, "w");
2250 if (state
->ekinstate
.bUpToDate
)
2252 flags_eks
= ((1 << eeksEKIN_N
) | (1 << eeksEKINH
) | (1 << eeksEKINF
) | (1 << eeksEKINO
)
2253 | (1 << eeksEKINSCALEF
) | (1 << eeksEKINSCALEH
) | (1 << eeksVSCALE
)
2254 | (1 << eeksDEKINDL
) | (1 << eeksMVCOS
));
2261 energyhistory_t
* enerhist
= observablesHistory
->energyHistory
.get();
2263 if (enerhist
!= nullptr && (enerhist
->nsum
> 0 || enerhist
->nsum_sim
> 0))
2265 flags_enh
|= (1 << eenhENERGY_N
) | (1 << eenhENERGY_NSTEPS
) | (1 << eenhENERGY_NSTEPS_SIM
);
2266 if (enerhist
->nsum
> 0)
2268 flags_enh
|= ((1 << eenhENERGY_AVER
) | (1 << eenhENERGY_SUM
) | (1 << eenhENERGY_NSUM
));
2270 if (enerhist
->nsum_sim
> 0)
2272 flags_enh
|= ((1 << eenhENERGY_SUM_SIM
) | (1 << eenhENERGY_NSUM_SIM
));
2274 if (enerhist
->deltaHForeignLambdas
!= nullptr)
2276 flags_enh
|= ((1 << eenhENERGY_DELTA_H_NN
) | (1 << eenhENERGY_DELTA_H_LIST
)
2277 | (1 << eenhENERGY_DELTA_H_STARTTIME
) | (1 << eenhENERGY_DELTA_H_STARTLAMBDA
));
2281 PullHistory
* pullHist
= observablesHistory
->pullHistory
.get();
2282 int flagsPullHistory
= 0;
2283 if (pullHist
!= nullptr && (pullHist
->numValuesInXSum
> 0 || pullHist
->numValuesInFSum
> 0))
2285 flagsPullHistory
|= (1 << epullhPULL_NUMCOORDINATES
);
2286 flagsPullHistory
|= ((1 << epullhPULL_NUMGROUPS
) | (1 << epullhPULL_NUMVALUESINXSUM
)
2287 | (1 << epullhPULL_NUMVALUESINFSUM
));
2293 flags_dfh
= ((1 << edfhBEQUIL
) | (1 << edfhNATLAMBDA
) | (1 << edfhSUMWEIGHTS
)
2294 | (1 << edfhSUMDG
) | (1 << edfhTIJ
) | (1 << edfhTIJEMP
));
2297 flags_dfh
|= ((1 << edfhWLDELTA
) | (1 << edfhWLHISTO
));
2299 if ((elamstats
== elamstatsMINVAR
) || (elamstats
== elamstatsBARKER
)
2300 || (elamstats
== elamstatsMETROPOLIS
))
2302 flags_dfh
|= ((1 << edfhACCUMP
) | (1 << edfhACCUMM
) | (1 << edfhACCUMP2
)
2303 | (1 << edfhACCUMM2
) | (1 << edfhSUMMINVAR
) | (1 << edfhSUMVAR
));
2312 if (state
->awhHistory
!= nullptr && !state
->awhHistory
->bias
.empty())
2314 flags_awhh
|= ((1 << eawhhIN_INITIAL
) | (1 << eawhhEQUILIBRATEHISTOGRAM
) | (1 << eawhhHISTSIZE
)
2315 | (1 << eawhhNPOINTS
) | (1 << eawhhCOORDPOINT
) | (1 << eawhhUMBRELLAGRIDPOINT
)
2316 | (1 << eawhhUPDATELIST
) | (1 << eawhhLOGSCALEDSAMPLEWEIGHT
)
2317 | (1 << eawhhNUMUPDATES
) | (1 << eawhhFORCECORRELATIONGRID
));
2320 /* We can check many more things now (CPU, acceleration, etc), but
2321 * it is highly unlikely to have two separate builds with exactly
2322 * the same version, user, time, and build host!
2325 int nlambda
= (state
->dfhist
? state
->dfhist
->nlambda
: 0);
2327 edsamhistory_t
* edsamhist
= observablesHistory
->edsamHistory
.get();
2328 int nED
= (edsamhist
? edsamhist
->nED
: 0);
2330 swaphistory_t
* swaphist
= observablesHistory
->swapHistory
.get();
2331 int eSwapCoords
= (swaphist
? swaphist
->eSwapCoords
: eswapNO
);
2333 CheckpointHeaderContents headerContents
= { 0,
2351 state
->nhchainlength
,
2361 std::strcpy(headerContents
.version
, gmx_version());
2362 std::strcpy(headerContents
.fprog
, gmx::getProgramContext().fullBinaryPath());
2363 std::strcpy(headerContents
.ftime
, timebuf
.c_str());
2364 if (DOMAINDECOMP(cr
))
2366 copy_ivec(domdecCells
, headerContents
.dd_nc
);
2369 do_cpt_header(gmx_fio_getxdr(fp
), FALSE
, nullptr, &headerContents
);
2371 if ((do_cpt_state(gmx_fio_getxdr(fp
), state
->flags
, state
, nullptr) < 0)
2372 || (do_cpt_ekinstate(gmx_fio_getxdr(fp
), flags_eks
, &state
->ekinstate
, nullptr) < 0)
2373 || (do_cpt_enerhist(gmx_fio_getxdr(fp
), FALSE
, flags_enh
, enerhist
, nullptr) < 0)
2374 || (doCptPullHist(gmx_fio_getxdr(fp
), FALSE
, flagsPullHistory
, pullHist
, StatePart::pullHistory
, nullptr)
2376 || (do_cpt_df_hist(gmx_fio_getxdr(fp
), flags_dfh
, nlambda
, &state
->dfhist
, nullptr) < 0)
2377 || (do_cpt_EDstate(gmx_fio_getxdr(fp
), FALSE
, nED
, edsamhist
, nullptr) < 0)
2378 || (do_cpt_awh(gmx_fio_getxdr(fp
), FALSE
, flags_awhh
, state
->awhHistory
.get(), nullptr) < 0)
2379 || (do_cpt_swapstate(gmx_fio_getxdr(fp
), FALSE
, eSwapCoords
, swaphist
, nullptr) < 0)
2380 || (do_cpt_files(gmx_fio_getxdr(fp
), FALSE
, &outputfiles
, nullptr, headerContents
.file_version
) < 0))
2382 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2385 // Checkpointing MdModules
2387 gmx::KeyValueTreeBuilder builder
;
2388 gmx::MdModulesWriteCheckpointData mdModulesWriteCheckpoint
= { builder
.rootObject(),
2389 headerContents
.file_version
};
2390 mdModulesNotifier
.notifier_
.notify(mdModulesWriteCheckpoint
);
2391 auto tree
= builder
.build();
2392 gmx::FileIOXdrSerializer
serializer(fp
);
2393 gmx::serializeKeyValueTree(tree
, &serializer
);
2396 do_cpt_footer(gmx_fio_getxdr(fp
), headerContents
.file_version
);
2398 /* we really, REALLY, want to make sure to physically write the checkpoint,
2399 and all the files it depends on, out to disk. Because we've
2400 opened the checkpoint with gmx_fio_open(), it's in our list
2402 ret
= gmx_fio_all_output_fsync();
2407 sprintf(buf
, "Cannot fsync '%s'; maybe you are out of disk space?", gmx_fio_getname(ret
));
2409 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV
) == nullptr)
2415 gmx_warning("%s", buf
);
2419 if (gmx_fio_close(fp
) != 0)
2421 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2424 /* we don't move the checkpoint if the user specified they didn't want it,
2425 or if the fsyncs failed */
2427 if (!bNumberAndKeep
&& !ret
)
2431 /* Rename the previous checkpoint file */
2432 mpiBarrierBeforeRename(applyMpiBarrierBeforeRename
, mpiBarrierCommunicator
);
2434 std::strcpy(buf
, fn
);
2435 buf
[std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1] = '\0';
2436 std::strcat(buf
, "_prev");
2437 std::strcat(buf
, fn
+ std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1);
2440 /* we copy here so that if something goes wrong between now and
2441 * the rename below, there's always a state.cpt.
2442 * If renames are atomic (such as in POSIX systems),
2443 * this copying should be unneccesary.
2445 gmx_file_copy(fn
, buf
, FALSE
);
2446 /* We don't really care if this fails:
2447 * there's already a new checkpoint.
2452 gmx_file_rename(fn
, buf
);
2456 /* Rename the checkpoint file from the temporary to the final name */
2457 mpiBarrierBeforeRename(applyMpiBarrierBeforeRename
, mpiBarrierCommunicator
);
2459 if (gmx_file_rename(fntemp
, fn
) != 0)
2461 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
2464 #endif /* GMX_NO_RENAME */
2469 /* Always FAH checkpoint immediately after a Gromacs checkpoint.
2471 * Note that it is critical that we save a FAH checkpoint directly
2472 * after writing a Gromacs checkpoint. If the program dies, either
2473 * by the machine powering off suddenly or the process being,
2474 * killed, FAH can recover files that have only appended data by
2475 * truncating them to the last recorded length. The Gromacs
2476 * checkpoint does not just append data, it is fully rewritten each
2477 * time so a crash between moving the new Gromacs checkpoint file in
2478 * to place and writing a FAH checkpoint is not recoverable. Thus
2479 * the time between these operations must be kept as short a
2486 static void check_int(FILE* fplog
, const char* type
, int p
, int f
, gmx_bool
* mm
)
2488 bool foundMismatch
= (p
!= f
);
2496 fprintf(fplog
, " %s mismatch,\n", type
);
2497 fprintf(fplog
, " current program: %d\n", p
);
2498 fprintf(fplog
, " checkpoint file: %d\n", f
);
2499 fprintf(fplog
, "\n");
2503 static void check_string(FILE* fplog
, const char* type
, const char* p
, const char* f
, gmx_bool
* mm
)
2505 bool foundMismatch
= (std::strcmp(p
, f
) != 0);
2513 fprintf(fplog
, " %s mismatch,\n", type
);
2514 fprintf(fplog
, " current program: %s\n", p
);
2515 fprintf(fplog
, " checkpoint file: %s\n", f
);
2516 fprintf(fplog
, "\n");
2520 static void check_match(FILE* fplog
,
2521 const t_commrec
* cr
,
2523 const CheckpointHeaderContents
& headerContents
,
2524 gmx_bool reproducibilityRequested
)
2526 /* Note that this check_string on the version will also print a message
2527 * when only the minor version differs. But we only print a warning
2528 * message further down with reproducibilityRequested=TRUE.
2530 gmx_bool versionDiffers
= FALSE
;
2531 check_string(fplog
, "Version", gmx_version(), headerContents
.version
, &versionDiffers
);
2533 gmx_bool precisionDiffers
= FALSE
;
2534 check_int(fplog
, "Double prec.", GMX_DOUBLE
, headerContents
.double_prec
, &precisionDiffers
);
2535 if (precisionDiffers
)
2537 const char msg_precision_difference
[] =
2538 "You are continuing a simulation with a different precision. Not matching\n"
2539 "single/double precision will lead to precision or performance loss.\n";
2542 fprintf(fplog
, "%s\n", msg_precision_difference
);
2546 gmx_bool mm
= (versionDiffers
|| precisionDiffers
);
2548 if (reproducibilityRequested
)
2550 check_string(fplog
, "Program name", gmx::getProgramContext().fullBinaryPath(),
2551 headerContents
.fprog
, &mm
);
2553 check_int(fplog
, "#ranks", cr
->nnodes
, headerContents
.nnodes
, &mm
);
2556 if (cr
->nnodes
> 1 && reproducibilityRequested
)
2558 check_int(fplog
, "#PME-ranks", cr
->npmenodes
, headerContents
.npme
, &mm
);
2560 int npp
= cr
->nnodes
;
2561 if (cr
->npmenodes
>= 0)
2563 npp
-= cr
->npmenodes
;
2565 int npp_f
= headerContents
.nnodes
;
2566 if (headerContents
.npme
>= 0)
2568 npp_f
-= headerContents
.npme
;
2572 check_int(fplog
, "#DD-cells[x]", dd_nc
[XX
], headerContents
.dd_nc
[XX
], &mm
);
2573 check_int(fplog
, "#DD-cells[y]", dd_nc
[YY
], headerContents
.dd_nc
[YY
], &mm
);
2574 check_int(fplog
, "#DD-cells[z]", dd_nc
[ZZ
], headerContents
.dd_nc
[ZZ
], &mm
);
2580 /* Gromacs should be able to continue from checkpoints between
2581 * different patch level versions, but we do not guarantee
2582 * compatibility between different major/minor versions - check this.
2586 sscanf(gmx_version(), "%5d", &gmx_major
);
2587 int ret
= sscanf(headerContents
.version
, "%5d", &cpt_major
);
2588 gmx_bool majorVersionDiffers
= (ret
< 1 || gmx_major
!= cpt_major
);
2590 const char msg_major_version_difference
[] =
2591 "The current GROMACS major version is not identical to the one that\n"
2592 "generated the checkpoint file. In principle GROMACS does not support\n"
2593 "continuation from checkpoints between different versions, so we advise\n"
2594 "against this. If you still want to try your luck we recommend that you use\n"
2595 "the -noappend flag to keep your output files from the two versions separate.\n"
2596 "This might also work around errors where the output fields in the energy\n"
2597 "file have changed between the different versions.\n";
2599 const char msg_mismatch_notice
[] =
2600 "GROMACS patchlevel, binary or parallel settings differ from previous run.\n"
2601 "Continuation is exact, but not guaranteed to be binary identical.\n";
2603 if (majorVersionDiffers
)
2607 fprintf(fplog
, "%s\n", msg_major_version_difference
);
2610 else if (reproducibilityRequested
)
2612 /* Major & minor versions match at least, but something is different. */
2615 fprintf(fplog
, "%s\n", msg_mismatch_notice
);
2621 static void read_checkpoint(const char* fn
,
2623 const t_commrec
* cr
,
2626 int* init_fep_state
,
2627 CheckpointHeaderContents
* headerContents
,
2629 ObservablesHistory
* observablesHistory
,
2630 gmx_bool reproducibilityRequested
,
2631 const gmx::MdModulesNotifier
& mdModulesNotifier
)
2634 char buf
[STEPSTRSIZE
];
2637 fp
= gmx_fio_open(fn
, "r");
2638 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, nullptr, headerContents
);
2640 // If we are appending, then we don't want write to the open log
2641 // file because we still need to compute a checksum for it. In
2642 // that case, the filehandle will be nullptr. Otherwise, we report
2643 // to the new log file about the checkpoint file that we are
2645 FILE* fplog
= gmx_fio_getfp(logfio
);
2648 fprintf(fplog
, "\n");
2649 fprintf(fplog
, "Reading checkpoint file %s\n", fn
);
2650 fprintf(fplog
, " file generated by: %s\n", headerContents
->fprog
);
2651 fprintf(fplog
, " file generated at: %s\n", headerContents
->ftime
);
2652 fprintf(fplog
, " GROMACS double prec.: %d\n", headerContents
->double_prec
);
2653 fprintf(fplog
, " simulation part #: %d\n", headerContents
->simulation_part
);
2654 fprintf(fplog
, " step: %s\n", gmx_step_str(headerContents
->step
, buf
));
2655 fprintf(fplog
, " time: %f\n", headerContents
->t
);
2656 fprintf(fplog
, "\n");
2659 if (headerContents
->natoms
!= state
->natoms
)
2662 "Checkpoint file is for a system of %d atoms, while the current system consists "
2664 headerContents
->natoms
, state
->natoms
);
2666 if (headerContents
->ngtc
!= state
->ngtc
)
2669 "Checkpoint file is for a system of %d T-coupling groups, while the current "
2670 "system consists of %d T-coupling groups",
2671 headerContents
->ngtc
, state
->ngtc
);
2673 if (headerContents
->nnhpres
!= state
->nnhpres
)
2676 "Checkpoint file is for a system of %d NH-pressure-coupling variables, while the "
2677 "current system consists of %d NH-pressure-coupling variables",
2678 headerContents
->nnhpres
, state
->nnhpres
);
2681 int nlambdaHistory
= (state
->dfhist
? state
->dfhist
->nlambda
: 0);
2682 if (headerContents
->nlambda
!= nlambdaHistory
)
2685 "Checkpoint file is for a system with %d lambda states, while the current system "
2686 "consists of %d lambda states",
2687 headerContents
->nlambda
, nlambdaHistory
);
2690 init_gtc_state(state
, state
->ngtc
, state
->nnhpres
,
2691 headerContents
->nhchainlength
); /* need to keep this here to keep the tpr format working */
2692 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
2694 if (headerContents
->eIntegrator
!= eIntegrator
)
2697 "Cannot change integrator during a checkpoint restart. Perhaps you should make a "
2698 "new .tpr with grompp -f new.mdp -t %s",
2702 if (headerContents
->flags_state
!= state
->flags
)
2705 "Cannot change a simulation algorithm during a checkpoint restart. Perhaps you "
2706 "should make a new .tpr with grompp -f new.mdp -t %s",
2712 check_match(fplog
, cr
, dd_nc
, *headerContents
, reproducibilityRequested
);
2715 ret
= do_cpt_state(gmx_fio_getxdr(fp
), headerContents
->flags_state
, state
, nullptr);
2716 *init_fep_state
= state
->fep_state
; /* there should be a better way to do this than setting it
2717 here. Investigate for 5.0. */
2722 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), headerContents
->flags_eks
, &state
->ekinstate
, nullptr);
2727 state
->ekinstate
.hasReadEkinState
= (((headerContents
->flags_eks
& (1 << eeksEKINH
)) != 0)
2728 || ((headerContents
->flags_eks
& (1 << eeksEKINF
)) != 0)
2729 || ((headerContents
->flags_eks
& (1 << eeksEKINO
)) != 0)
2730 || (((headerContents
->flags_eks
& (1 << eeksEKINSCALEF
))
2731 | (headerContents
->flags_eks
& (1 << eeksEKINSCALEH
))
2732 | (headerContents
->flags_eks
& (1 << eeksVSCALE
)))
2735 if (headerContents
->flags_enh
&& observablesHistory
->energyHistory
== nullptr)
2737 observablesHistory
->energyHistory
= std::make_unique
<energyhistory_t
>();
2739 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
, headerContents
->flags_enh
,
2740 observablesHistory
->energyHistory
.get(), nullptr);
2746 if (headerContents
->flagsPullHistory
)
2748 if (observablesHistory
->pullHistory
== nullptr)
2750 observablesHistory
->pullHistory
= std::make_unique
<PullHistory
>();
2752 ret
= doCptPullHist(gmx_fio_getxdr(fp
), TRUE
, headerContents
->flagsPullHistory
,
2753 observablesHistory
->pullHistory
.get(), StatePart::pullHistory
, nullptr);
2760 if (headerContents
->file_version
< 6)
2763 "Continuing from checkpoint files written before GROMACS 4.5 is not supported");
2766 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
), headerContents
->flags_dfh
, headerContents
->nlambda
,
2767 &state
->dfhist
, nullptr);
2773 if (headerContents
->nED
> 0 && observablesHistory
->edsamHistory
== nullptr)
2775 observablesHistory
->edsamHistory
= std::make_unique
<edsamhistory_t
>(edsamhistory_t
{});
2777 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
->nED
,
2778 observablesHistory
->edsamHistory
.get(), nullptr);
2784 if (headerContents
->flags_awhh
!= 0 && state
->awhHistory
== nullptr)
2786 state
->awhHistory
= std::make_shared
<gmx::AwhHistory
>();
2788 ret
= do_cpt_awh(gmx_fio_getxdr(fp
), TRUE
, headerContents
->flags_awhh
, state
->awhHistory
.get(), nullptr);
2794 if (headerContents
->eSwapCoords
!= eswapNO
&& observablesHistory
->swapHistory
== nullptr)
2796 observablesHistory
->swapHistory
= std::make_unique
<swaphistory_t
>(swaphistory_t
{});
2798 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
->eSwapCoords
,
2799 observablesHistory
->swapHistory
.get(), nullptr);
2805 std::vector
<gmx_file_position_t
> outputfiles
;
2806 ret
= do_cpt_files(gmx_fio_getxdr(fp
), TRUE
, &outputfiles
, nullptr, headerContents
->file_version
);
2811 do_cpt_mdmodules(headerContents
->file_version
, fp
, mdModulesNotifier
);
2812 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), headerContents
->file_version
);
2817 if (gmx_fio_close(fp
) != 0)
2819 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2824 void load_checkpoint(const char* fn
,
2826 const t_commrec
* cr
,
2830 ObservablesHistory
* observablesHistory
,
2831 gmx_bool reproducibilityRequested
,
2832 const gmx::MdModulesNotifier
& mdModulesNotifier
)
2834 CheckpointHeaderContents headerContents
;
2837 /* Read the state from the checkpoint file */
2838 read_checkpoint(fn
, logfio
, cr
, dd_nc
, ir
->eI
, &(ir
->fepvals
->init_fep_state
), &headerContents
,
2839 state
, observablesHistory
, reproducibilityRequested
, mdModulesNotifier
);
2843 gmx_bcast(sizeof(headerContents
.step
), &headerContents
.step
, cr
);
2844 gmx::MdModulesCheckpointReadingBroadcast broadcastCheckPointData
= { *cr
, headerContents
.file_version
};
2845 mdModulesNotifier
.notifier_
.notify(broadcastCheckPointData
);
2847 ir
->bContinuation
= TRUE
;
2848 if (ir
->nsteps
>= 0)
2850 // TODO Should the following condition be <=? Currently if you
2851 // pass a checkpoint written by an normal completion to a restart,
2852 // mdrun will read all input, does some work but no steps, and
2853 // write successful output. But perhaps that is not desirable.
2854 // Note that we do not intend to support the use of mdrun
2855 // -nsteps to circumvent this condition.
2856 if (ir
->nsteps
+ ir
->init_step
< headerContents
.step
)
2858 char buf
[STEPSTRSIZE
];
2859 std::string message
=
2860 gmx::formatString("The input requested %s steps, ", gmx_step_str(ir
->nsteps
, buf
));
2861 if (ir
->init_step
> 0)
2863 message
+= gmx::formatString("starting from step %s, ", gmx_step_str(ir
->init_step
, buf
));
2865 message
+= gmx::formatString(
2866 "however the checkpoint "
2867 "file has already reached step %s. The simulation will not "
2868 "proceed, because either your simulation is already complete, "
2869 "or your combination of input files don't match.",
2870 gmx_step_str(headerContents
.step
, buf
));
2871 gmx_fatal(FARGS
, "%s", message
.c_str());
2873 ir
->nsteps
+= ir
->init_step
- headerContents
.step
;
2875 ir
->init_step
= headerContents
.step
;
2876 ir
->simulation_part
= headerContents
.simulation_part
+ 1;
2879 void read_checkpoint_part_and_step(const char* filename
, int* simulation_part
, int64_t* step
)
2883 if (filename
== nullptr || !gmx_fexist(filename
) || ((fp
= gmx_fio_open(filename
, "r")) == nullptr))
2885 *simulation_part
= 0;
2890 CheckpointHeaderContents headerContents
;
2891 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, nullptr, &headerContents
);
2893 *simulation_part
= headerContents
.simulation_part
;
2894 *step
= headerContents
.step
;
2897 static CheckpointHeaderContents
read_checkpoint_data(t_fileio
* fp
,
2899 std::vector
<gmx_file_position_t
>* outputfiles
)
2901 CheckpointHeaderContents headerContents
;
2902 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, nullptr, &headerContents
);
2903 state
->natoms
= headerContents
.natoms
;
2904 state
->ngtc
= headerContents
.ngtc
;
2905 state
->nnhpres
= headerContents
.nnhpres
;
2906 state
->nhchainlength
= headerContents
.nhchainlength
;
2907 state
->flags
= headerContents
.flags_state
;
2908 int ret
= do_cpt_state(gmx_fio_getxdr(fp
), state
->flags
, state
, nullptr);
2913 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), headerContents
.flags_eks
, &state
->ekinstate
, nullptr);
2919 energyhistory_t enerhist
;
2920 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flags_enh
, &enerhist
, nullptr);
2925 PullHistory pullHist
= {};
2926 ret
= doCptPullHist(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flagsPullHistory
, &pullHist
,
2927 StatePart::pullHistory
, nullptr);
2933 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
), headerContents
.flags_dfh
, headerContents
.nlambda
,
2934 &state
->dfhist
, nullptr);
2940 edsamhistory_t edsamhist
= {};
2941 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
.nED
, &edsamhist
, nullptr);
2947 ret
= do_cpt_awh(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flags_awhh
, state
->awhHistory
.get(), nullptr);
2953 swaphistory_t swaphist
= {};
2954 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
.eSwapCoords
, &swaphist
, nullptr);
2960 ret
= do_cpt_files(gmx_fio_getxdr(fp
), TRUE
, outputfiles
, nullptr, headerContents
.file_version
);
2966 gmx::MdModulesNotifier mdModuleNotifier
;
2967 do_cpt_mdmodules(headerContents
.file_version
, fp
, mdModuleNotifier
);
2968 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), headerContents
.file_version
);
2973 return headerContents
;
2976 void read_checkpoint_trxframe(t_fileio
* fp
, t_trxframe
* fr
)
2979 std::vector
<gmx_file_position_t
> outputfiles
;
2980 CheckpointHeaderContents headerContents
= read_checkpoint_data(fp
, &state
, &outputfiles
);
2982 fr
->natoms
= state
.natoms
;
2984 fr
->step
= int64_to_int(headerContents
.step
, "conversion of checkpoint to trajectory");
2986 fr
->time
= headerContents
.t
;
2988 fr
->lambda
= state
.lambda
[efptFEP
];
2989 fr
->fep_state
= state
.fep_state
;
2991 fr
->bX
= ((state
.flags
& (1 << estX
)) != 0);
2994 fr
->x
= makeRvecArray(state
.x
, state
.natoms
);
2996 fr
->bV
= ((state
.flags
& (1 << estV
)) != 0);
2999 fr
->v
= makeRvecArray(state
.v
, state
.natoms
);
3002 fr
->bBox
= ((state
.flags
& (1 << estBOX
)) != 0);
3005 copy_mat(state
.box
, fr
->box
);
3009 void list_checkpoint(const char* fn
, FILE* out
)
3016 fp
= gmx_fio_open(fn
, "r");
3017 CheckpointHeaderContents headerContents
;
3018 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, out
, &headerContents
);
3019 state
.natoms
= headerContents
.natoms
;
3020 state
.ngtc
= headerContents
.ngtc
;
3021 state
.nnhpres
= headerContents
.nnhpres
;
3022 state
.nhchainlength
= headerContents
.nhchainlength
;
3023 state
.flags
= headerContents
.flags_state
;
3024 ret
= do_cpt_state(gmx_fio_getxdr(fp
), state
.flags
, &state
, out
);
3029 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), headerContents
.flags_eks
, &state
.ekinstate
, out
);
3035 energyhistory_t enerhist
;
3036 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flags_enh
, &enerhist
, out
);
3040 PullHistory pullHist
= {};
3041 ret
= doCptPullHist(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flagsPullHistory
, &pullHist
,
3042 StatePart::pullHistory
, out
);
3047 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
), headerContents
.flags_dfh
, headerContents
.nlambda
,
3048 &state
.dfhist
, out
);
3053 edsamhistory_t edsamhist
= {};
3054 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
.nED
, &edsamhist
, out
);
3059 ret
= do_cpt_awh(gmx_fio_getxdr(fp
), TRUE
, headerContents
.flags_awhh
, state
.awhHistory
.get(), out
);
3064 swaphistory_t swaphist
= {};
3065 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, headerContents
.eSwapCoords
, &swaphist
, out
);
3070 std::vector
<gmx_file_position_t
> outputfiles
;
3071 ret
= do_cpt_files(gmx_fio_getxdr(fp
), TRUE
, &outputfiles
, out
, headerContents
.file_version
);
3076 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), headerContents
.file_version
);
3083 if (gmx_fio_close(fp
) != 0)
3085 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
3089 /* This routine cannot print tons of data, since it is called before the log file is opened. */
3090 CheckpointHeaderContents
read_checkpoint_simulation_part_and_filenames(t_fileio
* fp
,
3091 std::vector
<gmx_file_position_t
>* outputfiles
)
3094 CheckpointHeaderContents headerContents
= read_checkpoint_data(fp
, &state
, outputfiles
);
3095 if (gmx_fio_close(fp
) != 0)
3097 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
3099 return headerContents
;