2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2008,2009,2010,2011,2012,2013,2014,2015, 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 "gromacs/legacyheaders/checkpoint.h"
49 #ifdef GMX_NATIVE_WINDOWS
51 #include <sys/locking.h>
54 #include "buildinfo.h"
55 #include "gromacs/fileio/filenm.h"
56 #include "gromacs/fileio/gmxfio.h"
57 #include "gromacs/fileio/gmxfio-xdr.h"
58 #include "gromacs/fileio/trx.h"
59 #include "gromacs/fileio/xdr_datatype.h"
60 #include "gromacs/fileio/xdrf.h"
61 #include "gromacs/gmxlib/df_history.h"
62 #include "gromacs/gmxlib/energyhistory.h"
63 #include "gromacs/legacyheaders/copyrite.h"
64 #include "gromacs/legacyheaders/names.h"
65 #include "gromacs/legacyheaders/network.h"
66 #include "gromacs/legacyheaders/txtdump.h"
67 #include "gromacs/legacyheaders/typedefs.h"
68 #include "gromacs/legacyheaders/types/commrec.h"
69 #include "gromacs/math/vec.h"
70 #include "gromacs/utility/baseversion.h"
71 #include "gromacs/utility/cstringutil.h"
72 #include "gromacs/utility/fatalerror.h"
73 #include "gromacs/utility/futil.h"
74 #include "gromacs/utility/smalloc.h"
75 #include "gromacs/utility/sysinfo.h"
81 #define CPT_MAGIC1 171817
82 #define CPT_MAGIC2 171819
83 #define CPTSTRLEN 1024
86 #define GMX_CPT_BUILD_DP 1
88 #define GMX_CPT_BUILD_DP 0
91 /* cpt_version should normally only be changed
92 * when the header of footer format changes.
93 * The state data format itself is backward and forward compatible.
94 * But old code can not read a new entry that is present in the file
95 * (but can read a new format when new entries are not present).
97 static const int cpt_version
= 16;
100 const char *est_names
[estNR
] =
103 "box", "box-rel", "box-v", "pres_prev",
104 "nosehoover-xi", "thermostat-integral",
105 "x", "v", "SDx", "CGp", "LD-rng", "LD-rng-i",
106 "disre_initf", "disre_rm3tav",
107 "orire_initf", "orire_Dtav",
108 "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng", "MC-rng-i"
112 eeksEKIN_N
, eeksEKINH
, eeksDEKINDL
, eeksMVCOS
, eeksEKINF
, eeksEKINO
, eeksEKINSCALEF
, eeksEKINSCALEH
, eeksVSCALE
, eeksEKINTOTAL
, eeksNR
115 const char *eeks_names
[eeksNR
] =
117 "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
118 "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
122 eenhENERGY_N
, eenhENERGY_AVER
, eenhENERGY_SUM
, eenhENERGY_NSUM
,
123 eenhENERGY_SUM_SIM
, eenhENERGY_NSUM_SIM
,
124 eenhENERGY_NSTEPS
, eenhENERGY_NSTEPS_SIM
,
125 eenhENERGY_DELTA_H_NN
,
126 eenhENERGY_DELTA_H_LIST
,
127 eenhENERGY_DELTA_H_STARTTIME
,
128 eenhENERGY_DELTA_H_STARTLAMBDA
,
132 const char *eenh_names
[eenhNR
] =
134 "energy_n", "energy_aver", "energy_sum", "energy_nsum",
135 "energy_sum_sim", "energy_nsum_sim",
136 "energy_nsteps", "energy_nsteps_sim",
138 "energy_delta_h_list",
139 "energy_delta_h_start_time",
140 "energy_delta_h_start_lambda"
143 /* free energy history variables -- need to be preserved over checkpoint */
145 edfhBEQUIL
, edfhNATLAMBDA
, edfhWLHISTO
, edfhWLDELTA
, edfhSUMWEIGHTS
, edfhSUMDG
, edfhSUMMINVAR
, edfhSUMVAR
,
146 edfhACCUMP
, edfhACCUMM
, edfhACCUMP2
, edfhACCUMM2
, edfhTIJ
, edfhTIJEMP
, edfhNR
148 /* free energy history variable names */
149 const char *edfh_names
[edfhNR
] =
151 "bEquilibrated", "N_at_state", "Wang-Landau Histogram", "Wang-Landau Delta", "Weights", "Free Energies", "minvar", "variance",
152 "accumulated_plus", "accumulated_minus", "accumulated_plus_2", "accumulated_minus_2", "Tij", "Tij_empirical"
156 ecprREAL
, ecprRVEC
, ecprMATRIX
160 cptpEST
, cptpEEKS
, cptpEENH
, cptpEDFH
162 /* enums for the different components of checkpoint variables, replacing the hard coded ones.
163 cptpEST - state variables.
164 cptpEEKS - Kinetic energy state variables.
165 cptpEENH - Energy history state variables.
166 cptpEDFH - free energy history variables.
170 static const char *st_names(int cptp
, int ecpt
)
174 case cptpEST
: return est_names
[ecpt
]; break;
175 case cptpEEKS
: return eeks_names
[ecpt
]; break;
176 case cptpEENH
: return eenh_names
[ecpt
]; break;
177 case cptpEDFH
: return edfh_names
[ecpt
]; break;
183 static void cp_warning(FILE *fp
)
185 fprintf(fp
, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
188 static void cp_error()
190 gmx_fatal(FARGS
, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
193 static void do_cpt_string_err(XDR
*xd
, gmx_bool bRead
, const char *desc
, char **s
, FILE *list
)
199 if (xdr_string(xd
, s
, CPTSTRLEN
) == 0)
205 fprintf(list
, "%s = %s\n", desc
, *s
);
210 static int do_cpt_int(XDR
*xd
, const char *desc
, int *i
, FILE *list
)
212 if (xdr_int(xd
, i
) == 0)
218 fprintf(list
, "%s = %d\n", desc
, *i
);
223 static int do_cpt_u_chars(XDR
*xd
, const char *desc
, int n
, unsigned char *i
, FILE *list
)
227 fprintf(list
, "%s = ", desc
);
230 for (int j
= 0; j
< n
&& res
; j
++)
232 res
&= xdr_u_char(xd
, &i
[j
]);
235 fprintf(list
, "%02x", i
[j
]);
250 static void do_cpt_int_err(XDR
*xd
, const char *desc
, int *i
, FILE *list
)
252 if (do_cpt_int(xd
, desc
, i
, list
) < 0)
258 static void do_cpt_step_err(XDR
*xd
, const char *desc
, gmx_int64_t
*i
, FILE *list
)
260 char buf
[STEPSTRSIZE
];
262 if (xdr_int64(xd
, i
) == 0)
268 fprintf(list
, "%s = %s\n", desc
, gmx_step_str(*i
, buf
));
272 static void do_cpt_double_err(XDR
*xd
, const char *desc
, double *f
, FILE *list
)
274 if (xdr_double(xd
, f
) == 0)
280 fprintf(list
, "%s = %f\n", desc
, *f
);
284 static void do_cpt_real_err(XDR
*xd
, real
*f
)
287 bool_t res
= xdr_double(xd
, f
);
289 bool_t res
= xdr_float(xd
, f
);
297 static void do_cpt_n_rvecs_err(XDR
*xd
, const char *desc
, int n
, rvec f
[], FILE *list
)
299 for (int i
= 0; i
< n
; i
++)
301 for (int j
= 0; j
< DIM
; j
++)
303 do_cpt_real_err(xd
, &f
[i
][j
]);
309 pr_rvecs(list
, 0, desc
, f
, n
);
313 /* If nval >= 0, nval is used; on read this should match the passed value.
314 * If nval n<0, *nptr is used; on read the value is stored in nptr
316 static int do_cpte_reals_low(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
317 int nval
, int *nptr
, real
**v
,
318 FILE *list
, int erealtype
)
322 int dtc
= xdr_datatype_float
;
324 int dtc
= xdr_datatype_double
;
326 real
*vp
, *va
= NULL
;
341 gmx_incons("*ntpr=NULL in do_cpte_reals_low");
346 res
= xdr_int(xd
, &nf
);
357 gmx_fatal(FARGS
, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp
, ecpt
), nval
, nf
);
366 res
= xdr_int(xd
, &dt
);
373 fprintf(stderr
, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
374 st_names(cptp
, ecpt
), xdr_datatype_names
[dtc
],
375 xdr_datatype_names
[dt
]);
377 if (list
|| !(sflags
& (1<<ecpt
)))
390 if (dt
== xdr_datatype_float
)
392 if (dtc
== xdr_datatype_float
)
394 vf
= reinterpret_cast<float *>(vp
);
400 res
= xdr_vector(xd
, reinterpret_cast<char *>(vf
), nf
,
401 static_cast<unsigned int>(sizeof(float)), (xdrproc_t
)xdr_float
);
406 if (dtc
!= xdr_datatype_float
)
408 for (i
= 0; i
< nf
; i
++)
417 if (dtc
== xdr_datatype_double
)
419 /* cppcheck-suppress invalidPointerCast
420 * Only executed if real is anyhow double */
427 res
= xdr_vector(xd
, reinterpret_cast<char *>(vd
), nf
,
428 static_cast<unsigned int>(sizeof(double)), (xdrproc_t
)xdr_double
);
433 if (dtc
!= xdr_datatype_double
)
435 for (i
= 0; i
< nf
; i
++)
448 pr_reals(list
, 0, st_names(cptp
, ecpt
), vp
, nf
);
451 pr_rvecs(list
, 0, st_names(cptp
, ecpt
), (rvec
*)vp
, nf
/3);
454 gmx_incons("Unknown checkpoint real type");
466 /* This function stores n along with the reals for reading,
467 * but on reading it assumes that n matches the value in the checkpoint file,
468 * a fatal error is generated when this is not the case.
470 static int do_cpte_reals(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
471 int n
, real
**v
, FILE *list
)
473 return do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
, n
, NULL
, v
, list
, ecprREAL
);
476 /* This function does the same as do_cpte_reals,
477 * except that on reading it ignores the passed value of *n
478 * and stored the value read from the checkpoint file in *n.
480 static int do_cpte_n_reals(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
481 int *n
, real
**v
, FILE *list
)
483 return do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
, -1, n
, v
, list
, ecprREAL
);
486 static int do_cpte_real(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
489 return do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
, 1, NULL
, &r
, list
, ecprREAL
);
492 static int do_cpte_ints(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
493 int n
, int **v
, FILE *list
)
496 int dtc
= xdr_datatype_int
;
501 res
= xdr_int(xd
, &nf
);
506 if (list
== NULL
&& v
!= NULL
&& nf
!= n
)
508 gmx_fatal(FARGS
, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp
, ecpt
), n
, nf
);
511 res
= xdr_int(xd
, &dt
);
518 gmx_fatal(FARGS
, "Type mismatch for state entry %s, code type is %s, file type is %s\n",
519 st_names(cptp
, ecpt
), xdr_datatype_names
[dtc
],
520 xdr_datatype_names
[dt
]);
522 if (list
|| !(sflags
& (1<<ecpt
)) || v
== NULL
)
535 res
= xdr_vector(xd
, reinterpret_cast<char *>(vp
), nf
,
536 static_cast<unsigned int>(sizeof(int)), (xdrproc_t
)xdr_int
);
543 pr_ivec(list
, 0, st_names(cptp
, ecpt
), vp
, nf
, TRUE
);
553 static int do_cpte_int(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
556 return do_cpte_ints(xd
, cptp
, ecpt
, sflags
, 1, &i
, list
);
559 static int do_cpte_doubles(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
560 int n
, double **v
, FILE *list
)
563 int dtc
= xdr_datatype_double
;
564 double *vp
, *va
= NULL
;
568 res
= xdr_int(xd
, &nf
);
573 if (list
== NULL
&& nf
!= n
)
575 gmx_fatal(FARGS
, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp
, ecpt
), n
, nf
);
578 res
= xdr_int(xd
, &dt
);
585 gmx_fatal(FARGS
, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
586 st_names(cptp
, ecpt
), xdr_datatype_names
[dtc
],
587 xdr_datatype_names
[dt
]);
589 if (list
|| !(sflags
& (1<<ecpt
)))
602 res
= xdr_vector(xd
, reinterpret_cast<char *>(vp
), nf
,
603 static_cast<unsigned int>(sizeof(double)), (xdrproc_t
)xdr_double
);
610 pr_doubles(list
, 0, st_names(cptp
, ecpt
), vp
, nf
);
620 static int do_cpte_double(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
621 double *r
, FILE *list
)
623 return do_cpte_doubles(xd
, cptp
, ecpt
, sflags
, 1, &r
, list
);
627 static int do_cpte_rvecs(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
628 int n
, rvec
**v
, FILE *list
)
630 return do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
,
631 n
*DIM
, NULL
, (real
**)v
, list
, ecprRVEC
);
634 static int do_cpte_matrix(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
635 matrix v
, FILE *list
)
641 ret
= do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
,
642 DIM
*DIM
, NULL
, &vr
, NULL
, ecprMATRIX
);
644 if (list
&& ret
== 0)
646 pr_rvecs(list
, 0, st_names(cptp
, ecpt
), v
, DIM
);
653 static int do_cpte_nmatrix(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
654 int n
, real
**v
, FILE *list
)
658 char name
[CPTSTRLEN
];
665 for (i
= 0; i
< n
; i
++)
667 reti
= do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
, n
, NULL
, &(v
[i
]), NULL
, ecprREAL
);
668 if (list
&& reti
== 0)
670 sprintf(name
, "%s[%d]", st_names(cptp
, ecpt
), i
);
671 pr_reals(list
, 0, name
, v
[i
], n
);
681 static int do_cpte_matrices(XDR
*xd
, int cptp
, int ecpt
, int sflags
,
682 int n
, matrix
**v
, FILE *list
)
685 matrix
*vp
, *va
= NULL
;
691 res
= xdr_int(xd
, &nf
);
696 if (list
== NULL
&& nf
!= n
)
698 gmx_fatal(FARGS
, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp
, ecpt
), n
, nf
);
700 if (list
|| !(sflags
& (1<<ecpt
)))
713 snew(vr
, nf
*DIM
*DIM
);
714 for (i
= 0; i
< nf
; i
++)
716 for (j
= 0; j
< DIM
; j
++)
718 for (k
= 0; k
< DIM
; k
++)
720 vr
[(i
*DIM
+j
)*DIM
+k
] = vp
[i
][j
][k
];
724 ret
= do_cpte_reals_low(xd
, cptp
, ecpt
, sflags
,
725 nf
*DIM
*DIM
, NULL
, &vr
, NULL
, ecprMATRIX
);
726 for (i
= 0; i
< nf
; i
++)
728 for (j
= 0; j
< DIM
; j
++)
730 for (k
= 0; k
< DIM
; k
++)
732 vp
[i
][j
][k
] = vr
[(i
*DIM
+j
)*DIM
+k
];
738 if (list
&& ret
== 0)
740 for (i
= 0; i
< nf
; i
++)
742 pr_rvecs(list
, 0, st_names(cptp
, ecpt
), vp
[i
], DIM
);
753 static void do_cpt_header(XDR
*xd
, gmx_bool bRead
, int *file_version
,
754 char **version
, char **btime
, char **buser
, char **bhost
,
756 char **fprog
, char **ftime
,
757 int *eIntegrator
, int *simulation_part
,
758 gmx_int64_t
*step
, double *t
,
759 int *nnodes
, int *dd_nc
, int *npme
,
760 int *natoms
, int *ngtc
, int *nnhpres
, int *nhchainlength
,
761 int *nlambda
, int *flags_state
,
762 int *flags_eks
, int *flags_enh
, int *flags_dfh
,
763 int *nED
, int *eSwapCoords
,
779 res
= xdr_int(xd
, &magic
);
782 gmx_fatal(FARGS
, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
784 if (magic
!= CPT_MAGIC1
)
786 gmx_fatal(FARGS
, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
787 "The checkpoint file is corrupted or not a checkpoint file",
793 gmx_gethostname(fhost
, 255);
795 do_cpt_string_err(xd
, bRead
, "GROMACS version", version
, list
);
796 do_cpt_string_err(xd
, bRead
, "GROMACS build time", btime
, list
);
797 do_cpt_string_err(xd
, bRead
, "GROMACS build user", buser
, list
);
798 do_cpt_string_err(xd
, bRead
, "GROMACS build host", bhost
, list
);
799 do_cpt_string_err(xd
, bRead
, "generating program", fprog
, list
);
800 do_cpt_string_err(xd
, bRead
, "generation time", ftime
, list
);
801 *file_version
= cpt_version
;
802 do_cpt_int_err(xd
, "checkpoint file version", file_version
, list
);
803 if (*file_version
> cpt_version
)
805 gmx_fatal(FARGS
, "Attempting to read a checkpoint file of version %d with code of version %d\n", *file_version
, cpt_version
);
807 if (*file_version
>= 13)
809 do_cpt_int_err(xd
, "GROMACS double precision", double_prec
, list
);
815 if (*file_version
>= 12)
817 do_cpt_string_err(xd
, bRead
, "generating host", &fhost
, list
);
823 do_cpt_int_err(xd
, "#atoms", natoms
, list
);
824 do_cpt_int_err(xd
, "#T-coupling groups", ngtc
, list
);
825 if (*file_version
>= 10)
827 do_cpt_int_err(xd
, "#Nose-Hoover T-chains", nhchainlength
, list
);
833 if (*file_version
>= 11)
835 do_cpt_int_err(xd
, "#Nose-Hoover T-chains for barostat ", nnhpres
, list
);
841 if (*file_version
>= 14)
843 do_cpt_int_err(xd
, "# of total lambda states ", nlambda
, list
);
849 do_cpt_int_err(xd
, "integrator", eIntegrator
, list
);
850 if (*file_version
>= 3)
852 do_cpt_int_err(xd
, "simulation part #", simulation_part
, list
);
856 *simulation_part
= 1;
858 if (*file_version
>= 5)
860 do_cpt_step_err(xd
, "step", step
, list
);
864 do_cpt_int_err(xd
, "step", &idum
, list
);
867 do_cpt_double_err(xd
, "t", t
, list
);
868 do_cpt_int_err(xd
, "#PP-ranks", nnodes
, list
);
870 do_cpt_int_err(xd
, "dd_nc[x]", dd_nc
? &(dd_nc
[0]) : &idum
, list
);
871 do_cpt_int_err(xd
, "dd_nc[y]", dd_nc
? &(dd_nc
[1]) : &idum
, list
);
872 do_cpt_int_err(xd
, "dd_nc[z]", dd_nc
? &(dd_nc
[2]) : &idum
, list
);
873 do_cpt_int_err(xd
, "#PME-only ranks", npme
, list
);
874 do_cpt_int_err(xd
, "state flags", flags_state
, list
);
875 if (*file_version
>= 4)
877 do_cpt_int_err(xd
, "ekin data flags", flags_eks
, list
);
878 do_cpt_int_err(xd
, "energy history flags", flags_enh
, list
);
883 *flags_enh
= (*flags_state
>> (estORIRE_DTAV
+1));
884 *flags_state
= (*flags_state
& ~((1<<(estORIRE_DTAV
+1)) |
885 (1<<(estORIRE_DTAV
+2)) |
886 (1<<(estORIRE_DTAV
+3))));
888 if (*file_version
>= 14)
890 do_cpt_int_err(xd
, "df history flags", flags_dfh
, list
);
897 if (*file_version
>= 15)
899 do_cpt_int_err(xd
, "ED data sets", nED
, list
);
905 if (*file_version
>= 16)
907 do_cpt_int_err(xd
, "swap", eSwapCoords
, list
);
911 static int do_cpt_footer(XDR
*xd
, int file_version
)
916 if (file_version
>= 2)
919 res
= xdr_int(xd
, &magic
);
924 if (magic
!= CPT_MAGIC2
)
933 static int do_cpt_state(XDR
*xd
, gmx_bool bRead
,
934 int fflags
, t_state
*state
,
944 nnht
= state
->nhchainlength
*state
->ngtc
;
945 nnhtp
= state
->nhchainlength
*state
->nnhpres
;
947 if (bRead
) /* we need to allocate space for dfhist if we are reading */
949 init_df_history(&state
->dfhist
, state
->dfhist
.nlambda
);
952 sflags
= state
->flags
;
953 for (i
= 0; (i
< estNR
&& ret
== 0); i
++)
959 case estLAMBDA
: ret
= do_cpte_reals(xd
, cptpEST
, i
, sflags
, efptNR
, &(state
->lambda
), list
); break;
960 case estFEPSTATE
: ret
= do_cpte_int (xd
, cptpEST
, i
, sflags
, &state
->fep_state
, list
); break;
961 case estBOX
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->box
, list
); break;
962 case estBOX_REL
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->box_rel
, list
); break;
963 case estBOXV
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->boxv
, list
); break;
964 case estPRES_PREV
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->pres_prev
, list
); break;
965 case estSVIR_PREV
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->svir_prev
, list
); break;
966 case estFVIR_PREV
: ret
= do_cpte_matrix(xd
, cptpEST
, i
, sflags
, state
->fvir_prev
, list
); break;
967 case estNH_XI
: ret
= do_cpte_doubles(xd
, cptpEST
, i
, sflags
, nnht
, &state
->nosehoover_xi
, list
); break;
968 case estNH_VXI
: ret
= do_cpte_doubles(xd
, cptpEST
, i
, sflags
, nnht
, &state
->nosehoover_vxi
, list
); break;
969 case estNHPRES_XI
: ret
= do_cpte_doubles(xd
, cptpEST
, i
, sflags
, nnhtp
, &state
->nhpres_xi
, list
); break;
970 case estNHPRES_VXI
: ret
= do_cpte_doubles(xd
, cptpEST
, i
, sflags
, nnhtp
, &state
->nhpres_vxi
, list
); break;
971 case estTC_INT
: ret
= do_cpte_doubles(xd
, cptpEST
, i
, sflags
, state
->ngtc
, &state
->therm_integral
, list
); break;
972 case estVETA
: ret
= do_cpte_real(xd
, cptpEST
, i
, sflags
, &state
->veta
, list
); break;
973 case estVOL0
: ret
= do_cpte_real(xd
, cptpEST
, i
, sflags
, &state
->vol0
, list
); break;
974 case estX
: ret
= do_cpte_rvecs(xd
, cptpEST
, i
, sflags
, state
->natoms
, &state
->x
, list
); break;
975 case estV
: ret
= do_cpte_rvecs(xd
, cptpEST
, i
, sflags
, state
->natoms
, &state
->v
, list
); break;
976 case estSDX
: ret
= do_cpte_rvecs(xd
, cptpEST
, i
, sflags
, state
->natoms
, &state
->sd_X
, list
); break;
977 /* The RNG entries are no longer written,
978 * the next 4 lines are only for reading old files.
980 case estLD_RNG
: ret
= do_cpte_ints(xd
, cptpEST
, i
, sflags
, 0, NULL
, list
); break;
981 case estLD_RNGI
: ret
= do_cpte_ints(xd
, cptpEST
, i
, sflags
, 0, NULL
, list
); break;
982 case estMC_RNG
: ret
= do_cpte_ints(xd
, cptpEST
, i
, sflags
, 0, NULL
, list
); break;
983 case estMC_RNGI
: ret
= do_cpte_ints(xd
, cptpEST
, i
, sflags
, 0, NULL
, list
); break;
984 case estDISRE_INITF
: ret
= do_cpte_real (xd
, cptpEST
, i
, sflags
, &state
->hist
.disre_initf
, list
); break;
985 case estDISRE_RM3TAV
: ret
= do_cpte_n_reals(xd
, cptpEST
, i
, sflags
, &state
->hist
.ndisrepairs
, &state
->hist
.disre_rm3tav
, list
); break;
986 case estORIRE_INITF
: ret
= do_cpte_real (xd
, cptpEST
, i
, sflags
, &state
->hist
.orire_initf
, list
); break;
987 case estORIRE_DTAV
: ret
= do_cpte_n_reals(xd
, cptpEST
, i
, sflags
, &state
->hist
.norire_Dtav
, &state
->hist
.orire_Dtav
, list
); break;
989 gmx_fatal(FARGS
, "Unknown state entry %d\n"
990 "You are probably reading a new checkpoint file with old code", i
);
998 static int do_cpt_ekinstate(XDR
*xd
, int fflags
, ekinstate_t
*ekins
,
1006 for (i
= 0; (i
< eeksNR
&& ret
== 0); i
++)
1008 if (fflags
& (1<<i
))
1013 case eeksEKIN_N
: ret
= do_cpte_int(xd
, cptpEEKS
, i
, fflags
, &ekins
->ekin_n
, list
); break;
1014 case eeksEKINH
: ret
= do_cpte_matrices(xd
, cptpEEKS
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinh
, list
); break;
1015 case eeksEKINF
: ret
= do_cpte_matrices(xd
, cptpEEKS
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinf
, list
); break;
1016 case eeksEKINO
: ret
= do_cpte_matrices(xd
, cptpEEKS
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinh_old
, list
); break;
1017 case eeksEKINTOTAL
: ret
= do_cpte_matrix(xd
, cptpEEKS
, i
, fflags
, ekins
->ekin_total
, list
); break;
1018 case eeksEKINSCALEF
: ret
= do_cpte_doubles(xd
, cptpEEKS
, i
, fflags
, ekins
->ekin_n
, &ekins
->ekinscalef_nhc
, list
); break;
1019 case eeksVSCALE
: ret
= do_cpte_doubles(xd
, 1, cptpEEKS
, fflags
, ekins
->ekin_n
, &ekins
->vscale_nhc
, list
); break;
1020 case eeksEKINSCALEH
: ret
= do_cpte_doubles(xd
, 1, cptpEEKS
, fflags
, ekins
->ekin_n
, &ekins
->ekinscaleh_nhc
, list
); break;
1021 case eeksDEKINDL
: ret
= do_cpte_real(xd
, 1, cptpEEKS
, fflags
, &ekins
->dekindl
, list
); break;
1022 case eeksMVCOS
: ret
= do_cpte_real(xd
, 1, cptpEEKS
, fflags
, &ekins
->mvcos
, list
); break;
1024 gmx_fatal(FARGS
, "Unknown ekin data state entry %d\n"
1025 "You are probably reading a new checkpoint file with old code", i
);
1034 static int do_cpt_swapstate(XDR
*xd
, gmx_bool bRead
, swapstate_t
*swapstate
, FILE *list
)
1038 int swap_cpt_version
= 1;
1041 if (eswapNO
== swapstate
->eSwapCoords
)
1046 swapstate
->bFromCpt
= bRead
;
1048 do_cpt_int_err(xd
, "swap checkpoint version", &swap_cpt_version
, list
);
1049 do_cpt_int_err(xd
, "swap coupling steps", &swapstate
->nAverage
, list
);
1051 /* When reading, init_swapcoords has not been called yet,
1052 * so we have to allocate memory first. */
1054 for (ic
= 0; ic
< eCompNR
; ic
++)
1056 for (ii
= 0; ii
< eIonNR
; ii
++)
1060 do_cpt_int_err(xd
, "swap requested atoms", &swapstate
->nat_req
[ic
][ii
], list
);
1064 do_cpt_int_err(xd
, "swap requested atoms p", swapstate
->nat_req_p
[ic
][ii
], list
);
1069 do_cpt_int_err(xd
, "swap influx netto", &swapstate
->inflow_netto
[ic
][ii
], list
);
1073 do_cpt_int_err(xd
, "swap influx netto p", swapstate
->inflow_netto_p
[ic
][ii
], list
);
1076 if (bRead
&& (NULL
== swapstate
->nat_past
[ic
][ii
]) )
1078 snew(swapstate
->nat_past
[ic
][ii
], swapstate
->nAverage
);
1081 for (j
= 0; j
< swapstate
->nAverage
; j
++)
1085 do_cpt_int_err(xd
, "swap past atom counts", &swapstate
->nat_past
[ic
][ii
][j
], list
);
1089 do_cpt_int_err(xd
, "swap past atom counts p", &swapstate
->nat_past_p
[ic
][ii
][j
], list
);
1095 /* Ion flux per channel */
1096 for (ic
= 0; ic
< eChanNR
; ic
++)
1098 for (ii
= 0; ii
< eIonNR
; ii
++)
1102 do_cpt_int_err(xd
, "channel flux", &swapstate
->fluxfromAtoB
[ic
][ii
], list
);
1106 do_cpt_int_err(xd
, "channel flux p", swapstate
->fluxfromAtoB_p
[ic
][ii
], list
);
1111 /* Ion flux leakage */
1114 snew(swapstate
->fluxleak
, 1);
1116 do_cpt_int_err(xd
, "flux leakage", swapstate
->fluxleak
, list
);
1119 do_cpt_int_err(xd
, "number of ions", &swapstate
->nions
, list
);
1123 snew(swapstate
->channel_label
, swapstate
->nions
);
1124 snew(swapstate
->comp_from
, swapstate
->nions
);
1127 do_cpt_u_chars(xd
, "channel history", swapstate
->nions
, swapstate
->channel_label
, list
);
1128 do_cpt_u_chars(xd
, "domain history", swapstate
->nions
, swapstate
->comp_from
, list
);
1130 /* Save the last known whole positions to checkpoint
1131 * file to be able to also make multimeric channels whole in PBC */
1132 do_cpt_int_err(xd
, "Ch0 atoms", &swapstate
->nat
[eChan0
], list
);
1133 do_cpt_int_err(xd
, "Ch1 atoms", &swapstate
->nat
[eChan1
], list
);
1136 snew(swapstate
->xc_old_whole
[eChan0
], swapstate
->nat
[eChan0
]);
1137 snew(swapstate
->xc_old_whole
[eChan1
], swapstate
->nat
[eChan1
]);
1138 do_cpt_n_rvecs_err(xd
, "Ch0 whole x", swapstate
->nat
[eChan0
], swapstate
->xc_old_whole
[eChan0
], list
);
1139 do_cpt_n_rvecs_err(xd
, "Ch1 whole x", swapstate
->nat
[eChan1
], swapstate
->xc_old_whole
[eChan1
], list
);
1143 do_cpt_n_rvecs_err(xd
, "Ch0 whole x", swapstate
->nat
[eChan0
], *swapstate
->xc_old_whole_p
[eChan0
], list
);
1144 do_cpt_n_rvecs_err(xd
, "Ch1 whole x", swapstate
->nat
[eChan1
], *swapstate
->xc_old_whole_p
[eChan1
], list
);
1151 static int do_cpt_enerhist(XDR
*xd
, gmx_bool bRead
,
1152 int fflags
, energyhistory_t
*enerhist
,
1163 enerhist
->nsteps
= 0;
1165 enerhist
->nsteps_sim
= 0;
1166 enerhist
->nsum_sim
= 0;
1167 enerhist
->dht
= NULL
;
1169 if (fflags
& (1<< eenhENERGY_DELTA_H_NN
) )
1171 snew(enerhist
->dht
, 1);
1172 enerhist
->dht
->ndh
= NULL
;
1173 enerhist
->dht
->dh
= NULL
;
1174 enerhist
->dht
->start_lambda_set
= FALSE
;
1178 for (i
= 0; (i
< eenhNR
&& ret
== 0); i
++)
1180 if (fflags
& (1<<i
))
1184 case eenhENERGY_N
: ret
= do_cpte_int(xd
, cptpEENH
, i
, fflags
, &enerhist
->nener
, list
); break;
1185 case eenhENERGY_AVER
: ret
= do_cpte_doubles(xd
, cptpEENH
, i
, fflags
, enerhist
->nener
, &enerhist
->ener_ave
, list
); break;
1186 case eenhENERGY_SUM
: ret
= do_cpte_doubles(xd
, cptpEENH
, i
, fflags
, enerhist
->nener
, &enerhist
->ener_sum
, list
); break;
1187 case eenhENERGY_NSUM
: do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsum
, list
); break;
1188 case eenhENERGY_SUM_SIM
: ret
= do_cpte_doubles(xd
, cptpEENH
, i
, fflags
, enerhist
->nener
, &enerhist
->ener_sum_sim
, list
); break;
1189 case eenhENERGY_NSUM_SIM
: do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsum_sim
, list
); break;
1190 case eenhENERGY_NSTEPS
: do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsteps
, list
); break;
1191 case eenhENERGY_NSTEPS_SIM
: do_cpt_step_err(xd
, eenh_names
[i
], &enerhist
->nsteps_sim
, list
); break;
1192 case eenhENERGY_DELTA_H_NN
: do_cpt_int_err(xd
, eenh_names
[i
], &(enerhist
->dht
->nndh
), list
);
1193 if (bRead
) /* now allocate memory for it */
1195 snew(enerhist
->dht
->dh
, enerhist
->dht
->nndh
);
1196 snew(enerhist
->dht
->ndh
, enerhist
->dht
->nndh
);
1197 for (j
= 0; j
< enerhist
->dht
->nndh
; j
++)
1199 enerhist
->dht
->ndh
[j
] = 0;
1200 enerhist
->dht
->dh
[j
] = NULL
;
1204 case eenhENERGY_DELTA_H_LIST
:
1205 for (j
= 0; j
< enerhist
->dht
->nndh
; j
++)
1207 ret
= do_cpte_n_reals(xd
, cptpEENH
, i
, fflags
, &enerhist
->dht
->ndh
[j
], &(enerhist
->dht
->dh
[j
]), list
);
1210 case eenhENERGY_DELTA_H_STARTTIME
:
1211 ret
= do_cpte_double(xd
, cptpEENH
, i
, fflags
, &(enerhist
->dht
->start_time
), list
); break;
1212 case eenhENERGY_DELTA_H_STARTLAMBDA
:
1213 ret
= do_cpte_double(xd
, cptpEENH
, i
, fflags
, &(enerhist
->dht
->start_lambda
), list
); break;
1215 gmx_fatal(FARGS
, "Unknown energy history entry %d\n"
1216 "You are probably reading a new checkpoint file with old code", i
);
1221 if ((fflags
& (1<<eenhENERGY_SUM
)) && !(fflags
& (1<<eenhENERGY_SUM_SIM
)))
1223 /* Assume we have an old file format and copy sum to sum_sim */
1224 srenew(enerhist
->ener_sum_sim
, enerhist
->nener
);
1225 for (i
= 0; i
< enerhist
->nener
; i
++)
1227 enerhist
->ener_sum_sim
[i
] = enerhist
->ener_sum
[i
];
1231 if ( (fflags
& (1<<eenhENERGY_NSUM
)) &&
1232 !(fflags
& (1<<eenhENERGY_NSTEPS
)))
1234 /* Assume we have an old file format and copy nsum to nsteps */
1235 enerhist
->nsteps
= enerhist
->nsum
;
1237 if ( (fflags
& (1<<eenhENERGY_NSUM_SIM
)) &&
1238 !(fflags
& (1<<eenhENERGY_NSTEPS_SIM
)))
1240 /* Assume we have an old file format and copy nsum to nsteps */
1241 enerhist
->nsteps_sim
= enerhist
->nsum_sim
;
1247 static int do_cpt_df_hist(XDR
*xd
, int fflags
, df_history_t
*dfhist
, FILE *list
)
1252 nlambda
= dfhist
->nlambda
;
1255 for (i
= 0; (i
< edfhNR
&& ret
== 0); i
++)
1257 if (fflags
& (1<<i
))
1261 case edfhBEQUIL
: ret
= do_cpte_int(xd
, cptpEDFH
, i
, fflags
, &dfhist
->bEquil
, list
); break;
1262 case edfhNATLAMBDA
: ret
= do_cpte_ints(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->n_at_lam
, list
); break;
1263 case edfhWLHISTO
: ret
= do_cpte_reals(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->wl_histo
, list
); break;
1264 case edfhWLDELTA
: ret
= do_cpte_real(xd
, cptpEDFH
, i
, fflags
, &dfhist
->wl_delta
, list
); break;
1265 case edfhSUMWEIGHTS
: ret
= do_cpte_reals(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->sum_weights
, list
); break;
1266 case edfhSUMDG
: ret
= do_cpte_reals(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->sum_dg
, list
); break;
1267 case edfhSUMMINVAR
: ret
= do_cpte_reals(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->sum_minvar
, list
); break;
1268 case edfhSUMVAR
: ret
= do_cpte_reals(xd
, cptpEDFH
, i
, fflags
, nlambda
, &dfhist
->sum_variance
, list
); break;
1269 case edfhACCUMP
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->accum_p
, list
); break;
1270 case edfhACCUMM
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->accum_m
, list
); break;
1271 case edfhACCUMP2
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->accum_p2
, list
); break;
1272 case edfhACCUMM2
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->accum_m2
, list
); break;
1273 case edfhTIJ
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->Tij
, list
); break;
1274 case edfhTIJEMP
: ret
= do_cpte_nmatrix(xd
, cptpEDFH
, i
, fflags
, nlambda
, dfhist
->Tij_empirical
, list
); break;
1277 gmx_fatal(FARGS
, "Unknown df history entry %d\n"
1278 "You are probably reading a new checkpoint file with old code", i
);
1287 /* This function stores the last whole configuration of the reference and
1288 * average structure in the .cpt file
1290 static int do_cpt_EDstate(XDR
*xd
, gmx_bool bRead
,
1291 edsamstate_t
*EDstate
, FILE *list
)
1298 EDstate
->bFromCpt
= bRead
;
1300 if (EDstate
->nED
<= 0)
1305 /* When reading, init_edsam has not been called yet,
1306 * so we have to allocate memory first. */
1309 snew(EDstate
->nref
, EDstate
->nED
);
1310 snew(EDstate
->old_sref
, EDstate
->nED
);
1311 snew(EDstate
->nav
, EDstate
->nED
);
1312 snew(EDstate
->old_sav
, EDstate
->nED
);
1315 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1316 for (i
= 0; i
< EDstate
->nED
; i
++)
1318 /* Reference structure SREF */
1319 sprintf(buf
, "ED%d # of atoms in reference structure", i
+1);
1320 do_cpt_int_err(xd
, buf
, &EDstate
->nref
[i
], list
);
1321 sprintf(buf
, "ED%d x_ref", i
+1);
1324 snew(EDstate
->old_sref
[i
], EDstate
->nref
[i
]);
1325 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nref
[i
], EDstate
->old_sref
[i
], list
);
1329 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nref
[i
], EDstate
->old_sref_p
[i
], list
);
1332 /* Average structure SAV */
1333 sprintf(buf
, "ED%d # of atoms in average structure", i
+1);
1334 do_cpt_int_err(xd
, buf
, &EDstate
->nav
[i
], list
);
1335 sprintf(buf
, "ED%d x_av", i
+1);
1338 snew(EDstate
->old_sav
[i
], EDstate
->nav
[i
]);
1339 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nav
[i
], EDstate
->old_sav
[i
], list
);
1343 do_cpt_n_rvecs_err(xd
, buf
, EDstate
->nav
[i
], EDstate
->old_sav_p
[i
], list
);
1351 static int do_cpt_files(XDR
*xd
, gmx_bool bRead
,
1352 gmx_file_position_t
**p_outputfiles
, int *nfiles
,
1353 FILE *list
, int file_version
)
1357 gmx_off_t mask
= 0xFFFFFFFFL
;
1358 int offset_high
, offset_low
;
1360 gmx_file_position_t
*outputfiles
;
1362 if (do_cpt_int(xd
, "number of output files", nfiles
, list
) != 0)
1369 snew(*p_outputfiles
, *nfiles
);
1372 outputfiles
= *p_outputfiles
;
1374 for (i
= 0; i
< *nfiles
; i
++)
1376 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1379 do_cpt_string_err(xd
, bRead
, "output filename", &buf
, list
);
1380 std::strncpy(outputfiles
[i
].filename
, buf
, CPTSTRLEN
-1);
1386 if (do_cpt_int(xd
, "file_offset_high", &offset_high
, list
) != 0)
1390 if (do_cpt_int(xd
, "file_offset_low", &offset_low
, list
) != 0)
1394 outputfiles
[i
].offset
= (static_cast<gmx_off_t
>(offset_high
) << 32 ) | ( static_cast<gmx_off_t
>(offset_low
) & mask
);
1398 buf
= outputfiles
[i
].filename
;
1399 do_cpt_string_err(xd
, bRead
, "output filename", &buf
, list
);
1401 offset
= outputfiles
[i
].offset
;
1409 offset_low
= static_cast<int>(offset
& mask
);
1410 offset_high
= static_cast<int>((offset
>> 32) & mask
);
1412 if (do_cpt_int(xd
, "file_offset_high", &offset_high
, list
) != 0)
1416 if (do_cpt_int(xd
, "file_offset_low", &offset_low
, list
) != 0)
1421 if (file_version
>= 8)
1423 if (do_cpt_int(xd
, "file_checksum_size", &(outputfiles
[i
].chksum_size
),
1428 if (do_cpt_u_chars(xd
, "file_checksum", 16, outputfiles
[i
].chksum
, list
) != 0)
1435 outputfiles
[i
].chksum_size
= -1;
1442 void write_checkpoint(const char *fn
, gmx_bool bNumberAndKeep
,
1443 FILE *fplog
, t_commrec
*cr
,
1444 int eIntegrator
, int simulation_part
,
1445 gmx_bool bExpanded
, int elamstats
,
1446 gmx_int64_t step
, double t
, t_state
*state
)
1456 char *fntemp
; /* the temporary checkpoint file name */
1457 char timebuf
[STRLEN
];
1458 int nppnodes
, npmenodes
;
1459 char buf
[1024], suffix
[5+STEPSTRSIZE
], sbuf
[STEPSTRSIZE
];
1460 gmx_file_position_t
*outputfiles
;
1463 int flags_eks
, flags_enh
, flags_dfh
;
1466 if (DOMAINDECOMP(cr
))
1468 nppnodes
= cr
->dd
->nnodes
;
1469 npmenodes
= cr
->npmenodes
;
1478 /* make the new temporary filename */
1479 snew(fntemp
, std::strlen(fn
)+5+STEPSTRSIZE
);
1480 std::strcpy(fntemp
, fn
);
1481 fntemp
[std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1] = '\0';
1482 sprintf(suffix
, "_%s%s", "step", gmx_step_str(step
, sbuf
));
1483 std::strcat(fntemp
, suffix
);
1484 std::strcat(fntemp
, fn
+std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1);
1486 /* if we can't rename, we just overwrite the cpt file.
1487 * dangerous if interrupted.
1489 snew(fntemp
, std::strlen(fn
));
1490 std::strcpy(fntemp
, fn
);
1492 gmx_format_current_time(timebuf
, STRLEN
);
1496 fprintf(fplog
, "Writing checkpoint, step %s at %s\n\n",
1497 gmx_step_str(step
, buf
), timebuf
);
1500 /* Get offsets for open files */
1501 gmx_fio_get_output_file_positions(&outputfiles
, &noutputfiles
);
1503 fp
= gmx_fio_open(fntemp
, "w");
1505 if (state
->ekinstate
.bUpToDate
)
1508 ((1<<eeksEKIN_N
) | (1<<eeksEKINH
) | (1<<eeksEKINF
) |
1509 (1<<eeksEKINO
) | (1<<eeksEKINSCALEF
) | (1<<eeksEKINSCALEH
) |
1510 (1<<eeksVSCALE
) | (1<<eeksDEKINDL
) | (1<<eeksMVCOS
));
1518 if (state
->enerhist
->nsum
> 0 || state
->enerhist
->nsum_sim
> 0)
1520 flags_enh
|= (1<<eenhENERGY_N
) | (1<<eenhENERGY_NSTEPS
) | (1<<eenhENERGY_NSTEPS_SIM
);
1521 if (state
->enerhist
->nsum
> 0)
1523 flags_enh
|= ((1<<eenhENERGY_AVER
) | (1<<eenhENERGY_SUM
) |
1524 (1<<eenhENERGY_NSUM
));
1526 if (state
->enerhist
->nsum_sim
> 0)
1528 flags_enh
|= ((1<<eenhENERGY_SUM_SIM
) | (1<<eenhENERGY_NSUM_SIM
));
1530 if (state
->enerhist
->dht
)
1532 flags_enh
|= ( (1<< eenhENERGY_DELTA_H_NN
) |
1533 (1<< eenhENERGY_DELTA_H_LIST
) |
1534 (1<< eenhENERGY_DELTA_H_STARTTIME
) |
1535 (1<< eenhENERGY_DELTA_H_STARTLAMBDA
) );
1541 flags_dfh
= ((1<<edfhBEQUIL
) | (1<<edfhNATLAMBDA
) | (1<<edfhSUMWEIGHTS
) | (1<<edfhSUMDG
) |
1542 (1<<edfhTIJ
) | (1<<edfhTIJEMP
));
1545 flags_dfh
|= ((1<<edfhWLDELTA
) | (1<<edfhWLHISTO
));
1547 if ((elamstats
== elamstatsMINVAR
) || (elamstats
== elamstatsBARKER
) || (elamstats
== elamstatsMETROPOLIS
))
1549 flags_dfh
|= ((1<<edfhACCUMP
) | (1<<edfhACCUMM
) | (1<<edfhACCUMP2
) | (1<<edfhACCUMM2
)
1550 | (1<<edfhSUMMINVAR
) | (1<<edfhSUMVAR
));
1558 /* We can check many more things now (CPU, acceleration, etc), but
1559 * it is highly unlikely to have two separate builds with exactly
1560 * the same version, user, time, and build host!
1563 version
= gmx_strdup(gmx_version());
1564 btime
= gmx_strdup(BUILD_TIME
);
1565 buser
= gmx_strdup(BUILD_USER
);
1566 bhost
= gmx_strdup(BUILD_HOST
);
1568 double_prec
= GMX_CPT_BUILD_DP
;
1569 fprog
= gmx_strdup(Program());
1571 ftime
= &(timebuf
[0]);
1573 do_cpt_header(gmx_fio_getxdr(fp
), FALSE
, &file_version
,
1574 &version
, &btime
, &buser
, &bhost
, &double_prec
, &fprog
, &ftime
,
1575 &eIntegrator
, &simulation_part
, &step
, &t
, &nppnodes
,
1576 DOMAINDECOMP(cr
) ? cr
->dd
->nc
: NULL
, &npmenodes
,
1577 &state
->natoms
, &state
->ngtc
, &state
->nnhpres
,
1578 &state
->nhchainlength
, &(state
->dfhist
.nlambda
), &state
->flags
, &flags_eks
, &flags_enh
, &flags_dfh
,
1579 &state
->edsamstate
.nED
, &state
->swapstate
.eSwapCoords
,
1588 if ((do_cpt_state(gmx_fio_getxdr(fp
), FALSE
, state
->flags
, state
, NULL
) < 0) ||
1589 (do_cpt_ekinstate(gmx_fio_getxdr(fp
), flags_eks
, &state
->ekinstate
, NULL
) < 0) ||
1590 (do_cpt_enerhist(gmx_fio_getxdr(fp
), FALSE
, flags_enh
, state
->enerhist
, NULL
) < 0) ||
1591 (do_cpt_df_hist(gmx_fio_getxdr(fp
), flags_dfh
, &state
->dfhist
, NULL
) < 0) ||
1592 (do_cpt_EDstate(gmx_fio_getxdr(fp
), FALSE
, &state
->edsamstate
, NULL
) < 0) ||
1593 (do_cpt_swapstate(gmx_fio_getxdr(fp
), FALSE
, &state
->swapstate
, NULL
) < 0) ||
1594 (do_cpt_files(gmx_fio_getxdr(fp
), FALSE
, &outputfiles
, &noutputfiles
, NULL
,
1597 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1600 do_cpt_footer(gmx_fio_getxdr(fp
), file_version
);
1602 /* we really, REALLY, want to make sure to physically write the checkpoint,
1603 and all the files it depends on, out to disk. Because we've
1604 opened the checkpoint with gmx_fio_open(), it's in our list
1606 ret
= gmx_fio_all_output_fsync();
1612 "Cannot fsync '%s'; maybe you are out of disk space?",
1613 gmx_fio_getname(ret
));
1615 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV
) == NULL
)
1625 if (gmx_fio_close(fp
) != 0)
1627 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1630 /* we don't move the checkpoint if the user specified they didn't want it,
1631 or if the fsyncs failed */
1633 if (!bNumberAndKeep
&& !ret
)
1637 /* Rename the previous checkpoint file */
1638 std::strcpy(buf
, fn
);
1639 buf
[std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1] = '\0';
1640 std::strcat(buf
, "_prev");
1641 std::strcat(buf
, fn
+std::strlen(fn
) - std::strlen(ftp2ext(fn2ftp(fn
))) - 1);
1643 /* we copy here so that if something goes wrong between now and
1644 * the rename below, there's always a state.cpt.
1645 * If renames are atomic (such as in POSIX systems),
1646 * this copying should be unneccesary.
1648 gmx_file_copy(fn
, buf
, FALSE
);
1649 /* We don't really care if this fails:
1650 * there's already a new checkpoint.
1653 gmx_file_rename(fn
, buf
);
1656 if (gmx_file_rename(fntemp
, fn
) != 0)
1658 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
1661 #endif /* GMX_NO_RENAME */
1667 /*code for alternate checkpointing scheme. moved from top of loop over
1669 fcRequestCheckPoint();
1670 if (fcCheckPointParallel( cr
->nodeid
, NULL
, 0) == 0)
1672 gmx_fatal( 3, __FILE__
, __LINE__
, "Checkpoint error on step %d\n", step
);
1674 #endif /* end GMX_FAHCORE block */
1677 static void print_flag_mismatch(FILE *fplog
, int sflags
, int fflags
)
1681 fprintf(fplog
, "\nState entry mismatch between the simulation and the checkpoint file\n");
1682 fprintf(fplog
, "Entries which are not present in the checkpoint file will not be updated\n");
1683 fprintf(fplog
, " %24s %11s %11s\n", "", "simulation", "checkpoint");
1684 for (i
= 0; i
< estNR
; i
++)
1686 if ((sflags
& (1<<i
)) || (fflags
& (1<<i
)))
1688 fprintf(fplog
, " %24s %11s %11s\n",
1690 (sflags
& (1<<i
)) ? " present " : "not present",
1691 (fflags
& (1<<i
)) ? " present " : "not present");
1696 static void check_int(FILE *fplog
, const char *type
, int p
, int f
, gmx_bool
*mm
)
1698 FILE *fp
= fplog
? fplog
: stderr
;
1702 fprintf(fp
, " %s mismatch,\n", type
);
1703 fprintf(fp
, " current program: %d\n", p
);
1704 fprintf(fp
, " checkpoint file: %d\n", f
);
1710 static void check_string(FILE *fplog
, const char *type
, const char *p
,
1711 const char *f
, gmx_bool
*mm
)
1713 FILE *fp
= fplog
? fplog
: stderr
;
1715 if (std::strcmp(p
, f
) != 0)
1717 fprintf(fp
, " %s mismatch,\n", type
);
1718 fprintf(fp
, " current program: %s\n", p
);
1719 fprintf(fp
, " checkpoint file: %s\n", f
);
1725 static void check_match(FILE *fplog
,
1727 char *btime
, char *buser
, char *bhost
, int double_prec
,
1729 t_commrec
*cr
, int npp_f
, int npme_f
,
1730 ivec dd_nc
, ivec dd_nc_f
)
1733 gmx_bool mm
= FALSE
;
1734 gmx_bool patchlevel_differs
= FALSE
;
1735 gmx_bool version_differs
= FALSE
;
1737 check_string(fplog
, "Version", gmx_version(), version
, &mm
);
1738 patchlevel_differs
= mm
;
1740 if (patchlevel_differs
)
1742 /* Gromacs should be able to continue from checkpoints between
1743 * different patch level versions, but we do not guarantee
1744 * compatibility between different major/minor versions - check this.
1746 int gmx_major
, gmx_minor
;
1747 int cpt_major
, cpt_minor
;
1748 sscanf(gmx_version(), "VERSION %5d.%5d", &gmx_major
, &gmx_minor
);
1749 int ret
= sscanf(version
, "VERSION %5d.%5d", &cpt_major
, &cpt_minor
);
1750 version_differs
= (ret
< 2 || gmx_major
!= cpt_major
||
1751 gmx_minor
!= cpt_minor
);
1754 check_string(fplog
, "Build time", BUILD_TIME
, btime
, &mm
);
1755 check_string(fplog
, "Build user", BUILD_USER
, buser
, &mm
);
1756 check_string(fplog
, "Build host", BUILD_HOST
, bhost
, &mm
);
1757 check_int (fplog
, "Double prec.", GMX_CPT_BUILD_DP
, double_prec
, &mm
);
1758 check_string(fplog
, "Program name", Program(), fprog
, &mm
);
1760 check_int (fplog
, "#ranks", cr
->nnodes
, npp_f
+npme_f
, &mm
);
1763 check_int (fplog
, "#PME-ranks", cr
->npmenodes
, npme_f
, &mm
);
1766 if (cr
->npmenodes
>= 0)
1768 npp
-= cr
->npmenodes
;
1772 check_int (fplog
, "#DD-cells[x]", dd_nc
[XX
], dd_nc_f
[XX
], &mm
);
1773 check_int (fplog
, "#DD-cells[y]", dd_nc
[YY
], dd_nc_f
[YY
], &mm
);
1774 check_int (fplog
, "#DD-cells[z]", dd_nc
[ZZ
], dd_nc_f
[ZZ
], &mm
);
1780 const char msg_version_difference
[] =
1781 "The current GROMACS major & minor version are not identical to those that\n"
1782 "generated the checkpoint file. In principle GROMACS does not support\n"
1783 "continuation from checkpoints between different versions, so we advise\n"
1784 "against this. If you still want to try your luck we recommend that you use\n"
1785 "the -noappend flag to keep your output files from the two versions separate.\n"
1786 "This might also work around errors where the output fields in the energy\n"
1787 "file have changed between the different major & minor versions.\n";
1789 const char msg_mismatch_notice
[] =
1790 "GROMACS patchlevel, binary or parallel settings differ from previous run.\n"
1791 "Continuation is exact, but not guaranteed to be binary identical.\n";
1793 const char msg_logdetails
[] =
1794 "See the log file for details.\n";
1796 if (version_differs
)
1798 fprintf(stderr
, "%s%s\n", msg_version_difference
, fplog
? msg_logdetails
: "");
1802 fprintf(fplog
, "%s\n", msg_version_difference
);
1807 /* Major & minor versions match at least, but something is different. */
1808 fprintf(stderr
, "%s%s\n", msg_mismatch_notice
, fplog
? msg_logdetails
: "");
1811 fprintf(fplog
, "%s\n", msg_mismatch_notice
);
1817 static void read_checkpoint(const char *fn
, FILE **pfplog
,
1818 t_commrec
*cr
, ivec dd_nc
,
1819 int eIntegrator
, int *init_fep_state
, gmx_int64_t
*step
, double *t
,
1820 t_state
*state
, gmx_bool
*bReadEkin
,
1821 int *simulation_part
,
1822 gmx_bool bAppendOutputFiles
, gmx_bool bForceAppend
)
1827 char *version
, *btime
, *buser
, *bhost
, *fprog
, *ftime
;
1829 char buf
[STEPSTRSIZE
];
1830 int eIntegrator_f
, nppnodes_f
, npmenodes_f
;
1832 int natoms
, ngtc
, nnhpres
, nhchainlength
, nlambda
, fflags
, flags_eks
, flags_enh
, flags_dfh
;
1835 gmx_file_position_t
*outputfiles
;
1837 t_fileio
*chksum_file
;
1838 FILE * fplog
= *pfplog
;
1839 unsigned char digest
[16];
1840 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1841 struct flock fl
; /* don't initialize here: the struct order is OS
1845 const char *int_warn
=
1846 "WARNING: The checkpoint file was generated with integrator %s,\n"
1847 " while the simulation uses integrator %s\n\n";
1849 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1850 fl
.l_type
= F_WRLCK
;
1851 fl
.l_whence
= SEEK_SET
;
1857 fp
= gmx_fio_open(fn
, "r");
1858 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, &file_version
,
1859 &version
, &btime
, &buser
, &bhost
, &double_prec
, &fprog
, &ftime
,
1860 &eIntegrator_f
, simulation_part
, step
, t
,
1861 &nppnodes_f
, dd_nc_f
, &npmenodes_f
,
1862 &natoms
, &ngtc
, &nnhpres
, &nhchainlength
, &nlambda
,
1863 &fflags
, &flags_eks
, &flags_enh
, &flags_dfh
,
1864 &state
->edsamstate
.nED
, &state
->swapstate
.eSwapCoords
, NULL
);
1866 if (bAppendOutputFiles
&&
1867 file_version
>= 13 && double_prec
!= GMX_CPT_BUILD_DP
)
1869 gmx_fatal(FARGS
, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
1872 if (cr
== NULL
|| MASTER(cr
))
1874 fprintf(stderr
, "\nReading checkpoint file %s generated: %s\n\n",
1878 /* This will not be written if we do appending, since fplog is still NULL then */
1881 fprintf(fplog
, "\n");
1882 fprintf(fplog
, "Reading checkpoint file %s\n", fn
);
1883 fprintf(fplog
, " file generated by: %s\n", fprog
);
1884 fprintf(fplog
, " file generated at: %s\n", ftime
);
1885 fprintf(fplog
, " GROMACS build time: %s\n", btime
);
1886 fprintf(fplog
, " GROMACS build user: %s\n", buser
);
1887 fprintf(fplog
, " GROMACS build host: %s\n", bhost
);
1888 fprintf(fplog
, " GROMACS double prec.: %d\n", double_prec
);
1889 fprintf(fplog
, " simulation part #: %d\n", *simulation_part
);
1890 fprintf(fplog
, " step: %s\n", gmx_step_str(*step
, buf
));
1891 fprintf(fplog
, " time: %f\n", *t
);
1892 fprintf(fplog
, "\n");
1895 if (natoms
!= state
->natoms
)
1897 gmx_fatal(FARGS
, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", natoms
, state
->natoms
);
1899 if (ngtc
!= state
->ngtc
)
1901 gmx_fatal(FARGS
, "Checkpoint file is for a system of %d T-coupling groups, while the current system consists of %d T-coupling groups", ngtc
, state
->ngtc
);
1903 if (nnhpres
!= state
->nnhpres
)
1905 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", nnhpres
, state
->nnhpres
);
1908 if (nlambda
!= state
->dfhist
.nlambda
)
1910 gmx_fatal(FARGS
, "Checkpoint file is for a system with %d lambda states, while the current system consists of %d lambda states", nlambda
, state
->dfhist
.nlambda
);
1913 init_gtc_state(state
, state
->ngtc
, state
->nnhpres
, nhchainlength
); /* need to keep this here to keep the tpr format working */
1914 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
1916 if (eIntegrator_f
!= eIntegrator
)
1920 fprintf(stderr
, int_warn
, EI(eIntegrator_f
), EI(eIntegrator
));
1922 if (bAppendOutputFiles
)
1925 "Output file appending requested, but input/checkpoint integrators do not match.\n"
1926 "Stopping the run to prevent you from ruining all your data...\n"
1927 "If you _really_ know what you are doing, try with the -noappend option.\n");
1931 fprintf(fplog
, int_warn
, EI(eIntegrator_f
), EI(eIntegrator
));
1939 else if (cr
->nnodes
== nppnodes_f
+ npmenodes_f
)
1941 if (cr
->npmenodes
< 0)
1943 cr
->npmenodes
= npmenodes_f
;
1945 int nppnodes
= cr
->nnodes
- cr
->npmenodes
;
1946 if (nppnodes
== nppnodes_f
)
1948 for (d
= 0; d
< DIM
; d
++)
1952 dd_nc
[d
] = dd_nc_f
[d
];
1958 if (fflags
!= state
->flags
)
1963 if (bAppendOutputFiles
)
1966 "Output file appending requested, but input and checkpoint states are not identical.\n"
1967 "Stopping the run to prevent you from ruining all your data...\n"
1968 "You can try with the -noappend option, and get more info in the log file.\n");
1971 if (getenv("GMX_ALLOW_CPT_MISMATCH") == NULL
)
1973 gmx_fatal(FARGS
, "You seem to have switched ensemble, integrator, T and/or P-coupling algorithm between the cpt and tpr file. The recommended way of doing this is passing the cpt file to grompp (with option -t) instead of to mdrun. If you know what you are doing, you can override this error by setting the env.var. GMX_ALLOW_CPT_MISMATCH");
1978 "WARNING: The checkpoint state entries do not match the simulation,\n"
1979 " see the log file for details\n\n");
1985 print_flag_mismatch(fplog
, state
->flags
, fflags
);
1992 check_match(fplog
, version
, btime
, buser
, bhost
, double_prec
, fprog
,
1993 cr
, nppnodes_f
, npmenodes_f
, dd_nc
, dd_nc_f
);
1996 ret
= do_cpt_state(gmx_fio_getxdr(fp
), TRUE
, fflags
, state
, NULL
);
1997 *init_fep_state
= state
->fep_state
; /* there should be a better way to do this than setting it here.
1998 Investigate for 5.0. */
2003 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), flags_eks
, &state
->ekinstate
, NULL
);
2008 *bReadEkin
= ((flags_eks
& (1<<eeksEKINH
)) || (flags_eks
& (1<<eeksEKINF
)) || (flags_eks
& (1<<eeksEKINO
)) ||
2009 ((flags_eks
& (1<<eeksEKINSCALEF
)) | (flags_eks
& (1<<eeksEKINSCALEH
)) | (flags_eks
& (1<<eeksVSCALE
))));
2011 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
,
2012 flags_enh
, state
->enerhist
, NULL
);
2018 if (file_version
< 6)
2020 const char *warn
= "Reading checkpoint file in old format, assuming that the run that generated this file started at step 0, if this is not the case the averages stored in the energy file will be incorrect.";
2022 fprintf(stderr
, "\nWARNING: %s\n\n", warn
);
2025 fprintf(fplog
, "\nWARNING: %s\n\n", warn
);
2027 state
->enerhist
->nsum
= *step
;
2028 state
->enerhist
->nsum_sim
= *step
;
2031 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
), flags_dfh
, &state
->dfhist
, NULL
);
2037 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, &state
->edsamstate
, NULL
);
2043 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, &state
->swapstate
, NULL
);
2049 ret
= do_cpt_files(gmx_fio_getxdr(fp
), TRUE
, &outputfiles
, &nfiles
, NULL
, file_version
);
2055 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), file_version
);
2060 if (gmx_fio_close(fp
) != 0)
2062 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2071 /* If the user wants to append to output files,
2072 * we use the file pointer positions of the output files stored
2073 * in the checkpoint file and truncate the files such that any frames
2074 * written after the checkpoint time are removed.
2075 * All files are md5sum checked such that we can be sure that
2076 * we do not truncate other (maybe imprortant) files.
2078 if (bAppendOutputFiles
)
2080 if (fn2ftp(outputfiles
[0].filename
) != efLOG
)
2082 /* make sure first file is log file so that it is OK to use it for
2085 gmx_fatal(FARGS
, "The first output file should always be the log "
2086 "file but instead is: %s. Cannot do appending because of this condition.", outputfiles
[0].filename
);
2088 for (i
= 0; i
< nfiles
; i
++)
2090 if (outputfiles
[i
].offset
< 0)
2092 gmx_fatal(FARGS
, "The original run wrote a file called '%s' which "
2093 "is larger than 2 GB, but mdrun did not support large file"
2094 " offsets. Can not append. Run mdrun with -noappend",
2095 outputfiles
[i
].filename
);
2098 chksum_file
= gmx_fio_open(outputfiles
[i
].filename
, "a");
2101 chksum_file
= gmx_fio_open(outputfiles
[i
].filename
, "r+");
2106 /* Note that there are systems where the lock operation
2107 * will succeed, but a second process can also lock the file.
2108 * We should probably try to detect this.
2110 #if defined __native_client__
2114 #elif defined GMX_NATIVE_WINDOWS
2115 if (_locking(fileno(gmx_fio_getfp(chksum_file
)), _LK_NBLCK
, LONG_MAX
) == -1)
2117 if (fcntl(fileno(gmx_fio_getfp(chksum_file
)), F_SETLK
, &fl
) == -1)
2120 if (errno
== ENOSYS
)
2124 gmx_fatal(FARGS
, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2128 fprintf(stderr
, "\nNOTE: File locking is not supported on this system, will not lock %s\n\n", outputfiles
[i
].filename
);
2131 fprintf(fplog
, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", outputfiles
[i
].filename
);
2135 else if (errno
== EACCES
|| errno
== EAGAIN
)
2137 gmx_fatal(FARGS
, "Failed to lock: %s. Already running "
2138 "simulation?", outputfiles
[i
].filename
);
2142 gmx_fatal(FARGS
, "Failed to lock: %s. %s.",
2143 outputfiles
[i
].filename
, std::strerror(errno
));
2148 /* compute md5 chksum */
2149 if (outputfiles
[i
].chksum_size
!= -1)
2151 if (gmx_fio_get_file_md5(chksum_file
, outputfiles
[i
].offset
,
2152 digest
) != outputfiles
[i
].chksum_size
) /*at the end of the call the file position is at the end of the file*/
2154 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.",
2155 outputfiles
[i
].chksum_size
,
2156 outputfiles
[i
].filename
);
2159 if (i
== 0) /*log file needs to be seeked in case we need to truncate (other files are truncated below)*/
2161 if (gmx_fio_seek(chksum_file
, outputfiles
[i
].offset
))
2163 gmx_fatal(FARGS
, "Seek error! Failed to truncate log-file: %s.", std::strerror(errno
));
2168 if (i
== 0) /*open log file here - so that lock is never lifted
2169 after chksum is calculated */
2171 *pfplog
= gmx_fio_getfp(chksum_file
);
2175 gmx_fio_close(chksum_file
);
2178 /* compare md5 chksum */
2179 if (outputfiles
[i
].chksum_size
!= -1 &&
2180 memcmp(digest
, outputfiles
[i
].chksum
, 16) != 0)
2184 fprintf(debug
, "chksum for %s: ", outputfiles
[i
].filename
);
2185 for (j
= 0; j
< 16; j
++)
2187 fprintf(debug
, "%02x", digest
[j
]);
2189 fprintf(debug
, "\n");
2191 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.",
2192 outputfiles
[i
].filename
);
2197 if (i
!= 0) /*log file is already seeked to correct position */
2199 #if !defined(GMX_NATIVE_WINDOWS) || !defined(GMX_FAHCORE)
2200 /* For FAHCORE, we do this elsewhere*/
2201 rc
= gmx_truncate(outputfiles
[i
].filename
, outputfiles
[i
].offset
);
2204 gmx_fatal(FARGS
, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfiles
[i
].filename
);
2215 void load_checkpoint(const char *fn
, FILE **fplog
,
2216 t_commrec
*cr
, ivec dd_nc
,
2217 t_inputrec
*ir
, t_state
*state
,
2218 gmx_bool
*bReadEkin
,
2219 gmx_bool bAppend
, gmx_bool bForceAppend
)
2226 /* Read the state from the checkpoint file */
2227 read_checkpoint(fn
, fplog
,
2229 ir
->eI
, &(ir
->fepvals
->init_fep_state
), &step
, &t
, state
, bReadEkin
,
2230 &ir
->simulation_part
, bAppend
, bForceAppend
);
2234 gmx_bcast(sizeof(cr
->npmenodes
), &cr
->npmenodes
, cr
);
2235 gmx_bcast(DIM
*sizeof(dd_nc
[0]), dd_nc
, cr
);
2236 gmx_bcast(sizeof(step
), &step
, cr
);
2237 gmx_bcast(sizeof(*bReadEkin
), bReadEkin
, cr
);
2239 ir
->bContinuation
= TRUE
;
2240 if (ir
->nsteps
>= 0)
2242 ir
->nsteps
+= ir
->init_step
- step
;
2244 ir
->init_step
= step
;
2245 ir
->simulation_part
+= 1;
2248 void read_checkpoint_part_and_step(const char *filename
,
2249 int *simulation_part
,
2253 char *version
, *btime
, *buser
, *bhost
, *fprog
, *ftime
;
2258 int flags_eks
, flags_enh
, flags_dfh
;
2263 if (filename
== NULL
||
2264 !gmx_fexist(filename
) ||
2265 (!(fp
= gmx_fio_open(filename
, "r"))))
2267 *simulation_part
= 0;
2272 /* Not calling initializing state before use is nasty, but all we
2273 do is read into its member variables and throw the struct away
2274 again immediately. */
2276 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, &file_version
,
2277 &version
, &btime
, &buser
, &bhost
, &double_prec
, &fprog
, &ftime
,
2278 &eIntegrator
, simulation_part
, step
, &t
, &nppnodes
, dd_nc
, &npme
,
2279 &state
.natoms
, &state
.ngtc
, &state
.nnhpres
, &state
.nhchainlength
,
2280 &(state
.dfhist
.nlambda
), &state
.flags
, &flags_eks
, &flags_enh
, &flags_dfh
,
2281 &state
.edsamstate
.nED
, &state
.swapstate
.eSwapCoords
, NULL
);
2286 static void read_checkpoint_data(t_fileio
*fp
, int *simulation_part
,
2287 gmx_int64_t
*step
, double *t
, t_state
*state
,
2288 int *nfiles
, gmx_file_position_t
**outputfiles
)
2291 char *version
, *btime
, *buser
, *bhost
, *fprog
, *ftime
;
2296 int flags_eks
, flags_enh
, flags_dfh
;
2298 gmx_file_position_t
*files_loc
= NULL
;
2301 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, &file_version
,
2302 &version
, &btime
, &buser
, &bhost
, &double_prec
, &fprog
, &ftime
,
2303 &eIntegrator
, simulation_part
, step
, t
, &nppnodes
, dd_nc
, &npme
,
2304 &state
->natoms
, &state
->ngtc
, &state
->nnhpres
, &state
->nhchainlength
,
2305 &(state
->dfhist
.nlambda
), &state
->flags
, &flags_eks
, &flags_enh
, &flags_dfh
,
2306 &state
->edsamstate
.nED
, &state
->swapstate
.eSwapCoords
, NULL
);
2308 do_cpt_state(gmx_fio_getxdr(fp
), TRUE
, state
->flags
, state
, NULL
);
2313 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), flags_eks
, &state
->ekinstate
, NULL
);
2318 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
,
2319 flags_enh
, state
->enerhist
, NULL
);
2324 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
), flags_dfh
, &state
->dfhist
, NULL
);
2330 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, &state
->edsamstate
, NULL
);
2336 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, &state
->swapstate
, NULL
);
2342 ret
= do_cpt_files(gmx_fio_getxdr(fp
), TRUE
,
2343 outputfiles
!= NULL
? outputfiles
: &files_loc
,
2344 outputfiles
!= NULL
? nfiles
: &nfiles_loc
,
2345 NULL
, file_version
);
2346 if (files_loc
!= NULL
)
2356 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), file_version
);
2370 read_checkpoint_state(const char *fn
, int *simulation_part
,
2371 gmx_int64_t
*step
, double *t
, t_state
*state
)
2375 fp
= gmx_fio_open(fn
, "r");
2376 read_checkpoint_data(fp
, simulation_part
, step
, t
, state
, NULL
, NULL
);
2377 if (gmx_fio_close(fp
) != 0)
2379 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2383 void read_checkpoint_trxframe(t_fileio
*fp
, t_trxframe
*fr
)
2385 /* This next line is nasty because the sub-structures of t_state
2386 * cannot be assumed to be zeroed (or even initialized in ways the
2387 * rest of the code might assume). Using snew would be better, but
2388 * this will all go away for 5.0. */
2390 int simulation_part
;
2394 init_state(&state
, 0, 0, 0, 0, 0);
2396 read_checkpoint_data(fp
, &simulation_part
, &step
, &t
, &state
, NULL
, NULL
);
2398 fr
->natoms
= state
.natoms
;
2401 fr
->step
= gmx_int64_to_int(step
,
2402 "conversion of checkpoint to trajectory");
2406 fr
->lambda
= state
.lambda
[efptFEP
];
2407 fr
->fep_state
= state
.fep_state
;
2409 fr
->bX
= (state
.flags
& (1<<estX
));
2415 fr
->bV
= (state
.flags
& (1<<estV
));
2422 fr
->bBox
= (state
.flags
& (1<<estBOX
));
2425 copy_mat(state
.box
, fr
->box
);
2430 void list_checkpoint(const char *fn
, FILE *out
)
2434 char *version
, *btime
, *buser
, *bhost
, *fprog
, *ftime
;
2436 int eIntegrator
, simulation_part
, nppnodes
, npme
;
2441 int flags_eks
, flags_enh
, flags_dfh
;
2443 gmx_file_position_t
*outputfiles
;
2446 init_state(&state
, -1, -1, -1, -1, 0);
2448 fp
= gmx_fio_open(fn
, "r");
2449 do_cpt_header(gmx_fio_getxdr(fp
), TRUE
, &file_version
,
2450 &version
, &btime
, &buser
, &bhost
, &double_prec
, &fprog
, &ftime
,
2451 &eIntegrator
, &simulation_part
, &step
, &t
, &nppnodes
, dd_nc
, &npme
,
2452 &state
.natoms
, &state
.ngtc
, &state
.nnhpres
, &state
.nhchainlength
,
2453 &(state
.dfhist
.nlambda
), &state
.flags
,
2454 &flags_eks
, &flags_enh
, &flags_dfh
, &state
.edsamstate
.nED
,
2455 &state
.swapstate
.eSwapCoords
, out
);
2456 ret
= do_cpt_state(gmx_fio_getxdr(fp
), TRUE
, state
.flags
, &state
, out
);
2461 ret
= do_cpt_ekinstate(gmx_fio_getxdr(fp
), flags_eks
, &state
.ekinstate
, out
);
2466 ret
= do_cpt_enerhist(gmx_fio_getxdr(fp
), TRUE
,
2467 flags_enh
, state
.enerhist
, out
);
2471 ret
= do_cpt_df_hist(gmx_fio_getxdr(fp
),
2472 flags_dfh
, &state
.dfhist
, out
);
2477 ret
= do_cpt_EDstate(gmx_fio_getxdr(fp
), TRUE
, &state
.edsamstate
, out
);
2482 ret
= do_cpt_swapstate(gmx_fio_getxdr(fp
), TRUE
, &state
.swapstate
, out
);
2487 do_cpt_files(gmx_fio_getxdr(fp
), TRUE
, &outputfiles
, &nfiles
, out
, file_version
);
2492 ret
= do_cpt_footer(gmx_fio_getxdr(fp
), file_version
);
2499 if (gmx_fio_close(fp
) != 0)
2501 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2507 /* This routine cannot print tons of data, since it is called before the log file is opened. */
2509 read_checkpoint_simulation_part_and_filenames(t_fileio
*fp
,
2510 int *simulation_part
,
2512 gmx_file_position_t
**outputfiles
)
2514 gmx_int64_t step
= 0;
2518 init_state(&state
, 0, 0, 0, 0, 0);
2520 read_checkpoint_data(fp
, simulation_part
, &step
, &t
, &state
,
2521 nfiles
, outputfiles
);
2522 if (gmx_fio_close(fp
) != 0)
2524 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");