Code beautification with uncrustify
[alexxy/gromacs.git] / src / gromacs / gmxlib / checkpoint.c
1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
2  *
3  *
4  * This file is part of Gromacs        Copyright (c) 1991-2008
5  * David van der Spoel, Erik Lindahl, Berk Hess, University of Groningen.
6  *
7  * This program is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU General Public License
9  * as published by the Free Software Foundation; either version 2
10  * of the License, or (at your option) any later version.
11  *
12  * To help us fund GROMACS development, we humbly ask that you cite
13  * the research papers on the package. Check out http://www.gromacs.org
14  *
15  * And Hey:
16  * Gnomes, ROck Monsters And Chili Sauce
17  */
18
19 /* The source code in this file should be thread-safe.
20    Please keep it that way. */
21
22
23 #ifdef HAVE_CONFIG_H
24 #include <config.h>
25 #endif
26
27 #include <string.h>
28 #include <time.h>
29
30 #ifdef HAVE_SYS_TIME_H
31 #include <sys/time.h>
32 #endif
33
34 #ifdef HAVE_UNISTD_H
35 #include <unistd.h>
36 #endif
37
38 #ifdef GMX_NATIVE_WINDOWS
39 /* _chsize_s */
40 #include <io.h>
41 #include <sys/locking.h>
42 #endif
43
44
45 #include "filenm.h"
46 #include "names.h"
47 #include "typedefs.h"
48 #include "smalloc.h"
49 #include "gmxfio.h"
50 #include "xdrf.h"
51 #include "statutil.h"
52 #include "txtdump.h"
53 #include "vec.h"
54 #include "network.h"
55 #include "gmx_random.h"
56 #include "checkpoint.h"
57 #include "futil.h"
58 #include "string2.h"
59 #include <fcntl.h>
60
61 #include "buildinfo.h"
62
63 #ifdef GMX_FAHCORE
64 #include "corewrap.h"
65 #endif
66
67
68 /* Portable version of ctime_r implemented in src/gmxlib/string2.c, but we do not want it declared in public installed headers */
69 char *
70 gmx_ctime_r(const time_t *clock, char *buf, int n);
71
72
73 #define CPT_MAGIC1 171817
74 #define CPT_MAGIC2 171819
75 #define CPTSTRLEN 1024
76
77 #ifdef GMX_DOUBLE
78 #define GMX_CPT_BUILD_DP 1
79 #else
80 #define GMX_CPT_BUILD_DP 0
81 #endif
82
83 /* cpt_version should normally only be changed
84  * when the header of footer format changes.
85  * The state data format itself is backward and forward compatible.
86  * But old code can not read a new entry that is present in the file
87  * (but can read a new format when new entries are not present).
88  */
89 static const int cpt_version = 15;
90
91
92 const char *est_names[estNR] =
93 {
94     "FE-lambda",
95     "box", "box-rel", "box-v", "pres_prev",
96     "nosehoover-xi", "thermostat-integral",
97     "x", "v", "SDx", "CGp", "LD-rng", "LD-rng-i",
98     "disre_initf", "disre_rm3tav",
99     "orire_initf", "orire_Dtav",
100     "svir_prev", "nosehoover-vxi", "v_eta", "vol0", "nhpres_xi", "nhpres_vxi", "fvir_prev", "fep_state", "MC-rng", "MC-rng-i"
101 };
102
103 enum {
104     eeksEKIN_N, eeksEKINH, eeksDEKINDL, eeksMVCOS, eeksEKINF, eeksEKINO, eeksEKINSCALEF, eeksEKINSCALEH, eeksVSCALE, eeksEKINTOTAL, eeksNR
105 };
106
107 const char *eeks_names[eeksNR] =
108 {
109     "Ekin_n", "Ekinh", "dEkindlambda", "mv_cos",
110     "Ekinf", "Ekinh_old", "EkinScaleF_NHC", "EkinScaleH_NHC", "Vscale_NHC", "Ekin_Total"
111 };
112
113 enum {
114     eenhENERGY_N, eenhENERGY_AVER, eenhENERGY_SUM, eenhENERGY_NSUM,
115     eenhENERGY_SUM_SIM, eenhENERGY_NSUM_SIM,
116     eenhENERGY_NSTEPS, eenhENERGY_NSTEPS_SIM,
117     eenhENERGY_DELTA_H_NN,
118     eenhENERGY_DELTA_H_LIST,
119     eenhENERGY_DELTA_H_STARTTIME,
120     eenhENERGY_DELTA_H_STARTLAMBDA,
121     eenhNR
122 };
123
124 const char *eenh_names[eenhNR] =
125 {
126     "energy_n", "energy_aver", "energy_sum", "energy_nsum",
127     "energy_sum_sim", "energy_nsum_sim",
128     "energy_nsteps", "energy_nsteps_sim",
129     "energy_delta_h_nn",
130     "energy_delta_h_list",
131     "energy_delta_h_start_time",
132     "energy_delta_h_start_lambda"
133 };
134
135 /* free energy history variables -- need to be preserved over checkpoint */
136 enum {
137     edfhBEQUIL, edfhNATLAMBDA, edfhWLHISTO, edfhWLDELTA, edfhSUMWEIGHTS, edfhSUMDG, edfhSUMMINVAR, edfhSUMVAR,
138     edfhACCUMP, edfhACCUMM, edfhACCUMP2, edfhACCUMM2, edfhTIJ, edfhTIJEMP, edfhNR
139 };
140 /* free energy history variable names  */
141 const char *edfh_names[edfhNR] =
142 {
143     "bEquilibrated", "N_at_state", "Wang-Landau_Histogram", "Wang-Landau-delta", "Weights", "Free Energies", "minvar", "variance",
144     "accumulated_plus", "accumulated_minus", "accumulated_plus_2",  "accumulated_minus_2", "Tij", "Tij_empirical"
145 };
146
147 #ifdef GMX_NATIVE_WINDOWS
148 static int
149 gmx_wintruncate(const char *filename, __int64 size)
150 {
151 #ifdef GMX_FAHCORE
152     /*we do this elsewhere*/
153     return 0;
154 #else
155     FILE *fp;
156     int   rc;
157
158     fp = fopen(filename, "rb+");
159
160     if (fp == NULL)
161     {
162         return -1;
163     }
164
165     return _chsize_s( fileno(fp), size);
166 #endif
167 }
168 #endif
169
170
171 enum {
172     ecprREAL, ecprRVEC, ecprMATRIX
173 };
174
175 enum {
176     cptpEST, cptpEEKS, cptpEENH, cptpEDFH
177 };
178 /* enums for the different components of checkpoint variables, replacing the hard coded ones.
179    cptpEST - state variables.
180    cptpEEKS - Kinetic energy state variables.
181    cptpEENH - Energy history state variables.
182    cptpEDFH - free energy history variables.
183  */
184
185
186 static const char *st_names(int cptp, int ecpt)
187 {
188     switch (cptp)
189     {
190         case cptpEST: return est_names [ecpt]; break;
191         case cptpEEKS: return eeks_names[ecpt]; break;
192         case cptpEENH: return eenh_names[ecpt]; break;
193         case cptpEDFH: return edfh_names[ecpt]; break;
194     }
195
196     return NULL;
197 }
198
199 static void cp_warning(FILE *fp)
200 {
201     fprintf(fp, "\nWARNING: Checkpoint file is corrupted or truncated\n\n");
202 }
203
204 static void cp_error()
205 {
206     gmx_fatal(FARGS, "Checkpoint file corrupted/truncated, or maybe you are out of disk space?");
207 }
208
209 static void do_cpt_string_err(XDR *xd, gmx_bool bRead, const char *desc, char **s, FILE *list)
210 {
211     bool_t res = 0;
212
213     if (bRead)
214     {
215         snew(*s, CPTSTRLEN);
216     }
217     res = xdr_string(xd, s, CPTSTRLEN);
218     if (res == 0)
219     {
220         cp_error();
221     }
222     if (list)
223     {
224         fprintf(list, "%s = %s\n", desc, *s);
225         sfree(*s);
226     }
227 }
228
229 static int do_cpt_int(XDR *xd, const char *desc, int *i, FILE *list)
230 {
231     bool_t res = 0;
232
233     res = xdr_int(xd, i);
234     if (res == 0)
235     {
236         return -1;
237     }
238     if (list)
239     {
240         fprintf(list, "%s = %d\n", desc, *i);
241     }
242     return 0;
243 }
244
245 static int do_cpt_u_chars(XDR *xd, const char *desc, int n, unsigned char *i, FILE *list)
246 {
247     bool_t res = 1;
248     int    j;
249     if (list)
250     {
251         fprintf(list, "%s = ", desc);
252     }
253     for (j = 0; j < n && res; j++)
254     {
255         res &= xdr_u_char(xd, &i[j]);
256         if (list)
257         {
258             fprintf(list, "%02x", i[j]);
259         }
260     }
261     if (list)
262     {
263         fprintf(list, "\n");
264     }
265     if (res == 0)
266     {
267         return -1;
268     }
269
270     return 0;
271 }
272
273 static void do_cpt_int_err(XDR *xd, const char *desc, int *i, FILE *list)
274 {
275     if (do_cpt_int(xd, desc, i, list) < 0)
276     {
277         cp_error();
278     }
279 }
280
281 static void do_cpt_step_err(XDR *xd, const char *desc, gmx_large_int_t *i, FILE *list)
282 {
283     bool_t res = 0;
284     char   buf[STEPSTRSIZE];
285
286     res = xdr_gmx_large_int(xd, i, "reading checkpoint file");
287     if (res == 0)
288     {
289         cp_error();
290     }
291     if (list)
292     {
293         fprintf(list, "%s = %s\n", desc, gmx_step_str(*i, buf));
294     }
295 }
296
297 static void do_cpt_double_err(XDR *xd, const char *desc, double *f, FILE *list)
298 {
299     bool_t res = 0;
300
301     res = xdr_double(xd, f);
302     if (res == 0)
303     {
304         cp_error();
305     }
306     if (list)
307     {
308         fprintf(list, "%s = %f\n", desc, *f);
309     }
310 }
311
312 static void do_cpt_real_err(XDR *xd, const char *desc, real *f)
313 {
314     bool_t res = 0;
315
316 #ifdef GMX_DOUBLE
317     res = xdr_double(xd, f);
318 #else
319     res = xdr_float(xd, f);
320 #endif
321     if (res == 0)
322     {
323         cp_error();
324     }
325 }
326
327 static void do_cpt_n_rvecs_err(XDR *xd, const char *desc, int n, rvec f[], FILE *list)
328 {
329     int i, j;
330
331     for (i = 0; i < n; i++)
332     {
333         for (j = 0; j < DIM; j++)
334         {
335             do_cpt_real_err(xd, desc, &f[i][j]);
336         }
337     }
338
339     if (list)
340     {
341         pr_rvecs(list, 0, desc, f, n);
342     }
343 }
344
345 /* If nval >= 0, nval is used; on read this should match the passed value.
346  * If nval n<0, *nptr is used; on read the value is stored in nptr
347  */
348 static int do_cpte_reals_low(XDR *xd, int cptp, int ecpt, int sflags,
349                              int nval, int *nptr, real **v,
350                              FILE *list, int erealtype)
351 {
352     bool_t  res = 0;
353 #ifndef GMX_DOUBLE
354     int     dtc = xdr_datatype_float;
355 #else
356     int     dtc = xdr_datatype_double;
357 #endif
358     real   *vp, *va = NULL;
359     float  *vf;
360     double *vd;
361     int     nf, dt, i;
362
363     if (list == NULL)
364     {
365         if (nval >= 0)
366         {
367             nf = nval;
368         }
369         else
370         {
371             if (nptr == NULL)
372             {
373                 gmx_incons("*ntpr=NULL in do_cpte_reals_low");
374             }
375             nf = *nptr;
376         }
377     }
378     res = xdr_int(xd, &nf);
379     if (res == 0)
380     {
381         return -1;
382     }
383     if (list == NULL)
384     {
385         if (nval >= 0)
386         {
387             if (nf != nval)
388             {
389                 gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), nval, nf);
390             }
391         }
392         else
393         {
394             *nptr = nf;
395         }
396     }
397     dt  = dtc;
398     res = xdr_int(xd, &dt);
399     if (res == 0)
400     {
401         return -1;
402     }
403     if (dt != dtc)
404     {
405         fprintf(stderr, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
406                 st_names(cptp, ecpt), xdr_datatype_names[dtc],
407                 xdr_datatype_names[dt]);
408     }
409     if (list || !(sflags & (1<<ecpt)))
410     {
411         snew(va, nf);
412         vp = va;
413     }
414     else
415     {
416         if (*v == NULL)
417         {
418             snew(*v, nf);
419         }
420         vp = *v;
421     }
422     if (dt == xdr_datatype_float)
423     {
424         if (dtc == xdr_datatype_float)
425         {
426             vf = (float *)vp;
427         }
428         else
429         {
430             snew(vf, nf);
431         }
432         res = xdr_vector(xd, (char *)vf, nf,
433                          (unsigned int)sizeof(float), (xdrproc_t)xdr_float);
434         if (res == 0)
435         {
436             return -1;
437         }
438         if (dtc != xdr_datatype_float)
439         {
440             for (i = 0; i < nf; i++)
441             {
442                 vp[i] = vf[i];
443             }
444             sfree(vf);
445         }
446     }
447     else
448     {
449         if (dtc == xdr_datatype_double)
450         {
451             vd = (double *)vp;
452         }
453         else
454         {
455             snew(vd, nf);
456         }
457         res = xdr_vector(xd, (char *)vd, nf,
458                          (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
459         if (res == 0)
460         {
461             return -1;
462         }
463         if (dtc != xdr_datatype_double)
464         {
465             for (i = 0; i < nf; i++)
466             {
467                 vp[i] = vd[i];
468             }
469             sfree(vd);
470         }
471     }
472
473     if (list)
474     {
475         switch (erealtype)
476         {
477             case ecprREAL:
478                 pr_reals(list, 0, st_names(cptp, ecpt), vp, nf);
479                 break;
480             case ecprRVEC:
481                 pr_rvecs(list, 0, st_names(cptp, ecpt), (rvec *)vp, nf/3);
482                 break;
483             default:
484                 gmx_incons("Unknown checkpoint real type");
485         }
486     }
487     if (va)
488     {
489         sfree(va);
490     }
491
492     return 0;
493 }
494
495
496 /* This function stores n along with the reals for reading,
497  * but on reading it assumes that n matches the value in the checkpoint file,
498  * a fatal error is generated when this is not the case.
499  */
500 static int do_cpte_reals(XDR *xd, int cptp, int ecpt, int sflags,
501                          int n, real **v, FILE *list)
502 {
503     return do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, v, list, ecprREAL);
504 }
505
506 /* This function does the same as do_cpte_reals,
507  * except that on reading it ignores the passed value of *n
508  * and stored the value read from the checkpoint file in *n.
509  */
510 static int do_cpte_n_reals(XDR *xd, int cptp, int ecpt, int sflags,
511                            int *n, real **v, FILE *list)
512 {
513     return do_cpte_reals_low(xd, cptp, ecpt, sflags, -1, n, v, list, ecprREAL);
514 }
515
516 static int do_cpte_real(XDR *xd, int cptp, int ecpt, int sflags,
517                         real *r, FILE *list)
518 {
519     int n;
520
521     return do_cpte_reals_low(xd, cptp, ecpt, sflags, 1, NULL, &r, list, ecprREAL);
522 }
523
524 static int do_cpte_ints(XDR *xd, int cptp, int ecpt, int sflags,
525                         int n, int **v, FILE *list)
526 {
527     bool_t res = 0;
528     int    dtc = xdr_datatype_int;
529     int   *vp, *va = NULL;
530     int    nf, dt, i;
531
532     nf  = n;
533     res = xdr_int(xd, &nf);
534     if (res == 0)
535     {
536         return -1;
537     }
538     if (list == NULL && v != NULL && nf != n)
539     {
540         gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
541     }
542     dt  = dtc;
543     res = xdr_int(xd, &dt);
544     if (res == 0)
545     {
546         return -1;
547     }
548     if (dt != dtc)
549     {
550         gmx_fatal(FARGS, "Type mismatch for state entry %s, code type is %s, file type is %s\n",
551                   st_names(cptp, ecpt), xdr_datatype_names[dtc],
552                   xdr_datatype_names[dt]);
553     }
554     if (list || !(sflags & (1<<ecpt)) || v == NULL)
555     {
556         snew(va, nf);
557         vp = va;
558     }
559     else
560     {
561         if (*v == NULL)
562         {
563             snew(*v, nf);
564         }
565         vp = *v;
566     }
567     res = xdr_vector(xd, (char *)vp, nf,
568                      (unsigned int)sizeof(int), (xdrproc_t)xdr_int);
569     if (res == 0)
570     {
571         return -1;
572     }
573     if (list)
574     {
575         pr_ivec(list, 0, st_names(cptp, ecpt), vp, nf, TRUE);
576     }
577     if (va)
578     {
579         sfree(va);
580     }
581
582     return 0;
583 }
584
585 static int do_cpte_int(XDR *xd, int cptp, int ecpt, int sflags,
586                        int *i, FILE *list)
587 {
588     return do_cpte_ints(xd, cptp, ecpt, sflags, 1, &i, list);
589 }
590
591 static int do_cpte_doubles(XDR *xd, int cptp, int ecpt, int sflags,
592                            int n, double **v, FILE *list)
593 {
594     bool_t  res = 0;
595     int     dtc = xdr_datatype_double;
596     double *vp, *va = NULL;
597     int     nf, dt, i;
598
599     nf  = n;
600     res = xdr_int(xd, &nf);
601     if (res == 0)
602     {
603         return -1;
604     }
605     if (list == NULL && nf != n)
606     {
607         gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
608     }
609     dt  = dtc;
610     res = xdr_int(xd, &dt);
611     if (res == 0)
612     {
613         return -1;
614     }
615     if (dt != dtc)
616     {
617         gmx_fatal(FARGS, "Precision mismatch for state entry %s, code precision is %s, file precision is %s\n",
618                   st_names(cptp, ecpt), xdr_datatype_names[dtc],
619                   xdr_datatype_names[dt]);
620     }
621     if (list || !(sflags & (1<<ecpt)))
622     {
623         snew(va, nf);
624         vp = va;
625     }
626     else
627     {
628         if (*v == NULL)
629         {
630             snew(*v, nf);
631         }
632         vp = *v;
633     }
634     res = xdr_vector(xd, (char *)vp, nf,
635                      (unsigned int)sizeof(double), (xdrproc_t)xdr_double);
636     if (res == 0)
637     {
638         return -1;
639     }
640     if (list)
641     {
642         pr_doubles(list, 0, st_names(cptp, ecpt), vp, nf);
643     }
644     if (va)
645     {
646         sfree(va);
647     }
648
649     return 0;
650 }
651
652 static int do_cpte_double(XDR *xd, int cptp, int ecpt, int sflags,
653                           double *r, FILE *list)
654 {
655     return do_cpte_doubles(xd, cptp, ecpt, sflags, 1, &r, list);
656 }
657
658
659 static int do_cpte_rvecs(XDR *xd, int cptp, int ecpt, int sflags,
660                          int n, rvec **v, FILE *list)
661 {
662     int n3;
663
664     return do_cpte_reals_low(xd, cptp, ecpt, sflags,
665                              n*DIM, NULL, (real **)v, list, ecprRVEC);
666 }
667
668 static int do_cpte_matrix(XDR *xd, int cptp, int ecpt, int sflags,
669                           matrix v, FILE *list)
670 {
671     real *vr;
672     real  ret;
673
674     vr  = (real *)&(v[0][0]);
675     ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
676                             DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
677
678     if (list && ret == 0)
679     {
680         pr_rvecs(list, 0, st_names(cptp, ecpt), v, DIM);
681     }
682
683     return ret;
684 }
685
686
687 static int do_cpte_nmatrix(XDR *xd, int cptp, int ecpt, int sflags,
688                            int n, real **v, FILE *list)
689 {
690     int   i;
691     real *vr;
692     real  ret, reti;
693     char  name[CPTSTRLEN];
694
695     ret = 0;
696     if (v == NULL)
697     {
698         snew(v, n);
699     }
700     for (i = 0; i < n; i++)
701     {
702         reti = 0;
703         vr   = v[i];
704         reti = do_cpte_reals_low(xd, cptp, ecpt, sflags, n, NULL, &(v[i]), NULL, ecprREAL);
705         if (list && reti == 0)
706         {
707             sprintf(name, "%s[%d]", st_names(cptp, ecpt), i);
708             pr_reals(list, 0, name, v[i], n);
709         }
710         if (reti == 0)
711         {
712             ret = 0;
713         }
714     }
715     return ret;
716 }
717
718 static int do_cpte_matrices(XDR *xd, int cptp, int ecpt, int sflags,
719                             int n, matrix **v, FILE *list)
720 {
721     bool_t  res = 0;
722     matrix *vp, *va = NULL;
723     real   *vr;
724     int     nf, i, j, k;
725     int     ret;
726
727     nf  = n;
728     res = xdr_int(xd, &nf);
729     if (res == 0)
730     {
731         return -1;
732     }
733     if (list == NULL && nf != n)
734     {
735         gmx_fatal(FARGS, "Count mismatch for state entry %s, code count is %d, file count is %d\n", st_names(cptp, ecpt), n, nf);
736     }
737     if (list || !(sflags & (1<<ecpt)))
738     {
739         snew(va, nf);
740         vp = va;
741     }
742     else
743     {
744         if (*v == NULL)
745         {
746             snew(*v, nf);
747         }
748         vp = *v;
749     }
750     snew(vr, nf*DIM*DIM);
751     for (i = 0; i < nf; i++)
752     {
753         for (j = 0; j < DIM; j++)
754         {
755             for (k = 0; k < DIM; k++)
756             {
757                 vr[(i*DIM+j)*DIM+k] = vp[i][j][k];
758             }
759         }
760     }
761     ret = do_cpte_reals_low(xd, cptp, ecpt, sflags,
762                             nf*DIM*DIM, NULL, &vr, NULL, ecprMATRIX);
763     for (i = 0; i < nf; i++)
764     {
765         for (j = 0; j < DIM; j++)
766         {
767             for (k = 0; k < DIM; k++)
768             {
769                 vp[i][j][k] = vr[(i*DIM+j)*DIM+k];
770             }
771         }
772     }
773     sfree(vr);
774
775     if (list && ret == 0)
776     {
777         for (i = 0; i < nf; i++)
778         {
779             pr_rvecs(list, 0, st_names(cptp, ecpt), vp[i], DIM);
780         }
781     }
782     if (va)
783     {
784         sfree(va);
785     }
786
787     return ret;
788 }
789
790 static void do_cpt_header(XDR *xd, gmx_bool bRead, int *file_version,
791                           char **version, char **btime, char **buser, char **bhost,
792                           int *double_prec,
793                           char **fprog, char **ftime,
794                           int *eIntegrator, int *simulation_part,
795                           gmx_large_int_t *step, double *t,
796                           int *nnodes, int *dd_nc, int *npme,
797                           int *natoms, int *ngtc, int *nnhpres, int *nhchainlength,
798                           int *nlambda, int *flags_state,
799                           int *flags_eks, int *flags_enh, int *flags_dfh,
800                           int *nED,
801                           FILE *list)
802 {
803     bool_t res = 0;
804     int    magic;
805     int    idum = 0;
806     int    i;
807     char  *fhost;
808
809     if (bRead)
810     {
811         magic = -1;
812     }
813     else
814     {
815         magic = CPT_MAGIC1;
816     }
817     res = xdr_int(xd, &magic);
818     if (res == 0)
819     {
820         gmx_fatal(FARGS, "The checkpoint file is empty/corrupted, or maybe you are out of disk space?");
821     }
822     if (magic != CPT_MAGIC1)
823     {
824         gmx_fatal(FARGS, "Start of file magic number mismatch, checkpoint file has %d, should be %d\n"
825                   "The checkpoint file is corrupted or not a checkpoint file",
826                   magic, CPT_MAGIC1);
827     }
828     if (!bRead)
829     {
830         snew(fhost, 255);
831 #ifdef HAVE_UNISTD_H
832         if (gethostname(fhost, 255) != 0)
833         {
834             sprintf(fhost, "unknown");
835         }
836 #else
837         sprintf(fhost, "unknown");
838 #endif
839     }
840     do_cpt_string_err(xd, bRead, "GROMACS version", version, list);
841     do_cpt_string_err(xd, bRead, "GROMACS build time", btime, list);
842     do_cpt_string_err(xd, bRead, "GROMACS build user", buser, list);
843     do_cpt_string_err(xd, bRead, "GROMACS build host", bhost, list);
844     do_cpt_string_err(xd, bRead, "generating program", fprog, list);
845     do_cpt_string_err(xd, bRead, "generation time", ftime, list);
846     *file_version = cpt_version;
847     do_cpt_int_err(xd, "checkpoint file version", file_version, list);
848     if (*file_version > cpt_version)
849     {
850         gmx_fatal(FARGS, "Attempting to read a checkpoint file of version %d with code of version %d\n", *file_version, cpt_version);
851     }
852     if (*file_version >= 13)
853     {
854         do_cpt_int_err(xd, "GROMACS double precision", double_prec, list);
855     }
856     else
857     {
858         *double_prec = -1;
859     }
860     if (*file_version >= 12)
861     {
862         do_cpt_string_err(xd, bRead, "generating host", &fhost, list);
863         if (list == NULL)
864         {
865             sfree(fhost);
866         }
867     }
868     do_cpt_int_err(xd, "#atoms", natoms, list);
869     do_cpt_int_err(xd, "#T-coupling groups", ngtc, list);
870     if (*file_version >= 10)
871     {
872         do_cpt_int_err(xd, "#Nose-Hoover T-chains", nhchainlength, list);
873     }
874     else
875     {
876         *nhchainlength = 1;
877     }
878     if (*file_version >= 11)
879     {
880         do_cpt_int_err(xd, "#Nose-Hoover T-chains for barostat ", nnhpres, list);
881     }
882     else
883     {
884         *nnhpres = 0;
885     }
886     if (*file_version >= 14)
887     {
888         do_cpt_int_err(xd, "# of total lambda states ", nlambda, list);
889     }
890     else
891     {
892         *nlambda = 0;
893     }
894     do_cpt_int_err(xd, "integrator", eIntegrator, list);
895     if (*file_version >= 3)
896     {
897         do_cpt_int_err(xd, "simulation part #", simulation_part, list);
898     }
899     else
900     {
901         *simulation_part = 1;
902     }
903     if (*file_version >= 5)
904     {
905         do_cpt_step_err(xd, "step", step, list);
906     }
907     else
908     {
909         do_cpt_int_err(xd, "step", &idum, list);
910         *step = idum;
911     }
912     do_cpt_double_err(xd, "t", t, list);
913     do_cpt_int_err(xd, "#PP-nodes", nnodes, list);
914     idum = 1;
915     do_cpt_int_err(xd, "dd_nc[x]", dd_nc ? &(dd_nc[0]) : &idum, list);
916     do_cpt_int_err(xd, "dd_nc[y]", dd_nc ? &(dd_nc[1]) : &idum, list);
917     do_cpt_int_err(xd, "dd_nc[z]", dd_nc ? &(dd_nc[2]) : &idum, list);
918     do_cpt_int_err(xd, "#PME-only nodes", npme, list);
919     do_cpt_int_err(xd, "state flags", flags_state, list);
920     if (*file_version >= 4)
921     {
922         do_cpt_int_err(xd, "ekin data flags", flags_eks, list);
923         do_cpt_int_err(xd, "energy history flags", flags_enh, list);
924     }
925     else
926     {
927         *flags_eks   = 0;
928         *flags_enh   = (*flags_state >> (estORIRE_DTAV+1));
929         *flags_state = (*flags_state & ~((1<<(estORIRE_DTAV+1)) |
930                                          (1<<(estORIRE_DTAV+2)) |
931                                          (1<<(estORIRE_DTAV+3))));
932     }
933     if (*file_version >= 14)
934     {
935         do_cpt_int_err(xd, "df history flags", flags_dfh, list);
936     }
937     else
938     {
939         *flags_dfh = 0;
940     }
941
942     if (*file_version >= 15)
943     {
944         do_cpt_int_err(xd, "ED data sets", nED, list);
945     }
946     else
947     {
948         *nED = 0;
949     }
950 }
951
952 static int do_cpt_footer(XDR *xd, gmx_bool bRead, int file_version)
953 {
954     bool_t res = 0;
955     int    magic;
956
957     if (file_version >= 2)
958     {
959         magic = CPT_MAGIC2;
960         res   = xdr_int(xd, &magic);
961         if (res == 0)
962         {
963             cp_error();
964         }
965         if (magic != CPT_MAGIC2)
966         {
967             return -1;
968         }
969     }
970
971     return 0;
972 }
973
974 static int do_cpt_state(XDR *xd, gmx_bool bRead,
975                         int fflags, t_state *state,
976                         gmx_bool bReadRNG, FILE *list)
977 {
978     int    sflags;
979     int  **rng_p, **rngi_p;
980     int    i;
981     int    ret;
982     int    nnht, nnhtp;
983
984     ret = 0;
985
986     nnht  = state->nhchainlength*state->ngtc;
987     nnhtp = state->nhchainlength*state->nnhpres;
988
989     if (bReadRNG)
990     {
991         rng_p  = (int **)&state->ld_rng;
992         rngi_p = &state->ld_rngi;
993     }
994     else
995     {
996         /* Do not read the RNG data */
997         rng_p  = NULL;
998         rngi_p = NULL;
999     }
1000     /* We want the MC_RNG the same across all the notes for now -- lambda MC is global */
1001
1002     sflags = state->flags;
1003     for (i = 0; (i < estNR && ret == 0); i++)
1004     {
1005         if (fflags & (1<<i))
1006         {
1007             switch (i)
1008             {
1009                 case estLAMBDA:  ret      = do_cpte_reals(xd, cptpEST, i, sflags, efptNR, &(state->lambda), list); break;
1010                 case estFEPSTATE: ret     = do_cpte_int (xd, cptpEST, i, sflags, &state->fep_state, list); break;
1011                 case estBOX:     ret      = do_cpte_matrix(xd, cptpEST, i, sflags, state->box, list); break;
1012                 case estBOX_REL: ret      = do_cpte_matrix(xd, cptpEST, i, sflags, state->box_rel, list); break;
1013                 case estBOXV:    ret      = do_cpte_matrix(xd, cptpEST, i, sflags, state->boxv, list); break;
1014                 case estPRES_PREV: ret    = do_cpte_matrix(xd, cptpEST, i, sflags, state->pres_prev, list); break;
1015                 case estSVIR_PREV:  ret   = do_cpte_matrix(xd, cptpEST, i, sflags, state->svir_prev, list); break;
1016                 case estFVIR_PREV:  ret   = do_cpte_matrix(xd, cptpEST, i, sflags, state->fvir_prev, list); break;
1017                 case estNH_XI:   ret      = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_xi, list); break;
1018                 case estNH_VXI:  ret      = do_cpte_doubles(xd, cptpEST, i, sflags, nnht, &state->nosehoover_vxi, list); break;
1019                 case estNHPRES_XI:   ret  = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_xi, list); break;
1020                 case estNHPRES_VXI:  ret  = do_cpte_doubles(xd, cptpEST, i, sflags, nnhtp, &state->nhpres_vxi, list); break;
1021                 case estTC_INT:  ret      = do_cpte_doubles(xd, cptpEST, i, sflags, state->ngtc, &state->therm_integral, list); break;
1022                 case estVETA:    ret      = do_cpte_real(xd, cptpEST, i, sflags, &state->veta, list); break;
1023                 case estVOL0:    ret      = do_cpte_real(xd, cptpEST, i, sflags, &state->vol0, list); break;
1024                 case estX:       ret      = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->x, list); break;
1025                 case estV:       ret      = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->v, list); break;
1026                 case estSDX:     ret      = do_cpte_rvecs(xd, cptpEST, i, sflags, state->natoms, &state->sd_X, list); break;
1027                 case estLD_RNG:  ret      = do_cpte_ints(xd, cptpEST, i, sflags, state->nrng, rng_p, list); break;
1028                 case estLD_RNGI: ret      = do_cpte_ints(xd, cptpEST, i, sflags, state->nrngi, rngi_p, list); break;
1029                 case estMC_RNG:  ret      = do_cpte_ints(xd, cptpEST, i, sflags, state->nmcrng, (int **)&state->mc_rng, list); break;
1030                 case estMC_RNGI: ret      = do_cpte_ints(xd, cptpEST, i, sflags, 1, &state->mc_rngi, list); break;
1031                 case estDISRE_INITF:  ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.disre_initf, list); break;
1032                 case estDISRE_RM3TAV: ret = do_cpte_reals(xd, cptpEST, i, sflags, state->hist.ndisrepairs, &state->hist.disre_rm3tav, list); break;
1033                 case estORIRE_INITF:  ret = do_cpte_real (xd, cptpEST, i, sflags, &state->hist.orire_initf, list); break;
1034                 case estORIRE_DTAV:   ret = do_cpte_reals(xd, cptpEST, i, sflags, state->hist.norire_Dtav, &state->hist.orire_Dtav, list); break;
1035                 default:
1036                     gmx_fatal(FARGS, "Unknown state entry %d\n"
1037                               "You are probably reading a new checkpoint file with old code", i);
1038             }
1039         }
1040     }
1041
1042     return ret;
1043 }
1044
1045 static int do_cpt_ekinstate(XDR *xd, gmx_bool bRead,
1046                             int fflags, ekinstate_t *ekins,
1047                             FILE *list)
1048 {
1049     int  i;
1050     int  ret;
1051
1052     ret = 0;
1053
1054     for (i = 0; (i < eeksNR && ret == 0); i++)
1055     {
1056         if (fflags & (1<<i))
1057         {
1058             switch (i)
1059             {
1060
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;
1071                 default:
1072                     gmx_fatal(FARGS, "Unknown ekin data state entry %d\n"
1073                               "You are probably reading a new checkpoint file with old code", i);
1074             }
1075         }
1076     }
1077
1078     return ret;
1079 }
1080
1081
1082 static int do_cpt_enerhist(XDR *xd, gmx_bool bRead,
1083                            int fflags, energyhistory_t *enerhist,
1084                            FILE *list)
1085 {
1086     int  i;
1087     int  j;
1088     int  ret;
1089
1090     ret = 0;
1091
1092     if (bRead)
1093     {
1094         enerhist->nsteps     = 0;
1095         enerhist->nsum       = 0;
1096         enerhist->nsteps_sim = 0;
1097         enerhist->nsum_sim   = 0;
1098         enerhist->dht        = NULL;
1099
1100         if (fflags & (1<< eenhENERGY_DELTA_H_NN) )
1101         {
1102             snew(enerhist->dht, 1);
1103             enerhist->dht->ndh              = NULL;
1104             enerhist->dht->dh               = NULL;
1105             enerhist->dht->start_lambda_set = FALSE;
1106         }
1107     }
1108
1109     for (i = 0; (i < eenhNR && ret == 0); i++)
1110     {
1111         if (fflags & (1<<i))
1112         {
1113             switch (i)
1114             {
1115                 case eenhENERGY_N:     ret = do_cpte_int(xd, cptpEENH, i, fflags, &enerhist->nener, list); break;
1116                 case eenhENERGY_AVER:  ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_ave, list); break;
1117                 case eenhENERGY_SUM:   ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum, list); break;
1118                 case eenhENERGY_NSUM:  do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum, list); break;
1119                 case eenhENERGY_SUM_SIM: ret = do_cpte_doubles(xd, cptpEENH, i, fflags, enerhist->nener, &enerhist->ener_sum_sim, list); break;
1120                 case eenhENERGY_NSUM_SIM:   do_cpt_step_err(xd, eenh_names[i], &enerhist->nsum_sim, list); break;
1121                 case eenhENERGY_NSTEPS:     do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps, list); break;
1122                 case eenhENERGY_NSTEPS_SIM: do_cpt_step_err(xd, eenh_names[i], &enerhist->nsteps_sim, list); break;
1123                 case eenhENERGY_DELTA_H_NN: do_cpt_int_err(xd, eenh_names[i], &(enerhist->dht->nndh), list);
1124                     if (bRead) /* now allocate memory for it */
1125                     {
1126                         snew(enerhist->dht->dh, enerhist->dht->nndh);
1127                         snew(enerhist->dht->ndh, enerhist->dht->nndh);
1128                         for (j = 0; j < enerhist->dht->nndh; j++)
1129                         {
1130                             enerhist->dht->ndh[j] = 0;
1131                             enerhist->dht->dh[j]  = NULL;
1132                         }
1133                     }
1134                     break;
1135                 case eenhENERGY_DELTA_H_LIST:
1136                     for (j = 0; j < enerhist->dht->nndh; j++)
1137                     {
1138                         ret = do_cpte_n_reals(xd, cptpEENH, i, fflags, &enerhist->dht->ndh[j], &(enerhist->dht->dh[j]), list);
1139                     }
1140                     break;
1141                 case eenhENERGY_DELTA_H_STARTTIME:
1142                     ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_time), list); break;
1143                 case eenhENERGY_DELTA_H_STARTLAMBDA:
1144                     ret = do_cpte_double(xd, cptpEENH, i, fflags, &(enerhist->dht->start_lambda), list); break;
1145                 default:
1146                     gmx_fatal(FARGS, "Unknown energy history entry %d\n"
1147                               "You are probably reading a new checkpoint file with old code", i);
1148             }
1149         }
1150     }
1151
1152     if ((fflags & (1<<eenhENERGY_SUM)) && !(fflags & (1<<eenhENERGY_SUM_SIM)))
1153     {
1154         /* Assume we have an old file format and copy sum to sum_sim */
1155         srenew(enerhist->ener_sum_sim, enerhist->nener);
1156         for (i = 0; i < enerhist->nener; i++)
1157         {
1158             enerhist->ener_sum_sim[i] = enerhist->ener_sum[i];
1159         }
1160         fflags |= (1<<eenhENERGY_SUM_SIM);
1161     }
1162
1163     if ( (fflags & (1<<eenhENERGY_NSUM)) &&
1164          !(fflags & (1<<eenhENERGY_NSTEPS)))
1165     {
1166         /* Assume we have an old file format and copy nsum to nsteps */
1167         enerhist->nsteps = enerhist->nsum;
1168         fflags          |= (1<<eenhENERGY_NSTEPS);
1169     }
1170     if ( (fflags & (1<<eenhENERGY_NSUM_SIM)) &&
1171          !(fflags & (1<<eenhENERGY_NSTEPS_SIM)))
1172     {
1173         /* Assume we have an old file format and copy nsum to nsteps */
1174         enerhist->nsteps_sim = enerhist->nsum_sim;
1175         fflags              |= (1<<eenhENERGY_NSTEPS_SIM);
1176     }
1177
1178     return ret;
1179 }
1180
1181 static int do_cpt_df_hist(XDR *xd, gmx_bool bRead, int fflags, df_history_t *dfhist, FILE *list)
1182 {
1183     int  i, nlambda;
1184     int  ret;
1185
1186     nlambda = dfhist->nlambda;
1187     ret     = 0;
1188
1189     for (i = 0; (i < edfhNR && ret == 0); i++)
1190     {
1191         if (fflags & (1<<i))
1192         {
1193             switch (i)
1194             {
1195                 case edfhBEQUIL:       ret = do_cpte_int(xd, cptpEDFH, i, fflags, &dfhist->bEquil, list); break;
1196                 case edfhNATLAMBDA:    ret = do_cpte_ints(xd, cptpEDFH, i, fflags, nlambda, &dfhist->n_at_lam, list); break;
1197                 case edfhWLHISTO:      ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->wl_histo, list); break;
1198                 case edfhWLDELTA:      ret = do_cpte_real(xd, cptpEDFH, i, fflags, &dfhist->wl_delta, list); break;
1199                 case edfhSUMWEIGHTS:   ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_weights, list); break;
1200                 case edfhSUMDG:        ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_dg, list); break;
1201                 case edfhSUMMINVAR:    ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_minvar, list); break;
1202                 case edfhSUMVAR:       ret = do_cpte_reals(xd, cptpEDFH, i, fflags, nlambda, &dfhist->sum_variance, list); break;
1203                 case edfhACCUMP:       ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p, list); break;
1204                 case edfhACCUMM:       ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m, list); break;
1205                 case edfhACCUMP2:      ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_p2, list); break;
1206                 case edfhACCUMM2:      ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->accum_m2, list); break;
1207                 case edfhTIJ:          ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij, list); break;
1208                 case edfhTIJEMP:       ret = do_cpte_nmatrix(xd, cptpEDFH, i, fflags, nlambda, dfhist->Tij_empirical, list); break;
1209
1210                 default:
1211                     gmx_fatal(FARGS, "Unknown df history entry %d\n"
1212                               "You are probably reading a new checkpoint file with old code", i);
1213             }
1214         }
1215     }
1216
1217     return ret;
1218 }
1219
1220
1221 /* This function stores the last whole configuration of the reference and
1222  * average structure in the .cpt file
1223  */
1224 static int do_cpt_EDstate(XDR *xd, gmx_bool bRead,
1225                           edsamstate_t *EDstate, FILE *list)
1226 {
1227     int  i, j;
1228     int  ret = 0;
1229     char buf[STRLEN];
1230
1231
1232     EDstate->bFromCpt = bRead;
1233
1234     if (EDstate->nED <= 0)
1235     {
1236         return ret;
1237     }
1238
1239     /* When reading, init_edsam has not been called yet,
1240      * so we have to allocate memory first. */
1241     if (bRead)
1242     {
1243         snew(EDstate->nref, EDstate->nED);
1244         snew(EDstate->old_sref, EDstate->nED);
1245         snew(EDstate->nav, EDstate->nED);
1246         snew(EDstate->old_sav, EDstate->nED);
1247     }
1248
1249     /* Read/write the last whole conformation of SREF and SAV for each ED dataset (usually only one) */
1250     for (i = 0; i < EDstate->nED; i++)
1251     {
1252         /* Reference structure SREF */
1253         sprintf(buf, "ED%d # of atoms in reference structure", i+1);
1254         do_cpt_int_err(xd, buf, &EDstate->nref[i], list);
1255         sprintf(buf, "ED%d x_ref", i+1);
1256         if (bRead)
1257         {
1258             snew(EDstate->old_sref[i], EDstate->nref[i]);
1259             do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref[i], list);
1260         }
1261         else
1262         {
1263             do_cpt_n_rvecs_err(xd, buf, EDstate->nref[i], EDstate->old_sref_p[i], list);
1264         }
1265
1266         /* Average structure SAV */
1267         sprintf(buf, "ED%d # of atoms in average structure", i+1);
1268         do_cpt_int_err(xd, buf, &EDstate->nav[i], list);
1269         sprintf(buf, "ED%d x_av", i+1);
1270         if (bRead)
1271         {
1272             snew(EDstate->old_sav[i], EDstate->nav[i]);
1273             do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav[i], list);
1274         }
1275         else
1276         {
1277             do_cpt_n_rvecs_err(xd, buf, EDstate->nav[i], EDstate->old_sav_p[i], list);
1278         }
1279     }
1280
1281     return ret;
1282 }
1283
1284
1285 static int do_cpt_files(XDR *xd, gmx_bool bRead,
1286                         gmx_file_position_t **p_outputfiles, int *nfiles,
1287                         FILE *list, int file_version)
1288 {
1289     int                  i, j;
1290     gmx_off_t            offset;
1291     gmx_off_t            mask = 0xFFFFFFFFL;
1292     int                  offset_high, offset_low;
1293     char                *buf;
1294     gmx_file_position_t *outputfiles;
1295
1296     if (do_cpt_int(xd, "number of output files", nfiles, list) != 0)
1297     {
1298         return -1;
1299     }
1300
1301     if (bRead)
1302     {
1303         snew(*p_outputfiles, *nfiles);
1304     }
1305
1306     outputfiles = *p_outputfiles;
1307
1308     for (i = 0; i < *nfiles; i++)
1309     {
1310         /* 64-bit XDR numbers are not portable, so it is stored as separate high/low fractions */
1311         if (bRead)
1312         {
1313             do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1314             strncpy(outputfiles[i].filename, buf, CPTSTRLEN-1);
1315             if (list == NULL)
1316             {
1317                 sfree(buf);
1318             }
1319
1320             if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1321             {
1322                 return -1;
1323             }
1324             if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1325             {
1326                 return -1;
1327             }
1328 #if (SIZEOF_GMX_OFF_T > 4)
1329             outputfiles[i].offset = ( ((gmx_off_t) offset_high) << 32 ) | ( (gmx_off_t) offset_low & mask );
1330 #else
1331             outputfiles[i].offset = offset_low;
1332 #endif
1333         }
1334         else
1335         {
1336             buf = outputfiles[i].filename;
1337             do_cpt_string_err(xd, bRead, "output filename", &buf, list);
1338             /* writing */
1339             offset      = outputfiles[i].offset;
1340             if (offset == -1)
1341             {
1342                 offset_low  = -1;
1343                 offset_high = -1;
1344             }
1345             else
1346             {
1347 #if (SIZEOF_GMX_OFF_T > 4)
1348                 offset_low  = (int) (offset & mask);
1349                 offset_high = (int) ((offset >> 32) & mask);
1350 #else
1351                 offset_low  = offset;
1352                 offset_high = 0;
1353 #endif
1354             }
1355             if (do_cpt_int(xd, "file_offset_high", &offset_high, list) != 0)
1356             {
1357                 return -1;
1358             }
1359             if (do_cpt_int(xd, "file_offset_low", &offset_low, list) != 0)
1360             {
1361                 return -1;
1362             }
1363         }
1364         if (file_version >= 8)
1365         {
1366             if (do_cpt_int(xd, "file_checksum_size", &(outputfiles[i].chksum_size),
1367                            list) != 0)
1368             {
1369                 return -1;
1370             }
1371             if (do_cpt_u_chars(xd, "file_checksum", 16, outputfiles[i].chksum, list) != 0)
1372             {
1373                 return -1;
1374             }
1375         }
1376         else
1377         {
1378             outputfiles[i].chksum_size = -1;
1379         }
1380     }
1381     return 0;
1382 }
1383
1384
1385 void write_checkpoint(const char *fn, gmx_bool bNumberAndKeep,
1386                       FILE *fplog, t_commrec *cr,
1387                       int eIntegrator, int simulation_part,
1388                       gmx_bool bExpanded, int elamstats,
1389                       gmx_large_int_t step, double t, t_state *state)
1390 {
1391     t_fileio            *fp;
1392     int                  file_version;
1393     char                *version;
1394     char                *btime;
1395     char                *buser;
1396     char                *bhost;
1397     int                  double_prec;
1398     char                *fprog;
1399     char                *fntemp; /* the temporary checkpoint file name */
1400     time_t               now;
1401     char                 timebuf[STRLEN];
1402     int                  nppnodes, npmenodes, flag_64bit;
1403     char                 buf[1024], suffix[5+STEPSTRSIZE], sbuf[STEPSTRSIZE];
1404     gmx_file_position_t *outputfiles;
1405     int                  noutputfiles;
1406     char                *ftime;
1407     int                  flags_eks, flags_enh, flags_dfh, i;
1408     t_fileio            *ret;
1409
1410     if (PAR(cr))
1411     {
1412         if (DOMAINDECOMP(cr))
1413         {
1414             nppnodes  = cr->dd->nnodes;
1415             npmenodes = cr->npmenodes;
1416         }
1417         else
1418         {
1419             nppnodes  = cr->nnodes;
1420             npmenodes = 0;
1421         }
1422     }
1423     else
1424     {
1425         nppnodes  = 1;
1426         npmenodes = 0;
1427     }
1428
1429     /* make the new temporary filename */
1430     snew(fntemp, strlen(fn)+5+STEPSTRSIZE);
1431     strcpy(fntemp, fn);
1432     fntemp[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1433     sprintf(suffix, "_%s%s", "step", gmx_step_str(step, sbuf));
1434     strcat(fntemp, suffix);
1435     strcat(fntemp, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1436
1437     time(&now);
1438     gmx_ctime_r(&now, timebuf, STRLEN);
1439
1440     if (fplog)
1441     {
1442         fprintf(fplog, "Writing checkpoint, step %s at %s\n\n",
1443                 gmx_step_str(step, buf), timebuf);
1444     }
1445
1446     /* Get offsets for open files */
1447     gmx_fio_get_output_file_positions(&outputfiles, &noutputfiles);
1448
1449     fp = gmx_fio_open(fntemp, "w");
1450
1451     if (state->ekinstate.bUpToDate)
1452     {
1453         flags_eks =
1454             ((1<<eeksEKIN_N) | (1<<eeksEKINH) | (1<<eeksEKINF) |
1455              (1<<eeksEKINO) | (1<<eeksEKINSCALEF) | (1<<eeksEKINSCALEH) |
1456              (1<<eeksVSCALE) | (1<<eeksDEKINDL) | (1<<eeksMVCOS));
1457     }
1458     else
1459     {
1460         flags_eks = 0;
1461     }
1462
1463     flags_enh = 0;
1464     if (state->enerhist.nsum > 0 || state->enerhist.nsum_sim > 0)
1465     {
1466         flags_enh |= (1<<eenhENERGY_N);
1467         if (state->enerhist.nsum > 0)
1468         {
1469             flags_enh |= ((1<<eenhENERGY_AVER) | (1<<eenhENERGY_SUM) |
1470                           (1<<eenhENERGY_NSTEPS) | (1<<eenhENERGY_NSUM));
1471         }
1472         if (state->enerhist.nsum_sim > 0)
1473         {
1474             flags_enh |= ((1<<eenhENERGY_SUM_SIM) | (1<<eenhENERGY_NSTEPS_SIM) |
1475                           (1<<eenhENERGY_NSUM_SIM));
1476         }
1477         if (state->enerhist.dht)
1478         {
1479             flags_enh |= ( (1<< eenhENERGY_DELTA_H_NN) |
1480                            (1<< eenhENERGY_DELTA_H_LIST) |
1481                            (1<< eenhENERGY_DELTA_H_STARTTIME) |
1482                            (1<< eenhENERGY_DELTA_H_STARTLAMBDA) );
1483         }
1484     }
1485
1486     if (bExpanded)
1487     {
1488         flags_dfh = ((1<<edfhBEQUIL) | (1<<edfhNATLAMBDA) | (1<<edfhSUMWEIGHTS) |  (1<<edfhSUMDG)  |
1489                      (1<<edfhTIJ) | (1<<edfhTIJEMP));
1490         if (EWL(elamstats))
1491         {
1492             flags_dfh |= ((1<<edfhWLDELTA) | (1<<edfhWLHISTO));
1493         }
1494         if ((elamstats == elamstatsMINVAR) || (elamstats == elamstatsBARKER) || (elamstats == elamstatsMETROPOLIS))
1495         {
1496             flags_dfh |= ((1<<edfhACCUMP) | (1<<edfhACCUMM) | (1<<edfhACCUMP2) | (1<<edfhACCUMM2)
1497                           | (1<<edfhSUMMINVAR) | (1<<edfhSUMVAR));
1498         }
1499     }
1500     else
1501     {
1502         flags_dfh = 0;
1503     }
1504
1505     /* We can check many more things now (CPU, acceleration, etc), but
1506      * it is highly unlikely to have two separate builds with exactly
1507      * the same version, user, time, and build host!
1508      */
1509
1510     version = gmx_strdup(VERSION);
1511     btime   = gmx_strdup(BUILD_TIME);
1512     buser   = gmx_strdup(BUILD_USER);
1513     bhost   = gmx_strdup(BUILD_HOST);
1514
1515     double_prec = GMX_CPT_BUILD_DP;
1516     fprog       = gmx_strdup(Program());
1517
1518     ftime   = &(timebuf[0]);
1519
1520     do_cpt_header(gmx_fio_getxdr(fp), FALSE, &file_version,
1521                   &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1522                   &eIntegrator, &simulation_part, &step, &t, &nppnodes,
1523                   DOMAINDECOMP(cr) ? cr->dd->nc : NULL, &npmenodes,
1524                   &state->natoms, &state->ngtc, &state->nnhpres,
1525                   &state->nhchainlength, &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
1526                   &state->edsamstate.nED,
1527                   NULL);
1528
1529     sfree(version);
1530     sfree(btime);
1531     sfree(buser);
1532     sfree(bhost);
1533     sfree(fprog);
1534
1535     if ((do_cpt_state(gmx_fio_getxdr(fp), FALSE, state->flags, state, TRUE, NULL) < 0)        ||
1536         (do_cpt_ekinstate(gmx_fio_getxdr(fp), FALSE, flags_eks, &state->ekinstate, NULL) < 0) ||
1537         (do_cpt_enerhist(gmx_fio_getxdr(fp), FALSE, flags_enh, &state->enerhist, NULL) < 0)  ||
1538         (do_cpt_df_hist(gmx_fio_getxdr(fp), FALSE, flags_dfh, &state->dfhist, NULL) < 0)  ||
1539         (do_cpt_EDstate(gmx_fio_getxdr(fp), FALSE, &state->edsamstate, NULL) < 0)      ||
1540         (do_cpt_files(gmx_fio_getxdr(fp), FALSE, &outputfiles, &noutputfiles, NULL,
1541                       file_version) < 0))
1542     {
1543         gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1544     }
1545
1546     do_cpt_footer(gmx_fio_getxdr(fp), FALSE, file_version);
1547
1548     /* we really, REALLY, want to make sure to physically write the checkpoint,
1549        and all the files it depends on, out to disk. Because we've
1550        opened the checkpoint with gmx_fio_open(), it's in our list
1551        of open files.  */
1552     ret = gmx_fio_all_output_fsync();
1553
1554     if (ret)
1555     {
1556         char buf[STRLEN];
1557         sprintf(buf,
1558                 "Cannot fsync '%s'; maybe you are out of disk space?",
1559                 gmx_fio_getname(ret));
1560
1561         if (getenv(GMX_IGNORE_FSYNC_FAILURE_ENV) == NULL)
1562         {
1563             gmx_file(buf);
1564         }
1565         else
1566         {
1567             gmx_warning(buf);
1568         }
1569     }
1570
1571     if (gmx_fio_close(fp) != 0)
1572     {
1573         gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
1574     }
1575
1576     /* we don't move the checkpoint if the user specified they didn't want it,
1577        or if the fsyncs failed */
1578     if (!bNumberAndKeep && !ret)
1579     {
1580         if (gmx_fexist(fn))
1581         {
1582             /* Rename the previous checkpoint file */
1583             strcpy(buf, fn);
1584             buf[strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1] = '\0';
1585             strcat(buf, "_prev");
1586             strcat(buf, fn+strlen(fn) - strlen(ftp2ext(fn2ftp(fn))) - 1);
1587 #ifndef GMX_FAHCORE
1588             /* we copy here so that if something goes wrong between now and
1589              * the rename below, there's always a state.cpt.
1590              * If renames are atomic (such as in POSIX systems),
1591              * this copying should be unneccesary.
1592              */
1593             gmx_file_copy(fn, buf, FALSE);
1594             /* We don't really care if this fails:
1595              * there's already a new checkpoint.
1596              */
1597 #else
1598             gmx_file_rename(fn, buf);
1599 #endif
1600         }
1601         if (gmx_file_rename(fntemp, fn) != 0)
1602         {
1603             gmx_file("Cannot rename checkpoint file; maybe you are out of disk space?");
1604         }
1605     }
1606
1607     sfree(outputfiles);
1608     sfree(fntemp);
1609
1610 #ifdef GMX_FAHCORE
1611     /*code for alternate checkpointing scheme.  moved from top of loop over
1612        steps */
1613     fcRequestCheckPoint();
1614     if (fcCheckPointParallel( cr->nodeid, NULL, 0) == 0)
1615     {
1616         gmx_fatal( 3, __FILE__, __LINE__, "Checkpoint error on step %d\n", step );
1617     }
1618 #endif /* end GMX_FAHCORE block */
1619 }
1620
1621 static void print_flag_mismatch(FILE *fplog, int sflags, int fflags)
1622 {
1623     int i;
1624
1625     fprintf(fplog, "\nState entry mismatch between the simulation and the checkpoint file\n");
1626     fprintf(fplog, "Entries which are not present in the checkpoint file will not be updated\n");
1627     fprintf(fplog, "  %24s    %11s    %11s\n", "", "simulation", "checkpoint");
1628     for (i = 0; i < estNR; i++)
1629     {
1630         if ((sflags & (1<<i)) || (fflags & (1<<i)))
1631         {
1632             fprintf(fplog, "  %24s    %11s    %11s\n",
1633                     est_names[i],
1634                     (sflags & (1<<i)) ? "  present  " : "not present",
1635                     (fflags & (1<<i)) ? "  present  " : "not present");
1636         }
1637     }
1638 }
1639
1640 static void check_int(FILE *fplog, const char *type, int p, int f, gmx_bool *mm)
1641 {
1642     FILE *fp = fplog ? fplog : stderr;
1643
1644     if (p != f)
1645     {
1646         fprintf(fp, "  %s mismatch,\n", type);
1647         fprintf(fp, "    current program: %d\n", p);
1648         fprintf(fp, "    checkpoint file: %d\n", f);
1649         fprintf(fp, "\n");
1650         *mm = TRUE;
1651     }
1652 }
1653
1654 static void check_string(FILE *fplog, const char *type, const char *p,
1655                          const char *f, gmx_bool *mm)
1656 {
1657     FILE *fp = fplog ? fplog : stderr;
1658
1659     if (strcmp(p, f) != 0)
1660     {
1661         fprintf(fp, "  %s mismatch,\n", type);
1662         fprintf(fp, "    current program: %s\n", p);
1663         fprintf(fp, "    checkpoint file: %s\n", f);
1664         fprintf(fp, "\n");
1665         *mm = TRUE;
1666     }
1667 }
1668
1669 static void check_match(FILE *fplog,
1670                         char *version,
1671                         char *btime, char *buser, char *bhost, int double_prec,
1672                         char *fprog,
1673                         t_commrec *cr, gmx_bool bPartDecomp, int npp_f, int npme_f,
1674                         ivec dd_nc, ivec dd_nc_f)
1675 {
1676     int      npp;
1677     gmx_bool mm;
1678
1679     mm = FALSE;
1680
1681     check_string(fplog, "Version", VERSION, version, &mm);
1682     check_string(fplog, "Build time", BUILD_TIME, btime, &mm);
1683     check_string(fplog, "Build user", BUILD_USER, buser, &mm);
1684     check_string(fplog, "Build host", BUILD_HOST, bhost, &mm);
1685     check_int   (fplog, "Double prec.", GMX_CPT_BUILD_DP, double_prec, &mm);
1686     check_string(fplog, "Program name", Program(), fprog, &mm);
1687
1688     check_int   (fplog, "#nodes", cr->nnodes, npp_f+npme_f, &mm);
1689     if (bPartDecomp)
1690     {
1691         dd_nc[XX] = 1;
1692         dd_nc[YY] = 1;
1693         dd_nc[ZZ] = 1;
1694     }
1695     if (cr->nnodes > 1)
1696     {
1697         check_int (fplog, "#PME-nodes", cr->npmenodes, npme_f, &mm);
1698
1699         npp = cr->nnodes;
1700         if (cr->npmenodes >= 0)
1701         {
1702             npp -= cr->npmenodes;
1703         }
1704         if (npp == npp_f)
1705         {
1706             check_int (fplog, "#DD-cells[x]", dd_nc[XX], dd_nc_f[XX], &mm);
1707             check_int (fplog, "#DD-cells[y]", dd_nc[YY], dd_nc_f[YY], &mm);
1708             check_int (fplog, "#DD-cells[z]", dd_nc[ZZ], dd_nc_f[ZZ], &mm);
1709         }
1710     }
1711
1712     if (mm)
1713     {
1714         fprintf(stderr,
1715                 "Gromacs binary or parallel settings not identical to previous run.\n"
1716                 "Continuation is exact, but is not guaranteed to be binary identical%s.\n\n",
1717                 fplog ? ",\n see the log file for details" : "");
1718
1719         if (fplog)
1720         {
1721             fprintf(fplog,
1722                     "Gromacs binary or parallel settings not identical to previous run.\n"
1723                     "Continuation is exact, but is not guaranteed to be binary identical.\n\n");
1724         }
1725     }
1726 }
1727
1728 static void read_checkpoint(const char *fn, FILE **pfplog,
1729                             t_commrec *cr, gmx_bool bPartDecomp, ivec dd_nc,
1730                             int eIntegrator, int *init_fep_state, gmx_large_int_t *step, double *t,
1731                             t_state *state, gmx_bool *bReadRNG, gmx_bool *bReadEkin,
1732                             int *simulation_part,
1733                             gmx_bool bAppendOutputFiles, gmx_bool bForceAppend)
1734 {
1735     t_fileio            *fp;
1736     int                  i, j, rc;
1737     int                  file_version;
1738     char                *version, *btime, *buser, *bhost, *fprog, *ftime;
1739     int                  double_prec;
1740     char                 filename[STRLEN], buf[STEPSTRSIZE];
1741     int                  nppnodes, eIntegrator_f, nppnodes_f, npmenodes_f;
1742     ivec                 dd_nc_f;
1743     int                  natoms, ngtc, nnhpres, nhchainlength, nlambda, fflags, flags_eks, flags_enh, flags_dfh;
1744     int                  d;
1745     int                  ret;
1746     gmx_file_position_t *outputfiles;
1747     int                  nfiles;
1748     t_fileio            *chksum_file;
1749     FILE               * fplog = *pfplog;
1750     unsigned char        digest[16];
1751 #ifndef GMX_NATIVE_WINDOWS
1752     struct flock         fl; /* don't initialize here: the struct order is OS
1753                                 dependent! */
1754 #endif
1755
1756     const char *int_warn =
1757         "WARNING: The checkpoint file was generated with integrator %s,\n"
1758         "         while the simulation uses integrator %s\n\n";
1759     const char *sd_note =
1760         "NOTE: The checkpoint file was for %d nodes doing SD or BD,\n"
1761         "      while the simulation uses %d SD or BD nodes,\n"
1762         "      continuation will be exact, except for the random state\n\n";
1763
1764 #ifndef GMX_NATIVE_WINDOWS
1765     fl.l_type   = F_WRLCK;
1766     fl.l_whence = SEEK_SET;
1767     fl.l_start  = 0;
1768     fl.l_len    = 0;
1769     fl.l_pid    = 0;
1770 #endif
1771
1772     if (PARTDECOMP(cr))
1773     {
1774         gmx_fatal(FARGS,
1775                   "read_checkpoint not (yet) supported with particle decomposition");
1776     }
1777
1778     fp = gmx_fio_open(fn, "r");
1779     do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
1780                   &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
1781                   &eIntegrator_f, simulation_part, step, t,
1782                   &nppnodes_f, dd_nc_f, &npmenodes_f,
1783                   &natoms, &ngtc, &nnhpres, &nhchainlength, &nlambda,
1784                   &fflags, &flags_eks, &flags_enh, &flags_dfh,
1785                   &state->edsamstate.nED, NULL);
1786
1787     if (bAppendOutputFiles &&
1788         file_version >= 13 && double_prec != GMX_CPT_BUILD_DP)
1789     {
1790         gmx_fatal(FARGS, "Output file appending requested, but the code and checkpoint file precision (single/double) don't match");
1791     }
1792
1793     if (cr == NULL || MASTER(cr))
1794     {
1795         fprintf(stderr, "\nReading checkpoint file %s generated: %s\n\n",
1796                 fn, ftime);
1797     }
1798
1799     /* This will not be written if we do appending, since fplog is still NULL then */
1800     if (fplog)
1801     {
1802         fprintf(fplog, "\n");
1803         fprintf(fplog, "Reading checkpoint file %s\n", fn);
1804         fprintf(fplog, "  file generated by:     %s\n", fprog);
1805         fprintf(fplog, "  file generated at:     %s\n", ftime);
1806         fprintf(fplog, "  GROMACS build time:    %s\n", btime);
1807         fprintf(fplog, "  GROMACS build user:    %s\n", buser);
1808         fprintf(fplog, "  GROMACS build host:    %s\n", bhost);
1809         fprintf(fplog, "  GROMACS double prec.:  %d\n", double_prec);
1810         fprintf(fplog, "  simulation part #:     %d\n", *simulation_part);
1811         fprintf(fplog, "  step:                  %s\n", gmx_step_str(*step, buf));
1812         fprintf(fplog, "  time:                  %f\n", *t);
1813         fprintf(fplog, "\n");
1814     }
1815
1816     if (natoms != state->natoms)
1817     {
1818         gmx_fatal(FARGS, "Checkpoint file is for a system of %d atoms, while the current system consists of %d atoms", natoms, state->natoms);
1819     }
1820     if (ngtc != state->ngtc)
1821     {
1822         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);
1823     }
1824     if (nnhpres != state->nnhpres)
1825     {
1826         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);
1827     }
1828
1829     if (nlambda != state->dfhist.nlambda)
1830     {
1831         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);
1832     }
1833
1834     init_gtc_state(state, state->ngtc, state->nnhpres, nhchainlength); /* need to keep this here to keep the tpr format working */
1835     /* write over whatever was read; we use the number of Nose-Hoover chains from the checkpoint */
1836
1837     if (eIntegrator_f != eIntegrator)
1838     {
1839         if (MASTER(cr))
1840         {
1841             fprintf(stderr, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1842         }
1843         if (bAppendOutputFiles)
1844         {
1845             gmx_fatal(FARGS,
1846                       "Output file appending requested, but input/checkpoint integrators do not match.\n"
1847                       "Stopping the run to prevent you from ruining all your data...\n"
1848                       "If you _really_ know what you are doing, try with the -noappend option.\n");
1849         }
1850         if (fplog)
1851         {
1852             fprintf(fplog, int_warn, EI(eIntegrator_f), EI(eIntegrator));
1853         }
1854     }
1855
1856     if (!PAR(cr))
1857     {
1858         nppnodes      = 1;
1859         cr->npmenodes = 0;
1860     }
1861     else if (bPartDecomp)
1862     {
1863         nppnodes      = cr->nnodes;
1864         cr->npmenodes = 0;
1865     }
1866     else if (cr->nnodes == nppnodes_f + npmenodes_f)
1867     {
1868         if (cr->npmenodes < 0)
1869         {
1870             cr->npmenodes = npmenodes_f;
1871         }
1872         nppnodes = cr->nnodes - cr->npmenodes;
1873         if (nppnodes == nppnodes_f)
1874         {
1875             for (d = 0; d < DIM; d++)
1876             {
1877                 if (dd_nc[d] == 0)
1878                 {
1879                     dd_nc[d] = dd_nc_f[d];
1880                 }
1881             }
1882         }
1883     }
1884     else
1885     {
1886         /* The number of PP nodes has not been set yet */
1887         nppnodes = -1;
1888     }
1889
1890     if ((EI_SD(eIntegrator) || eIntegrator == eiBD) && nppnodes > 0)
1891     {
1892         /* Correct the RNG state size for the number of PP nodes.
1893          * Such assignments should all be moved to one central function.
1894          */
1895         state->nrng  = nppnodes*gmx_rng_n();
1896         state->nrngi = nppnodes;
1897     }
1898
1899     *bReadRNG = TRUE;
1900     if (fflags != state->flags)
1901     {
1902
1903         if (MASTER(cr))
1904         {
1905             if (bAppendOutputFiles)
1906             {
1907                 gmx_fatal(FARGS,
1908                           "Output file appending requested, but input and checkpoint states are not identical.\n"
1909                           "Stopping the run to prevent you from ruining all your data...\n"
1910                           "You can try with the -noappend option, and get more info in the log file.\n");
1911             }
1912
1913             if (getenv("GMX_ALLOW_CPT_MISMATCH") == NULL)
1914             {
1915                 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");
1916             }
1917             else
1918             {
1919                 fprintf(stderr,
1920                         "WARNING: The checkpoint state entries do not match the simulation,\n"
1921                         "         see the log file for details\n\n");
1922             }
1923         }
1924
1925         if (fplog)
1926         {
1927             print_flag_mismatch(fplog, state->flags, fflags);
1928         }
1929     }
1930     else
1931     {
1932         if ((EI_SD(eIntegrator) || eIntegrator == eiBD) &&
1933             nppnodes != nppnodes_f)
1934         {
1935             *bReadRNG = FALSE;
1936             if (MASTER(cr))
1937             {
1938                 fprintf(stderr, sd_note, nppnodes_f, nppnodes);
1939             }
1940             if (fplog)
1941             {
1942                 fprintf(fplog, sd_note, nppnodes_f, nppnodes);
1943             }
1944         }
1945         if (MASTER(cr))
1946         {
1947             check_match(fplog, version, btime, buser, bhost, double_prec, fprog,
1948                         cr, bPartDecomp, nppnodes_f, npmenodes_f, dd_nc, dd_nc_f);
1949         }
1950     }
1951     ret             = do_cpt_state(gmx_fio_getxdr(fp), TRUE, fflags, state, *bReadRNG, NULL);
1952     *init_fep_state = state->fep_state;  /* there should be a better way to do this than setting it here.
1953                                             Investigate for 5.0. */
1954     if (ret)
1955     {
1956         cp_error();
1957     }
1958     ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), TRUE,
1959                            flags_eks, &state->ekinstate, NULL);
1960     if (ret)
1961     {
1962         cp_error();
1963     }
1964     *bReadEkin = ((flags_eks & (1<<eeksEKINH)) || (flags_eks & (1<<eeksEKINF)) || (flags_eks & (1<<eeksEKINO)) ||
1965                   ((flags_eks & (1<<eeksEKINSCALEF)) | (flags_eks & (1<<eeksEKINSCALEH)) | (flags_eks & (1<<eeksVSCALE))));
1966
1967     ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
1968                           flags_enh, &state->enerhist, NULL);
1969     if (ret)
1970     {
1971         cp_error();
1972     }
1973
1974     ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
1975     if (ret)
1976     {
1977         cp_error();
1978     }
1979
1980     if (file_version < 6)
1981     {
1982         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.";
1983
1984         fprintf(stderr, "\nWARNING: %s\n\n", warn);
1985         if (fplog)
1986         {
1987             fprintf(fplog, "\nWARNING: %s\n\n", warn);
1988         }
1989         state->enerhist.nsum     = *step;
1990         state->enerhist.nsum_sim = *step;
1991     }
1992
1993     ret = do_cpt_df_hist(gmx_fio_getxdr(fp), TRUE,
1994                          flags_dfh, &state->dfhist, NULL);
1995     if (ret)
1996     {
1997         cp_error();
1998     }
1999
2000     ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, NULL, file_version);
2001     if (ret)
2002     {
2003         cp_error();
2004     }
2005
2006     ret = do_cpt_footer(gmx_fio_getxdr(fp), TRUE, file_version);
2007     if (ret)
2008     {
2009         cp_error();
2010     }
2011     if (gmx_fio_close(fp) != 0)
2012     {
2013         gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2014     }
2015
2016     sfree(fprog);
2017     sfree(ftime);
2018     sfree(btime);
2019     sfree(buser);
2020     sfree(bhost);
2021
2022     /* If the user wants to append to output files,
2023      * we use the file pointer positions of the output files stored
2024      * in the checkpoint file and truncate the files such that any frames
2025      * written after the checkpoint time are removed.
2026      * All files are md5sum checked such that we can be sure that
2027      * we do not truncate other (maybe imprortant) files.
2028      */
2029     if (bAppendOutputFiles)
2030     {
2031         if (fn2ftp(outputfiles[0].filename) != efLOG)
2032         {
2033             /* make sure first file is log file so that it is OK to use it for
2034              * locking
2035              */
2036             gmx_fatal(FARGS, "The first output file should always be the log "
2037                       "file but instead is: %s. Cannot do appending because of this condition.", outputfiles[0].filename);
2038         }
2039         for (i = 0; i < nfiles; i++)
2040         {
2041             if (outputfiles[i].offset < 0)
2042             {
2043                 gmx_fatal(FARGS, "The original run wrote a file called '%s' which "
2044                           "is larger than 2 GB, but mdrun did not support large file"
2045                           " offsets. Can not append. Run mdrun with -noappend",
2046                           outputfiles[i].filename);
2047             }
2048 #ifdef GMX_FAHCORE
2049             chksum_file = gmx_fio_open(outputfiles[i].filename, "a");
2050
2051 #else
2052             chksum_file = gmx_fio_open(outputfiles[i].filename, "r+");
2053
2054             /* lock log file */
2055             if (i == 0)
2056             {
2057                 /* Note that there are systems where the lock operation
2058                  * will succeed, but a second process can also lock the file.
2059                  * We should probably try to detect this.
2060                  */
2061 #ifndef GMX_NATIVE_WINDOWS
2062                 if (fcntl(fileno(gmx_fio_getfp(chksum_file)), F_SETLK, &fl)
2063                     == -1)
2064 #else
2065                 if (_locking(fileno(gmx_fio_getfp(chksum_file)), _LK_NBLCK, LONG_MAX) == -1)
2066 #endif
2067                 {
2068                     if (errno == ENOSYS)
2069                     {
2070                         if (!bForceAppend)
2071                         {
2072                             gmx_fatal(FARGS, "File locking is not supported on this system. Use -noappend or specify -append explicitly to append anyhow.");
2073                         }
2074                         else
2075                         {
2076                             fprintf(stderr, "\nNOTE: File locking is not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2077                             if (fplog)
2078                             {
2079                                 fprintf(fplog, "\nNOTE: File locking not supported on this system, will not lock %s\n\n", outputfiles[i].filename);
2080                             }
2081                         }
2082                     }
2083                     else if (errno == EACCES || errno == EAGAIN)
2084                     {
2085                         gmx_fatal(FARGS, "Failed to lock: %s. Already running "
2086                                   "simulation?", outputfiles[i].filename);
2087                     }
2088                     else
2089                     {
2090                         gmx_fatal(FARGS, "Failed to lock: %s. %s.",
2091                                   outputfiles[i].filename, strerror(errno));
2092                     }
2093                 }
2094             }
2095
2096             /* compute md5 chksum */
2097             if (outputfiles[i].chksum_size != -1)
2098             {
2099                 if (gmx_fio_get_file_md5(chksum_file, outputfiles[i].offset,
2100                                          digest) != outputfiles[i].chksum_size) /*at the end of the call the file position is at the end of the file*/
2101                 {
2102                     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.",
2103                               outputfiles[i].chksum_size,
2104                               outputfiles[i].filename);
2105                 }
2106             }
2107             if (i == 0)  /*log file needs to be seeked in case we need to truncate (other files are truncated below)*/
2108             {
2109                 if (gmx_fio_seek(chksum_file, outputfiles[i].offset))
2110                 {
2111                     gmx_fatal(FARGS, "Seek error! Failed to truncate log-file: %s.", strerror(errno));
2112                 }
2113             }
2114 #endif
2115
2116             if (i == 0) /*open log file here - so that lock is never lifted
2117                            after chksum is calculated */
2118             {
2119                 *pfplog = gmx_fio_getfp(chksum_file);
2120             }
2121             else
2122             {
2123                 gmx_fio_close(chksum_file);
2124             }
2125 #ifndef GMX_FAHCORE
2126             /* compare md5 chksum */
2127             if (outputfiles[i].chksum_size != -1 &&
2128                 memcmp(digest, outputfiles[i].chksum, 16) != 0)
2129             {
2130                 if (debug)
2131                 {
2132                     fprintf(debug, "chksum for %s: ", outputfiles[i].filename);
2133                     for (j = 0; j < 16; j++)
2134                     {
2135                         fprintf(debug, "%02x", digest[j]);
2136                     }
2137                     fprintf(debug, "\n");
2138                 }
2139                 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.",
2140                           outputfiles[i].filename);
2141             }
2142 #endif
2143
2144
2145             if (i != 0) /*log file is already seeked to correct position */
2146             {
2147 #ifdef GMX_NATIVE_WINDOWS
2148                 rc = gmx_wintruncate(outputfiles[i].filename, outputfiles[i].offset);
2149 #else
2150                 rc = truncate(outputfiles[i].filename, outputfiles[i].offset);
2151 #endif
2152                 if (rc != 0)
2153                 {
2154                     gmx_fatal(FARGS, "Truncation of file %s failed. Cannot do appending because of this failure.", outputfiles[i].filename);
2155                 }
2156             }
2157         }
2158     }
2159
2160     sfree(outputfiles);
2161 }
2162
2163
2164 void load_checkpoint(const char *fn, FILE **fplog,
2165                      t_commrec *cr, gmx_bool bPartDecomp, ivec dd_nc,
2166                      t_inputrec *ir, t_state *state,
2167                      gmx_bool *bReadRNG, gmx_bool *bReadEkin,
2168                      gmx_bool bAppend, gmx_bool bForceAppend)
2169 {
2170     gmx_large_int_t step;
2171     double          t;
2172
2173     if (SIMMASTER(cr))
2174     {
2175         /* Read the state from the checkpoint file */
2176         read_checkpoint(fn, fplog,
2177                         cr, bPartDecomp, dd_nc,
2178                         ir->eI, &(ir->fepvals->init_fep_state), &step, &t, state, bReadRNG, bReadEkin,
2179                         &ir->simulation_part, bAppend, bForceAppend);
2180     }
2181     if (PAR(cr))
2182     {
2183         gmx_bcast(sizeof(cr->npmenodes), &cr->npmenodes, cr);
2184         gmx_bcast(DIM*sizeof(dd_nc[0]), dd_nc, cr);
2185         gmx_bcast(sizeof(step), &step, cr);
2186         gmx_bcast(sizeof(*bReadRNG), bReadRNG, cr);
2187         gmx_bcast(sizeof(*bReadEkin), bReadEkin, cr);
2188     }
2189     ir->bContinuation    = TRUE;
2190     if (ir->nsteps >= 0)
2191     {
2192         ir->nsteps          += ir->init_step - step;
2193     }
2194     ir->init_step        = step;
2195     ir->simulation_part += 1;
2196 }
2197
2198 static void read_checkpoint_data(t_fileio *fp, int *simulation_part,
2199                                  gmx_large_int_t *step, double *t, t_state *state,
2200                                  gmx_bool bReadRNG,
2201                                  int *nfiles, gmx_file_position_t **outputfiles)
2202 {
2203     int                  file_version;
2204     char                *version, *btime, *buser, *bhost, *fprog, *ftime;
2205     int                  double_prec;
2206     int                  eIntegrator;
2207     int                  nppnodes, npme;
2208     ivec                 dd_nc;
2209     int                  flags_eks, flags_enh, flags_dfh;
2210     int                  nfiles_loc;
2211     gmx_file_position_t *files_loc = NULL;
2212     int                  ret;
2213
2214     do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2215                   &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2216                   &eIntegrator, simulation_part, step, t, &nppnodes, dd_nc, &npme,
2217                   &state->natoms, &state->ngtc, &state->nnhpres, &state->nhchainlength,
2218                   &(state->dfhist.nlambda), &state->flags, &flags_eks, &flags_enh, &flags_dfh,
2219                   &state->edsamstate.nED, NULL);
2220     ret =
2221         do_cpt_state(gmx_fio_getxdr(fp), TRUE, state->flags, state, bReadRNG, NULL);
2222     if (ret)
2223     {
2224         cp_error();
2225     }
2226     ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), TRUE,
2227                            flags_eks, &state->ekinstate, NULL);
2228     if (ret)
2229     {
2230         cp_error();
2231     }
2232     ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2233                           flags_enh, &state->enerhist, NULL);
2234     if (ret)
2235     {
2236         cp_error();
2237     }
2238     ret = do_cpt_df_hist(gmx_fio_getxdr(fp), TRUE,
2239                          flags_dfh, &state->dfhist, NULL);
2240     if (ret)
2241     {
2242         cp_error();
2243     }
2244
2245     ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state->edsamstate, NULL);
2246     if (ret)
2247     {
2248         cp_error();
2249     }
2250
2251     ret = do_cpt_files(gmx_fio_getxdr(fp), TRUE,
2252                        outputfiles != NULL ? outputfiles : &files_loc,
2253                        outputfiles != NULL ? nfiles : &nfiles_loc,
2254                        NULL, file_version);
2255     if (files_loc != NULL)
2256     {
2257         sfree(files_loc);
2258     }
2259
2260     if (ret)
2261     {
2262         cp_error();
2263     }
2264
2265     ret = do_cpt_footer(gmx_fio_getxdr(fp), TRUE, file_version);
2266     if (ret)
2267     {
2268         cp_error();
2269     }
2270
2271     sfree(fprog);
2272     sfree(ftime);
2273     sfree(btime);
2274     sfree(buser);
2275     sfree(bhost);
2276 }
2277
2278 void
2279 read_checkpoint_state(const char *fn, int *simulation_part,
2280                       gmx_large_int_t *step, double *t, t_state *state)
2281 {
2282     t_fileio *fp;
2283
2284     fp = gmx_fio_open(fn, "r");
2285     read_checkpoint_data(fp, simulation_part, step, t, state, FALSE, NULL, NULL);
2286     if (gmx_fio_close(fp) != 0)
2287     {
2288         gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2289     }
2290 }
2291
2292 void read_checkpoint_trxframe(t_fileio *fp, t_trxframe *fr)
2293 {
2294     t_state         state;
2295     int             simulation_part;
2296     gmx_large_int_t step;
2297     double          t;
2298
2299     init_state(&state, 0, 0, 0, 0, 0);
2300
2301     read_checkpoint_data(fp, &simulation_part, &step, &t, &state, FALSE, NULL, NULL);
2302
2303     fr->natoms  = state.natoms;
2304     fr->bTitle  = FALSE;
2305     fr->bStep   = TRUE;
2306     fr->step    = gmx_large_int_to_int(step,
2307                                        "conversion of checkpoint to trajectory");
2308     fr->bTime      = TRUE;
2309     fr->time       = t;
2310     fr->bLambda    = TRUE;
2311     fr->lambda     = state.lambda[efptFEP];
2312     fr->fep_state  = state.fep_state;
2313     fr->bAtoms     = FALSE;
2314     fr->bX         = (state.flags & (1<<estX));
2315     if (fr->bX)
2316     {
2317         fr->x     = state.x;
2318         state.x   = NULL;
2319     }
2320     fr->bV      = (state.flags & (1<<estV));
2321     if (fr->bV)
2322     {
2323         fr->v     = state.v;
2324         state.v   = NULL;
2325     }
2326     fr->bF      = FALSE;
2327     fr->bBox    = (state.flags & (1<<estBOX));
2328     if (fr->bBox)
2329     {
2330         copy_mat(state.box, fr->box);
2331     }
2332     done_state(&state);
2333 }
2334
2335 void list_checkpoint(const char *fn, FILE *out)
2336 {
2337     t_fileio            *fp;
2338     int                  file_version;
2339     char                *version, *btime, *buser, *bhost, *fprog, *ftime;
2340     int                  double_prec;
2341     int                  eIntegrator, simulation_part, nppnodes, npme;
2342     gmx_large_int_t      step;
2343     double               t;
2344     ivec                 dd_nc;
2345     t_state              state;
2346     int                  flags_eks, flags_enh, flags_dfh;
2347     int                  indent;
2348     int                  i, j;
2349     int                  ret;
2350     gmx_file_position_t *outputfiles;
2351     int                  nfiles;
2352
2353     init_state(&state, -1, -1, -1, -1, 0);
2354
2355     fp = gmx_fio_open(fn, "r");
2356     do_cpt_header(gmx_fio_getxdr(fp), TRUE, &file_version,
2357                   &version, &btime, &buser, &bhost, &double_prec, &fprog, &ftime,
2358                   &eIntegrator, &simulation_part, &step, &t, &nppnodes, dd_nc, &npme,
2359                   &state.natoms, &state.ngtc, &state.nnhpres, &state.nhchainlength,
2360                   &(state.dfhist.nlambda), &state.flags,
2361                   &flags_eks, &flags_enh, &flags_dfh, &state.edsamstate.nED, out);
2362     ret = do_cpt_state(gmx_fio_getxdr(fp), TRUE, state.flags, &state, TRUE, out);
2363     if (ret)
2364     {
2365         cp_error();
2366     }
2367     ret = do_cpt_ekinstate(gmx_fio_getxdr(fp), TRUE,
2368                            flags_eks, &state.ekinstate, out);
2369     if (ret)
2370     {
2371         cp_error();
2372     }
2373     ret = do_cpt_enerhist(gmx_fio_getxdr(fp), TRUE,
2374                           flags_enh, &state.enerhist, out);
2375
2376     if (ret == 0)
2377     {
2378         init_df_history(&state.dfhist, state.dfhist.nlambda, 0); /* reinitialize state with correct sizes */
2379         ret = do_cpt_df_hist(gmx_fio_getxdr(fp), TRUE,
2380                              flags_dfh, &state.dfhist, out);
2381     }
2382
2383     if (ret == 0)
2384     {
2385         ret = do_cpt_EDstate(gmx_fio_getxdr(fp), TRUE, &state.edsamstate, out);
2386     }
2387
2388     if (ret == 0)
2389     {
2390         do_cpt_files(gmx_fio_getxdr(fp), TRUE, &outputfiles, &nfiles, out, file_version);
2391     }
2392
2393     if (ret == 0)
2394     {
2395         ret = do_cpt_footer(gmx_fio_getxdr(fp), TRUE, file_version);
2396     }
2397
2398     if (ret)
2399     {
2400         cp_warning(out);
2401     }
2402     if (gmx_fio_close(fp) != 0)
2403     {
2404         gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2405     }
2406
2407     done_state(&state);
2408 }
2409
2410
2411 static gmx_bool exist_output_file(const char *fnm_cp, int nfile, const t_filenm fnm[])
2412 {
2413     int i;
2414
2415     /* Check if the output file name stored in the checkpoint file
2416      * is one of the output file names of mdrun.
2417      */
2418     i = 0;
2419     while (i < nfile &&
2420            !(is_output(&fnm[i]) && strcmp(fnm_cp, fnm[i].fns[0]) == 0))
2421     {
2422         i++;
2423     }
2424
2425     return (i < nfile && gmx_fexist(fnm_cp));
2426 }
2427
2428 /* This routine cannot print tons of data, since it is called before the log file is opened. */
2429 gmx_bool read_checkpoint_simulation_part(const char *filename, int *simulation_part,
2430                                          gmx_large_int_t *cpt_step, t_commrec *cr,
2431                                          gmx_bool bAppendReq,
2432                                          int nfile, const t_filenm fnm[],
2433                                          const char *part_suffix, gmx_bool *bAddPart)
2434 {
2435     t_fileio            *fp;
2436     gmx_large_int_t      step = 0;
2437     double               t;
2438     t_state              state;
2439     int                  nfiles;
2440     gmx_file_position_t *outputfiles;
2441     int                  nexist, f;
2442     gmx_bool             bAppend;
2443     char                *fn, suf_up[STRLEN];
2444
2445     bAppend = FALSE;
2446
2447     if (SIMMASTER(cr))
2448     {
2449         if (!gmx_fexist(filename) || (!(fp = gmx_fio_open(filename, "r")) ))
2450         {
2451             *simulation_part = 0;
2452         }
2453         else
2454         {
2455             init_state(&state, 0, 0, 0, 0, 0);
2456
2457             read_checkpoint_data(fp, simulation_part, &step, &t, &state, FALSE,
2458                                  &nfiles, &outputfiles);
2459             if (gmx_fio_close(fp) != 0)
2460             {
2461                 gmx_file("Cannot read/write checkpoint; corrupt file, or maybe you are out of disk space?");
2462             }
2463             done_state(&state);
2464
2465             if (bAppendReq)
2466             {
2467                 nexist = 0;
2468                 for (f = 0; f < nfiles; f++)
2469                 {
2470                     if (exist_output_file(outputfiles[f].filename, nfile, fnm))
2471                     {
2472                         nexist++;
2473                     }
2474                 }
2475                 if (nexist == nfiles)
2476                 {
2477                     bAppend = bAppendReq;
2478                 }
2479                 else if (nexist > 0)
2480                 {
2481                     fprintf(stderr,
2482                             "Output file appending has been requested,\n"
2483                             "but some output files listed in the checkpoint file %s\n"
2484                             "are not present or are named differently by the current program:\n",
2485                             filename);
2486                     fprintf(stderr, "output files present:");
2487                     for (f = 0; f < nfiles; f++)
2488                     {
2489                         if (exist_output_file(outputfiles[f].filename,
2490                                               nfile, fnm))
2491                         {
2492                             fprintf(stderr, " %s", outputfiles[f].filename);
2493                         }
2494                     }
2495                     fprintf(stderr, "\n");
2496                     fprintf(stderr, "output files not present or named differently:");
2497                     for (f = 0; f < nfiles; f++)
2498                     {
2499                         if (!exist_output_file(outputfiles[f].filename,
2500                                                nfile, fnm))
2501                         {
2502                             fprintf(stderr, " %s", outputfiles[f].filename);
2503                         }
2504                     }
2505                     fprintf(stderr, "\n");
2506
2507                     gmx_fatal(FARGS, "File appending requested, but only %d of the %d output files are present", nexist, nfiles);
2508                 }
2509             }
2510
2511             if (bAppend)
2512             {
2513                 if (nfiles == 0)
2514                 {
2515                     gmx_fatal(FARGS, "File appending requested, but no output file information is stored in the checkpoint file");
2516                 }
2517                 fn = outputfiles[0].filename;
2518                 if (strlen(fn) < 4 ||
2519                     gmx_strcasecmp(fn+strlen(fn)-4, ftp2ext(efLOG)) == 0)
2520                 {
2521                     gmx_fatal(FARGS, "File appending requested, but the log file is not the first file listed in the checkpoint file");
2522                 }
2523                 /* Set bAddPart to whether the suffix string '.part' is present
2524                  * in the log file name.
2525                  */
2526                 strcpy(suf_up, part_suffix);
2527                 upstring(suf_up);
2528                 *bAddPart = (strstr(fn, part_suffix) != NULL ||
2529                              strstr(fn, suf_up) != NULL);
2530             }
2531
2532             sfree(outputfiles);
2533         }
2534     }
2535     if (PAR(cr))
2536     {
2537         gmx_bcast(sizeof(*simulation_part), simulation_part, cr);
2538
2539         if (*simulation_part > 0 && bAppendReq)
2540         {
2541             gmx_bcast(sizeof(bAppend), &bAppend, cr);
2542             gmx_bcast(sizeof(*bAddPart), bAddPart, cr);
2543         }
2544     }
2545     if (NULL != cpt_step)
2546     {
2547         *cpt_step = step;
2548     }
2549
2550     return bAppend;
2551 }