ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-2.0/src/UseTheForce/DUFF.cpp
Revision: 1670
Committed: Thu Oct 28 16:56:20 2004 UTC (19 years, 8 months ago) by gezelter
File size: 44243 byte(s)
Log Message:
fixed duplicate declaration foo

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