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. */
40 #include "gromacs/legacyheaders/checkpoint.h"
51 #ifdef HAVE_SYS_TIME_H
57 #ifdef GMX_NATIVE_WINDOWS
60 #include <sys/locking.h>
63 #include "buildinfo.h"
64 #include "gromacs/fileio/filenm.h"
65 #include "gromacs/fileio/gmxfio.h"
66 #include "gromacs/fileio/xdr_datatype.h"
67 #include "gromacs/fileio/xdrf.h"
68 #include "gromacs/legacyheaders/copyrite.h"
69 #include "gromacs/legacyheaders/names.h"
70 #include "gromacs/legacyheaders/network.h"
71 #include "gromacs/legacyheaders/txtdump.h"
72 #include "gromacs/legacyheaders/typedefs.h"
73 #include "gromacs/legacyheaders/types/commrec.h"
74 #include "gromacs/math/vec.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/futil.h"
80 #include "gromacs/utility/smalloc.h"
86 #define CPT_MAGIC1 171817
87 #define CPT_MAGIC2 171819
88 #define CPTSTRLEN 1024
91 #define GMX_CPT_BUILD_DP 1
93 #define GMX_CPT_BUILD_DP 0
96 /* cpt_version should normally only be changed
97 * when the header of footer format changes.
98 * The state data format itself is backward and forward compatible.
99 * But old code can not read a new entry that is present in the file
100 * (but can read a new format when new entries are not present).
102 static const int cpt_version = 16;
105 const char *est_names[estNR] =
108 "box", "box-rel", "box-v", "pres_prev",
109 "nosehoover-xi", "thermostat-integral",
110 "x", "v", "SDx", "CGp", "LD-rng", "LD-rng-i",
111 "disre_initf", "disre_rm3tav",
112 "orire_initf", "orire_Dtav",
113 "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng", "MC-rng-i"
117 eeksEKIN_N, eeksEKINH, eeksDEKINDL, eeksMVCOS, eeksEKINF, eeksEKINO, eeksEKINSCALEF, eeksEKINSCALEH, eeksVSCALE, eeksEKINTOTAL, eeksNR
120 const char *eeks_names[eeksNR] =
122 "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
123 "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
127 eenhENERGY_N, eenhENERGY_AVER, eenhENERGY_SUM, eenhENERGY_NSUM,
128 eenhENERGY_SUM_SIM, eenhENERGY_NSUM_SIM,
129 eenhENERGY_NSTEPS, eenhENERGY_NSTEPS_SIM,
130 eenhENERGY_DELTA_H_NN,
131 eenhENERGY_DELTA_H_LIST,
132 eenhENERGY_DELTA_H_STARTTIME,
133 eenhENERGY_DELTA_H_STARTLAMBDA,
137 const char *eenh_names[eenhNR] =
139 "energy_n", "energy_aver", "energy_sum", "energy_nsum",
140 "energy_sum_sim", "energy_nsum_sim",
141 "energy_nsteps", "energy_nsteps_sim",
143 "energy_delta_h_list",
144 "energy_delta_h_start_time",
145 "energy_delta_h_start_lambda"
148 /* free energy history variables -- need to be preserved over checkpoint */
150 edfhBEQUIL, edfhNATLAMBDA, edfhWLHISTO, edfhWLDELTA, edfhSUMWEIGHTS, edfhSUMDG, edfhSUMMINVAR, edfhSUMVAR,
151 edfhACCUMP, edfhACCUMM, edfhACCUMP2, edfhACCUMM2, edfhTIJ, edfhTIJEMP, edfhNR
153 /* free energy history variable names */
154 const char *edfh_names[edfhNR] =
156 "bEquilibrated", "N_at_state", "Wang-Landau Histogram", "Wang-Landau Delta", "Weights", "Free Energies", "minvar", "variance",
157 "accumulated_plus", "accumulated_minus", "accumulated_plus_2", "accumulated_minus_2", "Tij", "Tij_empirical"
160 #ifdef GMX_NATIVE_WINDOWS
162 gmx_wintruncate(const char *filename, __int64 size)
165 /*we do this elsewhere*/
170 fp = fopen(filename, "rb+");
178 return _chsize_s( fileno(fp), size);
180 return _chsize( fileno(fp), size);
188 ecprREAL, ecprRVEC, ecprMATRIX
192 cptpEST, cptpEEKS, cptpEENH, cptpEDFH
194 /* enums for the different components of checkpoint variables, replacing the hard coded ones.
195 cptpEST - state variables.
196 cptpEEKS - Kinetic energy state variables.
197 cptpEENH - Energy history state variables.
198 cptpEDFH - free energy history variables.
202 static const char *st_names(int cptp, int ecpt)
206 case cptpEST: return est_names [ecpt]; break;
207 case cptpEEKS: return eeks_names[ecpt]; break;
208 case cptpEENH: return eenh_names[ecpt]; break;
209 case cptpEDFH: return edfh_names[ecpt]; break;
215 static void cp_warning(FILE *fp)
217 fprintf(fp, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
220 static void cp_error()
222 gmx_fatal(FARGS, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
225 static void do_cpt_string_err(XDR *xd, gmx_bool bRead, const char *desc, char **s, FILE *list)
233 res = xdr_string(xd, s, CPTSTRLEN);
240 fprintf(list, "%s = %s\n", desc, *s);
245 static int do_cpt_int(XDR *xd, const char *desc, int *i, FILE *list)
249 res = xdr_int(xd, i);
256 fprintf(list, "%s = %d\n", desc, *i);
261 static int do_cpt_u_chars(XDR *xd, const char *desc, int n, unsigned char *i, FILE *list)
267 fprintf(list, "%s = ", desc);
269 for (j = 0; j < n && res; j++)
271 res &= xdr_u_char(xd, &i[j]);
274 fprintf(list, "%02x", i[j]);
289 static void do_cpt_int_err(XDR *xd, const char *desc, int *i, FILE *list)
291 if (do_cpt_int(xd, desc, i, list) < 0)
297 static void do_cpt_step_err(XDR *xd, const char *desc, gmx_int64_t *i, FILE *list)
300 char buf[STEPSTRSIZE];
302 res = xdr_int64(xd, i);
309 fprintf(list, "%s = %s\n", desc, gmx_step_str(*i, buf));
313 static void do_cpt_double_err(XDR *xd, const char *desc, double *f, FILE *list)
317 res = xdr_double(xd, f);
324 fprintf(list, "%s = %f\n", desc, *f);
328 static void do_cpt_real_err(XDR *xd, real *f)
333 res = xdr_double(xd, f);
335 res = xdr_float(xd, f);
343 static void do_cpt_n_rvecs_err(XDR *xd, const char *desc, int n, rvec f[], FILE *list)
347 for (i = 0; i < n; i++)
349 for (j = 0; j < DIM; j++)
351 do_cpt_real_err(xd, &f[i][j]);
357 pr_rvecs(list, 0, desc, f, n);
361 /* If nval >= 0, nval is used; on read this should match the passed value.
362 * If nval n<0, *nptr is used; on read the value is stored in nptr
364 static int do_cpte_reals_low(XDR *xd, int cptp, int ecpt, int sflags,
365 int nval, int *nptr, real **v,
366 FILE *list, int erealtype)
370 int dtc = xdr_datatype_float;
372 int dtc = xdr_datatype_double;
374 real *vp, *va = NULL;
389 gmx_incons("*ntpr=NULL in do_cpte_reals_low");
394 res = xdr_int(xd, &nf);
405 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), nval, nf);
414 res = xdr_int(xd, &dt);
421 fprintf(stderr, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
422 st_names(cptp, ecpt), xdr_datatype_names[dtc],
423 xdr_datatype_names[dt]);
425 if (list || !(sflags & (1<<ecpt)))
438 if (dt == xdr_datatype_float)
440 if (dtc == xdr_datatype_float)
448 res = xdr_vector(xd, (char *)vf, nf,
449 (unsigned int)sizeof(float), (xdrproc_t)xdr_float);
454 if (dtc != xdr_datatype_float)
456 for (i = 0; i < nf; i++)
465 if (dtc == xdr_datatype_double)
467 /* cppcheck-suppress invalidPointerCast
468 * Only executed if real is anyhow double */
475 res = xdr_vector(xd, (char *)vd, nf,
476 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
481 if (dtc != xdr_datatype_double)
483 for (i = 0; i < nf; i++)
496 pr_reals(list, 0, st_names(cptp, ecpt), vp, nf);
499 pr_rvecs(list, 0, st_names(cptp, ecpt), (rvec *)vp, nf/3);
502 gmx_incons("Unknown checkpoint real type");
514 /* This function stores n along with the reals for reading,
515 * but on reading it assumes that n matches the value in the checkpoint file,
516 * a fatal error is generated when this is not the case.
518 static int do_cpte_reals(XDR *xd, int cptp, int ecpt, int sflags,
519 int n, real **v, FILE *list)
521 return do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, v, list, ecprREAL);
524 /* This function does the same as do_cpte_reals,
525 * except that on reading it ignores the passed value of *n
526 * and stored the value read from the checkpoint file in *n.
528 static int do_cpte_n_reals(XDR *xd, int cptp, int ecpt, int sflags,
529 int *n, real **v, FILE *list)
531 return do_cpte_reals_low(xd, cptp, ecpt, sflags, -1, n, v, list, ecprREAL);
534 static int do_cpte_real(XDR *xd, int cptp, int ecpt, int sflags,
537 return do_cpte_reals_low(xd, cptp, ecpt, sflags, 1, NULL, &r, list, ecprREAL);
540 static int do_cpte_ints(XDR *xd, int cptp, int ecpt, int sflags,
541 int n, int **v, FILE *list)
544 int dtc = xdr_datatype_int;
549 res = xdr_int(xd, &nf);
554 if (list == NULL && v != NULL && nf != n)
556 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
559 res = xdr_int(xd, &dt);
566 gmx_fatal(FARGS, "Type mismatch for state entry %s, code type is %s, file type is %s\n",
567 st_names(cptp, ecpt), xdr_datatype_names[dtc],
568 xdr_datatype_names[dt]);
570 if (list || !(sflags & (1<<ecpt)) || v == NULL)
583 res = xdr_vector(xd, (char *)vp, nf,
584 (unsigned int)sizeof(int), (xdrproc_t)xdr_int);
591 pr_ivec(list, 0, st_names(cptp, ecpt), vp, nf, TRUE);
601 static int do_cpte_int(XDR *xd, int cptp, int ecpt, int sflags,
604 return do_cpte_ints(xd, cptp, ecpt, sflags, 1, &i, list);
607 static int do_cpte_doubles(XDR *xd, int cptp, int ecpt, int sflags,
608 int n, double **v, FILE *list)
611 int dtc = xdr_datatype_double;
612 double *vp, *va = NULL;
616 res = xdr_int(xd, &nf);
621 if (list == NULL && nf != n)
623 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
626 res = xdr_int(xd, &dt);
633 gmx_fatal(FARGS, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
634 st_names(cptp, ecpt), xdr_datatype_names[dtc],
635 xdr_datatype_names[dt]);
637 if (list || !(sflags & (1<<ecpt)))
650 res = xdr_vector(xd, (char *)vp, nf,
651 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
658 pr_doubles(list, 0, st_names(cptp, ecpt), vp, nf);
668 static int do_cpte_double(XDR *xd, int cptp, int ecpt, int sflags,
669 double *r, FILE *list)
671 return do_cpte_doubles(xd, cptp, ecpt, sflags, 1, &r, list);
675 static int do_cpte_rvecs(XDR *xd, int cptp, int ecpt, int sflags,
676 int n, rvec **v, FILE *list)
678 return do_cpte_reals_low(xd, cptp, ecpt, sflags,
679 n*DIM, NULL, (real **)v, list, ecprRVEC);
682 static int do_cpte_matrix(XDR *xd, int cptp, int ecpt, int sflags,
683 matrix v, FILE *list)
688 vr = (real *)&(v[0][0]);
689 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
690 DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
692 if (list && ret == 0)
694 pr_rvecs(list, 0, st_names(cptp, ecpt), v, DIM);
701 static int do_cpte_nmatrix(XDR *xd, int cptp, int ecpt, int sflags,
702 int n, real **v, FILE *list)
706 char name[CPTSTRLEN];
713 for (i = 0; i < n; i++)
715 reti = do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, &(v[i]), NULL, ecprREAL);
716 if (list && reti == 0)
718 sprintf(name, "%s[%d]", st_names(cptp, ecpt), i);
719 pr_reals(list, 0, name, v[i], n);
729 static int do_cpte_matrices(XDR *xd, int cptp, int ecpt, int sflags,
730 int n, matrix **v, FILE *list)
733 matrix *vp, *va = NULL;
739 res = xdr_int(xd, &nf);
744 if (list == NULL && nf != n)
746 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
748 if (list || !(sflags & (1<<ecpt)))
761 snew(vr, nf*DIM*DIM);
762 for (i = 0; i < nf; i++)
764 for (j = 0; j < DIM; j++)
766 for (k = 0; k < DIM; k++)
768 vr[(i*DIM+j)*DIM+k] = vp[i][j][k];
772 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
773 nf*DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
774 for (i = 0; i < nf; i++)
776 for (j = 0; j < DIM; j++)
778 for (k = 0; k < DIM; k++)
780 vp[i][j][k] = vr[(i*DIM+j)*DIM+k];
786 if (list && ret == 0)
788 for (i = 0; i < nf; i++)
790 pr_rvecs(list, 0, st_names(cptp, ecpt), vp[i], DIM);
801 static void do_cpt_header(XDR *xd, gmx_bool bRead, int *file_version,
802 char **version, char **btime, char **buser, char **bhost,
804 char **fprog, char **ftime,
805 int *eIntegrator, int *simulation_part,
806 gmx_int64_t *step, double *t,
807 int *nnodes, int *dd_nc, int *npme,
808 int *natoms, int *ngtc, int *nnhpres, int *nhchainlength,
809 int *nlambda, int *flags_state,
810 int *flags_eks, int *flags_enh, int *flags_dfh,
811 int *nED, int *eSwapCoords,
827 res = xdr_int(xd, &magic);
830 gmx_fatal(FARGS, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
832 if (magic != CPT_MAGIC1)
834 gmx_fatal(FARGS, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
835 "The checkpoint file is corrupted or not a checkpoint file",
841 gmx_gethostname(fhost, 255);
843 do_cpt_string_err(xd, bRead, "GROMACS version", version, list);
844 do_cpt_string_err(xd, bRead, "GROMACS build time", btime, list);
845 do_cpt_string_err(xd, bRead, "GROMACS build user", buser, list);
846 do_cpt_string_err(xd, bRead, "GROMACS build host", bhost, list);
847 do_cpt_string_err(xd, bRead, "generating program", fprog, list);
848 do_cpt_string_err(xd, bRead, "generation time", ftime, list);
849 *file_version = cpt_version;
850 do_cpt_int_err(xd, "checkpoint file version", file_version, list);
851 if (*file_version > cpt_version)
853 gmx_fatal(FARGS, "Attempting to read a checkpoint file of version %d with code of version %d\n", *file_version, cpt_version);
855 if (*file_version >= 13)
857 do_cpt_int_err(xd, "GROMACS double precision", double_prec, list);
863 if (*file_version >= 12)
865 do_cpt_string_err(xd, bRead, "generating host", &fhost, list);
871 do_cpt_int_err(xd, "#atoms", natoms, list);
872 do_cpt_int_err(xd, "#T-coupling groups", ngtc, list);
873 if (*file_version >= 10)
875 do_cpt_int_err(xd, "#Nose-Hoover T-chains", nhchainlength, list);
881 if (*file_version >= 11)
883 do_cpt_int_err(xd, "#Nose-Hoover T-chains for barostat ", nnhpres, list);
889 if (*file_version >= 14)
891 do_cpt_int_err(xd, "# of total lambda states ", nlambda, list);
897 do_cpt_int_err(xd, "integrator", eIntegrator, list);
898 if (*file_version >= 3)
900 do_cpt_int_err(xd, "simulation part #", simulation_part, list);
904 *simulation_part = 1;
906 if (*file_version >= 5)
908 do_cpt_step_err(xd, "step", step, list);
912 do_cpt_int_err(xd, "step", &idum, list);
915 do_cpt_double_err(xd, "t", t, list);
916 do_cpt_int_err(xd, "#PP-ranks", nnodes, list);
918 do_cpt_int_err(xd, "dd_nc[x]", dd_nc ? &(dd_nc[0]) : &idum, list);
919 do_cpt_int_err(xd, "dd_nc[y]", dd_nc ? &(dd_nc[1]) : &idum, list);
920 do_cpt_int_err(xd, "dd_nc[z]", dd_nc ? &(dd_nc[2]) : &idum, list);
921 do_cpt_int_err(xd, "#PME-only ranks", npme, list);
922 do_cpt_int_err(xd, "state flags", flags_state, list);
923 if (*file_version >= 4)
925 do_cpt_int_err(xd, "ekin data flags", flags_eks, list);
926 do_cpt_int_err(xd, "energy history flags", flags_enh, list);
931 *flags_enh = (*flags_state >> (estORIRE_DTAV+1));
932 *flags_state = (*flags_state & ~((1<<(estORIRE_DTAV+1)) |
933 (1<<(estORIRE_DTAV+2)) |
934 (1<<(estORIRE_DTAV+3))));
936 if (*file_version >= 14)
938 do_cpt_int_err(xd, "df history flags", flags_dfh, list);
945 if (*file_version >= 15)
947 do_cpt_int_err(xd, "ED data sets", nED, list);
953 if (*file_version >= 16)
955 do_cpt_int_err(xd, "swap", eSwapCoords, list);
959 static int do_cpt_footer(XDR *xd, int file_version)
964 if (file_version >= 2)
967 res = xdr_int(xd, &magic);
972 if (magic != CPT_MAGIC2)
981 static int do_cpt_state(XDR *xd, gmx_bool bRead,
982 int fflags, t_state *state,
992 nnht = state->nhchainlength*state->ngtc;
993 nnhtp = state->nhchainlength*state->nnhpres;
995 if (bRead) /* we need to allocate space for dfhist if we are reading */
997 init_df_history(&state->dfhist, state->dfhist.nlambda);
1000 sflags = state->flags;
1001 for (i = 0; (i < estNR && ret == 0); i++)
1003 if (fflags & (1<<i))
1007 case estLAMBDA: ret = do_cpte_reals(xd, cptpEST, i, sflags, efptNR, &(state->lambda), list); break;
1008 case estFEPSTATE: ret = do_cpte_int (xd, cptpEST, i, sflags, &state->fep_state, list); break;
1009 case estBOX: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box, list); break;
1010 case estBOX_REL: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box_rel, list); break;
1011 case estBOXV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->boxv, list); break;
1012 case estPRES_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->pres_prev, list); break;
1013 case estSVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->svir_prev, list); break;
1014 case estFVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->fvir_prev, list); break;
1015 case estNH_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_xi, list); break;
1016 case estNH_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_vxi, list); break;
1017 case estNHPRES_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_xi, list); break;
1018 case estNHPRES_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_vxi, list); break;
1019 case estTC_INT: ret = do_cpte_doubles(xd, cptpEST, i, sflags, state->ngtc, &state->therm_integral, list); break;
1020 case estVETA: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->veta, list); break;
1021 case estVOL0: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->vol0, list); break;
1022 case estX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->x, list); break;
1023 case estV: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->v, list); break;
1024 case estSDX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->sd_X, list); break;
1025 /* The RNG entries are no longer written,
1026 * the next 4 lines are only for reading old files.
1028 case estLD_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1029 case estLD_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1030 case estMC_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1031 case estMC_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1032 case estDISRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.disre_initf, list); break;
1033 case estDISRE_RM3TAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.ndisrepairs, &state->hist.disre_rm3tav, list); break;
1034 case estORIRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.orire_initf, list); break;
1035 case estORIRE_DTAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.norire_Dtav, &state->hist.orire_Dtav, list); break;
1037 gmx_fatal(FARGS, "Unknown state entry %d\n"
1038 "You are probably reading a new checkpoint file with old code", i);
1046 static int do_cpt_ekinstate(XDR *xd, int fflags, ekinstate_t *ekins,
1054 for (i = 0; (i < eeksNR && ret == 0); i++)
1056 if (fflags & (1<<i))
1061 case eeksEKIN_N: ret = do_cpte_int(xd, cptpEEKS, i, fflags, &ekins->ekin_n, list); break;
1062 case eeksEKINH: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh, list); break;
1063 case eeksEKINF: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinf, list); break;
1064 case eeksEKINO: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh_old, list); break;
1065 case eeksEKINTOTAL: ret = do_cpte_matrix(xd, cptpEEKS, i, fflags, ekins->ekin_total, list); break;
1066 case eeksEKINSCALEF: ret = do_cpte_doubles(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinscalef_nhc, list); break;
1067 case eeksVSCALE: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->vscale_nhc, list); break;
1068 case eeksEKINSCALEH: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->ekinscaleh_nhc, list); break;
1069 case eeksDEKINDL: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->dekindl, list); break;
1070 case eeksMVCOS: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->mvcos, list); break;
1072 gmx_fatal(FARGS, "Unknown ekin data state entry %d\n"
1073 "You are probably reading a new checkpoint file with old code", i);
1082 static int do_cpt_swapstate(XDR *xd, gmx_bool bRead, swapstate_t *swapstate, FILE *list)
1086 int swap_cpt_version = 1;
1089 if (eswapNO == swapstate->eSwapCoords)
1094 swapstate->bFromCpt = bRead;
1096 do_cpt_int_err(xd, "swap checkpoint version", &swap_cpt_version, list);
1097 do_cpt_int_err(xd, "swap coupling steps", &swapstate->nAverage, list);
1099 /* When reading, init_swapcoords has not been called yet,
1100 * so we have to allocate memory first. */
1102 for (ic = 0; ic < eCompNR; ic++)
1104 for (ii = 0; ii < eIonNR; ii++)
1108 do_cpt_int_err(xd, "swap requested atoms", &swapstate->nat_req[ic][ii], list);
1112 do_cpt_int_err(xd, "swap requested atoms p", swapstate->nat_req_p[ic][ii], list);
1117 do_cpt_int_err(xd, "swap influx netto", &swapstate->inflow_netto[ic][ii], list);
1121 do_cpt_int_err(xd, "swap influx netto p", swapstate->inflow_netto_p[ic][ii], list);
1124 if (bRead && (NULL == swapstate->nat_past[ic][ii]) )
1126 snew(swapstate->nat_past[ic][ii], swapstate->nAverage);
1129 for (j = 0; j < swapstate->nAverage; j++)
1133 do_cpt_int_err(xd, "swap past atom counts", &swapstate->nat_past[ic][ii][j], list);
1137 do_cpt_int_err(xd, "swap past atom counts p", &swapstate->nat_past_p[ic][ii][j], list);
1143 /* Ion flux per channel */
1144 for (ic = 0; ic < eChanNR; ic++)
1146 for (ii = 0; ii < eIonNR; ii++)
1150 do_cpt_int_err(xd, "channel flux", &swapstate->fluxfromAtoB[ic][ii], list);
1154 do_cpt_int_err(xd, "channel flux p", swapstate->fluxfromAtoB_p[ic][ii], list);
1159 /* Ion flux leakage */
1162 snew(swapstate->fluxleak, 1);
1164 do_cpt_int_err(xd, "flux leakage", swapstate->fluxleak, list);
1167 do_cpt_int_err(xd, "number of ions", &swapstate->nions, list);
1171 snew(swapstate->channel_label, swapstate->nions);
1172 snew(swapstate->comp_from, swapstate->nions);
1175 do_cpt_u_chars(xd, "channel history", swapstate->nions, swapstate->channel_label, list);
1176 do_cpt_u_chars(xd, "domain history", swapstate->nions, swapstate->comp_from, list);
1178 /* Save the last known whole positions to checkpoint
1179 * file to be able to also make multimeric channels whole in PBC */
1180 do_cpt_int_err(xd, "Ch0 atoms", &swapstate->nat[eChan0], list);
1181 do_cpt_int_err(xd, "Ch1 atoms", &swapstate->nat[eChan1], list);
1184 snew(swapstate->xc_old_whole[eChan0], swapstate->nat[eChan0]);
1185 snew(swapstate->xc_old_whole[eChan1], swapstate->nat[eChan1]);
1186 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], swapstate->xc_old_whole[eChan0], list);
1187 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], swapstate->xc_old_whole[eChan1], list);
1191 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], *swapstate->xc_old_whole_p[eChan0], list);
1192 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], *swapstate->xc_old_whole_p[eChan1], list);
1199 static int do_cpt_enerhist(XDR *xd, gmx_bool bRead,
1200 int fflags, energyhistory_t *enerhist,
1211 enerhist->nsteps = 0;
1213 enerhist->nsteps_sim = 0;
1214 enerhist->nsum_sim = 0;
1215 enerhist->dht = NULL;
1217 if (fflags & (1<< eenhENERGY_DELTA_H_NN) )
1219 snew(enerhist->dht, 1);
1220 enerhist->dht->ndh = NULL;
1221 enerhist->dht->dh = NULL;
1222 enerhist->dht->start_lambda_set = FALSE;
1226 for (i = 0; (i < eenhNR && ret == 0); i++)
1228 if (fflags & (1<<i))
1232 case eenhENERGY_N: ret = do_cpte_int(xd, cptpEENH, i, fflags, &enerhist->nener, list); break;
1233 case eenhENERGY_AVER: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_ave, list); break;
1234 case eenhENERGY_SUM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum, list); break;
1235 case eenhENERGY_NSUM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum, list); break;
1236 case eenhENERGY_SUM_SIM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum_sim, list); break;
1237 case eenhENERGY_NSUM_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum_sim, list); break;
1238 case eenhENERGY_NSTEPS: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps, list); break;
1239 case eenhENERGY_NSTEPS_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps_sim, list); break;
1240 case eenhENERGY_DELTA_H_NN: do_cpt_int_err(xd, eenh_names[i], &(enerhist->dht->nndh), list);
1241 if (bRead) /* now allocate memory for it */
1243 snew(enerhist->dht->dh, enerhist->dht->nndh);
1244 snew(enerhist->dht->ndh, enerhist->dht->nndh);
1245 for (j = 0; j < enerhist->dht->nndh; j++)
1247 enerhist->dht->ndh[j] = 0;
1248 enerhist->dht->dh[j] = NULL;
1252 case eenhENERGY_DELTA_H_LIST:
1253 for (j = 0; j < enerhist->dht->nndh; j++)
1255 ret = do_cpte_n_reals(xd, cptpEENH, i, fflags, &enerhist->dht->ndh[j], &(enerhist->dht->dh[j]), list);
1258 case eenhENERGY_DELTA_H_STARTTIME:
1259 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_time), list); break;
1260 case eenhENERGY_DELTA_H_STARTLAMBDA:
1261 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_lambda), list); break;
1263 gmx_fatal(FARGS, "Unknown energy history entry %d\n"
1264 "You are probably reading a new checkpoint file with old code", i);
1269 if ((fflags & (1<<eenhENERGY_SUM)) && !(fflags & (1<<eenhENERGY_SUM_SIM)))
1271 /* Assume we have an old file format and copy sum to sum_sim */
1272 srenew(enerhist->ener_sum_sim, enerhist->nener);
1273 for (i = 0; i < enerhist->nener; i++)
1275 enerhist->ener_sum_sim[i] = enerhist->ener_sum[i];
1279 if ( (fflags & (1<<eenhENERGY_NSUM)) &&
1280 !(fflags & (1<<eenhENERGY_NSTEPS)))
1282 /* Assume we have an old file format and copy nsum to nsteps */
1283 enerhist->nsteps = enerhist->nsum;
1285 if ( (fflags & (1<<eenhENERGY_NSUM_SIM)) &&
1286 !(fflags & (1<<eenhENERGY_NSTEPS_SIM)))
1288 /* Assume we have an old file format and copy nsum to nsteps */
1289 enerhist->nsteps_sim = enerhist->nsum_sim;
1295 static int do_cpt_df_hist(XDR *xd, int fflags, df_history_t *dfhist, FILE *list)
1300 nlambda = dfhist->nlambda;
1303 for (i = 0; (i < edfhNR && ret == 0); i++)
1305 if (fflags & (1<<i))
1309 case edfhBEQUIL: ret = do_cpte_int(xd, cptpEDFH, i, fflags, &dfhist->bEquil, list); break;
1310 case edfhNATLAMBDA: ret = do_cpte_ints(xd, cptpEDFH, i, fflags, nlambda, &dfhist->n_at_lam, list); break;
1311 case edfhWLHISTO: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->wl_histo, list); break;
1312 case edfhWLDELTA: ret = do_cpte_real(xd, cptpEDFH, i, fflags, &dfhist->wl_delta, list); break;
1313 case edfhSUMWEIGHTS: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_weights, list); break;
1314 case edfhSUMDG: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_dg, list); break;
1315 case edfhSUMMINVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_minvar, list); break;
1316 case edfhSUMVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_variance, list); break;
1317 case edfhACCUMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p, list); break;
1318 case edfhACCUMM: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m, list); break;
1319 case edfhACCUMP2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p2, list); break;
1320 case edfhACCUMM2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m2, list); break;
1321 case edfhTIJ: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij, list); break;
1322 case edfhTIJEMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij_empirical, list); break;
1325 gmx_fatal(FARGS, "Unknown df history entry %d\n"
1326 "You are probably reading a new checkpoint file with old code", i);
1335 /* This function stores the last whole configuration of the reference and
1336 * average structure in the .cpt file
1338 static int do_cpt_EDstate(XDR *xd, gmx_bool bRead,
1339 edsamstate_t *EDstate, FILE *list)
1346 EDstate->bFromCpt = bRead;
1348 if (EDstate->nED <= 0)
1353 /* When reading, init_edsam has not been called yet,
1354 * so we have to allocate memory first. */
1357 snew(EDstate->nref, EDstate->nED);
1358 snew(EDstate->old_sref, EDstate->nED);
1359 snew(EDstate->nav, EDstate->nED);
1360 snew(EDstate->old_sav, EDstate->nED);
1363 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1364 for (i = 0; i < EDstate->nED; i++)
1366 /* Reference structure SREF */
1367 sprintf(buf, "ED%d # of atoms in reference structure", i+1);
1368 do_cpt_int_err(xd, buf, &EDstate->nref[i], list);
1369 sprintf(buf, "ED%d x_ref", i+1);
1372 snew(EDstate->old_sref[i], EDstate->nref[i]);
1373 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref[i], list);
1377 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref_p[i], list);
1380 /* Average structure SAV */
1381 sprintf(buf, "ED%d # of atoms in average structure", i+1);
1382 do_cpt_int_err(xd, buf, &EDstate->nav[i], list);
1383 sprintf(buf, "ED%d x_av", i+1);
1386 snew(EDstate->old_sav[i], EDstate->nav[i]);
1387 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav[i], list);
1391 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav_p[i], list);
1399 static int do_cpt_files(XDR *xd, gmx_bool bRead,
1400 gmx_file_position_t **p_outputfiles, int *nfiles,
1401 FILE *list, int file_version)
1405 gmx_off_t mask = 0xFFFFFFFFL;
1406 int offset_high, offset_low;
1408 gmx_file_position_t *outputfiles;
1410 if (do_cpt_int(xd, "number of output files", nfiles, list) != 0)
1417 snew(*p_outputfiles, *nfiles);
1420 outputfiles = *p_outputfiles;
1422 for (i = 0; i < *nfiles; i++)
1424 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1427 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1428 strncpy(outputfiles[i].filename, buf, CPTSTRLEN-1);
1434 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1438 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1442 outputfiles[i].offset = ( ((gmx_off_t) offset_high) << 32 ) | ( (gmx_off_t) offset_low & mask );
1446 buf = outputfiles[i].filename;
1447 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1449 offset = outputfiles[i].offset;
1457 offset_low = (int) (offset & mask);
1458 offset_high = (int) ((offset >> 32) & mask);
1460 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1464 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1469 if (file_version >= 8)
1471 if (do_cpt_int(xd, "file_checksum_size", &(outputfiles[i].chksum_size),
1476 if (do_cpt_u_chars(xd, "file_checksum", 16, outputfiles[i].chksum, list) != 0)
1483 outputfiles[i].chksum_size = -1;
1490 void write_checkpoint(const char *fn, gmx_bool bNumberAndKeep,
1491 FILE *fplog, t_commrec *cr,
1492 int eIntegrator, int simulation_part,
1493 gmx_bool bExpanded, int elamstats,
1494 gmx_int64_t step, double t, t_state *state)
1504 char *fntemp; /* the temporary checkpoint file name */
1506 char timebuf[STRLEN];
1507 int nppnodes, npmenodes;
1508 char buf[1024], suffix[5+STEPSTRSIZE], sbuf[STEPSTRSIZE];
1509 gmx_file_position_t *outputfiles;
1512 int flags_eks, flags_enh, flags_dfh;
1515 if (DOMAINDECOMP(cr))
1517 nppnodes = cr->dd->nnodes;
1518 npmenodes = cr->npmenodes;
1526 #ifndef GMX_NO_RENAME
1527 /* make the new temporary filename */
1528 snew(fntemp, strlen(fn)+5+STEPSTRSIZE);
1530 fntemp[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1531 sprintf(suffix, "_%s%s", "step", gmx_step_str(step, sbuf));
1532 strcat(fntemp, suffix);
1533 strcat(fntemp, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1535 /* if we can't rename, we just overwrite the cpt file.
1536 * dangerous if interrupted.
1538 snew(fntemp, strlen(fn));
1542 gmx_ctime_r(&now, timebuf, STRLEN);
1546 fprintf(fplog, "Writing checkpoint, step %s at %s\n\n",
1547 gmx_step_str(step, buf), timebuf);
1550 /* Get offsets for open files */
1551 gmx_fio_get_output_file_positions(&outputfiles, &noutputfiles);
1553 fp = gmx_fio_open(fntemp, "w");
1555 if (state->ekinstate.bUpToDate)
1558 ((1<<eeksEKIN_N) | (1<<eeksEKINH) | (1<<eeksEKINF) |
1559 (1<<eeksEKINO) | (1<<eeksEKINSCALEF) | (1<<eeksEKINSCALEH) |
1560 (1<<eeksVSCALE) | (1<<eeksDEKINDL) | (1<<eeksMVCOS));
1568 if (state->enerhist.nsum > 0 || state->enerhist.nsum_sim > 0)
1570 flags_enh |= (1<<eenhENERGY_N);
1571 if (state->enerhist.nsum > 0)
1573 flags_enh |= ((1<<eenhENERGY_AVER) | (1<<eenhENERGY_SUM) |
1574 (1<<eenhENERGY_NSTEPS) | (1<<eenhENERGY_NSUM));
1576 if (state->enerhist.nsum_sim > 0)
1578 flags_enh |= ((1<<eenhENERGY_SUM_SIM) | (1<<eenhENERGY_NSTEPS_SIM) |
1579 (1<<eenhENERGY_NSUM_SIM));
1581 if (state->enerhist.dht)
1583 flags_enh |= ( (1<< eenhENERGY_DELTA_H_NN) |
1584 (1<< eenhENERGY_DELTA_H_LIST) |
1585 (1<< eenhENERGY_DELTA_H_STARTTIME) |
1586 (1<< eenhENERGY_DELTA_H_STARTLAMBDA) );
1592 flags_dfh = ((1<<edfhBEQUIL) | (1<<edfhNATLAMBDA) | (1<<edfhSUMWEIGHTS) | (1<<edfhSUMDG) |
1593 (1<<edfhTIJ) | (1<<edfhTIJEMP));
1596 flags_dfh |= ((1<<edfhWLDELTA) | (1<<edfhWLHISTO));
1598 if ((elamstats == elamstatsMINVAR) || (elamstats == elamstatsBARKER) || (elamstats == elamstatsMETROPOLIS))
1600 flags_dfh |= ((1<<edfhACCUMP) | (1<<edfhACCUMM) | (1<<edfhACCUMP2) | (1<<edfhACCUMM2)
1601 | (1<<edfhSUMMINVAR) | (1<<edfhSUMVAR));
1609 /* We can check many more things now (CPU, acceleration, etc), but
1610 * it is highly unlikely to have two separate builds with exactly
1611 * the same version, user, time, and build host!
1614 version = gmx_strdup(gmx_version());
1615 btime = gmx_strdup(BUILD_TIME);
1616 buser = gmx_strdup(BUILD_USER);
1617 bhost = gmx_strdup(BUILD_HOST);
1619 double_prec = GMX_CPT_BUILD_DP;
1620 fprog = gmx_strdup(Program());
1622 ftime = &(timebuf[0]);
1624 do_cpt_header(gmx_fio_getxdr(fp), FALSE, &file_version,
1625 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1626 &eIntegrator, &simulation_part, &step, &t, &nppnodes,
1627 DOMAINDECOMP(cr) ? cr->dd->nc : NULL, &npmenodes,
1628 &state->natoms, &state->ngtc, &state->nnhpres,
1629 &state->nhchainlength, &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
1630 &state->edsamstate.nED, &state->swapstate.eSwapCoords,
1639 if ((do_cpt_state(gmx_fio_getxdr(fp), FALSE, state->flags, state, NULL) < 0) ||
1640 (do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL) < 0) ||
1641 (do_cpt_enerhist(gmx_fio_getxdr(fp), FALSE, flags_enh, &state->enerhist, NULL) < 0) ||
1642 (do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL) < 0) ||
1643 (do_cpt_EDstate(gmx_fio_getxdr(fp), FALSE, &state->edsamstate, NULL) < 0) ||
1644 (do_cpt_swapstate(gmx_fio_getxdr(fp), FALSE, &state->swapstate, NULL) < 0) ||
1645 (do_cpt_files(gmx_fio_getxdr(fp), FALSE, &outputfiles, &noutputfiles, NULL,
1648 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1651 do_cpt_footer(gmx_fio_getxdr(fp), file_version);
1653 /* we really, REALLY, want to make sure to physically write the checkpoint,
1654 and all the files it depends on, out to disk. Because we've
1655 opened the checkpoint with gmx_fio_open(), it's in our list
1657 ret = gmx_fio_all_output_fsync();
1663 "Cannot fsync '%s'; maybe you are out of disk space?",
1664 gmx_fio_getname(ret));
1666 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV) == NULL)
1676 if (gmx_fio_close(fp) != 0)
1678 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1681 /* we don't move the checkpoint if the user specified they didn't want it,
1682 or if the fsyncs failed */
1683 #ifndef GMX_NO_RENAME
1684 if (!bNumberAndKeep && !ret)
1688 /* Rename the previous checkpoint file */
1690 buf[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1691 strcat(buf, "_prev");
1692 strcat(buf, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1694 /* we copy here so that if something goes wrong between now and
1695 * the rename below, there's always a state.cpt.
1696 * If renames are atomic (such as in POSIX systems),
1697 * this copying should be unneccesary.
1699 gmx_file_copy(fn, buf, FALSE);
1700 /* We don't really care if this fails:
1701 * there's already a new checkpoint.
1704 gmx_file_rename(fn, buf);
1707 if (gmx_file_rename(fntemp, fn) != 0)
1709 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
1712 #endif /* GMX_NO_RENAME */
1718 /*code for alternate checkpointing scheme. moved from top of loop over
1720 fcRequestCheckPoint();
1721 if (fcCheckPointParallel( cr->nodeid, NULL, 0) == 0)
1723 gmx_fatal( 3, __FILE__, __LINE__, "Checkpoint error on step %d\n", step );
1725 #endif /* end GMX_FAHCORE block */
1728 static void print_flag_mismatch(FILE *fplog, int sflags, int fflags)
1732 fprintf(fplog, "\nState entry mismatch between the simulation and the checkpoint file\n");
1733 fprintf(fplog, "Entries which are not present in the checkpoint file will not be updated\n");
1734 fprintf(fplog, " %24s %11s %11s\n", "", "simulation", "checkpoint");
1735 for (i = 0; i < estNR; i++)
1737 if ((sflags & (1<<i)) || (fflags & (1<<i)))
1739 fprintf(fplog, " %24s %11s %11s\n",
1741 (sflags & (1<<i)) ? " present " : "not present",
1742 (fflags & (1<<i)) ? " present " : "not present");
1747 static void check_int(FILE *fplog, const char *type, int p, int f, gmx_bool *mm)
1749 FILE *fp = fplog ? fplog : stderr;
1753 fprintf(fp, " %s mismatch,\n", type);
1754 fprintf(fp, " current program: %d\n", p);
1755 fprintf(fp, " checkpoint file: %d\n", f);
1761 static void check_string(FILE *fplog, const char *type, const char *p,
1762 const char *f, gmx_bool *mm)
1764 FILE *fp = fplog ? fplog : stderr;
1766 if (strcmp(p, f) != 0)
1768 fprintf(fp, " %s mismatch,\n", type);
1769 fprintf(fp, " current program: %s\n", p);
1770 fprintf(fp, " checkpoint file: %s\n", f);
1776 static void check_match(FILE *fplog,
1778 char *btime, char *buser, char *bhost, int double_prec,
1780 t_commrec *cr, int npp_f, int npme_f,
1781 ivec dd_nc, ivec dd_nc_f)
1784 gmx_bool mm = FALSE;
1785 gmx_bool patchlevel_differs = FALSE;
1786 gmx_bool version_differs = FALSE;
1788 check_string(fplog, "Version", gmx_version(), version, &mm);
1789 patchlevel_differs = mm;
1791 if (patchlevel_differs)
1793 /* Gromacs should be able to continue from checkpoints between
1794 * different patch level versions, but we do not guarantee
1795 * compatibility between different major/minor versions - check this.
1797 int gmx_major, gmx_minor;
1798 int cpt_major, cpt_minor;
1799 sscanf(gmx_version(), "VERSION %5d.%5d", &gmx_major, &gmx_minor);
1800 sscanf(version, "VERSION %5d.%5d", &cpt_major, &cpt_minor);
1801 version_differs = (gmx_major != cpt_major || gmx_minor != cpt_minor);
1804 check_string(fplog, "Build time", BUILD_TIME, btime, &mm);
1805 check_string(fplog, "Build user", BUILD_USER, buser, &mm);
1806 check_string(fplog, "Build host", BUILD_HOST, bhost, &mm);
1807 check_int (fplog, "Double prec.", GMX_CPT_BUILD_DP, double_prec, &mm);
1808 check_string(fplog, "Program name", Program(), fprog, &mm);
1810 check_int (fplog, "#ranks", cr->nnodes, npp_f+npme_f, &mm);
1813 check_int (fplog, "#PME-ranks", cr->npmenodes, npme_f, &mm);
1816 if (cr->npmenodes >= 0)
1818 npp -= cr->npmenodes;
1822 check_int (fplog, "#DD-cells[x]", dd_nc[XX], dd_nc_f[XX], &mm);
1823 check_int (fplog, "#DD-cells[y]", dd_nc[YY], dd_nc_f[YY], &mm);
1824 check_int (fplog, "#DD-cells[z]", dd_nc[ZZ], dd_nc_f[ZZ], &mm);
1830 const char msg_version_difference[] =
1831 "The current Gromacs major & minor version are not identical to those that\n"
1832 "generated the checkpoint file. In principle Gromacs does not support\n"
1833 "continuation from checkpoints between different versions, so we advise\n"
1834 "against this. If you still want to try your luck we recommend that you use\n"
1835 "the -noappend flag to keep your output files from the two versions separate.\n"
1836 "This might also work around errors where the output fields in the energy\n"
1837 "file have changed between the different major & minor versions.\n";
1839 const char msg_mismatch_notice[] =
1840 "Gromacs patchlevel, binary or parallel settings differ from previous run.\n"
1841 "Continuation is exact, but not guaranteed to be binary identical.\n";
1843 const char msg_logdetails[] =
1844 "See the log file for details.\n";
1846 if (version_differs)
1848 fprintf(stderr, "%s%s\n", msg_version_difference, fplog ? msg_logdetails : "");
1852 fprintf(fplog, "%s\n", msg_version_difference);
1857 /* Major & minor versions match at least, but something is different. */
1858 fprintf(stderr, "%s%s\n", msg_mismatch_notice, fplog ? msg_logdetails : "");
1861 fprintf(fplog, "%s\n", msg_mismatch_notice);
1867 static void read_checkpoint(const char *fn, FILE **pfplog,
1868 t_commrec *cr, ivec dd_nc,
1869 int eIntegrator, int *init_fep_state, gmx_int64_t *step, double *t,
1870 t_state *state, gmx_bool *bReadEkin,
1871 int *simulation_part,
1872 gmx_bool bAppendOutputFiles, gmx_bool bForceAppend)
1877 char *version, *btime, *buser, *bhost, *fprog, *ftime;
1879 char buf[STEPSTRSIZE];
1880 int eIntegrator_f, nppnodes_f, npmenodes_f;
1882 int natoms, ngtc, nnhpres, nhchainlength, nlambda, fflags, flags_eks, flags_enh, flags_dfh;
1885 gmx_file_position_t *outputfiles;
1887 t_fileio *chksum_file;
1888 FILE * fplog = *pfplog;
1889 unsigned char digest[16];
1890 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1891 struct flock fl; /* don't initialize here: the struct order is OS
1895 const char *int_warn =
1896 "WARNING: The checkpoint file was generated with integrator %s,\n"
1897 " while the simulation uses integrator %s\n\n";
1899 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1900 fl.l_type = F_WRLCK;
1901 fl.l_whence = SEEK_SET;
1907 fp = gmx_fio_open(fn, "r");
1908 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
1909 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1910 &eIntegrator_f, simulation_part, step, t,
1911 &nppnodes_f, dd_nc_f, &npmenodes_f,
1912 &natoms, &ngtc, &nnhpres, &nhchainlength, &nlambda,
1913 &fflags, &flags_eks, &flags_enh, &flags_dfh,
1914 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
1916 if (bAppendOutputFiles &&
1917 file_version >= 13 && double_prec != GMX_CPT_BUILD_DP)
1919 gmx_fatal(FARGS, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
1922 if (cr == NULL || MASTER(cr))
1924 fprintf(stderr, "\nReading checkpoint file %s generated: %s\n\n",
1928 /* This will not be written if we do appending, since fplog is still NULL then */
1931 fprintf(fplog, "\n");
1932 fprintf(fplog, "Reading checkpoint file %s\n", fn);
1933 fprintf(fplog, " file generated by: %s\n", fprog);
1934 fprintf(fplog, " file generated at: %s\n", ftime);
1935 fprintf(fplog, " GROMACS build time: %s\n", btime);
1936 fprintf(fplog, " GROMACS build user: %s\n", buser);
1937 fprintf(fplog, " GROMACS build host: %s\n", bhost);
1938 fprintf(fplog, " GROMACS double prec.: %d\n", double_prec);
1939 fprintf(fplog, " simulation part #: %d\n", *simulation_part);
1940 fprintf(fplog, " step: %s\n", gmx_step_str(*step, buf));
1941 fprintf(fplog, " time: %f\n", *t);
1942 fprintf(fplog, "\n");
1945 if (natoms != state->natoms)
1947 gmx_fatal(FARGS, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", natoms, state->natoms);
1949 if (ngtc != state->ngtc)
1951 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);
1953 if (nnhpres != state->nnhpres)
1955 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);
1958 if (nlambda != state->dfhist.nlambda)
1960 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);
1963 init_gtc_state(state, state->ngtc, state->nnhpres, nhchainlength); /* need to keep this here to keep the tpr format working */
1964 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
1966 if (eIntegrator_f != eIntegrator)
1970 fprintf(stderr, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1972 if (bAppendOutputFiles)
1975 "Output file appending requested, but input/checkpoint integrators do not match.\n"
1976 "Stopping the run to prevent you from ruining all your data...\n"
1977 "If you _really_ know what you are doing, try with the -noappend option.\n");
1981 fprintf(fplog, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1989 else if (cr->nnodes == nppnodes_f + npmenodes_f)
1991 if (cr->npmenodes < 0)
1993 cr->npmenodes = npmenodes_f;
1995 int nppnodes = cr->nnodes - cr->npmenodes;
1996 if (nppnodes == nppnodes_f)
1998 for (d = 0; d < DIM; d++)
2002 dd_nc[d] = dd_nc_f[d];
2008 if (fflags != state->flags)
2013 if (bAppendOutputFiles)
2016 "Output file appending requested, but input and checkpoint states are not identical.\n"
2017 "Stopping the run to prevent you from ruining all your data...\n"
2018 "You can try with the -noappend option, and get more info in the log file.\n");
2021 if (getenv("GMX_ALLOW_CPT_MISMATCH") == NULL)
2023 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");
2028 "WARNING: The checkpoint state entries do not match the simulation,\n"
2029 " see the log file for details\n\n");
2035 print_flag_mismatch(fplog, state->flags, fflags);
2042 check_match(fplog, version, btime, buser, bhost, double_prec, fprog,
2043 cr, nppnodes_f, npmenodes_f, dd_nc, dd_nc_f);
2046 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, fflags, state, NULL);
2047 *init_fep_state = state->fep_state; /* there should be a better way to do this than setting it here.
2048 Investigate for 5.0. */
2053 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2058 *bReadEkin = ((flags_eks & (1<<eeksEKINH)) || (flags_eks & (1<<eeksEKINF)) || (flags_eks & (1<<eeksEKINO)) ||
2059 ((flags_eks & (1<<eeksEKINSCALEF)) | (flags_eks & (1<<eeksEKINSCALEH)) | (flags_eks & (1<<eeksVSCALE))));
2061 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2062 flags_enh, &state->enerhist, NULL);
2068 if (file_version < 6)
2070 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.";
2072 fprintf(stderr, "\nWARNING: %s\n\n", warn);
2075 fprintf(fplog, "\nWARNING: %s\n\n", warn);
2077 state->enerhist.nsum = *step;
2078 state->enerhist.nsum_sim = *step;
2081 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2087 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2093 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2099 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, NULL, file_version);
2105 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2110 if (gmx_fio_close(fp) != 0)
2112 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2121 /* If the user wants to append to output files,
2122 * we use the file pointer positions of the output files stored
2123 * in the checkpoint file and truncate the files such that any frames
2124 * written after the checkpoint time are removed.
2125 * All files are md5sum checked such that we can be sure that
2126 * we do not truncate other (maybe imprortant) files.
2128 if (bAppendOutputFiles)
2130 if (fn2ftp(outputfiles[0].filename) != efLOG)
2132 /* make sure first file is log file so that it is OK to use it for
2135 gmx_fatal(FARGS, "The first output file should always be the log "
2136 "file but instead is: %s. Cannot do appending because of this condition.", outputfiles[0].filename);
2138 for (i = 0; i < nfiles; i++)
2140 if (outputfiles[i].offset < 0)
2142 gmx_fatal(FARGS, "The original run wrote a file called '%s' which "
2143 "is larger than 2 GB, but mdrun did not support large file"
2144 " offsets. Can not append. Run mdrun with -noappend",
2145 outputfiles[i].filename);
2148 chksum_file = gmx_fio_open(outputfiles[i].filename, "a");
2151 chksum_file = gmx_fio_open(outputfiles[i].filename, "r+");
2156 /* Note that there are systems where the lock operation
2157 * will succeed, but a second process can also lock the file.
2158 * We should probably try to detect this.
2160 #if defined __native_client__
2164 #elif defined GMX_NATIVE_WINDOWS
2165 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX) == -1)
2167 if (fcntl(fileno(gmx_fio_getfp(chksum_file)), F_SETLK, &fl) == -1)
2170 if (errno == ENOSYS)
2174 gmx_fatal(FARGS, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2178 fprintf(stderr, "\nNOTE: File locking is not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2181 fprintf(fplog, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2185 else if (errno == EACCES || errno == EAGAIN)
2187 gmx_fatal(FARGS, "Failed to lock: %s. Already running "
2188 "simulation?", outputfiles[i].filename);
2192 gmx_fatal(FARGS, "Failed to lock: %s. %s.",
2193 outputfiles[i].filename, strerror(errno));
2198 /* compute md5 chksum */
2199 if (outputfiles[i].chksum_size != -1)
2201 if (gmx_fio_get_file_md5(chksum_file, outputfiles[i].offset,
2202 digest) != outputfiles[i].chksum_size) /*at the end of the call the file position is at the end of the file*/
2204 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.",
2205 outputfiles[i].chksum_size,
2206 outputfiles[i].filename);
2209 if (i == 0) /*log file needs to be seeked in case we need to truncate (other files are truncated below)*/
2211 if (gmx_fio_seek(chksum_file, outputfiles[i].offset))
2213 gmx_fatal(FARGS, "Seek error! Failed to truncate log-file: %s.", strerror(errno));
2218 if (i == 0) /*open log file here - so that lock is never lifted
2219 after chksum is calculated */
2221 *pfplog = gmx_fio_getfp(chksum_file);
2225 gmx_fio_close(chksum_file);
2228 /* compare md5 chksum */
2229 if (outputfiles[i].chksum_size != -1 &&
2230 memcmp(digest, outputfiles[i].chksum, 16) != 0)
2234 fprintf(debug, "chksum for %s: ", outputfiles[i].filename);
2235 for (j = 0; j < 16; j++)
2237 fprintf(debug, "%02x", digest[j]);
2239 fprintf(debug, "\n");
2241 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.",
2242 outputfiles[i].filename);
2247 if (i != 0) /*log file is already seeked to correct position */
2249 #ifdef GMX_NATIVE_WINDOWS
2250 rc = gmx_wintruncate(outputfiles[i].filename, outputfiles[i].offset);
2252 rc = truncate(outputfiles[i].filename, outputfiles[i].offset);
2256 gmx_fatal(FARGS, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfiles[i].filename);
2266 void load_checkpoint(const char *fn, FILE **fplog,
2267 t_commrec *cr, ivec dd_nc,
2268 t_inputrec *ir, t_state *state,
2269 gmx_bool *bReadEkin,
2270 gmx_bool bAppend, gmx_bool bForceAppend)
2277 /* Read the state from the checkpoint file */
2278 read_checkpoint(fn, fplog,
2280 ir->eI, &(ir->fepvals->init_fep_state), &step, &t, state, bReadEkin,
2281 &ir->simulation_part, bAppend, bForceAppend);
2285 gmx_bcast(sizeof(cr->npmenodes), &cr->npmenodes, cr);
2286 gmx_bcast(DIM*sizeof(dd_nc[0]), dd_nc, cr);
2287 gmx_bcast(sizeof(step), &step, cr);
2288 gmx_bcast(sizeof(*bReadEkin), bReadEkin, cr);
2290 ir->bContinuation = TRUE;
2291 if (ir->nsteps >= 0)
2293 ir->nsteps += ir->init_step - step;
2295 ir->init_step = step;
2296 ir->simulation_part += 1;
2299 static void read_checkpoint_data(t_fileio *fp, int *simulation_part,
2300 gmx_int64_t *step, double *t, t_state *state,
2301 int *nfiles, gmx_file_position_t **outputfiles)
2304 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2309 int flags_eks, flags_enh, flags_dfh;
2311 gmx_file_position_t *files_loc = NULL;
2314 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2315 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2316 &eIntegrator, simulation_part, step, t, &nppnodes, dd_nc, &npme,
2317 &state->natoms, &state->ngtc, &state->nnhpres, &state->nhchainlength,
2318 &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
2319 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
2321 do_cpt_state(gmx_fio_getxdr(fp), TRUE, state->flags, state, NULL);
2326 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2331 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2332 flags_enh, &state->enerhist, NULL);
2337 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2343 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2349 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2355 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE,
2356 outputfiles != NULL ? outputfiles : &files_loc,
2357 outputfiles != NULL ? nfiles : &nfiles_loc,
2358 NULL, file_version);
2359 if (files_loc != NULL)
2369 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2383 read_checkpoint_state(const char *fn, int *simulation_part,
2384 gmx_int64_t *step, double *t, t_state *state)
2388 fp = gmx_fio_open(fn, "r");
2389 read_checkpoint_data(fp, simulation_part, step, t, state, NULL, NULL);
2390 if (gmx_fio_close(fp) != 0)
2392 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2396 void read_checkpoint_trxframe(t_fileio *fp, t_trxframe *fr)
2398 /* This next line is nasty because the sub-structures of t_state
2399 * cannot be assumed to be zeroed (or even initialized in ways the
2400 * rest of the code might assume). Using snew would be better, but
2401 * this will all go away for 5.0. */
2403 int simulation_part;
2407 init_state(&state, 0, 0, 0, 0, 0);
2409 read_checkpoint_data(fp, &simulation_part, &step, &t, &state, NULL, NULL);
2411 fr->natoms = state.natoms;
2414 fr->step = gmx_int64_to_int(step,
2415 "conversion of checkpoint to trajectory");
2419 fr->lambda = state.lambda[efptFEP];
2420 fr->fep_state = state.fep_state;
2422 fr->bX = (state.flags & (1<<estX));
2428 fr->bV = (state.flags & (1<<estV));
2435 fr->bBox = (state.flags & (1<<estBOX));
2438 copy_mat(state.box, fr->box);
2443 void list_checkpoint(const char *fn, FILE *out)
2447 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2449 int eIntegrator, simulation_part, nppnodes, npme;
2454 int flags_eks, flags_enh, flags_dfh;
2456 gmx_file_position_t *outputfiles;
2459 init_state(&state, -1, -1, -1, -1, 0);
2461 fp = gmx_fio_open(fn, "r");
2462 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2463 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2464 &eIntegrator, &simulation_part, &step, &t, &nppnodes, dd_nc, &npme,
2465 &state.natoms, &state.ngtc, &state.nnhpres, &state.nhchainlength,
2466 &(state.dfhist.nlambda), &state.flags,
2467 &flags_eks, &flags_enh, &flags_dfh, &state.edsamstate.nED,
2468 &state.swapstate.eSwapCoords, out);
2469 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, state.flags, &state, out);
2474 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state.ekinstate, out);
2479 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2480 flags_enh, &state.enerhist, out);
2484 ret = do_cpt_df_hist(gmx_fio_getxdr(fp),
2485 flags_dfh, &state.dfhist, out);
2490 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state.edsamstate, out);
2495 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state.swapstate, out);
2500 do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, out, file_version);
2505 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2512 if (gmx_fio_close(fp) != 0)
2514 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2521 static gmx_bool exist_output_file(const char *fnm_cp, int nfile, const t_filenm fnm[])
2525 /* Check if the output file name stored in the checkpoint file
2526 * is one of the output file names of mdrun.
2530 !(is_output(&fnm[i]) && strcmp(fnm_cp, fnm[i].fns[0]) == 0))
2535 return (i < nfile && gmx_fexist(fnm_cp));
2538 /* This routine cannot print tons of data, since it is called before the log file is opened. */
2539 gmx_bool read_checkpoint_simulation_part(const char *filename, int *simulation_part,
2540 gmx_int64_t *cpt_step, t_commrec *cr,
2541 gmx_bool bAppendReq,
2542 int nfile, const t_filenm fnm[],
2543 const char *part_suffix, gmx_bool *bAddPart)
2546 gmx_int64_t step = 0;
2548 /* This next line is nasty because the sub-structures of t_state
2549 * cannot be assumed to be zeroed (or even initialized in ways the
2550 * rest of the code might assume). Using snew would be better, but
2551 * this will all go away for 5.0. */
2554 gmx_file_position_t *outputfiles;
2557 char *fn, suf_up[STRLEN];
2563 if (!gmx_fexist(filename) || (!(fp = gmx_fio_open(filename, "r")) ))
2565 *simulation_part = 0;
2569 init_state(&state, 0, 0, 0, 0, 0);
2571 read_checkpoint_data(fp, simulation_part, &step, &t, &state,
2572 &nfiles, &outputfiles);
2573 if (gmx_fio_close(fp) != 0)
2575 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2582 for (f = 0; f < nfiles; f++)
2584 if (exist_output_file(outputfiles[f].filename, nfile, fnm))
2589 if (nexist == nfiles)
2591 bAppend = bAppendReq;
2593 else if (nexist > 0)
2596 "Output file appending has been requested,\n"
2597 "but some output files listed in the checkpoint file %s\n"
2598 "are not present or are named differently by the current program:\n",
2600 fprintf(stderr, "output files present:");
2601 for (f = 0; f < nfiles; f++)
2603 if (exist_output_file(outputfiles[f].filename,
2606 fprintf(stderr, " %s", outputfiles[f].filename);
2609 fprintf(stderr, "\n");
2610 fprintf(stderr, "output files not present or named differently:");
2611 for (f = 0; f < nfiles; f++)
2613 if (!exist_output_file(outputfiles[f].filename,
2616 fprintf(stderr, " %s", outputfiles[f].filename);
2619 fprintf(stderr, "\n");
2621 gmx_fatal(FARGS, "File appending requested, but %d of the %d output files are not present or are named differently", nfiles-nexist, nfiles);
2629 gmx_fatal(FARGS, "File appending requested, but no output file information is stored in the checkpoint file");
2631 fn = outputfiles[0].filename;
2632 if (strlen(fn) < 4 ||
2633 gmx_strcasecmp(fn+strlen(fn)-4, ftp2ext(efLOG)) == 0)
2635 gmx_fatal(FARGS, "File appending requested, but the log file is not the first file listed in the checkpoint file");
2637 /* Set bAddPart to whether the suffix string '.part' is present
2638 * in the log file name.
2640 strcpy(suf_up, part_suffix);
2642 *bAddPart = (strstr(fn, part_suffix) != NULL ||
2643 strstr(fn, suf_up) != NULL);
2651 gmx_bcast(sizeof(*simulation_part), simulation_part, cr);
2653 if (*simulation_part > 0 && bAppendReq)
2655 gmx_bcast(sizeof(bAppend), &bAppend, cr);
2656 gmx_bcast(sizeof(*bAddPart), bAddPart, cr);
2659 if (NULL != cpt_step)