ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 493
Committed: Mon Apr 14 19:03:38 2003 UTC (21 years, 2 months ago) by mmeineke
File size: 40525 byte(s)
Log Message:
added Ghost bends to the TraPPE_Ex forceField

File Contents

# Content
1 #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 #include "simError.h"
11
12 #include "fortranWrappers.hpp"
13
14 #ifdef IS_MPI
15 #include "mpiForceField.h"
16 #endif // is_mpi
17
18 namespace TPE { // restrict the access of the folowing to this file only.
19
20
21 // Declare the structures that will be passed by MPI
22
23 typedef struct{
24 char name[15];
25 double mass;
26 double epslon;
27 double sigma;
28 double dipole;
29 double w0;
30 double v0;
31 int isSSD;
32 int isDipole;
33 int ident;
34 int last; // 0 -> default
35 // 1 -> tells nodes to stop listening
36 } atomStruct;
37
38
39 typedef struct{
40 char nameA[15];
41 char nameB[15];
42 char type[30];
43 double d0;
44 int last; // 0 -> default
45 // 1 -> tells nodes to stop listening
46 } bondStruct;
47
48
49 typedef struct{
50 char nameA[15];
51 char nameB[15];
52 char nameC[15];
53 char type[30];
54 double k1, k2, k3, t0;
55 int last; // 0 -> default
56 // 1 -> tells nodes to stop listening
57 } bendStruct;
58
59
60 typedef struct{
61 char nameA[15];
62 char nameB[15];
63 char nameC[15];
64 char nameD[15];
65 char type[30];
66 double k1, k2, k3, k4;
67 int last; // 0 -> default
68 // 1 -> tells nodes to stop listening
69 } torsionStruct;
70
71
72 int parseAtom( char *lineBuffer, int lineNum, atomStruct &info );
73 int parseBond( char *lineBuffer, int lineNum, bondStruct &info );
74 int parseBend( char *lineBuffer, int lineNum, bendStruct &info );
75 int parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info );
76
77
78 #ifdef IS_MPI
79
80 MPI_Datatype mpiAtomStructType;
81 MPI_Datatype mpiBondStructType;
82 MPI_Datatype mpiBendStructType;
83 MPI_Datatype mpiTorsionStructType;
84
85 #endif
86
87 class LinkedAtomType {
88 public:
89 LinkedAtomType(){
90 next = NULL;
91 name[0] = '\0';
92 }
93 ~LinkedAtomType(){ if( next != NULL ) delete next; }
94
95 LinkedAtomType* find(char* key){
96 if( !strcmp(name, key) ) return this;
97 if( next != NULL ) return next->find(key);
98 return NULL;
99 }
100
101 void printMe( void ){
102
103 std::cerr << "LinkedAtype " << name << ": ident = " << ident << "\n";
104 if( next != NULL ) next->printMe();
105
106 }
107
108 void add( atomStruct &info ){
109
110 // check for duplicates
111
112 if( !strcmp( info.name, name ) ){
113 sprintf( painCave.errMsg,
114 "Duplicate TraPPE_Ex atom type \"%s\" found in "
115 "the TraPPE_ExFF param file./n",
116 name );
117 painCave.isFatal = 1;
118 simError();
119 }
120
121 if( next != NULL ) next->add(info);
122 else{
123 next = new LinkedAtomType();
124 strcpy(next->name, info.name);
125 next->isDipole = info.isDipole;
126 next->isSSD = info.isSSD;
127 next->mass = info.mass;
128 next->epslon = info.epslon;
129 next->sigma = info.sigma;
130 next->dipole = info.dipole;
131 next->w0 = info.w0;
132 next->v0 = info.v0;
133 next->ident = info.ident;
134 }
135 }
136
137 #ifdef IS_MPI
138
139 void duplicate( atomStruct &info ){
140 strcpy(info.name, name);
141 info.isDipole = isDipole;
142 info.isSSD = isSSD;
143 info.mass = mass;
144 info.epslon = epslon;
145 info.sigma = sigma;
146 info.dipole = dipole;
147 info.w0 = w0;
148 info.v0 = v0;
149 info.ident = ident;
150 info.last = 0;
151 }
152
153
154 #endif
155
156 char name[15];
157 int isDipole;
158 int isSSD;
159 double mass;
160 double epslon;
161 double sigma;
162 double dipole;
163 double w0;
164 double v0;
165 int ident;
166 LinkedAtomType* next;
167 };
168
169 class LinkedBondType {
170 public:
171 LinkedBondType(){
172 next = NULL;
173 nameA[0] = '\0';
174 nameB[0] = '\0';
175 type[0] = '\0';
176 }
177 ~LinkedBondType(){ if( next != NULL ) delete next; }
178
179 LinkedBondType* find(char* key1, char* key2){
180 if( !strcmp(nameA, key1 ) && !strcmp( nameB, key2 ) ) return this;
181 if( !strcmp(nameA, key2 ) && !strcmp( nameB, key1 ) ) return this;
182 if( next != NULL ) return next->find(key1, key2);
183 return NULL;
184 }
185
186
187 void add( bondStruct &info ){
188
189 // check for duplicates
190 int dup = 0;
191
192 if( !strcmp(nameA, info.nameA ) && !strcmp( nameB, info.nameB ) ) dup = 1;
193 if( !strcmp(nameA, info.nameB ) && !strcmp( nameB, info.nameA ) ) dup = 1;
194
195 if(dup){
196 sprintf( painCave.errMsg,
197 "Duplicate TraPPE_Ex bond type \"%s - %s\" found in "
198 "the TraPPE_ExFF param file./n",
199 nameA, nameB );
200 painCave.isFatal = 1;
201 simError();
202 }
203
204
205 if( next != NULL ) next->add(info);
206 else{
207 next = new LinkedBondType();
208 strcpy(next->nameA, info.nameA);
209 strcpy(next->nameB, info.nameB);
210 strcpy(next->type, info.type);
211 next->d0 = info.d0;
212 }
213 }
214
215 #ifdef IS_MPI
216 void duplicate( bondStruct &info ){
217 strcpy(info.nameA, nameA);
218 strcpy(info.nameB, nameB);
219 strcpy(info.type, type);
220 info.d0 = d0;
221 info.last = 0;
222 }
223
224
225 #endif
226
227 char nameA[15];
228 char nameB[15];
229 char type[30];
230 double d0;
231
232 LinkedBondType* next;
233 };
234
235 class LinkedBendType {
236 public:
237 LinkedBendType(){
238 next = NULL;
239 nameA[0] = '\0';
240 nameB[0] = '\0';
241 nameC[0] = '\0';
242 type[0] = '\0';
243 }
244 ~LinkedBendType(){ if( next != NULL ) delete next; }
245
246 LinkedBendType* find( char* key1, char* key2, char* key3 ){
247 if( !strcmp( nameA, key1 ) && !strcmp( nameB, key2 )
248 && !strcmp( nameC, key3 ) ) return this;
249 if( !strcmp( nameA, key3 ) && !strcmp( nameB, key2 )
250 && !strcmp( nameC, key1 ) ) return this;
251 if( next != NULL ) return next->find(key1, key2, key3);
252 return NULL;
253 }
254
255 void add( bendStruct &info ){
256
257 // check for duplicates
258 int dup = 0;
259
260 if( !strcmp( nameA, info.nameA ) && !strcmp( nameB, info.nameB )
261 && !strcmp( nameC, info.nameC ) ) dup = 1;
262 if( !strcmp( nameA, info.nameC ) && !strcmp( nameB, info.nameB )
263 && !strcmp( nameC, info.nameA ) ) dup = 1;
264
265 if(dup){
266 sprintf( painCave.errMsg,
267 "Duplicate TraPPE_Ex bend type \"%s - %s - %s\" found in "
268 "the TraPPE_ExFF param file./n",
269 nameA, nameB, nameC );
270 painCave.isFatal = 1;
271 simError();
272 }
273
274 if( next != NULL ) next->add(info);
275 else{
276 next = new LinkedBendType();
277 strcpy(next->nameA, info.nameA);
278 strcpy(next->nameB, info.nameB);
279 strcpy(next->nameC, info.nameC);
280 strcpy(next->type, info.type);
281 next->k1 = info.k1;
282 next->k2 = info.k2;
283 next->k3 = info.k3;
284 next->t0 = info.t0;
285 }
286 }
287
288 #ifdef IS_MPI
289
290 void duplicate( bendStruct &info ){
291 strcpy(info.nameA, nameA);
292 strcpy(info.nameB, nameB);
293 strcpy(info.nameC, nameC);
294 strcpy(info.type, type);
295 info.k1 = k1;
296 info.k2 = k2;
297 info.k3 = k3;
298 info.t0 = t0;
299 info.last = 0;
300 }
301
302 #endif // is_mpi
303
304 char nameA[15];
305 char nameB[15];
306 char nameC[15];
307 char type[30];
308 double k1, k2, k3, t0;
309
310 LinkedBendType* next;
311 };
312
313 class LinkedTorsionType {
314 public:
315 LinkedTorsionType(){
316 next = NULL;
317 nameA[0] = '\0';
318 nameB[0] = '\0';
319 nameC[0] = '\0';
320 type[0] = '\0';
321 }
322 ~LinkedTorsionType(){ if( next != NULL ) delete next; }
323
324 LinkedTorsionType* find( char* key1, char* key2, char* key3, char* key4 ){
325
326
327
328
329 if( !strcmp( nameA, key1 ) && !strcmp( nameB, key2 ) &&
330 !strcmp( nameC, key3 ) && !strcmp( nameD, key4 ) ) return this;
331
332 if( !strcmp( nameA, key4 ) && !strcmp( nameB, key3 ) &&
333 !strcmp( nameC, key2 ) && !strcmp( nameD, key1 ) ) return this;
334
335 if( next != NULL ) return next->find(key1, key2, key3, key4);
336 return NULL;
337 }
338
339 void add( torsionStruct &info ){
340
341 // check for duplicates
342 int dup = 0;
343
344 if( !strcmp( nameA, info.nameA ) && !strcmp( nameB, info.nameB ) &&
345 !strcmp( nameC, info.nameC ) && !strcmp( nameD, info.nameD ) ) dup = 1;
346
347 if( !strcmp( nameA, info.nameD ) && !strcmp( nameB, info.nameC ) &&
348 !strcmp( nameC, info.nameB ) && !strcmp( nameD, info.nameA ) ) dup = 1;
349
350 if(dup){
351 sprintf( painCave.errMsg,
352 "Duplicate TraPPE_Ex torsion type \"%s - %s - %s - %s\" found in "
353 "the TraPPE_ExFF param file./n", nameA, nameB, nameC, nameD );
354 painCave.isFatal = 1;
355 simError();
356 }
357
358 if( next != NULL ) next->add(info);
359 else{
360 next = new LinkedTorsionType();
361 strcpy(next->nameA, info.nameA);
362 strcpy(next->nameB, info.nameB);
363 strcpy(next->nameC, info.nameC);
364 strcpy(next->nameD, info.nameD);
365 strcpy(next->type, info.type);
366 next->k1 = info.k1;
367 next->k2 = info.k2;
368 next->k3 = info.k3;
369 next->k4 = info.k4;
370
371 }
372 }
373
374 #ifdef IS_MPI
375
376 void duplicate( torsionStruct &info ){
377 strcpy(info.nameA, nameA);
378 strcpy(info.nameB, nameB);
379 strcpy(info.nameC, nameC);
380 strcpy(info.nameD, nameD);
381 strcpy(info.type, type);
382 info.k1 = k1;
383 info.k2 = k2;
384 info.k3 = k3;
385 info.k4 = k4;
386 info.last = 0;
387 }
388
389 #endif
390
391 char nameA[15];
392 char nameB[15];
393 char nameC[15];
394 char nameD[15];
395 char type[30];
396 double k1, k2, k3, k4;
397
398 LinkedTorsionType* next;
399 };
400
401
402 LinkedAtomType* headAtomType;
403 LinkedAtomType* currentAtomType;
404 LinkedBondType* headBondType;
405 LinkedBondType* currentBondType;
406 LinkedBendType* headBendType;
407 LinkedBendType* currentBendType;
408 LinkedTorsionType* headTorsionType;
409 LinkedTorsionType* currentTorsionType;
410
411 } // namespace
412
413 using namespace TPE;
414
415
416 //****************************************************************
417 // begins the actual forcefield stuff.
418 //****************************************************************
419
420
421 TraPPE_ExFF::TraPPE_ExFF(){
422
423 char fileName[200];
424 char* ffPath_env = "FORCE_PARAM_PATH";
425 char* ffPath;
426 char temp[200];
427 char errMsg[1000];
428
429 headAtomType = NULL;
430 currentAtomType = NULL;
431 headBondType = NULL;
432 currentBondType = NULL;
433 headBendType = NULL;
434 currentBendType = NULL;
435 headTorsionType = NULL;
436 currentTorsionType = NULL;
437
438 // do the funtion wrapping
439 wrapMeFF( this );
440
441
442 #ifdef IS_MPI
443 int i;
444
445 // **********************************************************************
446 // Init the atomStruct mpi type
447
448 atomStruct atomProto; // mpiPrototype
449 int atomBC[3] = {15,6,4}; // block counts
450 MPI_Aint atomDspls[3]; // displacements
451 MPI_Datatype atomMbrTypes[3]; // member mpi types
452
453 MPI_Address(&atomProto.name, &atomDspls[0]);
454 MPI_Address(&atomProto.mass, &atomDspls[1]);
455 MPI_Address(&atomProto.isSSD, &atomDspls[2]);
456
457 atomMbrTypes[0] = MPI_CHAR;
458 atomMbrTypes[1] = MPI_DOUBLE;
459 atomMbrTypes[2] = MPI_INT;
460
461 for (i=2; i >= 0; i--) atomDspls[i] -= atomDspls[0];
462
463 MPI_Type_struct(3, atomBC, atomDspls, atomMbrTypes, &mpiAtomStructType);
464 MPI_Type_commit(&mpiAtomStructType);
465
466
467 // **********************************************************************
468 // Init the bondStruct mpi type
469
470 bondStruct bondProto; // mpiPrototype
471 int bondBC[3] = {60,1,1}; // block counts
472 MPI_Aint bondDspls[3]; // displacements
473 MPI_Datatype bondMbrTypes[3]; // member mpi types
474
475 MPI_Address(&bondProto.nameA, &bondDspls[0]);
476 MPI_Address(&bondProto.d0, &bondDspls[1]);
477 MPI_Address(&bondProto.last, &bondDspls[2]);
478
479 bondMbrTypes[0] = MPI_CHAR;
480 bondMbrTypes[1] = MPI_DOUBLE;
481 bondMbrTypes[2] = MPI_INT;
482
483 for (i=2; i >= 0; i--) bondDspls[i] -= bondDspls[0];
484
485 MPI_Type_struct(3, bondBC, bondDspls, bondMbrTypes, &mpiBondStructType);
486 MPI_Type_commit(&mpiBondStructType);
487
488
489 // **********************************************************************
490 // Init the bendStruct mpi type
491
492 bendStruct bendProto; // mpiPrototype
493 int bendBC[3] = {75,4,1}; // block counts
494 MPI_Aint bendDspls[3]; // displacements
495 MPI_Datatype bendMbrTypes[3]; // member mpi types
496
497 MPI_Address(&bendProto.nameA, &bendDspls[0]);
498 MPI_Address(&bendProto.k1, &bendDspls[1]);
499 MPI_Address(&bendProto.last, &bendDspls[2]);
500
501 bendMbrTypes[0] = MPI_CHAR;
502 bendMbrTypes[1] = MPI_DOUBLE;
503 bendMbrTypes[2] = MPI_INT;
504
505 for (i=2; i >= 0; i--) bendDspls[i] -= bendDspls[0];
506
507 MPI_Type_struct(3, bendBC, bendDspls, bendMbrTypes, &mpiBendStructType);
508 MPI_Type_commit(&mpiBendStructType);
509
510
511 // **********************************************************************
512 // Init the torsionStruct mpi type
513
514 torsionStruct torsionProto; // mpiPrototype
515 int torsionBC[3] = {90,4,1}; // block counts
516 MPI_Aint torsionDspls[3]; // displacements
517 MPI_Datatype torsionMbrTypes[3]; // member mpi types
518
519 MPI_Address(&torsionProto.nameA, &torsionDspls[0]);
520 MPI_Address(&torsionProto.k1, &torsionDspls[1]);
521 MPI_Address(&torsionProto.last, &torsionDspls[2]);
522
523 torsionMbrTypes[0] = MPI_CHAR;
524 torsionMbrTypes[1] = MPI_DOUBLE;
525 torsionMbrTypes[2] = MPI_INT;
526
527 for (i=2; i >= 0; i--) torsionDspls[i] -= torsionDspls[0];
528
529 MPI_Type_struct(3, torsionBC, torsionDspls, torsionMbrTypes,
530 &mpiTorsionStructType);
531 MPI_Type_commit(&mpiTorsionStructType);
532
533 // ***********************************************************************
534
535 if( worldRank == 0 ){
536 #endif
537
538 // generate the force file name
539
540 strcpy( fileName, "TraPPE_Ex.frc" );
541 // fprintf( stderr,"Trying to open %s\n", fileName );
542
543 // attempt to open the file in the current directory first.
544
545 frcFile = fopen( fileName, "r" );
546
547 if( frcFile == NULL ){
548
549 // next see if the force path enviorment variable is set
550
551 ffPath = getenv( ffPath_env );
552 if( ffPath == NULL ) {
553 STR_DEFINE(ffPath, FRC_PATH );
554 }
555
556
557 strcpy( temp, ffPath );
558 strcat( temp, "/" );
559 strcat( temp, fileName );
560 strcpy( fileName, temp );
561
562 frcFile = fopen( fileName, "r" );
563
564 if( frcFile == NULL ){
565
566 sprintf( painCave.errMsg,
567 "Error opening the force field parameter file: %s\n"
568 "Have you tried setting the FORCE_PARAM_PATH environment "
569 "vairable?\n",
570 fileName );
571 painCave.isFatal = 1;
572 simError();
573 }
574 }
575
576 #ifdef IS_MPI
577 }
578
579 sprintf( checkPointMsg, "TraPPE_ExFF file opened sucessfully." );
580 MPIcheckPoint();
581
582 #endif // is_mpi
583 }
584
585
586 TraPPE_ExFF::~TraPPE_ExFF(){
587
588 if( headAtomType != NULL ) delete headAtomType;
589 if( headBondType != NULL ) delete headBondType;
590 if( headBendType != NULL ) delete headBendType;
591 if( headTorsionType != NULL ) delete headTorsionType;
592
593 #ifdef IS_MPI
594 if( worldRank == 0 ){
595 #endif // is_mpi
596
597 fclose( frcFile );
598
599 #ifdef IS_MPI
600 }
601 #endif // is_mpi
602 }
603
604 void TraPPE_ExFF::cleanMe( void ){
605
606 #ifdef IS_MPI
607
608 // keep the linked lists in the mpi version
609
610 #else // is_mpi
611
612 // delete the linked lists in the single processor version
613
614 if( headAtomType != NULL ) delete headAtomType;
615 if( headBondType != NULL ) delete headBondType;
616 if( headBendType != NULL ) delete headBendType;
617 if( headTorsionType != NULL ) delete headTorsionType;
618
619 #endif // is_mpi
620 }
621
622
623 void TraPPE_ExFF::initForceField( int ljMixRule ){
624
625 initFortran( ljMixRule, entry_plug->useReactionField );
626 }
627
628
629 void TraPPE_ExFF::readParams( void ){
630
631 int i, a, b, c, d;
632 int identNum;
633 char* atomA;
634 char* atomB;
635 char* atomC;
636 char* atomD;
637
638 atomStruct atomInfo;
639 bondStruct bondInfo;
640 bendStruct bendInfo;
641 torsionStruct torsionInfo;
642
643 bigSigma = 0.0;
644
645 atomInfo.last = 1;
646 bondInfo.last = 1;
647 bendInfo.last = 1;
648 torsionInfo.last = 1;
649
650 // read in the atom info
651
652 #ifdef IS_MPI
653 if( worldRank == 0 ){
654 #endif
655
656 // read in the atom types.
657
658 headAtomType = new LinkedAtomType;
659
660 fastForward( "AtomTypes", "initializeAtoms" );
661
662 // we are now at the AtomTypes section.
663
664 eof_test = fgets( readLine, sizeof(readLine), frcFile );
665 lineNum++;
666
667
668 // read a line, and start parseing out the atom types
669
670 if( eof_test == NULL ){
671 sprintf( painCave.errMsg,
672 "Error in reading Atoms from force file at line %d.\n",
673 lineNum );
674 painCave.isFatal = 1;
675 simError();
676 }
677
678 identNum = 1;
679 // stop reading at end of file, or at next section
680 while( readLine[0] != '#' && eof_test != NULL ){
681
682 // toss comment lines
683 if( readLine[0] != '!' ){
684
685 // the parser returns 0 if the line was blank
686 if( parseAtom( readLine, lineNum, atomInfo ) ){
687 atomInfo.ident = identNum;
688 headAtomType->add( atomInfo );;
689 identNum++;
690 }
691 }
692 eof_test = fgets( readLine, sizeof(readLine), frcFile );
693 lineNum++;
694 }
695
696 #ifdef IS_MPI
697
698 // send out the linked list to all the other processes
699
700 sprintf( checkPointMsg,
701 "TraPPE_ExFF atom structures read successfully." );
702 MPIcheckPoint();
703
704 currentAtomType = headAtomType->next; //skip the first element who is a place holder.
705 while( currentAtomType != NULL ){
706 currentAtomType->duplicate( atomInfo );
707
708
709
710 sendFrcStruct( &atomInfo, mpiAtomStructType );
711
712 sprintf( checkPointMsg,
713 "successfully sent TraPPE_Ex force type: \"%s\"\n",
714 atomInfo.name );
715 MPIcheckPoint();
716
717 currentAtomType = currentAtomType->next;
718 }
719 atomInfo.last = 1;
720 sendFrcStruct( &atomInfo, mpiAtomStructType );
721
722 }
723
724 else{
725
726 // listen for node 0 to send out the force params
727
728 MPIcheckPoint();
729
730 headAtomType = new LinkedAtomType;
731 recieveFrcStruct( &atomInfo, mpiAtomStructType );
732
733 while( !atomInfo.last ){
734
735
736
737 headAtomType->add( atomInfo );
738
739 MPIcheckPoint();
740
741 recieveFrcStruct( &atomInfo, mpiAtomStructType );
742 }
743 }
744
745 #endif // is_mpi
746
747
748
749 // call new A_types in fortran
750
751 int isError;
752
753 // dummy variables
754
755 int isGB = 0;
756 int isLJ = 1;
757 double GB_dummy = 0.0;
758
759
760 currentAtomType = headAtomType->next;;
761 while( currentAtomType != NULL ){
762
763 if(currentAtomType->isDipole) entry_plug->useDipole = 1;
764 if(currentAtomType->isSSD) {
765 entry_plug->useSticky = 1;
766 set_sticky_params( &(currentAtomType->w0), &(currentAtomType->v0));
767 }
768
769 if( currentAtomType->name[0] != '\0' ){
770 isError = 0;
771 makeAtype( &(currentAtomType->ident),
772 &isLJ,
773 &(currentAtomType->isSSD),
774 &(currentAtomType->isDipole),
775 &isGB,
776 &(currentAtomType->epslon),
777 &(currentAtomType->sigma),
778 &(currentAtomType->dipole),
779 &isError );
780 if( isError ){
781 sprintf( painCave.errMsg,
782 "Error initializing the \"%s\" atom type in fortran\n",
783 currentAtomType->name );
784 painCave.isFatal = 1;
785 simError();
786 }
787 }
788 currentAtomType = currentAtomType->next;
789 }
790
791 #ifdef IS_MPI
792 sprintf( checkPointMsg,
793 "TraPPE_ExFF atom structures successfully sent to fortran\n" );
794 MPIcheckPoint();
795 #endif // is_mpi
796
797
798
799 // read in the bonds
800
801 #ifdef IS_MPI
802 if( worldRank == 0 ){
803 #endif
804
805 // read in the bond types.
806
807 headBondType = new LinkedBondType;
808
809 fastForward( "BondTypes", "initializeBonds" );
810
811 // we are now at the bondTypes section
812
813 eof_test = fgets( readLine, sizeof(readLine), frcFile );
814 lineNum++;
815
816
817 // read a line, and start parseing out the atom types
818
819 if( eof_test == NULL ){
820 sprintf( painCave.errMsg,
821 "Error in reading bonds from force file at line %d.\n",
822 lineNum );
823 painCave.isFatal = 1;
824 simError();
825 }
826
827 // stop reading at end of file, or at next section
828 while( readLine[0] != '#' && eof_test != NULL ){
829
830 // toss comment lines
831 if( readLine[0] != '!' ){
832
833 // the parser returns 0 if the line was blank
834 if( parseBond( readLine, lineNum, bondInfo ) ){
835 headBondType->add( bondInfo );
836 }
837 }
838 eof_test = fgets( readLine, sizeof(readLine), frcFile );
839 lineNum++;
840 }
841
842 #ifdef IS_MPI
843
844 // send out the linked list to all the other processes
845
846 sprintf( checkPointMsg,
847 "TraPPE_Ex bond structures read successfully." );
848 MPIcheckPoint();
849
850 currentBondType = headBondType->next;
851 while( currentBondType != NULL ){
852 currentBondType->duplicate( bondInfo );
853 sendFrcStruct( &bondInfo, mpiBondStructType );
854 currentBondType = currentBondType->next;
855 }
856 bondInfo.last = 1;
857 sendFrcStruct( &bondInfo, mpiBondStructType );
858
859 }
860
861 else{
862
863 // listen for node 0 to send out the force params
864
865 MPIcheckPoint();
866
867 headBondType = new LinkedBondType;
868 recieveFrcStruct( &bondInfo, mpiBondStructType );
869 while( !bondInfo.last ){
870
871 headBondType->add( bondInfo );
872 recieveFrcStruct( &bondInfo, mpiBondStructType );
873 }
874 }
875 #endif // is_mpi
876
877
878 // read in the bends
879
880 #ifdef IS_MPI
881 if( worldRank == 0 ){
882 #endif
883
884 // read in the bend types.
885
886 headBendType = new LinkedBendType;
887
888 fastForward( "BendTypes", "initializeBends" );
889
890 // we are now at the bendTypes section
891
892 eof_test = fgets( readLine, sizeof(readLine), frcFile );
893 lineNum++;
894
895 // read a line, and start parseing out the bend types
896
897 if( eof_test == NULL ){
898 sprintf( painCave.errMsg,
899 "Error in reading bends from force file at line %d.\n",
900 lineNum );
901 painCave.isFatal = 1;
902 simError();
903 }
904
905 // stop reading at end of file, or at next section
906 while( readLine[0] != '#' && eof_test != NULL ){
907
908 // toss comment lines
909 if( readLine[0] != '!' ){
910
911 // the parser returns 0 if the line was blank
912 if( parseBend( readLine, lineNum, bendInfo ) ){
913 headBendType->add( bendInfo );
914 }
915 }
916 eof_test = fgets( readLine, sizeof(readLine), frcFile );
917 lineNum++;
918 }
919
920 #ifdef IS_MPI
921
922 // send out the linked list to all the other processes
923
924 sprintf( checkPointMsg,
925 "TraPPE_Ex bend structures read successfully." );
926 MPIcheckPoint();
927
928 currentBendType = headBendType->next;
929 while( currentBendType != NULL ){
930 currentBendType->duplicate( bendInfo );
931 sendFrcStruct( &bendInfo, mpiBendStructType );
932 currentBendType = currentBendType->next;
933 }
934 bendInfo.last = 1;
935 sendFrcStruct( &bendInfo, mpiBendStructType );
936
937 }
938
939 else{
940
941 // listen for node 0 to send out the force params
942
943 MPIcheckPoint();
944
945 headBendType = new LinkedBendType;
946 recieveFrcStruct( &bendInfo, mpiBendStructType );
947 while( !bendInfo.last ){
948
949 headBendType->add( bendInfo );
950 recieveFrcStruct( &bendInfo, mpiBendStructType );
951 }
952 }
953 #endif // is_mpi
954
955
956 // read in the torsions
957
958 #ifdef IS_MPI
959 if( worldRank == 0 ){
960 #endif
961
962 // read in the torsion types.
963
964 headTorsionType = new LinkedTorsionType;
965
966 fastForward( "TorsionTypes", "initializeTorsions" );
967
968 // we are now at the torsionTypes section
969
970 eof_test = fgets( readLine, sizeof(readLine), frcFile );
971 lineNum++;
972
973
974 // read a line, and start parseing out the atom types
975
976 if( eof_test == NULL ){
977 sprintf( painCave.errMsg,
978 "Error in reading torsions from force file at line %d.\n",
979 lineNum );
980 painCave.isFatal = 1;
981 simError();
982 }
983
984 // stop reading at end of file, or at next section
985 while( readLine[0] != '#' && eof_test != NULL ){
986
987 // toss comment lines
988 if( readLine[0] != '!' ){
989
990 // the parser returns 0 if the line was blank
991 if( parseTorsion( readLine, lineNum, torsionInfo ) ){
992 headTorsionType->add( torsionInfo );
993
994 }
995 }
996 eof_test = fgets( readLine, sizeof(readLine), frcFile );
997 lineNum++;
998 }
999
1000 #ifdef IS_MPI
1001
1002 // send out the linked list to all the other processes
1003
1004 sprintf( checkPointMsg,
1005 "TraPPE_Ex torsion structures read successfully." );
1006 MPIcheckPoint();
1007
1008 currentTorsionType = headTorsionType->next;
1009 while( currentTorsionType != NULL ){
1010 currentTorsionType->duplicate( torsionInfo );
1011 sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1012 currentTorsionType = currentTorsionType->next;
1013 }
1014 torsionInfo.last = 1;
1015 sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1016
1017 }
1018
1019 else{
1020
1021 // listen for node 0 to send out the force params
1022
1023 MPIcheckPoint();
1024
1025 headTorsionType = new LinkedTorsionType;
1026 recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1027 while( !torsionInfo.last ){
1028
1029 headTorsionType->add( torsionInfo );
1030 recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1031 }
1032 }
1033 #endif // is_mpi
1034
1035 entry_plug->useLJ = 1;
1036 }
1037
1038
1039
1040 void TraPPE_ExFF::initializeAtoms( int nAtoms, Atom** the_atoms ){
1041
1042
1043 //////////////////////////////////////////////////
1044 // a quick water fix
1045
1046 double waterI[3][3];
1047 waterI[0][0] = 1.76958347772500;
1048 waterI[0][1] = 0.0;
1049 waterI[0][2] = 0.0;
1050
1051 waterI[1][0] = 0.0;
1052 waterI[1][1] = 0.614537057924513;
1053 waterI[1][2] = 0.0;
1054
1055 waterI[2][0] = 0.0;
1056 waterI[2][1] = 0.0;
1057 waterI[2][2] = 1.15504641980049;
1058
1059
1060 double headI[3][3];
1061 headI[0][0] = 1125;
1062 headI[0][1] = 0.0;
1063 headI[0][2] = 0.0;
1064
1065 headI[1][0] = 0.0;
1066 headI[1][1] = 1125;
1067 headI[1][2] = 0.0;
1068
1069 headI[2][0] = 0.0;
1070 headI[2][1] = 0.0;
1071 headI[2][2] = 250;
1072
1073 //////////////////////////////////////////////////
1074
1075
1076 // initialize the atoms
1077
1078 DirectionalAtom* dAtom;
1079
1080 for(int i=0; i<nAtoms; i++ ){
1081
1082 currentAtomType = headAtomType->find( the_atoms[i]->getType() );
1083 if( currentAtomType == NULL ){
1084 sprintf( painCave.errMsg,
1085 "AtomType error, %s not found in force file.\n",
1086 the_atoms[i]->getType() );
1087 painCave.isFatal = 1;
1088 simError();
1089 }
1090
1091 the_atoms[i]->setMass( currentAtomType->mass );
1092 the_atoms[i]->setEpslon( currentAtomType->epslon );
1093 the_atoms[i]->setSigma( currentAtomType->sigma );
1094 the_atoms[i]->setIdent( currentAtomType->ident );
1095 the_atoms[i]->setLJ();
1096
1097 if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
1098
1099 if( currentAtomType->isDipole ){
1100 if( the_atoms[i]->isDirectional() ){
1101
1102 dAtom = (DirectionalAtom *) the_atoms[i];
1103 dAtom->setMu( currentAtomType->dipole );
1104 dAtom->setHasDipole( 1 );
1105 dAtom->setJx( 0.0 );
1106 dAtom->setJy( 0.0 );
1107 dAtom->setJz( 0.0 );
1108
1109 if(!strcmp("SSD",the_atoms[i]->getType())){
1110 dAtom->setI( waterI );
1111 dAtom->setSSD( 1 );
1112 }
1113 else if(!strcmp("HEAD",the_atoms[i]->getType())){
1114 dAtom->setI( headI );
1115 dAtom->setSSD( 0 );
1116 }
1117 else{
1118 sprintf(painCave.errMsg,
1119 "AtmType error, %s does not have a moment of inertia set.\n",
1120 the_atoms[i]->getType() );
1121 painCave.isFatal = 1;
1122 simError();
1123 }
1124 entry_plug->n_dipoles++;
1125 }
1126 else{
1127
1128 sprintf( painCave.errMsg,
1129 "TraPPE_ExFF error: Atom \"%s\" is a dipole, yet no standard"
1130 " orientation was specifed in the BASS file.\n",
1131 currentAtomType->name );
1132 painCave.isFatal = 1;
1133 simError();
1134 }
1135 }
1136 else{
1137 if( the_atoms[i]->isDirectional() ){
1138 sprintf( painCave.errMsg,
1139 "TraPPE_ExFF error: Atom \"%s\" was given a standard"
1140 "orientation in the BASS file, yet it is not a dipole.\n",
1141 currentAtomType->name);
1142 painCave.isFatal = 1;
1143 simError();
1144 }
1145 }
1146 }
1147 }
1148
1149 void TraPPE_ExFF::initializeBonds( int nBonds, Bond** bondArray,
1150 bond_pair* the_bonds ){
1151 int i,a,b;
1152 char* atomA;
1153 char* atomB;
1154
1155 Atom** the_atoms;
1156 the_atoms = entry_plug->atoms;
1157
1158
1159 // initialize the Bonds
1160
1161 for( i=0; i<nBonds; i++ ){
1162
1163 a = the_bonds[i].a;
1164 b = the_bonds[i].b;
1165
1166 atomA = the_atoms[a]->getType();
1167 atomB = the_atoms[b]->getType();
1168 currentBondType = headBondType->find( atomA, atomB );
1169 if( currentBondType == NULL ){
1170 sprintf( painCave.errMsg,
1171 "BondType error, %s - %s not found in force file.\n",
1172 atomA, atomB );
1173 painCave.isFatal = 1;
1174 simError();
1175 }
1176
1177 if( !strcmp( currentBondType->type, "fixed" ) ){
1178
1179 bondArray[i] = new ConstrainedBond( *the_atoms[a],
1180 *the_atoms[b],
1181 currentBondType->d0 );
1182 entry_plug->n_constraints++;
1183 }
1184 }
1185 }
1186
1187 void TraPPE_ExFF::initializeBends( int nBends, Bend** bendArray,
1188 bend_set* the_bends ){
1189
1190 QuadraticBend* qBend;
1191 GhostBend* gBend;
1192 Atom** the_atoms;
1193 the_atoms = entry_plug->atoms;
1194
1195 int i, a, b, c;
1196 char* atomA;
1197 char* atomB;
1198 char* atomC;
1199
1200 // initialize the Bends
1201
1202 for( i=0; i<nBends; i++ ){
1203
1204 a = the_bends[i].a;
1205 b = the_bends[i].b;
1206 c = the_bends[i].c;
1207
1208 atomA = the_atoms[a]->getType();
1209 atomB = the_atoms[b]->getType();
1210
1211 if( the_bends[i].isGhost ) atomC = "GHOST";
1212 else atomC = the_atoms[c]->getType();
1213
1214 currentBendType = headBendType->find( atomA, atomB, atomC );
1215 if( currentBendType == NULL ){
1216 sprintf( painCave.errMsg, "BendType error, %s - %s - %s not found"
1217 " in force file.\n",
1218 atomA, atomB, atomC );
1219 painCave.isFatal = 1;
1220 simError();
1221 }
1222
1223 if( !strcmp( currentBendType->type, "quadratic" ) ){
1224
1225 if( the_bends[i].isGhost){
1226
1227 if( the_bends[i].ghost == b ){
1228 // do nothing
1229 }
1230 else if( the_bends[i].ghost == a ){
1231 c = a;
1232 a = b;
1233 b = c;
1234 }
1235 else{
1236 sprintf( painCave.errMsg,
1237 "BendType error, %s - %s - %s,\n"
1238 " --> central atom is not "
1239 "correctly identified with the "
1240 "\"ghostVectorSource = \" tag.\n",
1241 atomA, atomB, atomC );
1242 painCave.isFatal = 1;
1243 simError();
1244 }
1245
1246 gBend = new GhostBend( *the_atoms[a],
1247 *the_atoms[b] );
1248 gBend->setConstants( currentBendType->k1,
1249 currentBendType->k2,
1250 currentBendType->k3,
1251 currentBendType->t0 );
1252 bendArray[i] = gBend;
1253 }
1254 else{
1255 qBend = new QuadraticBend( *the_atoms[a],
1256 *the_atoms[b],
1257 *the_atoms[c] );
1258 qBend->setConstants( currentBendType->k1,
1259 currentBendType->k2,
1260 currentBendType->k3,
1261 currentBendType->t0 );
1262 bendArray[i] = qBend;
1263 }
1264 }
1265 }
1266 }
1267
1268 void TraPPE_ExFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
1269 torsion_set* the_torsions ){
1270
1271 int i, a, b, c, d;
1272 char* atomA;
1273 char* atomB;
1274 char* atomC;
1275 char* atomD;
1276
1277 CubicTorsion* cTors;
1278 Atom** the_atoms;
1279 the_atoms = entry_plug->atoms;
1280
1281 // initialize the Torsions
1282
1283 for( i=0; i<nTorsions; i++ ){
1284
1285 a = the_torsions[i].a;
1286 b = the_torsions[i].b;
1287 c = the_torsions[i].c;
1288 d = the_torsions[i].d;
1289
1290 atomA = the_atoms[a]->getType();
1291 atomB = the_atoms[b]->getType();
1292 atomC = the_atoms[c]->getType();
1293 atomD = the_atoms[d]->getType();
1294 currentTorsionType = headTorsionType->find( atomA, atomB, atomC, atomD );
1295 if( currentTorsionType == NULL ){
1296 sprintf( painCave.errMsg,
1297 "TorsionType error, %s - %s - %s - %s not found"
1298 " in force file.\n",
1299 atomA, atomB, atomC, atomD );
1300 painCave.isFatal = 1;
1301 simError();
1302 }
1303
1304 if( !strcmp( currentTorsionType->type, "cubic" ) ){
1305
1306 cTors = new CubicTorsion( *the_atoms[a], *the_atoms[b],
1307 *the_atoms[c], *the_atoms[d] );
1308 cTors->setConstants( currentTorsionType->k1, currentTorsionType->k2,
1309 currentTorsionType->k3, currentTorsionType->k4 );
1310 torsionArray[i] = cTors;
1311 }
1312 }
1313 }
1314
1315 void TraPPE_ExFF::fastForward( char* stopText, char* searchOwner ){
1316
1317 int foundText = 0;
1318 char* the_token;
1319
1320 rewind( frcFile );
1321 lineNum = 0;
1322
1323 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1324 lineNum++;
1325 if( eof_test == NULL ){
1326 sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
1327 " file is empty.\n",
1328 searchOwner );
1329 painCave.isFatal = 1;
1330 simError();
1331 }
1332
1333
1334 while( !foundText ){
1335 while( eof_test != NULL && readLine[0] != '#' ){
1336 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1337 lineNum++;
1338 }
1339 if( eof_test == NULL ){
1340 sprintf( painCave.errMsg,
1341 "Error fast forwarding force file for %s at "
1342 "line %d: file ended unexpectedly.\n",
1343 searchOwner,
1344 lineNum );
1345 painCave.isFatal = 1;
1346 simError();
1347 }
1348
1349 the_token = strtok( readLine, " ,;\t#\n" );
1350 foundText = !strcmp( stopText, the_token );
1351
1352 if( !foundText ){
1353 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1354 lineNum++;
1355
1356 if( eof_test == NULL ){
1357 sprintf( painCave.errMsg,
1358 "Error fast forwarding force file for %s at "
1359 "line %d: file ended unexpectedly.\n",
1360 searchOwner,
1361 lineNum );
1362 painCave.isFatal = 1;
1363 simError();
1364 }
1365 }
1366 }
1367 }
1368
1369
1370 int TPE::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1371
1372 char* the_token;
1373
1374 the_token = strtok( lineBuffer, " \n\t,;" );
1375 if( the_token != NULL ){
1376
1377 strcpy( info.name, the_token );
1378
1379 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1380 sprintf( painCave.errMsg,
1381 "Error parseing AtomTypes: line %d\n", lineNum );
1382 painCave.isFatal = 1;
1383 simError();
1384 }
1385
1386 info.isDipole = atoi( the_token );
1387
1388 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1389 sprintf( painCave.errMsg,
1390 "Error parseing AtomTypes: line %d\n", lineNum );
1391 painCave.isFatal = 1;
1392 simError();
1393 }
1394
1395 info.isSSD = atoi( the_token );
1396
1397 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1398 sprintf( painCave.errMsg,
1399 "Error parseing AtomTypes: line %d\n", lineNum );
1400 painCave.isFatal = 1;
1401 simError();
1402 }
1403
1404 info.mass = atof( the_token );
1405
1406 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1407 sprintf( painCave.errMsg,
1408 "Error parseing AtomTypes: line %d\n", lineNum );
1409 painCave.isFatal = 1;
1410 simError();
1411 }
1412
1413 info.epslon = atof( the_token );
1414
1415 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1416 sprintf( painCave.errMsg,
1417 "Error parseing AtomTypes: line %d\n", lineNum );
1418 painCave.isFatal = 1;
1419 simError();
1420 }
1421
1422 info.sigma = atof( the_token );
1423
1424 if( info.isDipole ){
1425
1426 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1427 sprintf( painCave.errMsg,
1428 "Error parseing AtomTypes: line %d\n", lineNum );
1429 painCave.isFatal = 1;
1430 simError();
1431 }
1432
1433 info.dipole = atof( the_token );
1434 }
1435 else info.dipole = 0.0;
1436
1437 if( info.isSSD ){
1438
1439 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1440 sprintf( painCave.errMsg,
1441 "Error parseing AtomTypes: line %d\n", lineNum );
1442 painCave.isFatal = 1;
1443 simError();
1444 }
1445
1446 info.w0 = atof( the_token );
1447
1448 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1449 sprintf( painCave.errMsg,
1450 "Error parseing AtomTypes: line %d\n", lineNum );
1451 painCave.isFatal = 1;
1452 simError();
1453 }
1454
1455 info.v0 = atof( the_token );
1456 }
1457 else info.v0 = info.w0 = 0.0;
1458
1459 return 1;
1460 }
1461 else return 0;
1462 }
1463
1464 int TPE::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1465
1466 char* the_token;
1467
1468 the_token = strtok( lineBuffer, " \n\t,;" );
1469 if( the_token != NULL ){
1470
1471 strcpy( info.nameA, the_token );
1472
1473 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1474 sprintf( painCave.errMsg,
1475 "Error parseing BondTypes: line %d\n", lineNum );
1476 painCave.isFatal = 1;
1477 simError();
1478 }
1479
1480 strcpy( info.nameB, the_token );
1481
1482 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1483 sprintf( painCave.errMsg,
1484 "Error parseing BondTypes: line %d\n", lineNum );
1485 painCave.isFatal = 1;
1486 simError();
1487 }
1488
1489 strcpy( info.type, the_token );
1490
1491 if( !strcmp( info.type, "fixed" ) ){
1492 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1493 sprintf( painCave.errMsg,
1494 "Error parseing BondTypes: line %d\n", lineNum );
1495 painCave.isFatal = 1;
1496 simError();
1497 }
1498
1499 info.d0 = atof( the_token );
1500 }
1501 else{
1502 sprintf( painCave.errMsg,
1503 "Unknown TraPPE_Ex bond type \"%s\" at line %d\n",
1504 info.type,
1505 lineNum );
1506 painCave.isFatal = 1;
1507 simError();
1508 }
1509
1510 return 1;
1511 }
1512 else return 0;
1513 }
1514
1515
1516 int TPE::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1517
1518 char* the_token;
1519
1520 the_token = strtok( lineBuffer, " \n\t,;" );
1521 if( the_token != NULL ){
1522
1523 strcpy( info.nameA, the_token );
1524
1525 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1526 sprintf( painCave.errMsg,
1527 "Error parseing BendTypes: line %d\n", lineNum );
1528 painCave.isFatal = 1;
1529 simError();
1530 }
1531
1532 strcpy( info.nameB, the_token );
1533
1534 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1535 sprintf( painCave.errMsg,
1536 "Error parseing BendTypes: line %d\n", lineNum );
1537 painCave.isFatal = 1;
1538 simError();
1539 }
1540
1541 strcpy( info.nameC, the_token );
1542
1543 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1544 sprintf( painCave.errMsg,
1545 "Error parseing BendTypes: line %d\n", lineNum );
1546 painCave.isFatal = 1;
1547 simError();
1548 }
1549
1550 strcpy( info.type, the_token );
1551
1552 if( !strcmp( info.type, "quadratic" ) ){
1553 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1554 sprintf( painCave.errMsg,
1555 "Error parseing BendTypes: line %d\n", lineNum );
1556 painCave.isFatal = 1;
1557 simError();
1558 }
1559
1560 info.k1 = atof( the_token );
1561
1562 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1563 sprintf( painCave.errMsg,
1564 "Error parseing BendTypes: line %d\n", lineNum );
1565 painCave.isFatal = 1;
1566 simError();
1567 }
1568
1569 info.k2 = atof( the_token );
1570
1571 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1572 sprintf( painCave.errMsg,
1573 "Error parseing BendTypes: line %d\n", lineNum );
1574 painCave.isFatal = 1;
1575 simError();
1576 }
1577
1578 info.k3 = atof( the_token );
1579
1580 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1581 sprintf( painCave.errMsg,
1582 "Error parseing BendTypes: line %d\n", lineNum );
1583 painCave.isFatal = 1;
1584 simError();
1585 }
1586
1587 info.t0 = atof( the_token );
1588 }
1589
1590 else{
1591 sprintf( painCave.errMsg,
1592 "Unknown TraPPE_Ex bend type \"%s\" at line %d\n",
1593 info.type,
1594 lineNum );
1595 painCave.isFatal = 1;
1596 simError();
1597 }
1598
1599 return 1;
1600 }
1601 else return 0;
1602 }
1603
1604 int TPE::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1605
1606 char* the_token;
1607
1608 the_token = strtok( lineBuffer, " \n\t,;" );
1609 if( the_token != NULL ){
1610
1611 strcpy( info.nameA, the_token );
1612
1613 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1614 sprintf( painCave.errMsg,
1615 "Error parseing TorsionTypes: line %d\n", lineNum );
1616 painCave.isFatal = 1;
1617 simError();
1618 }
1619
1620 strcpy( info.nameB, the_token );
1621
1622 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1623 sprintf( painCave.errMsg,
1624 "Error parseing TorsionTypes: line %d\n", lineNum );
1625 painCave.isFatal = 1;
1626 simError();
1627 }
1628
1629 strcpy( info.nameC, the_token );
1630
1631 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1632 sprintf( painCave.errMsg,
1633 "Error parseing TorsionTypes: line %d\n", lineNum );
1634 painCave.isFatal = 1;
1635 simError();
1636 }
1637
1638 strcpy( info.nameD, the_token );
1639
1640 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1641 sprintf( painCave.errMsg,
1642 "Error parseing TorsionTypes: line %d\n", lineNum );
1643 painCave.isFatal = 1;
1644 simError();
1645 }
1646
1647 strcpy( info.type, the_token );
1648
1649 if( !strcmp( info.type, "cubic" ) ){
1650 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1651 sprintf( painCave.errMsg,
1652 "Error parseing TorsionTypes: line %d\n", lineNum );
1653 painCave.isFatal = 1;
1654 simError();
1655 }
1656
1657 info.k1 = atof( the_token );
1658
1659 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1660 sprintf( painCave.errMsg,
1661 "Error parseing TorsionTypes: line %d\n", lineNum );
1662 painCave.isFatal = 1;
1663 simError();
1664 }
1665
1666 info.k2 = atof( the_token );
1667
1668 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1669 sprintf( painCave.errMsg,
1670 "Error parseing TorsionTypes: line %d\n", lineNum );
1671 painCave.isFatal = 1;
1672 simError();
1673 }
1674
1675 info.k3 = atof( the_token );
1676
1677 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1678 sprintf( painCave.errMsg,
1679 "Error parseing TorsionTypes: line %d\n", lineNum );
1680 painCave.isFatal = 1;
1681 simError();
1682 }
1683
1684 info.k4 = atof( the_token );
1685
1686 }
1687
1688 else{
1689 sprintf( painCave.errMsg,
1690 "Unknown TraPPE_Ex torsion type \"%s\" at line %d\n",
1691 info.type,
1692 lineNum );
1693 painCave.isFatal = 1;
1694 simError();
1695 }
1696
1697 return 1;
1698 }
1699
1700 else return 0;
1701 }