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"
50 #ifdef HAVE_SYS_TIME_H
58 #ifdef GMX_NATIVE_WINDOWS
61 #include <sys/locking.h>
64 #include "gromacs/legacyheaders/copyrite.h"
65 #include "gromacs/legacyheaders/names.h"
66 #include "gromacs/legacyheaders/typedefs.h"
67 #include "gromacs/legacyheaders/types/commrec.h"
68 #include "gromacs/legacyheaders/txtdump.h"
69 #include "gromacs/math/vec.h"
70 #include "gromacs/legacyheaders/network.h"
72 #include "gromacs/fileio/filenm.h"
73 #include "gromacs/utility/futil.h"
74 #include "gromacs/fileio/gmxfio.h"
75 #include "gromacs/fileio/xdrf.h"
76 #include "gromacs/fileio/xdr_datatype.h"
77 #include "gromacs/utility/basenetwork.h"
78 #include "gromacs/utility/baseversion.h"
79 #include "gromacs/utility/cstringutil.h"
80 #include "gromacs/utility/fatalerror.h"
81 #include "gromacs/utility/smalloc.h"
83 #include "buildinfo.h"
89 #define CPT_MAGIC1 171817
90 #define CPT_MAGIC2 171819
91 #define CPTSTRLEN 1024
94 #define GMX_CPT_BUILD_DP 1
96 #define GMX_CPT_BUILD_DP 0
99 /* cpt_version should normally only be changed
100 * when the header of footer format changes.
101 * The state data format itself is backward and forward compatible.
102 * But old code can not read a new entry that is present in the file
103 * (but can read a new format when new entries are not present).
105 static const int cpt_version = 16;
108 const char *est_names[estNR] =
111 "box", "box-rel", "box-v", "pres_prev",
112 "nosehoover-xi", "thermostat-integral",
113 "x", "v", "SDx", "CGp", "LD-rng", "LD-rng-i",
114 "disre_initf", "disre_rm3tav",
115 "orire_initf", "orire_Dtav",
116 "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng", "MC-rng-i"
120 eeksEKIN_N, eeksEKINH, eeksDEKINDL, eeksMVCOS, eeksEKINF, eeksEKINO, eeksEKINSCALEF, eeksEKINSCALEH, eeksVSCALE, eeksEKINTOTAL, eeksNR
123 const char *eeks_names[eeksNR] =
125 "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
126 "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
130 eenhENERGY_N, eenhENERGY_AVER, eenhENERGY_SUM, eenhENERGY_NSUM,
131 eenhENERGY_SUM_SIM, eenhENERGY_NSUM_SIM,
132 eenhENERGY_NSTEPS, eenhENERGY_NSTEPS_SIM,
133 eenhENERGY_DELTA_H_NN,
134 eenhENERGY_DELTA_H_LIST,
135 eenhENERGY_DELTA_H_STARTTIME,
136 eenhENERGY_DELTA_H_STARTLAMBDA,
140 const char *eenh_names[eenhNR] =
142 "energy_n", "energy_aver", "energy_sum", "energy_nsum",
143 "energy_sum_sim", "energy_nsum_sim",
144 "energy_nsteps", "energy_nsteps_sim",
146 "energy_delta_h_list",
147 "energy_delta_h_start_time",
148 "energy_delta_h_start_lambda"
151 /* free energy history variables -- need to be preserved over checkpoint */
153 edfhBEQUIL, edfhNATLAMBDA, edfhWLHISTO, edfhWLDELTA, edfhSUMWEIGHTS, edfhSUMDG, edfhSUMMINVAR, edfhSUMVAR,
154 edfhACCUMP, edfhACCUMM, edfhACCUMP2, edfhACCUMM2, edfhTIJ, edfhTIJEMP, edfhNR
156 /* free energy history variable names */
157 const char *edfh_names[edfhNR] =
159 "bEquilibrated", "N_at_state", "Wang-Landau Histogram", "Wang-Landau Delta", "Weights", "Free Energies", "minvar", "variance",
160 "accumulated_plus", "accumulated_minus", "accumulated_plus_2", "accumulated_minus_2", "Tij", "Tij_empirical"
163 #ifdef GMX_NATIVE_WINDOWS
165 gmx_wintruncate(const char *filename, __int64 size)
168 /*we do this elsewhere*/
173 fp = fopen(filename, "rb+");
181 return _chsize_s( fileno(fp), size);
183 return _chsize( fileno(fp), size);
191 ecprREAL, ecprRVEC, ecprMATRIX
195 cptpEST, cptpEEKS, cptpEENH, cptpEDFH
197 /* enums for the different components of checkpoint variables, replacing the hard coded ones.
198 cptpEST - state variables.
199 cptpEEKS - Kinetic energy state variables.
200 cptpEENH - Energy history state variables.
201 cptpEDFH - free energy history variables.
205 static const char *st_names(int cptp, int ecpt)
209 case cptpEST: return est_names [ecpt]; break;
210 case cptpEEKS: return eeks_names[ecpt]; break;
211 case cptpEENH: return eenh_names[ecpt]; break;
212 case cptpEDFH: return edfh_names[ecpt]; break;
218 static void cp_warning(FILE *fp)
220 fprintf(fp, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
223 static void cp_error()
225 gmx_fatal(FARGS, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
228 static void do_cpt_string_err(XDR *xd, gmx_bool bRead, const char *desc, char **s, FILE *list)
236 res = xdr_string(xd, s, CPTSTRLEN);
243 fprintf(list, "%s = %s\n", desc, *s);
248 static int do_cpt_int(XDR *xd, const char *desc, int *i, FILE *list)
252 res = xdr_int(xd, i);
259 fprintf(list, "%s = %d\n", desc, *i);
264 static int do_cpt_u_chars(XDR *xd, const char *desc, int n, unsigned char *i, FILE *list)
270 fprintf(list, "%s = ", desc);
272 for (j = 0; j < n && res; j++)
274 res &= xdr_u_char(xd, &i[j]);
277 fprintf(list, "%02x", i[j]);
292 static void do_cpt_int_err(XDR *xd, const char *desc, int *i, FILE *list)
294 if (do_cpt_int(xd, desc, i, list) < 0)
300 static void do_cpt_step_err(XDR *xd, const char *desc, gmx_int64_t *i, FILE *list)
303 char buf[STEPSTRSIZE];
305 res = xdr_int64(xd, i);
312 fprintf(list, "%s = %s\n", desc, gmx_step_str(*i, buf));
316 static void do_cpt_double_err(XDR *xd, const char *desc, double *f, FILE *list)
320 res = xdr_double(xd, f);
327 fprintf(list, "%s = %f\n", desc, *f);
331 static void do_cpt_real_err(XDR *xd, real *f)
336 res = xdr_double(xd, f);
338 res = xdr_float(xd, f);
346 static void do_cpt_n_rvecs_err(XDR *xd, const char *desc, int n, rvec f[], FILE *list)
350 for (i = 0; i < n; i++)
352 for (j = 0; j < DIM; j++)
354 do_cpt_real_err(xd, &f[i][j]);
360 pr_rvecs(list, 0, desc, f, n);
364 /* If nval >= 0, nval is used; on read this should match the passed value.
365 * If nval n<0, *nptr is used; on read the value is stored in nptr
367 static int do_cpte_reals_low(XDR *xd, int cptp, int ecpt, int sflags,
368 int nval, int *nptr, real **v,
369 FILE *list, int erealtype)
373 int dtc = xdr_datatype_float;
375 int dtc = xdr_datatype_double;
377 real *vp, *va = NULL;
392 gmx_incons("*ntpr=NULL in do_cpte_reals_low");
397 res = xdr_int(xd, &nf);
408 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), nval, nf);
417 res = xdr_int(xd, &dt);
424 fprintf(stderr, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
425 st_names(cptp, ecpt), xdr_datatype_names[dtc],
426 xdr_datatype_names[dt]);
428 if (list || !(sflags & (1<<ecpt)))
441 if (dt == xdr_datatype_float)
443 if (dtc == xdr_datatype_float)
451 res = xdr_vector(xd, (char *)vf, nf,
452 (unsigned int)sizeof(float), (xdrproc_t)xdr_float);
457 if (dtc != xdr_datatype_float)
459 for (i = 0; i < nf; i++)
468 if (dtc == xdr_datatype_double)
470 /* cppcheck-suppress invalidPointerCast
471 * Only executed if real is anyhow double */
478 res = xdr_vector(xd, (char *)vd, nf,
479 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
484 if (dtc != xdr_datatype_double)
486 for (i = 0; i < nf; i++)
499 pr_reals(list, 0, st_names(cptp, ecpt), vp, nf);
502 pr_rvecs(list, 0, st_names(cptp, ecpt), (rvec *)vp, nf/3);
505 gmx_incons("Unknown checkpoint real type");
517 /* This function stores n along with the reals for reading,
518 * but on reading it assumes that n matches the value in the checkpoint file,
519 * a fatal error is generated when this is not the case.
521 static int do_cpte_reals(XDR *xd, int cptp, int ecpt, int sflags,
522 int n, real **v, FILE *list)
524 return do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, v, list, ecprREAL);
527 /* This function does the same as do_cpte_reals,
528 * except that on reading it ignores the passed value of *n
529 * and stored the value read from the checkpoint file in *n.
531 static int do_cpte_n_reals(XDR *xd, int cptp, int ecpt, int sflags,
532 int *n, real **v, FILE *list)
534 return do_cpte_reals_low(xd, cptp, ecpt, sflags, -1, n, v, list, ecprREAL);
537 static int do_cpte_real(XDR *xd, int cptp, int ecpt, int sflags,
540 return do_cpte_reals_low(xd, cptp, ecpt, sflags, 1, NULL, &r, list, ecprREAL);
543 static int do_cpte_ints(XDR *xd, int cptp, int ecpt, int sflags,
544 int n, int **v, FILE *list)
547 int dtc = xdr_datatype_int;
552 res = xdr_int(xd, &nf);
557 if (list == NULL && v != NULL && nf != n)
559 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
562 res = xdr_int(xd, &dt);
569 gmx_fatal(FARGS, "Type mismatch for state entry %s, code type is %s, file type is %s\n",
570 st_names(cptp, ecpt), xdr_datatype_names[dtc],
571 xdr_datatype_names[dt]);
573 if (list || !(sflags & (1<<ecpt)) || v == NULL)
586 res = xdr_vector(xd, (char *)vp, nf,
587 (unsigned int)sizeof(int), (xdrproc_t)xdr_int);
594 pr_ivec(list, 0, st_names(cptp, ecpt), vp, nf, TRUE);
604 static int do_cpte_int(XDR *xd, int cptp, int ecpt, int sflags,
607 return do_cpte_ints(xd, cptp, ecpt, sflags, 1, &i, list);
610 static int do_cpte_doubles(XDR *xd, int cptp, int ecpt, int sflags,
611 int n, double **v, FILE *list)
614 int dtc = xdr_datatype_double;
615 double *vp, *va = NULL;
619 res = xdr_int(xd, &nf);
624 if (list == NULL && nf != n)
626 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
629 res = xdr_int(xd, &dt);
636 gmx_fatal(FARGS, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
637 st_names(cptp, ecpt), xdr_datatype_names[dtc],
638 xdr_datatype_names[dt]);
640 if (list || !(sflags & (1<<ecpt)))
653 res = xdr_vector(xd, (char *)vp, nf,
654 (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
661 pr_doubles(list, 0, st_names(cptp, ecpt), vp, nf);
671 static int do_cpte_double(XDR *xd, int cptp, int ecpt, int sflags,
672 double *r, FILE *list)
674 return do_cpte_doubles(xd, cptp, ecpt, sflags, 1, &r, list);
678 static int do_cpte_rvecs(XDR *xd, int cptp, int ecpt, int sflags,
679 int n, rvec **v, FILE *list)
681 return do_cpte_reals_low(xd, cptp, ecpt, sflags,
682 n*DIM, NULL, (real **)v, list, ecprRVEC);
685 static int do_cpte_matrix(XDR *xd, int cptp, int ecpt, int sflags,
686 matrix v, FILE *list)
691 vr = (real *)&(v[0][0]);
692 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
693 DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
695 if (list && ret == 0)
697 pr_rvecs(list, 0, st_names(cptp, ecpt), v, DIM);
704 static int do_cpte_nmatrix(XDR *xd, int cptp, int ecpt, int sflags,
705 int n, real **v, FILE *list)
709 char name[CPTSTRLEN];
716 for (i = 0; i < n; i++)
718 reti = do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, &(v[i]), NULL, ecprREAL);
719 if (list && reti == 0)
721 sprintf(name, "%s[%d]", st_names(cptp, ecpt), i);
722 pr_reals(list, 0, name, v[i], n);
732 static int do_cpte_matrices(XDR *xd, int cptp, int ecpt, int sflags,
733 int n, matrix **v, FILE *list)
736 matrix *vp, *va = NULL;
742 res = xdr_int(xd, &nf);
747 if (list == NULL && nf != n)
749 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
751 if (list || !(sflags & (1<<ecpt)))
764 snew(vr, nf*DIM*DIM);
765 for (i = 0; i < nf; i++)
767 for (j = 0; j < DIM; j++)
769 for (k = 0; k < DIM; k++)
771 vr[(i*DIM+j)*DIM+k] = vp[i][j][k];
775 ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
776 nf*DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
777 for (i = 0; i < nf; i++)
779 for (j = 0; j < DIM; j++)
781 for (k = 0; k < DIM; k++)
783 vp[i][j][k] = vr[(i*DIM+j)*DIM+k];
789 if (list && ret == 0)
791 for (i = 0; i < nf; i++)
793 pr_rvecs(list, 0, st_names(cptp, ecpt), vp[i], DIM);
804 static void do_cpt_header(XDR *xd, gmx_bool bRead, int *file_version,
805 char **version, char **btime, char **buser, char **bhost,
807 char **fprog, char **ftime,
808 int *eIntegrator, int *simulation_part,
809 gmx_int64_t *step, double *t,
810 int *nnodes, int *dd_nc, int *npme,
811 int *natoms, int *ngtc, int *nnhpres, int *nhchainlength,
812 int *nlambda, int *flags_state,
813 int *flags_eks, int *flags_enh, int *flags_dfh,
814 int *nED, int *eSwapCoords,
830 res = xdr_int(xd, &magic);
833 gmx_fatal(FARGS, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
835 if (magic != CPT_MAGIC1)
837 gmx_fatal(FARGS, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
838 "The checkpoint file is corrupted or not a checkpoint file",
844 gmx_gethostname(fhost, 255);
846 do_cpt_string_err(xd, bRead, "GROMACS version", version, list);
847 do_cpt_string_err(xd, bRead, "GROMACS build time", btime, list);
848 do_cpt_string_err(xd, bRead, "GROMACS build user", buser, list);
849 do_cpt_string_err(xd, bRead, "GROMACS build host", bhost, list);
850 do_cpt_string_err(xd, bRead, "generating program", fprog, list);
851 do_cpt_string_err(xd, bRead, "generation time", ftime, list);
852 *file_version = cpt_version;
853 do_cpt_int_err(xd, "checkpoint file version", file_version, list);
854 if (*file_version > cpt_version)
856 gmx_fatal(FARGS, "Attempting to read a checkpoint file of version %d with code of version %d\n", *file_version, cpt_version);
858 if (*file_version >= 13)
860 do_cpt_int_err(xd, "GROMACS double precision", double_prec, list);
866 if (*file_version >= 12)
868 do_cpt_string_err(xd, bRead, "generating host", &fhost, list);
874 do_cpt_int_err(xd, "#atoms", natoms, list);
875 do_cpt_int_err(xd, "#T-coupling groups", ngtc, list);
876 if (*file_version >= 10)
878 do_cpt_int_err(xd, "#Nose-Hoover T-chains", nhchainlength, list);
884 if (*file_version >= 11)
886 do_cpt_int_err(xd, "#Nose-Hoover T-chains for barostat ", nnhpres, list);
892 if (*file_version >= 14)
894 do_cpt_int_err(xd, "# of total lambda states ", nlambda, list);
900 do_cpt_int_err(xd, "integrator", eIntegrator, list);
901 if (*file_version >= 3)
903 do_cpt_int_err(xd, "simulation part #", simulation_part, list);
907 *simulation_part = 1;
909 if (*file_version >= 5)
911 do_cpt_step_err(xd, "step", step, list);
915 do_cpt_int_err(xd, "step", &idum, list);
918 do_cpt_double_err(xd, "t", t, list);
919 do_cpt_int_err(xd, "#PP-ranks", nnodes, list);
921 do_cpt_int_err(xd, "dd_nc[x]", dd_nc ? &(dd_nc[0]) : &idum, list);
922 do_cpt_int_err(xd, "dd_nc[y]", dd_nc ? &(dd_nc[1]) : &idum, list);
923 do_cpt_int_err(xd, "dd_nc[z]", dd_nc ? &(dd_nc[2]) : &idum, list);
924 do_cpt_int_err(xd, "#PME-only ranks", npme, list);
925 do_cpt_int_err(xd, "state flags", flags_state, list);
926 if (*file_version >= 4)
928 do_cpt_int_err(xd, "ekin data flags", flags_eks, list);
929 do_cpt_int_err(xd, "energy history flags", flags_enh, list);
934 *flags_enh = (*flags_state >> (estORIRE_DTAV+1));
935 *flags_state = (*flags_state & ~((1<<(estORIRE_DTAV+1)) |
936 (1<<(estORIRE_DTAV+2)) |
937 (1<<(estORIRE_DTAV+3))));
939 if (*file_version >= 14)
941 do_cpt_int_err(xd, "df history flags", flags_dfh, list);
948 if (*file_version >= 15)
950 do_cpt_int_err(xd, "ED data sets", nED, list);
956 if (*file_version >= 16)
958 do_cpt_int_err(xd, "swap", eSwapCoords, list);
962 static int do_cpt_footer(XDR *xd, int file_version)
967 if (file_version >= 2)
970 res = xdr_int(xd, &magic);
975 if (magic != CPT_MAGIC2)
984 static int do_cpt_state(XDR *xd, gmx_bool bRead,
985 int fflags, t_state *state,
995 nnht = state->nhchainlength*state->ngtc;
996 nnhtp = state->nhchainlength*state->nnhpres;
998 if (bRead) /* we need to allocate space for dfhist if we are reading */
1000 init_df_history(&state->dfhist, state->dfhist.nlambda);
1003 sflags = state->flags;
1004 for (i = 0; (i < estNR && ret == 0); i++)
1006 if (fflags & (1<<i))
1010 case estLAMBDA: ret = do_cpte_reals(xd, cptpEST, i, sflags, efptNR, &(state->lambda), list); break;
1011 case estFEPSTATE: ret = do_cpte_int (xd, cptpEST, i, sflags, &state->fep_state, list); break;
1012 case estBOX: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box, list); break;
1013 case estBOX_REL: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->box_rel, list); break;
1014 case estBOXV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->boxv, list); break;
1015 case estPRES_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->pres_prev, list); break;
1016 case estSVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->svir_prev, list); break;
1017 case estFVIR_PREV: ret = do_cpte_matrix(xd, cptpEST, i, sflags, state->fvir_prev, list); break;
1018 case estNH_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_xi, list); break;
1019 case estNH_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_vxi, list); break;
1020 case estNHPRES_XI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_xi, list); break;
1021 case estNHPRES_VXI: ret = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_vxi, list); break;
1022 case estTC_INT: ret = do_cpte_doubles(xd, cptpEST, i, sflags, state->ngtc, &state->therm_integral, list); break;
1023 case estVETA: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->veta, list); break;
1024 case estVOL0: ret = do_cpte_real(xd, cptpEST, i, sflags, &state->vol0, list); break;
1025 case estX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->x, list); break;
1026 case estV: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->v, list); break;
1027 case estSDX: ret = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->sd_X, list); break;
1028 /* The RNG entries are no longer written,
1029 * the next 4 lines are only for reading old files.
1031 case estLD_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1032 case estLD_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1033 case estMC_RNG: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1034 case estMC_RNGI: ret = do_cpte_ints(xd, cptpEST, i, sflags, 0, NULL, list); break;
1035 case estDISRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.disre_initf, list); break;
1036 case estDISRE_RM3TAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.ndisrepairs, &state->hist.disre_rm3tav, list); break;
1037 case estORIRE_INITF: ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.orire_initf, list); break;
1038 case estORIRE_DTAV: ret = do_cpte_n_reals(xd, cptpEST, i, sflags, &state->hist.norire_Dtav, &state->hist.orire_Dtav, list); break;
1040 gmx_fatal(FARGS, "Unknown state entry %d\n"
1041 "You are probably reading a new checkpoint file with old code", i);
1049 static int do_cpt_ekinstate(XDR *xd, int fflags, ekinstate_t *ekins,
1057 for (i = 0; (i < eeksNR && ret == 0); i++)
1059 if (fflags & (1<<i))
1064 case eeksEKIN_N: ret = do_cpte_int(xd, cptpEEKS, i, fflags, &ekins->ekin_n, list); break;
1065 case eeksEKINH: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh, list); break;
1066 case eeksEKINF: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinf, list); break;
1067 case eeksEKINO: ret = do_cpte_matrices(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinh_old, list); break;
1068 case eeksEKINTOTAL: ret = do_cpte_matrix(xd, cptpEEKS, i, fflags, ekins->ekin_total, list); break;
1069 case eeksEKINSCALEF: ret = do_cpte_doubles(xd, cptpEEKS, i, fflags, ekins->ekin_n, &ekins->ekinscalef_nhc, list); break;
1070 case eeksVSCALE: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->vscale_nhc, list); break;
1071 case eeksEKINSCALEH: ret = do_cpte_doubles(xd, 1, cptpEEKS, fflags, ekins->ekin_n, &ekins->ekinscaleh_nhc, list); break;
1072 case eeksDEKINDL: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->dekindl, list); break;
1073 case eeksMVCOS: ret = do_cpte_real(xd, 1, cptpEEKS, fflags, &ekins->mvcos, list); break;
1075 gmx_fatal(FARGS, "Unknown ekin data state entry %d\n"
1076 "You are probably reading a new checkpoint file with old code", i);
1085 static int do_cpt_swapstate(XDR *xd, gmx_bool bRead, swapstate_t *swapstate, FILE *list)
1089 int swap_cpt_version = 1;
1092 if (eswapNO == swapstate->eSwapCoords)
1097 swapstate->bFromCpt = bRead;
1099 do_cpt_int_err(xd, "swap checkpoint version", &swap_cpt_version, list);
1100 do_cpt_int_err(xd, "swap coupling steps", &swapstate->nAverage, list);
1102 /* When reading, init_swapcoords has not been called yet,
1103 * so we have to allocate memory first. */
1105 for (ic = 0; ic < eCompNR; ic++)
1107 for (ii = 0; ii < eIonNR; ii++)
1111 do_cpt_int_err(xd, "swap requested atoms", &swapstate->nat_req[ic][ii], list);
1115 do_cpt_int_err(xd, "swap requested atoms p", swapstate->nat_req_p[ic][ii], list);
1120 do_cpt_int_err(xd, "swap influx netto", &swapstate->inflow_netto[ic][ii], list);
1124 do_cpt_int_err(xd, "swap influx netto p", swapstate->inflow_netto_p[ic][ii], list);
1127 if (bRead && (NULL == swapstate->nat_past[ic][ii]) )
1129 snew(swapstate->nat_past[ic][ii], swapstate->nAverage);
1132 for (j = 0; j < swapstate->nAverage; j++)
1136 do_cpt_int_err(xd, "swap past atom counts", &swapstate->nat_past[ic][ii][j], list);
1140 do_cpt_int_err(xd, "swap past atom counts p", &swapstate->nat_past_p[ic][ii][j], list);
1146 /* Ion flux per channel */
1147 for (ic = 0; ic < eChanNR; ic++)
1149 for (ii = 0; ii < eIonNR; ii++)
1153 do_cpt_int_err(xd, "channel flux", &swapstate->fluxfromAtoB[ic][ii], list);
1157 do_cpt_int_err(xd, "channel flux p", swapstate->fluxfromAtoB_p[ic][ii], list);
1162 /* Ion flux leakage */
1165 snew(swapstate->fluxleak, 1);
1167 do_cpt_int_err(xd, "flux leakage", swapstate->fluxleak, list);
1170 do_cpt_int_err(xd, "number of ions", &swapstate->nions, list);
1174 snew(swapstate->channel_label, swapstate->nions);
1175 snew(swapstate->comp_from, swapstate->nions);
1178 do_cpt_u_chars(xd, "channel history", swapstate->nions, swapstate->channel_label, list);
1179 do_cpt_u_chars(xd, "domain history", swapstate->nions, swapstate->comp_from, list);
1181 /* Save the last known whole positions to checkpoint
1182 * file to be able to also make multimeric channels whole in PBC */
1183 do_cpt_int_err(xd, "Ch0 atoms", &swapstate->nat[eChan0], list);
1184 do_cpt_int_err(xd, "Ch1 atoms", &swapstate->nat[eChan1], list);
1187 snew(swapstate->xc_old_whole[eChan0], swapstate->nat[eChan0]);
1188 snew(swapstate->xc_old_whole[eChan1], swapstate->nat[eChan1]);
1189 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], swapstate->xc_old_whole[eChan0], list);
1190 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], swapstate->xc_old_whole[eChan1], list);
1194 do_cpt_n_rvecs_err(xd, "Ch0 whole x", swapstate->nat[eChan0], *swapstate->xc_old_whole_p[eChan0], list);
1195 do_cpt_n_rvecs_err(xd, "Ch1 whole x", swapstate->nat[eChan1], *swapstate->xc_old_whole_p[eChan1], list);
1202 static int do_cpt_enerhist(XDR *xd, gmx_bool bRead,
1203 int fflags, energyhistory_t *enerhist,
1214 enerhist->nsteps = 0;
1216 enerhist->nsteps_sim = 0;
1217 enerhist->nsum_sim = 0;
1218 enerhist->dht = NULL;
1220 if (fflags & (1<< eenhENERGY_DELTA_H_NN) )
1222 snew(enerhist->dht, 1);
1223 enerhist->dht->ndh = NULL;
1224 enerhist->dht->dh = NULL;
1225 enerhist->dht->start_lambda_set = FALSE;
1229 for (i = 0; (i < eenhNR && ret == 0); i++)
1231 if (fflags & (1<<i))
1235 case eenhENERGY_N: ret = do_cpte_int(xd, cptpEENH, i, fflags, &enerhist->nener, list); break;
1236 case eenhENERGY_AVER: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_ave, list); break;
1237 case eenhENERGY_SUM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum, list); break;
1238 case eenhENERGY_NSUM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum, list); break;
1239 case eenhENERGY_SUM_SIM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum_sim, list); break;
1240 case eenhENERGY_NSUM_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum_sim, list); break;
1241 case eenhENERGY_NSTEPS: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps, list); break;
1242 case eenhENERGY_NSTEPS_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps_sim, list); break;
1243 case eenhENERGY_DELTA_H_NN: do_cpt_int_err(xd, eenh_names[i], &(enerhist->dht->nndh), list);
1244 if (bRead) /* now allocate memory for it */
1246 snew(enerhist->dht->dh, enerhist->dht->nndh);
1247 snew(enerhist->dht->ndh, enerhist->dht->nndh);
1248 for (j = 0; j < enerhist->dht->nndh; j++)
1250 enerhist->dht->ndh[j] = 0;
1251 enerhist->dht->dh[j] = NULL;
1255 case eenhENERGY_DELTA_H_LIST:
1256 for (j = 0; j < enerhist->dht->nndh; j++)
1258 ret = do_cpte_n_reals(xd, cptpEENH, i, fflags, &enerhist->dht->ndh[j], &(enerhist->dht->dh[j]), list);
1261 case eenhENERGY_DELTA_H_STARTTIME:
1262 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_time), list); break;
1263 case eenhENERGY_DELTA_H_STARTLAMBDA:
1264 ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_lambda), list); break;
1266 gmx_fatal(FARGS, "Unknown energy history entry %d\n"
1267 "You are probably reading a new checkpoint file with old code", i);
1272 if ((fflags & (1<<eenhENERGY_SUM)) && !(fflags & (1<<eenhENERGY_SUM_SIM)))
1274 /* Assume we have an old file format and copy sum to sum_sim */
1275 srenew(enerhist->ener_sum_sim, enerhist->nener);
1276 for (i = 0; i < enerhist->nener; i++)
1278 enerhist->ener_sum_sim[i] = enerhist->ener_sum[i];
1282 if ( (fflags & (1<<eenhENERGY_NSUM)) &&
1283 !(fflags & (1<<eenhENERGY_NSTEPS)))
1285 /* Assume we have an old file format and copy nsum to nsteps */
1286 enerhist->nsteps = enerhist->nsum;
1288 if ( (fflags & (1<<eenhENERGY_NSUM_SIM)) &&
1289 !(fflags & (1<<eenhENERGY_NSTEPS_SIM)))
1291 /* Assume we have an old file format and copy nsum to nsteps */
1292 enerhist->nsteps_sim = enerhist->nsum_sim;
1298 static int do_cpt_df_hist(XDR *xd, int fflags, df_history_t *dfhist, FILE *list)
1303 nlambda = dfhist->nlambda;
1306 for (i = 0; (i < edfhNR && ret == 0); i++)
1308 if (fflags & (1<<i))
1312 case edfhBEQUIL: ret = do_cpte_int(xd, cptpEDFH, i, fflags, &dfhist->bEquil, list); break;
1313 case edfhNATLAMBDA: ret = do_cpte_ints(xd, cptpEDFH, i, fflags, nlambda, &dfhist->n_at_lam, list); break;
1314 case edfhWLHISTO: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->wl_histo, list); break;
1315 case edfhWLDELTA: ret = do_cpte_real(xd, cptpEDFH, i, fflags, &dfhist->wl_delta, list); break;
1316 case edfhSUMWEIGHTS: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_weights, list); break;
1317 case edfhSUMDG: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_dg, list); break;
1318 case edfhSUMMINVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_minvar, list); break;
1319 case edfhSUMVAR: ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_variance, list); break;
1320 case edfhACCUMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p, list); break;
1321 case edfhACCUMM: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m, list); break;
1322 case edfhACCUMP2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p2, list); break;
1323 case edfhACCUMM2: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m2, list); break;
1324 case edfhTIJ: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij, list); break;
1325 case edfhTIJEMP: ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij_empirical, list); break;
1328 gmx_fatal(FARGS, "Unknown df history entry %d\n"
1329 "You are probably reading a new checkpoint file with old code", i);
1338 /* This function stores the last whole configuration of the reference and
1339 * average structure in the .cpt file
1341 static int do_cpt_EDstate(XDR *xd, gmx_bool bRead,
1342 edsamstate_t *EDstate, FILE *list)
1349 EDstate->bFromCpt = bRead;
1351 if (EDstate->nED <= 0)
1356 /* When reading, init_edsam has not been called yet,
1357 * so we have to allocate memory first. */
1360 snew(EDstate->nref, EDstate->nED);
1361 snew(EDstate->old_sref, EDstate->nED);
1362 snew(EDstate->nav, EDstate->nED);
1363 snew(EDstate->old_sav, EDstate->nED);
1366 /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1367 for (i = 0; i < EDstate->nED; i++)
1369 /* Reference structure SREF */
1370 sprintf(buf, "ED%d # of atoms in reference structure", i+1);
1371 do_cpt_int_err(xd, buf, &EDstate->nref[i], list);
1372 sprintf(buf, "ED%d x_ref", i+1);
1375 snew(EDstate->old_sref[i], EDstate->nref[i]);
1376 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref[i], list);
1380 do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref_p[i], list);
1383 /* Average structure SAV */
1384 sprintf(buf, "ED%d # of atoms in average structure", i+1);
1385 do_cpt_int_err(xd, buf, &EDstate->nav[i], list);
1386 sprintf(buf, "ED%d x_av", i+1);
1389 snew(EDstate->old_sav[i], EDstate->nav[i]);
1390 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav[i], list);
1394 do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav_p[i], list);
1402 static int do_cpt_files(XDR *xd, gmx_bool bRead,
1403 gmx_file_position_t **p_outputfiles, int *nfiles,
1404 FILE *list, int file_version)
1408 gmx_off_t mask = 0xFFFFFFFFL;
1409 int offset_high, offset_low;
1411 gmx_file_position_t *outputfiles;
1413 if (do_cpt_int(xd, "number of output files", nfiles, list) != 0)
1420 snew(*p_outputfiles, *nfiles);
1423 outputfiles = *p_outputfiles;
1425 for (i = 0; i < *nfiles; i++)
1427 /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1430 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1431 strncpy(outputfiles[i].filename, buf, CPTSTRLEN-1);
1437 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1441 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1445 outputfiles[i].offset = ( ((gmx_off_t) offset_high) << 32 ) | ( (gmx_off_t) offset_low & mask );
1449 buf = outputfiles[i].filename;
1450 do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1452 offset = outputfiles[i].offset;
1460 offset_low = (int) (offset & mask);
1461 offset_high = (int) ((offset >> 32) & mask);
1463 if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1467 if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1472 if (file_version >= 8)
1474 if (do_cpt_int(xd, "file_checksum_size", &(outputfiles[i].chksum_size),
1479 if (do_cpt_u_chars(xd, "file_checksum", 16, outputfiles[i].chksum, list) != 0)
1486 outputfiles[i].chksum_size = -1;
1493 void write_checkpoint(const char *fn, gmx_bool bNumberAndKeep,
1494 FILE *fplog, t_commrec *cr,
1495 int eIntegrator, int simulation_part,
1496 gmx_bool bExpanded, int elamstats,
1497 gmx_int64_t step, double t, t_state *state)
1507 char *fntemp; /* the temporary checkpoint file name */
1509 char timebuf[STRLEN];
1510 int nppnodes, npmenodes;
1511 char buf[1024], suffix[5+STEPSTRSIZE], sbuf[STEPSTRSIZE];
1512 gmx_file_position_t *outputfiles;
1515 int flags_eks, flags_enh, flags_dfh;
1518 if (DOMAINDECOMP(cr))
1520 nppnodes = cr->dd->nnodes;
1521 npmenodes = cr->npmenodes;
1529 #ifndef GMX_NO_RENAME
1530 /* make the new temporary filename */
1531 snew(fntemp, strlen(fn)+5+STEPSTRSIZE);
1533 fntemp[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1534 sprintf(suffix, "_%s%s", "step", gmx_step_str(step, sbuf));
1535 strcat(fntemp, suffix);
1536 strcat(fntemp, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1538 /* if we can't rename, we just overwrite the cpt file.
1539 * dangerous if interrupted.
1541 snew(fntemp, strlen(fn));
1545 gmx_ctime_r(&now, timebuf, STRLEN);
1549 fprintf(fplog, "Writing checkpoint, step %s at %s\n\n",
1550 gmx_step_str(step, buf), timebuf);
1553 /* Get offsets for open files */
1554 gmx_fio_get_output_file_positions(&outputfiles, &noutputfiles);
1556 fp = gmx_fio_open(fntemp, "w");
1558 if (state->ekinstate.bUpToDate)
1561 ((1<<eeksEKIN_N) | (1<<eeksEKINH) | (1<<eeksEKINF) |
1562 (1<<eeksEKINO) | (1<<eeksEKINSCALEF) | (1<<eeksEKINSCALEH) |
1563 (1<<eeksVSCALE) | (1<<eeksDEKINDL) | (1<<eeksMVCOS));
1571 if (state->enerhist.nsum > 0 || state->enerhist.nsum_sim > 0)
1573 flags_enh |= (1<<eenhENERGY_N);
1574 if (state->enerhist.nsum > 0)
1576 flags_enh |= ((1<<eenhENERGY_AVER) | (1<<eenhENERGY_SUM) |
1577 (1<<eenhENERGY_NSTEPS) | (1<<eenhENERGY_NSUM));
1579 if (state->enerhist.nsum_sim > 0)
1581 flags_enh |= ((1<<eenhENERGY_SUM_SIM) | (1<<eenhENERGY_NSTEPS_SIM) |
1582 (1<<eenhENERGY_NSUM_SIM));
1584 if (state->enerhist.dht)
1586 flags_enh |= ( (1<< eenhENERGY_DELTA_H_NN) |
1587 (1<< eenhENERGY_DELTA_H_LIST) |
1588 (1<< eenhENERGY_DELTA_H_STARTTIME) |
1589 (1<< eenhENERGY_DELTA_H_STARTLAMBDA) );
1595 flags_dfh = ((1<<edfhBEQUIL) | (1<<edfhNATLAMBDA) | (1<<edfhSUMWEIGHTS) | (1<<edfhSUMDG) |
1596 (1<<edfhTIJ) | (1<<edfhTIJEMP));
1599 flags_dfh |= ((1<<edfhWLDELTA) | (1<<edfhWLHISTO));
1601 if ((elamstats == elamstatsMINVAR) || (elamstats == elamstatsBARKER) || (elamstats == elamstatsMETROPOLIS))
1603 flags_dfh |= ((1<<edfhACCUMP) | (1<<edfhACCUMM) | (1<<edfhACCUMP2) | (1<<edfhACCUMM2)
1604 | (1<<edfhSUMMINVAR) | (1<<edfhSUMVAR));
1612 /* We can check many more things now (CPU, acceleration, etc), but
1613 * it is highly unlikely to have two separate builds with exactly
1614 * the same version, user, time, and build host!
1617 version = gmx_strdup(gmx_version());
1618 btime = gmx_strdup(BUILD_TIME);
1619 buser = gmx_strdup(BUILD_USER);
1620 bhost = gmx_strdup(BUILD_HOST);
1622 double_prec = GMX_CPT_BUILD_DP;
1623 fprog = gmx_strdup(Program());
1625 ftime = &(timebuf[0]);
1627 do_cpt_header(gmx_fio_getxdr(fp), FALSE, &file_version,
1628 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1629 &eIntegrator, &simulation_part, &step, &t, &nppnodes,
1630 DOMAINDECOMP(cr) ? cr->dd->nc : NULL, &npmenodes,
1631 &state->natoms, &state->ngtc, &state->nnhpres,
1632 &state->nhchainlength, &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
1633 &state->edsamstate.nED, &state->swapstate.eSwapCoords,
1642 if ((do_cpt_state(gmx_fio_getxdr(fp), FALSE, state->flags, state, NULL) < 0) ||
1643 (do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL) < 0) ||
1644 (do_cpt_enerhist(gmx_fio_getxdr(fp), FALSE, flags_enh, &state->enerhist, NULL) < 0) ||
1645 (do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL) < 0) ||
1646 (do_cpt_EDstate(gmx_fio_getxdr(fp), FALSE, &state->edsamstate, NULL) < 0) ||
1647 (do_cpt_swapstate(gmx_fio_getxdr(fp), FALSE, &state->swapstate, NULL) < 0) ||
1648 (do_cpt_files(gmx_fio_getxdr(fp), FALSE, &outputfiles, &noutputfiles, NULL,
1651 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1654 do_cpt_footer(gmx_fio_getxdr(fp), file_version);
1656 /* we really, REALLY, want to make sure to physically write the checkpoint,
1657 and all the files it depends on, out to disk. Because we've
1658 opened the checkpoint with gmx_fio_open(), it's in our list
1660 ret = gmx_fio_all_output_fsync();
1666 "Cannot fsync '%s'; maybe you are out of disk space?",
1667 gmx_fio_getname(ret));
1669 if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV) == NULL)
1679 if (gmx_fio_close(fp) != 0)
1681 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1684 /* we don't move the checkpoint if the user specified they didn't want it,
1685 or if the fsyncs failed */
1686 #ifndef GMX_NO_RENAME
1687 if (!bNumberAndKeep && !ret)
1691 /* Rename the previous checkpoint file */
1693 buf[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1694 strcat(buf, "_prev");
1695 strcat(buf, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1697 /* we copy here so that if something goes wrong between now and
1698 * the rename below, there's always a state.cpt.
1699 * If renames are atomic (such as in POSIX systems),
1700 * this copying should be unneccesary.
1702 gmx_file_copy(fn, buf, FALSE);
1703 /* We don't really care if this fails:
1704 * there's already a new checkpoint.
1707 gmx_file_rename(fn, buf);
1710 if (gmx_file_rename(fntemp, fn) != 0)
1712 gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
1715 #endif /* GMX_NO_RENAME */
1721 /*code for alternate checkpointing scheme. moved from top of loop over
1723 fcRequestCheckPoint();
1724 if (fcCheckPointParallel( cr->nodeid, NULL, 0) == 0)
1726 gmx_fatal( 3, __FILE__, __LINE__, "Checkpoint error on step %d\n", step );
1728 #endif /* end GMX_FAHCORE block */
1731 static void print_flag_mismatch(FILE *fplog, int sflags, int fflags)
1735 fprintf(fplog, "\nState entry mismatch between the simulation and the checkpoint file\n");
1736 fprintf(fplog, "Entries which are not present in the checkpoint file will not be updated\n");
1737 fprintf(fplog, " %24s %11s %11s\n", "", "simulation", "checkpoint");
1738 for (i = 0; i < estNR; i++)
1740 if ((sflags & (1<<i)) || (fflags & (1<<i)))
1742 fprintf(fplog, " %24s %11s %11s\n",
1744 (sflags & (1<<i)) ? " present " : "not present",
1745 (fflags & (1<<i)) ? " present " : "not present");
1750 static void check_int(FILE *fplog, const char *type, int p, int f, gmx_bool *mm)
1752 FILE *fp = fplog ? fplog : stderr;
1756 fprintf(fp, " %s mismatch,\n", type);
1757 fprintf(fp, " current program: %d\n", p);
1758 fprintf(fp, " checkpoint file: %d\n", f);
1764 static void check_string(FILE *fplog, const char *type, const char *p,
1765 const char *f, gmx_bool *mm)
1767 FILE *fp = fplog ? fplog : stderr;
1769 if (strcmp(p, f) != 0)
1771 fprintf(fp, " %s mismatch,\n", type);
1772 fprintf(fp, " current program: %s\n", p);
1773 fprintf(fp, " checkpoint file: %s\n", f);
1779 static void check_match(FILE *fplog,
1781 char *btime, char *buser, char *bhost, int double_prec,
1783 t_commrec *cr, int npp_f, int npme_f,
1784 ivec dd_nc, ivec dd_nc_f)
1787 gmx_bool mm = FALSE;
1788 gmx_bool patchlevel_differs = FALSE;
1789 gmx_bool version_differs = FALSE;
1791 check_string(fplog, "Version", gmx_version(), version, &mm);
1792 patchlevel_differs = mm;
1794 if (patchlevel_differs)
1796 /* Gromacs should be able to continue from checkpoints between
1797 * different patch level versions, but we do not guarantee
1798 * compatibility between different major/minor versions - check this.
1800 int gmx_major, gmx_minor;
1801 int cpt_major, cpt_minor;
1802 sscanf(gmx_version(), "VERSION %5d.%5d", &gmx_major, &gmx_minor);
1803 sscanf(version, "VERSION %5d.%5d", &cpt_major, &cpt_minor);
1804 version_differs = (gmx_major != cpt_major || gmx_minor != cpt_minor);
1807 check_string(fplog, "Build time", BUILD_TIME, btime, &mm);
1808 check_string(fplog, "Build user", BUILD_USER, buser, &mm);
1809 check_string(fplog, "Build host", BUILD_HOST, bhost, &mm);
1810 check_int (fplog, "Double prec.", GMX_CPT_BUILD_DP, double_prec, &mm);
1811 check_string(fplog, "Program name", Program(), fprog, &mm);
1813 check_int (fplog, "#ranks", cr->nnodes, npp_f+npme_f, &mm);
1816 check_int (fplog, "#PME-ranks", cr->npmenodes, npme_f, &mm);
1819 if (cr->npmenodes >= 0)
1821 npp -= cr->npmenodes;
1825 check_int (fplog, "#DD-cells[x]", dd_nc[XX], dd_nc_f[XX], &mm);
1826 check_int (fplog, "#DD-cells[y]", dd_nc[YY], dd_nc_f[YY], &mm);
1827 check_int (fplog, "#DD-cells[z]", dd_nc[ZZ], dd_nc_f[ZZ], &mm);
1833 const char msg_version_difference[] =
1834 "The current Gromacs major & minor version are not identical to those that\n"
1835 "generated the checkpoint file. In principle Gromacs does not support\n"
1836 "continuation from checkpoints between different versions, so we advise\n"
1837 "against this. If you still want to try your luck we recommend that you use\n"
1838 "the -noappend flag to keep your output files from the two versions separate.\n"
1839 "This might also work around errors where the output fields in the energy\n"
1840 "file have changed between the different major & minor versions.\n";
1842 const char msg_mismatch_notice[] =
1843 "Gromacs patchlevel, binary or parallel settings differ from previous run.\n"
1844 "Continuation is exact, but not guaranteed to be binary identical.\n";
1846 const char msg_logdetails[] =
1847 "See the log file for details.\n";
1849 if (version_differs)
1851 fprintf(stderr, "%s%s\n", msg_version_difference, fplog ? msg_logdetails : "");
1855 fprintf(fplog, "%s\n", msg_version_difference);
1860 /* Major & minor versions match at least, but something is different. */
1861 fprintf(stderr, "%s%s\n", msg_mismatch_notice, fplog ? msg_logdetails : "");
1864 fprintf(fplog, "%s\n", msg_mismatch_notice);
1870 static void read_checkpoint(const char *fn, FILE **pfplog,
1871 t_commrec *cr, ivec dd_nc,
1872 int eIntegrator, int *init_fep_state, gmx_int64_t *step, double *t,
1873 t_state *state, gmx_bool *bReadEkin,
1874 int *simulation_part,
1875 gmx_bool bAppendOutputFiles, gmx_bool bForceAppend)
1880 char *version, *btime, *buser, *bhost, *fprog, *ftime;
1882 char buf[STEPSTRSIZE];
1883 int eIntegrator_f, nppnodes_f, npmenodes_f;
1885 int natoms, ngtc, nnhpres, nhchainlength, nlambda, fflags, flags_eks, flags_enh, flags_dfh;
1888 gmx_file_position_t *outputfiles;
1890 t_fileio *chksum_file;
1891 FILE * fplog = *pfplog;
1892 unsigned char digest[16];
1893 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1894 struct flock fl; /* don't initialize here: the struct order is OS
1898 const char *int_warn =
1899 "WARNING: The checkpoint file was generated with integrator %s,\n"
1900 " while the simulation uses integrator %s\n\n";
1902 #if !defined __native_client__ && !defined GMX_NATIVE_WINDOWS
1903 fl.l_type = F_WRLCK;
1904 fl.l_whence = SEEK_SET;
1910 fp = gmx_fio_open(fn, "r");
1911 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
1912 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1913 &eIntegrator_f, simulation_part, step, t,
1914 &nppnodes_f, dd_nc_f, &npmenodes_f,
1915 &natoms, &ngtc, &nnhpres, &nhchainlength, &nlambda,
1916 &fflags, &flags_eks, &flags_enh, &flags_dfh,
1917 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
1919 if (bAppendOutputFiles &&
1920 file_version >= 13 && double_prec != GMX_CPT_BUILD_DP)
1922 gmx_fatal(FARGS, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
1925 if (cr == NULL || MASTER(cr))
1927 fprintf(stderr, "\nReading checkpoint file %s generated: %s\n\n",
1931 /* This will not be written if we do appending, since fplog is still NULL then */
1934 fprintf(fplog, "\n");
1935 fprintf(fplog, "Reading checkpoint file %s\n", fn);
1936 fprintf(fplog, " file generated by: %s\n", fprog);
1937 fprintf(fplog, " file generated at: %s\n", ftime);
1938 fprintf(fplog, " GROMACS build time: %s\n", btime);
1939 fprintf(fplog, " GROMACS build user: %s\n", buser);
1940 fprintf(fplog, " GROMACS build host: %s\n", bhost);
1941 fprintf(fplog, " GROMACS double prec.: %d\n", double_prec);
1942 fprintf(fplog, " simulation part #: %d\n", *simulation_part);
1943 fprintf(fplog, " step: %s\n", gmx_step_str(*step, buf));
1944 fprintf(fplog, " time: %f\n", *t);
1945 fprintf(fplog, "\n");
1948 if (natoms != state->natoms)
1950 gmx_fatal(FARGS, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", natoms, state->natoms);
1952 if (ngtc != state->ngtc)
1954 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);
1956 if (nnhpres != state->nnhpres)
1958 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);
1961 if (nlambda != state->dfhist.nlambda)
1963 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);
1966 init_gtc_state(state, state->ngtc, state->nnhpres, nhchainlength); /* need to keep this here to keep the tpr format working */
1967 /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
1969 if (eIntegrator_f != eIntegrator)
1973 fprintf(stderr, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1975 if (bAppendOutputFiles)
1978 "Output file appending requested, but input/checkpoint integrators do not match.\n"
1979 "Stopping the run to prevent you from ruining all your data...\n"
1980 "If you _really_ know what you are doing, try with the -noappend option.\n");
1984 fprintf(fplog, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1992 else if (cr->nnodes == nppnodes_f + npmenodes_f)
1994 if (cr->npmenodes < 0)
1996 cr->npmenodes = npmenodes_f;
1998 int nppnodes = cr->nnodes - cr->npmenodes;
1999 if (nppnodes == nppnodes_f)
2001 for (d = 0; d < DIM; d++)
2005 dd_nc[d] = dd_nc_f[d];
2011 if (fflags != state->flags)
2016 if (bAppendOutputFiles)
2019 "Output file appending requested, but input and checkpoint states are not identical.\n"
2020 "Stopping the run to prevent you from ruining all your data...\n"
2021 "You can try with the -noappend option, and get more info in the log file.\n");
2024 if (getenv("GMX_ALLOW_CPT_MISMATCH") == NULL)
2026 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");
2031 "WARNING: The checkpoint state entries do not match the simulation,\n"
2032 " see the log file for details\n\n");
2038 print_flag_mismatch(fplog, state->flags, fflags);
2045 check_match(fplog, version, btime, buser, bhost, double_prec, fprog,
2046 cr, nppnodes_f, npmenodes_f, dd_nc, dd_nc_f);
2049 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, fflags, state, NULL);
2050 *init_fep_state = state->fep_state; /* there should be a better way to do this than setting it here.
2051 Investigate for 5.0. */
2056 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2061 *bReadEkin = ((flags_eks & (1<<eeksEKINH)) || (flags_eks & (1<<eeksEKINF)) || (flags_eks & (1<<eeksEKINO)) ||
2062 ((flags_eks & (1<<eeksEKINSCALEF)) | (flags_eks & (1<<eeksEKINSCALEH)) | (flags_eks & (1<<eeksVSCALE))));
2064 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2065 flags_enh, &state->enerhist, NULL);
2071 if (file_version < 6)
2073 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.";
2075 fprintf(stderr, "\nWARNING: %s\n\n", warn);
2078 fprintf(fplog, "\nWARNING: %s\n\n", warn);
2080 state->enerhist.nsum = *step;
2081 state->enerhist.nsum_sim = *step;
2084 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2090 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2096 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2102 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, NULL, file_version);
2108 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2113 if (gmx_fio_close(fp) != 0)
2115 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2124 /* If the user wants to append to output files,
2125 * we use the file pointer positions of the output files stored
2126 * in the checkpoint file and truncate the files such that any frames
2127 * written after the checkpoint time are removed.
2128 * All files are md5sum checked such that we can be sure that
2129 * we do not truncate other (maybe imprortant) files.
2131 if (bAppendOutputFiles)
2133 if (fn2ftp(outputfiles[0].filename) != efLOG)
2135 /* make sure first file is log file so that it is OK to use it for
2138 gmx_fatal(FARGS, "The first output file should always be the log "
2139 "file but instead is: %s. Cannot do appending because of this condition.", outputfiles[0].filename);
2141 for (i = 0; i < nfiles; i++)
2143 if (outputfiles[i].offset < 0)
2145 gmx_fatal(FARGS, "The original run wrote a file called '%s' which "
2146 "is larger than 2 GB, but mdrun did not support large file"
2147 " offsets. Can not append. Run mdrun with -noappend",
2148 outputfiles[i].filename);
2151 chksum_file = gmx_fio_open(outputfiles[i].filename, "a");
2154 chksum_file = gmx_fio_open(outputfiles[i].filename, "r+");
2159 /* Note that there are systems where the lock operation
2160 * will succeed, but a second process can also lock the file.
2161 * We should probably try to detect this.
2163 #if defined __native_client__
2167 #elif defined GMX_NATIVE_WINDOWS
2168 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX) == -1)
2170 if (fcntl(fileno(gmx_fio_getfp(chksum_file)), F_SETLK, &fl) == -1)
2173 if (errno == ENOSYS)
2177 gmx_fatal(FARGS, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2181 fprintf(stderr, "\nNOTE: File locking is not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2184 fprintf(fplog, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2188 else if (errno == EACCES || errno == EAGAIN)
2190 gmx_fatal(FARGS, "Failed to lock: %s. Already running "
2191 "simulation?", outputfiles[i].filename);
2195 gmx_fatal(FARGS, "Failed to lock: %s. %s.",
2196 outputfiles[i].filename, strerror(errno));
2201 /* compute md5 chksum */
2202 if (outputfiles[i].chksum_size != -1)
2204 if (gmx_fio_get_file_md5(chksum_file, outputfiles[i].offset,
2205 digest) != outputfiles[i].chksum_size) /*at the end of the call the file position is at the end of the file*/
2207 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.",
2208 outputfiles[i].chksum_size,
2209 outputfiles[i].filename);
2212 if (i == 0) /*log file needs to be seeked in case we need to truncate (other files are truncated below)*/
2214 if (gmx_fio_seek(chksum_file, outputfiles[i].offset))
2216 gmx_fatal(FARGS, "Seek error! Failed to truncate log-file: %s.", strerror(errno));
2221 if (i == 0) /*open log file here - so that lock is never lifted
2222 after chksum is calculated */
2224 *pfplog = gmx_fio_getfp(chksum_file);
2228 gmx_fio_close(chksum_file);
2231 /* compare md5 chksum */
2232 if (outputfiles[i].chksum_size != -1 &&
2233 memcmp(digest, outputfiles[i].chksum, 16) != 0)
2237 fprintf(debug, "chksum for %s: ", outputfiles[i].filename);
2238 for (j = 0; j < 16; j++)
2240 fprintf(debug, "%02x", digest[j]);
2242 fprintf(debug, "\n");
2244 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.",
2245 outputfiles[i].filename);
2250 if (i != 0) /*log file is already seeked to correct position */
2252 #ifdef GMX_NATIVE_WINDOWS
2253 rc = gmx_wintruncate(outputfiles[i].filename, outputfiles[i].offset);
2255 rc = truncate(outputfiles[i].filename, outputfiles[i].offset);
2259 gmx_fatal(FARGS, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfiles[i].filename);
2269 void load_checkpoint(const char *fn, FILE **fplog,
2270 t_commrec *cr, ivec dd_nc,
2271 t_inputrec *ir, t_state *state,
2272 gmx_bool *bReadEkin,
2273 gmx_bool bAppend, gmx_bool bForceAppend)
2280 /* Read the state from the checkpoint file */
2281 read_checkpoint(fn, fplog,
2283 ir->eI, &(ir->fepvals->init_fep_state), &step, &t, state, bReadEkin,
2284 &ir->simulation_part, bAppend, bForceAppend);
2288 gmx_bcast(sizeof(cr->npmenodes), &cr->npmenodes, cr);
2289 gmx_bcast(DIM*sizeof(dd_nc[0]), dd_nc, cr);
2290 gmx_bcast(sizeof(step), &step, cr);
2291 gmx_bcast(sizeof(*bReadEkin), bReadEkin, cr);
2293 ir->bContinuation = TRUE;
2294 if (ir->nsteps >= 0)
2296 ir->nsteps += ir->init_step - step;
2298 ir->init_step = step;
2299 ir->simulation_part += 1;
2302 static void read_checkpoint_data(t_fileio *fp, int *simulation_part,
2303 gmx_int64_t *step, double *t, t_state *state,
2304 int *nfiles, gmx_file_position_t **outputfiles)
2307 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2312 int flags_eks, flags_enh, flags_dfh;
2314 gmx_file_position_t *files_loc = NULL;
2317 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2318 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2319 &eIntegrator, simulation_part, step, t, &nppnodes, dd_nc, &npme,
2320 &state->natoms, &state->ngtc, &state->nnhpres, &state->nhchainlength,
2321 &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
2322 &state->edsamstate.nED, &state->swapstate.eSwapCoords, NULL);
2324 do_cpt_state(gmx_fio_getxdr(fp), TRUE, state->flags, state, NULL);
2329 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state->ekinstate, NULL);
2334 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2335 flags_enh, &state->enerhist, NULL);
2340 ret = do_cpt_df_hist(gmx_fio_getxdr(fp), flags_dfh, &state->dfhist, NULL);
2346 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2352 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state->swapstate, NULL);
2358 ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE,
2359 outputfiles != NULL ? outputfiles : &files_loc,
2360 outputfiles != NULL ? nfiles : &nfiles_loc,
2361 NULL, file_version);
2362 if (files_loc != NULL)
2372 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2386 read_checkpoint_state(const char *fn, int *simulation_part,
2387 gmx_int64_t *step, double *t, t_state *state)
2391 fp = gmx_fio_open(fn, "r");
2392 read_checkpoint_data(fp, simulation_part, step, t, state, NULL, NULL);
2393 if (gmx_fio_close(fp) != 0)
2395 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2399 void read_checkpoint_trxframe(t_fileio *fp, t_trxframe *fr)
2401 /* This next line is nasty because the sub-structures of t_state
2402 * cannot be assumed to be zeroed (or even initialized in ways the
2403 * rest of the code might assume). Using snew would be better, but
2404 * this will all go away for 5.0. */
2406 int simulation_part;
2410 init_state(&state, 0, 0, 0, 0, 0);
2412 read_checkpoint_data(fp, &simulation_part, &step, &t, &state, NULL, NULL);
2414 fr->natoms = state.natoms;
2417 fr->step = gmx_int64_to_int(step,
2418 "conversion of checkpoint to trajectory");
2422 fr->lambda = state.lambda[efptFEP];
2423 fr->fep_state = state.fep_state;
2425 fr->bX = (state.flags & (1<<estX));
2431 fr->bV = (state.flags & (1<<estV));
2438 fr->bBox = (state.flags & (1<<estBOX));
2441 copy_mat(state.box, fr->box);
2446 void list_checkpoint(const char *fn, FILE *out)
2450 char *version, *btime, *buser, *bhost, *fprog, *ftime;
2452 int eIntegrator, simulation_part, nppnodes, npme;
2457 int flags_eks, flags_enh, flags_dfh;
2459 gmx_file_position_t *outputfiles;
2462 init_state(&state, -1, -1, -1, -1, 0);
2464 fp = gmx_fio_open(fn, "r");
2465 do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2466 &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2467 &eIntegrator, &simulation_part, &step, &t, &nppnodes, dd_nc, &npme,
2468 &state.natoms, &state.ngtc, &state.nnhpres, &state.nhchainlength,
2469 &(state.dfhist.nlambda), &state.flags,
2470 &flags_eks, &flags_enh, &flags_dfh, &state.edsamstate.nED,
2471 &state.swapstate.eSwapCoords, out);
2472 ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, state.flags, &state, out);
2477 ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), flags_eks, &state.ekinstate, out);
2482 ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2483 flags_enh, &state.enerhist, out);
2487 ret = do_cpt_df_hist(gmx_fio_getxdr(fp),
2488 flags_dfh, &state.dfhist, out);
2493 ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state.edsamstate, out);
2498 ret = do_cpt_swapstate(gmx_fio_getxdr(fp), TRUE, &state.swapstate, out);
2503 do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, out, file_version);
2508 ret = do_cpt_footer(gmx_fio_getxdr(fp), file_version);
2515 if (gmx_fio_close(fp) != 0)
2517 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2524 static gmx_bool exist_output_file(const char *fnm_cp, int nfile, const t_filenm fnm[])
2528 /* Check if the output file name stored in the checkpoint file
2529 * is one of the output file names of mdrun.
2533 !(is_output(&fnm[i]) && strcmp(fnm_cp, fnm[i].fns[0]) == 0))
2538 return (i < nfile && gmx_fexist(fnm_cp));
2541 /* This routine cannot print tons of data, since it is called before the log file is opened. */
2542 gmx_bool read_checkpoint_simulation_part(const char *filename, int *simulation_part,
2543 gmx_int64_t *cpt_step, t_commrec *cr,
2544 gmx_bool bAppendReq,
2545 int nfile, const t_filenm fnm[],
2546 const char *part_suffix, gmx_bool *bAddPart)
2549 gmx_int64_t step = 0;
2551 /* This next line is nasty because the sub-structures of t_state
2552 * cannot be assumed to be zeroed (or even initialized in ways the
2553 * rest of the code might assume). Using snew would be better, but
2554 * this will all go away for 5.0. */
2557 gmx_file_position_t *outputfiles;
2560 char *fn, suf_up[STRLEN];
2566 if (!gmx_fexist(filename) || (!(fp = gmx_fio_open(filename, "r")) ))
2568 *simulation_part = 0;
2572 init_state(&state, 0, 0, 0, 0, 0);
2574 read_checkpoint_data(fp, simulation_part, &step, &t, &state,
2575 &nfiles, &outputfiles);
2576 if (gmx_fio_close(fp) != 0)
2578 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2585 for (f = 0; f < nfiles; f++)
2587 if (exist_output_file(outputfiles[f].filename, nfile, fnm))
2592 if (nexist == nfiles)
2594 bAppend = bAppendReq;
2596 else if (nexist > 0)
2599 "Output file appending has been requested,\n"
2600 "but some output files listed in the checkpoint file %s\n"
2601 "are not present or are named differently by the current program:\n",
2603 fprintf(stderr, "output files present:");
2604 for (f = 0; f < nfiles; f++)
2606 if (exist_output_file(outputfiles[f].filename,
2609 fprintf(stderr, " %s", outputfiles[f].filename);
2612 fprintf(stderr, "\n");
2613 fprintf(stderr, "output files not present or named differently:");
2614 for (f = 0; f < nfiles; f++)
2616 if (!exist_output_file(outputfiles[f].filename,
2619 fprintf(stderr, " %s", outputfiles[f].filename);
2622 fprintf(stderr, "\n");
2624 gmx_fatal(FARGS, "File appending requested, but %d of the %d output files are not present or are named differently", nfiles-nexist, nfiles);
2632 gmx_fatal(FARGS, "File appending requested, but no output file information is stored in the checkpoint file");
2634 fn = outputfiles[0].filename;
2635 if (strlen(fn) < 4 ||
2636 gmx_strcasecmp(fn+strlen(fn)-4, ftp2ext(efLOG)) == 0)
2638 gmx_fatal(FARGS, "File appending requested, but the log file is not the first file listed in the checkpoint file");
2640 /* Set bAddPart to whether the suffix string '.part' is present
2641 * in the log file name.
2643 strcpy(suf_up, part_suffix);
2645 *bAddPart = (strstr(fn, part_suffix) != NULL ||
2646 strstr(fn, suf_up) != NULL);
2654 gmx_bcast(sizeof(*simulation_part), simulation_part, cr);
2656 if (*simulation_part > 0 && bAppendReq)
2658 gmx_bcast(sizeof(bAppend), &bAppend, cr);
2659 gmx_bcast(sizeof(*bAddPart), bAddPart, cr);
2662 if (NULL != cpt_step)