src/gmxlib/typedefs.c

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