mdtools/interface_implementation/DipoleTestFF.cpp

#include <cstdlib>
#include <cstdio>
#include <cstring>

#include <iostream>
using namespace std;

#include "ForceFields.hpp"
#include "SRI.hpp"


DipoleTestFF::DipoleTestFF(){

  char fileName[200];
  char* ffPath_env = "FORCE_PARAM_PATH";
  char* ffPath;
  char temp[200];

  // generate the force file name

  strcpy( fileName, "DipoleTest.frc" ); 

  // attempt to open the file in the current directory first.
  
  frcFile = fopen( fileName, "r" );
  
  if( frcFile == NULL ){

    // next see if the force path enviorment variable is set
    
    ffPath = getenv( ffPath_env );
    strcpy( temp, ffPath );
    strcat( temp, "/" );
    strcat( temp, fileName );
    strcpy( fileName, temp );

    frcFile = fopen( fileName, "r" );
    
    if( frcFile == NULL ){
      
      fprintf( stderr,
               "Error opening the force field parameter file: %s\n"
               "Have you tried setting the FORCE_PARAM_PATH environment "
               "vairable?\n",
               fileName );
      exit( 8 );
    }
  }
}

DipoleTestFF::~DipoleTestFF(){
  
  fclose( frcFile );
}

void DipoleTestFF::initializeAtoms( void ){
  
  class LinkedType {
  public:
    LinkedType(){ 
      next = NULL;
      name[0] = '\0';
    }
    ~LinkedType(){ if( next != NULL ) delete next; }

    LinkedType* find(char* key){ 
      if( !strcmp(name, key) ) return this;
      if( next != NULL ) return next->find(key);
      return NULL;
    }
    
    char name[15];
    int isDipole;
    double mass;
    double epslon;
    double sigma;
    double dipole;
    LinkedType* next;
  };
  
  LinkedType* headAtomType; 
  LinkedType* currentAtomType;
  LinkedType* tempAtomType;

  char readLine[500];
  char* the_token;
  char* eof_test;
  int foundAtom = 0;
  int lineNum = 0;
  int i;
  
  //////////////////////////////////////////////////
  // a quick water fix 

  double testI[3][3]; 
  testI[0][0] = 1.76958347772500;
  testI[0][1] = 0.0;
  testI[0][2] = 0.0;

  testI[1][0] = 0.0;
  testI[1][1] = 0.614537057924513;
  testI[1][2] = 0.0;

  testI[2][0] = 0.0;
  testI[2][1] = 0.0;
  testI[2][2] = 1.15504641980049;

  //////////////////////////////////////////////////

  Atom** the_atoms;
  int nAtoms;
  the_atoms = entry_plug->atoms;
  nAtoms = entry_plug->n_atoms;
  
  // read in the atom types.

  rewind( frcFile );
  headAtomType = new LinkedType;
  currentAtomType = headAtomType;
  
  eof_test = fgets( readLine, sizeof(readLine), frcFile );
  lineNum++;
  if( eof_test == NULL ){
    fprintf( stderr, "Error in reading Atoms from force file.\n" );
    exit(8);
  }

  
  while( !foundAtom ){
    while( eof_test != NULL && readLine[0] != '#' ){
      eof_test = fgets( readLine, sizeof(readLine), frcFile );
      lineNum++;
    }
    if( eof_test == NULL ){
      fprintf( stderr, 
               "Error in reading Atoms from force file at line %d.\n",
               lineNum );
      exit(8);
    }

    the_token = strtok( readLine, " ,;\t#\n" );
    foundAtom = !strcmp( "AtomTypes", the_token );
    
    if( !foundAtom ){
      eof_test = fgets( readLine, sizeof(readLine), frcFile );
      lineNum++;

      if( eof_test == NULL ){
        fprintf( stderr, 
                 "Error in reading Atoms from force file at line %d.\n",
                 lineNum );
        exit(8);
      } 
    }
  }

  // we are now at the AtomTypes section.

  eof_test = fgets( readLine, sizeof(readLine), frcFile );
  lineNum++;
  
  if( eof_test == NULL ){
    fprintf( stderr, 
             "Error in reading Atoms from force file at line %d.\n",
             lineNum );
    exit(8);
  }

  while( readLine[0] != '#' && eof_test != NULL ){
    
    if( readLine[0] != '!' ){
      
      the_token = strtok( readLine, " \n\t,;" );
      if( the_token != NULL ){
        
        strcpy( currentAtomType->name, the_token );

        if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
          fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
          exit(8);
        }

        sscanf( the_token, "%d", &currentAtomType->isDipole );

        if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
          fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
          exit(8);
        }

        sscanf( the_token, "%lf", &currentAtomType->mass );

        if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
          fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
          exit(8);
        }
        
        sscanf( the_token, "%lf", &currentAtomType->epslon );
        
        if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
          fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
          exit(8);
        }

        sscanf( the_token, "%lf", &currentAtomType->sigma );
        
        if( currentAtomType->isDipole ){
          if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
            fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
            exit(8);
          }
          
          sscanf( the_token, "%lf", &currentAtomType->dipole );
        }
      }
    }
    
    tempAtomType = new LinkedType;
    currentAtomType->next = tempAtomType;
    currentAtomType = tempAtomType;
    
    eof_test = fgets( readLine, sizeof(readLine), frcFile );
    lineNum++;
  }

  
  // initialize the atoms
  
  DirectionalAtom* dAtom;

  for( i=0; i<nAtoms; i++ ){
    
    currentAtomType = headAtomType->find( the_atoms[i]->getType() );
    if( currentAtomType == NULL ){
      fprintf( stderr, "AtomType error, %s not found in force file.\n",
               the_atoms[i]->getType() );
      exit(8);
    }
    
    the_atoms[i]->setMass( currentAtomType->mass );
    the_atoms[i]->setEpslon( currentAtomType->epslon );
    the_atoms[i]->setSigma( currentAtomType->sigma );
    the_atoms[i]->setLJ();

    if( currentAtomType->isDipole ){
      if( the_atoms[i]->isDirectional() ){
        dAtom = (DirectionalAtom *) the_atoms[i];
        dAtom->setMu( currentAtomType->dipole );
        dAtom->setHasDipole( 1 );
        dAtom->setSSD( 1 );
        dAtom->setJx( 0.0 );
        dAtom->setJy( 0.0 );
        dAtom->setJz( 0.0 );
        
        dAtom->setI( testI );

        entry_plug->n_dipoles++;
      }
      else{
        std::cerr 
          << "DipoleTestFF error: Atom \"" 
          << currentAtomType->name << "\" is a dipole, yet no standard"
          << " orientation was specifed in the BASS file.\n";
        exit(8);
      }
    }
    else{
      if( the_atoms[i]->isDirectional() ){
        std::cerr
          << "DipoleTestFF error: Atom \""
          << currentAtomType->name << "\" was given a standard orientation"
          << " in the BASS file, yet it is not a dipole.\n";
        exit(8);
      }
    }
  }


  // clean up the memory
  
  delete headAtomType;
}

