2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2008,2009,2010,2011,2012,2013,2014, 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. */
38 #include "checkpoint.h"
48 #ifdef HAVE_SYS_TIME_H
56 #ifdef GMX_NATIVE_WINDOWS
59 #include <sys/locking.h>
65 #include "types/commrec.h"
67 #include "gromacs/math/vec.h"
70 #include "gromacs/fileio/filenm.h"
71 #include "gromacs/utility/futil.h"
72 #include "gromacs/fileio/gmxfio.h"
73 #include "gromacs/fileio/xdrf.h"
74 #include "gromacs/fileio/xdr_datatype.h"
75 #include "gromacs/utility/basenetwork.h"
76 #include "gromacs/utility/baseversion.h"
77 #include "gromacs/utility/cstringutil.h"
78 #include "gromacs/utility/fatalerror.h"
79 #include "gromacs/utility/smalloc.h"
81 #include "buildinfo.h"
87 #define CPT_MAGIC1 171817
88 #define CPT_MAGIC2 171819
89 #define CPTSTRLEN 1024
92 #define GMX_CPT_BUILD_DP 1
94 #define GMX_CPT_BUILD_DP 0
97 /* cpt_version should normally only be changed
98 * when the header of footer format changes.
99 * The state data format itself is backward and forward compatible.
100 * But old code can not read a new entry that is present in the file
101 * (but can read a new format when new entries are not present).
103 static const int cpt_version = 16;
106 const char *est_names[estNR] =
109 "box", "box-rel", "box-v", "pres_prev",
110 "nosehoover-xi", "thermostat-integral",
111 "x", "v", "SDx", "CGp", "LD-rng", "LD-rng-i",
112 "disre_initf", "disre_rm3tav",
113 "orire_initf", "orire_Dtav",
114 "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng", "MC-rng-i"
118 eeksEKIN_N, eeksEKINH, eeksDEKINDL, eeksMVCOS, eeksEKINF, eeksEKINO, eeksEKINSCALEF, eeksEKINSCALEH, eeksVSCALE, eeksEKINTOTAL, eeksNR
121 const char *eeks_names[eeksNR] =
123 "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
124 "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
128 eenhENERGY_N, eenhENERGY_AVER, eenhENERGY_SUM, eenhENERGY_NSUM,
129 eenhENERGY_SUM_SIM, eenhENERGY_NSUM_SIM,
130 eenhENERGY_NSTEPS, eenhENERGY_NSTEPS_SIM,
131 eenhENERGY_DELTA_H_NN,
132 eenhENERGY_DELTA_H_LIST,
133 eenhENERGY_DELTA_H_STARTTIME,
134 eenhENERGY_DELTA_H_STARTLAMBDA,
138 const char *eenh_names[eenhNR] =
140 "energy_n", "energy_aver", "energy_sum", "energy_nsum",
141 "energy_sum_sim", "energy_nsum_sim",
142 "energy_nsteps", "energy_nsteps_sim",
144 "energy_delta_h_list",
145 "energy_delta_h_start_time",
146 "energy_delta_h_start_lambda"
149 /* free energy history variables -- need to be preserved over checkpoint */
151 edfhBEQUIL, edfhNATLAMBDA, edfhWLHISTO, edfhWLDELTA, edfhSUMWEIGHTS, edfhSUMDG, edfhSUMMINVAR, edfhSUMVAR,
152 edfhACCUMP, edfhACCUMM, edfhACCUMP2, edfhACCUMM2, edfhTIJ, edfhTIJEMP, edfhNR
154 /* free energy history variable names */
155 const char *edfh_names[edfhNR] =
157 "bEquilibrated", "N_at_state", "Wang-Landau Histogram", "Wang-Landau Delta", "Weights", "Free Energies", "minvar", "variance",
158 "accumulated_plus", "accumulated_minus", "accumulated_plus_2", "accumulated_minus_2", "Tij", "Tij_empirical"
161 #ifdef GMX_NATIVE_WINDOWS
163 gmx_wintruncate(const char *filename, __int64 size)
166 /*we do this elsewhere*/
171 fp = fopen(filename, "rb+");
178 return _chsize_s( fileno(fp), size);
185 ecprREAL, ecprRVEC, ecprMATRIX
189 cptpEST, cptpEEKS, cptpEENH, cptpEDFH
191 /* enums for the different components of checkpoint variables, replacing the hard coded ones.
192 cptpEST - state variables.
193 cptpEEKS - Kinetic energy state variables.
194 cptpEENH - Energy history state variables.
195 cptpEDFH - free energy history variables.
199 static const char *st_names(int cptp, int ecpt)
203 case cptpEST: return est_names [ecpt]; break;
204 case cptpEEKS: return eeks_names[ecpt]; break;
205 case cptpEENH: return eenh_names[ecpt]; break;
206 case cptpEDFH: return edfh_names[ecpt]; break;
212 static void cp_warning(FILE *fp)
214 fprintf(fp, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
217 static void cp_error()
219 gmx_fatal(FARGS, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
222 static void do_cpt_string_err(XDR *xd, gmx_bool bRead, const char *desc, char **s, FILE *list)
230 res = xdr_string(xd, s, CPTSTRLEN);
237 fprintf(list, "%s = %s\n", desc, *s);
242 static int do_cpt_int(XDR *xd, const char *desc, int *i, FILE *list)
246 res = xdr_int(xd, i);
253 fprintf(list, "%s = %d\n", desc, *i);
258 static int do_cpt_u_chars(XDR *xd, const char *desc, int n, unsigned char *i, FILE *list)
264 fprintf(list, "%s = ", desc);
266 for (j = 0; j < n && res; j++)
268 res &= xdr_u_char(xd, &i[j]);
271 fprintf(list, "%02x", i[j]);
286 static void do_cpt_int_err(XDR *xd, const char *desc, int *i, FILE *list)
288 if (do_cpt_int(xd, desc, i, list) < 0)
294 static void do_cpt_step_err(XDR *xd, const char *desc, gmx_int64_t *i, FILE *list)
297 char buf[STEPSTRSIZE];
299 res = xdr_int64(xd, i);
306 fprintf(list, "%s = %s\n", desc, gmx_step_str(*i, buf));
310 static void do_cpt_double_err(XDR *xd, const char *desc, double *f, FILE *list)
314 res = xdr_double(xd, f);
321 fprintf(list, "%s = %f\n", desc, *f);
325 static void do_cpt_real_err(XDR *xd, real *f)
330 res = xdr_double(xd, f);
332 res = xdr_float(xd, f);
340 static void do_cpt_n_rvecs_err(XDR *xd, const char *desc, int n, rvec f[], FILE *list)
344 for (i = 0; i < n; i++)
346 for (j = 0; j < DIM; j++)
348 do_cpt_real_err(xd, &f[i][j]);
354 pr_rvecs(list, 0, desc, f, n);
358 /* If nval >= 0, nval is used; on read this should match the passed value.
359 * If nval n<0, *nptr is used; on read the value is stored in nptr
361 static int do_cpte_reals_low(XDR *xd, int cptp, int ecpt, int sflags,
362 int nval, int *nptr, real **v,
363 FILE *list, int erealtype)
367 int dtc = xdr_datatype_float;
369 int dtc = xdr_datatype_double;
371 real *vp, *va = NULL;
386 gmx_incons("*ntpr=NULL in do_cpte_reals_low");
391 res = xdr_int(xd, &nf);
402 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), nval, nf);
411 res = xdr_int(xd, &dt);
418 fprintf(stderr, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
419 st_names(cptp, ecpt), xdr_datatype_names[dtc],
420 xdr_datatype_names[dt]);
422 if (list || !(sflags & (1<<ecpt)))
435 if (dt == xdr_datatype_float)
437 if (dtc == xdr_datatype_float)
445 res = xdr_vector(xd, (char *)vf, nf,
446 (unsigned int)sizeof(float), (xdrproc_t)xdr_float);
451 if (dtc != xdr_datatype_float)
453 for (i = 0; i < nf; i++)
462 if (dtc == xdr_datatype_double)
464 /* cppcheck-suppress invalidPointerCast
465 * Only executed if real is anyhow double */
472 res = xdr_vector(xd, (char *)vd, nf,
473 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
478 if (dtc != xdr_datatype_double)
480 for (i = 0; i < nf; i++)
493 pr_reals(list, 0, st_names(cptp, ecpt), vp, nf);
496 pr_rvecs(list, 0, st_names(cptp, ecpt), (rvec *)vp, nf/3);
499 gmx_incons("Unknown checkpoint real type");
511 /* This function stores n along with the reals for reading,
512 * but on reading it assumes that n matches the value in the checkpoint file,
513 * a fatal error is generated when this is not the case.
515 static int do_cpte_reals(XDR *xd, int cptp, int ecpt, int sflags,
516 int n, real **v, FILE *list)
518 return do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, v, list, ecprREAL);
521 /* This function does the same as do_cpte_reals,
522 * except that on reading it ignores the passed value of *n
523 * and stored the value read from the checkpoint file in *n.
525 static int do_cpte_n_reals(XDR *xd, int cptp, int ecpt, int sflags,
526 int *n, real **v, FILE *list)
528 return do_cpte_reals_low(xd, cptp, ecpt, sflags, -1, n, v, list, ecprREAL);
531 static int do_cpte_real(XDR *xd, int cptp, int ecpt, int sflags,
534 return do_cpte_reals_low(xd, cptp, ecpt, sflags, 1, NULL, &r, list, ecprREAL);
537 static int do_cpte_ints(XDR *xd, int cptp, int ecpt, int sflags,
538 int n, int **v, FILE *list)
541 int dtc = xdr_datatype_int;
546 res = xdr_int(xd, &nf);
551 if (list == NULL && v != NULL && nf != n)
553 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
556 res = xdr_int(xd, &dt);
563 gmx_fatal(FARGS, "Type mismatch for state entry %s, code type is %s, file type is %s\n",
564 st_names(cptp, ecpt), xdr_datatype_names[dtc],
565 xdr_datatype_names[dt]);
567 if (list || !(sflags & (1<<ecpt)) || v == NULL)
580 res = xdr_vector(xd, (char *)vp, nf,
581 (unsigned int)sizeof(int), (xdrproc_t)xdr_int);
588 pr_ivec(list, 0, st_names(cptp, ecpt), vp, nf, TRUE);
598 static int do_cpte_int(XDR *xd, int cptp, int ecpt, int sflags,
601 return do_cpte_ints(xd, cptp, ecpt, sflags, 1, &i, list);
604 static int do_cpte_doubles(XDR *xd, int cptp, int ecpt, int sflags,
605 int n, double **v, FILE *list)
608 int dtc = xdr_datatype_double;
609 double *vp, *va = NULL;
613 res = xdr_int(xd, &nf);
618 if (list == NULL && nf != n)
620 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
623 res = xdr_int(xd, &dt);
630 gmx_fatal(FARGS, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
631 st_names(cptp, ecpt), xdr_datatype_names[dtc],
632 xdr_datatype_names[dt]);
634 if (list || !(sflags & (1<<ecpt)))
647 res = xdr_vector(xd, (char *)vp, nf,
648 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
655 pr_doubles(list, 0, st_names(cptp, ecpt), vp, nf);
665 static int do_cpte_double(XDR *xd, int cptp, int ecpt, int sflags,
666 double *r, FILE *list)
668 return do_cpte_doubles(xd, cptp, ecpt, sflags, 1, &r, list);
672 static int do_cpte_rvecs(XDR *xd, int cptp, int ecpt, int sflags,
673 int n, rvec **v, FILE *list)
675 return do_cpte_reals_low(xd, cptp, ecpt, sflags,
676 n*DIM, NULL, (real **)v, list, ecprRVEC);
679 static int do_cpte_matrix(XDR *xd, int cptp, int ecpt, int sflags,
680 matrix v, FILE *list)
685 vr = (real *)&(v[0][0]);
686 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
687 DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
689 if (list && ret == 0)
691 pr_rvecs(list, 0, st_names(cptp, ecpt), v, DIM);
698 static int do_cpte_nmatrix(XDR *xd, int cptp, int ecpt, int sflags,
699 int n, real **v, FILE *list)
703 char name[CPTSTRLEN];
710 for (i = 0; i < n; i++)
712 reti = do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, &(v[i]), NULL, ecprREAL);
713 if (list && reti == 0)
715 sprintf(name, "%s[%d]", st_names(cptp, ecpt), i);
716 pr_reals(list, 0, name, v[i], n);
726 static int do_cpte_matrices(XDR *xd, int cptp, int ecpt, int sflags,
727 int n, matrix **v, FILE *list)
730 matrix *vp, *va = NULL;
736 res = xdr_int(xd, &nf);
741 if (list == NULL && nf != n)
743 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
745 if (list || !(sflags & (1<<ecpt)))
758 snew(vr, nf*DIM*DIM);
759 for (i = 0; i < nf; i++)
761 for (j = 0; j < DIM; j++)
763 for (k = 0; k < DIM; k++)
765 vr[(i*DIM+j)*DIM+k] = vp[i][j][k];
769 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
770 nf*DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
771 for (i = 0; i < nf; i++)
773 for (j = 0; j < DIM; j++)
775 for (k = 0; k < DIM; k++)
777 vp[i][j][k] = vr[(i*DIM+j)*DIM+k];
783 if (list && ret == 0)
785 for (i = 0; i < nf; i++)
787 pr_rvecs(list, 0, st_names(cptp, ecpt), vp[i], DIM);
798 static void do_cpt_header(XDR *xd, gmx_bool bRead, int *file_version,
799 char **version, char **btime, char **buser, char **bhost,
801 char **fprog, char **ftime,
802 int *eIntegrator, int *simulation_part,
803 gmx_int64_t *step, double *t,
804 int *nnodes, int *dd_nc, int *npme,
805 int *natoms, int *ngtc, int *nnhpres, int *nhchainlength,
806 int *nlambda, int *flags_state,
807 int *flags_eks, int *flags_enh, int *flags_dfh,
808 int *nED, int *eSwapCoords,
824 res = xdr_int(xd, &magic);
827 gmx_fatal(FARGS, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
829 if (magic != CPT_MAGIC1)
831 gmx_fatal(FARGS, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
832 "The checkpoint file is corrupted or not a checkpoint file",
838 gmx_gethostname(fhost, 255);
840 do_cpt_string_err(xd, bRead, "GROMACS version", version, list);
841 do_cpt_string_err(xd, bRead, "GROMACS build time", btime, list);
842 do_cpt_string_err(xd, bRead, "GROMACS build user", buser, list);
843 do_cpt_string_err(xd, bRead, "GROMACS build host", bhost, list);
844 do_cpt_string_err(xd, bRead, "generating program", fprog, list);
845 do_cpt_string_err(xd, bRead, "generation time", ftime, list);
846 *file_version = cpt_version;
847 do_cpt_int_err(xd, "checkpoint file version", file_version, list);
848 if (*file_version > cpt_version)
850 gmx_fatal(FARGS, "Attempting to read a checkpoint file of version %d with code of version %d\n", *file_version, cpt_version);
852 if (*file_version >= 13)
854 do_cpt_int_err(xd, "GROMACS double precision", double_prec, list);
860 if (*file_version >= 12)
862 do_cpt_string_err(xd, bRead, "generating host", &fhost, list);
868 do_cpt_int_err(xd, "#atoms", natoms, list);
869 do_cpt_int_err(xd, "#T-coupling groups", ngtc, list);
870 if (*file_version >= 10)
872 do_cpt_int_err(xd, "#Nose-Hoover T-chains", nhchainlength, list);
878 if (*file_version >= 11)
880 do_cpt_int_err(xd, "#Nose-Hoover T-chains for barostat ", nnhpres, list);
886 if (*file_version >= 14)
888 do_cpt_int_err(xd, "# of total lambda states ", nlambda, list);
894 do_cpt_int_err(xd, "integrator", eIntegrator, list);
895 if (*file_version >= 3)
897 do_cpt_int_err(xd, "simulation part #", simulation_part, list);
901 *simulation_part = 1;
903 if (*file_version >= 5)
905 do_cpt_step_err(xd, "step", step, list);
909 do_cpt_int_err(xd, "step", &idum, list);
912 do_cpt_double_err(xd, "t", t, list);
913 do_cpt_int_err(xd, "#PP-ranks", nnodes, list);
915 do_cpt_int_err(xd, "dd_nc[x]", dd_nc ? &(dd_nc[0]) : &idum, list);
916 do_cpt_int_err(xd, "dd_nc[y]", dd_nc ? &(dd_nc[1]) : &idum, list);
917 do_cpt_int_err(xd, "dd_nc[z]", dd_nc ? &(dd_nc[2]) : &idum, list);
918 do_cpt_int_err(xd, "#PME-only ranks", npme, list);
919 do_cpt_int_err(xd, "state flags", flags_state, list);
920 if (*file_version >= 4)
922 do_cpt_int_err(xd, "ekin data flags", flags_eks, list);
923 do_cpt_int_err(xd, "energy history flags", flags_enh, list);
928 *flags_enh = (*flags_state >> (estORIRE_DTAV+1));
929 *flags_state = (*flags_state & ~((1<<(estORIRE_DTAV+1)) |
930 (1<<(estORIRE_DTAV+2)) |
931 (1<<(estORIRE_DTAV+3))));
933 if (*file_version >= 14)
935 do_cpt_int_err(xd, "df history flags", flags_dfh, list);
942 if (*file_version >= 15)
944 do_cpt_int_err(xd, "ED data sets", nED, list);
950 if (*file_version >= 16)
952 do_cpt_int_err(xd, "swap", eSwapCoords, list);
956 static int do_cpt_footer(XDR *xd, int file_version)
961 if (file_version >= 2)
964 res = xdr_int(xd, &magic);
969 if (magic != CPT_MAGIC2)
978 static int do_cpt_state(XDR *xd, gmx_bool bRead,
979 int fflags, t_state *state,
989 nnht = state->nhchainlength*state->ngtc;
990 nnhtp = state->nhchainlength*state->nnhpres;
992 if (bRead) /* we need to allocate space for dfhist if we are reading */
994 init_df_history(&state->dfhist, state->dfhist.nlambda);
997 sflags = state->flags;
998 for (i = 0; (i < estNR && ret == 0); i++)
1000 if (fflags & (1<<i))
1004 case estLAMBDA: ret = do_cpte_reals(xd, cptpEST, i, sflags, efptNR, &(state->lambda), list); break;
1005 case estFEPSTATE: ret = do_cpte_int (xd, cptpEST, i, sflags, &state->fep_state, list); break;
1006 case estBOX: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box, list); break;
1007 case estBOX_REL: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box_rel, list); break;
1008 case estBOXV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->boxv, list); break;
1009 case estPRES_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->pres_prev, list); break;
1010 case estSVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->svir_prev, list); break;
1011 case estFVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->fvir_prev, list); break;
1012 case estNH_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_xi, list); break;
1013 case estNH_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_vxi, list); break;
1014 case estNHPRES_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_xi, list); break;
1015 case estNHPRES_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_vxi, list); break;
1016 case estTC_INT: ret = do_cpte_doubles(xd, cptpEST, i, sflags, state->ngtc, &state->therm_integral, list); break;
1017 case estVETA: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->veta, list); break;
1018 case estVOL0: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->vol0, list); break;
1019 case estX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->x, list); break;
1020 case estV: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->v, list); break;
1021 case estSDX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->sd_X, list); break;
1022 /* The RNG entries are no longer written,
1023 * the next 4 lines are only for reading old files.
1025 case estLD_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1026 case estLD_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1027 case estMC_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1028 case estMC_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1029 case estDISRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.disre_initf, list); break;
1030 case estDISRE_RM3TAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.ndisrepairs, &state->hist.disre_rm3tav, list); break;
1031 case estORIRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.orire_initf, list); break;
1032 case estORIRE_DTAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.norire_Dtav, &state->hist.orire_Dtav, list); break;
1034 gmx_fatal(FARGS, "Unknown state entry %d\n"
1035 "You are probably reading a new checkpoint file with old code", i);
1043 static int do_cpt_ekinstate(XDR *xd, int fflags, ekinstate_t *ekins,
1051 for (i = 0; (i < eeksNR && ret == 0); i++)
1053 if (fflags & (1<<i))
1058 case eeksEKIN_N: ret = do_cpte_int(xd, cptpEEKS, i, fflags, &ekins->ekin_n, list); break;
1059 case eeksEKINH: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh, list); break;
1060 case eeksEKINF: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinf, list); break;
1061 case eeksEKINO: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh_old, list); break;
1062 case eeksEKINTOTAL: ret = do_cpte_matrix(xd, cptpEEKS, i, fflags, ekins->ekin_total, list); break;
1063 case eeksEKINSCALEF: ret = do_cpte_doubles(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinscalef_nhc, list); break;
1064 case eeksVSCALE: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->vscale_nhc, list); break;
1065 case eeksEKINSCALEH: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->ekinscaleh_nhc, list); break;
1066 case eeksDEKINDL: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->dekindl, list); break;
1067 case eeksMVCOS: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->mvcos, list); break;
1069 gmx_fatal(FARGS, "Unknown ekin data state entry %d\n"
1070 "You are probably reading a new checkpoint file with old code", i);
1079 static int do_cpt_swapstate(XDR *xd, gmx_bool bRead, swapstate_t *swapstate, FILE *list)
1083 int swap_cpt_version = 1;
1086 if (eswapNO == swapstate->eSwapCoords)
1091 swapstate->bFromCpt = bRead;
1093 do_cpt_int_err(xd, "swap checkpoint version", &swap_cpt_version, list);
1094 do_cpt_int_err(xd, "swap coupling steps", &swapstate->nAverage, list);
1096 /* When reading, init_swapcoords has not been called yet,
1097 * so we have to allocate memory first. */
1099 for (ic = 0; ic < eCompNR; ic++)
1101 for (ii = 0; ii < eIonNR; ii++)
1105 do_cpt_int_err(xd, "swap requested atoms", &swapstate->nat_req[ic][ii], list);
1109 do_cpt_int_err(xd, "swap requested atoms p", swapstate->nat_req_p[ic][ii], list);
1114 do_cpt_int_err(xd, "swap influx netto", &swapstate->inflow_netto[ic][ii], list);
1118 do_cpt_int_err(xd, "swap influx netto p", swapstate->inflow_netto_p[ic][ii], list);
1121 if (bRead && (NULL == swapstate->nat_past[ic][ii]) )
1123 snew(swapstate->nat_past[ic][ii], swapstate->nAverage);
1126 for (j = 0; j < swapstate->nAverage; j++)
1130 do_cpt_int_err(xd, "swap past atom counts", &swapstate->nat_past[ic][ii][j], list);
1134 do_cpt_int_err(xd, "swap past atom counts p", &swapstate->nat_past_p[ic][ii][j], list);
1140 /* Ion flux per channel */
1141 for (ic = 0; ic < eChanNR; ic++)
1143 for (ii = 0; ii < eIonNR; ii++)
1147 do_cpt_int_err(xd, "channel flux", &swapstate->fluxfromAtoB[ic][ii], list);
1151 do_cpt_int_err(xd, "channel flux p", swapstate->fluxfromAtoB_p[ic][ii], list);
1156 /* Ion flux leakage */
1159 snew(swapstate->fluxleak, 1);
1161 do_cpt_int_err(xd, "flux leakage", swapstate->fluxleak, list);
1164 do_cpt_int_err(xd, "number of ions", &swapstate->nions, list);
1168 snew(swapstate->channel_label, swapstate->nions);
1169 snew(swapstate->comp_from, swapstate->nions);
1172 do_cpt_u_chars(xd, "channel history", swapstate->nions, swapstate->channel_label, list);
1173 do_cpt_u_chars(xd, "domain history", swapstate->nions, swapstate->comp_from, list);
1175 /* Save the last known whole positions to checkpoint
1176 * file to be able to also make multimeric channels whole in PBC */
1177 do_cpt_int_err(xd, "Ch0 atoms", &swapstate->nat[eChan0], list);
1178 do_cpt_int_err(xd, "Ch1 atoms", &swapstate->nat[eChan1], list);
1181 snew(swapstate->xc_old_whole[eChan0], swapstate->nat[eChan0]);
1182 snew(swapstate->xc_old_whole[eChan1], swapstate->nat[eChan1]);
1183 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], swapstate->xc_old_whole[eChan0], list);
1184 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], swapstate->xc_old_whole[eChan1], list);
1188 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], *swapstate->xc_old_whole_p[eChan0], list);
1189 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], *swapstate->xc_old_whole_p[eChan1], list);
1196 static int do_cpt_enerhist(XDR *xd, gmx_bool bRead,
1197 int fflags, energyhistory_t *enerhist,
1208 enerhist->nsteps = 0;
1210 enerhist->nsteps_sim = 0;
1211 enerhist->nsum_sim = 0;
1212 enerhist->dht = NULL;
1214 if (fflags & (1<< eenhENERGY_DELTA_H_NN) )
1216 snew(enerhist->dht, 1);
1217 enerhist->dht->ndh = NULL;
1218 enerhist->dht->dh = NULL;
1219 enerhist->dht->start_lambda_set = FALSE;
1223 for (i = 0; (i < eenhNR && ret == 0); i++)
1225 if (fflags & (1<<i))
1229 case eenhENERGY_N: ret = do_cpte_int(xd, cptpEENH, i, fflags, &enerhist->nener, list); break;
1230 case eenhENERGY_AVER: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_ave, list); break;
1231 case eenhENERGY_SUM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum, list); break;
1232 case eenhENERGY_NSUM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum, list); break;
1233 case eenhENERGY_SUM_SIM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum_sim, list); break;
1234 case eenhENERGY_NSUM_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum_sim, list); break;
1235 case eenhENERGY_NSTEPS: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps, list); break;
1236 case eenhENERGY_NSTEPS_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps_sim, list); break;
1237 case eenhENERGY_DELTA_H_NN: do_cpt_int_err(xd, eenh_names[i], &(enerhist->dht->nndh), list);
1238 if (bRead) /* now allocate memory for it */
1240 snew(enerhist->dht->dh, enerhist->dht->nndh);
1241 snew(enerhist->dht->ndh, enerhist->dht->nndh);
1242 for (j = 0; j < enerhist->dht->nndh; j++)
1244 enerhist->dht->ndh[j] = 0;
1245 enerhist->dht->dh[j] = NULL;
1249 case eenhENERGY_DELTA_H_LIST:
1250 for (j = 0; j < enerhist->dht->nndh; j++)
1252 ret = do_cpte_n_reals(xd, cptpEENH, i, fflags, &enerhist->dht->ndh[j], &(enerhist->dht->dh[j]), list);
1255 case eenhENERGY_DELTA_H_STARTTIME:
1256 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_time), list); break;
1257 case eenhENERGY_DELTA_H_STARTLAMBDA:
1258 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_lambda), list); break;
1260 gmx_fatal(FARGS, "Unknown energy history entry %d\n"
1261 "You are probably reading a new checkpoint file with old code", i);
1266 if ((fflags & (1<<eenhENERGY_SUM)) && !(fflags & (1<<eenhENERGY_SUM_SIM)))
1268 /* Assume we have an old file format and copy sum to sum_sim */
1269 srenew(enerhist->ener_sum_sim, enerhist->nener);
1270 for (i = 0; i < enerhist->nener; i++)
1272 enerhist->ener_sum_sim[i] = enerhist->ener_sum[i];
1276 if ( (fflags & (1<<eenhENERGY_NSUM)) &&
1277 !(fflags & (1<<eenhENERGY_NSTEPS)))
1279 /* Assume we have an old file format and copy nsum to nsteps */
1280 enerhist->nsteps = enerhist->nsum;
1282 if ( (fflags & (1<<eenhENERGY_NSUM_SIM)) &&
1283 !(fflags & (1<<eenhENERGY_NSTEPS_SIM)))
1285 /* Assume we have an old file format and copy nsum to nsteps */
1286 enerhist->nsteps_sim = enerhist->nsum_sim;
1292 static int do_cpt_df_hist(XDR *xd, int fflags, df_history_t *dfhist, FILE *list)
1297 nlambda = dfhist->nlambda;
1300 for (i = 0; (i < edfhNR && ret == 0); i++)
1302 if (fflags & (1<<i))
1306 case edfhBEQUIL: ret = do_cpte_int(xd, cptpEDFH, i, fflags, &dfhist->bEquil, list); break;
1307 case edfhNATLAMBDA: ret = do_cpte_ints(xd, cptpEDFH, i, fflags, nlambda, &dfhist->n_at_lam, list); break;
1308 case edfhWLHISTO: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->wl_histo, list); break;
1309 case edfhWLDELTA: ret = do_cpte_real(xd, cptpEDFH, i, fflags, &dfhist->wl_delta, list); break;
1310 case edfhSUMWEIGHTS: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_weights, list); break;
1311 case edfhSUMDG: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_dg, list); break;
1312 case edfhSUMMINVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_minvar, list); break;
1313 case edfhSUMVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_variance, list); break;
1314 case edfhACCUMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p, list); break;
1315 case edfhACCUMM: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m, list); break;
1316 case edfhACCUMP2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p2, list); break;
1317 case edfhACCUMM2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m2, list); break;
1318 case edfhTIJ: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij, list); break;
1319 case edfhTIJEMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij_empirical, list); break;
1322 gmx_fatal(FARGS, "Unknown df history entry %d\n"
1323 "You are probably reading a new checkpoint file with old code", i);
1332 /* This function stores the last whole configuration of the reference and
1333 * average structure in the .cpt file
1335 static int do_cpt_EDstate(XDR *xd, gmx_bool bRead,
1336 edsamstate_t *EDstate, FILE *list)
1343 EDstate->bFromCpt = bRead;
1345 if (EDstate->nED <= 0)
1350 /* When reading, init_edsam has not been called yet,
1351 * so we have to allocate memory first. */
1354 snew(EDstate->nref, EDstate->nED);
1355 snew(EDstate->old_sref, EDstate->nED);
1356 snew(EDstate->nav, EDstate->nED);
1357 snew(EDstate->old_sav, EDstate->nED);
1360 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1361 for (i = 0; i < EDstate->nED; i++)
1363 /* Reference structure SREF */
1364 sprintf(buf, "ED%d # of atoms in reference structure", i+1);
1365 do_cpt_int_err(xd, buf, &EDstate->nref[i], list);
1366 sprintf(buf, "ED%d x_ref", i+1);
1369 snew(EDstate->old_sref[i], EDstate->nref[i]);
1370 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref[i], list);
1374 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref_p[i], list);
1377 /* Average structure SAV */
1378 sprintf(buf, "ED%d # of atoms in average structure", i+1);
1379 do_cpt_int_err(xd, buf, &EDstate->nav[i], list);
1380 sprintf(buf, "ED%d x_av", i+1);
1383 snew(EDstate->old_sav[i], EDstate->nav[i]);
1384 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav[i], list);
1388 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav_p[i], list);
1396 static int do_cpt_files(XDR *xd, gmx_bool bRead,
1397 gmx_file_position_t **p_outputfiles, int *nfiles,
1398 FILE *list, int file_version)
1402 gmx_off_t mask = 0xFFFFFFFFL;
1403 int offset_high, offset_low;
1405 gmx_file_position_t *outputfiles;
1407 if (do_cpt_int(xd, "number of output files", nfiles, list) != 0)
1414 snew(*p_outputfiles, *nfiles);
1417 outputfiles = *p_outputfiles;
1419 for (i = 0; i < *nfiles; i++)
1421 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1424 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1425 strncpy(outputfiles[i].filename, buf, CPTSTRLEN-1);
1431 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1435 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1439 outputfiles[i].offset = ( ((gmx_off_t) offset_high) << 32 ) | ( (gmx_off_t) offset_low & mask );
1443 buf = outputfiles[i].filename;
1444 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1446 offset = outputfiles[i].offset;
1454 offset_low = (int) (offset & mask);
1455 offset_high = (int) ((offset >> 32) & mask);
1457 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1461 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1466 if (file_version >= 8)
1468 if (do_cpt_int(xd, "file_checksum_size", &(outputfiles[i].chksum_size),
1473 if (do_cpt_u_chars(xd, "file_checksum", 16, outputfiles[i].chksum, list) != 0)
1480 outputfiles[i].chksum_size = -1;
1487 void write_checkpoint(const char *fn, gmx_bool bNumberAndKeep,
1488 FILE *fplog, t_commrec *cr,
1489 int eIntegrator, int simulation_part,
1490 gmx_bool bExpanded, int elamstats,
1491 gmx_int64_t step, double t, t_state *state)
1501 char *fntemp; /* the temporary checkpoint file name */
1503 char timebuf[STRLEN];
1504 int nppnodes, npmenodes;
1505 char buf[1024], suffix[5+STEPSTRSIZE], sbuf[STEPSTRSIZE];
1506 gmx_file_position_t *outputfiles;
1509 int flags_eks, flags_enh, flags_dfh;
1512 if (DOMAINDECOMP(cr))
1514 nppnodes = cr->dd->nnodes;
1515 npmenodes = cr->npmenodes;
1523 #ifndef GMX_NO_RENAME
1524 /* make the new temporary filename */
1525 snew(fntemp, strlen(fn)+5+STEPSTRSIZE);
1527 fntemp[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1528 sprintf(suffix, "_%s%s", "step", gmx_step_str(step, sbuf));
1529 strcat(fntemp, suffix);
1530 strcat(fntemp, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1532 /* if we can't rename, we just overwrite the cpt file.
1533 * dangerous if interrupted.
1535 snew(fntemp, strlen(fn));
1539 gmx_ctime_r(&now, timebuf, STRLEN);
1543 fprintf(fplog, "Writing checkpoint, step %s at %s\n\n",
1544 gmx_step_str(step, buf), timebuf);
1547 /* Get offsets for open files */
1548 gmx_fio_get_output_file_positions(&outputfiles, &noutputfiles);
1550 fp = gmx_fio_open(fntemp, "w");
1552 if (state->ekinstate.bUpToDate)
1555 ((1<<eeksEKIN_N) | (1<<eeksEKINH) | (1<<eeksEKINF) |
1556 (1<<eeksEKINO) | (1<<eeksEKINSCALEF) | (1<<eeksEKINSCALEH) |
1557 (1<<eeksVSCALE) | (1<<eeksDEKINDL) | (1<<eeksMVCOS));
1565 if (state->enerhist.nsum > 0 || state->enerhist.nsum_sim > 0)
1567 flags_enh |= (1<<eenhENERGY_N);
1568 if (state->enerhist.nsum > 0)
1570 flags_enh |= ((1<<eenhENERGY_AVER) | (1<<eenhENERGY_SUM) |
1571 (1<<eenhENERGY_NSTEPS) | (1<<eenhENERGY_NSUM));
1573 if (state->enerhist.nsum_sim > 0)
1575 flags_enh |= ((1<<eenhENERGY_SUM_SIM) | (1<<eenhENERGY_NSTEPS_SIM) |
1576 (1<<eenhENERGY_NSUM_SIM));
1578 if (state->enerhist.dht)
1580 flags_enh |= ( (1<< eenhENERGY_DELTA_H_NN) |
1581 (1<< eenhENERGY_DELTA_H_LIST) |
1582 (1<< eenhENERGY_DELTA_H_STARTTIME) |
1583 (1<< eenhENERGY_DELTA_H_STARTLAMBDA) );
1589 flags_dfh = ((1<<edfhBEQUIL) | (1<<edfhNATLAMBDA) | (1<<edfhSUMWEIGHTS) | (1<<edfhSUMDG) |
1590 (1<<edfhTIJ) | (1<<edfhTIJEMP));
1593 flags_dfh |= ((1<<edfhWLDELTA) | (1<<edfhWLHISTO));
1595 if ((elamstats == elamstatsMINVAR) || (elamstats == elamstatsBARKER) || (elamstats == elamstatsMETROPOLIS))
1597 flags_dfh |= ((1<<edfhACCUMP) | (1<<edfhACCUMM) | (1<<edfhACCUMP2) | (1<<edfhACCUMM2)
1598 | (1<<edfhSUMMINVAR) | (1<<edfhSUMVAR));
1606 /* We can check many more things now (CPU, acceleration, etc), but
1607 * it is highly unlikely to have two separate builds with exactly
1608 * the same version, user, time, and build host!
1611 version = gmx_strdup(gmx_version());
1612 btime = gmx_strdup(BUILD_TIME);
1613 buser = gmx_strdup(BUILD_USER);
1614 bhost = gmx_strdup(BUILD_HOST);
1616 double_prec = GMX_CPT_BUILD_DP;
1617 fprog = gmx_strdup(Program());
1619 ftime = &(timebuf[0]);
1621 do_cpt_header(gmx_fio_getxdr(fp), FALSE, &file_version,
1622 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1623 &eIntegrator, &simulation_part, &step, &t, &nppnodes,
1624 DOMAINDECOMP(cr) ? cr->dd->nc : NULL, &npmenodes,
1625 &state->natoms, &state->ngtc, &state->nnhpres,
1626 &state->nhchainlength, &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
1627 &state->edsamstate.nED, &state->swapstate.eSwapCoords,
1636 if ((do_cpt_state(gmx_fio_getxdr(fp), FALSE, state->flags, state, NULL) < 0) ||
1637 (do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL) < 0) ||
1638 (do_cpt_enerhist(gmx_fio_getxdr(fp), FALSE, flags_enh, &state->enerhist, NULL) < 0) ||
1639 (do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL) < 0) ||
1640 (do_cpt_EDstate(gmx_fio_getxdr(fp), FALSE, &state->edsamstate, NULL) < 0) ||
1641 (do_cpt_swapstate(gmx_fio_getxdr(fp), FALSE, &state->swapstate, NULL) < 0) ||
1642 (do_cpt_files(gmx_fio_getxdr(fp), FALSE, &outputfiles, &noutputfiles, NULL,
1645 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1648 do_cpt_footer(gmx_fio_getxdr(fp), file_version);
1650 /* we really, REALLY, want to make sure to physically write the checkpoint,
1651 and all the files it depends on, out to disk. Because we've
1652 opened the checkpoint with gmx_fio_open(), it's in our list
1654 ret = gmx_fio_all_output_fsync();
1660 "Cannot fsync '%s'; maybe you are out of disk space?",
1661 gmx_fio_getname(ret));
1663 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV) == NULL)
1673 if (gmx_fio_close(fp) != 0)
1675 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1678 /* we don't move the checkpoint if the user specified they didn't want it,
1679 or if the fsyncs failed */
1680 #ifndef GMX_NO_RENAME
1681 if (!bNumberAndKeep && !ret)
1685 /* Rename the previous checkpoint file */
1687 buf[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1688 strcat(buf, "_prev");
1689 strcat(buf, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1691 /* we copy here so that if something goes wrong between now and
1692 * the rename below, there's always a state.cpt.
1693 * If renames are atomic (such as in POSIX systems),
1694 * this copying should be unneccesary.
1696 gmx_file_copy(fn, buf, FALSE);
1697 /* We don't really care if this fails:
1698 * there's already a new checkpoint.
1701 gmx_file_rename(fn, buf);
1704 if (gmx_file_rename(fntemp, fn) != 0)
1706 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
1709 #endif /* GMX_NO_RENAME */
1715 /*code for alternate checkpointing scheme. moved from top of loop over
1717 fcRequestCheckPoint();
1718 if (fcCheckPointParallel( cr->nodeid, NULL, 0) == 0)
1720 gmx_fatal( 3, __FILE__, __LINE__, "Checkpoint error on step %d\n", step );
1722 #endif /* end GMX_FAHCORE block */
1725 static void print_flag_mismatch(FILE *fplog, int sflags, int fflags)
1729 fprintf(fplog, "\nState entry mismatch between the simulation and the checkpoint file\n");
1730 fprintf(fplog, "Entries which are not present in the checkpoint file will not be updated\n");
1731 fprintf(fplog, " %24s %11s %11s\n", "", "simulation", "checkpoint");
1732 for (i = 0; i < estNR; i++)
1734 if ((sflags & (1<<i)) || (fflags & (1<<i)))
1736 fprintf(fplog, " %24s %11s %11s\n",
1738 (sflags & (1<<i)) ? " present " : "not present",
1739 (fflags & (1<<i)) ? " present " : "not present");
1744 static void check_int(FILE *fplog, const char *type, int p, int f, gmx_bool *mm)
1746 FILE *fp = fplog ? fplog : stderr;
1750 fprintf(fp, " %s mismatch,\n", type);
1751 fprintf(fp, " current program: %d\n", p);
1752 fprintf(fp, " checkpoint file: %d\n", f);
1758 static void check_string(FILE *fplog, const char *type, const char *p,
1759 const char *f, gmx_bool *mm)
1761 FILE *fp = fplog ? fplog : stderr;
1763 if (strcmp(p, f) != 0)
1765 fprintf(fp, " %s mismatch,\n", type);
1766 fprintf(fp, " current program: %s\n", p);
1767 fprintf(fp, " checkpoint file: %s\n", f);
1773 static void check_match(FILE *fplog,
1775 char *btime, char *buser, char *bhost, int double_prec,
1777 t_commrec *cr, int npp_f, int npme_f,
1778 ivec dd_nc, ivec dd_nc_f)
1781 gmx_bool mm = FALSE;
1782 gmx_bool patchlevel_differs = FALSE;
1783 gmx_bool version_differs = FALSE;
1785 check_string(fplog, "Version", gmx_version(), version, &mm);
1786 patchlevel_differs = mm;
1788 if (patchlevel_differs)
1790 /* Gromacs should be able to continue from checkpoints between
1791 * different patch level versions, but we do not guarantee
1792 * compatibility between different major/minor versions - check this.
1794 int gmx_major, gmx_minor;
1795 int cpt_major, cpt_minor;
1796 sscanf(gmx_version(), "VERSION %5d.%5d", &gmx_major, &gmx_minor);
1797 sscanf(version, "VERSION %5d.%5d", &cpt_major, &cpt_minor);
1798 version_differs = (gmx_major != cpt_major || gmx_minor != cpt_minor);
1801 check_string(fplog, "Build time", BUILD_TIME, btime, &mm);
1802 check_string(fplog, "Build user", BUILD_USER, buser, &mm);
1803 check_string(fplog, "Build host", BUILD_HOST, bhost, &mm);
1804 check_int (fplog, "Double prec.", GMX_CPT_BUILD_DP, double_prec, &mm);
1805 check_string(fplog, "Program name", Program(), fprog, &mm);
1807 check_int (fplog, "#ranks", cr->nnodes, npp_f+npme_f, &mm);
1810 check_int (fplog, "#PME-ranks", cr->npmenodes, npme_f, &mm);
1813 if (cr->npmenodes >= 0)
1815 npp -= cr->npmenodes;
1819 check_int (fplog, "#DD-cells[x]", dd_nc[XX], dd_nc_f[XX], &mm);
1820 check_int (fplog, "#DD-cells[y]", dd_nc[YY], dd_nc_f[YY], &mm);
1821 check_int (fplog, "#DD-cells[z]", dd_nc[ZZ], dd_nc_f[ZZ], &mm);
1827 const char msg_version_difference[] =
1828 "The current Gromacs major & minor version are not identical to those that\n"
1829 "generated the checkpoint file. In principle Gromacs does not support\n"
1830 "continuation from checkpoints between different versions, so we advise\n"
1831 "against this. If you still want to try your luck we recommend that you use\n"
1832 "the -noappend flag to keep your output files from the two versions separate.\n"
1833 "This might also work around errors where the output fields in the energy\n"
1834 "file have changed between the different major & minor versions.\n";
1836 const char msg_mismatch_notice[] =
1837 "Gromacs patchlevel, binary or parallel settings differ from previous run.\n"
1838 "Continuation is exact, but not guaranteed to be binary identical.\n";
1840 const char msg_logdetails[] =
1841 "See the log file for details.\n";
1843 if (version_differs)
1845 fprintf(stderr, "%s%s\n", msg_version_difference, fplog ? msg_logdetails : "");
1849 fprintf(fplog, "%s\n", msg_version_difference);
1854 /* Major & minor versions match at least, but something is different. */
1855 fprintf(stderr, "%s%s\n", msg_mismatch_notice, fplog ? msg_logdetails : "");
1858 fprintf(fplog, "%s\n", msg_mismatch_notice);
1864 static void read_checkpoint(const char *fn, FILE **pfplog,
1865 t_commrec *cr, ivec dd_nc,
1866 int eIntegrator, int *init_fep_state, gmx_int64_t *step, double *t,
1867 t_state *state, gmx_bool *bReadEkin,
1868 int *simulation_part,
1869 gmx_bool bAppendOutputFiles, gmx_bool bForceAppend)
1874 char *version, *btime, *buser, *bhost, *fprog, *ftime;
1876 char buf[STEPSTRSIZE];
1877 int eIntegrator_f, nppnodes_f, npmenodes_f;
1879 int natoms, ngtc, nnhpres, nhchainlength, nlambda, fflags, flags_eks, flags_enh, flags_dfh;
1882 gmx_file_position_t *outputfiles;
1884 t_fileio *chksum_file;
1885 FILE * fplog = *pfplog;
1886 unsigned char digest[16];
1887 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1888 struct flock fl; /* don't initialize here: the struct order is OS
1892 const char *int_warn =
1893 "WARNING: The checkpoint file was generated with integrator %s,\n"
1894 " while the simulation uses integrator %s\n\n";
1896 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1897 fl.l_type = F_WRLCK;
1898 fl.l_whence = SEEK_SET;
1904 fp = gmx_fio_open(fn, "r");
1905 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
1906 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1907 &eIntegrator_f, simulation_part, step, t,
1908 &nppnodes_f, dd_nc_f, &npmenodes_f,
1909 &natoms, &ngtc, &nnhpres, &nhchainlength, &nlambda,
1910 &fflags, &flags_eks, &flags_enh, &flags_dfh,
1911 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
1913 if (bAppendOutputFiles &&
1914 file_version >= 13 && double_prec != GMX_CPT_BUILD_DP)
1916 gmx_fatal(FARGS, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
1919 if (cr == NULL || MASTER(cr))
1921 fprintf(stderr, "\nReading checkpoint file %s generated: %s\n\n",
1925 /* This will not be written if we do appending, since fplog is still NULL then */
1928 fprintf(fplog, "\n");
1929 fprintf(fplog, "Reading checkpoint file %s\n", fn);
1930 fprintf(fplog, " file generated by: %s\n", fprog);
1931 fprintf(fplog, " file generated at: %s\n", ftime);
1932 fprintf(fplog, " GROMACS build time: %s\n", btime);
1933 fprintf(fplog, " GROMACS build user: %s\n", buser);
1934 fprintf(fplog, " GROMACS build host: %s\n", bhost);
1935 fprintf(fplog, " GROMACS double prec.: %d\n", double_prec);
1936 fprintf(fplog, " simulation part #: %d\n", *simulation_part);
1937 fprintf(fplog, " step: %s\n", gmx_step_str(*step, buf));
1938 fprintf(fplog, " time: %f\n", *t);
1939 fprintf(fplog, "\n");
1942 if (natoms != state->natoms)
1944 gmx_fatal(FARGS, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", natoms, state->natoms);
1946 if (ngtc != state->ngtc)
1948 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);
1950 if (nnhpres != state->nnhpres)
1952 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);
1955 if (nlambda != state->dfhist.nlambda)
1957 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);
1960 init_gtc_state(state, state->ngtc, state->nnhpres, nhchainlength); /* need to keep this here to keep the tpr format working */
1961 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
1963 if (eIntegrator_f != eIntegrator)
1967 fprintf(stderr, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1969 if (bAppendOutputFiles)
1972 "Output file appending requested, but input/checkpoint integrators do not match.\n"
1973 "Stopping the run to prevent you from ruining all your data...\n"
1974 "If you _really_ know what you are doing, try with the -noappend option.\n");
1978 fprintf(fplog, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1986 else if (cr->nnodes == nppnodes_f + npmenodes_f)
1988 if (cr->npmenodes < 0)
1990 cr->npmenodes = npmenodes_f;
1992 int nppnodes = cr->nnodes - cr->npmenodes;
1993 if (nppnodes == nppnodes_f)
1995 for (d = 0; d < DIM; d++)
1999 dd_nc[d] = dd_nc_f[d];
2005 if (fflags != state->flags)
2010 if (bAppendOutputFiles)
2013 "Output file appending requested, but input and checkpoint states are not identical.\n"
2014 "Stopping the run to prevent you from ruining all your data...\n"
2015 "You can try with the -noappend option, and get more info in the log file.\n");
2018 if (getenv("GMX_ALLOW_CPT_MISMATCH") == NULL)
2020 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");
2025 "WARNING: The checkpoint state entries do not match the simulation,\n"
2026 " see the log file for details\n\n");
2032 print_flag_mismatch(fplog, state->flags, fflags);
2039 check_match(fplog, version, btime, buser, bhost, double_prec, fprog,
2040 cr, nppnodes_f, npmenodes_f, dd_nc, dd_nc_f);
2043 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, fflags, state, NULL);
2044 *init_fep_state = state->fep_state; /* there should be a better way to do this than setting it here.
2045 Investigate for 5.0. */
2050 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2055 *bReadEkin = ((flags_eks & (1<<eeksEKINH)) || (flags_eks & (1<<eeksEKINF)) || (flags_eks & (1<<eeksEKINO)) ||
2056 ((flags_eks & (1<<eeksEKINSCALEF)) | (flags_eks & (1<<eeksEKINSCALEH)) | (flags_eks & (1<<eeksVSCALE))));
2058 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2059 flags_enh, &state->enerhist, NULL);
2065 if (file_version < 6)
2067 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.";
2069 fprintf(stderr, "\nWARNING: %s\n\n", warn);
2072 fprintf(fplog, "\nWARNING: %s\n\n", warn);
2074 state->enerhist.nsum = *step;
2075 state->enerhist.nsum_sim = *step;
2078 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2084 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2090 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2096 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, NULL, file_version);
2102 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2107 if (gmx_fio_close(fp) != 0)
2109 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2118 /* If the user wants to append to output files,
2119 * we use the file pointer positions of the output files stored
2120 * in the checkpoint file and truncate the files such that any frames
2121 * written after the checkpoint time are removed.
2122 * All files are md5sum checked such that we can be sure that
2123 * we do not truncate other (maybe imprortant) files.
2125 if (bAppendOutputFiles)
2127 if (fn2ftp(outputfiles[0].filename) != efLOG)
2129 /* make sure first file is log file so that it is OK to use it for
2132 gmx_fatal(FARGS, "The first output file should always be the log "
2133 "file but instead is: %s. Cannot do appending because of this condition.", outputfiles[0].filename);
2135 for (i = 0; i < nfiles; i++)
2137 if (outputfiles[i].offset < 0)
2139 gmx_fatal(FARGS, "The original run wrote a file called '%s' which "
2140 "is larger than 2 GB, but mdrun did not support large file"
2141 " offsets. Can not append. Run mdrun with -noappend",
2142 outputfiles[i].filename);
2145 chksum_file = gmx_fio_open(outputfiles[i].filename, "a");
2148 chksum_file = gmx_fio_open(outputfiles[i].filename, "r+");
2153 /* Note that there are systems where the lock operation
2154 * will succeed, but a second process can also lock the file.
2155 * We should probably try to detect this.
2157 #if defined __native_client__
2161 #elif defined GMX_NATIVE_WINDOWS
2162 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX) == -1)
2164 if (fcntl(fileno(gmx_fio_getfp(chksum_file)), F_SETLK, &fl) == -1)
2167 if (errno == ENOSYS)
2171 gmx_fatal(FARGS, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2175 fprintf(stderr, "\nNOTE: File locking is not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2178 fprintf(fplog, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2182 else if (errno == EACCES || errno == EAGAIN)
2184 gmx_fatal(FARGS, "Failed to lock: %s. Already running "
2185 "simulation?", outputfiles[i].filename);
2189 gmx_fatal(FARGS, "Failed to lock: %s. %s.",
2190 outputfiles[i].filename, strerror(errno));
2195 /* compute md5 chksum */
2196 if (outputfiles[i].chksum_size != -1)
2198 if (gmx_fio_get_file_md5(chksum_file, outputfiles[i].offset,
2199 digest) != outputfiles[i].chksum_size) /*at the end of the call the file position is at the end of the file*/
2201 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.",
2202 outputfiles[i].chksum_size,
2203 outputfiles[i].filename);
2206 if (i == 0) /*log file needs to be seeked in case we need to truncate (other files are truncated below)*/
2208 if (gmx_fio_seek(chksum_file, outputfiles[i].offset))
2210 gmx_fatal(FARGS, "Seek error! Failed to truncate log-file: %s.", strerror(errno));
2215 if (i == 0) /*open log file here - so that lock is never lifted
2216 after chksum is calculated */
2218 *pfplog = gmx_fio_getfp(chksum_file);
2222 gmx_fio_close(chksum_file);
2225 /* compare md5 chksum */
2226 if (outputfiles[i].chksum_size != -1 &&
2227 memcmp(digest, outputfiles[i].chksum, 16) != 0)
2231 fprintf(debug, "chksum for %s: ", outputfiles[i].filename);
2232 for (j = 0; j < 16; j++)
2234 fprintf(debug, "%02x", digest[j]);
2236 fprintf(debug, "\n");
2238 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.",
2239 outputfiles[i].filename);
2244 if (i != 0) /*log file is already seeked to correct position */
2246 #ifdef GMX_NATIVE_WINDOWS
2247 rc = gmx_wintruncate(outputfiles[i].filename, outputfiles[i].offset);
2249 rc = truncate(outputfiles[i].filename, outputfiles[i].offset);
2253 gmx_fatal(FARGS, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfiles[i].filename);
2263 void load_checkpoint(const char *fn, FILE **fplog,
2264 t_commrec *cr, ivec dd_nc,
2265 t_inputrec *ir, t_state *state,
2266 gmx_bool *bReadEkin,
2267 gmx_bool bAppend, gmx_bool bForceAppend)
2274 /* Read the state from the checkpoint file */
2275 read_checkpoint(fn, fplog,
2277 ir->eI, &(ir->fepvals->init_fep_state), &step, &t, state, bReadEkin,
2278 &ir->simulation_part, bAppend, bForceAppend);
2282 gmx_bcast(sizeof(cr->npmenodes), &cr->npmenodes, cr);
2283 gmx_bcast(DIM*sizeof(dd_nc[0]), dd_nc, cr);
2284 gmx_bcast(sizeof(step), &step, cr);
2285 gmx_bcast(sizeof(*bReadEkin), bReadEkin, cr);
2287 ir->bContinuation = TRUE;
2288 if (ir->nsteps >= 0)
2290 ir->nsteps += ir->init_step - step;
2292 ir->init_step = step;
2293 ir->simulation_part += 1;
2296 static void read_checkpoint_data(t_fileio *fp, int *simulation_part,
2297 gmx_int64_t *step, double *t, t_state *state,
2298 int *nfiles, gmx_file_position_t **outputfiles)
2301 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2306 int flags_eks, flags_enh, flags_dfh;
2308 gmx_file_position_t *files_loc = NULL;
2311 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2312 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2313 &eIntegrator, simulation_part, step, t, &nppnodes, dd_nc, &npme,
2314 &state->natoms, &state->ngtc, &state->nnhpres, &state->nhchainlength,
2315 &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
2316 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
2318 do_cpt_state(gmx_fio_getxdr(fp), TRUE, state->flags, state, NULL);
2323 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2328 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2329 flags_enh, &state->enerhist, NULL);
2334 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2340 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2346 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2352 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE,
2353 outputfiles != NULL ? outputfiles : &files_loc,
2354 outputfiles != NULL ? nfiles : &nfiles_loc,
2355 NULL, file_version);
2356 if (files_loc != NULL)
2366 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2380 read_checkpoint_state(const char *fn, int *simulation_part,
2381 gmx_int64_t *step, double *t, t_state *state)
2385 fp = gmx_fio_open(fn, "r");
2386 read_checkpoint_data(fp, simulation_part, step, t, state, NULL, NULL);
2387 if (gmx_fio_close(fp) != 0)
2389 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2393 void read_checkpoint_trxframe(t_fileio *fp, t_trxframe *fr)
2395 /* This next line is nasty because the sub-structures of t_state
2396 * cannot be assumed to be zeroed (or even initialized in ways the
2397 * rest of the code might assume). Using snew would be better, but
2398 * this will all go away for 5.0. */
2400 int simulation_part;
2404 init_state(&state, 0, 0, 0, 0, 0);
2406 read_checkpoint_data(fp, &simulation_part, &step, &t, &state, NULL, NULL);
2408 fr->natoms = state.natoms;
2411 fr->step = gmx_int64_to_int(step,
2412 "conversion of checkpoint to trajectory");
2416 fr->lambda = state.lambda[efptFEP];
2417 fr->fep_state = state.fep_state;
2419 fr->bX = (state.flags & (1<<estX));
2425 fr->bV = (state.flags & (1<<estV));
2432 fr->bBox = (state.flags & (1<<estBOX));
2435 copy_mat(state.box, fr->box);
2440 void list_checkpoint(const char *fn, FILE *out)
2444 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2446 int eIntegrator, simulation_part, nppnodes, npme;
2451 int flags_eks, flags_enh, flags_dfh;
2453 gmx_file_position_t *outputfiles;
2456 init_state(&state, -1, -1, -1, -1, 0);
2458 fp = gmx_fio_open(fn, "r");
2459 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2460 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2461 &eIntegrator, &simulation_part, &step, &t, &nppnodes, dd_nc, &npme,
2462 &state.natoms, &state.ngtc, &state.nnhpres, &state.nhchainlength,
2463 &(state.dfhist.nlambda), &state.flags,
2464 &flags_eks, &flags_enh, &flags_dfh, &state.edsamstate.nED,
2465 &state.swapstate.eSwapCoords, out);
2466 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, state.flags, &state, out);
2471 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state.ekinstate, out);
2476 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2477 flags_enh, &state.enerhist, out);
2481 ret = do_cpt_df_hist(gmx_fio_getxdr(fp),
2482 flags_dfh, &state.dfhist, out);
2487 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state.edsamstate, out);
2492 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state.swapstate, out);
2497 do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, out, file_version);
2502 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2509 if (gmx_fio_close(fp) != 0)
2511 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2518 static gmx_bool exist_output_file(const char *fnm_cp, int nfile, const t_filenm fnm[])
2522 /* Check if the output file name stored in the checkpoint file
2523 * is one of the output file names of mdrun.
2527 !(is_output(&fnm[i]) && strcmp(fnm_cp, fnm[i].fns[0]) == 0))
2532 return (i < nfile && gmx_fexist(fnm_cp));
2535 /* This routine cannot print tons of data, since it is called before the log file is opened. */
2536 gmx_bool read_checkpoint_simulation_part(const char *filename, int *simulation_part,
2537 gmx_int64_t *cpt_step, t_commrec *cr,
2538 gmx_bool bAppendReq,
2539 int nfile, const t_filenm fnm[],
2540 const char *part_suffix, gmx_bool *bAddPart)
2543 gmx_int64_t step = 0;
2545 /* This next line is nasty because the sub-structures of t_state
2546 * cannot be assumed to be zeroed (or even initialized in ways the
2547 * rest of the code might assume). Using snew would be better, but
2548 * this will all go away for 5.0. */
2551 gmx_file_position_t *outputfiles;
2554 char *fn, suf_up[STRLEN];
2560 if (!gmx_fexist(filename) || (!(fp = gmx_fio_open(filename, "r")) ))
2562 *simulation_part = 0;
2566 init_state(&state, 0, 0, 0, 0, 0);
2568 read_checkpoint_data(fp, simulation_part, &step, &t, &state,
2569 &nfiles, &outputfiles);
2570 if (gmx_fio_close(fp) != 0)
2572 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2579 for (f = 0; f < nfiles; f++)
2581 if (exist_output_file(outputfiles[f].filename, nfile, fnm))
2586 if (nexist == nfiles)
2588 bAppend = bAppendReq;
2590 else if (nexist > 0)
2593 "Output file appending has been requested,\n"
2594 "but some output files listed in the checkpoint file %s\n"
2595 "are not present or are named differently by the current program:\n",
2597 fprintf(stderr, "output files present:");
2598 for (f = 0; f < nfiles; f++)
2600 if (exist_output_file(outputfiles[f].filename,
2603 fprintf(stderr, " %s", outputfiles[f].filename);
2606 fprintf(stderr, "\n");
2607 fprintf(stderr, "output files not present or named differently:");
2608 for (f = 0; f < nfiles; f++)
2610 if (!exist_output_file(outputfiles[f].filename,
2613 fprintf(stderr, " %s", outputfiles[f].filename);
2616 fprintf(stderr, "\n");
2618 gmx_fatal(FARGS, "File appending requested, but %d of the %d output files are not present or are named differently", nfiles-nexist, nfiles);
2626 gmx_fatal(FARGS, "File appending requested, but no output file information is stored in the checkpoint file");
2628 fn = outputfiles[0].filename;
2629 if (strlen(fn) < 4 ||
2630 gmx_strcasecmp(fn+strlen(fn)-4, ftp2ext(efLOG)) == 0)
2632 gmx_fatal(FARGS, "File appending requested, but the log file is not the first file listed in the checkpoint file");
2634 /* Set bAddPart to whether the suffix string '.part' is present
2635 * in the log file name.
2637 strcpy(suf_up, part_suffix);
2639 *bAddPart = (strstr(fn, part_suffix) != NULL ||
2640 strstr(fn, suf_up) != NULL);
2648 gmx_bcast(sizeof(*simulation_part), simulation_part, cr);
2650 if (*simulation_part > 0 && bAppendReq)
2652 gmx_bcast(sizeof(bAppend), &bAppend, cr);
2653 gmx_bcast(sizeof(*bAddPart), bAddPart, cr);
2656 if (NULL != cpt_step)