?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
t_atoms_set_resinfo now sets rtp to NULL
[gromacs/rigid-bodies.git] / src / gmxlib / typedefs.c
blob2fe950246a8da3e824bfe58edc83e0bde08ae870
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 <string.h>
46
47 #ifdef GMX_THREADS
48 #include "thread_mpi.h"
49 #endif
50
51 /* The source code in this file should be thread-safe.
52 Please keep it that way. */
53
54
55
56 static bool bOverAllocDD=FALSE;
57 #ifdef GMX_THREADS
58 static tMPI_Thread_mutex_t over_alloc_mutex=TMPI_THREAD_MUTEX_INITIALIZER;
59 #endif
60
61
62 void set_over_alloc_dd(bool set)
63 {
64 #ifdef GMX_THREADS
65 tMPI_Thread_mutex_lock(&over_alloc_mutex);
66 /* we just make sure that we don't set this at the same time.
67 We don't worry too much about reading this rarely-set variable */
68 #endif
69 bOverAllocDD = set;
70 #ifdef GMX_THREADS
71 tMPI_Thread_mutex_unlock(&over_alloc_mutex);
72 #endif
73 }
74
75 int over_alloc_dd(int n)
76 {
77 if (bOverAllocDD)
78 return OVER_ALLOC_FAC*n + 100;
79 else
80 return n;
81 }
82
83 int gmx_large_int_to_int(gmx_large_int_t step,const char *warn)
84 {
85 int i;
86
87 i = (int)step;
88
89 if (warn != NULL && (step < INT_MIN || step > INT_MAX)) {
90 fprintf(stderr,"\nWARNING during %s:\n",warn);
91 fprintf(stderr,"step value ");
92 fprintf(stderr,gmx_large_int_pfmt,step);
93 fprintf(stderr," does not fit in int, converted to %d\n\n",i);
94 }
95
96 return i;
97 }
98
99 char *gmx_step_str(gmx_large_int_t i,char *buf)
100 {
101 sprintf(buf,gmx_large_int_pfmt,i);
102
103 return buf;
104 }
105
106 void init_block(t_block *block)
107 {
108 int i;
109
110 block->nr = 0;
111 block->nalloc_index = 1;
112 snew(block->index,block->nalloc_index);
113 block->index[0] = 0;
114 }
115
116 void init_blocka(t_blocka *block)
117 {
118 int i;
119
120 block->nr = 0;
121 block->nra = 0;
122 block->nalloc_index = 1;
123 snew(block->index,block->nalloc_index);
124 block->index[0] = 0;
125 block->nalloc_a = 0;
126 block->a = NULL;
127 }
128
129 void init_atom(t_atoms *at)
130 {
131 int i;
132
133 at->nr = 0;
134 at->nres = 0;
135 at->atom = NULL;
136 at->resinfo = NULL;
137 at->atomname = NULL;
138 at->atomtype = NULL;
139 at->atomtypeB= NULL;
140 at->pdbinfo = NULL;
141 }
142
143 void init_atomtypes(t_atomtypes *at)
144 {
145 at->nr = 0;
146 at->radius = NULL;
147 at->vol = NULL;
148 at->atomnumber = NULL;
149 at->gb_radius = NULL;
150 at->S_hct = NULL;
151 }
152
153 void init_groups(gmx_groups_t *groups)
154 {
155 int g;
156
157 groups->ngrpname = 0;
158 groups->grpname = NULL;
159 for(g=0; (g<egcNR); g++) {
160 groups->grps[g].nm_ind = NULL;
161 groups->ngrpnr[g] = 0;
162 groups->grpnr[g] = NULL;
163 }
164
165 }
166
167 void init_mtop(gmx_mtop_t *mtop)
168 {
169 mtop->name = NULL;
170 mtop->nmoltype = 0;
171 mtop->moltype = NULL;
172 mtop->nmolblock = 0;
173 mtop->molblock = NULL;
174 mtop->maxres_renum = 0;
175 mtop->maxresnr = -1;
176 init_groups(&mtop->groups);
177 init_block(&mtop->mols);
178 open_symtab(&mtop->symtab);
179 }
180
181 void init_top (t_topology *top)
182 {
183 int i;
184
185 top->name = NULL;
186 init_atom (&(top->atoms));
187 init_atomtypes(&(top->atomtypes));
188 init_block(&top->cgs);
189 init_block(&top->mols);
190 init_blocka(&top->excls);
191 open_symtab(&top->symtab);
192 }
193
194 void init_inputrec(t_inputrec *ir)
195 {
196 memset(ir,0,(size_t)sizeof(*ir));
197 }
198
199 void stupid_fill_block(t_block *grp,int natom,bool bOneIndexGroup)
200 {
201 int i;
202
203 if (bOneIndexGroup) {
204 grp->nalloc_index = 2;
205 snew(grp->index,grp->nalloc_index);
206 grp->index[0]=0;
207 grp->index[1]=natom;
208 grp->nr=1;
209 }
210 else {
211 grp->nalloc_index = natom+1;
212 snew(grp->index,grp->nalloc_index);
213 snew(grp->index,natom+1);
214 for(i=0; (i<=natom); i++)
215 grp->index[i]=i;
216 grp->nr=natom;
217 }
218 }
219
220 void stupid_fill_blocka(t_blocka *grp,int natom)
221 {
222 int i;
223
224 grp->nalloc_a = natom;
225 snew(grp->a,grp->nalloc_a);
226 for(i=0; (i<natom); i++)
227 grp->a[i]=i;
228 grp->nra=natom;
229
230 grp->nalloc_index = natom + 1;
231 snew(grp->index,grp->nalloc_index);
232 for(i=0; (i<=natom); i++)
233 grp->index[i]=i;
234 grp->nr=natom;
235 }
236
237 void copy_blocka(const t_blocka *src,t_blocka *dest)
238 {
239 int i;
240
241 dest->nr = src->nr;
242 dest->nalloc_index = dest->nr + 1;
243 snew(dest->index,dest->nalloc_index);
244 for(i=0; i<dest->nr+1; i++) {
245 dest->index[i] = src->index[i];
246 }
247 dest->nra = src->nra;
248 dest->nalloc_a = dest->nra + 1;
249 snew(dest->a,dest->nalloc_a);
250 for(i=0; i<dest->nra+1; i++) {
251 dest->a[i] = src->a[i];
252 }
253 }
254
255 void done_block(t_block *block)
256 {
257 block->nr = 0;
258 sfree(block->index);
259 block->nalloc_index = 0;
260 }
261
262 void done_blocka(t_blocka *block)
263 {
264 block->nr = 0;
265 block->nra = 0;
266 sfree(block->index);
267 if (block->a)
268 sfree(block->a);
269 block->nalloc_index = 0;
270 block->nalloc_a = 0;
271 }
272
273 void done_atom (t_atoms *at)
274 {
275 at->nr = 0;
276 at->nres = 0;
277 sfree(at->atom);
278 sfree(at->resinfo);
279 sfree(at->atomname);
280 }
281
282 void done_atomtypes(t_atomtypes *atype)
283 {
284 atype->nr = 0;
285 sfree(atype->radius);
286 sfree(atype->vol);
287 sfree(atype->surftens);
288 sfree(atype->gb_radius);
289 sfree(atype->S_hct);
290 }
291
292 void done_moltype(gmx_moltype_t *molt)
293 {
294 int f;
295
296 done_atom(&molt->atoms);
297 done_block(&molt->cgs);
298 done_blocka(&molt->excls);
299
300 for(f=0; f<F_NRE; f++) {
301 sfree(molt->ilist[f].iatoms);
302 molt->ilist[f].nalloc = 0;
303 }
304 }
305
306 void done_molblock(gmx_molblock_t *molb)
307 {
308 if (molb->nposres_xA > 0) {
309 molb->nposres_xA = 0;
310 free(molb->posres_xA);
311 }
312 if (molb->nposres_xB > 0) {
313 molb->nposres_xB = 0;
314 free(molb->posres_xB);
315 }
316 }
317
318 void done_mtop(gmx_mtop_t *mtop,bool bDoneSymtab)
319 {
320 int i;
321
322 if (bDoneSymtab) {
323 done_symtab(&mtop->symtab);
324 }
325
326 sfree(mtop->ffparams.functype);
327 sfree(mtop->ffparams.iparams);
328
329 for(i=0; i<mtop->nmoltype; i++) {
330 done_moltype(&mtop->moltype[i]);
331 }
332 sfree(mtop->moltype);
333 for(i=0; i<mtop->nmolblock; i++) {
334 done_molblock(&mtop->molblock[i]);
335 }
336 sfree(mtop->molblock);
337 done_block(&mtop->mols);
338 }
339
340 void done_top(t_topology *top)
341 {
342 int i;
343
344 done_atom (&(top->atoms));
345
346 /* For GB */
347 done_atomtypes(&(top->atomtypes));
348
349 done_symtab(&(top->symtab));
350 done_block(&(top->cgs));
351 done_block(&(top->mols));
352 done_blocka(&(top->excls));
353 }
354
355 static void done_pullgrp(t_pullgrp *pgrp)
356 {
357 sfree(pgrp->ind);
358 sfree(pgrp->ind_loc);
359 sfree(pgrp->weight);
360 sfree(pgrp->weight_loc);
361 }
362
363 static void done_pull(t_pull *pull)
364 {
365 int i;
366
367 for(i=0; i<pull->ngrp+1; i++) {
368 done_pullgrp(pull->grp);
369 done_pullgrp(pull->dyna);
370 }
371 }
372
373 void done_inputrec(t_inputrec *ir)
374 {
375 int m;
376
377 for(m=0; (m<DIM); m++) {
378 if (ir->ex[m].a) sfree(ir->ex[m].a);
379 if (ir->ex[m].phi) sfree(ir->ex[m].phi);
380 if (ir->et[m].a) sfree(ir->et[m].a);
381 if (ir->et[m].phi) sfree(ir->et[m].phi);
382 }
383
384 sfree(ir->opts.nrdf);
385 sfree(ir->opts.ref_t);
386 sfree(ir->opts.annealing);
387 sfree(ir->opts.anneal_npoints);
388 sfree(ir->opts.anneal_time);
389 sfree(ir->opts.anneal_temp);
390 sfree(ir->opts.tau_t);
391 sfree(ir->opts.acc);
392 sfree(ir->opts.nFreeze);
393 sfree(ir->opts.QMmethod);
394 sfree(ir->opts.QMbasis);
395 sfree(ir->opts.QMcharge);
396 sfree(ir->opts.QMmult);
397 sfree(ir->opts.bSH);
398 sfree(ir->opts.CASorbitals);
399 sfree(ir->opts.CASelectrons);
400 sfree(ir->opts.SAon);
401 sfree(ir->opts.SAoff);
402 sfree(ir->opts.SAsteps);
403 sfree(ir->opts.bOPT);
404 sfree(ir->opts.bTS);
405
406 if (ir->pull) {
407 done_pull(ir->pull);
408 sfree(ir->pull);
409 }
410 }
411
412 static void init_ekinstate(ekinstate_t *eks)
413 {
414 eks->ekin_n = 0;
415 eks->ekinh = NULL;
416 eks->ekinf = NULL;
417 eks->ekinh_old = NULL;
418 eks->ekinscalef_nhc = NULL;
419 eks->ekinscaleh_nhc = NULL;
420 eks->vscale_nhc = NULL;
421 eks->dekindl = 0;
422 eks->mvcos = 0;
423 }
424
425 static void init_energyhistory(energyhistory_t *enh)
426 {
427 enh->ener_ave = NULL;
428 enh->ener_sum = NULL;
429 enh->ener_sum_sim = NULL;
430 enh->nener = 0;
431 }
432
433 void init_gtc_state(t_state *state, int ngtc, int nnhpres, int nhchainlength)
434 {
435 int i,j;
436
437 state->ngtc = ngtc;
438 state->nnhpres = nnhpres;
439 state->nhchainlength = nhchainlength;
440 if (state->ngtc > 0)
441 {
442 snew(state->nosehoover_xi,state->nhchainlength*state->ngtc);
443 snew(state->nosehoover_vxi,state->nhchainlength*state->ngtc);
444 snew(state->therm_integral,state->ngtc);
445 for(i=0; i<state->ngtc; i++)
446 {
447 for (j=0;j<state->nhchainlength;j++)
448 {
449 state->nosehoover_xi[i*state->nhchainlength + j] = 0.0;
450 state->nosehoover_vxi[i*state->nhchainlength + j] = 0.0;
451 }
452 }
453 for(i=0; i<state->ngtc; i++) {
454 state->therm_integral[i] = 0.0;
455 }
456 }
457 else
458 {
459 state->nosehoover_xi = NULL;
460 state->nosehoover_vxi = NULL;
461 state->therm_integral = NULL;
462 }
463
464 if (state->nnhpres > 0)
465 {
466 snew(state->nhpres_xi,state->nhchainlength*nnhpres);
467 snew(state->nhpres_vxi,state->nhchainlength*nnhpres);
468 for(i=0; i<nnhpres; i++)
469 {
470 for (j=0;j<state->nhchainlength;j++)
471 {
472 state->nhpres_xi[i*nhchainlength + j] = 0.0;
473 state->nhpres_vxi[i*nhchainlength + j] = 0.0;
474 }
475 }
476 }
477 else
478 {
479 state->nhpres_xi = NULL;
480 state->nhpres_vxi = NULL;
481 }
482 }
483
484
485 void init_state(t_state *state, int natoms, int ngtc, int nnhpres, int nhchainlength)
486 {
487 int i;
488
489 state->natoms = natoms;
490 state->nrng = 0;
491 state->flags = 0;
492 state->lambda = 0;
493 state->veta = 0;
494 clear_mat(state->box);
495 clear_mat(state->box_rel);
496 clear_mat(state->boxv);
497 clear_mat(state->pres_prev);
498 clear_mat(state->svir_prev);
499 clear_mat(state->fvir_prev);
500 init_gtc_state(state,ngtc,nnhpres,nhchainlength);
501 state->nalloc = state->natoms;
502 if (state->nalloc > 0) {
503 snew(state->x,state->nalloc);
504 snew(state->v,state->nalloc);
505 } else {
506 state->x = NULL;
507 state->v = NULL;
508 }
509 state->sd_X = NULL;
510 state->cg_p = NULL;
511
512 init_ekinstate(&state->ekinstate);
513
514 init_energyhistory(&state->enerhist);
515
516 state->ddp_count = 0;
517 state->ddp_count_cg_gl = 0;
518 state->cg_gl = NULL;
519 state->cg_gl_nalloc = 0;
520 }
521
522 void done_state(t_state *state)
523 {
524 if (state->nosehoover_xi) sfree(state->nosehoover_xi);
525 if (state->x) sfree(state->x);
526 if (state->v) sfree(state->v);
527 if (state->sd_X) sfree(state->sd_X);
528 if (state->cg_p) sfree(state->cg_p);
529 state->nalloc = 0;
530 if (state->cg_gl) sfree(state->cg_gl);
531 state->cg_gl_nalloc = 0;
532 }
533
534 static void do_box_rel(t_inputrec *ir,matrix box_rel,matrix b,bool bInit)
535 {
536 int d,d2;
537
538 for(d=YY; d<=ZZ; d++) {
539 for(d2=XX; d2<=(ir->epct==epctSEMIISOTROPIC ? YY : ZZ); d2++) {
540 /* We need to check if this box component is deformed
541 * or if deformation of another component might cause
542 * changes in this component due to box corrections.
543 */
544 if (ir->deform[d][d2] == 0 &&
545 !(d == ZZ && d2 == XX && ir->deform[d][YY] != 0 &&
546 (b[YY][d2] != 0 || ir->deform[YY][d2] != 0))) {
547 if (bInit) {
548 box_rel[d][d2] = b[d][d2]/b[XX][XX];
549 } else {
550 b[d][d2] = b[XX][XX]*box_rel[d][d2];
551 }
552 }
553 }
554 }
555 }
556
557 void set_box_rel(t_inputrec *ir,t_state *state)
558 {
559 /* Make sure the box obeys the restrictions before we fix the ratios */
560 correct_box(NULL,0,state->box,NULL);
561
562 clear_mat(state->box_rel);
563
564 if (PRESERVE_SHAPE(*ir))
565 do_box_rel(ir,state->box_rel,state->box,TRUE);
566 }
567
568 void preserve_box_shape(t_inputrec *ir,matrix box_rel,matrix b)
569 {
570 if (PRESERVE_SHAPE(*ir))
571 do_box_rel(ir,box_rel,b,FALSE);
572 }
573
574 void add_t_atoms(t_atoms *atoms,int natom_extra,int nres_extra)
575 {
576 int i;
577
578 if (natom_extra > 0)
579 {
580 srenew(atoms->atomname,atoms->nr+natom_extra);
581 srenew(atoms->atom,atoms->nr+natom_extra);
582 if (NULL != atoms->pdbinfo)
583 srenew(atoms->pdbinfo,atoms->nr+natom_extra);
584 if (NULL != atoms->atomtype)
585 srenew(atoms->atomtype,atoms->nr+natom_extra);
586 if (NULL != atoms->atomtypeB)
587 srenew(atoms->atomtypeB,atoms->nr+natom_extra);
588 for(i=atoms->nr; (i<atoms->nr+natom_extra); i++) {
589 atoms->atomname[i] = NULL;
590 memset(&atoms->atom[i],0,sizeof(atoms->atom[i]));
591 if (NULL != atoms->pdbinfo)
592 memset(&atoms->pdbinfo[i],0,sizeof(atoms->pdbinfo[i]));
593 if (NULL != atoms->atomtype)
594 atoms->atomtype[i] = NULL;
595 if (NULL != atoms->atomtypeB)
596 atoms->atomtypeB[i] = NULL;
597 }
598 atoms->nr += natom_extra;
599 }
600 if (nres_extra > 0)
601 {
602 srenew(atoms->resinfo,atoms->nres+nres_extra);
603 for(i=atoms->nres; (i<atoms->nres+nres_extra); i++) {
604 memset(&atoms->resinfo[i],0,sizeof(atoms->resinfo[i]));
605 }
606 atoms->nres += nres_extra;
607 }
608 }
609
610 void init_t_atoms(t_atoms *atoms, int natoms, bool bPdbinfo)
611 {
612 atoms->nr=natoms;
613 atoms->nres=0;
614 snew(atoms->atomname,natoms);
615 atoms->atomtype=NULL;
616 atoms->atomtypeB=NULL;
617 snew(atoms->resinfo,natoms);
618 snew(atoms->atom,natoms);
619 if (bPdbinfo)
620 snew(atoms->pdbinfo,natoms);
621 else
622 atoms->pdbinfo=NULL;
623 }
624
625 t_atoms *copy_t_atoms(t_atoms *src)
626 {
627 t_atoms *dst;
628 int i;
629
630 snew(dst,1);
631 init_t_atoms(dst,src->nr,(NULL != src->pdbinfo));
632 dst->nr = src->nr;
633 if (NULL != src->atomname)
634 snew(dst->atomname,src->nr);
635 if (NULL != src->atomtype)
636 snew(dst->atomtype,src->nr);
637 if (NULL != src->atomtypeB)
638 snew(dst->atomtypeB,src->nr);
639 for(i=0; (i<src->nr); i++) {
640 dst->atom[i] = src->atom[i];
641 if (NULL != src->pdbinfo)
642 dst->pdbinfo[i] = src->pdbinfo[i];
643 if (NULL != src->atomname)
644 dst->atomname[i] = src->atomname[i];
645 if (NULL != src->atomtype)
646 dst->atomtype[i] = src->atomtype[i];
647 if (NULL != src->atomtypeB)
648 dst->atomtypeB[i] = src->atomtypeB[i];
649 }
650 dst->nres = src->nres;
651 for(i=0; (i<src->nres); i++) {
652 dst->resinfo[i] = src->resinfo[i];
653 }
654 return dst;
655 }
656
657 void t_atoms_set_resinfo(t_atoms *atoms,int atom_ind,t_symtab *symtab,
658 const char *resname,int resnr,unsigned char ic,
659 unsigned char chain)
660 {
661 t_resinfo *ri;
662
663 ri = &atoms->resinfo[atoms->atom[atom_ind].resind];
664 ri->name = put_symtab(symtab,resname);
665 ri->rtp = NULL;
666 ri->nr = resnr;
667 ri->ic = ic;
668 ri->chain = chain;
669 }
670
671 void free_t_atoms(t_atoms *atoms,bool bFreeNames)
672 {
673 int i;
674
675 if (bFreeNames) {
676 for(i=0; i<atoms->nr; i++) {
677 sfree(*atoms->atomname[i]);
678 *atoms->atomname[i]=NULL;
679 }
680 for(i=0; i<atoms->nres; i++) {
681 sfree(*atoms->resinfo[i].name);
682 *atoms->resinfo[i].name=NULL;
683 }
684 }
685 sfree(atoms->atomname);
686 /* Do we need to free atomtype and atomtypeB as well ? */
687 sfree(atoms->resinfo);
688 sfree(atoms->atom);
689 if (atoms->pdbinfo)
690 sfree(atoms->pdbinfo);
691 atoms->nr=0;
692 atoms->nres=0;
693 }
694
695 real max_cutoff(real cutoff1,real cutoff2)
696 {
697 if (cutoff1 == 0 || cutoff2 == 0)
698 {
699 return 0;
700 }
701 else
702 {
703 return max(cutoff1,cutoff2);
704 }
705 }