src/gromacs/gmxpreprocess/add_par.c

   1 /*
   2  *
   3  *                This source code is part of
   4  *
   5  *                 G   R   O   M   A   C   S
   6  *
   7  *          GROningen MAchine for Chemical Simulations
   8  *
   9  *                        VERSION 3.2.0
  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.
  14
  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.
  19  *
  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.
  26  *
  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.
  29  *
  30  * For more info, check our website at http://www.gromacs.org
  31  *
  32  * And Hey:
  33  * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
  34  */
  35 /* This file is completely threadsafe - keep it that way! */
  36 #ifdef HAVE_CONFIG_H
  37 #include <config.h>
  38 #endif
  39
  40 #include <string.h>
  41 #include "typedefs.h"
  42 #include "smalloc.h"
  43 #include "grompp.h"
  44 #include "macros.h"
  45 #include "toputil.h"
  46 #include "hackblock.h"
  47 #include "string2.h"
  48 #include "gmx_fatal.h"
  49
  50 static void clear_atom_list(int i0, atom_id a[])
  51 {
  52   int i;
  53
  54   for (i=i0; i < MAXATOMLIST; i++)
  55     a[i]=-1;
  56 }
  57
  58 static void clear_force_param(int i0, real c[])
  59 {
  60   int i;
  61
  62   for (i=i0; i < MAXFORCEPARAM; i++)
  63     c[i]=NOTSET;
  64 }
  65
  66 void add_param(t_params *ps,int ai,int aj, real *c, char *s)
  67 {
  68   int i;
  69
  70   if ((ai < 0) || (aj < 0))
  71     gmx_fatal(FARGS,"Trying to add impossible atoms: ai=%d, aj=%d",ai,aj);
  72   pr_alloc(1,ps);
  73   ps->param[ps->nr].AI=ai;
  74   ps->param[ps->nr].AJ=aj;
  75   clear_atom_list(2, ps->param[ps->nr].a);
  76   if (c==NULL)
  77     clear_force_param(0, ps->param[ps->nr].c);
  78   else
  79     for(i=0; (i < MAXFORCEPARAM); i++)
  80       ps->param[ps->nr].c[i]=c[i];
  81   set_p_string(&(ps->param[ps->nr]),s);
  82   ps->nr++;
  83 }
  84
  85 void add_imp_param(t_params *ps,int ai,int aj,int ak,int al,real c0, real c1,
  86                    char *s)
  87 {
  88   pr_alloc(1,ps);
  89   ps->param[ps->nr].AI=ai;
  90   ps->param[ps->nr].AJ=aj;
  91   ps->param[ps->nr].AK=ak;
  92   ps->param[ps->nr].AL=al;
  93   clear_atom_list  (4, ps->param[ps->nr].a);
  94   ps->param[ps->nr].C0=c0;
  95   ps->param[ps->nr].C1=c1;
  96   clear_force_param(2, ps->param[ps->nr].c);
  97   set_p_string(&(ps->param[ps->nr]),s);
  98   ps->nr++;
  99 }
 100
 101 void add_dih_param(t_params *ps,int ai,int aj,int ak,int al,real c0, real c1,
 102                    real c2,char *s)
 103 {
 104   pr_alloc(1,ps);
 105   ps->param[ps->nr].AI=ai;
 106   ps->param[ps->nr].AJ=aj;
 107   ps->param[ps->nr].AK=ak;
 108   ps->param[ps->nr].AL=al;
 109   clear_atom_list  (4, ps->param[ps->nr].a);
 110   ps->param[ps->nr].C0=c0;
 111   ps->param[ps->nr].C1=c1;
 112   ps->param[ps->nr].C2=c2;
 113   clear_force_param(3, ps->param[ps->nr].c);
 114   set_p_string(&(ps->param[ps->nr]),s);
 115   ps->nr++;
 116 }
 117
 118 void add_cmap_param(t_params *ps, int ai, int aj, int ak, int al, int am, char *s)
 119 {
 120         pr_alloc(1,ps);
 121         ps->param[ps->nr].AI=ai;
 122         ps->param[ps->nr].AJ=aj;
 123         ps->param[ps->nr].AK=ak;
 124         ps->param[ps->nr].AL=al;
 125         ps->param[ps->nr].AM=am;
 126         clear_atom_list(5,ps->param[ps->nr].a);
 127         clear_force_param(0,ps->param[ps->nr].c);
 128         set_p_string(&(ps->param[ps->nr]),s);
 129         ps->nr++;
 130 }
 131
 132 void add_vsite2_atoms(t_params *ps,int ai,int aj,int ak)
 133 {
 134   pr_alloc(1,ps);
 135   ps->param[ps->nr].AI=ai;
 136   ps->param[ps->nr].AJ=aj;
 137   ps->param[ps->nr].AK=ak;
 138   clear_atom_list  (3, ps->param[ps->nr].a);
 139   clear_force_param(0, ps->param[ps->nr].c);
 140   set_p_string(&(ps->param[ps->nr]),"");
 141   ps->nr++;
 142 }
 143
 144 void add_vsite2_param(t_params *ps,int ai,int aj,int ak,real c0)
 145 {
 146   pr_alloc(1,ps);
 147   ps->param[ps->nr].AI=ai;
 148   ps->param[ps->nr].AJ=aj;
 149   ps->param[ps->nr].AK=ak;
 150   clear_atom_list  (3, ps->param[ps->nr].a);
 151   ps->param[ps->nr].C0=c0;
 152   clear_force_param(1, ps->param[ps->nr].c);
 153   set_p_string(&(ps->param[ps->nr]),"");
 154   ps->nr++;
 155 }
 156
 157 void add_vsite3_param(t_params *ps,int ai,int aj,int ak,int al,
 158                     real c0, real c1)
 159 {
 160   pr_alloc(1,ps);
 161   ps->param[ps->nr].AI=ai;
 162   ps->param[ps->nr].AJ=aj;
 163   ps->param[ps->nr].AK=ak;
 164   ps->param[ps->nr].AL=al;
 165   clear_atom_list  (4, ps->param[ps->nr].a);
 166   ps->param[ps->nr].C0=c0;
 167   ps->param[ps->nr].C1=c1;
 168   clear_force_param(2, ps->param[ps->nr].c);
 169   set_p_string(&(ps->param[ps->nr]),"");
 170   ps->nr++;
 171 }
 172
 173 void add_vsite3_atoms(t_params *ps,int ai,int aj,int ak,int al, gmx_bool bSwapParity)
 174 {
 175   pr_alloc(1,ps);
 176   ps->param[ps->nr].AI=ai;
 177   ps->param[ps->nr].AJ=aj;
 178   ps->param[ps->nr].AK=ak;
 179   ps->param[ps->nr].AL=al;
 180   clear_atom_list  (4, ps->param[ps->nr].a);
 181   clear_force_param(0, ps->param[ps->nr].c);
 182   if (bSwapParity)
 183     ps->param[ps->nr].C1=-1;
 184   set_p_string(&(ps->param[ps->nr]),"");
 185   ps->nr++;
 186 }
 187
 188 void add_vsite4_atoms(t_params *ps,int ai,int aj,int ak,int al,int am)
 189 {
 190   pr_alloc(1,ps);
 191   ps->param[ps->nr].AI=ai;
 192   ps->param[ps->nr].AJ=aj;
 193   ps->param[ps->nr].AK=ak;
 194   ps->param[ps->nr].AL=al;
 195   ps->param[ps->nr].AM=am;
 196   clear_atom_list  (5, ps->param[ps->nr].a);
 197   clear_force_param(0, ps->param[ps->nr].c);
 198   set_p_string(&(ps->param[ps->nr]),"");
 199   ps->nr++;
 200 }
 201
 202 int search_jtype(t_restp *rtp,char *name,gmx_bool bNterm)
 203 {
 204   int  niter,iter,j,k,kmax,jmax,minstrlen;
 205   char *rtpname,searchname[12];
 206
 207   strcpy(searchname,name);
 208
 209   /* Do a best match comparison */
 210   /* for protein N-terminus, allow renaming of H1, H2 and H3 to H */
 211   if ( bNterm && (strlen(searchname)==2) && (searchname[0] == 'H') &&
 212        ( (searchname[1] == '1') || (searchname[1] == '2') ||
 213          (searchname[1] == '3') ) ) {
 214     niter = 2;
 215   } else {
 216     niter = 1;
 217   }
 218   kmax=0;
 219   jmax=-1;
 220   for(iter=0; (iter<niter && jmax==-1); iter++) {
 221     if (iter == 1) {
 222       /* Try without the hydrogen number in the N-terminus */
 223       searchname[1] = '\0';
 224     }
 225     for(j=0; (j<rtp->natom); j++) {
 226       rtpname=*(rtp->atomname[j]);
 227       if (gmx_strcasecmp(searchname,rtpname) == 0) {
 228         jmax=j;
 229         kmax=strlen(searchname);
 230         break;
 231       }
 232       if (iter == niter - 1) {
 233         minstrlen = min(strlen(searchname), strlen(rtpname));
 234         for(k=0; k < minstrlen; k++)
 235           if (searchname[k] != rtpname[k])
 236             break;
 237         if (k > kmax) {
 238           kmax=k;
 239           jmax=j;
 240         }
 241       }
 242     }
 243   }
 244   if (jmax == -1)
 245     gmx_fatal(FARGS,"Atom %s not found in rtp database in residue %s",
 246               searchname,rtp->resname);
 247   if (kmax != strlen(searchname))
 248     gmx_fatal(FARGS,"Atom %s not found in rtp database in residue %s, "
 249               "it looks a bit like %s",
 250               searchname,rtp->resname,*(rtp->atomname[jmax]));
 251   return jmax;
 252 }
 253