src/gromacs/gmxlib/checkpoint.c

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