src/gromacs/gmxlib/typedefs.c

   1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
   2  *
   3  *
   4  *                This source code is part of
   5  *
   6  *                 G   R   O   M   A   C   S
   7  *
   8  *          GROningen MAchine for Chemical Simulations
   9  *
  10  *                        VERSION 3.2.0
  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.
  15
  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.
  20  *
  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.
  27  *
  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.
  30  *
  31  * For more info, check our website at http://www.gromacs.org
  32  *
  33  * And Hey:
  34  * GROningen Mixture of Alchemy and Childrens' Stories
  35  */
  36 /* This file is completely threadsafe - keep it that way! */
  37 #ifdef HAVE_CONFIG_H
  38 #include <config.h>
  39 #endif
  40
  41 #include "smalloc.h"
  42 #include "symtab.h"
  43 #include "vec.h"
  44 #include "pbc.h"
  45 #include "macros.h"
  46 #include <string.h>
  47
  48 #ifdef GMX_THREAD_MPI
  49 #include "thread_mpi.h"
  50 #endif
  51
  52 /* The source code in this file should be thread-safe.
  53       Please keep it that way. */
  54
  55
  56
  57 static gmx_bool bOverAllocDD=FALSE;
  58 #ifdef GMX_THREAD_MPI
  59 static tMPI_Thread_mutex_t over_alloc_mutex=TMPI_THREAD_MUTEX_INITIALIZER;
  60 #endif
  61
  62
  63 void set_over_alloc_dd(gmx_bool set)
  64 {
  65 #ifdef GMX_THREAD_MPI
  66     tMPI_Thread_mutex_lock(&over_alloc_mutex);
  67     /* we just make sure that we don't set this at the same time.
  68        We don't worry too much about reading this rarely-set variable */
  69 #endif
  70     bOverAllocDD = set;
  71 #ifdef GMX_THREAD_MPI
  72     tMPI_Thread_mutex_unlock(&over_alloc_mutex);
  73 #endif
  74 }
  75
  76 int over_alloc_dd(int n)
  77 {
  78   if (bOverAllocDD)
  79     return OVER_ALLOC_FAC*n + 100;
  80   else
  81     return n;
  82 }
  83
  84 int gmx_large_int_to_int(gmx_large_int_t step,const char *warn)
  85 {
  86   int i;
  87
  88   i = (int)step;
  89
  90   if (warn != NULL && (step < INT_MIN || step > INT_MAX)) {
  91     fprintf(stderr,"\nWARNING during %s:\n",warn);
  92     fprintf(stderr,"step value ");
  93     fprintf(stderr,gmx_large_int_pfmt,step);
  94     fprintf(stderr," does not fit in int, converted to %d\n\n",i);
  95   }
  96
  97   return i;
  98 }
  99
 100 char *gmx_step_str(gmx_large_int_t i,char *buf)
 101 {
 102   sprintf(buf,gmx_large_int_pfmt,i);
 103
 104   return buf;
 105 }
 106
 107 void init_block(t_block *block)
 108 {
 109   int i;
 110
 111   block->nr           = 0;
 112   block->nalloc_index = 1;
 113   snew(block->index,block->nalloc_index);
 114   block->index[0]     = 0;
 115 }
 116
 117 void init_blocka(t_blocka *block)
 118 {
 119   int i;
 120
 121   block->nr           = 0;
 122   block->nra          = 0;
 123   block->nalloc_index = 1;
 124   snew(block->index,block->nalloc_index);
 125   block->index[0]     = 0;
 126   block->nalloc_a     = 0;
 127   block->a            = NULL;
 128 }
 129
 130 void init_atom(t_atoms *at)
 131 {
 132   int i;
 133
 134   at->nr       = 0;
 135   at->nres     = 0;
 136   at->atom     = NULL;
 137   at->resinfo  = NULL;
 138   at->atomname = NULL;
 139   at->atomtype = NULL;
 140   at->atomtypeB= NULL;
 141   at->pdbinfo  = NULL;
 142 }
 143
 144 void init_atomtypes(t_atomtypes *at)
 145 {
 146   at->nr = 0;
 147   at->radius = NULL;
 148   at->vol = NULL;
 149   at->atomnumber = NULL;
 150   at->gb_radius = NULL;
 151   at->S_hct = NULL;
 152 }
 153
 154 void init_groups(gmx_groups_t *groups)
 155 {
 156   int g;
 157
 158   groups->ngrpname = 0;
 159   groups->grpname  = NULL;
 160   for(g=0; (g<egcNR); g++) {
 161     groups->grps[g].nm_ind = NULL;
 162     groups->ngrpnr[g] = 0;
 163     groups->grpnr[g]  = NULL;
 164   }
 165
 166 }
 167
 168 void init_mtop(gmx_mtop_t *mtop)
 169 {
 170   mtop->name = NULL;
 171   mtop->nmoltype = 0;
 172   mtop->moltype = NULL;
 173   mtop->nmolblock = 0;
 174   mtop->molblock = NULL;
 175   mtop->maxres_renum = 0;
 176   mtop->maxresnr = -1;
 177   init_groups(&mtop->groups);
 178   init_block(&mtop->mols);
 179   open_symtab(&mtop->symtab);
 180 }
 181
 182 void init_top (t_topology *top)
 183 {
 184   int i;
 185
 186   top->name = NULL;
 187   init_atom (&(top->atoms));
 188   init_atomtypes(&(top->atomtypes));
 189   init_block(&top->cgs);
 190   init_block(&top->mols);
 191   init_blocka(&top->excls);
 192   open_symtab(&top->symtab);
 193 }
 194
 195 void init_inputrec(t_inputrec *ir)
 196 {
 197     memset(ir,0,(size_t)sizeof(*ir));
 198     snew(ir->fepvals,1);
 199     snew(ir->expandedvals,1);
 200     snew(ir->simtempvals,1);
 201 }
 202
 203 void stupid_fill_block(t_block *grp,int natom,gmx_bool bOneIndexGroup)
 204 {
 205   int i;
 206
 207   if (bOneIndexGroup) {
 208     grp->nalloc_index = 2;
 209     snew(grp->index,grp->nalloc_index);
 210     grp->index[0]=0;
 211     grp->index[1]=natom;
 212     grp->nr=1;
 213   }
 214   else {
 215     grp->nalloc_index = natom+1;
 216     snew(grp->index,grp->nalloc_index);
 217     snew(grp->index,natom+1);
 218     for(i=0; (i<=natom); i++)
 219       grp->index[i]=i;
 220     grp->nr=natom;
 221   }
 222 }
 223
 224 void stupid_fill_blocka(t_blocka *grp,int natom)
 225 {
 226   int i;
 227
 228   grp->nalloc_a = natom;
 229   snew(grp->a,grp->nalloc_a);
 230   for(i=0; (i<natom); i++)
 231     grp->a[i]=i;
 232   grp->nra=natom;
 233
 234   grp->nalloc_index = natom + 1;
 235   snew(grp->index,grp->nalloc_index);
 236   for(i=0; (i<=natom); i++)
 237     grp->index[i]=i;
 238   grp->nr=natom;
 239 }
 240
 241 void copy_blocka(const t_blocka *src,t_blocka *dest)
 242 {
 243   int i;
 244
 245   dest->nr = src->nr;
 246   dest->nalloc_index = dest->nr + 1;
 247   snew(dest->index,dest->nalloc_index);
 248   for(i=0; i<dest->nr+1; i++) {
 249     dest->index[i] = src->index[i];
 250   }
 251   dest->nra = src->nra;
 252   dest->nalloc_a = dest->nra + 1;
 253   snew(dest->a,dest->nalloc_a);
 254   for(i=0; i<dest->nra+1; i++) {
 255     dest->a[i] = src->a[i];
 256   }
 257 }
 258
 259 void done_block(t_block *block)
 260 {
 261   block->nr    = 0;
 262   sfree(block->index);
 263   block->nalloc_index = 0;
 264 }
 265
 266 void done_blocka(t_blocka *block)
 267 {
 268   block->nr    = 0;
 269   block->nra   = 0;
 270   sfree(block->index);
 271   if (block->a)
 272     sfree(block->a);
 273   block->nalloc_index = 0;
 274   block->nalloc_a = 0;
 275 }
 276
 277 void done_atom (t_atoms *at)
 278 {
 279   at->nr       = 0;
 280   at->nres     = 0;
 281   sfree(at->atom);
 282   sfree(at->resinfo);
 283   sfree(at->atomname);
 284   sfree(at->atomtype);
 285   sfree(at->atomtypeB);
 286 }
 287
 288 void done_atomtypes(t_atomtypes *atype)
 289 {
 290   atype->nr = 0;
 291   sfree(atype->radius);
 292   sfree(atype->vol);
 293   sfree(atype->surftens);
 294   sfree(atype->atomnumber);
 295   sfree(atype->gb_radius);
 296   sfree(atype->S_hct);
 297 }
 298
 299 void done_moltype(gmx_moltype_t *molt)
 300 {
 301   int f;
 302
 303   done_atom(&molt->atoms);
 304   done_block(&molt->cgs);
 305   done_blocka(&molt->excls);
 306
 307   for(f=0; f<F_NRE; f++) {
 308     sfree(molt->ilist[f].iatoms);
 309     molt->ilist[f].nalloc = 0;
 310   }
 311 }
 312
 313 void done_molblock(gmx_molblock_t *molb)
 314 {
 315   if (molb->nposres_xA > 0) {
 316     molb->nposres_xA = 0;
 317     free(molb->posres_xA);
 318   }
 319   if (molb->nposres_xB > 0) {
 320     molb->nposres_xB = 0;
 321     free(molb->posres_xB);
 322   }
 323 }
 324
 325 void done_mtop(gmx_mtop_t *mtop,gmx_bool bDoneSymtab)
 326 {
 327   int i;
 328
 329   if (bDoneSymtab) {
 330     done_symtab(&mtop->symtab);
 331   }
 332
 333   sfree(mtop->ffparams.functype);
 334   sfree(mtop->ffparams.iparams);
 335
 336   for(i=0; i<mtop->nmoltype; i++) {
 337     done_moltype(&mtop->moltype[i]);
 338   }
 339   sfree(mtop->moltype);
 340   for(i=0; i<mtop->nmolblock; i++) {
 341     done_molblock(&mtop->molblock[i]);
 342   }
 343   sfree(mtop->molblock);
 344   done_block(&mtop->mols);
 345 }
 346
 347 void done_top(t_topology *top)
 348 {
 349   int f;
 350
 351   sfree(top->idef.functype);
 352   sfree(top->idef.iparams);
 353   for (f = 0; f < F_NRE; ++f)
 354   {
 355       sfree(top->idef.il[f].iatoms);
 356       top->idef.il[f].iatoms = NULL;
 357       top->idef.il[f].nalloc = 0;
 358   }
 359
 360   done_atom (&(top->atoms));
 361
 362   /* For GB */
 363   done_atomtypes(&(top->atomtypes));
 364
 365   done_symtab(&(top->symtab));
 366   done_block(&(top->cgs));
 367   done_block(&(top->mols));
 368   done_blocka(&(top->excls));
 369 }
 370
 371 static void done_pullgrp(t_pullgrp *pgrp)
 372 {
 373   sfree(pgrp->ind);
 374   sfree(pgrp->ind_loc);
 375   sfree(pgrp->weight);
 376   sfree(pgrp->weight_loc);
 377 }
 378
 379 static void done_pull(t_pull *pull)
 380 {
 381   int i;
 382
 383   for(i=0; i<pull->ngrp+1; i++) {
 384     done_pullgrp(pull->grp);
 385     done_pullgrp(pull->dyna);
 386   }
 387 }
 388
 389 void done_inputrec(t_inputrec *ir)
 390 {
 391   int m;
 392
 393   for(m=0; (m<DIM); m++) {
 394     if (ir->ex[m].a)   sfree(ir->ex[m].a);
 395     if (ir->ex[m].phi) sfree(ir->ex[m].phi);
 396     if (ir->et[m].a)   sfree(ir->et[m].a);
 397     if (ir->et[m].phi) sfree(ir->et[m].phi);
 398   }
 399
 400   sfree(ir->opts.nrdf);
 401   sfree(ir->opts.ref_t);
 402   sfree(ir->opts.annealing);
 403   sfree(ir->opts.anneal_npoints);
 404   sfree(ir->opts.anneal_time);
 405   sfree(ir->opts.anneal_temp);
 406   sfree(ir->opts.tau_t);
 407   sfree(ir->opts.acc);
 408   sfree(ir->opts.nFreeze);
 409   sfree(ir->opts.QMmethod);
 410   sfree(ir->opts.QMbasis);
 411   sfree(ir->opts.QMcharge);
 412   sfree(ir->opts.QMmult);
 413   sfree(ir->opts.bSH);
 414   sfree(ir->opts.CASorbitals);
 415   sfree(ir->opts.CASelectrons);
 416   sfree(ir->opts.SAon);
 417   sfree(ir->opts.SAoff);
 418   sfree(ir->opts.SAsteps);
 419   sfree(ir->opts.bOPT);
 420   sfree(ir->opts.bTS);
 421
 422   if (ir->pull) {
 423     done_pull(ir->pull);
 424     sfree(ir->pull);
 425   }
 426 }
 427
 428 static void zero_ekinstate(ekinstate_t *eks)
 429 {
 430   eks->ekin_n         = 0;
 431   eks->ekinh          = NULL;
 432   eks->ekinf          = NULL;
 433   eks->ekinh_old      = NULL;
 434   eks->ekinscalef_nhc = NULL;
 435   eks->ekinscaleh_nhc = NULL;
 436   eks->vscale_nhc     = NULL;
 437   eks->dekindl        = 0;
 438   eks->mvcos          = 0;
 439 }
 440
 441 void init_energyhistory(energyhistory_t * enerhist)
 442 {
 443     enerhist->nener = 0;
 444
 445     enerhist->ener_ave     = NULL;
 446     enerhist->ener_sum     = NULL;
 447     enerhist->ener_sum_sim = NULL;
 448     enerhist->dht          = NULL;
 449
 450     enerhist->nsteps     = 0;
 451     enerhist->nsum       = 0;
 452     enerhist->nsteps_sim = 0;
 453     enerhist->nsum_sim   = 0;
 454
 455     enerhist->dht = NULL;
 456 }
 457
 458 static void done_delta_h_history(delta_h_history_t *dht)
 459 {
 460     int i;
 461
 462     for(i=0; i<dht->nndh; i++)
 463     {
 464         sfree(dht->dh[i]);
 465     }
 466     sfree(dht->dh);
 467     sfree(dht->ndh);
 468 }
 469
 470 void done_energyhistory(energyhistory_t * enerhist)
 471 {
 472     sfree(enerhist->ener_ave);
 473     sfree(enerhist->ener_sum);
 474     sfree(enerhist->ener_sum_sim);
 475
 476     if (enerhist->dht != NULL)
 477     {
 478         done_delta_h_history(enerhist->dht);
 479         sfree(enerhist->dht);
 480     }
 481 }
 482
 483 void init_gtc_state(t_state *state, int ngtc, int nnhpres, int nhchainlength)
 484 {
 485     int i,j;
 486
 487     state->ngtc = ngtc;
 488     state->nnhpres = nnhpres;
 489     state->nhchainlength = nhchainlength;
 490     if (state->ngtc > 0)
 491     {
 492         snew(state->nosehoover_xi,state->nhchainlength*state->ngtc);
 493         snew(state->nosehoover_vxi,state->nhchainlength*state->ngtc);
 494         snew(state->therm_integral,state->ngtc);
 495         for(i=0; i<state->ngtc; i++)
 496         {
 497             for (j=0;j<state->nhchainlength;j++)
 498             {
 499                 state->nosehoover_xi[i*state->nhchainlength + j]  = 0.0;
 500                 state->nosehoover_vxi[i*state->nhchainlength + j]  = 0.0;
 501             }
 502         }
 503         for(i=0; i<state->ngtc; i++) {
 504             state->therm_integral[i]  = 0.0;
 505         }
 506     }
 507     else
 508     {
 509         state->nosehoover_xi  = NULL;
 510         state->nosehoover_vxi = NULL;
 511         state->therm_integral = NULL;
 512     }
 513
 514     if (state->nnhpres > 0)
 515     {
 516         snew(state->nhpres_xi,state->nhchainlength*nnhpres);
 517         snew(state->nhpres_vxi,state->nhchainlength*nnhpres);
 518         for(i=0; i<nnhpres; i++)
 519         {
 520             for (j=0;j<state->nhchainlength;j++)
 521             {
 522                 state->nhpres_xi[i*nhchainlength + j]  = 0.0;
 523                 state->nhpres_vxi[i*nhchainlength + j]  = 0.0;
 524             }
 525         }
 526     }
 527     else
 528     {
 529         state->nhpres_xi  = NULL;
 530         state->nhpres_vxi = NULL;
 531     }
 532 }
 533
 534
 535 void init_state(t_state *state, int natoms, int ngtc, int nnhpres, int nhchainlength, int nlambda)
 536 {
 537   int i;
 538
 539   state->natoms = natoms;
 540   state->nrng   = 0;
 541   state->flags  = 0;
 542   state->lambda = 0;
 543   snew(state->lambda,efptNR);
 544   for (i=0;i<efptNR;i++)
 545   {
 546       state->lambda[i] = 0;
 547   }
 548   state->veta   = 0;
 549   clear_mat(state->box);
 550   clear_mat(state->box_rel);
 551   clear_mat(state->boxv);
 552   clear_mat(state->pres_prev);
 553   clear_mat(state->svir_prev);
 554   clear_mat(state->fvir_prev);
 555   init_gtc_state(state,ngtc,nnhpres,nhchainlength);
 556   state->nalloc = state->natoms;
 557   if (state->nalloc > 0) {
 558     snew(state->x,state->nalloc);
 559     snew(state->v,state->nalloc);
 560   } else {
 561     state->x = NULL;
 562     state->v = NULL;
 563   }
 564   state->sd_X = NULL;
 565   state->cg_p = NULL;
 566
 567   zero_ekinstate(&state->ekinstate);
 568
 569   init_energyhistory(&state->enerhist);
 570
 571   init_df_history(&state->dfhist,nlambda,0);
 572
 573   state->ddp_count = 0;
 574   state->ddp_count_cg_gl = 0;
 575   state->cg_gl = NULL;
 576   state->cg_gl_nalloc = 0;
 577 }
 578
 579 void done_state(t_state *state)
 580 {
 581   if (state->nosehoover_xi) sfree(state->nosehoover_xi);
 582   if (state->x) sfree(state->x);
 583   if (state->v) sfree(state->v);
 584   if (state->sd_X) sfree(state->sd_X);
 585   if (state->cg_p) sfree(state->cg_p);
 586   state->nalloc = 0;
 587   if (state->cg_gl) sfree(state->cg_gl);
 588   state->cg_gl_nalloc = 0;
 589 }
 590
 591 static void do_box_rel(t_inputrec *ir,matrix box_rel,matrix b,gmx_bool bInit)
 592 {
 593   int d,d2;
 594
 595   for(d=YY; d<=ZZ; d++) {
 596     for(d2=XX; d2<=(ir->epct==epctSEMIISOTROPIC ? YY : ZZ); d2++) {
 597       /* We need to check if this box component is deformed
 598        * or if deformation of another component might cause
 599        * changes in this component due to box corrections.
 600        */
 601       if (ir->deform[d][d2] == 0 &&
 602           !(d == ZZ && d2 == XX && ir->deform[d][YY] != 0 &&
 603             (b[YY][d2] != 0 || ir->deform[YY][d2] != 0))) {
 604         if (bInit) {
 605           box_rel[d][d2] = b[d][d2]/b[XX][XX];
 606         } else {
 607           b[d][d2] = b[XX][XX]*box_rel[d][d2];
 608         }
 609       }
 610     }
 611   }
 612 }
 613
 614 void set_box_rel(t_inputrec *ir,t_state *state)
 615 {
 616   /* Make sure the box obeys the restrictions before we fix the ratios */
 617   correct_box(NULL,0,state->box,NULL);
 618
 619   clear_mat(state->box_rel);
 620
 621   if (PRESERVE_SHAPE(*ir))
 622     do_box_rel(ir,state->box_rel,state->box,TRUE);
 623 }
 624
 625 void preserve_box_shape(t_inputrec *ir,matrix box_rel,matrix b)
 626 {
 627   if (PRESERVE_SHAPE(*ir))
 628     do_box_rel(ir,box_rel,b,FALSE);
 629 }
 630
 631 void add_t_atoms(t_atoms *atoms,int natom_extra,int nres_extra)
 632 {
 633     int i;
 634
 635     if (natom_extra > 0)
 636     {
 637         srenew(atoms->atomname,atoms->nr+natom_extra);
 638         srenew(atoms->atom,atoms->nr+natom_extra);
 639         if (NULL != atoms->pdbinfo)
 640             srenew(atoms->pdbinfo,atoms->nr+natom_extra);
 641         if (NULL != atoms->atomtype)
 642             srenew(atoms->atomtype,atoms->nr+natom_extra);
 643         if (NULL != atoms->atomtypeB)
 644             srenew(atoms->atomtypeB,atoms->nr+natom_extra);
 645         for(i=atoms->nr; (i<atoms->nr+natom_extra); i++) {
 646             atoms->atomname[i] = NULL;
 647             memset(&atoms->atom[i],0,sizeof(atoms->atom[i]));
 648             if (NULL != atoms->pdbinfo)
 649                 memset(&atoms->pdbinfo[i],0,sizeof(atoms->pdbinfo[i]));
 650             if (NULL != atoms->atomtype)
 651                 atoms->atomtype[i] = NULL;
 652             if (NULL != atoms->atomtypeB)
 653                 atoms->atomtypeB[i] = NULL;
 654         }
 655         atoms->nr += natom_extra;
 656     }
 657     if (nres_extra > 0)
 658     {
 659         srenew(atoms->resinfo,atoms->nres+nres_extra);
 660         for(i=atoms->nres; (i<atoms->nres+nres_extra); i++) {
 661             memset(&atoms->resinfo[i],0,sizeof(atoms->resinfo[i]));
 662         }
 663         atoms->nres += nres_extra;
 664     }
 665 }
 666
 667 void init_t_atoms(t_atoms *atoms, int natoms, gmx_bool bPdbinfo)
 668 {
 669   atoms->nr=natoms;
 670   atoms->nres=0;
 671   snew(atoms->atomname,natoms);
 672   atoms->atomtype=NULL;
 673   atoms->atomtypeB=NULL;
 674   snew(atoms->resinfo,natoms);
 675   snew(atoms->atom,natoms);
 676   if (bPdbinfo)
 677     snew(atoms->pdbinfo,natoms);
 678   else
 679     atoms->pdbinfo=NULL;
 680 }
 681
 682 t_atoms *copy_t_atoms(t_atoms *src)
 683 {
 684   t_atoms *dst;
 685   int i;
 686
 687   snew(dst,1);
 688   init_t_atoms(dst,src->nr,(NULL != src->pdbinfo));
 689   dst->nr = src->nr;
 690   if (NULL != src->atomname)
 691       snew(dst->atomname,src->nr);
 692   if (NULL != src->atomtype)
 693       snew(dst->atomtype,src->nr);
 694   if (NULL != src->atomtypeB)
 695       snew(dst->atomtypeB,src->nr);
 696   for(i=0; (i<src->nr); i++) {
 697     dst->atom[i] = src->atom[i];
 698     if (NULL != src->pdbinfo)
 699       dst->pdbinfo[i] = src->pdbinfo[i];
 700     if (NULL != src->atomname)
 701         dst->atomname[i]  = src->atomname[i];
 702     if (NULL != src->atomtype)
 703         dst->atomtype[i] = src->atomtype[i];
 704     if (NULL != src->atomtypeB)
 705         dst->atomtypeB[i] = src->atomtypeB[i];
 706   }
 707   dst->nres = src->nres;
 708   for(i=0; (i<src->nres); i++) {
 709     dst->resinfo[i] = src->resinfo[i];
 710   }
 711   return dst;
 712 }
 713
 714 void t_atoms_set_resinfo(t_atoms *atoms,int atom_ind,t_symtab *symtab,
 715                          const char *resname,int resnr,unsigned char ic,
 716                          int chainnum, char chainid)
 717 {
 718   t_resinfo *ri;
 719
 720   ri = &atoms->resinfo[atoms->atom[atom_ind].resind];
 721   ri->name  = put_symtab(symtab,resname);
 722   ri->rtp   = NULL;
 723   ri->nr    = resnr;
 724   ri->ic    = ic;
 725   ri->chainnum = chainnum;
 726   ri->chainid = chainid;
 727 }
 728
 729 void free_t_atoms(t_atoms *atoms,gmx_bool bFreeNames)
 730 {
 731   int i;
 732
 733   if (bFreeNames) {
 734     for(i=0; i<atoms->nr; i++) {
 735       sfree(*atoms->atomname[i]);
 736       *atoms->atomname[i]=NULL;
 737     }
 738     for(i=0; i<atoms->nres; i++) {
 739       sfree(*atoms->resinfo[i].name);
 740       *atoms->resinfo[i].name=NULL;
 741     }
 742   }
 743   sfree(atoms->atomname);
 744   /* Do we need to free atomtype and atomtypeB as well ? */
 745   sfree(atoms->resinfo);
 746   sfree(atoms->atom);
 747   if (atoms->pdbinfo)
 748     sfree(atoms->pdbinfo);
 749   atoms->nr=0;
 750   atoms->nres=0;
 751   atoms->atomname = NULL;
 752   atoms->resinfo = NULL;
 753   atoms->atom = NULL;
 754   atoms->pdbinfo = NULL;
 755 }
 756
 757 real max_cutoff(real cutoff1,real cutoff2)
 758 {
 759     if (cutoff1 == 0 || cutoff2 == 0)
 760     {
 761         return 0;
 762     }
 763     else
 764     {
 765         return max(cutoff1,cutoff2);
 766     }
 767 }
 768
 769 extern void init_df_history(df_history_t *dfhist, int nlambda, real wl_delta)
 770 {
 771     int i;
 772
 773     dfhist->bEquil = 0;
 774     dfhist->nlambda = nlambda;
 775     dfhist->wl_delta = wl_delta;
 776     snew(dfhist->sum_weights,dfhist->nlambda);
 777     snew(dfhist->sum_dg,dfhist->nlambda);
 778     snew(dfhist->sum_minvar,dfhist->nlambda);
 779     snew(dfhist->sum_variance,dfhist->nlambda);
 780     snew(dfhist->n_at_lam,dfhist->nlambda);
 781     snew(dfhist->wl_histo,dfhist->nlambda);
 782
 783     /* allocate transition matrices here */
 784     snew(dfhist->Tij,dfhist->nlambda);
 785     snew(dfhist->Tij_empirical,dfhist->nlambda);
 786
 787     for (i=0;i<dfhist->nlambda;i++) {
 788         snew(dfhist->Tij[i],dfhist->nlambda);
 789         snew(dfhist->Tij_empirical[i],dfhist->nlambda);
 790     }
 791
 792     snew(dfhist->accum_p,dfhist->nlambda);
 793     snew(dfhist->accum_m,dfhist->nlambda);
 794     snew(dfhist->accum_p2,dfhist->nlambda);
 795     snew(dfhist->accum_m2,dfhist->nlambda);
 796
 797     for (i=0;i<dfhist->nlambda;i++) {
 798         snew((dfhist->accum_p)[i],dfhist->nlambda);
 799         snew((dfhist->accum_m)[i],dfhist->nlambda);
 800         snew((dfhist->accum_p2)[i],dfhist->nlambda);
 801         snew((dfhist->accum_m2)[i],dfhist->nlambda);
 802     }
 803 }
 804
 805 extern void copy_df_history(df_history_t *df_dest, df_history_t *df_source)
 806 {
 807     int i,j;
 808
 809     init_df_history(df_dest,df_source->nlambda,df_source->wl_delta);
 810     df_dest->nlambda = df_source->nlambda;
 811     df_dest->bEquil = df_source->bEquil;
 812     for (i=0;i<df_dest->nlambda;i++)
 813     {
 814         df_dest->sum_weights[i]  = df_source->sum_weights[i];
 815         df_dest->sum_dg[i]       = df_source->sum_dg[i];
 816         df_dest->sum_minvar[i]   = df_source->sum_minvar[i];
 817         df_dest->sum_variance[i] = df_source->sum_variance[i];
 818         df_dest->n_at_lam[i]     = df_source->n_at_lam[i];
 819         df_dest->wl_histo[i]     = df_source->wl_histo[i];
 820         df_dest->accum_p[i]      = df_source->accum_p[i];
 821         df_dest->accum_m[i]      = df_source->accum_m[i];
 822         df_dest->accum_p2[i]     = df_source->accum_p2[i];
 823         df_dest->accum_m2[i]     = df_source->accum_m2[i];
 824     }
 825
 826     for (i=0;i<df_dest->nlambda;i++)
 827     {
 828         for (j=0;j<df_dest->nlambda;j++)
 829         {
 830             df_dest->Tij[i][j]  = df_source->Tij[i][j];
 831             df_dest->Tij_empirical[i][j]  = df_source->Tij_empirical[i][j];
 832         }
 833     }
 834 }