ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 976
Committed: Thu Jan 22 17:34:20 2004 UTC (20 years, 5 months ago) by chrisfen
File size: 40994 byte(s)
Log Message: 
Corrected spelling in several directories, and stated WATER.cpp

File Contents

# Content
1 #include <stdlib.h>
2 #include <stdio.h>
3 #include <string.h>
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 receiveFrcStruct( &atomInfo, mpiAtomStructType );
732
733 while( !atomInfo.last ){
734
735
736
737 headAtomType->add( atomInfo );
738
739 MPIcheckPoint();
740
741 receiveFrcStruct( &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 int isEAM = 0;
758 int isCharge = 0;
759 double GB_dummy = 0.0;
760 double charge = 0.0;
761
762
763 currentAtomType = headAtomType->next;;
764 while( currentAtomType != NULL ){
765
766 if(currentAtomType->isDipole) entry_plug->useDipole = 1;
767 if(currentAtomType->isSSD) {
768 entry_plug->useSticky = 1;
769 set_sticky_params( &(currentAtomType->w0), &(currentAtomType->v0));
770 }
771
772 if( currentAtomType->name[0] != '\0' ){
773 isError = 0;
774 makeAtype( &(currentAtomType->ident),
775 &isLJ,
776 &(currentAtomType->isSSD),
777 &(currentAtomType->isDipole),
778 &isGB,
779 &isEAM,
780 &isCharge,
781 &(currentAtomType->epslon),
782 &(currentAtomType->sigma),
783 &charge,
784 &(currentAtomType->dipole),
785 &isError );
786 if( isError ){
787 sprintf( painCave.errMsg,
788 "Error initializing the \"%s\" atom type in fortran\n",
789 currentAtomType->name );
790 painCave.isFatal = 1;
791 simError();
792 }
793 }
794 currentAtomType = currentAtomType->next;
795 }
796
797 #ifdef IS_MPI
798 sprintf( checkPointMsg,
799 "TraPPE_ExFF atom structures successfully sent to fortran\n" );
800 MPIcheckPoint();
801 #endif // is_mpi
802
803
804
805 // read in the bonds
806
807 #ifdef IS_MPI
808 if( worldRank == 0 ){
809 #endif
810
811 // read in the bond types.
812
813 headBondType = new LinkedBondType;
814
815 fastForward( "BondTypes", "initializeBonds" );
816
817 // we are now at the bondTypes section
818
819 eof_test = fgets( readLine, sizeof(readLine), frcFile );
820 lineNum++;
821
822
823 // read a line, and start parseing out the atom types
824
825 if( eof_test == NULL ){
826 sprintf( painCave.errMsg,
827 "Error in reading bonds from force file at line %d.\n",
828 lineNum );
829 painCave.isFatal = 1;
830 simError();
831 }
832
833 // stop reading at end of file, or at next section
834 while( readLine[0] != '#' && eof_test != NULL ){
835
836 // toss comment lines
837 if( readLine[0] != '!' ){
838
839 // the parser returns 0 if the line was blank
840 if( parseBond( readLine, lineNum, bondInfo ) ){
841 headBondType->add( bondInfo );
842 }
843 }
844 eof_test = fgets( readLine, sizeof(readLine), frcFile );
845 lineNum++;
846 }
847
848 #ifdef IS_MPI
849
850 // send out the linked list to all the other processes
851
852 sprintf( checkPointMsg,
853 "TraPPE_Ex bond structures read successfully." );
854 MPIcheckPoint();
855
856 currentBondType = headBondType->next;
857 while( currentBondType != NULL ){
858 currentBondType->duplicate( bondInfo );
859 sendFrcStruct( &bondInfo, mpiBondStructType );
860 currentBondType = currentBondType->next;
861 }
862 bondInfo.last = 1;
863 sendFrcStruct( &bondInfo, mpiBondStructType );
864
865 }
866
867 else{
868
869 // listen for node 0 to send out the force params
870
871 MPIcheckPoint();
872
873 headBondType = new LinkedBondType;
874 receiveFrcStruct( &bondInfo, mpiBondStructType );
875 while( !bondInfo.last ){
876
877 headBondType->add( bondInfo );
878 receiveFrcStruct( &bondInfo, mpiBondStructType );
879 }
880 }
881
882 sprintf( checkPointMsg,
883 "TraPPE_ExFF bond structures broadcast successfully." );
884 MPIcheckPoint();
885
886 #endif // is_mpi
887
888
889 // read in the bends
890
891 #ifdef IS_MPI
892 if( worldRank == 0 ){
893 #endif
894
895 // read in the bend types.
896
897 headBendType = new LinkedBendType;
898
899 fastForward( "BendTypes", "initializeBends" );
900
901 // we are now at the bendTypes section
902
903 eof_test = fgets( readLine, sizeof(readLine), frcFile );
904 lineNum++;
905
906 // read a line, and start parseing out the bend types
907
908 if( eof_test == NULL ){
909 sprintf( painCave.errMsg,
910 "Error in reading bends from force file at line %d.\n",
911 lineNum );
912 painCave.isFatal = 1;
913 simError();
914 }
915
916 // stop reading at end of file, or at next section
917 while( readLine[0] != '#' && eof_test != NULL ){
918
919 // toss comment lines
920 if( readLine[0] != '!' ){
921
922 // the parser returns 0 if the line was blank
923 if( parseBend( readLine, lineNum, bendInfo ) ){
924 headBendType->add( bendInfo );
925 }
926 }
927 eof_test = fgets( readLine, sizeof(readLine), frcFile );
928 lineNum++;
929 }
930
931 #ifdef IS_MPI
932
933 // send out the linked list to all the other processes
934
935 sprintf( checkPointMsg,
936 "TraPPE_Ex bend structures read successfully." );
937 MPIcheckPoint();
938
939 currentBendType = headBendType->next;
940 while( currentBendType != NULL ){
941 currentBendType->duplicate( bendInfo );
942 sendFrcStruct( &bendInfo, mpiBendStructType );
943 currentBendType = currentBendType->next;
944 }
945 bendInfo.last = 1;
946 sendFrcStruct( &bendInfo, mpiBendStructType );
947
948 }
949
950 else{
951
952 // listen for node 0 to send out the force params
953
954 MPIcheckPoint();
955
956 headBendType = new LinkedBendType;
957 receiveFrcStruct( &bendInfo, mpiBendStructType );
958 while( !bendInfo.last ){
959
960 headBendType->add( bendInfo );
961 receiveFrcStruct( &bendInfo, mpiBendStructType );
962 }
963 }
964
965 sprintf( checkPointMsg,
966 "TraPPE_ExFF bend structures broadcast successfully." );
967 MPIcheckPoint();
968
969 #endif // is_mpi
970
971
972 // read in the torsions
973
974 #ifdef IS_MPI
975 if( worldRank == 0 ){
976 #endif
977
978 // read in the torsion types.
979
980 headTorsionType = new LinkedTorsionType;
981
982 fastForward( "TorsionTypes", "initializeTorsions" );
983
984 // we are now at the torsionTypes section
985
986 eof_test = fgets( readLine, sizeof(readLine), frcFile );
987 lineNum++;
988
989
990 // read a line, and start parseing out the atom types
991
992 if( eof_test == NULL ){
993 sprintf( painCave.errMsg,
994 "Error in reading torsions from force file at line %d.\n",
995 lineNum );
996 painCave.isFatal = 1;
997 simError();
998 }
999
1000 // stop reading at end of file, or at next section
1001 while( readLine[0] != '#' && eof_test != NULL ){
1002
1003 // toss comment lines
1004 if( readLine[0] != '!' ){
1005
1006 // the parser returns 0 if the line was blank
1007 if( parseTorsion( readLine, lineNum, torsionInfo ) ){
1008 headTorsionType->add( torsionInfo );
1009
1010 }
1011 }
1012 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1013 lineNum++;
1014 }
1015
1016 #ifdef IS_MPI
1017
1018 // send out the linked list to all the other processes
1019
1020 sprintf( checkPointMsg,
1021 "TraPPE_Ex torsion structures read successfully." );
1022 MPIcheckPoint();
1023
1024 currentTorsionType = headTorsionType->next;
1025 while( currentTorsionType != NULL ){
1026 currentTorsionType->duplicate( torsionInfo );
1027 sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1028 currentTorsionType = currentTorsionType->next;
1029 }
1030 torsionInfo.last = 1;
1031 sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1032
1033 }
1034
1035 else{
1036
1037 // listen for node 0 to send out the force params
1038
1039 MPIcheckPoint();
1040
1041 headTorsionType = new LinkedTorsionType;
1042 receiveFrcStruct( &torsionInfo, mpiTorsionStructType );
1043 while( !torsionInfo.last ){
1044
1045 headTorsionType->add( torsionInfo );
1046 receiveFrcStruct( &torsionInfo, mpiTorsionStructType );
1047 }
1048 }
1049
1050 sprintf( checkPointMsg,
1051 "TraPPE_ExFF torsion structures broadcast successfully." );
1052 MPIcheckPoint();
1053
1054 #endif // is_mpi
1055
1056 entry_plug->useLJ = 1;
1057 }
1058
1059
1060
1061 void TraPPE_ExFF::initializeAtoms( int nAtoms, Atom** the_atoms ){
1062
1063
1064 //////////////////////////////////////////////////
1065 // a quick water fix
1066
1067 double waterI[3][3];
1068 waterI[0][0] = 1.76958347772500;
1069 waterI[0][1] = 0.0;
1070 waterI[0][2] = 0.0;
1071
1072 waterI[1][0] = 0.0;
1073 waterI[1][1] = 0.614537057924513;
1074 waterI[1][2] = 0.0;
1075
1076 waterI[2][0] = 0.0;
1077 waterI[2][1] = 0.0;
1078 waterI[2][2] = 1.15504641980049;
1079
1080
1081 double headI[3][3];
1082 headI[0][0] = 1125;
1083 headI[0][1] = 0.0;
1084 headI[0][2] = 0.0;
1085
1086 headI[1][0] = 0.0;
1087 headI[1][1] = 1125;
1088 headI[1][2] = 0.0;
1089
1090 headI[2][0] = 0.0;
1091 headI[2][1] = 0.0;
1092 headI[2][2] = 250;
1093
1094 //////////////////////////////////////////////////
1095
1096
1097 // initialize the atoms
1098
1099 DirectionalAtom* dAtom;
1100
1101 for(int i=0; i<nAtoms; i++ ){
1102
1103 currentAtomType = headAtomType->find( the_atoms[i]->getType() );
1104 if( currentAtomType == NULL ){
1105 sprintf( painCave.errMsg,
1106 "AtomType error, %s not found in force file.\n",
1107 the_atoms[i]->getType() );
1108 painCave.isFatal = 1;
1109 simError();
1110 }
1111
1112 the_atoms[i]->setMass( currentAtomType->mass );
1113 the_atoms[i]->setEpslon( currentAtomType->epslon );
1114 the_atoms[i]->setSigma( currentAtomType->sigma );
1115 the_atoms[i]->setIdent( currentAtomType->ident );
1116 the_atoms[i]->setLJ();
1117
1118 if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
1119
1120 if( currentAtomType->isDipole ){
1121 if( the_atoms[i]->isDirectional() ){
1122
1123 dAtom = (DirectionalAtom *) the_atoms[i];
1124 dAtom->setMu( currentAtomType->dipole );
1125 dAtom->setHasDipole( 1 );
1126 dAtom->setJx( 0.0 );
1127 dAtom->setJy( 0.0 );
1128 dAtom->setJz( 0.0 );
1129
1130 if(!strcmp("SSD",the_atoms[i]->getType())){
1131 dAtom->setI( waterI );
1132 dAtom->setSSD( 1 );
1133 }
1134 else if(!strcmp("HEAD",the_atoms[i]->getType())){
1135 dAtom->setI( headI );
1136 dAtom->setSSD( 0 );
1137 }
1138 else{
1139 sprintf(painCave.errMsg,
1140 "AtmType error, %s does not have a moment of inertia set.\n",
1141 the_atoms[i]->getType() );
1142 painCave.isFatal = 1;
1143 simError();
1144 }
1145 entry_plug->n_dipoles++;
1146 }
1147 else{
1148
1149 sprintf( painCave.errMsg,
1150 "TraPPE_ExFF error: Atom \"%s\" is a dipole, yet no standard"
1151 " orientation was specifed in the BASS file.\n",
1152 currentAtomType->name );
1153 painCave.isFatal = 1;
1154 simError();
1155 }
1156 }
1157 else{
1158 if( the_atoms[i]->isDirectional() ){
1159 sprintf( painCave.errMsg,
1160 "TraPPE_ExFF error: Atom \"%s\" was given a standard"
1161 "orientation in the BASS file, yet it is not a dipole.\n",
1162 currentAtomType->name);
1163 painCave.isFatal = 1;
1164 simError();
1165 }
1166 }
1167 }
1168 }
1169
1170 void TraPPE_ExFF::initializeBonds( int nBonds, Bond** bondArray,
1171 bond_pair* the_bonds ){
1172 int i,a,b;
1173 char* atomA;
1174 char* atomB;
1175
1176 Atom** the_atoms;
1177 the_atoms = entry_plug->atoms;
1178
1179
1180 // initialize the Bonds
1181
1182 for( i=0; i<nBonds; i++ ){
1183
1184 a = the_bonds[i].a;
1185 b = the_bonds[i].b;
1186
1187 atomA = the_atoms[a]->getType();
1188 atomB = the_atoms[b]->getType();
1189 currentBondType = headBondType->find( atomA, atomB );
1190 if( currentBondType == NULL ){
1191 sprintf( painCave.errMsg,
1192 "BondType error, %s - %s not found in force file.\n",
1193 atomA, atomB );
1194 painCave.isFatal = 1;
1195 simError();
1196 }
1197
1198 if( !strcmp( currentBondType->type, "fixed" ) ){
1199
1200 bondArray[i] = new ConstrainedBond( *the_atoms[a],
1201 *the_atoms[b],
1202 currentBondType->d0 );
1203 entry_plug->n_constraints++;
1204 }
1205 }
1206 }
1207
1208 void TraPPE_ExFF::initializeBends( int nBends, Bend** bendArray,
1209 bend_set* the_bends ){
1210
1211 QuadraticBend* qBend;
1212 GhostBend* gBend;
1213 Atom** the_atoms;
1214 the_atoms = entry_plug->atoms;
1215
1216 int i, a, b, c;
1217 char* atomA;
1218 char* atomB;
1219 char* atomC;
1220
1221 // initialize the Bends
1222
1223 for( i=0; i<nBends; i++ ){
1224
1225 a = the_bends[i].a;
1226 b = the_bends[i].b;
1227 c = the_bends[i].c;
1228
1229 atomA = the_atoms[a]->getType();
1230 atomB = the_atoms[b]->getType();
1231
1232 if( the_bends[i].isGhost ) atomC = "GHOST";
1233 else atomC = the_atoms[c]->getType();
1234
1235 currentBendType = headBendType->find( atomA, atomB, atomC );
1236 if( currentBendType == NULL ){
1237 sprintf( painCave.errMsg, "BendType error, %s - %s - %s not found"
1238 " in force file.\n",
1239 atomA, atomB, atomC );
1240 painCave.isFatal = 1;
1241 simError();
1242 }
1243
1244 if( !strcmp( currentBendType->type, "quadratic" ) ){
1245
1246 if( the_bends[i].isGhost){
1247
1248 if( the_bends[i].ghost == b ){
1249 // do nothing
1250 }
1251 else if( the_bends[i].ghost == a ){
1252 c = a;
1253 a = b;
1254 b = c;
1255 }
1256 else{
1257 sprintf( painCave.errMsg,
1258 "BendType error, %s - %s - %s,\n"
1259 " --> central atom is not "
1260 "correctly identified with the "
1261 "\"ghostVectorSource = \" tag.\n",
1262 atomA, atomB, atomC );
1263 painCave.isFatal = 1;
1264 simError();
1265 }
1266
1267 gBend = new GhostBend( *the_atoms[a],
1268 *the_atoms[b] );
1269 gBend->setConstants( currentBendType->k1,
1270 currentBendType->k2,
1271 currentBendType->k3,
1272 currentBendType->t0 );
1273 bendArray[i] = gBend;
1274 }
1275 else{
1276 qBend = new QuadraticBend( *the_atoms[a],
1277 *the_atoms[b],
1278 *the_atoms[c] );
1279 qBend->setConstants( currentBendType->k1,
1280 currentBendType->k2,
1281 currentBendType->k3,
1282 currentBendType->t0 );
1283 bendArray[i] = qBend;
1284 }
1285 }
1286 }
1287 }
1288
1289 void TraPPE_ExFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
1290 torsion_set* the_torsions ){
1291
1292 int i, a, b, c, d;
1293 char* atomA;
1294 char* atomB;
1295 char* atomC;
1296 char* atomD;
1297
1298 CubicTorsion* cTors;
1299 Atom** the_atoms;
1300 the_atoms = entry_plug->atoms;
1301
1302 // initialize the Torsions
1303
1304 for( i=0; i<nTorsions; i++ ){
1305
1306 a = the_torsions[i].a;
1307 b = the_torsions[i].b;
1308 c = the_torsions[i].c;
1309 d = the_torsions[i].d;
1310
1311 atomA = the_atoms[a]->getType();
1312 atomB = the_atoms[b]->getType();
1313 atomC = the_atoms[c]->getType();
1314 atomD = the_atoms[d]->getType();
1315 currentTorsionType = headTorsionType->find( atomA, atomB, atomC, atomD );
1316 if( currentTorsionType == NULL ){
1317 sprintf( painCave.errMsg,
1318 "TorsionType error, %s - %s - %s - %s not found"
1319 " in force file.\n",
1320 atomA, atomB, atomC, atomD );
1321 painCave.isFatal = 1;
1322 simError();
1323 }
1324
1325 if( !strcmp( currentTorsionType->type, "cubic" ) ){
1326
1327 cTors = new CubicTorsion( *the_atoms[a], *the_atoms[b],
1328 *the_atoms[c], *the_atoms[d] );
1329 cTors->setConstants( currentTorsionType->k1, currentTorsionType->k2,
1330 currentTorsionType->k3, currentTorsionType->k4 );
1331 torsionArray[i] = cTors;
1332 }
1333 }
1334 }
1335
1336 void TraPPE_ExFF::fastForward( char* stopText, char* searchOwner ){
1337
1338 int foundText = 0;
1339 char* the_token;
1340
1341 rewind( frcFile );
1342 lineNum = 0;
1343
1344 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1345 lineNum++;
1346 if( eof_test == NULL ){
1347 sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
1348 " file is empty.\n",
1349 searchOwner );
1350 painCave.isFatal = 1;
1351 simError();
1352 }
1353
1354
1355 while( !foundText ){
1356 while( eof_test != NULL && readLine[0] != '#' ){
1357 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1358 lineNum++;
1359 }
1360 if( eof_test == NULL ){
1361 sprintf( painCave.errMsg,
1362 "Error fast forwarding force file for %s at "
1363 "line %d: file ended unexpectedly.\n",
1364 searchOwner,
1365 lineNum );
1366 painCave.isFatal = 1;
1367 simError();
1368 }
1369
1370 the_token = strtok( readLine, " ,;\t#\n" );
1371 foundText = !strcmp( stopText, the_token );
1372
1373 if( !foundText ){
1374 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1375 lineNum++;
1376
1377 if( eof_test == NULL ){
1378 sprintf( painCave.errMsg,
1379 "Error fast forwarding force file for %s at "
1380 "line %d: file ended unexpectedly.\n",
1381 searchOwner,
1382 lineNum );
1383 painCave.isFatal = 1;
1384 simError();
1385 }
1386 }
1387 }
1388 }
1389
1390
1391 int TPE::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1392
1393 char* the_token;
1394
1395 the_token = strtok( lineBuffer, " \n\t,;" );
1396 if( the_token != NULL ){
1397
1398 strcpy( info.name, the_token );
1399
1400 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1401 sprintf( painCave.errMsg,
1402 "Error parseing AtomTypes: line %d\n", lineNum );
1403 painCave.isFatal = 1;
1404 simError();
1405 }
1406
1407 info.isDipole = atoi( the_token );
1408
1409 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1410 sprintf( painCave.errMsg,
1411 "Error parseing AtomTypes: line %d\n", lineNum );
1412 painCave.isFatal = 1;
1413 simError();
1414 }
1415
1416 info.isSSD = atoi( the_token );
1417
1418 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1419 sprintf( painCave.errMsg,
1420 "Error parseing AtomTypes: line %d\n", lineNum );
1421 painCave.isFatal = 1;
1422 simError();
1423 }
1424
1425 info.mass = atof( the_token );
1426
1427 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1428 sprintf( painCave.errMsg,
1429 "Error parseing AtomTypes: line %d\n", lineNum );
1430 painCave.isFatal = 1;
1431 simError();
1432 }
1433
1434 info.epslon = atof( the_token );
1435
1436 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1437 sprintf( painCave.errMsg,
1438 "Error parseing AtomTypes: line %d\n", lineNum );
1439 painCave.isFatal = 1;
1440 simError();
1441 }
1442
1443 info.sigma = atof( the_token );
1444
1445 if( info.isDipole ){
1446
1447 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1448 sprintf( painCave.errMsg,
1449 "Error parseing AtomTypes: line %d\n", lineNum );
1450 painCave.isFatal = 1;
1451 simError();
1452 }
1453
1454 info.dipole = atof( the_token );
1455 }
1456 else info.dipole = 0.0;
1457
1458 if( info.isSSD ){
1459
1460 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1461 sprintf( painCave.errMsg,
1462 "Error parseing AtomTypes: line %d\n", lineNum );
1463 painCave.isFatal = 1;
1464 simError();
1465 }
1466
1467 info.w0 = atof( the_token );
1468
1469 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1470 sprintf( painCave.errMsg,
1471 "Error parseing AtomTypes: line %d\n", lineNum );
1472 painCave.isFatal = 1;
1473 simError();
1474 }
1475
1476 info.v0 = atof( the_token );
1477 }
1478 else info.v0 = info.w0 = 0.0;
1479
1480 return 1;
1481 }
1482 else return 0;
1483 }
1484
1485 int TPE::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1486
1487 char* the_token;
1488
1489 the_token = strtok( lineBuffer, " \n\t,;" );
1490 if( the_token != NULL ){
1491
1492 strcpy( info.nameA, the_token );
1493
1494 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1495 sprintf( painCave.errMsg,
1496 "Error parseing BondTypes: line %d\n", lineNum );
1497 painCave.isFatal = 1;
1498 simError();
1499 }
1500
1501 strcpy( info.nameB, the_token );
1502
1503 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1504 sprintf( painCave.errMsg,
1505 "Error parseing BondTypes: line %d\n", lineNum );
1506 painCave.isFatal = 1;
1507 simError();
1508 }
1509
1510 strcpy( info.type, the_token );
1511
1512 if( !strcmp( info.type, "fixed" ) ){
1513 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1514 sprintf( painCave.errMsg,
1515 "Error parseing BondTypes: line %d\n", lineNum );
1516 painCave.isFatal = 1;
1517 simError();
1518 }
1519
1520 info.d0 = atof( the_token );
1521 }
1522 else{
1523 sprintf( painCave.errMsg,
1524 "Unknown TraPPE_Ex bond type \"%s\" at line %d\n",
1525 info.type,
1526 lineNum );
1527 painCave.isFatal = 1;
1528 simError();
1529 }
1530
1531 return 1;
1532 }
1533 else return 0;
1534 }
1535
1536
1537 int TPE::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1538
1539 char* the_token;
1540
1541 the_token = strtok( lineBuffer, " \n\t,;" );
1542 if( the_token != NULL ){
1543
1544 strcpy( info.nameA, the_token );
1545
1546 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1547 sprintf( painCave.errMsg,
1548 "Error parseing BendTypes: line %d\n", lineNum );
1549 painCave.isFatal = 1;
1550 simError();
1551 }
1552
1553 strcpy( info.nameB, the_token );
1554
1555 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1556 sprintf( painCave.errMsg,
1557 "Error parseing BendTypes: line %d\n", lineNum );
1558 painCave.isFatal = 1;
1559 simError();
1560 }
1561
1562 strcpy( info.nameC, the_token );
1563
1564 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1565 sprintf( painCave.errMsg,
1566 "Error parseing BendTypes: line %d\n", lineNum );
1567 painCave.isFatal = 1;
1568 simError();
1569 }
1570
1571 strcpy( info.type, the_token );
1572
1573 if( !strcmp( info.type, "quadratic" ) ){
1574 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1575 sprintf( painCave.errMsg,
1576 "Error parseing BendTypes: line %d\n", lineNum );
1577 painCave.isFatal = 1;
1578 simError();
1579 }
1580
1581 info.k1 = atof( the_token );
1582
1583 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1584 sprintf( painCave.errMsg,
1585 "Error parseing BendTypes: line %d\n", lineNum );
1586 painCave.isFatal = 1;
1587 simError();
1588 }
1589
1590 info.k2 = atof( the_token );
1591
1592 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1593 sprintf( painCave.errMsg,
1594 "Error parseing BendTypes: line %d\n", lineNum );
1595 painCave.isFatal = 1;
1596 simError();
1597 }
1598
1599 info.k3 = atof( the_token );
1600
1601 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1602 sprintf( painCave.errMsg,
1603 "Error parseing BendTypes: line %d\n", lineNum );
1604 painCave.isFatal = 1;
1605 simError();
1606 }
1607
1608 info.t0 = atof( the_token );
1609 }
1610
1611 else{
1612 sprintf( painCave.errMsg,
1613 "Unknown TraPPE_Ex bend type \"%s\" at line %d\n",
1614 info.type,
1615 lineNum );
1616 painCave.isFatal = 1;
1617 simError();
1618 }
1619
1620 return 1;
1621 }
1622 else return 0;
1623 }
1624
1625 int TPE::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1626
1627 char* the_token;
1628
1629 the_token = strtok( lineBuffer, " \n\t,;" );
1630 if( the_token != NULL ){
1631
1632 strcpy( info.nameA, the_token );
1633
1634 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1635 sprintf( painCave.errMsg,
1636 "Error parseing TorsionTypes: line %d\n", lineNum );
1637 painCave.isFatal = 1;
1638 simError();
1639 }
1640
1641 strcpy( info.nameB, the_token );
1642
1643 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1644 sprintf( painCave.errMsg,
1645 "Error parseing TorsionTypes: line %d\n", lineNum );
1646 painCave.isFatal = 1;
1647 simError();
1648 }
1649
1650 strcpy( info.nameC, the_token );
1651
1652 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1653 sprintf( painCave.errMsg,
1654 "Error parseing TorsionTypes: line %d\n", lineNum );
1655 painCave.isFatal = 1;
1656 simError();
1657 }
1658
1659 strcpy( info.nameD, the_token );
1660
1661 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1662 sprintf( painCave.errMsg,
1663 "Error parseing TorsionTypes: line %d\n", lineNum );
1664 painCave.isFatal = 1;
1665 simError();
1666 }
1667
1668 strcpy( info.type, the_token );
1669
1670 if( !strcmp( info.type, "cubic" ) ){
1671 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1672 sprintf( painCave.errMsg,
1673 "Error parseing TorsionTypes: line %d\n", lineNum );
1674 painCave.isFatal = 1;
1675 simError();
1676 }
1677
1678 info.k1 = atof( the_token );
1679
1680 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1681 sprintf( painCave.errMsg,
1682 "Error parseing TorsionTypes: line %d\n", lineNum );
1683 painCave.isFatal = 1;
1684 simError();
1685 }
1686
1687 info.k2 = atof( the_token );
1688
1689 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1690 sprintf( painCave.errMsg,
1691 "Error parseing TorsionTypes: line %d\n", lineNum );
1692 painCave.isFatal = 1;
1693 simError();
1694 }
1695
1696 info.k3 = atof( the_token );
1697
1698 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1699 sprintf( painCave.errMsg,
1700 "Error parseing TorsionTypes: line %d\n", lineNum );
1701 painCave.isFatal = 1;
1702 simError();
1703 }
1704
1705 info.k4 = atof( the_token );
1706
1707 }
1708
1709 else{
1710 sprintf( painCave.errMsg,
1711 "Unknown TraPPE_Ex torsion type \"%s\" at line %d\n",
1712 info.type,
1713 lineNum );
1714 painCave.isFatal = 1;
1715 simError();
1716 }
1717
1718 return 1;
1719 }
1720
1721 else return 0;
1722 }