3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
35 /* This file is completely threadsafe - keep it that way! */
50 #include "gmx_fatal.h"
54 #include "mtop_util.h"
56 static void cmp_int(FILE *fp
,const char *s
,int index
,int i1
,int i2
)
60 fprintf(fp
,"%s[%d] (%d - %d)\n",s
,index
,i1
,i2
);
62 fprintf(fp
,"%s (%d - %d)\n",s
,i1
,i2
);
66 static void cmp_gmx_large_int(FILE *fp
,const char *s
,gmx_large_int_t i1
,gmx_large_int_t i2
)
70 fprintf(fp
,gmx_large_int_pfmt
,i1
);
72 fprintf(fp
,gmx_large_int_pfmt
,i2
);
77 static void cmp_us(FILE *fp
,const char *s
,int index
,unsigned short i1
,unsigned short i2
)
81 fprintf(fp
,"%s[%d] (%d - %d)\n",s
,index
,i1
,i2
);
83 fprintf(fp
,"%s (%d - %d)\n",s
,i1
,i2
);
87 static void cmp_uc(FILE *fp
,const char *s
,int index
,unsigned char i1
,unsigned char i2
)
91 fprintf(fp
,"%s[%d] (%d - %d)\n",s
,index
,i1
,i2
);
93 fprintf(fp
,"%s (%d - %d)\n",s
,i1
,i2
);
97 static bool cmp_bool(FILE *fp
, const char *s
, int index
, bool b1
, bool b2
)
111 fprintf(fp
,"%s[%d] (%s - %s)\n",s
,index
,
112 bool_names
[b1
],bool_names
[b2
]);
114 fprintf(fp
,"%s (%s - %s)\n",s
,
115 bool_names
[b1
],bool_names
[b2
]);
120 static void cmp_str(FILE *fp
, const char *s
, int index
,
121 const char *s1
, const char *s2
)
123 if (strcmp(s1
,s2
) != 0) {
125 fprintf(fp
,"%s[%d] (%s - %s)\n",s
,index
,s1
,s2
);
127 fprintf(fp
,"%s (%s - %s)\n",s
,s1
,s2
);
131 static bool equal_real(real i1
,real i2
,real ftol
,real abstol
)
133 return ( ( 2*fabs(i1
- i2
) <= (fabs(i1
) + fabs(i2
))*ftol
) || fabs(i1
-i2
)<=abstol
);
136 static bool equal_double(double i1
,double i2
,real ftol
,real abstol
)
138 return ( ( 2*fabs(i1
- i2
) <= (fabs(i1
) + fabs(i2
))*ftol
) || fabs(i1
-i2
)<=abstol
);
142 cmp_real(FILE *fp
,const char *s
,int index
,real i1
,real i2
,real ftol
,real abstol
)
144 if (!equal_real(i1
,i2
,ftol
,abstol
)) {
146 fprintf(fp
,"%s[%2d] (%e - %e)\n",s
,index
,i1
,i2
);
148 fprintf(fp
,"%s (%e - %e)\n",s
,i1
,i2
);
153 cmp_double(FILE *fp
,const char *s
,int index
,double i1
,double i2
,double ftol
,double abstol
)
155 if (!equal_double(i1
,i2
,ftol
,abstol
)) {
157 fprintf(fp
,"%s[%2d] (%16.9e - %16.9e)\n",s
,index
,i1
,i2
);
159 fprintf(fp
,"%s (%16.9e - %16.9e)\n",s
,i1
,i2
);
163 static void cmp_rvec(FILE *fp
,const char *s
,int index
,rvec i1
,rvec i2
,real ftol
,real abstol
)
165 if(!equal_real(i1
[XX
],i2
[XX
],ftol
,abstol
) ||
166 !equal_real(i1
[YY
],i2
[YY
],ftol
,abstol
) ||
167 !equal_real(i1
[ZZ
],i2
[ZZ
],ftol
,abstol
))
170 fprintf(fp
,"%s[%5d] (%12.5e %12.5e %12.5e) - (%12.5e %12.5e %12.5e)\n",
171 s
,index
,i1
[XX
],i1
[YY
],i1
[ZZ
],i2
[XX
],i2
[YY
],i2
[ZZ
]);
173 fprintf(fp
,"%s (%12.5e %12.5e %12.5e) - (%12.5e %12.5e %12.5e)\n",
174 s
,i1
[XX
],i1
[YY
],i1
[ZZ
],i2
[XX
],i2
[YY
],i2
[ZZ
]);
178 static void cmp_ivec(FILE *fp
,const char *s
,int index
,ivec i1
,ivec i2
)
180 if ((i1
[XX
] != i2
[XX
]) || (i1
[YY
] != i2
[YY
]) || (i1
[ZZ
] != i2
[ZZ
])) {
182 fprintf(fp
,"%s[%5d] (%8d,%8d,%8d - %8d,%8d,%8d)\n",s
,index
,
183 i1
[XX
],i1
[YY
],i1
[ZZ
],i2
[XX
],i2
[YY
],i2
[ZZ
]);
185 fprintf(fp
,"%s (%8d,%8d,%8d - %8d,%8d,%8d)\n",s
,
186 i1
[XX
],i1
[YY
],i1
[ZZ
],i2
[XX
],i2
[YY
],i2
[ZZ
]);
190 static void cmp_ilist(FILE *fp
,int ftype
,t_ilist
*il1
,t_ilist
*il2
)
195 fprintf(fp
,"comparing ilist %s\n",interaction_function
[ftype
].name
);
196 sprintf(buf
,"%s->nr",interaction_function
[ftype
].name
);
197 cmp_int(fp
,buf
,-1,il1
->nr
,il2
->nr
);
198 sprintf(buf
,"%s->iatoms",interaction_function
[ftype
].name
);
199 if (((il1
->nr
> 0) && (!il1
->iatoms
)) ||
200 ((il2
->nr
> 0) && (!il2
->iatoms
)) ||
201 ((il1
->nr
!= il2
->nr
)))
202 fprintf(fp
,"Comparing radically different topologies - %s is different\n",
205 for(i
=0; (i
<il1
->nr
); i
++)
206 cmp_int(fp
,buf
,i
,il1
->iatoms
[i
],il2
->iatoms
[i
]);
209 void cmp_iparm(FILE *fp
,const char *s
,t_functype ft
,
210 t_iparams ip1
,t_iparams ip2
,real ftol
,real abstol
)
216 for(i
=0; i
<MAXFORCEPARAM
&& !bDiff
; i
++)
217 bDiff
= !equal_real(ip1
.generic
.buf
[i
],ip2
.generic
.buf
[i
],ftol
,abstol
);
219 fprintf(fp
,"%s1: ",s
);
220 pr_iparams(fp
,ft
,&ip1
);
221 fprintf(fp
,"%s2: ",s
);
222 pr_iparams(fp
,ft
,&ip2
);
226 void cmp_iparm_AB(FILE *fp
,const char *s
,t_functype ft
,t_iparams ip1
,real ftol
,real abstol
)
228 int nrfpA
,nrfpB
,p0
,i
;
231 /* Normally the first parameter is perturbable */
233 nrfpA
= interaction_function
[ft
].nrfpA
;
234 nrfpB
= interaction_function
[ft
].nrfpB
;
237 } else if (interaction_function
[ft
].flags
& IF_TABULATED
) {
238 /* For tabulated interactions only the second parameter is perturbable */
243 for(i
=0; i
<nrfpB
&& !bDiff
; i
++) {
244 bDiff
= !equal_real(ip1
.generic
.buf
[p0
+i
],ip1
.generic
.buf
[nrfpA
+i
],ftol
,abstol
);
247 fprintf(fp
,"%s: ",s
);
248 pr_iparams(fp
,ft
,&ip1
);
252 static void cmp_idef(FILE *fp
,t_idef
*id1
,t_idef
*id2
,real ftol
,real abstol
)
255 char buf1
[64],buf2
[64];
257 fprintf(fp
,"comparing idef\n");
259 cmp_int(fp
,"idef->ntypes",-1,id1
->ntypes
,id2
->ntypes
);
260 cmp_int(fp
,"idef->atnr", -1,id1
->atnr
,id2
->atnr
);
261 for(i
=0; (i
<id1
->ntypes
); i
++) {
262 sprintf(buf1
,"idef->functype[%d]",i
);
263 sprintf(buf2
,"idef->iparam[%d]",i
);
264 cmp_int(fp
,buf1
,i
,(int)id1
->functype
[i
],(int)id2
->functype
[i
]);
265 cmp_iparm(fp
,buf2
,id1
->functype
[i
],
266 id1
->iparams
[i
],id2
->iparams
[i
],ftol
,abstol
);
268 cmp_real(fp
,"fudgeQQ",-1,id1
->fudgeQQ
,id2
->fudgeQQ
,ftol
,abstol
);
269 for(i
=0; (i
<F_NRE
); i
++)
270 cmp_ilist(fp
,i
,&(id1
->il
[i
]),&(id2
->il
[i
]));
272 for(i
=0; (i
<id1
->ntypes
); i
++)
273 cmp_iparm_AB(fp
,"idef->iparam",id1
->functype
[i
],id1
->iparams
[i
],ftol
,abstol
);
277 static void cmp_block(FILE *fp
,t_block
*b1
,t_block
*b2
,const char *s
)
282 fprintf(fp
,"comparing block %s\n",s
);
283 sprintf(buf
,"%s.nr",s
);
284 cmp_int(fp
,buf
,-1,b1
->nr
,b2
->nr
);
287 static void cmp_blocka(FILE *fp
,t_blocka
*b1
,t_blocka
*b2
,const char *s
)
292 fprintf(fp
,"comparing blocka %s\n",s
);
293 sprintf(buf
,"%s.nr",s
);
294 cmp_int(fp
,buf
,-1,b1
->nr
,b2
->nr
);
295 sprintf(buf
,"%s.nra",s
);
296 cmp_int(fp
,buf
,-1,b1
->nra
,b2
->nra
);
299 static void cmp_atom(FILE *fp
,int index
,t_atom
*a1
,t_atom
*a2
,real ftol
,real abstol
)
305 cmp_us(fp
,"atom.type",index
,a1
->type
,a2
->type
);
306 cmp_us(fp
,"atom.ptype",index
,a1
->ptype
,a2
->ptype
);
307 cmp_int(fp
,"atom.resind",index
,a1
->resind
,a2
->resind
);
308 cmp_int(fp
,"atom.atomnumber",index
,a1
->atomnumber
,a2
->atomnumber
);
309 cmp_real(fp
,"atom.m",index
,a1
->m
,a2
->m
,ftol
,abstol
);
310 cmp_real(fp
,"atom.q",index
,a1
->q
,a2
->q
,ftol
,abstol
);
311 cmp_us(fp
,"atom.typeB",index
,a1
->typeB
,a2
->typeB
);
312 cmp_real(fp
,"atom.mB",index
,a1
->mB
,a2
->mB
,ftol
,abstol
);
313 cmp_real(fp
,"atom.qB",index
,a1
->qB
,a2
->qB
,ftol
,abstol
);
315 cmp_us(fp
,"atom.type",index
,a1
->type
,a1
->typeB
);
316 cmp_real(fp
,"atom.m",index
,a1
->m
,a1
->mB
,ftol
,abstol
);
317 cmp_real(fp
,"atom.q",index
,a1
->q
,a1
->qB
,ftol
,abstol
);
321 static void cmp_atoms(FILE *fp
,t_atoms
*a1
,t_atoms
*a2
,real ftol
, real abstol
)
325 fprintf(fp
,"comparing atoms\n");
328 cmp_int(fp
,"atoms->nr",-1,a1
->nr
,a2
->nr
);
329 for(i
=0; (i
<a1
->nr
); i
++)
330 cmp_atom(fp
,i
,&(a1
->atom
[i
]),&(a2
->atom
[i
]),ftol
,abstol
);
332 for(i
=0; (i
<a1
->nr
); i
++)
333 cmp_atom(fp
,i
,&(a1
->atom
[i
]),NULL
,ftol
,abstol
);
337 static void cmp_top(FILE *fp
,t_topology
*t1
,t_topology
*t2
,real ftol
, real abstol
)
341 fprintf(fp
,"comparing top\n");
343 cmp_idef(fp
,&(t1
->idef
),&(t2
->idef
),ftol
,abstol
);
344 cmp_atoms(fp
,&(t1
->atoms
),&(t2
->atoms
),ftol
,abstol
);
345 cmp_block(fp
,&t1
->cgs
,&t2
->cgs
,"cgs");
346 cmp_block(fp
,&t1
->mols
,&t2
->mols
,"mols");
347 cmp_blocka(fp
,&t1
->excls
,&t2
->excls
,"excls");
349 cmp_idef(fp
,&(t1
->idef
),NULL
,ftol
,abstol
);
350 cmp_atoms(fp
,&(t1
->atoms
),NULL
,ftol
,abstol
);
354 static void cmp_rvecs(FILE *fp
,const char *title
,int n
,rvec x1
[],rvec x2
[],
355 bool bRMSD
,real ftol
,real abstol
)
362 for(i
=0; (i
<n
); i
++) {
363 for(m
=0; m
<DIM
; m
++) {
364 d
= x1
[i
][m
] - x2
[i
][m
];
368 fprintf(fp
,"%s RMSD %f\n",title
,sqrt(ssd
/n
));
370 for(i
=0; (i
<n
); i
++) {
371 cmp_rvec(fp
,title
,i
,x1
[i
],x2
[i
],ftol
,abstol
);
376 static void cmp_grpopts(FILE *fp
,t_grpopts
*opt1
,t_grpopts
*opt2
,real ftol
, real abstol
)
379 char buf1
[256],buf2
[256];
381 cmp_int(fp
,"inputrec->grpopts.ngtc",-1, opt1
->ngtc
,opt2
->ngtc
);
382 cmp_int(fp
,"inputrec->grpopts.ngacc",-1, opt1
->ngacc
,opt2
->ngacc
);
383 cmp_int(fp
,"inputrec->grpopts.ngfrz",-1, opt1
->ngfrz
,opt2
->ngfrz
);
384 cmp_int(fp
,"inputrec->grpopts.ngener",-1,opt1
->ngener
,opt2
->ngener
);
385 for(i
=0; (i
<min(opt1
->ngtc
,opt2
->ngtc
)); i
++) {
386 cmp_real(fp
,"inputrec->grpopts.nrdf",i
,opt1
->nrdf
[i
],opt2
->nrdf
[i
],ftol
,abstol
);
387 cmp_real(fp
,"inputrec->grpopts.ref_t",i
,opt1
->ref_t
[i
],opt2
->ref_t
[i
],ftol
,abstol
);
388 cmp_real(fp
,"inputrec->grpopts.tau_t",i
,opt1
->tau_t
[i
],opt2
->tau_t
[i
],ftol
,abstol
);
389 cmp_int(fp
,"inputrec->grpopts.annealing",i
,opt1
->annealing
[i
],opt2
->annealing
[i
]);
390 cmp_int(fp
,"inputrec->grpopts.anneal_npoints",i
,
391 opt1
->anneal_npoints
[i
],opt2
->anneal_npoints
[i
]);
392 if(opt1
->anneal_npoints
[i
]==opt2
->anneal_npoints
[i
]) {
393 sprintf(buf1
,"inputrec->grpopts.anneal_time[%d]",i
);
394 sprintf(buf2
,"inputrec->grpopts.anneal_temp[%d]",i
);
395 for(j
=0;j
<opt1
->anneal_npoints
[i
];j
++) {
396 cmp_real(fp
,buf1
,j
,opt1
->anneal_time
[i
][j
],opt2
->anneal_time
[i
][j
],ftol
,abstol
);
397 cmp_real(fp
,buf2
,j
,opt1
->anneal_temp
[i
][j
],opt2
->anneal_temp
[i
][j
],ftol
,abstol
);
401 if (opt1
->ngener
== opt2
->ngener
) {
402 for(i
=0; i
<opt1
->ngener
; i
++)
403 for(j
=i
; j
<opt1
->ngener
; j
++) {
404 sprintf(buf1
,"inputrec->grpopts.egp_flags[%d]",i
);
406 opt1
->egp_flags
[opt1
->ngener
*i
+j
],
407 opt2
->egp_flags
[opt1
->ngener
*i
+j
]);
410 for(i
=0; (i
<min(opt1
->ngacc
,opt2
->ngacc
)); i
++)
411 cmp_rvec(fp
,"inputrec->grpopts.acc",i
,opt1
->acc
[i
],opt2
->acc
[i
],ftol
,abstol
);
412 for(i
=0; (i
<min(opt1
->ngfrz
,opt2
->ngfrz
)); i
++)
413 cmp_ivec(fp
,"inputrec->grpopts.nFreeze",i
,opt1
->nFreeze
[i
],opt2
->nFreeze
[i
]);
416 static void cmp_cosines(FILE *fp
,const char *s
,t_cosines c1
[DIM
],t_cosines c2
[DIM
],real ftol
, real abstol
)
421 for(m
=0; (m
<DIM
); m
++) {
422 sprintf(buf
,"inputrec->%s[%d]",s
,m
);
423 cmp_int(fp
,buf
,0,c1
->n
,c2
->n
);
424 for(i
=0; (i
<min(c1
->n
,c2
->n
)); i
++) {
425 cmp_real(fp
,buf
,i
,c1
->a
[i
],c2
->a
[i
],ftol
,abstol
);
426 cmp_real(fp
,buf
,i
,c1
->phi
[i
],c2
->phi
[i
],ftol
,abstol
);
431 static void cmp_pull(FILE *fp
,t_pull
*pull1
,t_pull
*pull2
,real ftol
, real abstol
)
433 fprintf(fp
,"WARNING: Both files use COM pulling, but comparing of the pull struct is not implemented (yet). The pull parameters could be the same or different.\n");
436 static void cmp_inputrec(FILE *fp
,t_inputrec
*ir1
,t_inputrec
*ir2
,real ftol
, real abstol
)
440 fprintf(fp
,"comparing inputrec\n");
442 /* gcc 2.96 doesnt like these defines at all, but issues a huge list
443 * of warnings. Maybe it will change in future versions, but for the
444 * moment I've spelled them out instead. /EL 000820
445 * #define CIB(s) cmp_int(fp,"inputrec->"#s,0,ir1->##s,ir2->##s)
446 * #define CII(s) cmp_int(fp,"inputrec->"#s,0,ir1->##s,ir2->##s)
447 * #define CIR(s) cmp_real(fp,"inputrec->"#s,0,ir1->##s,ir2->##s,ftol)
449 cmp_int(fp
,"inputrec->eI",-1,ir1
->eI
,ir2
->eI
);
450 cmp_gmx_large_int(fp
,"inputrec->nsteps",ir1
->nsteps
,ir2
->nsteps
);
451 cmp_gmx_large_int(fp
,"inputrec->init_step",ir1
->init_step
,ir2
->init_step
);
452 cmp_int(fp
,"inputrec->simulation_part",-1,ir1
->simulation_part
,ir2
->simulation_part
);
453 cmp_int(fp
,"inputrec->nstcalcenergy",-1,ir1
->nstcalcenergy
,ir2
->nstcalcenergy
);
454 cmp_int(fp
,"inputrec->ePBC",-1,ir1
->ePBC
,ir2
->ePBC
);
455 cmp_int(fp
,"inputrec->bPeriodicMols",-1,ir1
->bPeriodicMols
,ir2
->bPeriodicMols
);
456 cmp_int(fp
,"inputrec->ns_type",-1,ir1
->ns_type
,ir2
->ns_type
);
457 cmp_int(fp
,"inputrec->nstlist",-1,ir1
->nstlist
,ir2
->nstlist
);
458 cmp_int(fp
,"inputrec->ndelta",-1,ir1
->ndelta
,ir2
->ndelta
);
459 cmp_int(fp
,"inputrec->nstcomm",-1,ir1
->nstcomm
,ir2
->nstcomm
);
460 cmp_int(fp
,"inputrec->comm_mode",-1,ir1
->comm_mode
,ir2
->comm_mode
);
461 cmp_int(fp
,"inputrec->nstcheckpoint",-1,ir1
->nstcheckpoint
,ir2
->nstcheckpoint
);
462 cmp_int(fp
,"inputrec->nstlog",-1,ir1
->nstlog
,ir2
->nstlog
);
463 cmp_int(fp
,"inputrec->nstxout",-1,ir1
->nstxout
,ir2
->nstxout
);
464 cmp_int(fp
,"inputrec->nstvout",-1,ir1
->nstvout
,ir2
->nstvout
);
465 cmp_int(fp
,"inputrec->nstfout",-1,ir1
->nstfout
,ir2
->nstfout
);
466 cmp_int(fp
,"inputrec->nstenergy",-1,ir1
->nstenergy
,ir2
->nstenergy
);
467 cmp_int(fp
,"inputrec->nstxtcout",-1,ir1
->nstxtcout
,ir2
->nstxtcout
);
468 cmp_double(fp
,"inputrec->init_t",-1,ir1
->init_t
,ir2
->init_t
,ftol
,abstol
);
469 cmp_double(fp
,"inputrec->delta_t",-1,ir1
->delta_t
,ir2
->delta_t
,ftol
,abstol
);
470 cmp_real(fp
,"inputrec->xtcprec",-1,ir1
->xtcprec
,ir2
->xtcprec
,ftol
,abstol
);
471 cmp_int(fp
,"inputrec->nkx",-1,ir1
->nkx
,ir2
->nkx
);
472 cmp_int(fp
,"inputrec->nky",-1,ir1
->nky
,ir2
->nky
);
473 cmp_int(fp
,"inputrec->nkz",-1,ir1
->nkz
,ir2
->nkz
);
474 cmp_int(fp
,"inputrec->pme_order",-1,ir1
->pme_order
,ir2
->pme_order
);
475 cmp_real(fp
,"inputrec->ewald_rtol",-1,ir1
->ewald_rtol
,ir2
->ewald_rtol
,ftol
,abstol
);
476 cmp_int(fp
,"inputrec->ewald_geometry",-1,ir1
->ewald_geometry
,ir2
->ewald_geometry
);
477 cmp_real(fp
,"inputrec->epsilon_surface",-1,ir1
->epsilon_surface
,ir2
->epsilon_surface
,ftol
,abstol
);
478 cmp_int(fp
,"inputrec->bOptFFT",-1,ir1
->bOptFFT
,ir2
->bOptFFT
);
479 cmp_int(fp
,"inputrec->bContinuation",-1,ir1
->bContinuation
,ir2
->bContinuation
);
480 cmp_int(fp
,"inputrec->bShakeSOR",-1,ir1
->bShakeSOR
,ir2
->bShakeSOR
);
481 cmp_int(fp
,"inputrec->etc",-1,ir1
->etc
,ir2
->etc
);
482 cmp_int(fp
,"inputrec->epc",-1,ir1
->epc
,ir2
->epc
);
483 cmp_int(fp
,"inputrec->epct",-1,ir1
->epct
,ir2
->epct
);
484 cmp_real(fp
,"inputrec->tau_p",-1,ir1
->tau_p
,ir2
->tau_p
,ftol
,abstol
);
485 cmp_rvec(fp
,"inputrec->ref_p(x)",-1,ir1
->ref_p
[XX
],ir2
->ref_p
[XX
],ftol
,abstol
);
486 cmp_rvec(fp
,"inputrec->ref_p(y)",-1,ir1
->ref_p
[YY
],ir2
->ref_p
[YY
],ftol
,abstol
);
487 cmp_rvec(fp
,"inputrec->ref_p(z)",-1,ir1
->ref_p
[ZZ
],ir2
->ref_p
[ZZ
],ftol
,abstol
);
488 cmp_rvec(fp
,"inputrec->compress(x)",-1,ir1
->compress
[XX
],ir2
->compress
[XX
],ftol
,abstol
);
489 cmp_rvec(fp
,"inputrec->compress(y)",-1,ir1
->compress
[YY
],ir2
->compress
[YY
],ftol
,abstol
);
490 cmp_rvec(fp
,"inputrec->compress(z)",-1,ir1
->compress
[ZZ
],ir2
->compress
[ZZ
],ftol
,abstol
);
491 cmp_int(fp
,"refcoord_scaling",-1,ir1
->refcoord_scaling
,ir2
->refcoord_scaling
);
492 cmp_rvec(fp
,"inputrec->posres_com",-1,ir1
->posres_com
,ir2
->posres_com
,ftol
,abstol
);
493 cmp_rvec(fp
,"inputrec->posres_comB",-1,ir1
->posres_comB
,ir2
->posres_comB
,ftol
,abstol
);
494 cmp_int(fp
,"inputrec->andersen_seed",-1,ir1
->andersen_seed
,ir2
->andersen_seed
);
495 cmp_real(fp
,"inputrec->rlist",-1,ir1
->rlist
,ir2
->rlist
,ftol
,abstol
);
496 cmp_real(fp
,"inputrec->rlistlong",-1,ir1
->rlistlong
,ir2
->rlistlong
,ftol
,abstol
);
497 cmp_real(fp
,"inputrec->rtpi",-1,ir1
->rtpi
,ir2
->rtpi
,ftol
,abstol
);
498 cmp_int(fp
,"inputrec->coulombtype",-1,ir1
->coulombtype
,ir2
->coulombtype
);
499 cmp_real(fp
,"inputrec->rcoulomb_switch",-1,ir1
->rcoulomb_switch
,ir2
->rcoulomb_switch
,ftol
,abstol
);
500 cmp_real(fp
,"inputrec->rcoulomb",-1,ir1
->rcoulomb
,ir2
->rcoulomb
,ftol
,abstol
);
501 cmp_int(fp
,"inputrec->vdwtype",-1,ir1
->vdwtype
,ir2
->vdwtype
);
502 cmp_real(fp
,"inputrec->rvdw_switch",-1,ir1
->rvdw_switch
,ir2
->rvdw_switch
,ftol
,abstol
);
503 cmp_real(fp
,"inputrec->rvdw",-1,ir1
->rvdw
,ir2
->rvdw
,ftol
,abstol
);
504 cmp_real(fp
,"inputrec->epsilon_r",-1,ir1
->epsilon_r
,ir2
->epsilon_r
,ftol
,abstol
);
505 cmp_real(fp
,"inputrec->epsilon_rf",-1,ir1
->epsilon_rf
,ir2
->epsilon_rf
,ftol
,abstol
);
506 cmp_real(fp
,"inputrec->tabext",-1,ir1
->tabext
,ir2
->tabext
,ftol
,abstol
);
507 cmp_int(fp
,"inputrec->implicit_solvent",-1,ir1
->implicit_solvent
,ir2
->implicit_solvent
);
508 cmp_int(fp
,"inputrec->gb_algorithm",-1,ir1
->gb_algorithm
,ir2
->gb_algorithm
);
509 cmp_int(fp
,"inputrec->nstgbradii",-1,ir1
->nstgbradii
,ir2
->nstgbradii
);
510 cmp_real(fp
,"inputrec->rgbradii",-1,ir1
->rgbradii
,ir2
->rgbradii
,ftol
,abstol
);
511 cmp_real(fp
,"inputrec->gb_saltconc",-1,ir1
->gb_saltconc
,ir2
->gb_saltconc
,ftol
,abstol
);
512 cmp_real(fp
,"inputrec->gb_epsilon_solvent",-1,ir1
->gb_epsilon_solvent
,ir2
->gb_epsilon_solvent
,ftol
,abstol
);
513 cmp_real(fp
,"inputrec->gb_obc_alpha",-1,ir1
->gb_obc_alpha
,ir2
->gb_obc_alpha
,ftol
,abstol
);
514 cmp_real(fp
,"inputrec->gb_obc_beta",-1,ir1
->gb_obc_beta
,ir2
->gb_obc_beta
,ftol
,abstol
);
515 cmp_real(fp
,"inputrec->gb_obc_gamma",-1,ir1
->gb_obc_gamma
,ir2
->gb_obc_gamma
,ftol
,abstol
);
516 cmp_real(fp
,"inputrec->gb_dielectric_offset",-1,ir1
->gb_dielectric_offset
,ir2
->gb_dielectric_offset
,ftol
,abstol
);
517 cmp_int(fp
,"inputrec->sa_algorithm",-1,ir1
->sa_algorithm
,ir2
->sa_algorithm
);
518 cmp_real(fp
,"inputrec->sa_surface_tension",-1,ir1
->sa_surface_tension
,ir2
->sa_surface_tension
,ftol
,abstol
);
520 cmp_int(fp
,"inputrec->eDispCorr",-1,ir1
->eDispCorr
,ir2
->eDispCorr
);
521 cmp_real(fp
,"inputrec->shake_tol",-1,ir1
->shake_tol
,ir2
->shake_tol
,ftol
,abstol
);
522 cmp_int(fp
,"inputrec->efep",-1,ir1
->efep
,ir2
->efep
);
523 cmp_double(fp
,"inputrec->init_lambda",-1,ir1
->init_lambda
,ir2
->init_lambda
,ftol
,abstol
);
524 cmp_double(fp
,"inputrec->delta_lambda",-1,ir1
->delta_lambda
,ir2
->delta_lambda
,ftol
,abstol
);
525 cmp_int(fp
,"inputrec->n_foreign_lambda",-1,ir1
->n_flambda
,ir2
->n_flambda
);
526 for(i
=0; i
<min(ir1
->n_flambda
,ir2
->n_flambda
); i
++) {
527 cmp_double(fp
,"inputrec->foreign_lambda",-1,ir1
->flambda
[i
],ir2
->flambda
[i
],ftol
,abstol
);
529 cmp_real(fp
,"inputrec->sc_alpha",-1,ir1
->sc_alpha
,ir2
->sc_alpha
,ftol
,abstol
);
530 cmp_int(fp
,"inputrec->sc_power",-1,ir1
->sc_power
,ir2
->sc_power
);
531 cmp_real(fp
,"inputrec->sc_sigma",-1,ir1
->sc_sigma
,ir2
->sc_sigma
,ftol
,abstol
);
532 cmp_int(fp
,"inputrec->nstdhdl",-1,ir1
->nstdhdl
,ir2
->nstdhdl
);
534 cmp_int(fp
,"inputrec->nwall",-1,ir1
->nwall
,ir2
->nwall
);
535 cmp_int(fp
,"inputrec->wall_type",-1,ir1
->wall_type
,ir2
->wall_type
);
536 cmp_int(fp
,"inputrec->wall_atomtype[0]",-1,ir1
->wall_atomtype
[0],ir2
->wall_atomtype
[0]);
537 cmp_int(fp
,"inputrec->wall_atomtype[1]",-1,ir1
->wall_atomtype
[1],ir2
->wall_atomtype
[1]);
538 cmp_real(fp
,"inputrec->wall_density[0]",-1,ir1
->wall_density
[0],ir2
->wall_density
[0],ftol
,abstol
);
539 cmp_real(fp
,"inputrec->wall_density[1]",-1,ir1
->wall_density
[1],ir2
->wall_density
[1],ftol
,abstol
);
540 cmp_real(fp
,"inputrec->wall_ewald_zfac",-1,ir1
->wall_ewald_zfac
,ir2
->wall_ewald_zfac
,ftol
,abstol
);
542 cmp_int(fp
,"inputrec->ePull",-1,ir1
->ePull
,ir2
->ePull
);
543 if (ir1
->ePull
== ir2
->ePull
&& ir1
->ePull
!= epullNO
)
544 cmp_pull(fp
,ir1
->pull
,ir2
->pull
,ftol
,abstol
);
546 cmp_int(fp
,"inputrec->eDisre",-1,ir1
->eDisre
,ir2
->eDisre
);
547 cmp_real(fp
,"inputrec->dr_fc",-1,ir1
->dr_fc
,ir2
->dr_fc
,ftol
,abstol
);
548 cmp_int(fp
,"inputrec->eDisreWeighting",-1,ir1
->eDisreWeighting
,ir2
->eDisreWeighting
);
549 cmp_int(fp
,"inputrec->bDisreMixed",-1,ir1
->bDisreMixed
,ir2
->bDisreMixed
);
550 cmp_int(fp
,"inputrec->nstdisreout",-1,ir1
->nstdisreout
,ir2
->nstdisreout
);
551 cmp_real(fp
,"inputrec->dr_tau",-1,ir1
->dr_tau
,ir2
->dr_tau
,ftol
,abstol
);
552 cmp_real(fp
,"inputrec->orires_fc",-1,ir1
->orires_fc
,ir2
->orires_fc
,ftol
,abstol
);
553 cmp_real(fp
,"inputrec->orires_tau",-1,ir1
->orires_tau
,ir2
->orires_tau
,ftol
,abstol
);
554 cmp_int(fp
,"inputrec->nstorireout",-1,ir1
->nstorireout
,ir2
->nstorireout
);
555 cmp_real(fp
,"inputrec->dihre_fc",-1,ir1
->dihre_fc
,ir2
->dihre_fc
,ftol
,abstol
);
556 cmp_real(fp
,"inputrec->em_stepsize",-1,ir1
->em_stepsize
,ir2
->em_stepsize
,ftol
,abstol
);
557 cmp_real(fp
,"inputrec->em_tol",-1,ir1
->em_tol
,ir2
->em_tol
,ftol
,abstol
);
558 cmp_int(fp
,"inputrec->niter",-1,ir1
->niter
,ir2
->niter
);
559 cmp_real(fp
,"inputrec->fc_stepsize",-1,ir1
->fc_stepsize
,ir2
->fc_stepsize
,ftol
,abstol
);
560 cmp_int(fp
,"inputrec->nstcgsteep",-1,ir1
->nstcgsteep
,ir2
->nstcgsteep
);
561 cmp_int(fp
,"inputrec->nbfgscorr",0,ir1
->nbfgscorr
,ir2
->nbfgscorr
);
562 cmp_int(fp
,"inputrec->eConstrAlg",-1,ir1
->eConstrAlg
,ir2
->eConstrAlg
);
563 cmp_int(fp
,"inputrec->nProjOrder",-1,ir1
->nProjOrder
,ir2
->nProjOrder
);
564 cmp_real(fp
,"inputrec->LincsWarnAngle",-1,ir1
->LincsWarnAngle
,ir2
->LincsWarnAngle
,ftol
,abstol
);
565 cmp_int(fp
,"inputrec->nLincsIter",-1,ir1
->nLincsIter
,ir2
->nLincsIter
);
566 cmp_real(fp
,"inputrec->bd_fric",-1,ir1
->bd_fric
,ir2
->bd_fric
,ftol
,abstol
);
567 cmp_int(fp
,"inputrec->ld_seed",-1,ir1
->ld_seed
,ir2
->ld_seed
);
568 cmp_real(fp
,"inputrec->cos_accel",-1,ir1
->cos_accel
,ir2
->cos_accel
,ftol
,abstol
);
569 cmp_rvec(fp
,"inputrec->deform(a)",-1,ir1
->deform
[XX
],ir2
->deform
[XX
],ftol
,abstol
);
570 cmp_rvec(fp
,"inputrec->deform(b)",-1,ir1
->deform
[YY
],ir2
->deform
[YY
],ftol
,abstol
);
571 cmp_rvec(fp
,"inputrec->deform(c)",-1,ir1
->deform
[ZZ
],ir2
->deform
[ZZ
],ftol
,abstol
);
573 cmp_int(fp
,"ir->adress_type" ,-1,ir1
->adress_type
,ir2
->adress_type
);
574 cmp_int(fp
,"ir->badress_new_wf" ,-1,ir1
->badress_new_wf
,ir2
->badress_new_wf
);
575 cmp_real(fp
,"ir->adress_const_wf" ,-1,ir1
->adress_const_wf
,ir2
->adress_const_wf
,ftol
,abstol
);
576 cmp_real(fp
,"ir->adress_ex_width" ,-1,ir1
->adress_ex_width
,ir2
->adress_ex_width
,ftol
,abstol
);
577 cmp_real(fp
,"ir->adress_hy_width" ,-1,ir1
->adress_hy_width
,ir2
->adress_hy_width
,ftol
,abstol
);
578 cmp_int(fp
,"ir->adress_icor" ,-1,ir1
->adress_icor
,ir2
->adress_icor
);
579 cmp_int(fp
,"ir->adress_ivdw" ,-1,ir1
->adress_ivdw
,ir2
->adress_ivdw
);
580 cmp_int(fp
,"ir->adress_site" ,-1,ir1
->adress_site
,ir2
->adress_site
);
581 cmp_rvec(fp
,"ir->adress_refs" ,-1,ir1
->adress_refs
,ir2
->adress_refs
,ftol
,abstol
);
583 cmp_int(fp
,"inputrec->userint1",-1,ir1
->userint1
,ir2
->userint1
);
584 cmp_int(fp
,"inputrec->userint2",-1,ir1
->userint2
,ir2
->userint2
);
585 cmp_int(fp
,"inputrec->userint3",-1,ir1
->userint3
,ir2
->userint3
);
586 cmp_int(fp
,"inputrec->userint4",-1,ir1
->userint4
,ir2
->userint4
);
587 cmp_real(fp
,"inputrec->userreal1",-1,ir1
->userreal1
,ir2
->userreal1
,ftol
,abstol
);
588 cmp_real(fp
,"inputrec->userreal2",-1,ir1
->userreal2
,ir2
->userreal2
,ftol
,abstol
);
589 cmp_real(fp
,"inputrec->userreal3",-1,ir1
->userreal3
,ir2
->userreal3
,ftol
,abstol
);
590 cmp_real(fp
,"inputrec->userreal4",-1,ir1
->userreal4
,ir2
->userreal4
,ftol
,abstol
);
591 cmp_grpopts(fp
,&(ir1
->opts
),&(ir2
->opts
),ftol
,abstol
);
592 cmp_cosines(fp
,"ex",ir1
->ex
,ir2
->ex
,ftol
,abstol
);
593 cmp_cosines(fp
,"et",ir1
->et
,ir2
->et
,ftol
,abstol
);
596 static void comp_pull_AB(FILE *fp
,t_pull
*pull
,real ftol
,real abstol
)
600 for(i
=0; i
<pull
->ngrp
+1; i
++) {
601 fprintf(fp
,"comparing pull group %d\n",i
);
602 cmp_real(fp
,"pullgrp->k",-1,pull
->grp
[i
].k
,pull
->grp
[i
].kB
,ftol
,abstol
);
606 static void comp_state(t_state
*st1
, t_state
*st2
,
607 bool bRMSD
,real ftol
,real abstol
)
611 fprintf(stdout
,"comparing flags\n");
612 cmp_int(stdout
,"flags",-1,st1
->flags
,st2
->flags
);
613 fprintf(stdout
,"comparing box\n");
614 cmp_rvecs(stdout
,"box",DIM
,st1
->box
,st2
->box
,FALSE
,ftol
,abstol
);
615 fprintf(stdout
,"comparing box_rel\n");
616 cmp_rvecs(stdout
,"box_rel",DIM
,st1
->box_rel
,st2
->box_rel
,FALSE
,ftol
,abstol
);
617 fprintf(stdout
,"comparing boxv\n");
618 cmp_rvecs(stdout
,"boxv",DIM
,st1
->boxv
,st2
->boxv
,FALSE
,ftol
,abstol
);
619 if (st1
->flags
& (1<<estPRES_PREV
)) {
620 fprintf(stdout
,"comparing prev_pres\n");
621 cmp_rvecs(stdout
,"pres_prev",DIM
,st1
->pres_prev
,st2
->pres_prev
,FALSE
,ftol
,abstol
);
623 cmp_int(stdout
,"ngtc",-1,st1
->ngtc
,st2
->ngtc
);
624 if (st1
->ngtc
== st2
->ngtc
) {
625 for(i
=0; i
<st1
->ngtc
; i
++)
626 cmp_real(stdout
,"nosehoover_xi",
627 i
,st1
->nosehoover_xi
[i
],st2
->nosehoover_xi
[i
],ftol
,abstol
);
629 cmp_int(stdout
,"natoms",-1,st1
->natoms
,st2
->natoms
);
630 if (st1
->natoms
== st2
->natoms
) {
631 if ((st1
->flags
& (1<<estX
)) && (st2
->flags
& (1<<estX
))) {
632 fprintf(stdout
,"comparing x\n");
633 cmp_rvecs(stdout
,"x",st1
->natoms
,st1
->x
,st2
->x
,bRMSD
,ftol
,abstol
);
635 if ((st1
->flags
& (1<<estV
)) && (st2
->flags
& (1<<estV
))) {
636 fprintf(stdout
,"comparing v\n");
637 cmp_rvecs(stdout
,"v",st1
->natoms
,st1
->v
,st2
->v
,bRMSD
,ftol
,abstol
);
642 void comp_tpx(const char *fn1
,const char *fn2
,
643 bool bRMSD
,real ftol
,real abstol
)
655 for(i
=0; i
<(fn2
? 2 : 1); i
++) {
656 read_tpx_state(ff
[i
],&(ir
[i
]),&state
[i
],NULL
,&(mtop
[i
]));
659 cmp_inputrec(stdout
,&ir
[0],&ir
[1],ftol
,abstol
);
660 /* Convert gmx_mtop_t to t_topology.
661 * We should implement direct mtop comparison,
662 * but it might be useful to keep t_topology comparison as an option.
664 top
[0] = gmx_mtop_t_to_t_topology(&mtop
[0]);
665 top
[1] = gmx_mtop_t_to_t_topology(&mtop
[1]);
666 cmp_top(stdout
,&top
[0],&top
[1],ftol
,abstol
);
667 comp_state(&state
[0],&state
[1],bRMSD
,ftol
,abstol
);
669 if (ir
[0].efep
== efepNO
) {
670 fprintf(stdout
,"inputrec->efep = %s\n",efep_names
[ir
[0].efep
]);
672 if (ir
[0].ePull
!= epullNO
) {
673 comp_pull_AB(stdout
,ir
->pull
,ftol
,abstol
);
675 /* Convert gmx_mtop_t to t_topology.
676 * We should implement direct mtop comparison,
677 * but it might be useful to keep t_topology comparison as an option.
679 top
[0] = gmx_mtop_t_to_t_topology(&mtop
[0]);
680 cmp_top(stdout
,&top
[0],NULL
,ftol
,abstol
);
685 void comp_frame(FILE *fp
, t_trxframe
*fr1
, t_trxframe
*fr2
,
686 bool bRMSD
, real ftol
,real abstol
)
689 cmp_int(fp
,"flags",-1,fr1
->flags
,fr2
->flags
);
690 cmp_int(fp
,"not_ok",-1,fr1
->not_ok
,fr2
->not_ok
);
691 cmp_int(fp
,"natoms",-1,fr1
->natoms
,fr2
->natoms
);
692 cmp_real(fp
,"t0",-1,fr1
->t0
,fr2
->t0
,ftol
,abstol
);
693 if (cmp_bool(fp
,"bTitle",-1,fr1
->bTitle
,fr2
->bTitle
))
694 cmp_str(fp
,"title", -1, fr1
->title
, fr2
->title
);
695 if (cmp_bool(fp
,"bStep",-1,fr1
->bStep
,fr2
->bStep
))
696 cmp_int(fp
,"step",-1,fr1
->step
,fr2
->step
);
697 cmp_int(fp
,"step",-1,fr1
->step
,fr2
->step
);
698 if (cmp_bool(fp
,"bTime",-1,fr1
->bTime
,fr2
->bTime
))
699 cmp_real(fp
,"time",-1,fr1
->time
,fr2
->time
,ftol
,abstol
);
700 if (cmp_bool(fp
,"bLambda",-1,fr1
->bLambda
,fr2
->bLambda
))
701 cmp_real(fp
,"lambda",-1,fr1
->lambda
,fr2
->lambda
,ftol
,abstol
);
702 if (cmp_bool(fp
,"bAtoms",-1,fr1
->bAtoms
,fr2
->bAtoms
))
703 cmp_atoms(fp
,fr1
->atoms
,fr2
->atoms
,ftol
,abstol
);
704 if (cmp_bool(fp
,"bPrec",-1,fr1
->bPrec
,fr2
->bPrec
))
705 cmp_real(fp
,"prec",-1,fr1
->prec
,fr2
->prec
,ftol
,abstol
);
706 if (cmp_bool(fp
,"bX",-1,fr1
->bX
,fr2
->bX
))
707 cmp_rvecs(fp
,"x",min(fr1
->natoms
,fr2
->natoms
),fr1
->x
,fr2
->x
,bRMSD
,ftol
,abstol
);
708 if (cmp_bool(fp
,"bV",-1,fr1
->bV
,fr2
->bV
))
709 cmp_rvecs(fp
,"v",min(fr1
->natoms
,fr2
->natoms
),fr1
->v
,fr2
->v
,bRMSD
,ftol
,abstol
);
710 if (cmp_bool(fp
,"bF",-1,fr1
->bF
,fr2
->bF
))
711 cmp_rvecs(fp
,"f",min(fr1
->natoms
,fr2
->natoms
),fr1
->f
,fr2
->f
,bRMSD
,ftol
,abstol
);
712 if (cmp_bool(fp
,"bBox",-1,fr1
->bBox
,fr2
->bBox
))
713 cmp_rvecs(fp
,"box",3,fr1
->box
,fr2
->box
,FALSE
,ftol
,abstol
);
716 void comp_trx(const output_env_t oenv
,const char *fn1
, const char *fn2
,
717 bool bRMSD
,real ftol
,real abstol
)
727 fprintf(stderr
,"Comparing trajectory files %s and %s\n",fn1
,fn2
);
729 b
[i
] = read_first_frame(oenv
,&status
[i
],fn
[i
],&fr
[i
],TRX_READ_X
|TRX_READ_V
|TRX_READ_F
);
733 comp_frame(stdout
, &(fr
[0]), &(fr
[1]), bRMSD
, ftol
, abstol
);
736 b
[i
] = read_next_frame(oenv
,status
[i
],&fr
[i
]);
737 } while (b
[0] && b
[1]);
739 for (i
=0; i
<2; i
++) {
741 fprintf(stdout
,"\nEnd of file on %s but not on %s\n",fn
[i
],fn
[1-i
]);
742 close_trj(status
[i
]);
746 fprintf(stdout
,"\nBoth files read correctly\n");
749 static void cmp_energies(FILE *fp
,int step1
,int step2
,int nre
,
750 t_energy e1
[],t_energy e2
[],
751 gmx_enxnm_t
*enm1
,gmx_enxnm_t
*enm2
,
752 real ftol
,real abstol
,
757 for(i
=0; (i
<maxener
); i
++) {
758 if (!equal_real(e1
[i
].e
,e2
[i
].e
,ftol
,abstol
))
759 fprintf(fp
,"%-15s step %3d: %12g, %s step %3d: %12g\n",
760 enm1
[i
].name
,step1
,e1
[i
].e
,
761 strcmp(enm1
[i
].name
,enm2
[i
].name
)!=0 ? enm2
[i
].name
:"",
766 static void cmp_disres(t_enxframe
*fr1
,t_enxframe
*fr2
,real ftol
, real abstol
)
769 char bav
[64],bt
[64],bs
[22];
771 cmp_int(stdout
,"ndisre",-1,fr1
->ndisre
,fr2
->ndisre
);
772 if ((fr1
->ndisre
== fr2
->ndisre
) && (fr1
->ndisre
> 0)) {
773 sprintf(bav
,"step %s: disre rav",gmx_step_str(fr1
->step
,bs
));
774 sprintf(bt
, "step %s: disre rt",gmx_step_str(fr1
->step
,bs
));
775 for(i
=0; (i
<fr1
->ndisre
); i
++) {
776 cmp_real(stdout
,bav
,i
,fr1
->disre_rm3tav
[i
],fr2
->disre_rm3tav
[i
],ftol
,abstol
);
777 cmp_real(stdout
,bt
,i
,fr1
->disre_rt
[i
] ,fr2
->disre_rt
[i
] ,ftol
,abstol
);
782 static void cmp_eblocks(t_enxframe
*fr1
,t_enxframe
*fr2
,real ftol
, real abstol
)
787 cmp_int(stdout
,"nblock",-1,fr1
->nblock
,fr2
->nblock
);
788 if ((fr1
->nblock
== fr2
->nblock
) && (fr1
->nblock
> 0)) {
789 for(j
=0; (j
<fr1
->nblock
); j
++) {
790 sprintf(buf
,"step %s: block[%d]",gmx_step_str(fr1
->step
,bs
),j
);
791 cmp_int(stdout
,buf
,-1,fr1
->nr
[j
],fr2
->nr
[j
]);
792 if ((fr1
->nr
[j
] == fr2
->nr
[j
]) && (fr1
->nr
[j
] > 0)) {
793 for(i
=0; (i
<fr1
->nr
[j
]); i
++) {
794 cmp_real(stdout
,buf
,i
,fr1
->block
[i
][j
],fr2
->block
[i
][j
],ftol
,abstol
);
801 void comp_enx(const char *fn1
,const char *fn2
,real ftol
,real abstol
,const char *lastener
)
803 int nre
,nre1
,nre2
,block
;
804 ener_file_t in1
, in2
;
807 gmx_enxnm_t
*enm1
=NULL
,*enm2
=NULL
;
808 t_enxframe
*fr1
,*fr2
;
811 fprintf(stdout
,"comparing energy file %s and %s\n\n",fn1
,fn2
);
813 in1
= open_enx(fn1
,"r");
814 in2
= open_enx(fn2
,"r");
815 do_enxnms(in1
,&nre1
,&enm1
);
816 do_enxnms(in2
,&nre2
,&enm2
);
818 fprintf(stdout
,"%s: nre=%d, %s: nre=%d\n",fn1
,nre1
,fn2
,nre2
);
822 fprintf(stdout
,"There are %d terms in the energy files\n\n",nre
);
825 for(i
=0; (i
<nre
); i
++) {
826 cmp_str(stdout
,"enm",i
,enm1
[i
].name
,enm2
[i
].name
);
827 cmp_str(stdout
,"unit",i
,enm1
[i
].unit
,enm2
[i
].unit
);
828 if ((lastener
!= NULL
) && (strstr(enm1
[i
].name
,lastener
) != NULL
)) {
833 fprintf(stdout
,"There are %d terms to compare in the energy files\n\n",
839 b1
= do_enx(in1
,fr1
);
840 b2
= do_enx(in2
,fr2
);
842 fprintf(stdout
,"\nEnd of file on %s but not on %s\n",fn2
,fn1
);
844 fprintf(stdout
,"\nEnd of file on %s but not on %s\n",fn1
,fn2
);
846 fprintf(stdout
,"\nFiles read succesfully\n");
848 cmp_real(stdout
,"t",-1,fr1
->t
,fr2
->t
,ftol
,abstol
);
849 cmp_int(stdout
,"step",-1,fr1
->step
,fr2
->step
);
850 cmp_int(stdout
,"nre",-1,fr1
->nre
,fr2
->nre
);
851 if ((fr1
->nre
== nre
) && (fr2
->nre
== nre
))
852 cmp_energies(stdout
,fr1
->step
,fr1
->step
,nre
,fr1
->ener
,fr2
->ener
,
853 enm1
,enm2
,ftol
,abstol
,maxener
);
854 cmp_disres(fr1
,fr2
,ftol
,abstol
);
855 cmp_eblocks(fr1
,fr2
,ftol
,abstol
);