2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
39 /* This file is completely threadsafe - keep it that way! */
48 #include "gromacs/fileio/filenm.h"
49 #include "gromacs/fileio/gmxfio.h"
50 #include "gromacs/fileio/gmxfio-xdr.h"
51 #include "gromacs/legacyheaders/names.h"
52 #include "gromacs/legacyheaders/types/ifunc.h"
53 #include "gromacs/math/vec.h"
54 #include "gromacs/mdtypes/inputrec.h"
55 #include "gromacs/mdtypes/state.h"
56 #include "gromacs/pbcutil/boxutilities.h"
57 #include "gromacs/pbcutil/pbc.h"
58 #include "gromacs/topology/block.h"
59 #include "gromacs/topology/mtop_util.h"
60 #include "gromacs/topology/symtab.h"
61 #include "gromacs/topology/topology.h"
62 #include "gromacs/utility/arraysize.h"
63 #include "gromacs/utility/baseversion.h"
64 #include "gromacs/utility/cstringutil.h"
65 #include "gromacs/utility/fatalerror.h"
66 #include "gromacs/utility/futil.h"
67 #include "gromacs/utility/smalloc.h"
69 #define TPX_TAG_RELEASE "release"
71 /*! \brief Tag string for the file format written to run input files
72 * written by this version of the code.
74 * Change this if you want to change the run input format in a feature
75 * branch. This ensures that there will not be different run input
76 * formats around which cannot be distinguished, while not causing
77 * problems rebasing the feature branch onto upstream changes. When
78 * merging with mainstream GROMACS, set this tag string back to
79 * TPX_TAG_RELEASE, and instead add an element to tpxv.
81 static const char *tpx_tag
= TPX_TAG_RELEASE
;
83 /*! \brief Enum of values that describe the contents of a tpr file
84 * whose format matches a version number
86 * The enum helps the code be more self-documenting and ensure merges
87 * do not silently resolve when two patches make the same bump. When
88 * adding new functionality, add a new element just above tpxv_Count
89 * in this enumeration, and write code below that does the right thing
90 * according to the value of file_version.
93 tpxv_ComputationalElectrophysiology
= 96, /**< support for ion/water position swaps (computational electrophysiology) */
94 tpxv_Use64BitRandomSeed
, /**< change ld_seed from int to gmx_int64_t */
95 tpxv_RestrictedBendingAndCombinedAngleTorsionPotentials
, /**< potentials for supporting coarse-grained force fields */
96 tpxv_InteractiveMolecularDynamics
, /**< interactive molecular dynamics (IMD) */
97 tpxv_RemoveObsoleteParameters1
, /**< remove optimize_fft, dihre_fc, nstcheckpoint */
98 tpxv_PullCoordTypeGeom
, /**< add pull type and geometry per group and flat-bottom */
99 tpxv_PullGeomDirRel
, /**< add pull geometry direction-relative */
100 tpxv_IntermolecularBondeds
, /**< permit inter-molecular bonded interactions in the topology */
101 tpxv_CompElWithSwapLayerOffset
, /**< added parameters for improved CompEl setups */
102 tpxv_Count
/**< the total number of tpxv versions */
105 /*! \brief Version number of the file format written to run input
106 * files by this version of the code.
108 * The tpx_version increases whenever the file format in the main
109 * development branch changes, due to an extension of the tpxv enum above.
110 * Backward compatibility for reading old run input files is maintained
111 * by checking this version number against that of the file and then using
112 * the correct code path.
114 * When developing a feature branch that needs to change the run input
115 * file format, change tpx_tag instead. */
116 static const int tpx_version
= tpxv_Count
- 1;
119 /* This number should only be increased when you edit the TOPOLOGY section
120 * or the HEADER of the tpx format.
121 * This way we can maintain forward compatibility too for all analysis tools
122 * and/or external programs that only need to know the atom/residue names,
123 * charges, and bond connectivity.
125 * It first appeared in tpx version 26, when I also moved the inputrecord
126 * to the end of the tpx file, so we can just skip it if we only
129 * In particular, it must be increased when adding new elements to
130 * ftupd, so that old code can read new .tpr files.
132 static const int tpx_generation
= 26;
134 /* This number should be the most recent backwards incompatible version
135 * I.e., if this number is 9, we cannot read tpx version 9 with this code.
137 static const int tpx_incompatible_version
= 9;
141 /* Struct used to maintain tpx compatibility when function types are added */
143 int fvnr
; /* file version number in which the function type first appeared */
144 int ftype
; /* function type */
148 * The entries should be ordered in:
149 * 1. ascending function type number
150 * 2. ascending file version number
152 static const t_ftupd ftupd
[] = {
153 { 20, F_CUBICBONDS
},
158 { 43, F_TABBONDSNC
},
159 { 70, F_RESTRBONDS
},
160 { tpxv_RestrictedBendingAndCombinedAngleTorsionPotentials
, F_RESTRANGLES
},
161 { 76, F_LINEAR_ANGLES
},
162 { 30, F_CROSS_BOND_BONDS
},
163 { 30, F_CROSS_BOND_ANGLES
},
164 { 30, F_UREY_BRADLEY
},
165 { 34, F_QUARTIC_ANGLES
},
167 { tpxv_RestrictedBendingAndCombinedAngleTorsionPotentials
, F_RESTRDIHS
},
168 { tpxv_RestrictedBendingAndCombinedAngleTorsionPotentials
, F_CBTDIHS
},
177 { 72, F_NPSOLVATION
},
179 { 41, F_LJC_PAIRS_NB
},
182 { 32, F_COUL_RECIP
},
185 { 30, F_POLARIZATION
},
188 { 22, F_DISRESVIOL
},
192 { 26, F_DIHRESVIOL
},
197 { 46, F_ECONSERVED
},
198 { 69, F_VTEMP_NOLONGERUSED
},
200 { 54, F_DVDL_CONSTR
},
201 { 76, F_ANHARM_POL
},
204 { 79, F_DVDL_BONDED
, },
205 { 79, F_DVDL_RESTRAINT
},
206 { 79, F_DVDL_TEMPERATURE
},
208 #define NFTUPD asize(ftupd)
210 /* Needed for backward compatibility */
213 /**************************************************************
215 * Now the higer level routines that do io of the structures and arrays
217 **************************************************************/
218 static void do_pullgrp_tpx_pre95(t_fileio
*fio
,
226 gmx_fio_do_int(fio
, pgrp
->nat
);
229 snew(pgrp
->ind
, pgrp
->nat
);
231 gmx_fio_ndo_int(fio
, pgrp
->ind
, pgrp
->nat
);
232 gmx_fio_do_int(fio
, pgrp
->nweight
);
235 snew(pgrp
->weight
, pgrp
->nweight
);
237 gmx_fio_ndo_real(fio
, pgrp
->weight
, pgrp
->nweight
);
238 gmx_fio_do_int(fio
, pgrp
->pbcatom
);
239 gmx_fio_do_rvec(fio
, pcrd
->vec
);
240 clear_rvec(pcrd
->origin
);
241 gmx_fio_do_rvec(fio
, tmp
);
243 gmx_fio_do_real(fio
, pcrd
->rate
);
244 gmx_fio_do_real(fio
, pcrd
->k
);
245 if (file_version
>= 56)
247 gmx_fio_do_real(fio
, pcrd
->kB
);
255 static void do_pull_group(t_fileio
*fio
, t_pull_group
*pgrp
, gmx_bool bRead
)
257 gmx_fio_do_int(fio
, pgrp
->nat
);
260 snew(pgrp
->ind
, pgrp
->nat
);
262 gmx_fio_ndo_int(fio
, pgrp
->ind
, pgrp
->nat
);
263 gmx_fio_do_int(fio
, pgrp
->nweight
);
266 snew(pgrp
->weight
, pgrp
->nweight
);
268 gmx_fio_ndo_real(fio
, pgrp
->weight
, pgrp
->nweight
);
269 gmx_fio_do_int(fio
, pgrp
->pbcatom
);
272 static void do_pull_coord(t_fileio
*fio
, t_pull_coord
*pcrd
, int file_version
,
273 int ePullOld
, int eGeomOld
, ivec dimOld
)
275 gmx_fio_do_int(fio
, pcrd
->group
[0]);
276 gmx_fio_do_int(fio
, pcrd
->group
[1]);
277 if (file_version
>= tpxv_PullCoordTypeGeom
)
279 gmx_fio_do_int(fio
, pcrd
->eType
);
280 gmx_fio_do_int(fio
, pcrd
->eGeom
);
281 if (pcrd
->eGeom
== epullgDIRRELATIVE
)
283 gmx_fio_do_int(fio
, pcrd
->group
[2]);
284 gmx_fio_do_int(fio
, pcrd
->group
[3]);
286 gmx_fio_do_ivec(fio
, pcrd
->dim
);
290 pcrd
->eType
= ePullOld
;
291 pcrd
->eGeom
= eGeomOld
;
292 copy_ivec(dimOld
, pcrd
->dim
);
294 gmx_fio_do_rvec(fio
, pcrd
->origin
);
295 gmx_fio_do_rvec(fio
, pcrd
->vec
);
296 if (file_version
>= tpxv_PullCoordTypeGeom
)
298 gmx_fio_do_gmx_bool(fio
, pcrd
->bStart
);
302 /* This parameter is only printed, but not actually used by mdrun */
303 pcrd
->bStart
= FALSE
;
305 gmx_fio_do_real(fio
, pcrd
->init
);
306 gmx_fio_do_real(fio
, pcrd
->rate
);
307 gmx_fio_do_real(fio
, pcrd
->k
);
308 gmx_fio_do_real(fio
, pcrd
->kB
);
311 static void do_expandedvals(t_fileio
*fio
, t_expanded
*expand
, t_lambda
*fepvals
, gmx_bool bRead
, int file_version
)
313 int n_lambda
= fepvals
->n_lambda
;
315 /* reset the lambda calculation window */
316 fepvals
->lambda_start_n
= 0;
317 fepvals
->lambda_stop_n
= n_lambda
;
318 if (file_version
>= 79)
324 snew(expand
->init_lambda_weights
, n_lambda
);
326 gmx_fio_ndo_real(fio
, expand
->init_lambda_weights
, n_lambda
);
327 gmx_fio_do_gmx_bool(fio
, expand
->bInit_weights
);
330 gmx_fio_do_int(fio
, expand
->nstexpanded
);
331 gmx_fio_do_int(fio
, expand
->elmcmove
);
332 gmx_fio_do_int(fio
, expand
->elamstats
);
333 gmx_fio_do_int(fio
, expand
->lmc_repeats
);
334 gmx_fio_do_int(fio
, expand
->gibbsdeltalam
);
335 gmx_fio_do_int(fio
, expand
->lmc_forced_nstart
);
336 gmx_fio_do_int(fio
, expand
->lmc_seed
);
337 gmx_fio_do_real(fio
, expand
->mc_temp
);
338 gmx_fio_do_int(fio
, expand
->bSymmetrizedTMatrix
);
339 gmx_fio_do_int(fio
, expand
->nstTij
);
340 gmx_fio_do_int(fio
, expand
->minvarmin
);
341 gmx_fio_do_int(fio
, expand
->c_range
);
342 gmx_fio_do_real(fio
, expand
->wl_scale
);
343 gmx_fio_do_real(fio
, expand
->wl_ratio
);
344 gmx_fio_do_real(fio
, expand
->init_wl_delta
);
345 gmx_fio_do_gmx_bool(fio
, expand
->bWLoneovert
);
346 gmx_fio_do_int(fio
, expand
->elmceq
);
347 gmx_fio_do_int(fio
, expand
->equil_steps
);
348 gmx_fio_do_int(fio
, expand
->equil_samples
);
349 gmx_fio_do_int(fio
, expand
->equil_n_at_lam
);
350 gmx_fio_do_real(fio
, expand
->equil_wl_delta
);
351 gmx_fio_do_real(fio
, expand
->equil_ratio
);
355 static void do_simtempvals(t_fileio
*fio
, t_simtemp
*simtemp
, int n_lambda
, gmx_bool bRead
,
358 if (file_version
>= 79)
360 gmx_fio_do_int(fio
, simtemp
->eSimTempScale
);
361 gmx_fio_do_real(fio
, simtemp
->simtemp_high
);
362 gmx_fio_do_real(fio
, simtemp
->simtemp_low
);
367 snew(simtemp
->temperatures
, n_lambda
);
369 gmx_fio_ndo_real(fio
, simtemp
->temperatures
, n_lambda
);
374 static void do_imd(t_fileio
*fio
, t_IMD
*imd
, gmx_bool bRead
)
376 gmx_fio_do_int(fio
, imd
->nat
);
379 snew(imd
->ind
, imd
->nat
);
381 gmx_fio_ndo_int(fio
, imd
->ind
, imd
->nat
);
384 static void do_fepvals(t_fileio
*fio
, t_lambda
*fepvals
, gmx_bool bRead
, int file_version
)
386 /* i is defined in the ndo_double macro; use g to iterate. */
390 /* free energy values */
392 if (file_version
>= 79)
394 gmx_fio_do_int(fio
, fepvals
->init_fep_state
);
395 gmx_fio_do_double(fio
, fepvals
->init_lambda
);
396 gmx_fio_do_double(fio
, fepvals
->delta_lambda
);
398 else if (file_version
>= 59)
400 gmx_fio_do_double(fio
, fepvals
->init_lambda
);
401 gmx_fio_do_double(fio
, fepvals
->delta_lambda
);
405 gmx_fio_do_real(fio
, rdum
);
406 fepvals
->init_lambda
= rdum
;
407 gmx_fio_do_real(fio
, rdum
);
408 fepvals
->delta_lambda
= rdum
;
410 if (file_version
>= 79)
412 gmx_fio_do_int(fio
, fepvals
->n_lambda
);
415 snew(fepvals
->all_lambda
, efptNR
);
417 for (g
= 0; g
< efptNR
; g
++)
419 if (fepvals
->n_lambda
> 0)
423 snew(fepvals
->all_lambda
[g
], fepvals
->n_lambda
);
425 gmx_fio_ndo_double(fio
, fepvals
->all_lambda
[g
], fepvals
->n_lambda
);
426 gmx_fio_ndo_int(fio
, fepvals
->separate_dvdl
, efptNR
);
428 else if (fepvals
->init_lambda
>= 0)
430 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
434 else if (file_version
>= 64)
436 gmx_fio_do_int(fio
, fepvals
->n_lambda
);
441 snew(fepvals
->all_lambda
, efptNR
);
442 /* still allocate the all_lambda array's contents. */
443 for (g
= 0; g
< efptNR
; g
++)
445 if (fepvals
->n_lambda
> 0)
447 snew(fepvals
->all_lambda
[g
], fepvals
->n_lambda
);
451 gmx_fio_ndo_double(fio
, fepvals
->all_lambda
[efptFEP
],
453 if (fepvals
->init_lambda
>= 0)
457 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
461 /* copy the contents of the efptFEP lambda component to all
462 the other components */
463 for (g
= 0; g
< efptNR
; g
++)
465 for (h
= 0; h
< fepvals
->n_lambda
; h
++)
469 fepvals
->all_lambda
[g
][h
] =
470 fepvals
->all_lambda
[efptFEP
][h
];
479 fepvals
->n_lambda
= 0;
480 fepvals
->all_lambda
= NULL
;
481 if (fepvals
->init_lambda
>= 0)
483 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
486 if (file_version
>= 13)
488 gmx_fio_do_real(fio
, fepvals
->sc_alpha
);
492 fepvals
->sc_alpha
= 0;
494 if (file_version
>= 38)
496 gmx_fio_do_int(fio
, fepvals
->sc_power
);
500 fepvals
->sc_power
= 2;
502 if (file_version
>= 79)
504 gmx_fio_do_real(fio
, fepvals
->sc_r_power
);
508 fepvals
->sc_r_power
= 6.0;
510 if (file_version
>= 15)
512 gmx_fio_do_real(fio
, fepvals
->sc_sigma
);
516 fepvals
->sc_sigma
= 0.3;
520 if (file_version
>= 71)
522 fepvals
->sc_sigma_min
= fepvals
->sc_sigma
;
526 fepvals
->sc_sigma_min
= 0;
529 if (file_version
>= 79)
531 gmx_fio_do_int(fio
, fepvals
->bScCoul
);
535 fepvals
->bScCoul
= TRUE
;
537 if (file_version
>= 64)
539 gmx_fio_do_int(fio
, fepvals
->nstdhdl
);
543 fepvals
->nstdhdl
= 1;
546 if (file_version
>= 73)
548 gmx_fio_do_int(fio
, fepvals
->separate_dhdl_file
);
549 gmx_fio_do_int(fio
, fepvals
->dhdl_derivatives
);
553 fepvals
->separate_dhdl_file
= esepdhdlfileYES
;
554 fepvals
->dhdl_derivatives
= edhdlderivativesYES
;
556 if (file_version
>= 71)
558 gmx_fio_do_int(fio
, fepvals
->dh_hist_size
);
559 gmx_fio_do_double(fio
, fepvals
->dh_hist_spacing
);
563 fepvals
->dh_hist_size
= 0;
564 fepvals
->dh_hist_spacing
= 0.1;
566 if (file_version
>= 79)
568 gmx_fio_do_int(fio
, fepvals
->edHdLPrintEnergy
);
572 fepvals
->edHdLPrintEnergy
= edHdLPrintEnergyNO
;
575 /* handle lambda_neighbors */
576 if ((file_version
>= 83 && file_version
< 90) || file_version
>= 92)
578 gmx_fio_do_int(fio
, fepvals
->lambda_neighbors
);
579 if ( (fepvals
->lambda_neighbors
>= 0) && (fepvals
->init_fep_state
>= 0) &&
580 (fepvals
->init_lambda
< 0) )
582 fepvals
->lambda_start_n
= (fepvals
->init_fep_state
-
583 fepvals
->lambda_neighbors
);
584 fepvals
->lambda_stop_n
= (fepvals
->init_fep_state
+
585 fepvals
->lambda_neighbors
+ 1);
586 if (fepvals
->lambda_start_n
< 0)
588 fepvals
->lambda_start_n
= 0;;
590 if (fepvals
->lambda_stop_n
>= fepvals
->n_lambda
)
592 fepvals
->lambda_stop_n
= fepvals
->n_lambda
;
597 fepvals
->lambda_start_n
= 0;
598 fepvals
->lambda_stop_n
= fepvals
->n_lambda
;
603 fepvals
->lambda_start_n
= 0;
604 fepvals
->lambda_stop_n
= fepvals
->n_lambda
;
608 static void do_pull(t_fileio
*fio
, pull_params_t
*pull
, gmx_bool bRead
,
609 int file_version
, int ePullOld
)
615 if (file_version
>= 95)
617 gmx_fio_do_int(fio
, pull
->ngroup
);
619 gmx_fio_do_int(fio
, pull
->ncoord
);
620 if (file_version
< 95)
622 pull
->ngroup
= pull
->ncoord
+ 1;
624 if (file_version
< tpxv_PullCoordTypeGeom
)
628 gmx_fio_do_int(fio
, eGeomOld
);
629 gmx_fio_do_ivec(fio
, dimOld
);
630 /* The inner cylinder radius, now removed */
631 gmx_fio_do_real(fio
, dum
);
633 gmx_fio_do_real(fio
, pull
->cylinder_r
);
634 gmx_fio_do_real(fio
, pull
->constr_tol
);
635 if (file_version
>= 95)
637 gmx_fio_do_int(fio
, pull
->bPrintCOM1
);
638 /* With file_version < 95 this value is set below */
640 if (file_version
>= tpxv_PullCoordTypeGeom
)
642 gmx_fio_do_int(fio
, pull
->bPrintCOM2
);
643 gmx_fio_do_int(fio
, pull
->bPrintRefValue
);
644 gmx_fio_do_int(fio
, pull
->bPrintComp
);
648 pull
->bPrintCOM2
= FALSE
;
649 pull
->bPrintRefValue
= FALSE
;
650 pull
->bPrintComp
= TRUE
;
652 gmx_fio_do_int(fio
, pull
->nstxout
);
653 gmx_fio_do_int(fio
, pull
->nstfout
);
656 snew(pull
->group
, pull
->ngroup
);
657 snew(pull
->coord
, pull
->ncoord
);
659 if (file_version
< 95)
661 /* epullgPOS for position pulling, before epullgDIRPBC was removed */
662 if (eGeomOld
== epullgDIRPBC
)
664 gmx_fatal(FARGS
, "pull-geometry=position is no longer supported");
666 if (eGeomOld
> epullgDIRPBC
)
671 for (g
= 0; g
< pull
->ngroup
; g
++)
673 /* We read and ignore a pull coordinate for group 0 */
674 do_pullgrp_tpx_pre95(fio
, &pull
->group
[g
], &pull
->coord
[std::max(g
-1, 0)],
675 bRead
, file_version
);
678 pull
->coord
[g
-1].group
[0] = 0;
679 pull
->coord
[g
-1].group
[1] = g
;
683 pull
->bPrintCOM1
= (pull
->group
[0].nat
> 0);
687 for (g
= 0; g
< pull
->ngroup
; g
++)
689 do_pull_group(fio
, &pull
->group
[g
], bRead
);
691 for (g
= 0; g
< pull
->ncoord
; g
++)
693 do_pull_coord(fio
, &pull
->coord
[g
],
694 file_version
, ePullOld
, eGeomOld
, dimOld
);
700 static void do_rotgrp(t_fileio
*fio
, t_rotgrp
*rotg
, gmx_bool bRead
)
702 gmx_fio_do_int(fio
, rotg
->eType
);
703 gmx_fio_do_int(fio
, rotg
->bMassW
);
704 gmx_fio_do_int(fio
, rotg
->nat
);
707 snew(rotg
->ind
, rotg
->nat
);
709 gmx_fio_ndo_int(fio
, rotg
->ind
, rotg
->nat
);
712 snew(rotg
->x_ref
, rotg
->nat
);
714 gmx_fio_ndo_rvec(fio
, rotg
->x_ref
, rotg
->nat
);
715 gmx_fio_do_rvec(fio
, rotg
->vec
);
716 gmx_fio_do_rvec(fio
, rotg
->pivot
);
717 gmx_fio_do_real(fio
, rotg
->rate
);
718 gmx_fio_do_real(fio
, rotg
->k
);
719 gmx_fio_do_real(fio
, rotg
->slab_dist
);
720 gmx_fio_do_real(fio
, rotg
->min_gaussian
);
721 gmx_fio_do_real(fio
, rotg
->eps
);
722 gmx_fio_do_int(fio
, rotg
->eFittype
);
723 gmx_fio_do_int(fio
, rotg
->PotAngle_nstep
);
724 gmx_fio_do_real(fio
, rotg
->PotAngle_step
);
727 static void do_rot(t_fileio
*fio
, t_rot
*rot
, gmx_bool bRead
)
731 gmx_fio_do_int(fio
, rot
->ngrp
);
732 gmx_fio_do_int(fio
, rot
->nstrout
);
733 gmx_fio_do_int(fio
, rot
->nstsout
);
736 snew(rot
->grp
, rot
->ngrp
);
738 for (g
= 0; g
< rot
->ngrp
; g
++)
740 do_rotgrp(fio
, &rot
->grp
[g
], bRead
);
745 static void do_swapcoords(t_fileio
*fio
, t_swapcoords
*swap
, gmx_bool bRead
, int file_version
)
749 gmx_fio_do_int(fio
, swap
->nat
);
750 gmx_fio_do_int(fio
, swap
->nat_sol
);
751 for (j
= 0; j
< 2; j
++)
753 gmx_fio_do_int(fio
, swap
->nat_split
[j
]);
754 gmx_fio_do_int(fio
, swap
->massw_split
[j
]);
756 gmx_fio_do_int(fio
, swap
->nstswap
);
757 gmx_fio_do_int(fio
, swap
->nAverage
);
758 gmx_fio_do_real(fio
, swap
->threshold
);
759 gmx_fio_do_real(fio
, swap
->cyl0r
);
760 gmx_fio_do_real(fio
, swap
->cyl0u
);
761 gmx_fio_do_real(fio
, swap
->cyl0l
);
762 gmx_fio_do_real(fio
, swap
->cyl1r
);
763 gmx_fio_do_real(fio
, swap
->cyl1u
);
764 gmx_fio_do_real(fio
, swap
->cyl1l
);
768 snew(swap
->ind
, swap
->nat
);
769 snew(swap
->ind_sol
, swap
->nat_sol
);
770 for (j
= 0; j
< 2; j
++)
772 snew(swap
->ind_split
[j
], swap
->nat_split
[j
]);
776 gmx_fio_ndo_int(fio
, swap
->ind
, swap
->nat
);
777 gmx_fio_ndo_int(fio
, swap
->ind_sol
, swap
->nat_sol
);
778 for (j
= 0; j
< 2; j
++)
780 gmx_fio_ndo_int(fio
, swap
->ind_split
[j
], swap
->nat_split
[j
]);
783 for (j
= 0; j
< eCompNR
; j
++)
785 gmx_fio_do_int(fio
, swap
->nanions
[j
]);
786 gmx_fio_do_int(fio
, swap
->ncations
[j
]);
789 if (file_version
>= tpxv_CompElWithSwapLayerOffset
)
791 for (j
= 0; j
< eCompNR
; j
++)
793 gmx_fio_do_real(fio
, swap
->bulkOffset
[j
]);
800 static void do_inputrec(t_fileio
*fio
, t_inputrec
*ir
, gmx_bool bRead
,
801 int file_version
, real
*fudgeQQ
)
803 int i
, j
, k
, *tmp
, idum
= 0;
806 gmx_bool bSimAnn
, bdum
= 0;
807 real zerotemptime
, finish_t
, init_temp
, finish_temp
;
809 if (file_version
!= tpx_version
)
811 /* Give a warning about features that are not accessible */
812 fprintf(stderr
, "Note: file tpx version %d, software tpx version %d\n",
813 file_version
, tpx_version
);
821 if (file_version
== 0)
826 /* Basic inputrec stuff */
827 gmx_fio_do_int(fio
, ir
->eI
);
828 if (file_version
>= 62)
830 gmx_fio_do_int64(fio
, ir
->nsteps
);
834 gmx_fio_do_int(fio
, idum
);
838 if (file_version
> 25)
840 if (file_version
>= 62)
842 gmx_fio_do_int64(fio
, ir
->init_step
);
846 gmx_fio_do_int(fio
, idum
);
847 ir
->init_step
= idum
;
855 if (file_version
>= 58)
857 gmx_fio_do_int(fio
, ir
->simulation_part
);
861 ir
->simulation_part
= 1;
864 if (file_version
>= 67)
866 gmx_fio_do_int(fio
, ir
->nstcalcenergy
);
870 ir
->nstcalcenergy
= 1;
872 if (file_version
< 53)
874 /* The pbc info has been moved out of do_inputrec,
875 * since we always want it, also without reading the inputrec.
877 gmx_fio_do_int(fio
, ir
->ePBC
);
878 if ((file_version
<= 15) && (ir
->ePBC
== 2))
882 if (file_version
>= 45)
884 gmx_fio_do_int(fio
, ir
->bPeriodicMols
);
891 ir
->bPeriodicMols
= TRUE
;
895 ir
->bPeriodicMols
= FALSE
;
899 if (file_version
>= 81)
901 gmx_fio_do_int(fio
, ir
->cutoff_scheme
);
902 if (file_version
< 94)
904 ir
->cutoff_scheme
= 1 - ir
->cutoff_scheme
;
909 ir
->cutoff_scheme
= ecutsGROUP
;
911 gmx_fio_do_int(fio
, ir
->ns_type
);
912 gmx_fio_do_int(fio
, ir
->nstlist
);
913 gmx_fio_do_int(fio
, idum
); /* used to be ndelta; not used anymore */
914 if (file_version
< 41)
916 gmx_fio_do_int(fio
, idum
);
917 gmx_fio_do_int(fio
, idum
);
919 if (file_version
>= 45)
921 gmx_fio_do_real(fio
, ir
->rtpi
);
927 gmx_fio_do_int(fio
, ir
->nstcomm
);
928 if (file_version
> 34)
930 gmx_fio_do_int(fio
, ir
->comm_mode
);
932 else if (ir
->nstcomm
< 0)
934 ir
->comm_mode
= ecmANGULAR
;
938 ir
->comm_mode
= ecmLINEAR
;
940 ir
->nstcomm
= abs(ir
->nstcomm
);
942 /* ignore nstcheckpoint */
943 if (file_version
> 25 && file_version
< tpxv_RemoveObsoleteParameters1
)
945 gmx_fio_do_int(fio
, idum
);
948 gmx_fio_do_int(fio
, ir
->nstcgsteep
);
950 if (file_version
>= 30)
952 gmx_fio_do_int(fio
, ir
->nbfgscorr
);
959 gmx_fio_do_int(fio
, ir
->nstlog
);
960 gmx_fio_do_int(fio
, ir
->nstxout
);
961 gmx_fio_do_int(fio
, ir
->nstvout
);
962 gmx_fio_do_int(fio
, ir
->nstfout
);
963 gmx_fio_do_int(fio
, ir
->nstenergy
);
964 gmx_fio_do_int(fio
, ir
->nstxout_compressed
);
965 if (file_version
>= 59)
967 gmx_fio_do_double(fio
, ir
->init_t
);
968 gmx_fio_do_double(fio
, ir
->delta_t
);
972 gmx_fio_do_real(fio
, rdum
);
974 gmx_fio_do_real(fio
, rdum
);
977 gmx_fio_do_real(fio
, ir
->x_compression_precision
);
978 if (file_version
< 19)
980 gmx_fio_do_int(fio
, idum
);
981 gmx_fio_do_int(fio
, idum
);
983 if (file_version
< 18)
985 gmx_fio_do_int(fio
, idum
);
987 if (file_version
>= 81)
989 gmx_fio_do_real(fio
, ir
->verletbuf_tol
);
993 ir
->verletbuf_tol
= 0;
995 gmx_fio_do_real(fio
, ir
->rlist
);
996 if (file_version
>= 67)
998 gmx_fio_do_real(fio
, ir
->rlistlong
);
1000 if (file_version
>= 82 && file_version
!= 90)
1002 gmx_fio_do_int(fio
, ir
->nstcalclr
);
1006 /* Calculate at NS steps */
1007 ir
->nstcalclr
= ir
->nstlist
;
1009 gmx_fio_do_int(fio
, ir
->coulombtype
);
1010 if (file_version
< 32 && ir
->coulombtype
== eelRF
)
1012 ir
->coulombtype
= eelRF_NEC
;
1014 if (file_version
>= 81)
1016 gmx_fio_do_int(fio
, ir
->coulomb_modifier
);
1020 ir
->coulomb_modifier
= (ir
->cutoff_scheme
== ecutsVERLET
? eintmodPOTSHIFT
: eintmodNONE
);
1022 gmx_fio_do_real(fio
, ir
->rcoulomb_switch
);
1023 gmx_fio_do_real(fio
, ir
->rcoulomb
);
1024 gmx_fio_do_int(fio
, ir
->vdwtype
);
1025 if (file_version
>= 81)
1027 gmx_fio_do_int(fio
, ir
->vdw_modifier
);
1031 ir
->vdw_modifier
= (ir
->cutoff_scheme
== ecutsVERLET
? eintmodPOTSHIFT
: eintmodNONE
);
1033 gmx_fio_do_real(fio
, ir
->rvdw_switch
);
1034 gmx_fio_do_real(fio
, ir
->rvdw
);
1035 if (file_version
< 67)
1037 ir
->rlistlong
= max_cutoff(ir
->rlist
, max_cutoff(ir
->rvdw
, ir
->rcoulomb
));
1039 gmx_fio_do_int(fio
, ir
->eDispCorr
);
1040 gmx_fio_do_real(fio
, ir
->epsilon_r
);
1041 if (file_version
>= 37)
1043 gmx_fio_do_real(fio
, ir
->epsilon_rf
);
1047 if (EEL_RF(ir
->coulombtype
))
1049 ir
->epsilon_rf
= ir
->epsilon_r
;
1050 ir
->epsilon_r
= 1.0;
1054 ir
->epsilon_rf
= 1.0;
1057 if (file_version
>= 29)
1059 gmx_fio_do_real(fio
, ir
->tabext
);
1066 if (file_version
> 25)
1068 gmx_fio_do_int(fio
, ir
->gb_algorithm
);
1069 gmx_fio_do_int(fio
, ir
->nstgbradii
);
1070 gmx_fio_do_real(fio
, ir
->rgbradii
);
1071 gmx_fio_do_real(fio
, ir
->gb_saltconc
);
1072 gmx_fio_do_int(fio
, ir
->implicit_solvent
);
1076 ir
->gb_algorithm
= egbSTILL
;
1079 ir
->gb_saltconc
= 0;
1080 ir
->implicit_solvent
= eisNO
;
1082 if (file_version
>= 55)
1084 gmx_fio_do_real(fio
, ir
->gb_epsilon_solvent
);
1085 gmx_fio_do_real(fio
, ir
->gb_obc_alpha
);
1086 gmx_fio_do_real(fio
, ir
->gb_obc_beta
);
1087 gmx_fio_do_real(fio
, ir
->gb_obc_gamma
);
1088 if (file_version
>= 60)
1090 gmx_fio_do_real(fio
, ir
->gb_dielectric_offset
);
1091 gmx_fio_do_int(fio
, ir
->sa_algorithm
);
1095 ir
->gb_dielectric_offset
= 0.009;
1096 ir
->sa_algorithm
= esaAPPROX
;
1098 gmx_fio_do_real(fio
, ir
->sa_surface_tension
);
1100 /* Override sa_surface_tension if it is not changed in the mpd-file */
1101 if (ir
->sa_surface_tension
< 0)
1103 if (ir
->gb_algorithm
== egbSTILL
)
1105 ir
->sa_surface_tension
= 0.0049 * 100 * CAL2JOULE
;
1107 else if (ir
->gb_algorithm
== egbHCT
|| ir
->gb_algorithm
== egbOBC
)
1109 ir
->sa_surface_tension
= 0.0054 * 100 * CAL2JOULE
;
1116 /* Better use sensible values than insane (0.0) ones... */
1117 ir
->gb_epsilon_solvent
= 80;
1118 ir
->gb_obc_alpha
= 1.0;
1119 ir
->gb_obc_beta
= 0.8;
1120 ir
->gb_obc_gamma
= 4.85;
1121 ir
->sa_surface_tension
= 2.092;
1125 if (file_version
>= 81)
1127 gmx_fio_do_real(fio
, ir
->fourier_spacing
);
1131 ir
->fourier_spacing
= 0.0;
1133 gmx_fio_do_int(fio
, ir
->nkx
);
1134 gmx_fio_do_int(fio
, ir
->nky
);
1135 gmx_fio_do_int(fio
, ir
->nkz
);
1136 gmx_fio_do_int(fio
, ir
->pme_order
);
1137 gmx_fio_do_real(fio
, ir
->ewald_rtol
);
1139 if (file_version
>= 93)
1141 gmx_fio_do_real(fio
, ir
->ewald_rtol_lj
);
1145 ir
->ewald_rtol_lj
= ir
->ewald_rtol
;
1148 if (file_version
>= 24)
1150 gmx_fio_do_int(fio
, ir
->ewald_geometry
);
1154 ir
->ewald_geometry
= eewg3D
;
1157 if (file_version
<= 17)
1159 ir
->epsilon_surface
= 0;
1160 if (file_version
== 17)
1162 gmx_fio_do_int(fio
, idum
);
1167 gmx_fio_do_real(fio
, ir
->epsilon_surface
);
1170 /* ignore bOptFFT */
1171 if (file_version
< tpxv_RemoveObsoleteParameters1
)
1173 gmx_fio_do_gmx_bool(fio
, bdum
);
1176 if (file_version
>= 93)
1178 gmx_fio_do_int(fio
, ir
->ljpme_combination_rule
);
1180 gmx_fio_do_gmx_bool(fio
, ir
->bContinuation
);
1181 gmx_fio_do_int(fio
, ir
->etc
);
1182 /* before version 18, ir->etc was a gmx_bool (ir->btc),
1183 * but the values 0 and 1 still mean no and
1184 * berendsen temperature coupling, respectively.
1186 if (file_version
>= 79)
1188 gmx_fio_do_gmx_bool(fio
, ir
->bPrintNHChains
);
1190 if (file_version
>= 71)
1192 gmx_fio_do_int(fio
, ir
->nsttcouple
);
1196 ir
->nsttcouple
= ir
->nstcalcenergy
;
1198 if (file_version
<= 15)
1200 gmx_fio_do_int(fio
, idum
);
1202 if (file_version
<= 17)
1204 gmx_fio_do_int(fio
, ir
->epct
);
1205 if (file_version
<= 15)
1209 ir
->epct
= epctSURFACETENSION
;
1211 gmx_fio_do_int(fio
, idum
);
1214 /* we have removed the NO alternative at the beginning */
1218 ir
->epct
= epctISOTROPIC
;
1222 ir
->epc
= epcBERENDSEN
;
1227 gmx_fio_do_int(fio
, ir
->epc
);
1228 gmx_fio_do_int(fio
, ir
->epct
);
1230 if (file_version
>= 71)
1232 gmx_fio_do_int(fio
, ir
->nstpcouple
);
1236 ir
->nstpcouple
= ir
->nstcalcenergy
;
1238 gmx_fio_do_real(fio
, ir
->tau_p
);
1239 if (file_version
<= 15)
1241 gmx_fio_do_rvec(fio
, vdum
);
1242 clear_mat(ir
->ref_p
);
1243 for (i
= 0; i
< DIM
; i
++)
1245 ir
->ref_p
[i
][i
] = vdum
[i
];
1250 gmx_fio_do_rvec(fio
, ir
->ref_p
[XX
]);
1251 gmx_fio_do_rvec(fio
, ir
->ref_p
[YY
]);
1252 gmx_fio_do_rvec(fio
, ir
->ref_p
[ZZ
]);
1254 if (file_version
<= 15)
1256 gmx_fio_do_rvec(fio
, vdum
);
1257 clear_mat(ir
->compress
);
1258 for (i
= 0; i
< DIM
; i
++)
1260 ir
->compress
[i
][i
] = vdum
[i
];
1265 gmx_fio_do_rvec(fio
, ir
->compress
[XX
]);
1266 gmx_fio_do_rvec(fio
, ir
->compress
[YY
]);
1267 gmx_fio_do_rvec(fio
, ir
->compress
[ZZ
]);
1269 if (file_version
>= 47)
1271 gmx_fio_do_int(fio
, ir
->refcoord_scaling
);
1272 gmx_fio_do_rvec(fio
, ir
->posres_com
);
1273 gmx_fio_do_rvec(fio
, ir
->posres_comB
);
1277 ir
->refcoord_scaling
= erscNO
;
1278 clear_rvec(ir
->posres_com
);
1279 clear_rvec(ir
->posres_comB
);
1281 if ((file_version
> 25) && (file_version
< 79))
1283 gmx_fio_do_int(fio
, ir
->andersen_seed
);
1287 ir
->andersen_seed
= 0;
1289 if (file_version
< 26)
1291 gmx_fio_do_gmx_bool(fio
, bSimAnn
);
1292 gmx_fio_do_real(fio
, zerotemptime
);
1295 if (file_version
< 37)
1297 gmx_fio_do_real(fio
, rdum
);
1300 gmx_fio_do_real(fio
, ir
->shake_tol
);
1301 if (file_version
< 54)
1303 gmx_fio_do_real(fio
, *fudgeQQ
);
1306 gmx_fio_do_int(fio
, ir
->efep
);
1307 if (file_version
<= 14 && ir
->efep
!= efepNO
)
1311 do_fepvals(fio
, ir
->fepvals
, bRead
, file_version
);
1313 if (file_version
>= 79)
1315 gmx_fio_do_gmx_bool(fio
, ir
->bSimTemp
);
1318 ir
->bSimTemp
= TRUE
;
1323 ir
->bSimTemp
= FALSE
;
1327 do_simtempvals(fio
, ir
->simtempvals
, ir
->fepvals
->n_lambda
, bRead
, file_version
);
1330 if (file_version
>= 79)
1332 gmx_fio_do_gmx_bool(fio
, ir
->bExpanded
);
1335 ir
->bExpanded
= TRUE
;
1339 ir
->bExpanded
= FALSE
;
1344 do_expandedvals(fio
, ir
->expandedvals
, ir
->fepvals
, bRead
, file_version
);
1346 if (file_version
>= 57)
1348 gmx_fio_do_int(fio
, ir
->eDisre
);
1350 gmx_fio_do_int(fio
, ir
->eDisreWeighting
);
1351 if (file_version
< 22)
1353 if (ir
->eDisreWeighting
== 0)
1355 ir
->eDisreWeighting
= edrwEqual
;
1359 ir
->eDisreWeighting
= edrwConservative
;
1362 gmx_fio_do_gmx_bool(fio
, ir
->bDisreMixed
);
1363 gmx_fio_do_real(fio
, ir
->dr_fc
);
1364 gmx_fio_do_real(fio
, ir
->dr_tau
);
1365 gmx_fio_do_int(fio
, ir
->nstdisreout
);
1366 if (file_version
>= 22)
1368 gmx_fio_do_real(fio
, ir
->orires_fc
);
1369 gmx_fio_do_real(fio
, ir
->orires_tau
);
1370 gmx_fio_do_int(fio
, ir
->nstorireout
);
1376 ir
->nstorireout
= 0;
1379 /* ignore dihre_fc */
1380 if (file_version
>= 26 && file_version
< 79)
1382 gmx_fio_do_real(fio
, rdum
);
1383 if (file_version
< 56)
1385 gmx_fio_do_real(fio
, rdum
);
1386 gmx_fio_do_int(fio
, idum
);
1390 gmx_fio_do_real(fio
, ir
->em_stepsize
);
1391 gmx_fio_do_real(fio
, ir
->em_tol
);
1392 if (file_version
>= 22)
1394 gmx_fio_do_gmx_bool(fio
, ir
->bShakeSOR
);
1398 ir
->bShakeSOR
= TRUE
;
1400 if (file_version
>= 11)
1402 gmx_fio_do_int(fio
, ir
->niter
);
1407 fprintf(stderr
, "Note: niter not in run input file, setting it to %d\n",
1410 if (file_version
>= 21)
1412 gmx_fio_do_real(fio
, ir
->fc_stepsize
);
1416 ir
->fc_stepsize
= 0;
1418 gmx_fio_do_int(fio
, ir
->eConstrAlg
);
1419 gmx_fio_do_int(fio
, ir
->nProjOrder
);
1420 gmx_fio_do_real(fio
, ir
->LincsWarnAngle
);
1421 if (file_version
<= 14)
1423 gmx_fio_do_int(fio
, idum
);
1425 if (file_version
>= 26)
1427 gmx_fio_do_int(fio
, ir
->nLincsIter
);
1432 fprintf(stderr
, "Note: nLincsIter not in run input file, setting it to %d\n",
1435 if (file_version
< 33)
1437 gmx_fio_do_real(fio
, bd_temp
);
1439 gmx_fio_do_real(fio
, ir
->bd_fric
);
1440 if (file_version
>= tpxv_Use64BitRandomSeed
)
1442 gmx_fio_do_int64(fio
, ir
->ld_seed
);
1446 gmx_fio_do_int(fio
, idum
);
1449 if (file_version
>= 33)
1451 for (i
= 0; i
< DIM
; i
++)
1453 gmx_fio_do_rvec(fio
, ir
->deform
[i
]);
1458 for (i
= 0; i
< DIM
; i
++)
1460 clear_rvec(ir
->deform
[i
]);
1463 if (file_version
>= 14)
1465 gmx_fio_do_real(fio
, ir
->cos_accel
);
1471 gmx_fio_do_int(fio
, ir
->userint1
);
1472 gmx_fio_do_int(fio
, ir
->userint2
);
1473 gmx_fio_do_int(fio
, ir
->userint3
);
1474 gmx_fio_do_int(fio
, ir
->userint4
);
1475 gmx_fio_do_real(fio
, ir
->userreal1
);
1476 gmx_fio_do_real(fio
, ir
->userreal2
);
1477 gmx_fio_do_real(fio
, ir
->userreal3
);
1478 gmx_fio_do_real(fio
, ir
->userreal4
);
1481 if (file_version
>= 77)
1483 gmx_fio_do_gmx_bool(fio
, ir
->bAdress
);
1488 snew(ir
->adress
, 1);
1490 gmx_fio_do_int(fio
, ir
->adress
->type
);
1491 gmx_fio_do_real(fio
, ir
->adress
->const_wf
);
1492 gmx_fio_do_real(fio
, ir
->adress
->ex_width
);
1493 gmx_fio_do_real(fio
, ir
->adress
->hy_width
);
1494 gmx_fio_do_int(fio
, ir
->adress
->icor
);
1495 gmx_fio_do_int(fio
, ir
->adress
->site
);
1496 gmx_fio_do_rvec(fio
, ir
->adress
->refs
);
1497 gmx_fio_do_int(fio
, ir
->adress
->n_tf_grps
);
1498 gmx_fio_do_real(fio
, ir
->adress
->ex_forcecap
);
1499 gmx_fio_do_int(fio
, ir
->adress
->n_energy_grps
);
1500 gmx_fio_do_int(fio
, ir
->adress
->do_hybridpairs
);
1504 snew(ir
->adress
->tf_table_index
, ir
->adress
->n_tf_grps
);
1506 if (ir
->adress
->n_tf_grps
> 0)
1508 gmx_fio_ndo_int(fio
, ir
->adress
->tf_table_index
, ir
->adress
->n_tf_grps
);
1512 snew(ir
->adress
->group_explicit
, ir
->adress
->n_energy_grps
);
1514 if (ir
->adress
->n_energy_grps
> 0)
1516 gmx_fio_ndo_int(fio
, ir
->adress
->group_explicit
, ir
->adress
->n_energy_grps
);
1522 ir
->bAdress
= FALSE
;
1526 if (file_version
>= 48)
1530 if (file_version
>= tpxv_PullCoordTypeGeom
)
1532 gmx_fio_do_gmx_bool(fio
, ir
->bPull
);
1536 gmx_fio_do_int(fio
, ePullOld
);
1537 ir
->bPull
= (ePullOld
> 0);
1538 /* We removed the first ePull=ePullNo for the enum */
1547 do_pull(fio
, ir
->pull
, bRead
, file_version
, ePullOld
);
1555 /* Enforced rotation */
1556 if (file_version
>= 74)
1558 gmx_fio_do_int(fio
, ir
->bRot
);
1559 if (ir
->bRot
== TRUE
)
1565 do_rot(fio
, ir
->rot
, bRead
);
1573 /* Interactive molecular dynamics */
1574 if (file_version
>= tpxv_InteractiveMolecularDynamics
)
1576 gmx_fio_do_int(fio
, ir
->bIMD
);
1577 if (TRUE
== ir
->bIMD
)
1583 do_imd(fio
, ir
->imd
, bRead
);
1588 /* We don't support IMD sessions for old .tpr files */
1593 gmx_fio_do_int(fio
, ir
->opts
.ngtc
);
1594 if (file_version
>= 69)
1596 gmx_fio_do_int(fio
, ir
->opts
.nhchainlength
);
1600 ir
->opts
.nhchainlength
= 1;
1602 gmx_fio_do_int(fio
, ir
->opts
.ngacc
);
1603 gmx_fio_do_int(fio
, ir
->opts
.ngfrz
);
1604 gmx_fio_do_int(fio
, ir
->opts
.ngener
);
1608 snew(ir
->opts
.nrdf
, ir
->opts
.ngtc
);
1609 snew(ir
->opts
.ref_t
, ir
->opts
.ngtc
);
1610 snew(ir
->opts
.annealing
, ir
->opts
.ngtc
);
1611 snew(ir
->opts
.anneal_npoints
, ir
->opts
.ngtc
);
1612 snew(ir
->opts
.anneal_time
, ir
->opts
.ngtc
);
1613 snew(ir
->opts
.anneal_temp
, ir
->opts
.ngtc
);
1614 snew(ir
->opts
.tau_t
, ir
->opts
.ngtc
);
1615 snew(ir
->opts
.nFreeze
, ir
->opts
.ngfrz
);
1616 snew(ir
->opts
.acc
, ir
->opts
.ngacc
);
1617 snew(ir
->opts
.egp_flags
, ir
->opts
.ngener
*ir
->opts
.ngener
);
1619 if (ir
->opts
.ngtc
> 0)
1621 if (bRead
&& file_version
< 13)
1623 snew(tmp
, ir
->opts
.ngtc
);
1624 gmx_fio_ndo_int(fio
, tmp
, ir
->opts
.ngtc
);
1625 for (i
= 0; i
< ir
->opts
.ngtc
; i
++)
1627 ir
->opts
.nrdf
[i
] = tmp
[i
];
1633 gmx_fio_ndo_real(fio
, ir
->opts
.nrdf
, ir
->opts
.ngtc
);
1635 gmx_fio_ndo_real(fio
, ir
->opts
.ref_t
, ir
->opts
.ngtc
);
1636 gmx_fio_ndo_real(fio
, ir
->opts
.tau_t
, ir
->opts
.ngtc
);
1637 if (file_version
< 33 && ir
->eI
== eiBD
)
1639 for (i
= 0; i
< ir
->opts
.ngtc
; i
++)
1641 ir
->opts
.tau_t
[i
] = bd_temp
;
1645 if (ir
->opts
.ngfrz
> 0)
1647 gmx_fio_ndo_ivec(fio
, ir
->opts
.nFreeze
, ir
->opts
.ngfrz
);
1649 if (ir
->opts
.ngacc
> 0)
1651 gmx_fio_ndo_rvec(fio
, ir
->opts
.acc
, ir
->opts
.ngacc
);
1653 if (file_version
>= 12)
1655 gmx_fio_ndo_int(fio
, ir
->opts
.egp_flags
,
1656 ir
->opts
.ngener
*ir
->opts
.ngener
);
1659 if (bRead
&& file_version
< 26)
1661 for (i
= 0; i
< ir
->opts
.ngtc
; i
++)
1665 ir
->opts
.annealing
[i
] = eannSINGLE
;
1666 ir
->opts
.anneal_npoints
[i
] = 2;
1667 snew(ir
->opts
.anneal_time
[i
], 2);
1668 snew(ir
->opts
.anneal_temp
[i
], 2);
1669 /* calculate the starting/ending temperatures from reft, zerotemptime, and nsteps */
1670 finish_t
= ir
->init_t
+ ir
->nsteps
* ir
->delta_t
;
1671 init_temp
= ir
->opts
.ref_t
[i
]*(1-ir
->init_t
/zerotemptime
);
1672 finish_temp
= ir
->opts
.ref_t
[i
]*(1-finish_t
/zerotemptime
);
1673 ir
->opts
.anneal_time
[i
][0] = ir
->init_t
;
1674 ir
->opts
.anneal_time
[i
][1] = finish_t
;
1675 ir
->opts
.anneal_temp
[i
][0] = init_temp
;
1676 ir
->opts
.anneal_temp
[i
][1] = finish_temp
;
1680 ir
->opts
.annealing
[i
] = eannNO
;
1681 ir
->opts
.anneal_npoints
[i
] = 0;
1687 /* file version 26 or later */
1688 /* First read the lists with annealing and npoints for each group */
1689 gmx_fio_ndo_int(fio
, ir
->opts
.annealing
, ir
->opts
.ngtc
);
1690 gmx_fio_ndo_int(fio
, ir
->opts
.anneal_npoints
, ir
->opts
.ngtc
);
1691 for (j
= 0; j
< (ir
->opts
.ngtc
); j
++)
1693 k
= ir
->opts
.anneal_npoints
[j
];
1696 snew(ir
->opts
.anneal_time
[j
], k
);
1697 snew(ir
->opts
.anneal_temp
[j
], k
);
1699 gmx_fio_ndo_real(fio
, ir
->opts
.anneal_time
[j
], k
);
1700 gmx_fio_ndo_real(fio
, ir
->opts
.anneal_temp
[j
], k
);
1704 if (file_version
>= 45)
1706 gmx_fio_do_int(fio
, ir
->nwall
);
1707 gmx_fio_do_int(fio
, ir
->wall_type
);
1708 if (file_version
>= 50)
1710 gmx_fio_do_real(fio
, ir
->wall_r_linpot
);
1714 ir
->wall_r_linpot
= -1;
1716 gmx_fio_do_int(fio
, ir
->wall_atomtype
[0]);
1717 gmx_fio_do_int(fio
, ir
->wall_atomtype
[1]);
1718 gmx_fio_do_real(fio
, ir
->wall_density
[0]);
1719 gmx_fio_do_real(fio
, ir
->wall_density
[1]);
1720 gmx_fio_do_real(fio
, ir
->wall_ewald_zfac
);
1726 ir
->wall_atomtype
[0] = -1;
1727 ir
->wall_atomtype
[1] = -1;
1728 ir
->wall_density
[0] = 0;
1729 ir
->wall_density
[1] = 0;
1730 ir
->wall_ewald_zfac
= 3;
1732 /* Cosine stuff for electric fields */
1733 for (j
= 0; (j
< DIM
); j
++)
1735 gmx_fio_do_int(fio
, ir
->ex
[j
].n
);
1736 gmx_fio_do_int(fio
, ir
->et
[j
].n
);
1739 snew(ir
->ex
[j
].a
, ir
->ex
[j
].n
);
1740 snew(ir
->ex
[j
].phi
, ir
->ex
[j
].n
);
1741 snew(ir
->et
[j
].a
, ir
->et
[j
].n
);
1742 snew(ir
->et
[j
].phi
, ir
->et
[j
].n
);
1744 gmx_fio_ndo_real(fio
, ir
->ex
[j
].a
, ir
->ex
[j
].n
);
1745 gmx_fio_ndo_real(fio
, ir
->ex
[j
].phi
, ir
->ex
[j
].n
);
1746 gmx_fio_ndo_real(fio
, ir
->et
[j
].a
, ir
->et
[j
].n
);
1747 gmx_fio_ndo_real(fio
, ir
->et
[j
].phi
, ir
->et
[j
].n
);
1751 if (file_version
>= tpxv_ComputationalElectrophysiology
)
1753 gmx_fio_do_int(fio
, ir
->eSwapCoords
);
1754 if (ir
->eSwapCoords
!= eswapNO
)
1760 do_swapcoords(fio
, ir
->swap
, bRead
, file_version
);
1765 if (file_version
>= 39)
1767 gmx_fio_do_gmx_bool(fio
, ir
->bQMMM
);
1768 gmx_fio_do_int(fio
, ir
->QMMMscheme
);
1769 gmx_fio_do_real(fio
, ir
->scalefactor
);
1770 gmx_fio_do_int(fio
, ir
->opts
.ngQM
);
1773 snew(ir
->opts
.QMmethod
, ir
->opts
.ngQM
);
1774 snew(ir
->opts
.QMbasis
, ir
->opts
.ngQM
);
1775 snew(ir
->opts
.QMcharge
, ir
->opts
.ngQM
);
1776 snew(ir
->opts
.QMmult
, ir
->opts
.ngQM
);
1777 snew(ir
->opts
.bSH
, ir
->opts
.ngQM
);
1778 snew(ir
->opts
.CASorbitals
, ir
->opts
.ngQM
);
1779 snew(ir
->opts
.CASelectrons
, ir
->opts
.ngQM
);
1780 snew(ir
->opts
.SAon
, ir
->opts
.ngQM
);
1781 snew(ir
->opts
.SAoff
, ir
->opts
.ngQM
);
1782 snew(ir
->opts
.SAsteps
, ir
->opts
.ngQM
);
1783 snew(ir
->opts
.bOPT
, ir
->opts
.ngQM
);
1784 snew(ir
->opts
.bTS
, ir
->opts
.ngQM
);
1786 if (ir
->opts
.ngQM
> 0)
1788 gmx_fio_ndo_int(fio
, ir
->opts
.QMmethod
, ir
->opts
.ngQM
);
1789 gmx_fio_ndo_int(fio
, ir
->opts
.QMbasis
, ir
->opts
.ngQM
);
1790 gmx_fio_ndo_int(fio
, ir
->opts
.QMcharge
, ir
->opts
.ngQM
);
1791 gmx_fio_ndo_int(fio
, ir
->opts
.QMmult
, ir
->opts
.ngQM
);
1792 gmx_fio_ndo_gmx_bool(fio
, ir
->opts
.bSH
, ir
->opts
.ngQM
);
1793 gmx_fio_ndo_int(fio
, ir
->opts
.CASorbitals
, ir
->opts
.ngQM
);
1794 gmx_fio_ndo_int(fio
, ir
->opts
.CASelectrons
, ir
->opts
.ngQM
);
1795 gmx_fio_ndo_real(fio
, ir
->opts
.SAon
, ir
->opts
.ngQM
);
1796 gmx_fio_ndo_real(fio
, ir
->opts
.SAoff
, ir
->opts
.ngQM
);
1797 gmx_fio_ndo_int(fio
, ir
->opts
.SAsteps
, ir
->opts
.ngQM
);
1798 gmx_fio_ndo_gmx_bool(fio
, ir
->opts
.bOPT
, ir
->opts
.ngQM
);
1799 gmx_fio_ndo_gmx_bool(fio
, ir
->opts
.bTS
, ir
->opts
.ngQM
);
1801 /* end of QMMM stuff */
1806 static void do_harm(t_fileio
*fio
, t_iparams
*iparams
)
1808 gmx_fio_do_real(fio
, iparams
->harmonic
.rA
);
1809 gmx_fio_do_real(fio
, iparams
->harmonic
.krA
);
1810 gmx_fio_do_real(fio
, iparams
->harmonic
.rB
);
1811 gmx_fio_do_real(fio
, iparams
->harmonic
.krB
);
1814 void do_iparams(t_fileio
*fio
, t_functype ftype
, t_iparams
*iparams
,
1815 gmx_bool bRead
, int file_version
)
1828 do_harm(fio
, iparams
);
1829 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && bRead
)
1831 /* Correct incorrect storage of parameters */
1832 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
1833 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
1837 gmx_fio_do_real(fio
, iparams
->harmonic
.rA
);
1838 gmx_fio_do_real(fio
, iparams
->harmonic
.krA
);
1840 case F_LINEAR_ANGLES
:
1841 gmx_fio_do_real(fio
, iparams
->linangle
.klinA
);
1842 gmx_fio_do_real(fio
, iparams
->linangle
.aA
);
1843 gmx_fio_do_real(fio
, iparams
->linangle
.klinB
);
1844 gmx_fio_do_real(fio
, iparams
->linangle
.aB
);
1847 gmx_fio_do_real(fio
, iparams
->fene
.bm
);
1848 gmx_fio_do_real(fio
, iparams
->fene
.kb
);
1852 gmx_fio_do_real(fio
, iparams
->restraint
.lowA
);
1853 gmx_fio_do_real(fio
, iparams
->restraint
.up1A
);
1854 gmx_fio_do_real(fio
, iparams
->restraint
.up2A
);
1855 gmx_fio_do_real(fio
, iparams
->restraint
.kA
);
1856 gmx_fio_do_real(fio
, iparams
->restraint
.lowB
);
1857 gmx_fio_do_real(fio
, iparams
->restraint
.up1B
);
1858 gmx_fio_do_real(fio
, iparams
->restraint
.up2B
);
1859 gmx_fio_do_real(fio
, iparams
->restraint
.kB
);
1865 gmx_fio_do_real(fio
, iparams
->tab
.kA
);
1866 gmx_fio_do_int(fio
, iparams
->tab
.table
);
1867 gmx_fio_do_real(fio
, iparams
->tab
.kB
);
1869 case F_CROSS_BOND_BONDS
:
1870 gmx_fio_do_real(fio
, iparams
->cross_bb
.r1e
);
1871 gmx_fio_do_real(fio
, iparams
->cross_bb
.r2e
);
1872 gmx_fio_do_real(fio
, iparams
->cross_bb
.krr
);
1874 case F_CROSS_BOND_ANGLES
:
1875 gmx_fio_do_real(fio
, iparams
->cross_ba
.r1e
);
1876 gmx_fio_do_real(fio
, iparams
->cross_ba
.r2e
);
1877 gmx_fio_do_real(fio
, iparams
->cross_ba
.r3e
);
1878 gmx_fio_do_real(fio
, iparams
->cross_ba
.krt
);
1880 case F_UREY_BRADLEY
:
1881 gmx_fio_do_real(fio
, iparams
->u_b
.thetaA
);
1882 gmx_fio_do_real(fio
, iparams
->u_b
.kthetaA
);
1883 gmx_fio_do_real(fio
, iparams
->u_b
.r13A
);
1884 gmx_fio_do_real(fio
, iparams
->u_b
.kUBA
);
1885 if (file_version
>= 79)
1887 gmx_fio_do_real(fio
, iparams
->u_b
.thetaB
);
1888 gmx_fio_do_real(fio
, iparams
->u_b
.kthetaB
);
1889 gmx_fio_do_real(fio
, iparams
->u_b
.r13B
);
1890 gmx_fio_do_real(fio
, iparams
->u_b
.kUBB
);
1894 iparams
->u_b
.thetaB
= iparams
->u_b
.thetaA
;
1895 iparams
->u_b
.kthetaB
= iparams
->u_b
.kthetaA
;
1896 iparams
->u_b
.r13B
= iparams
->u_b
.r13A
;
1897 iparams
->u_b
.kUBB
= iparams
->u_b
.kUBA
;
1900 case F_QUARTIC_ANGLES
:
1901 gmx_fio_do_real(fio
, iparams
->qangle
.theta
);
1902 gmx_fio_ndo_real(fio
, iparams
->qangle
.c
, 5);
1905 gmx_fio_do_real(fio
, iparams
->bham
.a
);
1906 gmx_fio_do_real(fio
, iparams
->bham
.b
);
1907 gmx_fio_do_real(fio
, iparams
->bham
.c
);
1910 gmx_fio_do_real(fio
, iparams
->morse
.b0A
);
1911 gmx_fio_do_real(fio
, iparams
->morse
.cbA
);
1912 gmx_fio_do_real(fio
, iparams
->morse
.betaA
);
1913 if (file_version
>= 79)
1915 gmx_fio_do_real(fio
, iparams
->morse
.b0B
);
1916 gmx_fio_do_real(fio
, iparams
->morse
.cbB
);
1917 gmx_fio_do_real(fio
, iparams
->morse
.betaB
);
1921 iparams
->morse
.b0B
= iparams
->morse
.b0A
;
1922 iparams
->morse
.cbB
= iparams
->morse
.cbA
;
1923 iparams
->morse
.betaB
= iparams
->morse
.betaA
;
1927 gmx_fio_do_real(fio
, iparams
->cubic
.b0
);
1928 gmx_fio_do_real(fio
, iparams
->cubic
.kb
);
1929 gmx_fio_do_real(fio
, iparams
->cubic
.kcub
);
1933 case F_POLARIZATION
:
1934 gmx_fio_do_real(fio
, iparams
->polarize
.alpha
);
1937 gmx_fio_do_real(fio
, iparams
->anharm_polarize
.alpha
);
1938 gmx_fio_do_real(fio
, iparams
->anharm_polarize
.drcut
);
1939 gmx_fio_do_real(fio
, iparams
->anharm_polarize
.khyp
);
1942 if (file_version
< 31)
1944 gmx_fatal(FARGS
, "Old tpr files with water_polarization not supported. Make a new.");
1946 gmx_fio_do_real(fio
, iparams
->wpol
.al_x
);
1947 gmx_fio_do_real(fio
, iparams
->wpol
.al_y
);
1948 gmx_fio_do_real(fio
, iparams
->wpol
.al_z
);
1949 gmx_fio_do_real(fio
, iparams
->wpol
.rOH
);
1950 gmx_fio_do_real(fio
, iparams
->wpol
.rHH
);
1951 gmx_fio_do_real(fio
, iparams
->wpol
.rOD
);
1954 gmx_fio_do_real(fio
, iparams
->thole
.a
);
1955 gmx_fio_do_real(fio
, iparams
->thole
.alpha1
);
1956 gmx_fio_do_real(fio
, iparams
->thole
.alpha2
);
1957 gmx_fio_do_real(fio
, iparams
->thole
.rfac
);
1960 gmx_fio_do_real(fio
, iparams
->lj
.c6
);
1961 gmx_fio_do_real(fio
, iparams
->lj
.c12
);
1964 gmx_fio_do_real(fio
, iparams
->lj14
.c6A
);
1965 gmx_fio_do_real(fio
, iparams
->lj14
.c12A
);
1966 gmx_fio_do_real(fio
, iparams
->lj14
.c6B
);
1967 gmx_fio_do_real(fio
, iparams
->lj14
.c12B
);
1970 gmx_fio_do_real(fio
, iparams
->ljc14
.fqq
);
1971 gmx_fio_do_real(fio
, iparams
->ljc14
.qi
);
1972 gmx_fio_do_real(fio
, iparams
->ljc14
.qj
);
1973 gmx_fio_do_real(fio
, iparams
->ljc14
.c6
);
1974 gmx_fio_do_real(fio
, iparams
->ljc14
.c12
);
1976 case F_LJC_PAIRS_NB
:
1977 gmx_fio_do_real(fio
, iparams
->ljcnb
.qi
);
1978 gmx_fio_do_real(fio
, iparams
->ljcnb
.qj
);
1979 gmx_fio_do_real(fio
, iparams
->ljcnb
.c6
);
1980 gmx_fio_do_real(fio
, iparams
->ljcnb
.c12
);
1986 gmx_fio_do_real(fio
, iparams
->pdihs
.phiA
);
1987 gmx_fio_do_real(fio
, iparams
->pdihs
.cpA
);
1988 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && file_version
< 42)
1990 /* Read the incorrectly stored multiplicity */
1991 gmx_fio_do_real(fio
, iparams
->harmonic
.rB
);
1992 gmx_fio_do_real(fio
, iparams
->harmonic
.krB
);
1993 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
1994 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
1998 gmx_fio_do_real(fio
, iparams
->pdihs
.phiB
);
1999 gmx_fio_do_real(fio
, iparams
->pdihs
.cpB
);
2000 gmx_fio_do_int(fio
, iparams
->pdihs
.mult
);
2004 gmx_fio_do_real(fio
, iparams
->pdihs
.phiA
);
2005 gmx_fio_do_real(fio
, iparams
->pdihs
.cpA
);
2008 gmx_fio_do_int(fio
, iparams
->disres
.label
);
2009 gmx_fio_do_int(fio
, iparams
->disres
.type
);
2010 gmx_fio_do_real(fio
, iparams
->disres
.low
);
2011 gmx_fio_do_real(fio
, iparams
->disres
.up1
);
2012 gmx_fio_do_real(fio
, iparams
->disres
.up2
);
2013 gmx_fio_do_real(fio
, iparams
->disres
.kfac
);
2016 gmx_fio_do_int(fio
, iparams
->orires
.ex
);
2017 gmx_fio_do_int(fio
, iparams
->orires
.label
);
2018 gmx_fio_do_int(fio
, iparams
->orires
.power
);
2019 gmx_fio_do_real(fio
, iparams
->orires
.c
);
2020 gmx_fio_do_real(fio
, iparams
->orires
.obs
);
2021 gmx_fio_do_real(fio
, iparams
->orires
.kfac
);
2024 if (file_version
< 82)
2026 gmx_fio_do_int(fio
, idum
);
2027 gmx_fio_do_int(fio
, idum
);
2029 gmx_fio_do_real(fio
, iparams
->dihres
.phiA
);
2030 gmx_fio_do_real(fio
, iparams
->dihres
.dphiA
);
2031 gmx_fio_do_real(fio
, iparams
->dihres
.kfacA
);
2032 if (file_version
>= 82)
2034 gmx_fio_do_real(fio
, iparams
->dihres
.phiB
);
2035 gmx_fio_do_real(fio
, iparams
->dihres
.dphiB
);
2036 gmx_fio_do_real(fio
, iparams
->dihres
.kfacB
);
2040 iparams
->dihres
.phiB
= iparams
->dihres
.phiA
;
2041 iparams
->dihres
.dphiB
= iparams
->dihres
.dphiA
;
2042 iparams
->dihres
.kfacB
= iparams
->dihres
.kfacA
;
2046 gmx_fio_do_rvec(fio
, iparams
->posres
.pos0A
);
2047 gmx_fio_do_rvec(fio
, iparams
->posres
.fcA
);
2048 if (bRead
&& file_version
< 27)
2050 copy_rvec(iparams
->posres
.pos0A
, iparams
->posres
.pos0B
);
2051 copy_rvec(iparams
->posres
.fcA
, iparams
->posres
.fcB
);
2055 gmx_fio_do_rvec(fio
, iparams
->posres
.pos0B
);
2056 gmx_fio_do_rvec(fio
, iparams
->posres
.fcB
);
2060 gmx_fio_do_int(fio
, iparams
->fbposres
.geom
);
2061 gmx_fio_do_rvec(fio
, iparams
->fbposres
.pos0
);
2062 gmx_fio_do_real(fio
, iparams
->fbposres
.r
);
2063 gmx_fio_do_real(fio
, iparams
->fbposres
.k
);
2066 gmx_fio_ndo_real(fio
, iparams
->cbtdihs
.cbtcA
, NR_CBTDIHS
);
2069 gmx_fio_ndo_real(fio
, iparams
->rbdihs
.rbcA
, NR_RBDIHS
);
2070 if (file_version
>= 25)
2072 gmx_fio_ndo_real(fio
, iparams
->rbdihs
.rbcB
, NR_RBDIHS
);
2076 /* Fourier dihedrals are internally represented
2077 * as Ryckaert-Bellemans since those are faster to compute.
2079 gmx_fio_ndo_real(fio
, iparams
->rbdihs
.rbcA
, NR_RBDIHS
);
2080 gmx_fio_ndo_real(fio
, iparams
->rbdihs
.rbcB
, NR_RBDIHS
);
2084 gmx_fio_do_real(fio
, iparams
->constr
.dA
);
2085 gmx_fio_do_real(fio
, iparams
->constr
.dB
);
2088 gmx_fio_do_real(fio
, iparams
->settle
.doh
);
2089 gmx_fio_do_real(fio
, iparams
->settle
.dhh
);
2092 gmx_fio_do_real(fio
, iparams
->vsite
.a
);
2097 gmx_fio_do_real(fio
, iparams
->vsite
.a
);
2098 gmx_fio_do_real(fio
, iparams
->vsite
.b
);
2103 gmx_fio_do_real(fio
, iparams
->vsite
.a
);
2104 gmx_fio_do_real(fio
, iparams
->vsite
.b
);
2105 gmx_fio_do_real(fio
, iparams
->vsite
.c
);
2108 gmx_fio_do_int(fio
, iparams
->vsiten
.n
);
2109 gmx_fio_do_real(fio
, iparams
->vsiten
.a
);
2114 /* We got rid of some parameters in version 68 */
2115 if (bRead
&& file_version
< 68)
2117 gmx_fio_do_real(fio
, rdum
);
2118 gmx_fio_do_real(fio
, rdum
);
2119 gmx_fio_do_real(fio
, rdum
);
2120 gmx_fio_do_real(fio
, rdum
);
2122 gmx_fio_do_real(fio
, iparams
->gb
.sar
);
2123 gmx_fio_do_real(fio
, iparams
->gb
.st
);
2124 gmx_fio_do_real(fio
, iparams
->gb
.pi
);
2125 gmx_fio_do_real(fio
, iparams
->gb
.gbr
);
2126 gmx_fio_do_real(fio
, iparams
->gb
.bmlt
);
2129 gmx_fio_do_int(fio
, iparams
->cmap
.cmapA
);
2130 gmx_fio_do_int(fio
, iparams
->cmap
.cmapB
);
2133 gmx_fatal(FARGS
, "unknown function type %d (%s) in %s line %d",
2134 ftype
, interaction_function
[ftype
].name
, __FILE__
, __LINE__
);
2138 static void do_ilist(t_fileio
*fio
, t_ilist
*ilist
, gmx_bool bRead
, int file_version
)
2142 if (file_version
< 44)
2144 for (i
= 0; i
< MAXNODES
; i
++)
2146 gmx_fio_do_int(fio
, idum
);
2149 gmx_fio_do_int(fio
, ilist
->nr
);
2152 snew(ilist
->iatoms
, ilist
->nr
);
2154 gmx_fio_ndo_int(fio
, ilist
->iatoms
, ilist
->nr
);
2157 static void do_ffparams(t_fileio
*fio
, gmx_ffparams_t
*ffparams
,
2158 gmx_bool bRead
, int file_version
)
2163 gmx_fio_do_int(fio
, ffparams
->atnr
);
2164 if (file_version
< 57)
2166 gmx_fio_do_int(fio
, idum
);
2168 gmx_fio_do_int(fio
, ffparams
->ntypes
);
2171 snew(ffparams
->functype
, ffparams
->ntypes
);
2172 snew(ffparams
->iparams
, ffparams
->ntypes
);
2174 /* Read/write all the function types */
2175 gmx_fio_ndo_int(fio
, ffparams
->functype
, ffparams
->ntypes
);
2177 if (file_version
>= 66)
2179 gmx_fio_do_double(fio
, ffparams
->reppow
);
2183 ffparams
->reppow
= 12.0;
2186 if (file_version
>= 57)
2188 gmx_fio_do_real(fio
, ffparams
->fudgeQQ
);
2191 /* Check whether all these function types are supported by the code.
2192 * In practice the code is backwards compatible, which means that the
2193 * numbering may have to be altered from old numbering to new numbering
2195 for (i
= 0; (i
< ffparams
->ntypes
); i
++)
2199 /* Loop over file versions */
2200 for (k
= 0; (k
< NFTUPD
); k
++)
2202 /* Compare the read file_version to the update table */
2203 if ((file_version
< ftupd
[k
].fvnr
) &&
2204 (ffparams
->functype
[i
] >= ftupd
[k
].ftype
))
2206 ffparams
->functype
[i
] += 1;
2211 do_iparams(fio
, ffparams
->functype
[i
], &ffparams
->iparams
[i
], bRead
,
2216 static void add_settle_atoms(t_ilist
*ilist
)
2220 /* Settle used to only store the first atom: add the other two */
2221 srenew(ilist
->iatoms
, 2*ilist
->nr
);
2222 for (i
= ilist
->nr
/2-1; i
>= 0; i
--)
2224 ilist
->iatoms
[4*i
+0] = ilist
->iatoms
[2*i
+0];
2225 ilist
->iatoms
[4*i
+1] = ilist
->iatoms
[2*i
+1];
2226 ilist
->iatoms
[4*i
+2] = ilist
->iatoms
[2*i
+1] + 1;
2227 ilist
->iatoms
[4*i
+3] = ilist
->iatoms
[2*i
+1] + 2;
2229 ilist
->nr
= 2*ilist
->nr
;
2232 static void do_ilists(t_fileio
*fio
, t_ilist
*ilist
, gmx_bool bRead
,
2239 for (j
= 0; (j
< F_NRE
); j
++)
2244 for (k
= 0; k
< NFTUPD
; k
++)
2246 if ((file_version
< ftupd
[k
].fvnr
) && (j
== ftupd
[k
].ftype
))
2255 ilist
[j
].iatoms
= NULL
;
2259 do_ilist(fio
, &ilist
[j
], bRead
, file_version
);
2260 if (file_version
< 78 && j
== F_SETTLE
&& ilist
[j
].nr
> 0)
2262 add_settle_atoms(&ilist
[j
]);
2268 static void do_idef(t_fileio
*fio
, gmx_ffparams_t
*ffparams
, gmx_moltype_t
*molt
,
2269 gmx_bool bRead
, int file_version
)
2271 do_ffparams(fio
, ffparams
, bRead
, file_version
);
2273 if (file_version
>= 54)
2275 gmx_fio_do_real(fio
, ffparams
->fudgeQQ
);
2278 do_ilists(fio
, molt
->ilist
, bRead
, file_version
);
2281 static void do_block(t_fileio
*fio
, t_block
*block
, gmx_bool bRead
, int file_version
)
2283 int i
, idum
, dum_nra
, *dum_a
;
2285 if (file_version
< 44)
2287 for (i
= 0; i
< MAXNODES
; i
++)
2289 gmx_fio_do_int(fio
, idum
);
2292 gmx_fio_do_int(fio
, block
->nr
);
2293 if (file_version
< 51)
2295 gmx_fio_do_int(fio
, dum_nra
);
2299 if ((block
->nalloc_index
> 0) && (NULL
!= block
->index
))
2301 sfree(block
->index
);
2303 block
->nalloc_index
= block
->nr
+1;
2304 snew(block
->index
, block
->nalloc_index
);
2306 gmx_fio_ndo_int(fio
, block
->index
, block
->nr
+1);
2308 if (file_version
< 51 && dum_nra
> 0)
2310 snew(dum_a
, dum_nra
);
2311 gmx_fio_ndo_int(fio
, dum_a
, dum_nra
);
2316 static void do_blocka(t_fileio
*fio
, t_blocka
*block
, gmx_bool bRead
,
2321 if (file_version
< 44)
2323 for (i
= 0; i
< MAXNODES
; i
++)
2325 gmx_fio_do_int(fio
, idum
);
2328 gmx_fio_do_int(fio
, block
->nr
);
2329 gmx_fio_do_int(fio
, block
->nra
);
2332 block
->nalloc_index
= block
->nr
+1;
2333 snew(block
->index
, block
->nalloc_index
);
2334 block
->nalloc_a
= block
->nra
;
2335 snew(block
->a
, block
->nalloc_a
);
2337 gmx_fio_ndo_int(fio
, block
->index
, block
->nr
+1);
2338 gmx_fio_ndo_int(fio
, block
->a
, block
->nra
);
2341 /* This is a primitive routine to make it possible to translate atomic numbers
2342 * to element names when reading TPR files, without making the Gromacs library
2343 * directory a dependency on mdrun (which is the case if we need elements.dat).
2346 atomicnumber_to_element(int atomicnumber
)
2350 /* This does not have to be complete, so we only include elements likely
2351 * to occur in PDB files.
2353 switch (atomicnumber
)
2355 case 1: p
= "H"; break;
2356 case 5: p
= "B"; break;
2357 case 6: p
= "C"; break;
2358 case 7: p
= "N"; break;
2359 case 8: p
= "O"; break;
2360 case 9: p
= "F"; break;
2361 case 11: p
= "Na"; break;
2362 case 12: p
= "Mg"; break;
2363 case 15: p
= "P"; break;
2364 case 16: p
= "S"; break;
2365 case 17: p
= "Cl"; break;
2366 case 18: p
= "Ar"; break;
2367 case 19: p
= "K"; break;
2368 case 20: p
= "Ca"; break;
2369 case 25: p
= "Mn"; break;
2370 case 26: p
= "Fe"; break;
2371 case 28: p
= "Ni"; break;
2372 case 29: p
= "Cu"; break;
2373 case 30: p
= "Zn"; break;
2374 case 35: p
= "Br"; break;
2375 case 47: p
= "Ag"; break;
2376 default: p
= ""; break;
2382 static void do_atom(t_fileio
*fio
, t_atom
*atom
, int ngrp
, gmx_bool bRead
,
2383 int file_version
, gmx_groups_t
*groups
, int atnr
)
2387 gmx_fio_do_real(fio
, atom
->m
);
2388 gmx_fio_do_real(fio
, atom
->q
);
2389 gmx_fio_do_real(fio
, atom
->mB
);
2390 gmx_fio_do_real(fio
, atom
->qB
);
2391 gmx_fio_do_ushort(fio
, atom
->type
);
2392 gmx_fio_do_ushort(fio
, atom
->typeB
);
2393 gmx_fio_do_int(fio
, atom
->ptype
);
2394 gmx_fio_do_int(fio
, atom
->resind
);
2395 if (file_version
>= 52)
2397 gmx_fio_do_int(fio
, atom
->atomnumber
);
2400 /* Set element string from atomic number if present.
2401 * This routine returns an empty string if the name is not found.
2403 std::strncpy(atom
->elem
, atomicnumber_to_element(atom
->atomnumber
), 4);
2404 /* avoid warnings about potentially unterminated string */
2405 atom
->elem
[3] = '\0';
2410 atom
->atomnumber
= 0;
2412 if (file_version
< 23)
2416 else if (file_version
< 39)
2425 if (file_version
< 57)
2427 unsigned char uchar
[egcNR
];
2428 gmx_fio_ndo_uchar(fio
, uchar
, myngrp
);
2429 for (i
= myngrp
; (i
< ngrp
); i
++)
2433 /* Copy the old data format to the groups struct */
2434 for (i
= 0; i
< ngrp
; i
++)
2436 groups
->grpnr
[i
][atnr
] = uchar
[i
];
2441 static void do_grps(t_fileio
*fio
, int ngrp
, t_grps grps
[], gmx_bool bRead
,
2446 if (file_version
< 23)
2450 else if (file_version
< 39)
2459 for (j
= 0; (j
< ngrp
); j
++)
2463 gmx_fio_do_int(fio
, grps
[j
].nr
);
2466 snew(grps
[j
].nm_ind
, grps
[j
].nr
);
2468 gmx_fio_ndo_int(fio
, grps
[j
].nm_ind
, grps
[j
].nr
);
2473 snew(grps
[j
].nm_ind
, grps
[j
].nr
);
2478 static void do_symstr(t_fileio
*fio
, char ***nm
, gmx_bool bRead
, t_symtab
*symtab
)
2484 gmx_fio_do_int(fio
, ls
);
2485 *nm
= get_symtab_handle(symtab
, ls
);
2489 ls
= lookup_symtab(symtab
, *nm
);
2490 gmx_fio_do_int(fio
, ls
);
2494 static void do_strstr(t_fileio
*fio
, int nstr
, char ***nm
, gmx_bool bRead
,
2499 for (j
= 0; (j
< nstr
); j
++)
2501 do_symstr(fio
, &(nm
[j
]), bRead
, symtab
);
2505 static void do_resinfo(t_fileio
*fio
, int n
, t_resinfo
*ri
, gmx_bool bRead
,
2506 t_symtab
*symtab
, int file_version
)
2510 for (j
= 0; (j
< n
); j
++)
2512 do_symstr(fio
, &(ri
[j
].name
), bRead
, symtab
);
2513 if (file_version
>= 63)
2515 gmx_fio_do_int(fio
, ri
[j
].nr
);
2516 gmx_fio_do_uchar(fio
, ri
[j
].ic
);
2526 static void do_atoms(t_fileio
*fio
, t_atoms
*atoms
, gmx_bool bRead
, t_symtab
*symtab
,
2528 gmx_groups_t
*groups
)
2532 gmx_fio_do_int(fio
, atoms
->nr
);
2533 gmx_fio_do_int(fio
, atoms
->nres
);
2534 if (file_version
< 57)
2536 gmx_fio_do_int(fio
, groups
->ngrpname
);
2537 for (i
= 0; i
< egcNR
; i
++)
2539 groups
->ngrpnr
[i
] = atoms
->nr
;
2540 snew(groups
->grpnr
[i
], groups
->ngrpnr
[i
]);
2545 snew(atoms
->atom
, atoms
->nr
);
2546 snew(atoms
->atomname
, atoms
->nr
);
2547 snew(atoms
->atomtype
, atoms
->nr
);
2548 snew(atoms
->atomtypeB
, atoms
->nr
);
2549 snew(atoms
->resinfo
, atoms
->nres
);
2550 if (file_version
< 57)
2552 snew(groups
->grpname
, groups
->ngrpname
);
2554 atoms
->pdbinfo
= NULL
;
2556 for (i
= 0; (i
< atoms
->nr
); i
++)
2558 do_atom(fio
, &atoms
->atom
[i
], egcNR
, bRead
, file_version
, groups
, i
);
2560 do_strstr(fio
, atoms
->nr
, atoms
->atomname
, bRead
, symtab
);
2561 if (bRead
&& (file_version
<= 20))
2563 for (i
= 0; i
< atoms
->nr
; i
++)
2565 atoms
->atomtype
[i
] = put_symtab(symtab
, "?");
2566 atoms
->atomtypeB
[i
] = put_symtab(symtab
, "?");
2571 do_strstr(fio
, atoms
->nr
, atoms
->atomtype
, bRead
, symtab
);
2572 do_strstr(fio
, atoms
->nr
, atoms
->atomtypeB
, bRead
, symtab
);
2574 do_resinfo(fio
, atoms
->nres
, atoms
->resinfo
, bRead
, symtab
, file_version
);
2576 if (file_version
< 57)
2578 do_strstr(fio
, groups
->ngrpname
, groups
->grpname
, bRead
, symtab
);
2580 do_grps(fio
, egcNR
, groups
->grps
, bRead
, file_version
);
2584 static void do_groups(t_fileio
*fio
, gmx_groups_t
*groups
,
2585 gmx_bool bRead
, t_symtab
*symtab
,
2590 do_grps(fio
, egcNR
, groups
->grps
, bRead
, file_version
);
2591 gmx_fio_do_int(fio
, groups
->ngrpname
);
2594 snew(groups
->grpname
, groups
->ngrpname
);
2596 do_strstr(fio
, groups
->ngrpname
, groups
->grpname
, bRead
, symtab
);
2597 for (g
= 0; g
< egcNR
; g
++)
2599 gmx_fio_do_int(fio
, groups
->ngrpnr
[g
]);
2600 if (groups
->ngrpnr
[g
] == 0)
2604 groups
->grpnr
[g
] = NULL
;
2611 snew(groups
->grpnr
[g
], groups
->ngrpnr
[g
]);
2613 gmx_fio_ndo_uchar(fio
, groups
->grpnr
[g
], groups
->ngrpnr
[g
]);
2618 static void do_atomtypes(t_fileio
*fio
, t_atomtypes
*atomtypes
, gmx_bool bRead
,
2623 if (file_version
> 25)
2625 gmx_fio_do_int(fio
, atomtypes
->nr
);
2629 snew(atomtypes
->radius
, j
);
2630 snew(atomtypes
->vol
, j
);
2631 snew(atomtypes
->surftens
, j
);
2632 snew(atomtypes
->atomnumber
, j
);
2633 snew(atomtypes
->gb_radius
, j
);
2634 snew(atomtypes
->S_hct
, j
);
2636 gmx_fio_ndo_real(fio
, atomtypes
->radius
, j
);
2637 gmx_fio_ndo_real(fio
, atomtypes
->vol
, j
);
2638 gmx_fio_ndo_real(fio
, atomtypes
->surftens
, j
);
2639 if (file_version
>= 40)
2641 gmx_fio_ndo_int(fio
, atomtypes
->atomnumber
, j
);
2643 if (file_version
>= 60)
2645 gmx_fio_ndo_real(fio
, atomtypes
->gb_radius
, j
);
2646 gmx_fio_ndo_real(fio
, atomtypes
->S_hct
, j
);
2651 /* File versions prior to 26 cannot do GBSA,
2652 * so they dont use this structure
2655 atomtypes
->radius
= NULL
;
2656 atomtypes
->vol
= NULL
;
2657 atomtypes
->surftens
= NULL
;
2658 atomtypes
->atomnumber
= NULL
;
2659 atomtypes
->gb_radius
= NULL
;
2660 atomtypes
->S_hct
= NULL
;
2664 static void do_symtab(t_fileio
*fio
, t_symtab
*symtab
, gmx_bool bRead
)
2670 gmx_fio_do_int(fio
, symtab
->nr
);
2674 snew(symtab
->symbuf
, 1);
2675 symbuf
= symtab
->symbuf
;
2676 symbuf
->bufsize
= nr
;
2677 snew(symbuf
->buf
, nr
);
2678 for (i
= 0; (i
< nr
); i
++)
2680 gmx_fio_do_string(fio
, buf
);
2681 symbuf
->buf
[i
] = gmx_strdup(buf
);
2686 symbuf
= symtab
->symbuf
;
2687 while (symbuf
!= NULL
)
2689 for (i
= 0; (i
< symbuf
->bufsize
) && (i
< nr
); i
++)
2691 gmx_fio_do_string(fio
, symbuf
->buf
[i
]);
2694 symbuf
= symbuf
->next
;
2698 gmx_fatal(FARGS
, "nr of symtab strings left: %d", nr
);
2703 static void do_cmap(t_fileio
*fio
, gmx_cmap_t
*cmap_grid
, gmx_bool bRead
)
2705 int i
, j
, ngrid
, gs
, nelem
;
2707 gmx_fio_do_int(fio
, cmap_grid
->ngrid
);
2708 gmx_fio_do_int(fio
, cmap_grid
->grid_spacing
);
2710 ngrid
= cmap_grid
->ngrid
;
2711 gs
= cmap_grid
->grid_spacing
;
2716 snew(cmap_grid
->cmapdata
, ngrid
);
2718 for (i
= 0; i
< cmap_grid
->ngrid
; i
++)
2720 snew(cmap_grid
->cmapdata
[i
].cmap
, 4*nelem
);
2724 for (i
= 0; i
< cmap_grid
->ngrid
; i
++)
2726 for (j
= 0; j
< nelem
; j
++)
2728 gmx_fio_do_real(fio
, cmap_grid
->cmapdata
[i
].cmap
[j
*4]);
2729 gmx_fio_do_real(fio
, cmap_grid
->cmapdata
[i
].cmap
[j
*4+1]);
2730 gmx_fio_do_real(fio
, cmap_grid
->cmapdata
[i
].cmap
[j
*4+2]);
2731 gmx_fio_do_real(fio
, cmap_grid
->cmapdata
[i
].cmap
[j
*4+3]);
2737 static void do_moltype(t_fileio
*fio
, gmx_moltype_t
*molt
, gmx_bool bRead
,
2738 t_symtab
*symtab
, int file_version
,
2739 gmx_groups_t
*groups
)
2741 if (file_version
>= 57)
2743 do_symstr(fio
, &(molt
->name
), bRead
, symtab
);
2746 do_atoms(fio
, &molt
->atoms
, bRead
, symtab
, file_version
, groups
);
2748 if (file_version
>= 57)
2750 do_ilists(fio
, molt
->ilist
, bRead
, file_version
);
2752 do_block(fio
, &molt
->cgs
, bRead
, file_version
);
2755 /* This used to be in the atoms struct */
2756 do_blocka(fio
, &molt
->excls
, bRead
, file_version
);
2759 static void do_molblock(t_fileio
*fio
, gmx_molblock_t
*molb
, gmx_bool bRead
)
2761 gmx_fio_do_int(fio
, molb
->type
);
2762 gmx_fio_do_int(fio
, molb
->nmol
);
2763 gmx_fio_do_int(fio
, molb
->natoms_mol
);
2764 /* Position restraint coordinates */
2765 gmx_fio_do_int(fio
, molb
->nposres_xA
);
2766 if (molb
->nposres_xA
> 0)
2770 snew(molb
->posres_xA
, molb
->nposres_xA
);
2772 gmx_fio_ndo_rvec(fio
, molb
->posres_xA
, molb
->nposres_xA
);
2774 gmx_fio_do_int(fio
, molb
->nposres_xB
);
2775 if (molb
->nposres_xB
> 0)
2779 snew(molb
->posres_xB
, molb
->nposres_xB
);
2781 gmx_fio_ndo_rvec(fio
, molb
->posres_xB
, molb
->nposres_xB
);
2786 static t_block
mtop_mols(gmx_mtop_t
*mtop
)
2792 for (mb
= 0; mb
< mtop
->nmolblock
; mb
++)
2794 mols
.nr
+= mtop
->molblock
[mb
].nmol
;
2796 mols
.nalloc_index
= mols
.nr
+ 1;
2797 snew(mols
.index
, mols
.nalloc_index
);
2802 for (mb
= 0; mb
< mtop
->nmolblock
; mb
++)
2804 for (mol
= 0; mol
< mtop
->molblock
[mb
].nmol
; mol
++)
2806 a
+= mtop
->molblock
[mb
].natoms_mol
;
2815 static void add_posres_molblock(gmx_mtop_t
*mtop
)
2820 gmx_molblock_t
*molb
;
2823 /* posres reference positions are stored in ip->posres (if present) and
2824 in ip->fbposres (if present). If normal and flat-bottomed posres are present,
2825 posres.pos0A are identical to fbposres.pos0. */
2826 il
= &mtop
->moltype
[0].ilist
[F_POSRES
];
2827 ilfb
= &mtop
->moltype
[0].ilist
[F_FBPOSRES
];
2828 if (il
->nr
== 0 && ilfb
->nr
== 0)
2834 for (i
= 0; i
< il
->nr
; i
+= 2)
2836 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
2837 am
= std::max(am
, il
->iatoms
[i
+1]);
2838 if (ip
->posres
.pos0B
[XX
] != ip
->posres
.pos0A
[XX
] ||
2839 ip
->posres
.pos0B
[YY
] != ip
->posres
.pos0A
[YY
] ||
2840 ip
->posres
.pos0B
[ZZ
] != ip
->posres
.pos0A
[ZZ
])
2845 /* This loop is required if we have only flat-bottomed posres:
2847 - bFE == FALSE (no B-state for flat-bottomed posres) */
2850 for (i
= 0; i
< ilfb
->nr
; i
+= 2)
2852 am
= std::max(am
, ilfb
->iatoms
[i
+1]);
2855 /* Make the posres coordinate block end at a molecule end */
2857 while (am
>= mtop
->mols
.index
[mol
+1])
2861 molb
= &mtop
->molblock
[0];
2862 molb
->nposres_xA
= mtop
->mols
.index
[mol
+1];
2863 snew(molb
->posres_xA
, molb
->nposres_xA
);
2866 molb
->nposres_xB
= molb
->nposres_xA
;
2867 snew(molb
->posres_xB
, molb
->nposres_xB
);
2871 molb
->nposres_xB
= 0;
2873 for (i
= 0; i
< il
->nr
; i
+= 2)
2875 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
2876 a
= il
->iatoms
[i
+1];
2877 molb
->posres_xA
[a
][XX
] = ip
->posres
.pos0A
[XX
];
2878 molb
->posres_xA
[a
][YY
] = ip
->posres
.pos0A
[YY
];
2879 molb
->posres_xA
[a
][ZZ
] = ip
->posres
.pos0A
[ZZ
];
2882 molb
->posres_xB
[a
][XX
] = ip
->posres
.pos0B
[XX
];
2883 molb
->posres_xB
[a
][YY
] = ip
->posres
.pos0B
[YY
];
2884 molb
->posres_xB
[a
][ZZ
] = ip
->posres
.pos0B
[ZZ
];
2889 /* If only flat-bottomed posres are present, take reference pos from them.
2890 Here: bFE == FALSE */
2891 for (i
= 0; i
< ilfb
->nr
; i
+= 2)
2893 ip
= &mtop
->ffparams
.iparams
[ilfb
->iatoms
[i
]];
2894 a
= ilfb
->iatoms
[i
+1];
2895 molb
->posres_xA
[a
][XX
] = ip
->fbposres
.pos0
[XX
];
2896 molb
->posres_xA
[a
][YY
] = ip
->fbposres
.pos0
[YY
];
2897 molb
->posres_xA
[a
][ZZ
] = ip
->fbposres
.pos0
[ZZ
];
2902 static void set_disres_npair(gmx_mtop_t
*mtop
)
2905 gmx_mtop_ilistloop_t iloop
;
2906 t_ilist
*ilist
, *il
;
2910 ip
= mtop
->ffparams
.iparams
;
2912 iloop
= gmx_mtop_ilistloop_init(mtop
);
2913 while (gmx_mtop_ilistloop_next(iloop
, &ilist
, &nmol
))
2915 il
= &ilist
[F_DISRES
];
2921 for (i
= 0; i
< il
->nr
; i
+= 3)
2924 if (i
+3 == il
->nr
|| ip
[a
[i
]].disres
.label
!= ip
[a
[i
+3]].disres
.label
)
2926 ip
[a
[i
]].disres
.npair
= npair
;
2934 static void do_mtop(t_fileio
*fio
, gmx_mtop_t
*mtop
, gmx_bool bRead
,
2944 do_symtab(fio
, &(mtop
->symtab
), bRead
);
2946 do_symstr(fio
, &(mtop
->name
), bRead
, &(mtop
->symtab
));
2948 if (file_version
>= 57)
2950 do_ffparams(fio
, &mtop
->ffparams
, bRead
, file_version
);
2952 gmx_fio_do_int(fio
, mtop
->nmoltype
);
2960 snew(mtop
->moltype
, mtop
->nmoltype
);
2961 if (file_version
< 57)
2963 mtop
->moltype
[0].name
= mtop
->name
;
2966 for (mt
= 0; mt
< mtop
->nmoltype
; mt
++)
2968 do_moltype(fio
, &mtop
->moltype
[mt
], bRead
, &mtop
->symtab
, file_version
,
2972 if (file_version
>= 57)
2974 gmx_fio_do_int(fio
, mtop
->nmolblock
);
2978 mtop
->nmolblock
= 1;
2982 snew(mtop
->molblock
, mtop
->nmolblock
);
2984 if (file_version
>= 57)
2986 for (mb
= 0; mb
< mtop
->nmolblock
; mb
++)
2988 do_molblock(fio
, &mtop
->molblock
[mb
], bRead
);
2990 gmx_fio_do_int(fio
, mtop
->natoms
);
2994 mtop
->molblock
[0].type
= 0;
2995 mtop
->molblock
[0].nmol
= 1;
2996 mtop
->molblock
[0].natoms_mol
= mtop
->moltype
[0].atoms
.nr
;
2997 mtop
->molblock
[0].nposres_xA
= 0;
2998 mtop
->molblock
[0].nposres_xB
= 0;
3001 if (file_version
>= tpxv_IntermolecularBondeds
)
3003 gmx_fio_do_gmx_bool(fio
, mtop
->bIntermolecularInteractions
);
3004 if (mtop
->bIntermolecularInteractions
)
3008 snew(mtop
->intermolecular_ilist
, F_NRE
);
3010 do_ilists(fio
, mtop
->intermolecular_ilist
, bRead
, file_version
);
3015 mtop
->bIntermolecularInteractions
= FALSE
;
3018 do_atomtypes (fio
, &(mtop
->atomtypes
), bRead
, file_version
);
3020 if (file_version
< 57)
3022 do_idef (fio
, &mtop
->ffparams
, &mtop
->moltype
[0], bRead
, file_version
);
3023 mtop
->natoms
= mtop
->moltype
[0].atoms
.nr
;
3026 if (file_version
>= 65)
3028 do_cmap(fio
, &mtop
->ffparams
.cmap_grid
, bRead
);
3032 mtop
->ffparams
.cmap_grid
.ngrid
= 0;
3033 mtop
->ffparams
.cmap_grid
.grid_spacing
= 0;
3034 mtop
->ffparams
.cmap_grid
.cmapdata
= NULL
;
3037 if (file_version
>= 57)
3039 do_groups(fio
, &mtop
->groups
, bRead
, &(mtop
->symtab
), file_version
);
3042 if (file_version
< 57)
3044 do_block(fio
, &mtop
->moltype
[0].cgs
, bRead
, file_version
);
3045 do_block(fio
, &mtop
->mols
, bRead
, file_version
);
3046 /* Add the posres coordinates to the molblock */
3047 add_posres_molblock(mtop
);
3051 if (file_version
>= 57)
3053 done_block(&mtop
->mols
);
3054 mtop
->mols
= mtop_mols(mtop
);
3058 if (file_version
< 51)
3060 /* Here used to be the shake blocks */
3061 do_blocka(fio
, &dumb
, bRead
, file_version
);
3074 close_symtab(&(mtop
->symtab
));
3078 /* If TopOnlyOK is TRUE then we can read even future versions
3079 * of tpx files, provided the file_generation hasn't changed.
3080 * If it is FALSE, we need the inputrecord too, and bail out
3081 * if the file is newer than the program.
3083 * The version and generation if the topology (see top of this file)
3084 * are returned in the two last arguments.
3086 * If possible, we will read the inputrec even when TopOnlyOK is TRUE.
3088 static void do_tpxheader(t_fileio
*fio
, gmx_bool bRead
, t_tpxheader
*tpx
,
3089 gmx_bool TopOnlyOK
, int *file_version
,
3090 int *file_generation
)
3093 char file_tag
[STRLEN
];
3100 /* XDR binary topology file */
3101 precision
= sizeof(real
);
3104 gmx_fio_do_string(fio
, buf
);
3105 if (std::strncmp(buf
, "VERSION", 7))
3107 gmx_fatal(FARGS
, "Can not read file %s,\n"
3108 " this file is from a GROMACS version which is older than 2.0\n"
3109 " Make a new one with grompp or use a gro or pdb file, if possible",
3110 gmx_fio_getname(fio
));
3112 gmx_fio_do_int(fio
, precision
);
3113 bDouble
= (precision
== sizeof(double));
3114 if ((precision
!= sizeof(float)) && !bDouble
)
3116 gmx_fatal(FARGS
, "Unknown precision in file %s: real is %d bytes "
3117 "instead of %d or %d",
3118 gmx_fio_getname(fio
), precision
, sizeof(float), sizeof(double));
3120 gmx_fio_setprecision(fio
, bDouble
);
3121 fprintf(stderr
, "Reading file %s, %s (%s precision)\n",
3122 gmx_fio_getname(fio
), buf
, bDouble
? "double" : "single");
3126 gmx_fio_write_string(fio
, gmx_version());
3127 bDouble
= (precision
== sizeof(double));
3128 gmx_fio_setprecision(fio
, bDouble
);
3129 gmx_fio_do_int(fio
, precision
);
3131 sprintf(file_tag
, "%s", tpx_tag
);
3132 fgen
= tpx_generation
;
3135 /* Check versions! */
3136 gmx_fio_do_int(fio
, fver
);
3138 /* This is for backward compatibility with development versions 77-79
3139 * where the tag was, mistakenly, placed before the generation,
3140 * which would cause a segv instead of a proper error message
3141 * when reading the topology only from tpx with <77 code.
3143 if (fver
>= 77 && fver
<= 79)
3145 gmx_fio_do_string(fio
, file_tag
);
3150 gmx_fio_do_int(fio
, fgen
);
3159 gmx_fio_do_string(fio
, file_tag
);
3165 /* Versions before 77 don't have the tag, set it to release */
3166 sprintf(file_tag
, "%s", TPX_TAG_RELEASE
);
3169 if (std::strcmp(file_tag
, tpx_tag
) != 0)
3171 fprintf(stderr
, "Note: file tpx tag '%s', software tpx tag '%s'\n",
3174 /* We only support reading tpx files with the same tag as the code
3175 * or tpx files with the release tag and with lower version number.
3177 if (std::strcmp(file_tag
, TPX_TAG_RELEASE
) != 0 && fver
< tpx_version
)
3179 gmx_fatal(FARGS
, "tpx tag/version mismatch: reading tpx file (%s) version %d, tag '%s' with program for tpx version %d, tag '%s'",
3180 gmx_fio_getname(fio
), fver
, file_tag
,
3181 tpx_version
, tpx_tag
);
3186 if (file_version
!= NULL
)
3188 *file_version
= fver
;
3190 if (file_generation
!= NULL
)
3192 *file_generation
= fgen
;
3196 if ((fver
<= tpx_incompatible_version
) ||
3197 ((fver
> tpx_version
) && !TopOnlyOK
) ||
3198 (fgen
> tpx_generation
) ||
3199 tpx_version
== 80) /*80 was used by both 5.0-dev and 4.6-dev*/
3201 gmx_fatal(FARGS
, "reading tpx file (%s) version %d with version %d program",
3202 gmx_fio_getname(fio
), fver
, tpx_version
);
3205 gmx_fio_do_int(fio
, tpx
->natoms
);
3208 gmx_fio_do_int(fio
, tpx
->ngtc
);
3216 gmx_fio_do_int(fio
, idum
);
3217 gmx_fio_do_real(fio
, rdum
);
3221 gmx_fio_do_int(fio
, tpx
->fep_state
);
3223 gmx_fio_do_real(fio
, tpx
->lambda
);
3224 gmx_fio_do_int(fio
, tpx
->bIr
);
3225 gmx_fio_do_int(fio
, tpx
->bTop
);
3226 gmx_fio_do_int(fio
, tpx
->bX
);
3227 gmx_fio_do_int(fio
, tpx
->bV
);
3228 gmx_fio_do_int(fio
, tpx
->bF
);
3229 gmx_fio_do_int(fio
, tpx
->bBox
);
3231 if ((fgen
> tpx_generation
))
3233 /* This can only happen if TopOnlyOK=TRUE */
3238 static int do_tpx(t_fileio
*fio
, gmx_bool bRead
,
3239 t_inputrec
*ir
, t_state
*state
, gmx_mtop_t
*mtop
,
3240 gmx_bool bXVallocated
)
3246 int file_version
, file_generation
;
3249 gmx_bool bPeriodicMols
;
3253 tpx
.natoms
= state
->natoms
;
3254 tpx
.ngtc
= state
->ngtc
;
3255 tpx
.fep_state
= state
->fep_state
;
3256 tpx
.lambda
= state
->lambda
[efptFEP
];
3257 tpx
.bIr
= (ir
!= NULL
);
3258 tpx
.bTop
= (mtop
!= NULL
);
3259 tpx
.bX
= (state
->x
!= NULL
);
3260 tpx
.bV
= (state
->v
!= NULL
);
3265 TopOnlyOK
= (ir
== NULL
);
3267 do_tpxheader(fio
, bRead
, &tpx
, TopOnlyOK
, &file_version
, &file_generation
);
3276 init_state(state
, 0, tpx
.ngtc
, 0, 0, 0);
3277 state
->natoms
= tpx
.natoms
;
3278 state
->nalloc
= tpx
.natoms
;
3284 init_state(state
, tpx
.natoms
, tpx
.ngtc
, 0, 0, 0);
3288 #define do_test(fio, b, p) if (bRead && (p != NULL) && !b) gmx_fatal(FARGS, "No %s in %s",#p, gmx_fio_getname(fio))
3290 do_test(fio
, tpx
.bBox
, state
->box
);
3293 gmx_fio_ndo_rvec(fio
, state
->box
, DIM
);
3294 if (file_version
>= 51)
3296 gmx_fio_ndo_rvec(fio
, state
->box_rel
, DIM
);
3300 /* We initialize box_rel after reading the inputrec */
3301 clear_mat(state
->box_rel
);
3303 if (file_version
>= 28)
3305 gmx_fio_ndo_rvec(fio
, state
->boxv
, DIM
);
3306 if (file_version
< 56)
3309 gmx_fio_ndo_rvec(fio
, mdum
, DIM
);
3314 if (state
->ngtc
> 0 && file_version
>= 28)
3317 snew(dumv
, state
->ngtc
);
3318 if (file_version
< 69)
3320 gmx_fio_ndo_real(fio
, dumv
, state
->ngtc
);
3322 /* These used to be the Berendsen tcoupl_lambda's */
3323 gmx_fio_ndo_real(fio
, dumv
, state
->ngtc
);
3327 /* Prior to tpx version 26, the inputrec was here.
3328 * I moved it to enable partial forward-compatibility
3329 * for analysis/viewer programs.
3331 if (file_version
< 26)
3333 do_test(fio
, tpx
.bIr
, ir
);
3338 do_inputrec(fio
, ir
, bRead
, file_version
,
3339 mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
3343 do_inputrec(fio
, &dum_ir
, bRead
, file_version
,
3344 mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
3345 done_inputrec(&dum_ir
);
3351 do_test(fio
, tpx
.bTop
, mtop
);
3356 do_mtop(fio
, mtop
, bRead
, file_version
);
3360 do_mtop(fio
, &dum_top
, bRead
, file_version
);
3361 done_mtop(&dum_top
, TRUE
);
3364 do_test(fio
, tpx
.bX
, state
->x
);
3369 state
->flags
|= (1<<estX
);
3371 gmx_fio_ndo_rvec(fio
, state
->x
, state
->natoms
);
3374 do_test(fio
, tpx
.bV
, state
->v
);
3379 state
->flags
|= (1<<estV
);
3381 gmx_fio_ndo_rvec(fio
, state
->v
, state
->natoms
);
3384 // No need to run do_test when the last argument is NULL
3388 snew(dummyForces
, state
->natoms
);
3389 gmx_fio_ndo_rvec(fio
, dummyForces
, state
->natoms
);
3393 /* Starting with tpx version 26, we have the inputrec
3394 * at the end of the file, so we can ignore it
3395 * if the file is never than the software (but still the
3396 * same generation - see comments at the top of this file.
3401 bPeriodicMols
= FALSE
;
3402 if (file_version
>= 26)
3404 do_test(fio
, tpx
.bIr
, ir
);
3407 if (file_version
>= 53)
3409 /* Removed the pbc info from do_inputrec, since we always want it */
3413 bPeriodicMols
= ir
->bPeriodicMols
;
3415 gmx_fio_do_int(fio
, ePBC
);
3416 gmx_fio_do_gmx_bool(fio
, bPeriodicMols
);
3418 if (file_generation
<= tpx_generation
&& ir
)
3420 do_inputrec(fio
, ir
, bRead
, file_version
, mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
3421 if (file_version
< 51)
3423 set_box_rel(ir
, state
);
3425 if (file_version
< 53)
3428 bPeriodicMols
= ir
->bPeriodicMols
;
3431 if (bRead
&& ir
&& file_version
>= 53)
3433 /* We need to do this after do_inputrec, since that initializes ir */
3435 ir
->bPeriodicMols
= bPeriodicMols
;
3444 if (state
->ngtc
== 0)
3446 /* Reading old version without tcoupl state data: set it */
3447 init_gtc_state(state
, ir
->opts
.ngtc
, 0, ir
->opts
.nhchainlength
);
3449 if (tpx
.bTop
&& mtop
)
3451 if (file_version
< 57)
3453 if (mtop
->moltype
[0].ilist
[F_DISRES
].nr
> 0)
3455 ir
->eDisre
= edrSimple
;
3459 ir
->eDisre
= edrNone
;
3462 set_disres_npair(mtop
);
3466 if (tpx
.bTop
&& mtop
)
3468 gmx_mtop_finalize(mtop
);
3475 static t_fileio
*open_tpx(const char *fn
, const char *mode
)
3477 return gmx_fio_open(fn
, mode
);
3480 static void close_tpx(t_fileio
*fio
)
3485 /************************************************************
3487 * The following routines are the exported ones
3489 ************************************************************/
3491 void read_tpxheader(const char *fn
, t_tpxheader
*tpx
, gmx_bool TopOnlyOK
,
3492 int *file_version
, int *file_generation
)
3496 fio
= open_tpx(fn
, "r");
3497 do_tpxheader(fio
, TRUE
, tpx
, TopOnlyOK
, file_version
, file_generation
);
3501 void write_tpx_state(const char *fn
,
3502 t_inputrec
*ir
, t_state
*state
, gmx_mtop_t
*mtop
)
3506 fio
= open_tpx(fn
, "w");
3507 do_tpx(fio
, FALSE
, ir
, state
, mtop
, FALSE
);
3511 void read_tpx_state(const char *fn
,
3512 t_inputrec
*ir
, t_state
*state
, gmx_mtop_t
*mtop
)
3516 fio
= open_tpx(fn
, "r");
3517 do_tpx(fio
, TRUE
, ir
, state
, mtop
, FALSE
);
3521 int read_tpx(const char *fn
,
3522 t_inputrec
*ir
, matrix box
, int *natoms
,
3523 rvec
*x
, rvec
*v
, gmx_mtop_t
*mtop
)
3531 fio
= open_tpx(fn
, "r");
3532 ePBC
= do_tpx(fio
, TRUE
, ir
, &state
, mtop
, TRUE
);
3534 *natoms
= state
.natoms
;
3537 copy_mat(state
.box
, box
);
3546 int read_tpx_top(const char *fn
,
3547 t_inputrec
*ir
, matrix box
, int *natoms
,
3548 rvec
*x
, rvec
*v
, t_topology
*top
)
3553 ePBC
= read_tpx(fn
, ir
, box
, natoms
, x
, v
, &mtop
);
3555 *top
= gmx_mtop_t_to_t_topology(&mtop
);
3560 gmx_bool
fn2bTPX(const char *file
)
3562 return (efTPR
== fn2ftp(file
));