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