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 /* This number should be increased whenever the file format changes! */
67 static const int tpx_version
= 73;
69 /* This number should only be increased when you edit the TOPOLOGY section
70 * of the tpx format. This way we can maintain forward compatibility too
71 * for all analysis tools and/or external programs that only need to
72 * know the atom/residue names, charges, and bond connectivity.
74 * It first appeared in tpx version 26, when I also moved the inputrecord
75 * to the end of the tpx file, so we can just skip it if we only
78 static const int tpx_generation
= 23;
80 /* This number should be the most recent backwards incompatible version
81 * I.e., if this number is 9, we cannot read tpx version 9 with this code.
83 static const int tpx_incompatible_version
= 9;
87 /* Struct used to maintain tpx compatibility when function types are added */
89 int fvnr
; /* file version number in which the function type first appeared */
90 int ftype
; /* function type */
94 *The entries should be ordered in:
95 * 1. ascending file version number
96 * 2. ascending function type number
98 /*static const t_ftupd ftupd[] = {
103 { 22, F_DISRESVIOL },
109 { 26, F_DIHRESVIOL },
110 { 30, F_CROSS_BOND_BONDS },
111 { 30, F_CROSS_BOND_ANGLES },
112 { 30, F_UREY_BRADLEY },
113 { 30, F_POLARIZATION }
117 *The entries should be ordered in:
118 * 1. ascending function type number
119 * 2. ascending file version number
121 static const t_ftupd ftupd
[] = {
122 { 20, F_CUBICBONDS
},
127 { 43, F_TABBONDSNC
},
128 { 70, F_RESTRBONDS
},
129 { 30, F_CROSS_BOND_BONDS
},
130 { 30, F_CROSS_BOND_ANGLES
},
131 { 30, F_UREY_BRADLEY
},
132 { 34, F_QUARTIC_ANGLES
},
142 { 72, F_NPSOLVATION
},
144 { 41, F_LJC_PAIRS_NB
},
147 { 32, F_COUL_RECIP
},
149 { 30, F_POLARIZATION
},
151 { 22, F_DISRESVIOL
},
155 { 26, F_DIHRESVIOL
},
160 { 46, F_ECONSERVED
},
165 #define NFTUPD asize(ftupd)
167 /* Needed for backward compatibility */
170 static void _do_section(t_fileio
*fio
,int key
,gmx_bool bRead
,const char *src
,
176 if (gmx_fio_getftp(fio
) == efTPA
) {
178 gmx_fio_write_string(fio
,itemstr
[key
]);
179 bDbg
= gmx_fio_getdebug(fio
);
180 gmx_fio_setdebug(fio
,FALSE
);
181 gmx_fio_write_string(fio
,comment_str
[key
]);
182 gmx_fio_setdebug(fio
,bDbg
);
185 if (gmx_fio_getdebug(fio
))
186 fprintf(stderr
,"Looking for section %s (%s, %d)",
187 itemstr
[key
],src
,line
);
190 gmx_fio_do_string(fio
,buf
);
191 } while ((gmx_strcasecmp(buf
,itemstr
[key
]) != 0));
193 if (gmx_strcasecmp(buf
,itemstr
[key
]) != 0)
194 gmx_fatal(FARGS
,"\nCould not find section heading %s",itemstr
[key
]);
195 else if (gmx_fio_getdebug(fio
))
196 fprintf(stderr
," and found it\n");
201 #define do_section(fio,key,bRead) _do_section(fio,key,bRead,__FILE__,__LINE__)
203 /**************************************************************
205 * Now the higer level routines that do io of the structures and arrays
207 **************************************************************/
208 static void do_pullgrp(t_fileio
*fio
, t_pullgrp
*pgrp
, gmx_bool bRead
,
214 gmx_fio_do_int(fio
,pgrp
->nat
);
216 snew(pgrp
->ind
,pgrp
->nat
);
217 bDum
=gmx_fio_ndo_int(fio
,pgrp
->ind
,pgrp
->nat
);
218 gmx_fio_do_int(fio
,pgrp
->nweight
);
220 snew(pgrp
->weight
,pgrp
->nweight
);
221 bDum
=gmx_fio_ndo_real(fio
,pgrp
->weight
,pgrp
->nweight
);
222 gmx_fio_do_int(fio
,pgrp
->pbcatom
);
223 gmx_fio_do_rvec(fio
,pgrp
->vec
);
224 gmx_fio_do_rvec(fio
,pgrp
->init
);
225 gmx_fio_do_real(fio
,pgrp
->rate
);
226 gmx_fio_do_real(fio
,pgrp
->k
);
227 if (file_version
>= 56) {
228 gmx_fio_do_real(fio
,pgrp
->kB
);
234 static void do_pull(t_fileio
*fio
, t_pull
*pull
,gmx_bool bRead
, int file_version
)
238 gmx_fio_do_int(fio
,pull
->ngrp
);
239 gmx_fio_do_int(fio
,pull
->eGeom
);
240 gmx_fio_do_ivec(fio
,pull
->dim
);
241 gmx_fio_do_real(fio
,pull
->cyl_r1
);
242 gmx_fio_do_real(fio
,pull
->cyl_r0
);
243 gmx_fio_do_real(fio
,pull
->constr_tol
);
244 gmx_fio_do_int(fio
,pull
->nstxout
);
245 gmx_fio_do_int(fio
,pull
->nstfout
);
247 snew(pull
->grp
,pull
->ngrp
+1);
248 for(g
=0; g
<pull
->ngrp
+1; g
++)
249 do_pullgrp(fio
,&pull
->grp
[g
],bRead
,file_version
);
252 static void do_inputrec(t_fileio
*fio
, t_inputrec
*ir
,gmx_bool bRead
,
253 int file_version
, real
*fudgeQQ
)
255 int i
,j
,k
,*tmp
,idum
=0;
260 real zerotemptime
,finish_t
,init_temp
,finish_temp
;
262 if (file_version
!= tpx_version
)
264 /* Give a warning about features that are not accessible */
265 fprintf(stderr
,"Note: tpx file_version %d, software version %d\n",
266 file_version
,tpx_version
);
274 if (file_version
== 0)
279 /* Basic inputrec stuff */
280 gmx_fio_do_int(fio
,ir
->eI
);
281 if (file_version
>= 62) {
282 gmx_fio_do_gmx_large_int(fio
, ir
->nsteps
);
284 gmx_fio_do_int(fio
,idum
);
287 if(file_version
> 25) {
288 if (file_version
>= 62) {
289 gmx_fio_do_gmx_large_int(fio
, ir
->init_step
);
291 gmx_fio_do_int(fio
,idum
);
292 ir
->init_step
= idum
;
298 if(file_version
>= 58)
299 gmx_fio_do_int(fio
,ir
->simulation_part
);
301 ir
->simulation_part
=1;
303 if (file_version
>= 67) {
304 gmx_fio_do_int(fio
,ir
->nstcalcenergy
);
306 ir
->nstcalcenergy
= 1;
308 if (file_version
< 53) {
309 /* The pbc info has been moved out of do_inputrec,
310 * since we always want it, also without reading the inputrec.
312 gmx_fio_do_int(fio
,ir
->ePBC
);
313 if ((file_version
<= 15) && (ir
->ePBC
== 2))
315 if (file_version
>= 45) {
316 gmx_fio_do_int(fio
,ir
->bPeriodicMols
);
320 ir
->bPeriodicMols
= TRUE
;
322 ir
->bPeriodicMols
= FALSE
;
326 gmx_fio_do_int(fio
,ir
->ns_type
);
327 gmx_fio_do_int(fio
,ir
->nstlist
);
328 gmx_fio_do_int(fio
,ir
->ndelta
);
329 if (file_version
< 41) {
330 gmx_fio_do_int(fio
,idum
);
331 gmx_fio_do_int(fio
,idum
);
333 if (file_version
>= 45)
334 gmx_fio_do_real(fio
,ir
->rtpi
);
337 gmx_fio_do_int(fio
,ir
->nstcomm
);
338 if (file_version
> 34)
339 gmx_fio_do_int(fio
,ir
->comm_mode
);
340 else if (ir
->nstcomm
< 0)
341 ir
->comm_mode
= ecmANGULAR
;
343 ir
->comm_mode
= ecmLINEAR
;
344 ir
->nstcomm
= abs(ir
->nstcomm
);
346 if(file_version
> 25)
347 gmx_fio_do_int(fio
,ir
->nstcheckpoint
);
351 gmx_fio_do_int(fio
,ir
->nstcgsteep
);
354 gmx_fio_do_int(fio
,ir
->nbfgscorr
);
358 gmx_fio_do_int(fio
,ir
->nstlog
);
359 gmx_fio_do_int(fio
,ir
->nstxout
);
360 gmx_fio_do_int(fio
,ir
->nstvout
);
361 gmx_fio_do_int(fio
,ir
->nstfout
);
362 gmx_fio_do_int(fio
,ir
->nstenergy
);
363 gmx_fio_do_int(fio
,ir
->nstxtcout
);
364 if (file_version
>= 59) {
365 gmx_fio_do_double(fio
,ir
->init_t
);
366 gmx_fio_do_double(fio
,ir
->delta_t
);
368 gmx_fio_do_real(fio
,rdum
);
370 gmx_fio_do_real(fio
,rdum
);
373 gmx_fio_do_real(fio
,ir
->xtcprec
);
374 if (file_version
< 19) {
375 gmx_fio_do_int(fio
,idum
);
376 gmx_fio_do_int(fio
,idum
);
378 if(file_version
< 18)
379 gmx_fio_do_int(fio
,idum
);
380 gmx_fio_do_real(fio
,ir
->rlist
);
381 if (file_version
>= 67) {
382 gmx_fio_do_real(fio
,ir
->rlistlong
);
384 gmx_fio_do_int(fio
,ir
->coulombtype
);
385 if (file_version
< 32 && ir
->coulombtype
== eelRF
)
386 ir
->coulombtype
= eelRF_NEC
;
387 gmx_fio_do_real(fio
,ir
->rcoulomb_switch
);
388 gmx_fio_do_real(fio
,ir
->rcoulomb
);
389 gmx_fio_do_int(fio
,ir
->vdwtype
);
390 gmx_fio_do_real(fio
,ir
->rvdw_switch
);
391 gmx_fio_do_real(fio
,ir
->rvdw
);
392 if (file_version
< 67) {
393 ir
->rlistlong
= max_cutoff(ir
->rlist
,max_cutoff(ir
->rvdw
,ir
->rcoulomb
));
395 gmx_fio_do_int(fio
,ir
->eDispCorr
);
396 gmx_fio_do_real(fio
,ir
->epsilon_r
);
397 if (file_version
>= 37) {
398 gmx_fio_do_real(fio
,ir
->epsilon_rf
);
400 if (EEL_RF(ir
->coulombtype
)) {
401 ir
->epsilon_rf
= ir
->epsilon_r
;
404 ir
->epsilon_rf
= 1.0;
407 if (file_version
>= 29)
408 gmx_fio_do_real(fio
,ir
->tabext
);
412 if(file_version
> 25) {
413 gmx_fio_do_int(fio
,ir
->gb_algorithm
);
414 gmx_fio_do_int(fio
,ir
->nstgbradii
);
415 gmx_fio_do_real(fio
,ir
->rgbradii
);
416 gmx_fio_do_real(fio
,ir
->gb_saltconc
);
417 gmx_fio_do_int(fio
,ir
->implicit_solvent
);
419 ir
->gb_algorithm
=egbSTILL
;
423 ir
->implicit_solvent
=eisNO
;
427 gmx_fio_do_real(fio
,ir
->gb_epsilon_solvent
);
428 gmx_fio_do_real(fio
,ir
->gb_obc_alpha
);
429 gmx_fio_do_real(fio
,ir
->gb_obc_beta
);
430 gmx_fio_do_real(fio
,ir
->gb_obc_gamma
);
433 gmx_fio_do_real(fio
,ir
->gb_dielectric_offset
);
434 gmx_fio_do_int(fio
,ir
->sa_algorithm
);
438 ir
->gb_dielectric_offset
= 0.009;
439 ir
->sa_algorithm
= esaAPPROX
;
441 gmx_fio_do_real(fio
,ir
->sa_surface_tension
);
443 /* Override sa_surface_tension if it is not changed in the mpd-file */
444 if(ir
->sa_surface_tension
<0)
446 if(ir
->gb_algorithm
==egbSTILL
)
448 ir
->sa_surface_tension
= 0.0049 * 100 * CAL2JOULE
;
450 else if(ir
->gb_algorithm
==egbHCT
|| ir
->gb_algorithm
==egbOBC
)
452 ir
->sa_surface_tension
= 0.0054 * 100 * CAL2JOULE
;
459 /* Better use sensible values than insane (0.0) ones... */
460 ir
->gb_epsilon_solvent
= 80;
461 ir
->gb_obc_alpha
= 1.0;
462 ir
->gb_obc_beta
= 0.8;
463 ir
->gb_obc_gamma
= 4.85;
464 ir
->sa_surface_tension
= 2.092;
468 gmx_fio_do_int(fio
,ir
->nkx
);
469 gmx_fio_do_int(fio
,ir
->nky
);
470 gmx_fio_do_int(fio
,ir
->nkz
);
471 gmx_fio_do_int(fio
,ir
->pme_order
);
472 gmx_fio_do_real(fio
,ir
->ewald_rtol
);
474 if (file_version
>=24)
475 gmx_fio_do_int(fio
,ir
->ewald_geometry
);
477 ir
->ewald_geometry
=eewg3D
;
479 if (file_version
<=17) {
480 ir
->epsilon_surface
=0;
481 if (file_version
==17)
482 gmx_fio_do_int(fio
,idum
);
485 gmx_fio_do_real(fio
,ir
->epsilon_surface
);
487 gmx_fio_do_gmx_bool(fio
,ir
->bOptFFT
);
489 gmx_fio_do_gmx_bool(fio
,ir
->bContinuation
);
490 gmx_fio_do_int(fio
,ir
->etc
);
491 /* before version 18, ir->etc was a gmx_bool (ir->btc),
492 * but the values 0 and 1 still mean no and
493 * berendsen temperature coupling, respectively.
495 if (file_version
>= 71)
497 gmx_fio_do_int(fio
,ir
->nsttcouple
);
501 ir
->nsttcouple
= ir
->nstcalcenergy
;
503 if (file_version
<= 15)
505 gmx_fio_do_int(fio
,idum
);
507 if (file_version
<=17)
509 gmx_fio_do_int(fio
,ir
->epct
);
510 if (file_version
<= 15)
514 ir
->epct
= epctSURFACETENSION
;
516 gmx_fio_do_int(fio
,idum
);
519 /* we have removed the NO alternative at the beginning */
523 ir
->epct
=epctISOTROPIC
;
527 ir
->epc
=epcBERENDSEN
;
532 gmx_fio_do_int(fio
,ir
->epc
);
533 gmx_fio_do_int(fio
,ir
->epct
);
535 if (file_version
>= 71)
537 gmx_fio_do_int(fio
,ir
->nstpcouple
);
541 ir
->nstpcouple
= ir
->nstcalcenergy
;
543 gmx_fio_do_real(fio
,ir
->tau_p
);
544 if (file_version
<= 15) {
545 gmx_fio_do_rvec(fio
,vdum
);
546 clear_mat(ir
->ref_p
);
548 ir
->ref_p
[i
][i
] = vdum
[i
];
550 gmx_fio_do_rvec(fio
,ir
->ref_p
[XX
]);
551 gmx_fio_do_rvec(fio
,ir
->ref_p
[YY
]);
552 gmx_fio_do_rvec(fio
,ir
->ref_p
[ZZ
]);
554 if (file_version
<= 15) {
555 gmx_fio_do_rvec(fio
,vdum
);
556 clear_mat(ir
->compress
);
558 ir
->compress
[i
][i
] = vdum
[i
];
561 gmx_fio_do_rvec(fio
,ir
->compress
[XX
]);
562 gmx_fio_do_rvec(fio
,ir
->compress
[YY
]);
563 gmx_fio_do_rvec(fio
,ir
->compress
[ZZ
]);
565 if (file_version
>= 47) {
566 gmx_fio_do_int(fio
,ir
->refcoord_scaling
);
567 gmx_fio_do_rvec(fio
,ir
->posres_com
);
568 gmx_fio_do_rvec(fio
,ir
->posres_comB
);
570 ir
->refcoord_scaling
= erscNO
;
571 clear_rvec(ir
->posres_com
);
572 clear_rvec(ir
->posres_comB
);
574 if(file_version
> 25)
575 gmx_fio_do_int(fio
,ir
->andersen_seed
);
579 if(file_version
< 26) {
580 gmx_fio_do_gmx_bool(fio
,bSimAnn
);
581 gmx_fio_do_real(fio
,zerotemptime
);
584 if (file_version
< 37)
585 gmx_fio_do_real(fio
,rdum
);
587 gmx_fio_do_real(fio
,ir
->shake_tol
);
588 if (file_version
< 54)
589 gmx_fio_do_real(fio
,*fudgeQQ
);
590 gmx_fio_do_int(fio
,ir
->efep
);
591 if (file_version
<= 14 && ir
->efep
> efepNO
)
593 if (file_version
>= 59) {
594 gmx_fio_do_double(fio
, ir
->init_lambda
);
595 gmx_fio_do_double(fio
, ir
->delta_lambda
);
597 gmx_fio_do_real(fio
,rdum
);
598 ir
->init_lambda
= rdum
;
599 gmx_fio_do_real(fio
,rdum
);
600 ir
->delta_lambda
= rdum
;
602 if (file_version
>= 64) {
603 gmx_fio_do_int(fio
,ir
->n_flambda
);
605 snew(ir
->flambda
,ir
->n_flambda
);
607 bDum
=gmx_fio_ndo_double(fio
,ir
->flambda
,ir
->n_flambda
);
612 if (file_version
>= 13)
613 gmx_fio_do_real(fio
,ir
->sc_alpha
);
616 if (file_version
>= 38)
617 gmx_fio_do_int(fio
,ir
->sc_power
);
620 if (file_version
>= 15)
621 gmx_fio_do_real(fio
,ir
->sc_sigma
);
626 if (file_version
>= 71)
628 ir
->sc_sigma_min
= ir
->sc_sigma
;
632 ir
->sc_sigma_min
= 0;
635 if (file_version
>= 64) {
636 gmx_fio_do_int(fio
,ir
->nstdhdl
);
641 if (file_version
>= 73)
643 gmx_fio_do_int(fio
, ir
->separate_dhdl_file
);
644 gmx_fio_do_int(fio
, ir
->dhdl_derivatives
);
648 ir
->separate_dhdl_file
= sepdhdlfileYES
;
649 ir
->dhdl_derivatives
= dhdlderivativesYES
;
652 if (file_version
>= 71)
654 gmx_fio_do_int(fio
,ir
->dh_hist_size
);
655 gmx_fio_do_double(fio
,ir
->dh_hist_spacing
);
659 ir
->dh_hist_size
= 0;
660 ir
->dh_hist_spacing
= 0.1;
662 if (file_version
>= 57) {
663 gmx_fio_do_int(fio
,ir
->eDisre
);
665 gmx_fio_do_int(fio
,ir
->eDisreWeighting
);
666 if (file_version
< 22) {
667 if (ir
->eDisreWeighting
== 0)
668 ir
->eDisreWeighting
= edrwEqual
;
670 ir
->eDisreWeighting
= edrwConservative
;
672 gmx_fio_do_gmx_bool(fio
,ir
->bDisreMixed
);
673 gmx_fio_do_real(fio
,ir
->dr_fc
);
674 gmx_fio_do_real(fio
,ir
->dr_tau
);
675 gmx_fio_do_int(fio
,ir
->nstdisreout
);
676 if (file_version
>= 22) {
677 gmx_fio_do_real(fio
,ir
->orires_fc
);
678 gmx_fio_do_real(fio
,ir
->orires_tau
);
679 gmx_fio_do_int(fio
,ir
->nstorireout
);
685 if(file_version
>= 26) {
686 gmx_fio_do_real(fio
,ir
->dihre_fc
);
687 if (file_version
< 56) {
688 gmx_fio_do_real(fio
,rdum
);
689 gmx_fio_do_int(fio
,idum
);
695 gmx_fio_do_real(fio
,ir
->em_stepsize
);
696 gmx_fio_do_real(fio
,ir
->em_tol
);
697 if (file_version
>= 22)
698 gmx_fio_do_gmx_bool(fio
,ir
->bShakeSOR
);
700 ir
->bShakeSOR
= TRUE
;
701 if (file_version
>= 11)
702 gmx_fio_do_int(fio
,ir
->niter
);
705 fprintf(stderr
,"Note: niter not in run input file, setting it to %d\n",
708 if (file_version
>= 21)
709 gmx_fio_do_real(fio
,ir
->fc_stepsize
);
712 gmx_fio_do_int(fio
,ir
->eConstrAlg
);
713 gmx_fio_do_int(fio
,ir
->nProjOrder
);
714 gmx_fio_do_real(fio
,ir
->LincsWarnAngle
);
715 if (file_version
<= 14)
716 gmx_fio_do_int(fio
,idum
);
717 if (file_version
>=26)
718 gmx_fio_do_int(fio
,ir
->nLincsIter
);
721 fprintf(stderr
,"Note: nLincsIter not in run input file, setting it to %d\n",
724 if (file_version
< 33)
725 gmx_fio_do_real(fio
,bd_temp
);
726 gmx_fio_do_real(fio
,ir
->bd_fric
);
727 gmx_fio_do_int(fio
,ir
->ld_seed
);
728 if (file_version
>= 33) {
730 gmx_fio_do_rvec(fio
,ir
->deform
[i
]);
733 clear_rvec(ir
->deform
[i
]);
735 if (file_version
>= 14)
736 gmx_fio_do_real(fio
,ir
->cos_accel
);
739 gmx_fio_do_int(fio
,ir
->userint1
);
740 gmx_fio_do_int(fio
,ir
->userint2
);
741 gmx_fio_do_int(fio
,ir
->userint3
);
742 gmx_fio_do_int(fio
,ir
->userint4
);
743 gmx_fio_do_real(fio
,ir
->userreal1
);
744 gmx_fio_do_real(fio
,ir
->userreal2
);
745 gmx_fio_do_real(fio
,ir
->userreal3
);
746 gmx_fio_do_real(fio
,ir
->userreal4
);
749 if (file_version
>= 48) {
750 gmx_fio_do_int(fio
,ir
->ePull
);
751 if (ir
->ePull
!= epullNO
) {
754 do_pull(fio
, ir
->pull
,bRead
,file_version
);
761 gmx_fio_do_int(fio
,ir
->opts
.ngtc
);
762 if (file_version
>= 69) {
763 gmx_fio_do_int(fio
,ir
->opts
.nhchainlength
);
765 ir
->opts
.nhchainlength
= 1;
767 gmx_fio_do_int(fio
,ir
->opts
.ngacc
);
768 gmx_fio_do_int(fio
,ir
->opts
.ngfrz
);
769 gmx_fio_do_int(fio
,ir
->opts
.ngener
);
772 snew(ir
->opts
.nrdf
, ir
->opts
.ngtc
);
773 snew(ir
->opts
.ref_t
, ir
->opts
.ngtc
);
774 snew(ir
->opts
.annealing
, ir
->opts
.ngtc
);
775 snew(ir
->opts
.anneal_npoints
, ir
->opts
.ngtc
);
776 snew(ir
->opts
.anneal_time
, ir
->opts
.ngtc
);
777 snew(ir
->opts
.anneal_temp
, ir
->opts
.ngtc
);
778 snew(ir
->opts
.tau_t
, ir
->opts
.ngtc
);
779 snew(ir
->opts
.nFreeze
,ir
->opts
.ngfrz
);
780 snew(ir
->opts
.acc
, ir
->opts
.ngacc
);
781 snew(ir
->opts
.egp_flags
,ir
->opts
.ngener
*ir
->opts
.ngener
);
783 if (ir
->opts
.ngtc
> 0) {
784 if (bRead
&& file_version
<13) {
785 snew(tmp
,ir
->opts
.ngtc
);
786 bDum
=gmx_fio_ndo_int(fio
,tmp
, ir
->opts
.ngtc
);
787 for(i
=0; i
<ir
->opts
.ngtc
; i
++)
788 ir
->opts
.nrdf
[i
] = tmp
[i
];
791 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.nrdf
, ir
->opts
.ngtc
);
793 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.ref_t
,ir
->opts
.ngtc
);
794 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.tau_t
,ir
->opts
.ngtc
);
795 if (file_version
<33 && ir
->eI
==eiBD
) {
796 for(i
=0; i
<ir
->opts
.ngtc
; i
++)
797 ir
->opts
.tau_t
[i
] = bd_temp
;
800 if (ir
->opts
.ngfrz
> 0)
801 bDum
=gmx_fio_ndo_ivec(fio
,ir
->opts
.nFreeze
,ir
->opts
.ngfrz
);
802 if (ir
->opts
.ngacc
> 0)
803 gmx_fio_ndo_rvec(fio
,ir
->opts
.acc
,ir
->opts
.ngacc
);
804 if (file_version
>= 12)
805 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.egp_flags
,
806 ir
->opts
.ngener
*ir
->opts
.ngener
);
808 if(bRead
&& file_version
< 26) {
809 for(i
=0;i
<ir
->opts
.ngtc
;i
++) {
811 ir
->opts
.annealing
[i
] = eannSINGLE
;
812 ir
->opts
.anneal_npoints
[i
] = 2;
813 snew(ir
->opts
.anneal_time
[i
],2);
814 snew(ir
->opts
.anneal_temp
[i
],2);
815 /* calculate the starting/ending temperatures from reft, zerotemptime, and nsteps */
816 finish_t
= ir
->init_t
+ ir
->nsteps
* ir
->delta_t
;
817 init_temp
= ir
->opts
.ref_t
[i
]*(1-ir
->init_t
/zerotemptime
);
818 finish_temp
= ir
->opts
.ref_t
[i
]*(1-finish_t
/zerotemptime
);
819 ir
->opts
.anneal_time
[i
][0] = ir
->init_t
;
820 ir
->opts
.anneal_time
[i
][1] = finish_t
;
821 ir
->opts
.anneal_temp
[i
][0] = init_temp
;
822 ir
->opts
.anneal_temp
[i
][1] = finish_temp
;
824 ir
->opts
.annealing
[i
] = eannNO
;
825 ir
->opts
.anneal_npoints
[i
] = 0;
829 /* file version 26 or later */
830 /* First read the lists with annealing and npoints for each group */
831 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.annealing
,ir
->opts
.ngtc
);
832 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.anneal_npoints
,ir
->opts
.ngtc
);
833 for(j
=0;j
<(ir
->opts
.ngtc
);j
++) {
834 k
=ir
->opts
.anneal_npoints
[j
];
836 snew(ir
->opts
.anneal_time
[j
],k
);
837 snew(ir
->opts
.anneal_temp
[j
],k
);
839 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.anneal_time
[j
],k
);
840 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.anneal_temp
[j
],k
);
844 if (file_version
>= 45) {
845 gmx_fio_do_int(fio
,ir
->nwall
);
846 gmx_fio_do_int(fio
,ir
->wall_type
);
847 if (file_version
>= 50)
848 gmx_fio_do_real(fio
,ir
->wall_r_linpot
);
850 ir
->wall_r_linpot
= -1;
851 gmx_fio_do_int(fio
,ir
->wall_atomtype
[0]);
852 gmx_fio_do_int(fio
,ir
->wall_atomtype
[1]);
853 gmx_fio_do_real(fio
,ir
->wall_density
[0]);
854 gmx_fio_do_real(fio
,ir
->wall_density
[1]);
855 gmx_fio_do_real(fio
,ir
->wall_ewald_zfac
);
859 ir
->wall_atomtype
[0] = -1;
860 ir
->wall_atomtype
[1] = -1;
861 ir
->wall_density
[0] = 0;
862 ir
->wall_density
[1] = 0;
863 ir
->wall_ewald_zfac
= 3;
865 /* Cosine stuff for electric fields */
866 for(j
=0; (j
<DIM
); j
++) {
867 gmx_fio_do_int(fio
,ir
->ex
[j
].n
);
868 gmx_fio_do_int(fio
,ir
->et
[j
].n
);
870 snew(ir
->ex
[j
].a
, ir
->ex
[j
].n
);
871 snew(ir
->ex
[j
].phi
,ir
->ex
[j
].n
);
872 snew(ir
->et
[j
].a
, ir
->et
[j
].n
);
873 snew(ir
->et
[j
].phi
,ir
->et
[j
].n
);
875 bDum
=gmx_fio_ndo_real(fio
,ir
->ex
[j
].a
, ir
->ex
[j
].n
);
876 bDum
=gmx_fio_ndo_real(fio
,ir
->ex
[j
].phi
,ir
->ex
[j
].n
);
877 bDum
=gmx_fio_ndo_real(fio
,ir
->et
[j
].a
, ir
->et
[j
].n
);
878 bDum
=gmx_fio_ndo_real(fio
,ir
->et
[j
].phi
,ir
->et
[j
].n
);
882 if(file_version
>=39){
883 gmx_fio_do_gmx_bool(fio
,ir
->bQMMM
);
884 gmx_fio_do_int(fio
,ir
->QMMMscheme
);
885 gmx_fio_do_real(fio
,ir
->scalefactor
);
886 gmx_fio_do_int(fio
,ir
->opts
.ngQM
);
888 snew(ir
->opts
.QMmethod
, ir
->opts
.ngQM
);
889 snew(ir
->opts
.QMbasis
, ir
->opts
.ngQM
);
890 snew(ir
->opts
.QMcharge
, ir
->opts
.ngQM
);
891 snew(ir
->opts
.QMmult
, ir
->opts
.ngQM
);
892 snew(ir
->opts
.bSH
, ir
->opts
.ngQM
);
893 snew(ir
->opts
.CASorbitals
, ir
->opts
.ngQM
);
894 snew(ir
->opts
.CASelectrons
,ir
->opts
.ngQM
);
895 snew(ir
->opts
.SAon
, ir
->opts
.ngQM
);
896 snew(ir
->opts
.SAoff
, ir
->opts
.ngQM
);
897 snew(ir
->opts
.SAsteps
, ir
->opts
.ngQM
);
898 snew(ir
->opts
.bOPT
, ir
->opts
.ngQM
);
899 snew(ir
->opts
.bTS
, ir
->opts
.ngQM
);
901 if (ir
->opts
.ngQM
> 0) {
902 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMmethod
,ir
->opts
.ngQM
);
903 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMbasis
,ir
->opts
.ngQM
);
904 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMcharge
,ir
->opts
.ngQM
);
905 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.QMmult
,ir
->opts
.ngQM
);
906 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bSH
,ir
->opts
.ngQM
);
907 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.CASorbitals
,ir
->opts
.ngQM
);
908 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.CASelectrons
,ir
->opts
.ngQM
);
909 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.SAon
,ir
->opts
.ngQM
);
910 bDum
=gmx_fio_ndo_real(fio
,ir
->opts
.SAoff
,ir
->opts
.ngQM
);
911 bDum
=gmx_fio_ndo_int(fio
,ir
->opts
.SAsteps
,ir
->opts
.ngQM
);
912 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bOPT
,ir
->opts
.ngQM
);
913 bDum
=gmx_fio_ndo_gmx_bool(fio
,ir
->opts
.bTS
,ir
->opts
.ngQM
);
915 /* end of QMMM stuff */
920 static void do_harm(t_fileio
*fio
, t_iparams
*iparams
,gmx_bool bRead
)
922 gmx_fio_do_real(fio
,iparams
->harmonic
.rA
);
923 gmx_fio_do_real(fio
,iparams
->harmonic
.krA
);
924 gmx_fio_do_real(fio
,iparams
->harmonic
.rB
);
925 gmx_fio_do_real(fio
,iparams
->harmonic
.krB
);
928 void do_iparams(t_fileio
*fio
, t_functype ftype
,t_iparams
*iparams
,
929 gmx_bool bRead
, int file_version
)
936 gmx_fio_set_comment(fio
, interaction_function
[ftype
].name
);
944 do_harm(fio
, iparams
,bRead
);
945 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && bRead
) {
946 /* Correct incorrect storage of parameters */
947 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
948 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
952 gmx_fio_do_real(fio
,iparams
->fene
.bm
);
953 gmx_fio_do_real(fio
,iparams
->fene
.kb
);
956 gmx_fio_do_real(fio
,iparams
->restraint
.lowA
);
957 gmx_fio_do_real(fio
,iparams
->restraint
.up1A
);
958 gmx_fio_do_real(fio
,iparams
->restraint
.up2A
);
959 gmx_fio_do_real(fio
,iparams
->restraint
.kA
);
960 gmx_fio_do_real(fio
,iparams
->restraint
.lowB
);
961 gmx_fio_do_real(fio
,iparams
->restraint
.up1B
);
962 gmx_fio_do_real(fio
,iparams
->restraint
.up2B
);
963 gmx_fio_do_real(fio
,iparams
->restraint
.kB
);
969 gmx_fio_do_real(fio
,iparams
->tab
.kA
);
970 gmx_fio_do_int(fio
,iparams
->tab
.table
);
971 gmx_fio_do_real(fio
,iparams
->tab
.kB
);
973 case F_CROSS_BOND_BONDS
:
974 gmx_fio_do_real(fio
,iparams
->cross_bb
.r1e
);
975 gmx_fio_do_real(fio
,iparams
->cross_bb
.r2e
);
976 gmx_fio_do_real(fio
,iparams
->cross_bb
.krr
);
978 case F_CROSS_BOND_ANGLES
:
979 gmx_fio_do_real(fio
,iparams
->cross_ba
.r1e
);
980 gmx_fio_do_real(fio
,iparams
->cross_ba
.r2e
);
981 gmx_fio_do_real(fio
,iparams
->cross_ba
.r3e
);
982 gmx_fio_do_real(fio
,iparams
->cross_ba
.krt
);
985 gmx_fio_do_real(fio
,iparams
->u_b
.theta
);
986 gmx_fio_do_real(fio
,iparams
->u_b
.ktheta
);
987 gmx_fio_do_real(fio
,iparams
->u_b
.r13
);
988 gmx_fio_do_real(fio
,iparams
->u_b
.kUB
);
990 case F_QUARTIC_ANGLES
:
991 gmx_fio_do_real(fio
,iparams
->qangle
.theta
);
992 bDum
=gmx_fio_ndo_real(fio
,iparams
->qangle
.c
,5);
995 gmx_fio_do_real(fio
,iparams
->bham
.a
);
996 gmx_fio_do_real(fio
,iparams
->bham
.b
);
997 gmx_fio_do_real(fio
,iparams
->bham
.c
);
1000 gmx_fio_do_real(fio
,iparams
->morse
.b0
);
1001 gmx_fio_do_real(fio
,iparams
->morse
.cb
);
1002 gmx_fio_do_real(fio
,iparams
->morse
.beta
);
1005 gmx_fio_do_real(fio
,iparams
->cubic
.b0
);
1006 gmx_fio_do_real(fio
,iparams
->cubic
.kb
);
1007 gmx_fio_do_real(fio
,iparams
->cubic
.kcub
);
1011 case F_POLARIZATION
:
1012 gmx_fio_do_real(fio
,iparams
->polarize
.alpha
);
1015 if (file_version
< 31)
1016 gmx_fatal(FARGS
,"Old tpr files with water_polarization not supported. Make a new.");
1017 gmx_fio_do_real(fio
,iparams
->wpol
.al_x
);
1018 gmx_fio_do_real(fio
,iparams
->wpol
.al_y
);
1019 gmx_fio_do_real(fio
,iparams
->wpol
.al_z
);
1020 gmx_fio_do_real(fio
,iparams
->wpol
.rOH
);
1021 gmx_fio_do_real(fio
,iparams
->wpol
.rHH
);
1022 gmx_fio_do_real(fio
,iparams
->wpol
.rOD
);
1025 gmx_fio_do_real(fio
,iparams
->thole
.a
);
1026 gmx_fio_do_real(fio
,iparams
->thole
.alpha1
);
1027 gmx_fio_do_real(fio
,iparams
->thole
.alpha2
);
1028 gmx_fio_do_real(fio
,iparams
->thole
.rfac
);
1031 gmx_fio_do_real(fio
,iparams
->lj
.c6
);
1032 gmx_fio_do_real(fio
,iparams
->lj
.c12
);
1035 gmx_fio_do_real(fio
,iparams
->lj14
.c6A
);
1036 gmx_fio_do_real(fio
,iparams
->lj14
.c12A
);
1037 gmx_fio_do_real(fio
,iparams
->lj14
.c6B
);
1038 gmx_fio_do_real(fio
,iparams
->lj14
.c12B
);
1041 gmx_fio_do_real(fio
,iparams
->ljc14
.fqq
);
1042 gmx_fio_do_real(fio
,iparams
->ljc14
.qi
);
1043 gmx_fio_do_real(fio
,iparams
->ljc14
.qj
);
1044 gmx_fio_do_real(fio
,iparams
->ljc14
.c6
);
1045 gmx_fio_do_real(fio
,iparams
->ljc14
.c12
);
1047 case F_LJC_PAIRS_NB
:
1048 gmx_fio_do_real(fio
,iparams
->ljcnb
.qi
);
1049 gmx_fio_do_real(fio
,iparams
->ljcnb
.qj
);
1050 gmx_fio_do_real(fio
,iparams
->ljcnb
.c6
);
1051 gmx_fio_do_real(fio
,iparams
->ljcnb
.c12
);
1057 gmx_fio_do_real(fio
,iparams
->pdihs
.phiA
);
1058 gmx_fio_do_real(fio
,iparams
->pdihs
.cpA
);
1059 if ((ftype
== F_ANGRES
|| ftype
== F_ANGRESZ
) && file_version
< 42) {
1060 /* Read the incorrectly stored multiplicity */
1061 gmx_fio_do_real(fio
,iparams
->harmonic
.rB
);
1062 gmx_fio_do_real(fio
,iparams
->harmonic
.krB
);
1063 iparams
->pdihs
.phiB
= iparams
->pdihs
.phiA
;
1064 iparams
->pdihs
.cpB
= iparams
->pdihs
.cpA
;
1066 gmx_fio_do_real(fio
,iparams
->pdihs
.phiB
);
1067 gmx_fio_do_real(fio
,iparams
->pdihs
.cpB
);
1068 gmx_fio_do_int(fio
,iparams
->pdihs
.mult
);
1072 gmx_fio_do_int(fio
,iparams
->disres
.label
);
1073 gmx_fio_do_int(fio
,iparams
->disres
.type
);
1074 gmx_fio_do_real(fio
,iparams
->disres
.low
);
1075 gmx_fio_do_real(fio
,iparams
->disres
.up1
);
1076 gmx_fio_do_real(fio
,iparams
->disres
.up2
);
1077 gmx_fio_do_real(fio
,iparams
->disres
.kfac
);
1080 gmx_fio_do_int(fio
,iparams
->orires
.ex
);
1081 gmx_fio_do_int(fio
,iparams
->orires
.label
);
1082 gmx_fio_do_int(fio
,iparams
->orires
.power
);
1083 gmx_fio_do_real(fio
,iparams
->orires
.c
);
1084 gmx_fio_do_real(fio
,iparams
->orires
.obs
);
1085 gmx_fio_do_real(fio
,iparams
->orires
.kfac
);
1088 gmx_fio_do_int(fio
,iparams
->dihres
.power
);
1089 gmx_fio_do_int(fio
,iparams
->dihres
.label
);
1090 gmx_fio_do_real(fio
,iparams
->dihres
.phi
);
1091 gmx_fio_do_real(fio
,iparams
->dihres
.dphi
);
1092 gmx_fio_do_real(fio
,iparams
->dihres
.kfac
);
1095 gmx_fio_do_rvec(fio
,iparams
->posres
.pos0A
);
1096 gmx_fio_do_rvec(fio
,iparams
->posres
.fcA
);
1097 if (bRead
&& file_version
< 27) {
1098 copy_rvec(iparams
->posres
.pos0A
,iparams
->posres
.pos0B
);
1099 copy_rvec(iparams
->posres
.fcA
,iparams
->posres
.fcB
);
1101 gmx_fio_do_rvec(fio
,iparams
->posres
.pos0B
);
1102 gmx_fio_do_rvec(fio
,iparams
->posres
.fcB
);
1106 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcA
,NR_RBDIHS
);
1107 if(file_version
>=25)
1108 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcB
,NR_RBDIHS
);
1111 /* Fourier dihedrals are internally represented
1112 * as Ryckaert-Bellemans since those are faster to compute.
1114 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcA
, NR_RBDIHS
);
1115 bDum
=gmx_fio_ndo_real(fio
,iparams
->rbdihs
.rbcB
, NR_RBDIHS
);
1119 gmx_fio_do_real(fio
,iparams
->constr
.dA
);
1120 gmx_fio_do_real(fio
,iparams
->constr
.dB
);
1123 gmx_fio_do_real(fio
,iparams
->settle
.doh
);
1124 gmx_fio_do_real(fio
,iparams
->settle
.dhh
);
1127 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1132 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1133 gmx_fio_do_real(fio
,iparams
->vsite
.b
);
1138 gmx_fio_do_real(fio
,iparams
->vsite
.a
);
1139 gmx_fio_do_real(fio
,iparams
->vsite
.b
);
1140 gmx_fio_do_real(fio
,iparams
->vsite
.c
);
1143 gmx_fio_do_int(fio
,iparams
->vsiten
.n
);
1144 gmx_fio_do_real(fio
,iparams
->vsiten
.a
);
1149 /* We got rid of some parameters in version 68 */
1150 if(bRead
&& file_version
<68)
1152 gmx_fio_do_real(fio
,rdum
);
1153 gmx_fio_do_real(fio
,rdum
);
1154 gmx_fio_do_real(fio
,rdum
);
1155 gmx_fio_do_real(fio
,rdum
);
1157 gmx_fio_do_real(fio
,iparams
->gb
.sar
);
1158 gmx_fio_do_real(fio
,iparams
->gb
.st
);
1159 gmx_fio_do_real(fio
,iparams
->gb
.pi
);
1160 gmx_fio_do_real(fio
,iparams
->gb
.gbr
);
1161 gmx_fio_do_real(fio
,iparams
->gb
.bmlt
);
1164 gmx_fio_do_int(fio
,iparams
->cmap
.cmapA
);
1165 gmx_fio_do_int(fio
,iparams
->cmap
.cmapB
);
1168 gmx_fatal(FARGS
,"unknown function type %d (%s) in %s line %d",
1170 ftype
,interaction_function
[ftype
].name
,__FILE__
,__LINE__
);
1173 gmx_fio_unset_comment(fio
);
1176 static void do_ilist(t_fileio
*fio
, t_ilist
*ilist
,gmx_bool bRead
,int file_version
,
1183 gmx_fio_set_comment(fio
, interaction_function
[ftype
].name
);
1185 if (file_version
< 44) {
1186 for(i
=0; i
<MAXNODES
; i
++)
1187 gmx_fio_do_int(fio
,idum
);
1189 gmx_fio_do_int(fio
,ilist
->nr
);
1191 snew(ilist
->iatoms
,ilist
->nr
);
1192 bDum
=gmx_fio_ndo_int(fio
,ilist
->iatoms
,ilist
->nr
);
1194 gmx_fio_unset_comment(fio
);
1197 static void do_ffparams(t_fileio
*fio
, gmx_ffparams_t
*ffparams
,
1198 gmx_bool bRead
, int file_version
)
1204 gmx_fio_do_int(fio
,ffparams
->atnr
);
1205 if (file_version
< 57) {
1206 gmx_fio_do_int(fio
,idum
);
1208 gmx_fio_do_int(fio
,ffparams
->ntypes
);
1210 fprintf(debug
,"ffparams->atnr = %d, ntypes = %d\n",
1211 ffparams
->atnr
,ffparams
->ntypes
);
1213 snew(ffparams
->functype
,ffparams
->ntypes
);
1214 snew(ffparams
->iparams
,ffparams
->ntypes
);
1216 /* Read/write all the function types */
1217 bDum
=gmx_fio_ndo_int(fio
,ffparams
->functype
,ffparams
->ntypes
);
1219 pr_ivec(debug
,0,"functype",ffparams
->functype
,ffparams
->ntypes
,TRUE
);
1221 if (file_version
>= 66) {
1222 gmx_fio_do_double(fio
,ffparams
->reppow
);
1224 ffparams
->reppow
= 12.0;
1227 if (file_version
>= 57) {
1228 gmx_fio_do_real(fio
,ffparams
->fudgeQQ
);
1231 /* Check whether all these function types are supported by the code.
1232 * In practice the code is backwards compatible, which means that the
1233 * numbering may have to be altered from old numbering to new numbering
1235 for (i
=0; (i
<ffparams
->ntypes
); i
++) {
1237 /* Loop over file versions */
1238 for (k
=0; (k
<NFTUPD
); k
++)
1239 /* Compare the read file_version to the update table */
1240 if ((file_version
< ftupd
[k
].fvnr
) &&
1241 (ffparams
->functype
[i
] >= ftupd
[k
].ftype
)) {
1242 ffparams
->functype
[i
] += 1;
1244 fprintf(debug
,"Incrementing function type %d to %d (due to %s)\n",
1245 i
,ffparams
->functype
[i
],
1246 interaction_function
[ftupd
[k
].ftype
].longname
);
1251 do_iparams(fio
, ffparams
->functype
[i
],&ffparams
->iparams
[i
],bRead
,
1254 pr_iparams(debug
,ffparams
->functype
[i
],&ffparams
->iparams
[i
]);
1258 static void do_ilists(t_fileio
*fio
, t_ilist
*ilist
,gmx_bool bRead
,
1261 int i
,j
,renum
[F_NRE
];
1262 gmx_bool bDum
=TRUE
,bClear
;
1265 for(j
=0; (j
<F_NRE
); j
++) {
1268 for (k
=0; k
<NFTUPD
; k
++)
1269 if ((file_version
< ftupd
[k
].fvnr
) && (j
== ftupd
[k
].ftype
))
1273 ilist
[j
].iatoms
= NULL
;
1275 do_ilist(fio
, &ilist
[j
],bRead
,file_version
,j
);
1278 if (bRead && gmx_debug_at)
1279 pr_ilist(debug,0,interaction_function[j].longname,
1280 functype,&ilist[j],TRUE);
1285 static void do_idef(t_fileio
*fio
, gmx_ffparams_t
*ffparams
,gmx_moltype_t
*molt
,
1286 gmx_bool bRead
, int file_version
)
1288 do_ffparams(fio
, ffparams
,bRead
,file_version
);
1290 if (file_version
>= 54) {
1291 gmx_fio_do_real(fio
,ffparams
->fudgeQQ
);
1294 do_ilists(fio
, molt
->ilist
,bRead
,file_version
);
1297 static void do_block(t_fileio
*fio
, t_block
*block
,gmx_bool bRead
,int file_version
)
1299 int i
,idum
,dum_nra
,*dum_a
;
1302 if (file_version
< 44)
1303 for(i
=0; i
<MAXNODES
; i
++)
1304 gmx_fio_do_int(fio
,idum
);
1305 gmx_fio_do_int(fio
,block
->nr
);
1306 if (file_version
< 51)
1307 gmx_fio_do_int(fio
,dum_nra
);
1309 block
->nalloc_index
= block
->nr
+1;
1310 snew(block
->index
,block
->nalloc_index
);
1312 bDum
=gmx_fio_ndo_int(fio
,block
->index
,block
->nr
+1);
1314 if (file_version
< 51 && dum_nra
> 0) {
1315 snew(dum_a
,dum_nra
);
1316 bDum
=gmx_fio_ndo_int(fio
,dum_a
,dum_nra
);
1321 static void do_blocka(t_fileio
*fio
, t_blocka
*block
,gmx_bool bRead
,
1327 if (file_version
< 44)
1328 for(i
=0; i
<MAXNODES
; i
++)
1329 gmx_fio_do_int(fio
,idum
);
1330 gmx_fio_do_int(fio
,block
->nr
);
1331 gmx_fio_do_int(fio
,block
->nra
);
1333 block
->nalloc_index
= block
->nr
+1;
1334 snew(block
->index
,block
->nalloc_index
);
1335 block
->nalloc_a
= block
->nra
;
1336 snew(block
->a
,block
->nalloc_a
);
1338 bDum
=gmx_fio_ndo_int(fio
,block
->index
,block
->nr
+1);
1339 bDum
=gmx_fio_ndo_int(fio
,block
->a
,block
->nra
);
1342 static void do_atom(t_fileio
*fio
, t_atom
*atom
,int ngrp
,gmx_bool bRead
,
1343 int file_version
, gmx_groups_t
*groups
,int atnr
)
1347 gmx_fio_do_real(fio
,atom
->m
);
1348 gmx_fio_do_real(fio
,atom
->q
);
1349 gmx_fio_do_real(fio
,atom
->mB
);
1350 gmx_fio_do_real(fio
,atom
->qB
);
1351 gmx_fio_do_ushort(fio
, atom
->type
);
1352 gmx_fio_do_ushort(fio
, atom
->typeB
);
1353 gmx_fio_do_int(fio
,atom
->ptype
);
1354 gmx_fio_do_int(fio
,atom
->resind
);
1355 if (file_version
>= 52)
1356 gmx_fio_do_int(fio
,atom
->atomnumber
);
1358 atom
->atomnumber
= NOTSET
;
1359 if (file_version
< 23)
1361 else if (file_version
< 39)
1366 if (file_version
< 57) {
1367 unsigned char uchar
[egcNR
];
1368 gmx_fio_ndo_uchar(fio
,uchar
,myngrp
);
1369 for(i
=myngrp
; (i
<ngrp
); i
++) {
1372 /* Copy the old data format to the groups struct */
1373 for(i
=0; i
<ngrp
; i
++) {
1374 groups
->grpnr
[i
][atnr
] = uchar
[i
];
1379 static void do_grps(t_fileio
*fio
, int ngrp
,t_grps grps
[],gmx_bool bRead
,
1385 if (file_version
< 23)
1387 else if (file_version
< 39)
1392 for(j
=0; (j
<ngrp
); j
++) {
1394 gmx_fio_do_int(fio
,grps
[j
].nr
);
1396 snew(grps
[j
].nm_ind
,grps
[j
].nr
);
1397 bDum
=gmx_fio_ndo_int(fio
,grps
[j
].nm_ind
,grps
[j
].nr
);
1401 snew(grps
[j
].nm_ind
,grps
[j
].nr
);
1406 static void do_symstr(t_fileio
*fio
, char ***nm
,gmx_bool bRead
,t_symtab
*symtab
)
1411 gmx_fio_do_int(fio
,ls
);
1412 *nm
= get_symtab_handle(symtab
,ls
);
1415 ls
= lookup_symtab(symtab
,*nm
);
1416 gmx_fio_do_int(fio
,ls
);
1420 static void do_strstr(t_fileio
*fio
, int nstr
,char ***nm
,gmx_bool bRead
,
1425 for (j
=0; (j
<nstr
); j
++)
1426 do_symstr(fio
, &(nm
[j
]),bRead
,symtab
);
1429 static void do_resinfo(t_fileio
*fio
, int n
,t_resinfo
*ri
,gmx_bool bRead
,
1430 t_symtab
*symtab
, int file_version
)
1434 for (j
=0; (j
<n
); j
++) {
1435 do_symstr(fio
, &(ri
[j
].name
),bRead
,symtab
);
1436 if (file_version
>= 63) {
1437 gmx_fio_do_int(fio
,ri
[j
].nr
);
1438 gmx_fio_do_uchar(fio
, ri
[j
].ic
);
1446 static void do_atoms(t_fileio
*fio
, t_atoms
*atoms
,gmx_bool bRead
,t_symtab
*symtab
,
1448 gmx_groups_t
*groups
)
1452 gmx_fio_do_int(fio
,atoms
->nr
);
1453 gmx_fio_do_int(fio
,atoms
->nres
);
1454 if (file_version
< 57) {
1455 gmx_fio_do_int(fio
,groups
->ngrpname
);
1456 for(i
=0; i
<egcNR
; i
++) {
1457 groups
->ngrpnr
[i
] = atoms
->nr
;
1458 snew(groups
->grpnr
[i
],groups
->ngrpnr
[i
]);
1462 snew(atoms
->atom
,atoms
->nr
);
1463 snew(atoms
->atomname
,atoms
->nr
);
1464 snew(atoms
->atomtype
,atoms
->nr
);
1465 snew(atoms
->atomtypeB
,atoms
->nr
);
1466 snew(atoms
->resinfo
,atoms
->nres
);
1467 if (file_version
< 57) {
1468 snew(groups
->grpname
,groups
->ngrpname
);
1470 atoms
->pdbinfo
= NULL
;
1472 for(i
=0; (i
<atoms
->nr
); i
++) {
1473 do_atom(fio
, &atoms
->atom
[i
],egcNR
,bRead
, file_version
,groups
,i
);
1475 do_strstr(fio
, atoms
->nr
,atoms
->atomname
,bRead
,symtab
);
1476 if (bRead
&& (file_version
<= 20)) {
1477 for(i
=0; i
<atoms
->nr
; i
++) {
1478 atoms
->atomtype
[i
] = put_symtab(symtab
,"?");
1479 atoms
->atomtypeB
[i
] = put_symtab(symtab
,"?");
1482 do_strstr(fio
, atoms
->nr
,atoms
->atomtype
,bRead
,symtab
);
1483 do_strstr(fio
, atoms
->nr
,atoms
->atomtypeB
,bRead
,symtab
);
1485 do_resinfo(fio
, atoms
->nres
,atoms
->resinfo
,bRead
,symtab
,file_version
);
1487 if (file_version
< 57) {
1488 do_strstr(fio
, groups
->ngrpname
,groups
->grpname
,bRead
,symtab
);
1490 do_grps(fio
, egcNR
,groups
->grps
,bRead
,file_version
);
1494 static void do_groups(t_fileio
*fio
, gmx_groups_t
*groups
,
1495 gmx_bool bRead
,t_symtab
*symtab
,
1501 do_grps(fio
, egcNR
,groups
->grps
,bRead
,file_version
);
1502 gmx_fio_do_int(fio
,groups
->ngrpname
);
1504 snew(groups
->grpname
,groups
->ngrpname
);
1506 do_strstr(fio
, groups
->ngrpname
,groups
->grpname
,bRead
,symtab
);
1507 for(g
=0; g
<egcNR
; g
++) {
1508 gmx_fio_do_int(fio
,groups
->ngrpnr
[g
]);
1509 if (groups
->ngrpnr
[g
] == 0) {
1511 groups
->grpnr
[g
] = NULL
;
1515 snew(groups
->grpnr
[g
],groups
->ngrpnr
[g
]);
1517 bDum
=gmx_fio_ndo_uchar(fio
, groups
->grpnr
[g
],groups
->ngrpnr
[g
]);
1522 static void do_atomtypes(t_fileio
*fio
, t_atomtypes
*atomtypes
,gmx_bool bRead
,
1523 t_symtab
*symtab
,int file_version
)
1526 gmx_bool bDum
= TRUE
;
1528 if (file_version
> 25) {
1529 gmx_fio_do_int(fio
,atomtypes
->nr
);
1532 snew(atomtypes
->radius
,j
);
1533 snew(atomtypes
->vol
,j
);
1534 snew(atomtypes
->surftens
,j
);
1535 snew(atomtypes
->atomnumber
,j
);
1536 snew(atomtypes
->gb_radius
,j
);
1537 snew(atomtypes
->S_hct
,j
);
1539 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->radius
,j
);
1540 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->vol
,j
);
1541 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->surftens
,j
);
1542 if(file_version
>= 40)
1544 bDum
=gmx_fio_ndo_int(fio
,atomtypes
->atomnumber
,j
);
1546 if(file_version
>= 60)
1548 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->gb_radius
,j
);
1549 bDum
=gmx_fio_ndo_real(fio
,atomtypes
->S_hct
,j
);
1552 /* File versions prior to 26 cannot do GBSA,
1553 * so they dont use this structure
1556 atomtypes
->radius
= NULL
;
1557 atomtypes
->vol
= NULL
;
1558 atomtypes
->surftens
= NULL
;
1559 atomtypes
->atomnumber
= NULL
;
1560 atomtypes
->gb_radius
= NULL
;
1561 atomtypes
->S_hct
= NULL
;
1565 static void do_symtab(t_fileio
*fio
, t_symtab
*symtab
,gmx_bool bRead
)
1571 gmx_fio_do_int(fio
,symtab
->nr
);
1574 snew(symtab
->symbuf
,1);
1575 symbuf
= symtab
->symbuf
;
1576 symbuf
->bufsize
= nr
;
1577 snew(symbuf
->buf
,nr
);
1578 for (i
=0; (i
<nr
); i
++) {
1579 gmx_fio_do_string(fio
,buf
);
1580 symbuf
->buf
[i
]=strdup(buf
);
1584 symbuf
= symtab
->symbuf
;
1585 while (symbuf
!=NULL
) {
1586 for (i
=0; (i
<symbuf
->bufsize
) && (i
<nr
); i
++)
1587 gmx_fio_do_string(fio
,symbuf
->buf
[i
]);
1589 symbuf
=symbuf
->next
;
1592 gmx_fatal(FARGS
,"nr of symtab strings left: %d",nr
);
1596 static void do_cmap(t_fileio
*fio
, gmx_cmap_t
*cmap_grid
, gmx_bool bRead
)
1598 int i
,j
,ngrid
,gs
,nelem
;
1600 gmx_fio_do_int(fio
,cmap_grid
->ngrid
);
1601 gmx_fio_do_int(fio
,cmap_grid
->grid_spacing
);
1603 ngrid
= cmap_grid
->ngrid
;
1604 gs
= cmap_grid
->grid_spacing
;
1609 snew(cmap_grid
->cmapdata
,ngrid
);
1611 for(i
=0;i
<cmap_grid
->ngrid
;i
++)
1613 snew(cmap_grid
->cmapdata
[i
].cmap
,4*nelem
);
1617 for(i
=0;i
<cmap_grid
->ngrid
;i
++)
1619 for(j
=0;j
<nelem
;j
++)
1621 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4]);
1622 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+1]);
1623 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+2]);
1624 gmx_fio_do_real(fio
,cmap_grid
->cmapdata
[i
].cmap
[j
*4+3]);
1630 void tpx_make_chain_identifiers(t_atoms
*atoms
,t_block
*mols
)
1636 /* We always assign a new chain number, but save the chain id characters
1637 * for larger molecules.
1639 #define CHAIN_MIN_ATOMS 15
1643 for(m
=0; m
<mols
->nr
; m
++)
1646 a1
=mols
->index
[m
+1];
1647 if ((a1
-a0
>= CHAIN_MIN_ATOMS
) && (chainid
<= 'Z'))
1656 for(a
=a0
; a
<a1
; a
++)
1658 atoms
->resinfo
[atoms
->atom
[a
].resind
].chainnum
= chainnum
;
1659 atoms
->resinfo
[atoms
->atom
[a
].resind
].chainid
= c
;
1664 /* Blank out the chain id if there was only one chain */
1667 for(r
=0; r
<atoms
->nres
; r
++)
1669 atoms
->resinfo
[r
].chainid
= ' ';
1674 static void do_moltype(t_fileio
*fio
, gmx_moltype_t
*molt
,gmx_bool bRead
,
1675 t_symtab
*symtab
, int file_version
,
1676 gmx_groups_t
*groups
)
1680 if (file_version
>= 57) {
1681 do_symstr(fio
, &(molt
->name
),bRead
,symtab
);
1684 do_atoms(fio
, &molt
->atoms
, bRead
, symtab
, file_version
, groups
);
1686 if (bRead
&& gmx_debug_at
) {
1687 pr_atoms(debug
,0,"atoms",&molt
->atoms
,TRUE
);
1690 if (file_version
>= 57) {
1691 do_ilists(fio
, molt
->ilist
,bRead
,file_version
);
1693 do_block(fio
, &molt
->cgs
,bRead
,file_version
);
1694 if (bRead
&& gmx_debug_at
) {
1695 pr_block(debug
,0,"cgs",&molt
->cgs
,TRUE
);
1699 /* This used to be in the atoms struct */
1700 do_blocka(fio
, &molt
->excls
, bRead
, file_version
);
1703 static void do_molblock(t_fileio
*fio
, gmx_molblock_t
*molb
,gmx_bool bRead
,
1708 gmx_fio_do_int(fio
,molb
->type
);
1709 gmx_fio_do_int(fio
,molb
->nmol
);
1710 gmx_fio_do_int(fio
,molb
->natoms_mol
);
1711 /* Position restraint coordinates */
1712 gmx_fio_do_int(fio
,molb
->nposres_xA
);
1713 if (molb
->nposres_xA
> 0) {
1715 snew(molb
->posres_xA
,molb
->nposres_xA
);
1717 gmx_fio_ndo_rvec(fio
,molb
->posres_xA
,molb
->nposres_xA
);
1719 gmx_fio_do_int(fio
,molb
->nposres_xB
);
1720 if (molb
->nposres_xB
> 0) {
1722 snew(molb
->posres_xB
,molb
->nposres_xB
);
1724 gmx_fio_ndo_rvec(fio
,molb
->posres_xB
,molb
->nposres_xB
);
1729 static t_block
mtop_mols(gmx_mtop_t
*mtop
)
1735 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
1736 mols
.nr
+= mtop
->molblock
[mb
].nmol
;
1738 mols
.nalloc_index
= mols
.nr
+ 1;
1739 snew(mols
.index
,mols
.nalloc_index
);
1744 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
1745 for(mol
=0; mol
<mtop
->molblock
[mb
].nmol
; mol
++) {
1746 a
+= mtop
->molblock
[mb
].natoms_mol
;
1755 static void add_posres_molblock(gmx_mtop_t
*mtop
)
1760 gmx_molblock_t
*molb
;
1763 il
= &mtop
->moltype
[0].ilist
[F_POSRES
];
1769 for(i
=0; i
<il
->nr
; i
+=2) {
1770 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
1771 am
= max(am
,il
->iatoms
[i
+1]);
1772 if (ip
->posres
.pos0B
[XX
] != ip
->posres
.pos0A
[XX
] ||
1773 ip
->posres
.pos0B
[YY
] != ip
->posres
.pos0A
[YY
] ||
1774 ip
->posres
.pos0B
[ZZ
] != ip
->posres
.pos0A
[ZZ
]) {
1778 /* Make the posres coordinate block end at a molecule end */
1780 while(am
>= mtop
->mols
.index
[mol
+1]) {
1783 molb
= &mtop
->molblock
[0];
1784 molb
->nposres_xA
= mtop
->mols
.index
[mol
+1];
1785 snew(molb
->posres_xA
,molb
->nposres_xA
);
1787 molb
->nposres_xB
= molb
->nposres_xA
;
1788 snew(molb
->posres_xB
,molb
->nposres_xB
);
1790 molb
->nposres_xB
= 0;
1792 for(i
=0; i
<il
->nr
; i
+=2) {
1793 ip
= &mtop
->ffparams
.iparams
[il
->iatoms
[i
]];
1794 a
= il
->iatoms
[i
+1];
1795 molb
->posres_xA
[a
][XX
] = ip
->posres
.pos0A
[XX
];
1796 molb
->posres_xA
[a
][YY
] = ip
->posres
.pos0A
[YY
];
1797 molb
->posres_xA
[a
][ZZ
] = ip
->posres
.pos0A
[ZZ
];
1799 molb
->posres_xB
[a
][XX
] = ip
->posres
.pos0B
[XX
];
1800 molb
->posres_xB
[a
][YY
] = ip
->posres
.pos0B
[YY
];
1801 molb
->posres_xB
[a
][ZZ
] = ip
->posres
.pos0B
[ZZ
];
1806 static void set_disres_npair(gmx_mtop_t
*mtop
)
1813 ip
= mtop
->ffparams
.iparams
;
1815 for(mt
=0; mt
<mtop
->nmoltype
; mt
++) {
1816 il
= &mtop
->moltype
[mt
].ilist
[F_DISRES
];
1820 for(i
=0; i
<il
->nr
; i
+=3) {
1822 if (i
+3 == il
->nr
|| ip
[a
[i
]].disres
.label
!= ip
[a
[i
+3]].disres
.label
) {
1823 ip
[a
[i
]].disres
.npair
= npair
;
1831 static void do_mtop(t_fileio
*fio
, gmx_mtop_t
*mtop
,gmx_bool bRead
,
1839 do_symtab(fio
, &(mtop
->symtab
),bRead
);
1841 pr_symtab(debug
,0,"symtab",&mtop
->symtab
);
1843 do_symstr(fio
, &(mtop
->name
),bRead
,&(mtop
->symtab
));
1845 if (file_version
>= 57) {
1846 do_ffparams(fio
, &mtop
->ffparams
,bRead
,file_version
);
1848 gmx_fio_do_int(fio
,mtop
->nmoltype
);
1853 snew(mtop
->moltype
,mtop
->nmoltype
);
1854 if (file_version
< 57) {
1855 mtop
->moltype
[0].name
= mtop
->name
;
1858 for(mt
=0; mt
<mtop
->nmoltype
; mt
++) {
1859 do_moltype(fio
, &mtop
->moltype
[mt
],bRead
,&mtop
->symtab
,file_version
,
1863 if (file_version
>= 57) {
1864 gmx_fio_do_int(fio
,mtop
->nmolblock
);
1866 mtop
->nmolblock
= 1;
1869 snew(mtop
->molblock
,mtop
->nmolblock
);
1871 if (file_version
>= 57) {
1872 for(mb
=0; mb
<mtop
->nmolblock
; mb
++) {
1873 do_molblock(fio
, &mtop
->molblock
[mb
],bRead
,file_version
);
1875 gmx_fio_do_int(fio
,mtop
->natoms
);
1877 mtop
->molblock
[0].type
= 0;
1878 mtop
->molblock
[0].nmol
= 1;
1879 mtop
->molblock
[0].natoms_mol
= mtop
->moltype
[0].atoms
.nr
;
1880 mtop
->molblock
[0].nposres_xA
= 0;
1881 mtop
->molblock
[0].nposres_xB
= 0;
1884 do_atomtypes (fio
, &(mtop
->atomtypes
),bRead
,&(mtop
->symtab
), file_version
);
1886 pr_atomtypes(debug
,0,"atomtypes",&mtop
->atomtypes
,TRUE
);
1888 if (file_version
< 57) {
1889 /* Debug statements are inside do_idef */
1890 do_idef (fio
, &mtop
->ffparams
,&mtop
->moltype
[0],bRead
,file_version
);
1891 mtop
->natoms
= mtop
->moltype
[0].atoms
.nr
;
1894 if(file_version
>= 65)
1896 do_cmap(fio
, &mtop
->ffparams
.cmap_grid
,bRead
);
1900 mtop
->ffparams
.cmap_grid
.ngrid
= 0;
1901 mtop
->ffparams
.cmap_grid
.grid_spacing
= 0.1;
1902 mtop
->ffparams
.cmap_grid
.cmapdata
= NULL
;
1905 if (file_version
>= 57) {
1906 do_groups(fio
, &mtop
->groups
,bRead
,&(mtop
->symtab
),file_version
);
1909 if (file_version
< 57) {
1910 do_block(fio
, &mtop
->moltype
[0].cgs
,bRead
,file_version
);
1911 if (bRead
&& gmx_debug_at
) {
1912 pr_block(debug
,0,"cgs",&mtop
->moltype
[0].cgs
,TRUE
);
1914 do_block(fio
, &mtop
->mols
,bRead
,file_version
);
1915 /* Add the posres coordinates to the molblock */
1916 add_posres_molblock(mtop
);
1919 if (file_version
>= 57) {
1920 mtop
->mols
= mtop_mols(mtop
);
1923 pr_block(debug
,0,"mols",&mtop
->mols
,TRUE
);
1927 if (file_version
< 51) {
1928 /* Here used to be the shake blocks */
1929 do_blocka(fio
, &dumb
,bRead
,file_version
);
1937 close_symtab(&(mtop
->symtab
));
1941 /* If TopOnlyOK is TRUE then we can read even future versions
1942 * of tpx files, provided the file_generation hasn't changed.
1943 * If it is FALSE, we need the inputrecord too, and bail out
1944 * if the file is newer than the program.
1946 * The version and generation if the topology (see top of this file)
1947 * are returned in the two last arguments.
1949 * If possible, we will read the inputrec even when TopOnlyOK is TRUE.
1951 static void do_tpxheader(t_fileio
*fio
,gmx_bool bRead
,t_tpxheader
*tpx
,
1952 gmx_bool TopOnlyOK
, int *file_version
,
1953 int *file_generation
)
1962 gmx_fio_checktype(fio
);
1963 gmx_fio_setdebug(fio
,bDebugMode());
1965 /* NEW! XDR tpb file */
1966 precision
= sizeof(real
);
1968 gmx_fio_do_string(fio
,buf
);
1969 if (strncmp(buf
,"VERSION",7))
1970 gmx_fatal(FARGS
,"Can not read file %s,\n"
1971 " this file is from a Gromacs version which is older than 2.0\n"
1972 " Make a new one with grompp or use a gro or pdb file, if possible",
1973 gmx_fio_getname(fio
));
1974 gmx_fio_do_int(fio
,precision
);
1975 bDouble
= (precision
== sizeof(double));
1976 if ((precision
!= sizeof(float)) && !bDouble
)
1977 gmx_fatal(FARGS
,"Unknown precision in file %s: real is %d bytes "
1978 "instead of %d or %d",
1979 gmx_fio_getname(fio
),precision
,sizeof(float),sizeof(double));
1980 gmx_fio_setprecision(fio
,bDouble
);
1981 fprintf(stderr
,"Reading file %s, %s (%s precision)\n",
1982 gmx_fio_getname(fio
),buf
,bDouble
? "double" : "single");
1985 gmx_fio_write_string(fio
,GromacsVersion());
1986 bDouble
= (precision
== sizeof(double));
1987 gmx_fio_setprecision(fio
,bDouble
);
1988 gmx_fio_do_int(fio
,precision
);
1990 fgen
= tpx_generation
;
1993 /* Check versions! */
1994 gmx_fio_do_int(fio
,fver
);
1997 gmx_fio_do_int(fio
,fgen
);
2001 if(file_version
!=NULL
)
2002 *file_version
= fver
;
2003 if(file_generation
!=NULL
)
2004 *file_generation
= fgen
;
2007 if ((fver
<= tpx_incompatible_version
) ||
2008 ((fver
> tpx_version
) && !TopOnlyOK
) ||
2009 (fgen
> tpx_generation
))
2010 gmx_fatal(FARGS
,"reading tpx file (%s) version %d with version %d program",
2011 gmx_fio_getname(fio
),fver
,tpx_version
);
2013 do_section(fio
,eitemHEADER
,bRead
);
2014 gmx_fio_do_int(fio
,tpx
->natoms
);
2016 gmx_fio_do_int(fio
,tpx
->ngtc
);
2020 gmx_fio_do_int(fio
,idum
);
2021 gmx_fio_do_real(fio
,rdum
);
2023 gmx_fio_do_real(fio
,tpx
->lambda
);
2024 gmx_fio_do_int(fio
,tpx
->bIr
);
2025 gmx_fio_do_int(fio
,tpx
->bTop
);
2026 gmx_fio_do_int(fio
,tpx
->bX
);
2027 gmx_fio_do_int(fio
,tpx
->bV
);
2028 gmx_fio_do_int(fio
,tpx
->bF
);
2029 gmx_fio_do_int(fio
,tpx
->bBox
);
2031 if((fgen
> tpx_generation
)) {
2032 /* This can only happen if TopOnlyOK=TRUE */
2037 static int do_tpx(t_fileio
*fio
, gmx_bool bRead
,
2038 t_inputrec
*ir
,t_state
*state
,rvec
*f
,gmx_mtop_t
*mtop
,
2039 gmx_bool bXVallocated
)
2044 gmx_bool TopOnlyOK
,bDum
=TRUE
;
2045 int file_version
,file_generation
;
2049 gmx_bool bPeriodicMols
;
2052 tpx
.natoms
= state
->natoms
;
2053 tpx
.ngtc
= state
->ngtc
;
2054 tpx
.lambda
= state
->lambda
;
2055 tpx
.bIr
= (ir
!= NULL
);
2056 tpx
.bTop
= (mtop
!= NULL
);
2057 tpx
.bX
= (state
->x
!= NULL
);
2058 tpx
.bV
= (state
->v
!= NULL
);
2059 tpx
.bF
= (f
!= NULL
);
2063 TopOnlyOK
= (ir
==NULL
);
2065 do_tpxheader(fio
,bRead
,&tpx
,TopOnlyOK
,&file_version
,&file_generation
);
2069 state
->lambda
= tpx
.lambda
;
2070 /* The init_state calls initialize the Nose-Hoover xi integrals to zero */
2074 init_state(state
,0,tpx
.ngtc
,0,0); /* nose-hoover chains */ /* eventually, need to add nnhpres here? */
2075 state
->natoms
= tpx
.natoms
;
2076 state
->nalloc
= tpx
.natoms
;
2080 init_state(state
,tpx
.natoms
,tpx
.ngtc
,0,0); /* nose-hoover chains */
2084 #define do_test(fio,b,p) if (bRead && (p!=NULL) && !b) gmx_fatal(FARGS,"No %s in %s",#p,gmx_fio_getname(fio))
2086 do_test(fio
,tpx
.bBox
,state
->box
);
2087 do_section(fio
,eitemBOX
,bRead
);
2089 gmx_fio_ndo_rvec(fio
,state
->box
,DIM
);
2090 if (file_version
>= 51) {
2091 gmx_fio_ndo_rvec(fio
,state
->box_rel
,DIM
);
2093 /* We initialize box_rel after reading the inputrec */
2094 clear_mat(state
->box_rel
);
2096 if (file_version
>= 28) {
2097 gmx_fio_ndo_rvec(fio
,state
->boxv
,DIM
);
2098 if (file_version
< 56) {
2100 gmx_fio_ndo_rvec(fio
,mdum
,DIM
);
2105 if (state
->ngtc
> 0 && file_version
>= 28) {
2107 /*ndo_double(state->nosehoover_xi,state->ngtc,bDum);*/
2108 /*ndo_double(state->nosehoover_vxi,state->ngtc,bDum);*/
2109 /*ndo_double(state->therm_integral,state->ngtc,bDum);*/
2110 snew(dumv
,state
->ngtc
);
2111 if (file_version
< 69) {
2112 bDum
=gmx_fio_ndo_real(fio
,dumv
,state
->ngtc
);
2114 /* These used to be the Berendsen tcoupl_lambda's */
2115 bDum
=gmx_fio_ndo_real(fio
,dumv
,state
->ngtc
);
2119 /* Prior to tpx version 26, the inputrec was here.
2120 * I moved it to enable partial forward-compatibility
2121 * for analysis/viewer programs.
2123 if(file_version
<26) {
2124 do_test(fio
,tpx
.bIr
,ir
);
2125 do_section(fio
,eitemIR
,bRead
);
2128 do_inputrec(fio
, ir
,bRead
,file_version
,
2129 mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
2131 pr_inputrec(debug
,0,"inputrec",ir
,FALSE
);
2134 do_inputrec(fio
, &dum_ir
,bRead
,file_version
,
2135 mtop
? &mtop
->ffparams
.fudgeQQ
:NULL
);
2137 pr_inputrec(debug
,0,"inputrec",&dum_ir
,FALSE
);
2138 done_inputrec(&dum_ir
);
2144 do_test(fio
,tpx
.bTop
,mtop
);
2145 do_section(fio
,eitemTOP
,bRead
);
2148 do_mtop(fio
,mtop
,bRead
, file_version
);
2150 do_mtop(fio
,&dum_top
,bRead
,file_version
);
2151 done_mtop(&dum_top
,TRUE
);
2154 do_test(fio
,tpx
.bX
,state
->x
);
2155 do_section(fio
,eitemX
,bRead
);
2158 state
->flags
|= (1<<estX
);
2160 gmx_fio_ndo_rvec(fio
,state
->x
,state
->natoms
);
2163 do_test(fio
,tpx
.bV
,state
->v
);
2164 do_section(fio
,eitemV
,bRead
);
2167 state
->flags
|= (1<<estV
);
2169 gmx_fio_ndo_rvec(fio
,state
->v
,state
->natoms
);
2172 do_test(fio
,tpx
.bF
,f
);
2173 do_section(fio
,eitemF
,bRead
);
2174 if (tpx
.bF
) gmx_fio_ndo_rvec(fio
,f
,state
->natoms
);
2176 /* Starting with tpx version 26, we have the inputrec
2177 * at the end of the file, so we can ignore it
2178 * if the file is never than the software (but still the
2179 * same generation - see comments at the top of this file.
2184 bPeriodicMols
= FALSE
;
2185 if (file_version
>= 26) {
2186 do_test(fio
,tpx
.bIr
,ir
);
2187 do_section(fio
,eitemIR
,bRead
);
2189 if (file_version
>= 53) {
2190 /* Removed the pbc info from do_inputrec, since we always want it */
2193 bPeriodicMols
= ir
->bPeriodicMols
;
2195 gmx_fio_do_int(fio
,ePBC
);
2196 gmx_fio_do_gmx_bool(fio
,bPeriodicMols
);
2198 if (file_generation
<= tpx_generation
&& ir
) {
2199 do_inputrec(fio
, ir
,bRead
,file_version
,mtop
? &mtop
->ffparams
.fudgeQQ
: NULL
);
2201 pr_inputrec(debug
,0,"inputrec",ir
,FALSE
);
2202 if (file_version
< 51)
2203 set_box_rel(ir
,state
);
2204 if (file_version
< 53) {
2206 bPeriodicMols
= ir
->bPeriodicMols
;
2209 if (bRead
&& ir
&& file_version
>= 53) {
2210 /* We need to do this after do_inputrec, since that initializes ir */
2212 ir
->bPeriodicMols
= bPeriodicMols
;
2221 if (state
->ngtc
== 0)
2223 /* Reading old version without tcoupl state data: set it */
2224 init_gtc_state(state
,ir
->opts
.ngtc
,0,ir
->opts
.nhchainlength
);
2226 if (tpx
.bTop
&& mtop
)
2228 if (file_version
< 57)
2230 if (mtop
->moltype
[0].ilist
[F_DISRES
].nr
> 0)
2232 ir
->eDisre
= edrSimple
;
2236 ir
->eDisre
= edrNone
;
2239 set_disres_npair(mtop
);
2243 if (tpx
.bTop
&& mtop
)
2245 gmx_mtop_finalize(mtop
);
2248 if (file_version
>= 57)
2252 env
= getenv("GMX_NOCHARGEGROUPS");
2255 sscanf(env
,"%d",&ienv
);
2256 fprintf(stderr
,"\nFound env.var. GMX_NOCHARGEGROUPS = %d\n",
2261 "Will make single atomic charge groups in non-solvent%s\n",
2262 ienv
> 1 ? " and solvent" : "");
2263 gmx_mtop_make_atomic_charge_groups(mtop
,ienv
==1);
2265 fprintf(stderr
,"\n");
2273 /************************************************************
2275 * The following routines are the exported ones
2277 ************************************************************/
2279 t_fileio
*open_tpx(const char *fn
,const char *mode
)
2281 return gmx_fio_open(fn
,mode
);
2284 void close_tpx(t_fileio
*fio
)
2289 void read_tpxheader(const char *fn
, t_tpxheader
*tpx
, gmx_bool TopOnlyOK
,
2290 int *file_version
, int *file_generation
)
2294 fio
= open_tpx(fn
,"r");
2295 do_tpxheader(fio
,TRUE
,tpx
,TopOnlyOK
,file_version
,file_generation
);
2299 void write_tpx_state(const char *fn
,
2300 t_inputrec
*ir
,t_state
*state
,gmx_mtop_t
*mtop
)
2304 fio
= open_tpx(fn
,"w");
2305 do_tpx(fio
,FALSE
,ir
,state
,NULL
,mtop
,FALSE
);
2309 void read_tpx_state(const char *fn
,
2310 t_inputrec
*ir
,t_state
*state
,rvec
*f
,gmx_mtop_t
*mtop
)
2314 fio
= open_tpx(fn
,"r");
2315 do_tpx(fio
,TRUE
,ir
,state
,f
,mtop
,FALSE
);
2319 int read_tpx(const char *fn
,
2320 t_inputrec
*ir
, matrix box
,int *natoms
,
2321 rvec
*x
,rvec
*v
,rvec
*f
,gmx_mtop_t
*mtop
)
2329 fio
= open_tpx(fn
,"r");
2330 ePBC
= do_tpx(fio
,TRUE
,ir
,&state
,f
,mtop
,TRUE
);
2332 *natoms
= state
.natoms
;
2334 copy_mat(state
.box
,box
);
2342 int read_tpx_top(const char *fn
,
2343 t_inputrec
*ir
, matrix box
,int *natoms
,
2344 rvec
*x
,rvec
*v
,rvec
*f
,t_topology
*top
)
2350 ePBC
= read_tpx(fn
,ir
,box
,natoms
,x
,v
,f
,&mtop
);
2352 *top
= gmx_mtop_t_to_t_topology(&mtop
);
2357 gmx_bool
fn2bTPX(const char *file
)
2359 switch (fn2ftp(file
)) {
2369 gmx_bool
read_tps_conf(const char *infile
,char *title
,t_topology
*top
,int *ePBC
,
2370 rvec
**x
,rvec
**v
,matrix box
,gmx_bool bMass
)
2373 int natoms
,i
,version
,generation
;
2374 gmx_bool bTop
,bXNULL
;
2376 t_topology
*topconv
;
2379 bTop
= fn2bTPX(infile
);
2382 read_tpxheader(infile
,&header
,TRUE
,&version
,&generation
);
2384 snew(*x
,header
.natoms
);
2386 snew(*v
,header
.natoms
);
2388 *ePBC
= read_tpx(infile
,NULL
,box
,&natoms
,
2389 (x
==NULL
) ? NULL
: *x
,(v
==NULL
) ? NULL
: *v
,NULL
,mtop
);
2390 *top
= gmx_mtop_t_to_t_topology(mtop
);
2392 strcpy(title
,*top
->name
);
2393 tpx_make_chain_identifiers(&top
->atoms
,&top
->mols
);
2396 get_stx_coordnum(infile
,&natoms
);
2397 init_t_atoms(&top
->atoms
,natoms
,FALSE
);
2398 bXNULL
= (x
== NULL
);
2402 read_stx_conf(infile
,title
,&top
->atoms
,*x
,(v
==NULL
) ? NULL
: *v
,ePBC
,box
);
2408 aps
= gmx_atomprop_init();
2409 for(i
=0; (i
<natoms
); i
++)
2410 if (!gmx_atomprop_query(aps
,epropMass
,
2411 *top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].name
,
2412 *top
->atoms
.atomname
[i
],
2413 &(top
->atoms
.atom
[i
].m
))) {
2415 fprintf(debug
,"Can not find mass for atom %s %d %s, setting to 1\n",
2416 *top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].name
,
2417 top
->atoms
.resinfo
[top
->atoms
.atom
[i
].resind
].nr
,
2418 *top
->atoms
.atomname
[i
]);
2420 gmx_atomprop_destroy(aps
);
2422 top
->idef
.ntypes
=-1;