1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
34 * GROningen Mixture of Alchemy and Childrens' Stories
40 /* This file is completely threadsafe - keep it that way! */
42 #include <thread_mpi.h>
50 #include "gmx_fatal.h"
64 #include "mtop_util.h"
66 #define TPX_TAG_RELEASE "release"
68 /* This is the tag string which is stored in the tpx file.
69 * Change this if you want to change the tpx format in a feature branch.
70 * This ensures that there will not be different tpx formats around which
71 * can not be distinguished.
73 static const char *tpx_tag
= TPX_TAG_RELEASE
;
75 /* This number should be increased whenever the file format changes! */
76 static const int tpx_version
= 80;
78 /* This number should only be increased when you edit the TOPOLOGY section
79 * or the HEADER of the tpx format.
80 * This way we can maintain forward compatibility too for all analysis tools
81 * and/or external programs that only need to know the atom/residue names,
82 * charges, and bond connectivity.
84 * It first appeared in tpx version 26, when I also moved the inputrecord
85 * to the end of the tpx file, so we can just skip it if we only
88 static const int tpx_generation
= 24;
90 /* This number should be the most recent backwards incompatible version
91 * I.e., if this number is 9, we cannot read tpx version 9 with this code.
93 static const int tpx_incompatible_version
= 9;
97 /* Struct used to maintain tpx compatibility when function types are added */
99 int fvnr
; /* file version number in which the function type first appeared */
100 int ftype
; /* function type */
104 *The entries should be ordered in:
105 * 1. ascending file version number
106 * 2. ascending function type number
108 /*static const t_ftupd ftupd[] = {
109 { 20, F_CUBICBONDS },
113 { 22, F_DISRESVIOL },
119 { 26, F_DIHRESVIOL },
120 { 30, F_CROSS_BOND_BONDS },
121 { 30, F_CROSS_BOND_ANGLES },
122 { 30, F_UREY_BRADLEY },
123 { 30, F_POLARIZATION },
127 *The entries should be ordered in:
128 * 1. ascending function type number
129 * 2. ascending file version number
131 /* question; what is the purpose of the commented code above? */
132 static const t_ftupd ftupd
[] = {
133 { 20, F_CUBICBONDS
},
138 { 43, F_TABBONDSNC
},
139 { 70, F_RESTRBONDS
},
140 { 76, F_LINEAR_ANGLES
},
141 { 30, F_CROSS_BOND_BONDS
},
142 { 30, F_CROSS_BOND_ANGLES
},
143 { 30, F_UREY_BRADLEY
},
144 { 34, F_QUARTIC_ANGLES
},
154 { 72, F_NPSOLVATION
},
156 { 41, F_LJC_PAIRS_NB
},
159 { 32, F_COUL_RECIP
},
161 { 30, F_POLARIZATION
},
163 { 22, F_DISRESVIOL
},
167 { 26, F_DIHRESVIOL
},
172 { 46, F_ECONSERVED
},
176 { 76, F_ANHARM_POL
},
179 { 79, F_DVDL_BONDED
, },
180 { 79, F_DVDL_RESTRAINT
},
181 { 79, F_DVDL_TEMPERATURE
},
184 #define NFTUPD asize(ftupd)
186 /* Needed for backward compatibility */
189 static void _do_section(t_fileio
*fio
,int key
,gmx_bool bRead
,const char *src
,
195 if (gmx_fio_getftp(fio
) == efTPA
) {
197 gmx_fio_write_string(fio
,itemstr
[key
]);
198 bDbg
= gmx_fio_getdebug(fio
);
199 gmx_fio_setdebug(fio
,FALSE
);
200 gmx_fio_write_string(fio
,comment_str
[key
]);
201 gmx_fio_setdebug(fio
,bDbg
);
204 if (gmx_fio_getdebug(fio
))
205 fprintf(stderr
,"Looking for section %s (%s, %d)",
206 itemstr
[key
],src
,line
);
209 gmx_fio_do_string(fio
,buf
);
210 } while ((gmx_strcasecmp(buf
,itemstr
[key
]) != 0));
212 if (gmx_strcasecmp(buf
,itemstr
[key
]) != 0)
213 gmx_fatal(FARGS
,"\nCould not find section heading %s",itemstr
[key
]);
214 else if (gmx_fio_getdebug(fio
))
215 fprintf(stderr
," and found it\n");
220 #define do_section(fio,key,bRead) _do_section(fio,key,bRead,__FILE__,__LINE__)
222 /**************************************************************
224 * Now the higer level routines that do io of the structures and arrays
226 **************************************************************/
227 static void do_pullgrp(t_fileio
*fio
, t_pullgrp
*pgrp
, gmx_bool bRead
,
233 gmx_fio_do_int(fio
,pgrp
->nat
);
235 snew(pgrp
->ind
,pgrp
->nat
);
236 bDum
=gmx_fio_ndo_int(fio
,pgrp
->ind
,pgrp
->nat
);
237 gmx_fio_do_int(fio
,pgrp
->nweight
);
239 snew(pgrp
->weight
,pgrp
->nweight
);
240 bDum
=gmx_fio_ndo_real(fio
,pgrp
->weight
,pgrp
->nweight
);
241 gmx_fio_do_int(fio
,pgrp
->pbcatom
);
242 gmx_fio_do_rvec(fio
,pgrp
->vec
);
243 gmx_fio_do_rvec(fio
,pgrp
->init
);
244 gmx_fio_do_real(fio
,pgrp
->rate
);
245 gmx_fio_do_real(fio
,pgrp
->k
);
246 if (file_version
>= 56) {
247 gmx_fio_do_real(fio
,pgrp
->kB
);
253 static void do_expandedvals(t_fileio
*fio
,t_expanded
*expand
,int n_lambda
, gmx_bool bRead
, int file_version
)
255 /* i is used in the ndo_double macro*/
261 if (file_version
>= 79)
267 snew(expand
->init_lambda_weights
,n_lambda
);
269 bDum
=gmx_fio_ndo_real(fio
,expand
->init_lambda_weights
,n_lambda
);
270 gmx_fio_do_gmx_bool(fio
,expand
->bInit_weights
);
273 gmx_fio_do_int(fio
,expand
->nstexpanded
);
274 gmx_fio_do_int(fio
,expand
->elmcmove
);
275 gmx_fio_do_int(fio
,expand
->elamstats
);
276 gmx_fio_do_int(fio
,expand
->lmc_repeats
);
277 gmx_fio_do_int(fio
,expand
->gibbsdeltalam
);
278 gmx_fio_do_int(fio
,expand
->lmc_forced_nstart
);
279 gmx_fio_do_int(fio
,expand
->lmc_seed
);
280 gmx_fio_do_real(fio
,expand
->mc_temp
);
281 gmx_fio_do_int(fio
,expand
->bSymmetrizedTMatrix
);
282 gmx_fio_do_int(fio
,expand
->nstTij
);
283 gmx_fio_do_int(fio
,expand
->minvarmin
);
284 gmx_fio_do_int(fio
,expand
->c_range
);
285 gmx_fio_do_real(fio
,expand
->wl_scale
);
286 gmx_fio_do_real(fio
,expand
->wl_ratio
);
287 gmx_fio_do_real(fio
,expand
->init_wl_delta
);
288 gmx_fio_do_gmx_bool(fio
,expand
->bWLoneovert
);
289 gmx_fio_do_int(fio
,expand
->elmceq
);
290 gmx_fio_do_int(fio
,expand
->equil_steps
);
291 gmx_fio_do_int(fio
,expand
->equil_samples
);
292 gmx_fio_do_int(fio
,expand
->equil_n_at_lam
);
293 gmx_fio_do_real(fio
,expand
->equil_wl_delta
);
294 gmx_fio_do_real(fio
,expand
->equil_ratio
);
298 static void do_simtempvals(t_fileio
*fio
,t_simtemp
*simtemp
, int n_lambda
, gmx_bool bRead
,
303 if (file_version
>= 79)
305 gmx_fio_do_int(fio
,simtemp
->eSimTempScale
);
306 gmx_fio_do_real(fio
,simtemp
->simtemp_high
);
307 gmx_fio_do_real(fio
,simtemp
->simtemp_low
);
312 snew(simtemp
->temperatures
,n_lambda
);
314 bDum
=gmx_fio_ndo_real(fio
,simtemp
->temperatures
,n_lambda
);
319 static void do_fepvals(t_fileio
*fio
,t_lambda
*fepvals
,gmx_bool bRead
, int file_version
)
321 /* i is defined in the ndo_double macro; use g to iterate. */
327 /* free energy values */
328 if (file_version
>= 79)
330 gmx_fio_do_int(fio
,fepvals
->init_fep_state
);
331 gmx_fio_do_double(fio
,fepvals
->init_lambda
);
332 gmx_fio_do_double(fio
,fepvals
->delta_lambda
);
334 else if (file_version
>= 59) {
335 gmx_fio_do_double(fio
,fepvals
->init_lambda
);
336 gmx_fio_do_double(fio
,fepvals
->delta_lambda
);
338 gmx_fio_do_real(fio
,rdum
);
339 fepvals
->init_lambda
= rdum
;
340 gmx_fio_do_real(fio
,rdum
);
341 fepvals
->delta_lambda
= rdum
;
343 if (file_version
>= 79)
345 gmx_fio_do_int(fio
,fepvals
->n_lambda
);
348 snew(fepvals
->all_lambda
,efptNR
);
350 for (g
=0;g
<efptNR
;g
++)
352 if (fepvals
->n_lambda
> 0) {
355 snew(fepvals
->all_lambda
[g
],fepvals
->n_lambda
);
357 bDum
=gmx_fio_ndo_double(fio
,fepvals
->all_lambda
[g
],fepvals
->n_lambda
);
358 bDum
=gmx_fio_ndo_int(fio
,fepvals
->separate_dvdl
,efptNR
);
360 else if (fepvals
->init_lambda
>= 0)
362 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
366 else if (file_version
>= 64)
368 gmx_fio_do_int(fio
,fepvals
->n_lambda
);
369 snew(fepvals
->all_lambda
,efptNR
);
372 snew(fepvals
->all_lambda
[efptFEP
],fepvals
->n_lambda
);
374 bDum
=gmx_fio_ndo_double(fio
,fepvals
->all_lambda
[efptFEP
],fepvals
->n_lambda
);
375 if (fepvals
->init_lambda
>= 0)
377 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
379 /* still allocate the all_lambda array's contents. */
380 for (g
=0;g
<efptNR
;g
++)
382 if (fepvals
->n_lambda
> 0) {
385 snew(fepvals
->all_lambda
[g
],fepvals
->n_lambda
);
392 fepvals
->n_lambda
= 0;
393 fepvals
->all_lambda
= NULL
;
394 if (fepvals
->init_lambda
>= 0)
396 fepvals
->separate_dvdl
[efptFEP
] = TRUE
;
399 if (file_version
>= 13)
401 gmx_fio_do_real(fio
,fepvals
->sc_alpha
);
405 fepvals
->sc_alpha
= 0;
407 if (file_version
>= 38)
409 gmx_fio_do_int(fio
,fepvals
->sc_power
);
413 fepvals
->sc_power
= 2;
415 if (file_version
>= 79)
417 gmx_fio_do_real(fio
,fepvals
->sc_r_power
);
421 fepvals
->sc_r_power
= 6.0;
423 if (file_version
>= 15)
425 gmx_fio_do_real(fio
,fepvals
->sc_sigma
);
429 fepvals
->sc_sigma
= 0.3;
433 if (file_version
>= 71)
435 fepvals
->sc_sigma_min
= fepvals
->sc_sigma
;
439 fepvals
->sc_sigma_min
= 0;
442 if (file_version
>= 79)
444 gmx_fio_do_int(fio
,fepvals
->bScCoul
);
448 fepvals
->bScCoul
= TRUE
;
450 if (file_version
>= 64) {
451 gmx_fio_do_int(fio
,fepvals
->nstdhdl
);
453 fepvals
->nstdhdl
= 1;
456 if (file_version
>= 73)
458 gmx_fio_do_int(fio
, fepvals
->separate_dhdl_file
);
459 gmx_fio_do_int(fio
, fepvals
->dhdl_derivatives
);
463 fepvals
->separate_dhdl_file
= esepdhdlfileYES
;
464 fepvals
->dhdl_derivatives
= edhdlderivativesYES
;
466 if (file_version
>= 71)
468 gmx_fio_do_int(fio
,fepvals
->dh_hist_size
);
469 gmx_fio_do_double(fio
,fepvals
->dh_hist_spacing
);
473 fepvals
->dh_hist_size
= 0;
474 fepvals
->dh_hist_spacing
= 0.1;
476 if (file_version
>= 79)
478 gmx_fio_do_int(fio
,fepvals
->bPrintEnergy
);
482 fepvals
->bPrintEnergy
= FALSE
;
486 static void do_pull(t_fileio
*fio
, t_pull
*pull
,gmx_bool bRead
, int file_version
)
490 gmx_fio_do_int(fio
,pull
->ngrp
);
491 gmx_fio_do_int(fio
,pull
->eGeom
);
492 gmx_fio_do_ivec(fio
,pull
->dim
);
493 gmx_fio_do_real(fio
,pull
->cyl_r1
);
494 gmx_fio_do_real(fio
,pull
->cyl_r0
);
495 gmx_fio_do_real(fio
,pull
->constr_tol
);
496 gmx_fio_do_int(fio
,pull
->nstxout
);
497 gmx_fio_do_int(fio
,pull
->nstfout
);
499 snew(pull
->grp
,pull
->ngrp
+1);
500 for(g
=0; g
<pull
->ngrp
+1; g
++)
501 do_pullgrp(fio
,&pull
->grp
[g
],bRead
,file_version
);
505 static void do_rotgrp(t_fileio
*fio
, t_rotgrp
*rotg
,gmx_bool bRead
, int file_version
)
510 gmx_fio_do_int(fio
,rotg
->eType
);
511 gmx_fio_do_int(fio
,rotg
->bMassW
);
512 gmx_fio_do_int(fio
,rotg
->nat
);
514 snew(rotg
->ind
,rotg
->nat
);
515 gmx_fio_ndo_int(fio
,rotg
->ind
,rotg
->nat
);
517 snew(rotg
->x_ref
,rotg
->nat
);
518 gmx_fio_ndo_rvec(fio
,rotg
->x_ref
,rotg
->nat
);
519 gmx_fio_do_rvec(fio
,rotg
->vec
);
520 gmx_fio_do_rvec(fio
,rotg
->pivot
);
521 gmx_fio_do_real(fio
,rotg
->rate
);
522 gmx_fio_do_real(fio
,rotg
->k
);
523 gmx_fio_do_real(fio
,rotg
->slab_dist
);
524 gmx_fio_do_real(fio
,rotg
->min_gaussian
);
525 gmx_fio_do_real(fio
,rotg
->eps
);
526 gmx_fio_do_int(fio
,rotg
->eFittype
);
527 gmx_fio_do_int(fio
,rotg
->PotAngle_nstep
);
528 gmx_fio_do_real(fio
,rotg
->PotAngle_step
);
531 static void do_rot(t_fileio
*fio
, t_rot
*rot
,gmx_bool bRead
, int file_version
)
535 gmx_fio_do_int(fio
,rot
->ngrp
);
536 gmx_fio_do_int(fio
,rot
->nstrout
);
537 gmx_fio_do_int(fio
,rot
->nstsout
);
539 snew(rot
->grp
,rot
->ngrp
);
540 for(g
=0; g
<rot
->ngrp
; g
++)
541 do_rotgrp(fio
, &rot
->grp
[g
],bRead
,file_version
);
545 static void do_inputrec(t_fileio
*fio
, t_inputrec
*ir
,gmx_bool bRead
,
546 int file_version
, real
*fudgeQQ
)
548 int i
,j
,k
,*tmp
,idum
=0;
553 real zerotemptime
,finish_t
,init_temp
,finish_temp
;
555 if (file_version
!= tpx_version
)
557 /* Give a warning about features that are not accessible */
558 fprintf(stderr
,"Note: file tpx version %d, software tpx version %d\n",
559 file_version
,tpx_version
);
567 if (file_version
== 0)
572 /* Basic inputrec stuff */
573 gmx_fio_do_int(fio
,ir
->eI
);
574 if (file_version
>= 62) {
575 gmx_fio_do_gmx_large_int(fio
, ir
->nsteps
);
577 gmx_fio_do_int(fio
,idum
);
580 if(file_version
> 25) {
581 if (file_version
>= 62) {
582 gmx_fio_do_gmx_large_int(fio
, ir
->init_step
);
584 gmx_fio_do_int(fio
,idum
);
585 ir
->init_step
= idum
;
591 if(file_version
>= 58)
592 gmx_fio_do_int(fio
,ir
->simulation_part
);
594 ir
->simulation_part
=1;
596 if (file_version
>= 67) {
597 gmx_fio_do_int(fio
,ir
->nstcalcenergy
);
599 ir
->nstcalcenergy
= 1;
601 if (file_version
< 53) {
602 /* The pbc info has been moved out of do_inputrec,
603 * since we always want it, also without reading the inputrec.
605 gmx_fio_do_int(fio
,ir
->ePBC
);
606 if ((file_version
<= 15) && (ir
->ePBC
== 2))
608 if (file_version
>= 45) {
609 gmx_fio_do_int(fio
,ir
->bPeriodicMols
);
613 ir
->bPeriodicMols
= TRUE
;
615 ir
->bPeriodicMols
= FALSE
;
619 if (file_version
>= 80)
621 gmx_fio_do_int(fio
,ir
->cutoff_scheme
);
625 ir
->cutoff_scheme
= ecutsGROUP
;
627 gmx_fio_do_int(fio
,ir
->ns_type
);
628 gmx_fio_do_int(fio
,ir
->nstlist
);
629 gmx_fio_do_int(fio
,ir
->ndelta
);
630 if (file_version
< 41) {
631 gmx_fio_do_int(fio
,idum
);
632 gmx_fio_do_int(fio
,idum
);
634 if (file_version
>= 45)
635 gmx_fio_do_real(fio
,ir
->rtpi
);
638 gmx_fio_do_int(fio
,ir
->nstcomm
);
639 if (file_version
> 34)
640 gmx_fio_do_int(fio
,ir
->comm_mode
);
641 else if (ir
->nstcomm
< 0)
642 ir
->comm_mode
= ecmANGULAR
;
644 ir
->comm_mode
= ecmLINEAR
;
645 ir
->nstcomm
= abs(ir
->nstcomm
);
647 if(file_version
> 25)
648 gmx_fio_do_int(fio
,ir
->nstcheckpoint
);
652 gmx_fio_do_int(fio
,ir
->nstcgsteep
);
655 gmx_fio_do_int(fio
,ir
->nbfgscorr
);
659 gmx_fio_do_int(fio
,ir
->nstlog
);
660 gmx_fio_do_int(fio
,ir
->nstxout
);
661 gmx_fio_do_int(fio
,ir
->nstvout
);
662 gmx_fio_do_int(fio
,ir
->nstfout
);
663 gmx_fio_do_int(fio
,ir
->nstenergy
);
664 gmx_fio_do_int(fio
,ir
->nstxtcout
);
665 if (file_version
>= 59) {
666 gmx_fio_do_double(fio
,ir
->init_t
);
667 gmx_fio_do_double(fio
,ir
->delta_t
);
669 gmx_fio_do_real(fio
,rdum
);
671 gmx_fio_do_real(fio
,rdum
);
674 gmx_fio_do_real(fio
,ir
->xtcprec
);
675 if (file_version
< 19) {
676 gmx_fio_do_int(fio
,idum
);
677 gmx_fio_do_int(fio
,idum
);
679 if(file_version
< 18)
680 gmx_fio_do_int(fio
,idum
);
681 if (file_version
>= 80) {
682 gmx_fio_do_real(fio
,ir
->verletbuf_drift
);
684 ir
->verletbuf_drift
= 0;
686 gmx_fio_do_real(fio
,ir
->rlist
);
687 if (file_version
>= 67) {
688 gmx_fio_do_real(fio
,ir
->rlistlong
);
690 gmx_fio_do_int(fio
,ir
->coulombtype
);
691 if (file_version
< 32 && ir
->coulombtype
== eelRF
)
692 ir
->coulombtype
= eelRF_NEC
;
693 gmx_fio_do_real(fio
,ir
->rcoulomb_switch
);
694 gmx_fio_do_real(fio
,ir
->rcoulomb
);
695 gmx_fio_do_int(fio
,ir
->vdwtype
);
696 gmx_fio_do_real(fio
,ir
->rvdw_switch
);
697 gmx_fio_do_real(fio
,ir
->rvdw
);
698 if (file_version
< 67) {
699 ir
->rlistlong
= max_cutoff(ir
->rlist
,max_cutoff(ir
->rvdw
,ir
->rcoulomb
));
701 gmx_fio_do_int(fio
,ir
->eDispCorr
);
702 gmx_fio_do_real(fio
,ir
->epsilon_r
);
703 if (file_version
>= 37) {
704 gmx_fio_do_real(fio
,ir
->epsilon_rf
);
706 if (EEL_RF(ir
->coulombtype
)) {
707 ir
->epsilon_rf
= ir
->epsilon_r
;
710 ir
->epsilon_rf
= 1.0;
713 if (file_version
>= 29)
714 gmx_fio_do_real(fio
,ir
->tabext
);
718 if(file_version
> 25) {
719 gmx_fio_do_int(fio
,ir
->gb_algorithm
);
720 gmx_fio_do_int(fio
,ir
->nstgbradii
);
721 gmx_fio_do_real(fio
,ir
->rgbradii
);
722 gmx_fio_do_real(fio
,ir
->gb_saltconc
);
723 gmx_fio_do_int(fio
,ir
->implicit_solvent
);
725 ir
->gb_algorithm
=egbSTILL
;
729 ir
->implicit_solvent
=eisNO
;
733 gmx_fio_do_real(fio
,ir
->gb_epsilon_solvent
);
734 gmx_fio_do_real(fio
,ir
->gb_obc_alpha
);
735 gmx_fio_do_real(fio
,ir
->gb_obc_beta
);
736 gmx_fio_do_real(fio
,ir
->gb_obc_gamma
);
739 gmx_fio_do_real(fio
,ir
->gb_dielectric_offset
);
740 gmx_fio_do_int(fio
,ir
->sa_algorithm
);
744 ir
->gb_dielectric_offset
= 0.009;
745 ir
->sa_algorithm
= esaAPPROX
;
747 gmx_fio_do_real(fio
,ir
->sa_surface_tension
);
749 /* Override sa_surface_tension if it is not changed in the mpd-file */
750 if(ir
->sa_surface_tension
<0)
752 if(ir
->gb_algorithm
==egbSTILL
)
754 ir
->sa_surface_tension
= 0.0049 * 100 * CAL2JOULE
;
756 else if(ir
->gb_algorithm
==egbHCT
|| ir
->gb_algorithm
==egbOBC
)
758 ir
->sa_surface_tension
= 0.0054 * 100 * CAL2JOULE
;
765 /* Better use sensible values than insane (0.0) ones... */
766 ir
->gb_epsilon_solvent
= 80;
767 ir
->gb_obc_alpha
= 1.0;
768 ir
->gb_obc_beta
= 0.8;
769 ir
->gb_obc_gamma
= 4.85;
770 ir
->sa_surface_tension
= 2.092;
774 if (file_version
>= 80)
776 gmx_fio_do_real(fio
,ir
->fourier_spacing
);
780 ir
->fourier_spacing
= 0.0;
782 gmx_fio_do_int(fio
,ir
->nkx
);
783 gmx_fio_do_int(fio
,ir
->nky
);
784 gmx_fio_do_int(fio
,ir
->nkz
);
785 gmx_fio_do_int(fio
,ir
->pme_order
);
786 gmx_fio_do_real(fio
,ir
->ewald_rtol
);
788 if (file_version
>=24)
789 gmx_fio_do_int(fio
,ir
->ewald_geometry
);
791 ir
->ewald_geometry
=eewg3D
;
793 if (file_version
<=17) {
794 ir
->epsilon_surface
=0;
795 if (file_version
==17)
796 gmx_fio_do_int(fio
,idum
);
799 gmx_fio_do_real(fio
,ir
->epsilon_surface
);
801 gmx_fio_do_gmx_bool(fio
,ir
->bOptFFT
);
803 gmx_fio_do_gmx_bool(fio
,ir
->bContinuation
);
804 gmx_fio_do_int(fio
,ir
->etc
);
805 /* before version 18, ir->etc was a gmx_bool (ir->btc),
806 * but the values 0 and 1 still mean no and
807 * berendsen temperature coupling, respectively.
809 if (file_version
>= 79) {
810 gmx_fio_do_gmx_bool(fio
,ir
->bPrintNHChains
);
812 if (file_version
>= 71)
814 gmx_fio_do_int(fio
,ir
->nsttcouple
);
818 ir
->nsttcouple
= ir
->nstcalcenergy
;
820 if (file_version
<= 15)
822 gmx_fio_do_int(fio
,idum
);
824 if (file_version
<=17)
826 gmx_fio_do_int(fio
,ir
->epct
);
827 if (file_version
<= 15)
831 ir
->epct
= epctSURFACETENSION
;
833 gmx_fio_do_int(fio
,idum
);
836 /* we have removed the NO alternative at the beginning */
840 ir
->epct
=epctISOTROPIC
;
844 ir
->epc
=epcBERENDSEN
;
849 gmx_fio_do_int(fio
,ir
->epc
);
850 gmx_fio_do_int(fio
,ir
->epct
);
852 if (file_version
>= 71)
854 gmx_fio_do_int(fio
,ir
->nstpcouple
);
858 ir
->nstpcouple
= ir
->nstcalcenergy
;
860 gmx_fio_do_real(fio
,ir
->tau_p
);
861 if (file_version
<= 15) {
862 gmx_fio_do_rvec(fio
,vdum
);
863 clear_mat(ir
->ref_p
);
865 ir
->ref_p
[i
][i
] = vdum
[i
];
867 gmx_fio_do_rvec(fio
,ir
->ref_p
[XX
]);
868 gmx_fio_do_rvec(fio
,ir
->ref_p
[YY
]);
869 gmx_fio_do_rvec(fio
,ir
->ref_p
[ZZ
]);
871 if (file_version
<= 15) {
872 gmx_fio_do_rvec(fio
,vdum
);
873 clear_mat(ir
->compress
);
875 ir
->compress
[i
][i
] = vdum
[i
];
878 gmx_fio_do_rvec(fio
,ir
->compress
[XX
]);
879 gmx_fio_do_rvec(fio
,ir
->compress
[YY
]);
880 gmx_fio_do_rvec(fio
,ir
->compress
[ZZ
]);
882 if (file_version
>= 47) {
883 gmx_fio_do_int(fio
,ir
->refcoord_scaling
);
884 gmx_fio_do_rvec(fio
,ir
->posres_com
);
885 gmx_fio_do_rvec(fio
,ir
->posres_comB
);
887 ir
->refcoord_scaling
= erscNO
;
888 clear_rvec(ir
->posres_com
);
889 clear_rvec(ir
->posres_comB
);
891 if((file_version
> 25) && (file_version
< 79))
892 gmx_fio_do_int(fio
,ir
->andersen_seed
);
895 if(file_version
< 26) {
896 gmx_fio_do_gmx_bool(fio
,bSimAnn
);
897 gmx_fio_do_real(fio
,zerotemptime
);
900 if (file_version
< 37)
901 gmx_fio_do_real(fio
,rdum
);
903 gmx_fio_do_real(fio
,ir
->shake_tol
);
904 if (file_version
< 54)
905 gmx_fio_do_real(fio
,*fudgeQQ
);
907 gmx_fio_do_int(fio
,ir
->efep
);
908 if (file_version
<= 14 && ir
->efep
!= efepNO
)
912 do_fepvals(fio
,ir
->fepvals
,bRead
,file_version
);
914 if (file_version
>= 79)
916 gmx_fio_do_gmx_bool(fio
,ir
->bSimTemp
);
924 ir
->bSimTemp
= FALSE
;
928 do_simtempvals(fio
,ir
->simtempvals
,ir
->fepvals
->n_lambda
,bRead
,file_version
);
931 if (file_version
>= 79)
933 gmx_fio_do_gmx_bool(fio
,ir
->bExpanded
);
936 ir
->bExpanded
= TRUE
;
940 ir
->bExpanded
= FALSE
;
945 do_expandedvals(fio
,ir
->expandedvals
,ir
->fepvals
->n_lambda
,bRead
,file_version
);
947 if (file_version
>= 57) {
948 gmx_fio_do_int(fio
,ir
->eDisre
);
950 gmx_fio_do_int(fio
,ir
->eDisreWeighting
);
951 if (file_version
< 22) {
952 if (ir
->eDisreWeighting
== 0)
953 ir
->eDisreWeighting
= edrwEqual
;
955 ir
->eDisreWeighting
= edrwConservative
;
957 gmx_fio_do_gmx_bool(fio
,ir
->bDisreMixed
);
958 gmx_fio_do_real(fio
,ir
->dr_fc
);
959 gmx_fio_do_real(fio
,ir
->dr_tau
);
960 gmx_fio_do_int(fio
,ir
->nstdisreout
);
961 if (file_version
>= 22) {
962 gmx_fio_do_real(fio
,ir
->orires_fc
);
963 gmx_fio_do_real(fio
,ir
->orires_tau
);
964 gmx_fio_do_int(fio
,ir
->nstorireout
);
970 if(file_version
>= 26 && file_version
< 79) {
971 gmx_fio_do_real(fio
,ir
->dihre_fc
);
972 if (file_version
< 56)
974 gmx_fio_do_real(fio
,rdum
);
975 gmx_fio_do_int(fio
,idum
);
981 gmx_fio_do_real(fio
,ir
->em_stepsize
);
982 gmx_fio_do_real(fio
,ir
->em_tol
);
983 if (file_version
>= 22)
984 gmx_fio_do_gmx_bool(fio
,ir
->bShakeSOR
);
986 ir
->bShakeSOR
= TRUE
;
987 if (file_version
>= 11)
988 gmx_fio_do_int(fio
,ir
->niter
);
991 fprintf(stderr
,"Note: niter not in run input file, setting it to %d\n",
994 if (file_version
>= 21)
995 gmx_fio_do_real(fio
,ir
->fc_stepsize
);
998 gmx_fio_do_int(fio
,ir
->eConstrAlg
);
999 gmx_fio_do_int(fio
,ir
->nProjOrder
);
1000 gmx_fio_do_real(fio
,ir
->LincsWarnAngle
);
1001 if (file_version
<= 14)
1002 gmx_fio_do_int(fio
,idum
);
1003 if (file_version
>=26)
1004 gmx_fio_do_int(fio
,ir
->nLincsIter
);
1007 fprintf(stderr
,"Note: nLincsIter not in run input file, setting it to %d\n",
1010 if (file_version
< 33)
1011 gmx_fio_do_real(fio
,bd_temp
);
1012 gmx_fio_do_real(fio
,ir
->bd_fric
);
1013 gmx_fio_do_int(fio
,ir
->ld_seed
);
1014 if (file_version
>= 33) {
1015 for(i
=0; i
<DIM
; i
++)
1016 gmx_fio_do_rvec(fio
,ir
->deform
[i
]);
1018 for(i
=0; i
<DIM
; i
++)
1019 clear_rvec(ir
->deform
[i
]);
1021 if (file_version
>= 14)
1022 gmx_fio_do_real(fio
,ir
->cos_accel
);
1025 gmx_fio_do_int(fio
,ir
->userint1
);
1026 gmx_fio_do_int(fio
,ir
->userint2
);
1027 gmx_fio_do_int(fio
,ir
->userint3
);
1028 gmx_fio_do_int(fio
,ir
->userint4
);
1029 gmx_fio_do_real(fio
,ir
->userreal1
);
1030 gmx_fio_do_real(fio
,ir
->userreal2
);
1031 gmx_fio_do_real(fio
,ir
->userreal3
);
1032 gmx_fio_do_real(fio
,ir
->userreal4
);
1035 if (file_version
>= 77) {
1036 gmx_fio_do_gmx_bool(fio
,ir
->bAdress
);
1038 if (bRead
) snew(ir
->adress
, 1);
1039 gmx_fio_do_int(fio
,ir
->adress
->type
);
1040 gmx_fio_do_real(fio
,ir
->adress
->const_wf
);
1041 gmx_fio_do_real(fio
,ir
->adress
->ex_width
);
1042 gmx_fio_do_real(fio
,ir
->adress
->hy_width
);
1043 gmx_fio_do_int(fio
,ir
->adress
->icor
);
1044 gmx_fio_do_int(fio
,ir
->adress
->site
);
1045 gmx_fio_do_rvec(fio
,ir
->adress
->refs
);
1046 gmx_fio_do_int(fio
,ir
->adress
->n_tf_grps
);
1047 gmx_fio_do_real(fio
, ir
->adress
->ex_forcecap
);
1048 gmx_fio_do_int(fio
, ir
->adress
->n_energy_grps
);
1049 gmx_fio_do_int(fio
,ir
->adress
->do_hybridpairs
);
1051 if (bRead
)snew(ir
->adress
->tf_table_index
,ir
->adress
->n_tf_grps
);
1052 if (ir
->adress
->n_tf_grps
> 0) {
1053 bDum
=gmx_fio_ndo_int(fio
,ir
->adress
->tf_table_index
,ir
->adress
->n_tf_grps
);
1055 if (bRead
)snew(ir
->adress
->group_explicit
,ir
->adress
->n_energy_grps
);
1056 if (ir
->adress
->n_energy_grps
> 0) {
1057 bDum
=gmx_fio_ndo_int(fio
, ir
->adress
->group_explicit
,ir
->adress
->n_energy_grps
);
1061 ir
->bAdress
= FALSE
;
1065 if (file_version
>= 48) {
1066 gmx_fio_do_int(fio
,ir
->ePull
);
1067 if (ir
->ePull
!= epullNO
) {
1070 do_pull(fio
, ir
->pull
,bRead
,file_version
);
1073 ir
->ePull
= epullNO
;
1076 /* Enforced rotation */
1077 if (file_version
>= 74) {
1078 gmx_fio_do_int(fio
,ir
->bRot
);
1079 if (ir
->bRot
== TRUE
) {
1082 do_rot(fio
, ir
->rot
,bRead
,file_version
);
1089 gmx_fio_do_int(fio
,ir
->opts
.ngtc
);
1090 if (file_version
>= 69) {
1091 gmx_fio_do_int(fio
,ir
->opts
.nhchainlength
);
1093 ir
->opts
.nhchainlength
= 1;
1095 gmx_fio_do_int(fio
,ir
->opts
.ngacc
);
1096 gmx_fio_do_int(fio
,ir
->opts
.ngfrz
);
1097 gmx_fio_do_int(fio
,ir
->opts
.ngener
);
1100 snew(ir
->opts
.nrdf
, ir
->opts
.ngtc
);
1101 snew(ir
->opts
.ref_t
, ir
->opts
.ngtc
);
1102 snew(ir
->opts
.annealing
, ir
->opts
.ngtc
);
1103 snew(ir
->opts
.anneal_npoints
, ir
->opts
.ngtc
);
1104 snew(ir
->opts
.anneal_time
, ir
->opts
.ngtc
);
1105 snew(ir
->opts
.anneal_temp
, ir
->opts
.ngtc
);
1106 snew(ir
->opts
.tau_t
, ir
->opts
.ngtc
);
1107 snew(ir
->opts
.nFreeze
,ir
->opts
.ngfrz
);
1108 snew(ir
->opts
.acc
, ir
->opts
.ngacc
);
1109 snew(ir
->opts
.egp_flags
,ir
->opts
.ngener
*ir
->opts
.ngener
);
1111 if (ir
->opts
.ngtc
> 0) {
1112 if (bRead
&& file_version
<13) {
1113 snew(tmp
,ir
->opts
.ngtc
);
1114 bDum
=gmx_fio_ndo_int(fio
,tmp
, ir
->opts
.ngtc
);
1115 for(i
=0; i
<ir
->opts
.ngtc
; i
++)
1116 ir
->opts
.nrdf
[i
] = tmp
[i
];
1119 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.nrdf
, ir
->opts
.ngtc
);
1121 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.ref_t
,ir
->opts
.ngtc
);
1122 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.tau_t
,ir
->opts
.ngtc
);
1123 if (file_version
<33 && ir
->eI
==eiBD
) {
1124 for(i
=0; i
<ir
->opts
.ngtc
; i
++)
1125 ir
->opts
.tau_t
[i
] = bd_temp
;
1128 if (ir
->opts
.ngfrz
> 0)
1129 bDum
=gmx_fio_ndo_ivec(fio
,ir
->opts
.nFreeze
,ir
->opts
.ngfrz
);
1130 if (ir
->opts
.ngacc
> 0)
1131 gmx_fio_ndo_rvec(fio
,ir
->opts
.acc
,ir
->opts
.ngacc
);
1132 if (file_version
>= 12)
1133 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.egp_flags
,
1134 ir
->opts
.ngener
*ir
->opts
.ngener
);
1136 if(bRead
&& file_version
< 26) {
1137 for(i
=0;i
<ir
->opts
.ngtc
;i
++) {
1139 ir
->opts
.annealing
[i
] = eannSINGLE
;
1140 ir
->opts
.anneal_npoints
[i
] = 2;
1141 snew(ir
->opts
.anneal_time
[i
],2);
1142 snew(ir
->opts
.anneal_temp
[i
],2);
1143 /* calculate the starting/ending temperatures from reft, zerotemptime, and nsteps */
1144 finish_t
= ir
->init_t
+ ir
->nsteps
* ir
->delta_t
;
1145 init_temp
= ir
->opts
.ref_t
[i
]*(1-ir
->init_t
/zerotemptime
);
1146 finish_temp
= ir
->opts
.ref_t
[i
]*(1-finish_t
/zerotemptime
);
1147 ir
->opts
.anneal_time
[i
][0] = ir
->init_t
;
1148 ir
->opts
.anneal_time
[i
][1] = finish_t
;
1149 ir
->opts
.anneal_temp
[i
][0] = init_temp
;
1150 ir
->opts
.anneal_temp
[i
][1] = finish_temp
;
1152 ir
->opts
.annealing
[i
] = eannNO
;
1153 ir
->opts
.anneal_npoints
[i
] = 0;
1157 /* file version 26 or later */
1158 /* First read the lists with annealing and npoints for each group */
1159 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.annealing
,ir
->opts
.ngtc
);
1160 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.anneal_npoints
,ir
->opts
.ngtc
);
1161 for(j
=0;j
<(ir
->opts
.ngtc
);j
++) {
1162 k
=ir
->opts
.anneal_npoints
[j
];
1164 snew(ir
->opts
.anneal_time
[j
],k
);
1165 snew(ir
->opts
.anneal_temp
[j
],k
);
1167 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.anneal_time
[j
],k
);
1168 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.anneal_temp
[j
],k
);
1172 if (file_version
>= 45) {
1173 gmx_fio_do_int(fio
,ir
->nwall
);
1174 gmx_fio_do_int(fio
,ir
->wall_type
);
1175 if (file_version
>= 50)
1176 gmx_fio_do_real(fio
,ir
->wall_r_linpot
);
1178 ir
->wall_r_linpot
= -1;
1179 gmx_fio_do_int(fio
,ir
->wall_atomtype
[0]);
1180 gmx_fio_do_int(fio
,ir
->wall_atomtype
[1]);
1181 gmx_fio_do_real(fio
,ir
->wall_density
[0]);
1182 gmx_fio_do_real(fio
,ir
->wall_density
[1]);
1183 gmx_fio_do_real(fio
,ir
->wall_ewald_zfac
);
1187 ir
->wall_atomtype
[0] = -1;
1188 ir
->wall_atomtype
[1] = -1;
1189 ir
->wall_density
[0] = 0;
1190 ir
->wall_density
[1] = 0;
1191 ir
->wall_ewald_zfac
= 3;
1193 /* Cosine stuff for electric fields */
1194 for(j
=0; (j
<DIM
); j
++) {
1195 gmx_fio_do_int(fio
,ir
->ex
[j
].n
);
1196 gmx_fio_do_int(fio
,ir
->et
[j
].n
);
1198 snew(ir
->ex
[j
].a
, ir
->ex
[j
].n
);
1199 snew(ir
->ex
[j
].phi
,ir
->ex
[j
].n
);
1200 snew(ir
->et
[j
].a
, ir
->et
[j
].n
);
1201 snew(ir
->et
[j
].phi
,ir
->et
[j
].n
);
1203 bDum
=gmx_fio_ndo_real(fio
,ir
->ex
[j
].a
, ir
->ex
[j
].n
);
1204 bDum
=gmx_fio_ndo_real(fio
,ir
->ex
[j
].phi
,ir
->ex
[j
].n
);
1205 bDum
=gmx_fio_ndo_real(fio
,ir
->et
[j
].a
, ir
->et
[j
].n
);
1206 bDum
=gmx_fio_ndo_real(fio
,ir
->et
[j
].phi
,ir
->et
[j
].n
);
1210 if(file_version
>=39){
1211 gmx_fio_do_gmx_bool(fio
,ir
->bQMMM
);
1212 gmx_fio_do_int(fio
,ir
->QMMMscheme
);
1213 gmx_fio_do_real(fio
,ir
->scalefactor
);
1214 gmx_fio_do_int(fio
,ir
->opts
.ngQM
);
1216 snew(ir
->opts
.QMmethod
, ir
->opts
.ngQM
);
1217 snew(ir
->opts
.QMbasis
, ir
->opts
.ngQM
);
1218 snew(ir
->opts
.QMcharge
, ir
->opts
.ngQM
);
1219 snew(ir
->opts
.QMmult
, ir
->opts
.ngQM
);
1220 snew(ir
->opts
.bSH
, ir
->opts
.ngQM
);
1221 snew(ir
->opts
.CASorbitals
, ir
->opts
.ngQM
);
1222 snew(ir
->opts
.CASelectrons
,ir
->opts
.ngQM
);
1223 snew(ir
->opts
.SAon
, ir
->opts
.ngQM
);
1224 snew(ir
->opts
.SAoff
, ir
->opts
.ngQM
);
1225 snew(ir
->opts
.SAsteps
, ir
->opts
.ngQM
);
1226 snew(ir
->opts
.bOPT
, ir
->opts
.ngQM
);
1227 snew(ir
->opts
.bTS
, ir
->opts
.ngQM
);
1229 if (ir
->opts
.ngQM
> 0) {
1230 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMmethod
,ir
->opts
.ngQM
);
1231 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMbasis
,ir
->opts
.ngQM
);
1232 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMcharge
,ir
->opts
.ngQM
);
1233 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMmult
,ir
->opts
.ngQM
);
1234 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bSH
,ir
->opts
.ngQM
);
1235 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.CASorbitals
,ir
->opts
.ngQM
);
1236 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.CASelectrons
,ir
->opts
.ngQM
);
1237 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.SAon
,ir
->opts
.ngQM
);
1238 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.SAoff
,ir
->opts
.ngQM
);
1239 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.SAsteps
,ir
->opts
.ngQM
);
1240 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bOPT
,ir
->opts
.ngQM
);
1241 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bTS
,ir
->opts
.ngQM
);
1243 /* end of QMMM stuff */
1248 static void do_harm(t_fileio
*fio
, t_iparams
*iparams
,gmx_bool bRead
)
1250 gmx_fio_do_real(fio
,iparams
->harmonic
.rA
);
1251 gmx_fio_do_real(fio
,iparams
->harmonic
.krA
);
1252 gmx_fio_do_real(fio
,iparams
->harmonic
.rB
);
1253 gmx_fio_do_real(fio
,iparams
->harmonic
.krB
);
1256 void do_iparams(t_fileio
*fio
, t_functype ftype
,t_iparams
*iparams
,
1257 gmx_bool bRead
, int file_version
)
1264 gmx_fio_set_comment(fio
, interaction_function
[ftype
].name
);
1272 do_harm(fio
, iparams
,bRead
);
1273 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && bRead
) {
1274 /* Correct incorrect storage of parameters */
1275 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
1276 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
1279 case F_LINEAR_ANGLES
:
1280 gmx_fio_do_real(fio
,iparams
->linangle
.klinA
);
1281 gmx_fio_do_real(fio
,iparams
->linangle
.aA
);
1282 gmx_fio_do_real(fio
,iparams
->linangle
.klinB
);
1283 gmx_fio_do_real(fio
,iparams
->linangle
.aB
);
1286 gmx_fio_do_real(fio
,iparams
->fene
.bm
);
1287 gmx_fio_do_real(fio
,iparams
->fene
.kb
);
1290 gmx_fio_do_real(fio
,iparams
->restraint
.lowA
);
1291 gmx_fio_do_real(fio
,iparams
->restraint
.up1A
);
1292 gmx_fio_do_real(fio
,iparams
->restraint
.up2A
);
1293 gmx_fio_do_real(fio
,iparams
->restraint
.kA
);
1294 gmx_fio_do_real(fio
,iparams
->restraint
.lowB
);
1295 gmx_fio_do_real(fio
,iparams
->restraint
.up1B
);
1296 gmx_fio_do_real(fio
,iparams
->restraint
.up2B
);
1297 gmx_fio_do_real(fio
,iparams
->restraint
.kB
);
1303 gmx_fio_do_real(fio
,iparams
->tab
.kA
);
1304 gmx_fio_do_int(fio
,iparams
->tab
.table
);
1305 gmx_fio_do_real(fio
,iparams
->tab
.kB
);
1307 case F_CROSS_BOND_BONDS
:
1308 gmx_fio_do_real(fio
,iparams
->cross_bb
.r1e
);
1309 gmx_fio_do_real(fio
,iparams
->cross_bb
.r2e
);
1310 gmx_fio_do_real(fio
,iparams
->cross_bb
.krr
);
1312 case F_CROSS_BOND_ANGLES
:
1313 gmx_fio_do_real(fio
,iparams
->cross_ba
.r1e
);
1314 gmx_fio_do_real(fio
,iparams
->cross_ba
.r2e
);
1315 gmx_fio_do_real(fio
,iparams
->cross_ba
.r3e
);
1316 gmx_fio_do_real(fio
,iparams
->cross_ba
.krt
);
1318 case F_UREY_BRADLEY
:
1319 gmx_fio_do_real(fio
,iparams
->u_b
.thetaA
);
1320 gmx_fio_do_real(fio
,iparams
->u_b
.kthetaA
);
1321 gmx_fio_do_real(fio
,iparams
->u_b
.r13A
);
1322 gmx_fio_do_real(fio
,iparams
->u_b
.kUBA
);
1323 if (file_version
>= 79) {
1324 gmx_fio_do_real(fio
,iparams
->u_b
.thetaB
);
1325 gmx_fio_do_real(fio
,iparams
->u_b
.kthetaB
);
1326 gmx_fio_do_real(fio
,iparams
->u_b
.r13B
);
1327 gmx_fio_do_real(fio
,iparams
->u_b
.kUBB
);
1329 iparams
->u_b
.thetaB
=iparams
->u_b
.thetaA
;
1330 iparams
->u_b
.kthetaB
=iparams
->u_b
.kthetaA
;
1331 iparams
->u_b
.r13B
=iparams
->u_b
.r13A
;
1332 iparams
->u_b
.kUBB
=iparams
->u_b
.kUBA
;
1335 case F_QUARTIC_ANGLES
:
1336 gmx_fio_do_real(fio
,iparams
->qangle
.theta
);
1337 bDum
=gmx_fio_ndo_real(fio
,iparams
->qangle
.c
,5);
1340 gmx_fio_do_real(fio
,iparams
->bham
.a
);
1341 gmx_fio_do_real(fio
,iparams
->bham
.b
);
1342 gmx_fio_do_real(fio
,iparams
->bham
.c
);
1345 gmx_fio_do_real(fio
,iparams
->morse
.b0A
);
1346 gmx_fio_do_real(fio
,iparams
->morse
.cbA
);
1347 gmx_fio_do_real(fio
,iparams
->morse
.betaA
);
1348 if (file_version
>= 79) {
1349 gmx_fio_do_real(fio
,iparams
->morse
.b0B
);
1350 gmx_fio_do_real(fio
,iparams
->morse
.cbB
);
1351 gmx_fio_do_real(fio
,iparams
->morse
.betaB
);
1353 iparams
->morse
.b0B
= iparams
->morse
.b0A
;
1354 iparams
->morse
.cbB
= iparams
->morse
.cbA
;
1355 iparams
->morse
.betaB
= iparams
->morse
.betaA
;
1359 gmx_fio_do_real(fio
,iparams
->cubic
.b0
);
1360 gmx_fio_do_real(fio
,iparams
->cubic
.kb
);
1361 gmx_fio_do_real(fio
,iparams
->cubic
.kcub
);
1365 case F_POLARIZATION
:
1366 gmx_fio_do_real(fio
,iparams
->polarize
.alpha
);
1369 gmx_fio_do_real(fio
,iparams
->anharm_polarize
.alpha
);
1370 gmx_fio_do_real(fio
,iparams
->anharm_polarize
.drcut
);
1371 gmx_fio_do_real(fio
,iparams
->anharm_polarize
.khyp
);
1374 if (file_version
< 31)
1375 gmx_fatal(FARGS
,"Old tpr files with water_polarization not supported. Make a new.");
1376 gmx_fio_do_real(fio
,iparams
->wpol
.al_x
);
1377 gmx_fio_do_real(fio
,iparams
->wpol
.al_y
);
1378 gmx_fio_do_real(fio
,iparams
->wpol
.al_z
);
1379 gmx_fio_do_real(fio
,iparams
->wpol
.rOH
);
1380 gmx_fio_do_real(fio
,iparams
->wpol
.rHH
);
1381 gmx_fio_do_real(fio
,iparams
->wpol
.rOD
);
1384 gmx_fio_do_real(fio
,iparams
->thole
.a
);
1385 gmx_fio_do_real(fio
,iparams
->thole
.alpha1
);
1386 gmx_fio_do_real(fio
,iparams
->thole
.alpha2
);
1387 gmx_fio_do_real(fio
,iparams
->thole
.rfac
);
1390 gmx_fio_do_real(fio
,iparams
->lj
.c6
);
1391 gmx_fio_do_real(fio
,iparams
->lj
.c12
);
1394 gmx_fio_do_real(fio
,iparams
->lj14
.c6A
);
1395 gmx_fio_do_real(fio
,iparams
->lj14
.c12A
);
1396 gmx_fio_do_real(fio
,iparams
->lj14
.c6B
);
1397 gmx_fio_do_real(fio
,iparams
->lj14
.c12B
);
1400 gmx_fio_do_real(fio
,iparams
->ljc14
.fqq
);
1401 gmx_fio_do_real(fio
,iparams
->ljc14
.qi
);
1402 gmx_fio_do_real(fio
,iparams
->ljc14
.qj
);
1403 gmx_fio_do_real(fio
,iparams
->ljc14
.c6
);
1404 gmx_fio_do_real(fio
,iparams
->ljc14
.c12
);
1406 case F_LJC_PAIRS_NB
:
1407 gmx_fio_do_real(fio
,iparams
->ljcnb
.qi
);
1408 gmx_fio_do_real(fio
,iparams
->ljcnb
.qj
);
1409 gmx_fio_do_real(fio
,iparams
->ljcnb
.c6
);
1410 gmx_fio_do_real(fio
,iparams
->ljcnb
.c12
);
1416 gmx_fio_do_real(fio
,iparams
->pdihs
.phiA
);
1417 gmx_fio_do_real(fio
,iparams
->pdihs
.cpA
);
1418 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && file_version
< 42) {
1419 /* Read the incorrectly stored multiplicity */
1420 gmx_fio_do_real(fio
,iparams
->harmonic
.rB
);
1421 gmx_fio_do_real(fio
,iparams
->harmonic
.krB
);
1422 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
1423 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
1425 gmx_fio_do_real(fio
,iparams
->pdihs
.phiB
);
1426 gmx_fio_do_real(fio
,iparams
->pdihs
.cpB
);
1427 gmx_fio_do_int(fio
,iparams
->pdihs
.mult
);
1431 gmx_fio_do_int(fio
,iparams
->disres
.label
);
1432 gmx_fio_do_int(fio
,iparams
->disres
.type
);
1433 gmx_fio_do_real(fio
,iparams
->disres
.low
);
1434 gmx_fio_do_real(fio
,iparams
->disres
.up1
);
1435 gmx_fio_do_real(fio
,iparams
->disres
.up2
);
1436 gmx_fio_do_real(fio
,iparams
->disres
.kfac
);
1439 gmx_fio_do_int(fio
,iparams
->orires
.ex
);
1440 gmx_fio_do_int(fio
,iparams
->orires
.label
);
1441 gmx_fio_do_int(fio
,iparams
->orires
.power
);
1442 gmx_fio_do_real(fio
,iparams
->orires
.c
);
1443 gmx_fio_do_real(fio
,iparams
->orires
.obs
);
1444 gmx_fio_do_real(fio
,iparams
->orires
.kfac
);
1447 gmx_fio_do_real(fio
,iparams
->dihres
.phiA
);
1448 gmx_fio_do_real(fio
,iparams
->dihres
.dphiA
);
1449 gmx_fio_do_real(fio
,iparams
->dihres
.kfacA
);
1450 if (file_version
>= 72) {
1451 gmx_fio_do_real(fio
,iparams
->dihres
.phiB
);
1452 gmx_fio_do_real(fio
,iparams
->dihres
.dphiB
);
1453 gmx_fio_do_real(fio
,iparams
->dihres
.kfacB
);
1455 iparams
->dihres
.phiB
=iparams
->dihres
.phiA
;
1456 iparams
->dihres
.dphiB
=iparams
->dihres
.dphiA
;
1457 iparams
->dihres
.kfacB
=iparams
->dihres
.kfacA
;
1461 gmx_fio_do_rvec(fio
,iparams
->posres
.pos0A
);
1462 gmx_fio_do_rvec(fio
,iparams
->posres
.fcA
);
1463 if (bRead
&& file_version
< 27) {
1464 copy_rvec(iparams
->posres
.pos0A
,iparams
->posres
.pos0B
);
1465 copy_rvec(iparams
->posres
.fcA
,iparams
->posres
.fcB
);
1467 gmx_fio_do_rvec(fio
,iparams
->posres
.pos0B
);
1468 gmx_fio_do_rvec(fio
,iparams
->posres
.fcB
);
1472 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcA
,NR_RBDIHS
);
1473 if(file_version
>=25)
1474 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcB
,NR_RBDIHS
);
1477 /* Fourier dihedrals are internally represented
1478 * as Ryckaert-Bellemans since those are faster to compute.
1480 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcA
, NR_RBDIHS
);
1481 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcB
, NR_RBDIHS
);
1485 gmx_fio_do_real(fio
,iparams
->constr
.dA
);
1486 gmx_fio_do_real(fio
,iparams
->constr
.dB
);
1489 gmx_fio_do_real(fio
,iparams
->settle
.doh
);
1490 gmx_fio_do_real(fio
,iparams
->settle
.dhh
);
1493 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1498 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1499 gmx_fio_do_real(fio
,iparams
->vsite
.b
);
1504 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1505 gmx_fio_do_real(fio
,iparams
->vsite
.b
);
1506 gmx_fio_do_real(fio
,iparams
->vsite
.c
);
1509 gmx_fio_do_int(fio
,iparams
->vsiten
.n
);
1510 gmx_fio_do_real(fio
,iparams
->vsiten
.a
);
1515 /* We got rid of some parameters in version 68 */
1516 if(bRead
&& file_version
<68)
1518 gmx_fio_do_real(fio
,rdum
);
1519 gmx_fio_do_real(fio
,rdum
);
1520 gmx_fio_do_real(fio
,rdum
);
1521 gmx_fio_do_real(fio
,rdum
);
1523 gmx_fio_do_real(fio
,iparams
->gb
.sar
);
1524 gmx_fio_do_real(fio
,iparams
->gb
.st
);
1525 gmx_fio_do_real(fio
,iparams
->gb
.pi
);
1526 gmx_fio_do_real(fio
,iparams
->gb
.gbr
);
1527 gmx_fio_do_real(fio
,iparams
->gb
.bmlt
);
1530 gmx_fio_do_int(fio
,iparams
->cmap
.cmapA
);
1531 gmx_fio_do_int(fio
,iparams
->cmap
.cmapB
);
1534 gmx_fatal(FARGS
,"unknown function type %d (%s) in %s line %d",
1535 ftype
,interaction_function
[ftype
].name
,__FILE__
,__LINE__
);
1538 gmx_fio_unset_comment(fio
);
1541 static void do_ilist(t_fileio
*fio
, t_ilist
*ilist
,gmx_bool bRead
,int file_version
,
1548 gmx_fio_set_comment(fio
, interaction_function
[ftype
].name
);
1550 if (file_version
< 44) {
1551 for(i
=0; i
<MAXNODES
; i
++)
1552 gmx_fio_do_int(fio
,idum
);
1554 gmx_fio_do_int(fio
,ilist
->nr
);
1556 snew(ilist
->iatoms
,ilist
->nr
);
1557 bDum
=gmx_fio_ndo_int(fio
,ilist
->iatoms
,ilist
->nr
);
1559 gmx_fio_unset_comment(fio
);
1562 static void do_ffparams(t_fileio
*fio
, gmx_ffparams_t
*ffparams
,
1563 gmx_bool bRead
, int file_version
)
1569 gmx_fio_do_int(fio
,ffparams
->atnr
);
1570 if (file_version
< 57) {
1571 gmx_fio_do_int(fio
,idum
);
1573 gmx_fio_do_int(fio
,ffparams
->ntypes
);
1575 fprintf(debug
,"ffparams->atnr = %d, ntypes = %d\n",
1576 ffparams
->atnr
,ffparams
->ntypes
);
1578 snew(ffparams
->functype
,ffparams
->ntypes
);
1579 snew(ffparams
->iparams
,ffparams
->ntypes
);
1581 /* Read/write all the function types */
1582 bDum
=gmx_fio_ndo_int(fio
,ffparams
->functype
,ffparams
->ntypes
);
1584 pr_ivec(debug
,0,"functype",ffparams
->functype
,ffparams
->ntypes
,TRUE
);
1586 if (file_version
>= 66) {
1587 gmx_fio_do_double(fio
,ffparams
->reppow
);
1589 ffparams
->reppow
= 12.0;
1592 if (file_version
>= 57) {
1593 gmx_fio_do_real(fio
,ffparams
->fudgeQQ
);
1596 /* Check whether all these function types are supported by the code.
1597 * In practice the code is backwards compatible, which means that the
1598 * numbering may have to be altered from old numbering to new numbering
1600 for (i
=0; (i
<ffparams
->ntypes
); i
++) {
1602 /* Loop over file versions */
1603 for (k
=0; (k
<NFTUPD
); k
++)
1604 /* Compare the read file_version to the update table */
1605 if ((file_version
< ftupd
[k
].fvnr
) &&
1606 (ffparams
->functype
[i
] >= ftupd
[k
].ftype
)) {
1607 ffparams
->functype
[i
] += 1;
1609 fprintf(debug
,"Incrementing function type %d to %d (due to %s)\n",
1610 i
,ffparams
->functype
[i
],
1611 interaction_function
[ftupd
[k
].ftype
].longname
);
1616 do_iparams(fio
, ffparams
->functype
[i
],&ffparams
->iparams
[i
],bRead
,
1619 pr_iparams(debug
,ffparams
->functype
[i
],&ffparams
->iparams
[i
]);
1623 static void add_settle_atoms(t_ilist
*ilist
)
1627 /* Settle used to only store the first atom: add the other two */
1628 srenew(ilist
->iatoms
,2*ilist
->nr
);
1629 for(i
=ilist
->nr
/2-1; i
>=0; i
--)
1631 ilist
->iatoms
[4*i
+0] = ilist
->iatoms
[2*i
+0];
1632 ilist
->iatoms
[4*i
+1] = ilist
->iatoms
[2*i
+1];
1633 ilist
->iatoms
[4*i
+2] = ilist
->iatoms
[2*i
+1] + 1;
1634 ilist
->iatoms
[4*i
+3] = ilist
->iatoms
[2*i
+1] + 2;
1636 ilist
->nr
= 2*ilist
->nr
;
1639 static void do_ilists(t_fileio
*fio
, t_ilist
*ilist
,gmx_bool bRead
,
1642 int i
,j
,renum
[F_NRE
];
1643 gmx_bool bDum
=TRUE
,bClear
;
1646 for(j
=0; (j
<F_NRE
); j
++) {
1649 for (k
=0; k
<NFTUPD
; k
++)
1650 if ((file_version
< ftupd
[k
].fvnr
) && (j
== ftupd
[k
].ftype
))
1654 ilist
[j
].iatoms
= NULL
;
1656 do_ilist(fio
, &ilist
[j
],bRead
,file_version
,j
);
1657 if (file_version
< 78 && j
== F_SETTLE
&& ilist
[j
].nr
> 0)
1659 add_settle_atoms(&ilist
[j
]);
1663 if (bRead && gmx_debug_at)
1664 pr_ilist(debug,0,interaction_function[j].longname,
1665 functype,&ilist[j],TRUE);
1670 static void do_idef(t_fileio
*fio
, gmx_ffparams_t
*ffparams
,gmx_moltype_t
*molt
,
1671 gmx_bool bRead
, int file_version
)
1673 do_ffparams(fio
, ffparams
,bRead
,file_version
);
1675 if (file_version
>= 54) {
1676 gmx_fio_do_real(fio
,ffparams
->fudgeQQ
);
1679 do_ilists(fio
, molt
->ilist
,bRead
,file_version
);
1682 static void do_block(t_fileio
*fio
, t_block
*block
,gmx_bool bRead
,int file_version
)
1684 int i
,idum
,dum_nra
,*dum_a
;
1687 if (file_version
< 44)
1688 for(i
=0; i
<MAXNODES
; i
++)
1689 gmx_fio_do_int(fio
,idum
);
1690 gmx_fio_do_int(fio
,block
->nr
);
1691 if (file_version
< 51)
1692 gmx_fio_do_int(fio
,dum_nra
);
1694 block
->nalloc_index
= block
->nr
+1;
1695 snew(block
->index
,block
->nalloc_index
);
1697 bDum
=gmx_fio_ndo_int(fio
,block
->index
,block
->nr
+1);
1699 if (file_version
< 51 && dum_nra
> 0) {
1700 snew(dum_a
,dum_nra
);
1701 bDum
=gmx_fio_ndo_int(fio
,dum_a
,dum_nra
);
1706 static void do_blocka(t_fileio
*fio
, t_blocka
*block
,gmx_bool bRead
,
1712 if (file_version
< 44)
1713 for(i
=0; i
<MAXNODES
; i
++)
1714 gmx_fio_do_int(fio
,idum
);
1715 gmx_fio_do_int(fio
,block
->nr
);
1716 gmx_fio_do_int(fio
,block
->nra
);
1718 block
->nalloc_index
= block
->nr
+1;
1719 snew(block
->index
,block
->nalloc_index
);
1720 block
->nalloc_a
= block
->nra
;
1721 snew(block
->a
,block
->nalloc_a
);
1723 bDum
=gmx_fio_ndo_int(fio
,block
->index
,block
->nr
+1);
1724 bDum
=gmx_fio_ndo_int(fio
,block
->a
,block
->nra
);
1727 static void do_atom(t_fileio
*fio
, t_atom
*atom
,int ngrp
,gmx_bool bRead
,
1728 int file_version
, gmx_groups_t
*groups
,int atnr
)
1732 gmx_fio_do_real(fio
,atom
->m
);
1733 gmx_fio_do_real(fio
,atom
->q
);
1734 gmx_fio_do_real(fio
,atom
->mB
);
1735 gmx_fio_do_real(fio
,atom
->qB
);
1736 gmx_fio_do_ushort(fio
, atom
->type
);
1737 gmx_fio_do_ushort(fio
, atom
->typeB
);
1738 gmx_fio_do_int(fio
,atom
->ptype
);
1739 gmx_fio_do_int(fio
,atom
->resind
);
1740 if (file_version
>= 52)
1741 gmx_fio_do_int(fio
,atom
->atomnumber
);
1743 atom
->atomnumber
= NOTSET
;
1744 if (file_version
< 23)
1746 else if (file_version
< 39)
1751 if (file_version
< 57) {
1752 unsigned char uchar
[egcNR
];
1753 gmx_fio_ndo_uchar(fio
,uchar
,myngrp
);
1754 for(i
=myngrp
; (i
<ngrp
); i
++) {
1757 /* Copy the old data format to the groups struct */
1758 for(i
=0; i
<ngrp
; i
++) {
1759 groups
->grpnr
[i
][atnr
] = uchar
[i
];
1764 static void do_grps(t_fileio
*fio
, int ngrp
,t_grps grps
[],gmx_bool bRead
,
1770 if (file_version
< 23)
1772 else if (file_version
< 39)
1777 for(j
=0; (j
<ngrp
); j
++) {
1779 gmx_fio_do_int(fio
,grps
[j
].nr
);
1781 snew(grps
[j
].nm_ind
,grps
[j
].nr
);
1782 bDum
=gmx_fio_ndo_int(fio
,grps
[j
].nm_ind
,grps
[j
].nr
);
1786 snew(grps
[j
].nm_ind
,grps
[j
].nr
);
1791 static void do_symstr(t_fileio
*fio
, char ***nm
,gmx_bool bRead
,t_symtab
*symtab
)
1796 gmx_fio_do_int(fio
,ls
);
1797 *nm
= get_symtab_handle(symtab
,ls
);
1800 ls
= lookup_symtab(symtab
,*nm
);
1801 gmx_fio_do_int(fio
,ls
);
1805 static void do_strstr(t_fileio
*fio
, int nstr
,char ***nm
,gmx_bool bRead
,
1810 for (j
=0; (j
<nstr
); j
++)
1811 do_symstr(fio
, &(nm
[j
]),bRead
,symtab
);
1814 static void do_resinfo(t_fileio
*fio
, int n
,t_resinfo
*ri
,gmx_bool bRead
,
1815 t_symtab
*symtab
, int file_version
)
1819 for (j
=0; (j
<n
); j
++) {
1820 do_symstr(fio
, &(ri
[j
].name
),bRead
,symtab
);
1821 if (file_version
>= 63) {
1822 gmx_fio_do_int(fio
,ri
[j
].nr
);
1823 gmx_fio_do_uchar(fio
, ri
[j
].ic
);
1831 static void do_atoms(t_fileio
*fio
, t_atoms
*atoms
,gmx_bool bRead
,t_symtab
*symtab
,
1833 gmx_groups_t
*groups
)
1837 gmx_fio_do_int(fio
,atoms
->nr
);
1838 gmx_fio_do_int(fio
,atoms
->nres
);
1839 if (file_version
< 57) {
1840 gmx_fio_do_int(fio
,groups
->ngrpname
);
1841 for(i
=0; i
<egcNR
; i
++) {
1842 groups
->ngrpnr
[i
] = atoms
->nr
;
1843 snew(groups
->grpnr
[i
],groups
->ngrpnr
[i
]);
1847 snew(atoms
->atom
,atoms
->nr
);
1848 snew(atoms
->atomname
,atoms
->nr
);
1849 snew(atoms
->atomtype
,atoms
->nr
);
1850 snew(atoms
->atomtypeB
,atoms
->nr
);
1851 snew(atoms
->resinfo
,atoms
->nres
);
1852 if (file_version
< 57) {
1853 snew(groups
->grpname
,groups
->ngrpname
);
1855 atoms
->pdbinfo
= NULL
;
1857 for(i
=0; (i
<atoms
->nr
); i
++) {
1858 do_atom(fio
, &atoms
->atom
[i
],egcNR
,bRead
, file_version
,groups
,i
);
1860 do_strstr(fio
, atoms
->nr
,atoms
->atomname
,bRead
,symtab
);
1861 if (bRead
&& (file_version
<= 20)) {
1862 for(i
=0; i
<atoms
->nr
; i
++) {
1863 atoms
->atomtype
[i
] = put_symtab(symtab
,"?");
1864 atoms
->atomtypeB
[i
] = put_symtab(symtab
,"?");
1867 do_strstr(fio
, atoms
->nr
,atoms
->atomtype
,bRead
,symtab
);
1868 do_strstr(fio
, atoms
->nr
,atoms
->atomtypeB
,bRead
,symtab
);
1870 do_resinfo(fio
, atoms
->nres
,atoms
->resinfo
,bRead
,symtab
,file_version
);
1872 if (file_version
< 57) {
1873 do_strstr(fio
, groups
->ngrpname
,groups
->grpname
,bRead
,symtab
);
1875 do_grps(fio
, egcNR
,groups
->grps
,bRead
,file_version
);
1879 static void do_groups(t_fileio
*fio
, gmx_groups_t
*groups
,
1880 gmx_bool bRead
,t_symtab
*symtab
,
1886 do_grps(fio
, egcNR
,groups
->grps
,bRead
,file_version
);
1887 gmx_fio_do_int(fio
,groups
->ngrpname
);
1889 snew(groups
->grpname
,groups
->ngrpname
);
1891 do_strstr(fio
, groups
->ngrpname
,groups
->grpname
,bRead
,symtab
);
1892 for(g
=0; g
<egcNR
; g
++) {
1893 gmx_fio_do_int(fio
,groups
->ngrpnr
[g
]);
1894 if (groups
->ngrpnr
[g
] == 0) {
1896 groups
->grpnr
[g
] = NULL
;
1900 snew(groups
->grpnr
[g
],groups
->ngrpnr
[g
]);
1902 bDum
=gmx_fio_ndo_uchar(fio
, groups
->grpnr
[g
],groups
->ngrpnr
[g
]);
1907 static void do_atomtypes(t_fileio
*fio
, t_atomtypes
*atomtypes
,gmx_bool bRead
,
1908 t_symtab
*symtab
,int file_version
)
1911 gmx_bool bDum
= TRUE
;
1913 if (file_version
> 25) {
1914 gmx_fio_do_int(fio
,atomtypes
->nr
);
1917 snew(atomtypes
->radius
,j
);
1918 snew(atomtypes
->vol
,j
);
1919 snew(atomtypes
->surftens
,j
);
1920 snew(atomtypes
->atomnumber
,j
);
1921 snew(atomtypes
->gb_radius
,j
);
1922 snew(atomtypes
->S_hct
,j
);
1924 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->radius
,j
);
1925 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->vol
,j
);
1926 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->surftens
,j
);
1927 if(file_version
>= 40)
1929 bDum
=gmx_fio_ndo_int(fio
,atomtypes
->atomnumber
,j
);
1931 if(file_version
>= 60)
1933 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->gb_radius
,j
);
1934 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->S_hct
,j
);
1937 /* File versions prior to 26 cannot do GBSA,
1938 * so they dont use this structure
1941 atomtypes
->radius
= NULL
;
1942 atomtypes
->vol
= NULL
;
1943 atomtypes
->surftens
= NULL
;
1944 atomtypes
->atomnumber
= NULL
;
1945 atomtypes
->gb_radius
= NULL
;
1946 atomtypes
->S_hct
= NULL
;
1950 static void do_symtab(t_fileio
*fio
, t_symtab
*symtab
,gmx_bool bRead
)
1956 gmx_fio_do_int(fio
,symtab
->nr
);
1959 snew(symtab
->symbuf
,1);
1960 symbuf
= symtab
->symbuf
;
1961 symbuf
->bufsize
= nr
;
1962 snew(symbuf
->buf
,nr
);
1963 for (i
=0; (i
<nr
); i
++) {
1964 gmx_fio_do_string(fio
,buf
);
1965 symbuf
->buf
[i
]=strdup(buf
);
1969 symbuf
= symtab
->symbuf
;
1970 while (symbuf
!=NULL
) {
1971 for (i
=0; (i
<symbuf
->bufsize
) && (i
<nr
); i
++)
1972 gmx_fio_do_string(fio
,symbuf
->buf
[i
]);
1974 symbuf
=symbuf
->next
;
1977 gmx_fatal(FARGS
,"nr of symtab strings left: %d",nr
);
1981 static void do_cmap(t_fileio
*fio
, gmx_cmap_t
*cmap_grid
, gmx_bool bRead
)
1983 int i
,j
,ngrid
,gs
,nelem
;
1985 gmx_fio_do_int(fio
,cmap_grid
->ngrid
);
1986 gmx_fio_do_int(fio
,cmap_grid
->grid_spacing
);
1988 ngrid
= cmap_grid
->ngrid
;
1989 gs
= cmap_grid
->grid_spacing
;
1994 snew(cmap_grid
->cmapdata
,ngrid
);
1996 for(i
=0;i
<cmap_grid
->ngrid
;i
++)
1998 snew(cmap_grid
->cmapdata
[i
].cmap
,4*nelem
);
2002 for(i
=0;i
<cmap_grid
->ngrid
;i
++)
2004 for(j
=0;j
<nelem
;j
++)
2006 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4]);
2007 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+1]);
2008 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+2]);
2009 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+3]);
2015 void tpx_make_chain_identifiers(t_atoms
*atoms
,t_block
*mols
)
2021 /* We always assign a new chain number, but save the chain id characters
2022 * for larger molecules.
2024 #define CHAIN_MIN_ATOMS 15
2028 for(m
=0; m
<mols
->nr
; m
++)
2031 a1
=mols
->index
[m
+1];
2032 if ((a1
-a0
>= CHAIN_MIN_ATOMS
) && (chainid
<= 'Z'))
2041 for(a
=a0
; a
<a1
; a
++)
2043 atoms
->resinfo
[atoms
->atom
[a
].resind
].chainnum
= chainnum
;
2044 atoms
->resinfo
[atoms
->atom
[a
].resind
].chainid
= c
;
2049 /* Blank out the chain id if there was only one chain */
2052 for(r
=0; r
<atoms
->nres
; r
++)
2054 atoms
->resinfo
[r
].chainid
= ' ';
2059 static void do_moltype(t_fileio
*fio
, gmx_moltype_t
*molt
,gmx_bool bRead
,
2060 t_symtab
*symtab
, int file_version
,
2061 gmx_groups_t
*groups
)
2065 if (file_version
>= 57) {
2066 do_symstr(fio
, &(molt
->name
),bRead
,symtab
);
2069 do_atoms(fio
, &molt
->atoms
, bRead
, symtab
, file_version
, groups
);
2071 if (bRead
&& gmx_debug_at
) {
2072 pr_atoms(debug
,0,"atoms",&molt
->atoms
,TRUE
);
2075 if (file_version
>= 57) {
2076 do_ilists(fio
, molt
->ilist
,bRead
,file_version
);
2078 do_block(fio
, &molt
->cgs
,bRead
,file_version
);
2079 if (bRead
&& gmx_debug_at
) {
2080 pr_block(debug
,0,"cgs",&molt
->cgs
,TRUE
);
2084 /* This used to be in the atoms struct */
2085 do_blocka(fio
, &molt
->excls
, bRead
, file_version
);
2088 static void do_molblock(t_fileio
*fio
, gmx_molblock_t
*molb
,gmx_bool bRead
,
2093 gmx_fio_do_int(fio
,molb
->type
);
2094 gmx_fio_do_int(fio
,molb
->nmol
);
2095 gmx_fio_do_int(fio
,molb
->natoms_mol
);
2096 /* Position restraint coordinates */
2097 gmx_fio_do_int(fio
,molb
->nposres_xA
);
2098 if (molb
->nposres_xA
> 0) {
2100 snew(molb
->posres_xA
,molb
->nposres_xA
);
2102 gmx_fio_ndo_rvec(fio
,molb
->posres_xA
,molb
->nposres_xA
);
2104 gmx_fio_do_int(fio
,molb
->nposres_xB
);
2105 if (molb
->nposres_xB
> 0) {
2107 snew(molb
->posres_xB
,molb
->nposres_xB
);
2109 gmx_fio_ndo_rvec(fio
,molb
->posres_xB
,molb
->nposres_xB
);
2114 static t_block
mtop_mols(gmx_mtop_t
*mtop
)
2120 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
2121 mols
.nr
+= mtop
->molblock
[mb
].nmol
;
2123 mols
.nalloc_index
= mols
.nr
+ 1;
2124 snew(mols
.index
,mols
.nalloc_index
);
2129 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
2130 for(mol
=0; mol
<mtop
->molblock
[mb
].nmol
; mol
++) {
2131 a
+= mtop
->molblock
[mb
].natoms_mol
;
2140 static void add_posres_molblock(gmx_mtop_t
*mtop
)
2145 gmx_molblock_t
*molb
;
2148 il
= &mtop
->moltype
[0].ilist
[F_POSRES
];
2154 for(i
=0; i
<il
->nr
; i
+=2) {
2155 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
2156 am
= max(am
,il
->iatoms
[i
+1]);
2157 if (ip
->posres
.pos0B
[XX
] != ip
->posres
.pos0A
[XX
] ||
2158 ip
->posres
.pos0B
[YY
] != ip
->posres
.pos0A
[YY
] ||
2159 ip
->posres
.pos0B
[ZZ
] != ip
->posres
.pos0A
[ZZ
]) {
2163 /* Make the posres coordinate block end at a molecule end */
2165 while(am
>= mtop
->mols
.index
[mol
+1]) {
2168 molb
= &mtop
->molblock
[0];
2169 molb
->nposres_xA
= mtop
->mols
.index
[mol
+1];
2170 snew(molb
->posres_xA
,molb
->nposres_xA
);
2172 molb
->nposres_xB
= molb
->nposres_xA
;
2173 snew(molb
->posres_xB
,molb
->nposres_xB
);
2175 molb
->nposres_xB
= 0;
2177 for(i
=0; i
<il
->nr
; i
+=2) {
2178 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
2179 a
= il
->iatoms
[i
+1];
2180 molb
->posres_xA
[a
][XX
] = ip
->posres
.pos0A
[XX
];
2181 molb
->posres_xA
[a
][YY
] = ip
->posres
.pos0A
[YY
];
2182 molb
->posres_xA
[a
][ZZ
] = ip
->posres
.pos0A
[ZZ
];
2184 molb
->posres_xB
[a
][XX
] = ip
->posres
.pos0B
[XX
];
2185 molb
->posres_xB
[a
][YY
] = ip
->posres
.pos0B
[YY
];
2186 molb
->posres_xB
[a
][ZZ
] = ip
->posres
.pos0B
[ZZ
];
2191 static void set_disres_npair(gmx_mtop_t
*mtop
)
2198 ip
= mtop
->ffparams
.iparams
;
2200 for(mt
=0; mt
<mtop
->nmoltype
; mt
++) {
2201 il
= &mtop
->moltype
[mt
].ilist
[F_DISRES
];
2205 for(i
=0; i
<il
->nr
; i
+=3) {
2207 if (i
+3 == il
->nr
|| ip
[a
[i
]].disres
.label
!= ip
[a
[i
+3]].disres
.label
) {
2208 ip
[a
[i
]].disres
.npair
= npair
;
2216 static void do_mtop(t_fileio
*fio
, gmx_mtop_t
*mtop
,gmx_bool bRead
,
2224 do_symtab(fio
, &(mtop
->symtab
),bRead
);
2226 pr_symtab(debug
,0,"symtab",&mtop
->symtab
);
2228 do_symstr(fio
, &(mtop
->name
),bRead
,&(mtop
->symtab
));
2230 if (file_version
>= 57) {
2231 do_ffparams(fio
, &mtop
->ffparams
,bRead
,file_version
);
2233 gmx_fio_do_int(fio
,mtop
->nmoltype
);
2238 snew(mtop
->moltype
,mtop
->nmoltype
);
2239 if (file_version
< 57) {
2240 mtop
->moltype
[0].name
= mtop
->name
;
2243 for(mt
=0; mt
<mtop
->nmoltype
; mt
++) {
2244 do_moltype(fio
, &mtop
->moltype
[mt
],bRead
,&mtop
->symtab
,file_version
,
2248 if (file_version
>= 57) {
2249 gmx_fio_do_int(fio
,mtop
->nmolblock
);
2251 mtop
->nmolblock
= 1;
2254 snew(mtop
->molblock
,mtop
->nmolblock
);
2256 if (file_version
>= 57) {
2257 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
2258 do_molblock(fio
, &mtop
->molblock
[mb
],bRead
,file_version
);
2260 gmx_fio_do_int(fio
,mtop
->natoms
);
2262 mtop
->molblock
[0].type
= 0;
2263 mtop
->molblock
[0].nmol
= 1;
2264 mtop
->molblock
[0].natoms_mol
= mtop
->moltype
[0].atoms
.nr
;
2265 mtop
->molblock
[0].nposres_xA
= 0;
2266 mtop
->molblock
[0].nposres_xB
= 0;
2269 do_atomtypes (fio
, &(mtop
->atomtypes
),bRead
,&(mtop
->symtab
), file_version
);
2271 pr_atomtypes(debug
,0,"atomtypes",&mtop
->atomtypes
,TRUE
);
2273 if (file_version
< 57) {
2274 /* Debug statements are inside do_idef */
2275 do_idef (fio
, &mtop
->ffparams
,&mtop
->moltype
[0],bRead
,file_version
);
2276 mtop
->natoms
= mtop
->moltype
[0].atoms
.nr
;
2279 if(file_version
>= 65)
2281 do_cmap(fio
, &mtop
->ffparams
.cmap_grid
,bRead
);
2285 mtop
->ffparams
.cmap_grid
.ngrid
= 0;
2286 mtop
->ffparams
.cmap_grid
.grid_spacing
= 0;
2287 mtop
->ffparams
.cmap_grid
.cmapdata
= NULL
;
2290 if (file_version
>= 57) {
2291 do_groups(fio
, &mtop
->groups
,bRead
,&(mtop
->symtab
),file_version
);
2294 if (file_version
< 57) {
2295 do_block(fio
, &mtop
->moltype
[0].cgs
,bRead
,file_version
);
2296 if (bRead
&& gmx_debug_at
) {
2297 pr_block(debug
,0,"cgs",&mtop
->moltype
[0].cgs
,TRUE
);
2299 do_block(fio
, &mtop
->mols
,bRead
,file_version
);
2300 /* Add the posres coordinates to the molblock */
2301 add_posres_molblock(mtop
);
2304 if (file_version
>= 57) {
2305 mtop
->mols
= mtop_mols(mtop
);
2308 pr_block(debug
,0,"mols",&mtop
->mols
,TRUE
);
2312 if (file_version
< 51) {
2313 /* Here used to be the shake blocks */
2314 do_blocka(fio
, &dumb
,bRead
,file_version
);
2322 close_symtab(&(mtop
->symtab
));
2326 /* If TopOnlyOK is TRUE then we can read even future versions
2327 * of tpx files, provided the file_generation hasn't changed.
2328 * If it is FALSE, we need the inputrecord too, and bail out
2329 * if the file is newer than the program.
2331 * The version and generation if the topology (see top of this file)
2332 * are returned in the two last arguments.
2334 * If possible, we will read the inputrec even when TopOnlyOK is TRUE.
2336 static void do_tpxheader(t_fileio
*fio
,gmx_bool bRead
,t_tpxheader
*tpx
,
2337 gmx_bool TopOnlyOK
, int *file_version
,
2338 int *file_generation
)
2341 char file_tag
[STRLEN
];
2348 gmx_fio_checktype(fio
);
2349 gmx_fio_setdebug(fio
,bDebugMode());
2351 /* NEW! XDR tpb file */
2352 precision
= sizeof(real
);
2354 gmx_fio_do_string(fio
,buf
);
2355 if (strncmp(buf
,"VERSION",7))
2356 gmx_fatal(FARGS
,"Can not read file %s,\n"
2357 " this file is from a Gromacs version which is older than 2.0\n"
2358 " Make a new one with grompp or use a gro or pdb file, if possible",
2359 gmx_fio_getname(fio
));
2360 gmx_fio_do_int(fio
,precision
);
2361 bDouble
= (precision
== sizeof(double));
2362 if ((precision
!= sizeof(float)) && !bDouble
)
2363 gmx_fatal(FARGS
,"Unknown precision in file %s: real is %d bytes "
2364 "instead of %d or %d",
2365 gmx_fio_getname(fio
),precision
,sizeof(float),sizeof(double));
2366 gmx_fio_setprecision(fio
,bDouble
);
2367 fprintf(stderr
,"Reading file %s, %s (%s precision)\n",
2368 gmx_fio_getname(fio
),buf
,bDouble
? "double" : "single");
2371 gmx_fio_write_string(fio
,GromacsVersion());
2372 bDouble
= (precision
== sizeof(double));
2373 gmx_fio_setprecision(fio
,bDouble
);
2374 gmx_fio_do_int(fio
,precision
);
2376 sprintf(file_tag
,"%s",tpx_tag
);
2377 fgen
= tpx_generation
;
2380 /* Check versions! */
2381 gmx_fio_do_int(fio
,fver
);
2383 /* This is for backward compatibility with development versions 77-79
2384 * where the tag was, mistakenly, placed before the generation,
2385 * which would cause a segv instead of a proper error message
2386 * when reading the topology only from tpx with <77 code.
2388 if (fver
>= 77 && fver
<= 79)
2390 gmx_fio_do_string(fio
,file_tag
);
2395 gmx_fio_do_int(fio
,fgen
);
2404 gmx_fio_do_string(fio
,file_tag
);
2410 /* Versions before 77 don't have the tag, set it to release */
2411 sprintf(file_tag
,"%s",TPX_TAG_RELEASE
);
2414 if (strcmp(file_tag
,tpx_tag
) != 0)
2416 fprintf(stderr
,"Note: file tpx tag '%s', software tpx tag '%s'\n",
2419 /* We only support reading tpx files with the same tag as the code
2420 * or tpx files with the release tag and with lower version number.
2422 if (!strcmp(file_tag
,TPX_TAG_RELEASE
) == 0 && fver
< tpx_version
)
2424 gmx_fatal(FARGS
,"tpx tag/version mismatch: reading tpx file (%s) version %d, tag '%s' with program for tpx version %d, tag '%s'",
2425 gmx_fio_getname(fio
),fver
,file_tag
,
2426 tpx_version
,tpx_tag
);
2431 if (file_version
!= NULL
)
2433 *file_version
= fver
;
2435 if (file_generation
!= NULL
)
2437 *file_generation
= fgen
;
2441 if ((fver
<= tpx_incompatible_version
) ||
2442 ((fver
> tpx_version
) && !TopOnlyOK
) ||
2443 (fgen
> tpx_generation
))
2444 gmx_fatal(FARGS
,"reading tpx file (%s) version %d with version %d program",
2445 gmx_fio_getname(fio
),fver
,tpx_version
);
2447 do_section(fio
,eitemHEADER
,bRead
);
2448 gmx_fio_do_int(fio
,tpx
->natoms
);
2450 gmx_fio_do_int(fio
,tpx
->ngtc
);
2454 gmx_fio_do_int(fio
,idum
);
2455 gmx_fio_do_real(fio
,rdum
);
2457 /*a better decision will eventually (5.0 or later) need to be made
2458 on how to treat the alchemical state of the system, which can now
2459 vary through a simulation, and cannot be completely described
2460 though a single lambda variable, or even a single state
2461 index. Eventually, should probably be a vector. MRS*/
2464 gmx_fio_do_int(fio
,tpx
->fep_state
);
2466 gmx_fio_do_real(fio
,tpx
->lambda
);
2467 gmx_fio_do_int(fio
,tpx
->bIr
);
2468 gmx_fio_do_int(fio
,tpx
->bTop
);
2469 gmx_fio_do_int(fio
,tpx
->bX
);
2470 gmx_fio_do_int(fio
,tpx
->bV
);
2471 gmx_fio_do_int(fio
,tpx
->bF
);
2472 gmx_fio_do_int(fio
,tpx
->bBox
);
2474 if((fgen
> tpx_generation
)) {
2475 /* This can only happen if TopOnlyOK=TRUE */
2480 static int do_tpx(t_fileio
*fio
, gmx_bool bRead
,
2481 t_inputrec
*ir
,t_state
*state
,rvec
*f
,gmx_mtop_t
*mtop
,
2482 gmx_bool bXVallocated
)
2487 gmx_bool TopOnlyOK
,bDum
=TRUE
;
2488 int file_version
,file_generation
;
2492 gmx_bool bPeriodicMols
;
2495 tpx
.natoms
= state
->natoms
;
2496 tpx
.ngtc
= state
->ngtc
; /* need to add nnhpres here? */
2497 tpx
.fep_state
= state
->fep_state
;
2498 tpx
.lambda
= state
->lambda
[efptFEP
];
2499 tpx
.bIr
= (ir
!= NULL
);
2500 tpx
.bTop
= (mtop
!= NULL
);
2501 tpx
.bX
= (state
->x
!= NULL
);
2502 tpx
.bV
= (state
->v
!= NULL
);
2503 tpx
.bF
= (f
!= NULL
);
2507 TopOnlyOK
= (ir
==NULL
);
2509 do_tpxheader(fio
,bRead
,&tpx
,TopOnlyOK
,&file_version
,&file_generation
);
2513 /* state->lambda = tpx.lambda;*/ /*remove this eventually? */
2514 /* The init_state calls initialize the Nose-Hoover xi integrals to zero */
2518 init_state(state
,0,tpx
.ngtc
,0,0,0); /* nose-hoover chains */ /* eventually, need to add nnhpres here? */
2519 state
->natoms
= tpx
.natoms
;
2520 state
->nalloc
= tpx
.natoms
;
2524 init_state(state
,tpx
.natoms
,tpx
.ngtc
,0,0,0); /* nose-hoover chains */
2528 #define do_test(fio,b,p) if (bRead && (p!=NULL) && !b) gmx_fatal(FARGS,"No %s in %s",#p,gmx_fio_getname(fio))
2530 do_test(fio
,tpx
.bBox
,state
->box
);
2531 do_section(fio
,eitemBOX
,bRead
);
2533 gmx_fio_ndo_rvec(fio
,state
->box
,DIM
);
2534 if (file_version
>= 51) {
2535 gmx_fio_ndo_rvec(fio
,state
->box_rel
,DIM
);
2537 /* We initialize box_rel after reading the inputrec */
2538 clear_mat(state
->box_rel
);
2540 if (file_version
>= 28) {
2541 gmx_fio_ndo_rvec(fio
,state
->boxv
,DIM
);
2542 if (file_version
< 56) {
2544 gmx_fio_ndo_rvec(fio
,mdum
,DIM
);
2549 if (state
->ngtc
> 0 && file_version
>= 28) {
2551 /*ndo_double(state->nosehoover_xi,state->ngtc,bDum);*/
2552 /*ndo_double(state->nosehoover_vxi,state->ngtc,bDum);*/
2553 /*ndo_double(state->therm_integral,state->ngtc,bDum);*/
2554 snew(dumv
,state
->ngtc
);
2555 if (file_version
< 69) {
2556 bDum
=gmx_fio_ndo_real(fio
,dumv
,state
->ngtc
);
2558 /* These used to be the Berendsen tcoupl_lambda's */
2559 bDum
=gmx_fio_ndo_real(fio
,dumv
,state
->ngtc
);
2563 /* Prior to tpx version 26, the inputrec was here.
2564 * I moved it to enable partial forward-compatibility
2565 * for analysis/viewer programs.
2567 if(file_version
<26) {
2568 do_test(fio
,tpx
.bIr
,ir
);
2569 do_section(fio
,eitemIR
,bRead
);
2572 do_inputrec(fio
, ir
,bRead
,file_version
,
2573 mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
2575 pr_inputrec(debug
,0,"inputrec",ir
,FALSE
);
2578 do_inputrec(fio
, &dum_ir
,bRead
,file_version
,
2579 mtop
? &mtop
->ffparams
.fudgeQQ
:NULL
);
2581 pr_inputrec(debug
,0,"inputrec",&dum_ir
,FALSE
);
2582 done_inputrec(&dum_ir
);
2588 do_test(fio
,tpx
.bTop
,mtop
);
2589 do_section(fio
,eitemTOP
,bRead
);
2592 do_mtop(fio
,mtop
,bRead
, file_version
);
2594 do_mtop(fio
,&dum_top
,bRead
,file_version
);
2595 done_mtop(&dum_top
,TRUE
);
2598 do_test(fio
,tpx
.bX
,state
->x
);
2599 do_section(fio
,eitemX
,bRead
);
2602 state
->flags
|= (1<<estX
);
2604 gmx_fio_ndo_rvec(fio
,state
->x
,state
->natoms
);
2607 do_test(fio
,tpx
.bV
,state
->v
);
2608 do_section(fio
,eitemV
,bRead
);
2611 state
->flags
|= (1<<estV
);
2613 gmx_fio_ndo_rvec(fio
,state
->v
,state
->natoms
);
2616 do_test(fio
,tpx
.bF
,f
);
2617 do_section(fio
,eitemF
,bRead
);
2618 if (tpx
.bF
) gmx_fio_ndo_rvec(fio
,f
,state
->natoms
);
2620 /* Starting with tpx version 26, we have the inputrec
2621 * at the end of the file, so we can ignore it
2622 * if the file is never than the software (but still the
2623 * same generation - see comments at the top of this file.
2628 bPeriodicMols
= FALSE
;
2629 if (file_version
>= 26) {
2630 do_test(fio
,tpx
.bIr
,ir
);
2631 do_section(fio
,eitemIR
,bRead
);
2633 if (file_version
>= 53) {
2634 /* Removed the pbc info from do_inputrec, since we always want it */
2637 bPeriodicMols
= ir
->bPeriodicMols
;
2639 gmx_fio_do_int(fio
,ePBC
);
2640 gmx_fio_do_gmx_bool(fio
,bPeriodicMols
);
2642 if (file_generation
<= tpx_generation
&& ir
) {
2643 do_inputrec(fio
, ir
,bRead
,file_version
,mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
2645 pr_inputrec(debug
,0,"inputrec",ir
,FALSE
);
2646 if (file_version
< 51)
2647 set_box_rel(ir
,state
);
2648 if (file_version
< 53) {
2650 bPeriodicMols
= ir
->bPeriodicMols
;
2653 if (bRead
&& ir
&& file_version
>= 53) {
2654 /* We need to do this after do_inputrec, since that initializes ir */
2656 ir
->bPeriodicMols
= bPeriodicMols
;
2665 if (state
->ngtc
== 0)
2667 /* Reading old version without tcoupl state data: set it */
2668 init_gtc_state(state
,ir
->opts
.ngtc
,0,ir
->opts
.nhchainlength
);
2670 if (tpx
.bTop
&& mtop
)
2672 if (file_version
< 57)
2674 if (mtop
->moltype
[0].ilist
[F_DISRES
].nr
> 0)
2676 ir
->eDisre
= edrSimple
;
2680 ir
->eDisre
= edrNone
;
2683 set_disres_npair(mtop
);
2687 if (tpx
.bTop
&& mtop
)
2689 gmx_mtop_finalize(mtop
);
2692 if (file_version
>= 57)
2696 env
= getenv("GMX_NOCHARGEGROUPS");
2699 sscanf(env
,"%d",&ienv
);
2700 fprintf(stderr
,"\nFound env.var. GMX_NOCHARGEGROUPS = %d\n",
2705 "Will make single atomic charge groups in non-solvent%s\n",
2706 ienv
> 1 ? " and solvent" : "");
2707 gmx_mtop_make_atomic_charge_groups(mtop
,ienv
==1);
2709 fprintf(stderr
,"\n");
2717 /************************************************************
2719 * The following routines are the exported ones
2721 ************************************************************/
2723 t_fileio
*open_tpx(const char *fn
,const char *mode
)
2725 return gmx_fio_open(fn
,mode
);
2728 void close_tpx(t_fileio
*fio
)
2733 void read_tpxheader(const char *fn
, t_tpxheader
*tpx
, gmx_bool TopOnlyOK
,
2734 int *file_version
, int *file_generation
)
2738 fio
= open_tpx(fn
,"r");
2739 do_tpxheader(fio
,TRUE
,tpx
,TopOnlyOK
,file_version
,file_generation
);
2743 void write_tpx_state(const char *fn
,
2744 t_inputrec
*ir
,t_state
*state
,gmx_mtop_t
*mtop
)
2748 fio
= open_tpx(fn
,"w");
2749 do_tpx(fio
,FALSE
,ir
,state
,NULL
,mtop
,FALSE
);
2753 void read_tpx_state(const char *fn
,
2754 t_inputrec
*ir
,t_state
*state
,rvec
*f
,gmx_mtop_t
*mtop
)
2758 fio
= open_tpx(fn
,"r");
2759 do_tpx(fio
,TRUE
,ir
,state
,f
,mtop
,FALSE
);
2763 int read_tpx(const char *fn
,
2764 t_inputrec
*ir
, matrix box
,int *natoms
,
2765 rvec
*x
,rvec
*v
,rvec
*f
,gmx_mtop_t
*mtop
)
2773 fio
= open_tpx(fn
,"r");
2774 ePBC
= do_tpx(fio
,TRUE
,ir
,&state
,f
,mtop
,TRUE
);
2776 *natoms
= state
.natoms
;
2778 copy_mat(state
.box
,box
);
2786 int read_tpx_top(const char *fn
,
2787 t_inputrec
*ir
, matrix box
,int *natoms
,
2788 rvec
*x
,rvec
*v
,rvec
*f
,t_topology
*top
)
2794 ePBC
= read_tpx(fn
,ir
,box
,natoms
,x
,v
,f
,&mtop
);
2796 *top
= gmx_mtop_t_to_t_topology(&mtop
);
2801 gmx_bool
fn2bTPX(const char *file
)
2803 switch (fn2ftp(file
)) {
2813 gmx_bool
read_tps_conf(const char *infile
,char *title
,t_topology
*top
,int *ePBC
,
2814 rvec
**x
,rvec
**v
,matrix box
,gmx_bool bMass
)
2817 int natoms
,i
,version
,generation
;
2818 gmx_bool bTop
,bXNULL
=FALSE
;
2820 t_topology
*topconv
;
2823 bTop
= fn2bTPX(infile
);
2826 read_tpxheader(infile
,&header
,TRUE
,&version
,&generation
);
2828 snew(*x
,header
.natoms
);
2830 snew(*v
,header
.natoms
);
2832 *ePBC
= read_tpx(infile
,NULL
,box
,&natoms
,
2833 (x
==NULL
) ? NULL
: *x
,(v
==NULL
) ? NULL
: *v
,NULL
,mtop
);
2834 *top
= gmx_mtop_t_to_t_topology(mtop
);
2836 strcpy(title
,*top
->name
);
2837 tpx_make_chain_identifiers(&top
->atoms
,&top
->mols
);
2840 get_stx_coordnum(infile
,&natoms
);
2841 init_t_atoms(&top
->atoms
,natoms
,(fn2ftp(infile
) == efPDB
));
2850 read_stx_conf(infile
,title
,&top
->atoms
,*x
,(v
==NULL
) ? NULL
: *v
,ePBC
,box
);
2857 aps
= gmx_atomprop_init();
2858 for(i
=0; (i
<natoms
); i
++)
2859 if (!gmx_atomprop_query(aps
,epropMass
,
2860 *top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].name
,
2861 *top
->atoms
.atomname
[i
],
2862 &(top
->atoms
.atom
[i
].m
))) {
2864 fprintf(debug
,"Can not find mass for atom %s %d %s, setting to 1\n",
2865 *top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].name
,
2866 top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].nr
,
2867 *top
->atoms
.atomname
[i
]);
2869 gmx_atomprop_destroy(aps
);
2871 top
->idef
.ntypes
=-1;