src/gmxlib/checkpoint.c

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