void DipoleTestFF::initializeBonds( bond_pair* the_bonds ){
  
  if( entry_plug->n_bonds ){
    std::cerr << "DipoleTest does not support bonds!\n";
    exit(8);
  }
}

void DipoleTestFF::initializeBends( bend_set* the_bends ){
  
  if( entry_plug->n_bends ){
    std::cerr << "DipoleTest does not support bends!\n";
    exit(8);
  }
}

void DipoleTestFF::initializeTorsions( torsion_set* the_torsions ){

  if( entry_plug->n_torsions ){  
    std::cerr << "DipoleTest does not support torsions!\n";
    exit(8);
  }
}
Revision:	11
Committed:	Tue Jul 9 18:40:59 2002 UTC (22 years ago) by mmeineke
File size:	6869 byte(s)
Log Message:	This commit was generated by cvs2svn to compensate for changes in r10, which included commits to RCS files with non-trunk default branches.
#	User	Rev	Content
1	mmeineke	10	#include <cstdlib>
2			#include <cstdio>
3			#include <cstring>
4
5			#include <iostream>
6			using namespace std;
7
8			#include "ForceFields.hpp"
9			#include "SRI.hpp"
10
11
12			DipoleTestFF::DipoleTestFF(){
13
14			char fileName[200];
15			char* ffPath_env = "FORCE_PARAM_PATH";
16			char* ffPath;
17			char temp[200];
18
19			// generate the force file name
20
21			strcpy( fileName, "DipoleTest.frc" );
22
23			// attempt to open the file in the current directory first.
24
25			frcFile = fopen( fileName, "r" );
26
27			if( frcFile == NULL ){
28
29			// next see if the force path enviorment variable is set
30
31			ffPath = getenv( ffPath_env );
32			strcpy( temp, ffPath );
33			strcat( temp, "/" );
34			strcat( temp, fileName );
35			strcpy( fileName, temp );
36
37			frcFile = fopen( fileName, "r" );
38
39			if( frcFile == NULL ){
40
41			fprintf( stderr,
42			"Error opening the force field parameter file: %s\n"
43			"Have you tried setting the FORCE_PARAM_PATH environment "
44			"vairable?\n",
45			fileName );
46			exit( 8 );
47			}
48			}
49			}
50
51			DipoleTestFF::~DipoleTestFF(){
52
53			fclose( frcFile );
54			}
55
56			void DipoleTestFF::initializeAtoms( void ){
57
58			class LinkedType {
59			public:
60			LinkedType(){
61			next = NULL;
62			name[0] = '\0';
63			}
64			~LinkedType(){ if( next != NULL ) delete next; }
65
66			LinkedType* find(char* key){
67			if( !strcmp(name, key) ) return this;
68			if( next != NULL ) return next->find(key);
69			return NULL;
70			}
71
72			char name[15];
73			int isDipole;
74			double mass;
75			double epslon;
76			double sigma;
77			double dipole;
78			LinkedType* next;
79			};
80
81			LinkedType* headAtomType;
82			LinkedType* currentAtomType;
83			LinkedType* tempAtomType;
84
85			char readLine[500];
86			char* the_token;
87			char* eof_test;
88			int foundAtom = 0;
89			int lineNum = 0;
90			int i;
91
92			//////////////////////////////////////////////////
93			// a quick water fix
94
95			double testI[3][3];
96			testI[0][0] = 1.76958347772500;
97			testI[0][1] = 0.0;
98			testI[0][2] = 0.0;
99
100			testI[1][0] = 0.0;
101			testI[1][1] = 0.614537057924513;
102			testI[1][2] = 0.0;
103
104			testI[2][0] = 0.0;
105			testI[2][1] = 0.0;
106			testI[2][2] = 1.15504641980049;
107
108			//////////////////////////////////////////////////
109
110			Atom** the_atoms;
111			int nAtoms;
112			the_atoms = entry_plug->atoms;
113			nAtoms = entry_plug->n_atoms;
114
115			// read in the atom types.
116
117			rewind( frcFile );
118			headAtomType = new LinkedType;
119			currentAtomType = headAtomType;
120
121			eof_test = fgets( readLine, sizeof(readLine), frcFile );
122			lineNum++;
123			if( eof_test == NULL ){
124			fprintf( stderr, "Error in reading Atoms from force file.\n" );
125			exit(8);
126			}
127
128
129			while( !foundAtom ){
130			while( eof_test != NULL && readLine[0] != '#' ){
131			eof_test = fgets( readLine, sizeof(readLine), frcFile );
132			lineNum++;
133			}
134			if( eof_test == NULL ){
135			fprintf( stderr,
136			"Error in reading Atoms from force file at line %d.\n",
137			lineNum );
138			exit(8);
139			}
140
141			the_token = strtok( readLine, " ,;\t#\n" );
142			foundAtom = !strcmp( "AtomTypes", the_token );
143
144			if( !foundAtom ){
145			eof_test = fgets( readLine, sizeof(readLine), frcFile );
146			lineNum++;
147
148			if( eof_test == NULL ){
149			fprintf( stderr,
150			"Error in reading Atoms from force file at line %d.\n",
151			lineNum );
152			exit(8);
153			}
154			}
155			}
156
157			// we are now at the AtomTypes section.
158
159			eof_test = fgets( readLine, sizeof(readLine), frcFile );
160			lineNum++;
161
162			if( eof_test == NULL ){
163			fprintf( stderr,
164			"Error in reading Atoms from force file at line %d.\n",
165			lineNum );
166			exit(8);
167			}
168
169			while( readLine[0] != '#' && eof_test != NULL ){
170
171			if( readLine[0] != '!' ){
172
173			the_token = strtok( readLine, " \n\t,;" );
174			if( the_token != NULL ){
175
176			strcpy( currentAtomType->name, the_token );
177
178			if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
179			fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
180			exit(8);
181			}
182
183			sscanf( the_token, "%d", &currentAtomType->isDipole );
184
185			if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
186			fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
187			exit(8);
188			}
189
190			sscanf( the_token, "%lf", &currentAtomType->mass );
191
192			if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
193			fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
194			exit(8);
195			}
196
197			sscanf( the_token, "%lf", &currentAtomType->epslon );
198
199			if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
200			fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
201			exit(8);
202			}
203
204			sscanf( the_token, "%lf", &currentAtomType->sigma );
205
206			if( currentAtomType->isDipole ){
207			if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
208			fprintf( stderr, "Error parseing AtomTypes: line %d\n", lineNum );
209			exit(8);
210			}
211
212			sscanf( the_token, "%lf", &currentAtomType->dipole );
213			}
214			}
215			}
216
217			tempAtomType = new LinkedType;
218			currentAtomType->next = tempAtomType;
219			currentAtomType = tempAtomType;
220
221			eof_test = fgets( readLine, sizeof(readLine), frcFile );
222			lineNum++;
223			}
224
225
226			// initialize the atoms
227
228			DirectionalAtom* dAtom;
229
230			for( i=0; i<nAtoms; i++ ){
231
232			currentAtomType = headAtomType->find( the_atoms[i]->getType() );
233			if( currentAtomType == NULL ){
234			fprintf( stderr, "AtomType error, %s not found in force file.\n",
235			the_atoms[i]->getType() );
236			exit(8);
237			}
238
239			the_atoms[i]->setMass( currentAtomType->mass );
240			the_atoms[i]->setEpslon( currentAtomType->epslon );
241			the_atoms[i]->setSigma( currentAtomType->sigma );
242			the_atoms[i]->setLJ();
243
244			if( currentAtomType->isDipole ){
245			if( the_atoms[i]->isDirectional() ){
246			dAtom = (DirectionalAtom *) the_atoms[i];
247			dAtom->setMu( currentAtomType->dipole );
248			dAtom->setHasDipole( 1 );
249			dAtom->setSSD( 1 );
250			dAtom->setJx( 0.0 );
251			dAtom->setJy( 0.0 );
252			dAtom->setJz( 0.0 );
253
254			dAtom->setI( testI );
255
256			entry_plug->n_dipoles++;
257			}
258			else{
259			std::cerr
260			<< "DipoleTestFF error: Atom \""
261			<< currentAtomType->name << "\" is a dipole, yet no standard"
262			<< " orientation was specifed in the BASS file.\n";
263			exit(8);
264			}
265			}
266			else{
267			if( the_atoms[i]->isDirectional() ){
268			std::cerr
269			<< "DipoleTestFF error: Atom \""
270			<< currentAtomType->name << "\" was given a standard orientation"
271			<< " in the BASS file, yet it is not a dipole.\n";
272			exit(8);
273			}
274			}
275			}
276
277
278			// clean up the memory
279
280			delete headAtomType;
281			}
282
283			void DipoleTestFF::initializeBonds( bond_pair* the_bonds ){
284
285			if( entry_plug->n_bonds ){
286			std::cerr << "DipoleTest does not support bonds!\n";
287			exit(8);
288			}
289			}
290
291			void DipoleTestFF::initializeBends( bend_set* the_bends ){
292
293			if( entry_plug->n_bends ){
294			std::cerr << "DipoleTest does not support bends!\n";
295			exit(8);
296			}
297			}
298
299			void DipoleTestFF::initializeTorsions( torsion_set* the_torsions ){
300
301			if( entry_plug->n_torsions ){
302			std::cerr << "DipoleTest does not support torsions!\n";
303			exit(8);
304			}
305			}