ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/DUFF.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/DUFF.cpp (file contents):
Revision 561 by mmeineke, Fri Jun 20 20:29:36 2003 UTC vs.
Revision 1097 by gezelter, Mon Apr 12 20:32:20 2004 UTC

# Line 1 | Line 1
1 < #include <cstdlib>
2 < #include <cstdio>
3 < #include <cstring>
1 > #include <stdlib.h>
2 > #include <stdio.h>
3 > #include <string.h>
4  
5   #include <iostream>
6   using namespace std;
# Line 15 | Line 15 | namespace TPE {  // restrict the access of the folowin
15   #include "mpiForceField.h"
16   #endif // is_mpi
17  
18 – namespace TPE {  // restrict the access of the folowing to this file only.
18  
19 + // define some bond Types
20  
21 + #define FIXED_BOND    0
22 + #define HARMONIC_BOND 1
23 +
24 +
25 + namespace DUFF_NS {  // restrict the access of the folowing to this file only.
26 +
27 +
28    // Declare the structures that will be passed by MPI
29    
30    typedef struct{
# Line 25 | Line 32 | namespace TPE {  // restrict the access of the folowin
32      double mass;
33      double epslon;
34      double sigma;
35 +    double charge;
36      double dipole;
37      double w0;
38      double v0;
39 +    double v0p;
40 +    double rl;
41 +    double ru;
42 +    double rlp;
43 +    double rup;
44      int isSSD;
45 +    int isCharge;
46      int isDipole;
47      int ident;
48      int last;      //  0  -> default
# Line 39 | Line 53 | namespace TPE {  // restrict the access of the folowin
53    typedef struct{
54      char nameA[15];
55      char nameB[15];
42 –    char type[30];
56      double d0;
57 +    double k0;
58      int last;      //  0  -> default
59                     //  1  -> tells nodes to stop listening
60 +    int type;
61    } bondStruct;
62    
63    
# Line 101 | Line 116 | namespace TPE {  // restrict the access of the folowin
116      void printMe( void ){
117        
118        std::cerr << "LinkedAtype " << name << ": ident = " << ident << "\n";
119 <      if( next != NULL ) next->printMe();
119 >      //      if( next != NULL ) next->printMe();
120  
121      }
122  
# Line 130 | Line 145 | namespace TPE {  // restrict the access of the folowin
145          next->dipole   = info.dipole;
146          next->w0       = info.w0;
147          next->v0       = info.v0;
148 +        next->v0p      = info.v0p;
149 +        next->rl       = info.rl;
150 +        next->ru       = info.ru;
151 +        next->rlp      = info.rlp;
152 +        next->rup      = info.rup;
153          next->ident    = info.ident;
154        }
155      }
# Line 146 | Line 166 | namespace TPE {  // restrict the access of the folowin
166        info.dipole   = dipole;
167        info.w0       = w0;
168        info.v0       = v0;
169 +      info.v0p      = v0p;
170 +      info.rl       = rl;
171 +      info.ru       = ru;
172 +      info.rlp      = rlp;
173 +      info.rup      = rup;
174        info.ident    = ident;
175        info.last     = 0;
176      }
# Line 162 | Line 187 | namespace TPE {  // restrict the access of the folowin
187      double dipole;
188      double w0;
189      double v0;
190 +    double v0p;
191 +    double rl;
192 +    double ru;
193 +    double rlp;
194 +    double rup;
195      int ident;
196      LinkedAtomType* next;
197    };
# Line 172 | Line 202 | namespace TPE {  // restrict the access of the folowin
202        next = NULL;
203        nameA[0] = '\0';
204        nameB[0] = '\0';
175 –      type[0] = '\0';
205      }
206      ~LinkedBondType(){ if( next != NULL ) delete next; }
207  
# Line 207 | Line 236 | namespace TPE {  // restrict the access of the folowin
236          next = new LinkedBondType();
237          strcpy(next->nameA, info.nameA);
238          strcpy(next->nameB, info.nameB);
239 <        strcpy(next->type,  info.type);
239 >        next->type = info.type;
240          next->d0 = info.d0;
241 +        next->k0 = info.k0;
242        }
243      }
244      
# Line 216 | Line 246 | namespace TPE {  // restrict the access of the folowin
246      void duplicate( bondStruct &info ){
247        strcpy(info.nameA, nameA);
248        strcpy(info.nameB, nameB);
249 <      strcpy(info.type,  type);
249 >      info.type = type;
250        info.d0   = d0;
251 +      info.k0   = k0;
252        info.last = 0;
253      }
254  
# Line 226 | Line 257 | namespace TPE {  // restrict the access of the folowin
257  
258      char nameA[15];
259      char nameB[15];
260 <    char type[30];
260 >    int type;
261      double d0;
262 +    double k0;
263  
264      LinkedBondType* next;
265    };
# Line 410 | Line 442 | using namespace TPE;
442  
443   } // namespace
444  
445 < using namespace TPE;
445 > using namespace DUFF_NS;
446  
447  
448   //****************************************************************
# Line 424 | Line 456 | DUFF::DUFF(){
456    char* ffPath_env = "FORCE_PARAM_PATH";
457    char* ffPath;
458    char temp[200];
427 –  char errMsg[1000];
459  
460    headAtomType       = NULL;
461    currentAtomType    = NULL;
# Line 446 | Line 477 | DUFF::DUFF(){
477    // Init the atomStruct mpi type
478  
479    atomStruct atomProto; // mpiPrototype
480 <  int atomBC[3] = {15,6,4};  // block counts
480 >  int atomBC[3] = {15,12,5};  // block counts
481    MPI_Aint atomDspls[3];           // displacements
482    MPI_Datatype atomMbrTypes[3];    // member mpi types
483  
# Line 468 | Line 499 | DUFF::DUFF(){
499    // Init the bondStruct mpi type
500    
501    bondStruct bondProto; // mpiPrototype
502 <  int bondBC[3] = {60,1,1};  // block counts
502 >  int bondBC[3] = {30,2,2};  // block counts
503    MPI_Aint bondDspls[3];           // displacements
504    MPI_Datatype bondMbrTypes[3];    // member mpi types
505    
# Line 628 | Line 659 | void DUFF::readParams( void ){
659  
660   void DUFF::readParams( void ){
661  
631 –  int i, a, b, c, d;
662    int identNum;
633 –  char* atomA;
634 –  char* atomB;
635 –  char* atomC;
636 –  char* atomD;
663    
664    atomStruct atomInfo;
665    bondStruct bondInfo;
# Line 705 | Line 731 | void DUFF::readParams( void ){
731      while( currentAtomType != NULL ){
732        currentAtomType->duplicate( atomInfo );
733  
708 –
709 –
734        sendFrcStruct( &atomInfo, mpiAtomStructType );
735  
736        sprintf( checkPointMsg,
# Line 724 | Line 748 | void DUFF::readParams( void ){
748    else{
749      
750      // listen for node 0 to send out the force params
751 <    
751 >
752      MPIcheckPoint();
753  
754      headAtomType = new LinkedAtomType;
755 <    recieveFrcStruct( &atomInfo, mpiAtomStructType );
755 >    receiveFrcStruct( &atomInfo, mpiAtomStructType );
756      
757      while( !atomInfo.last ){
758  
735 –
736 –
759        headAtomType->add( atomInfo );
760        
761        MPIcheckPoint();
762  
763 <      recieveFrcStruct( &atomInfo, mpiAtomStructType );
763 >      receiveFrcStruct( &atomInfo, mpiAtomStructType );
764      }
765    }
766  
# Line 754 | Line 776 | void DUFF::readParams( void ){
776    
777    int isGB = 0;
778    int isLJ = 1;
779 <  double GB_dummy = 0.0;
780 <  
781 <  
779 >  int isEAM =0;
780 >  int isCharge = 0;
781 >  double charge=0.0;
782 >    
783    currentAtomType = headAtomType->next;;
784    while( currentAtomType != NULL ){
785      
786 <    if(currentAtomType->isDipole) entry_plug->useDipole = 1;
786 >    if(currentAtomType->isDipole) entry_plug->useDipoles = 1;
787      if(currentAtomType->isSSD) {
788        entry_plug->useSticky = 1;
789 <      set_sticky_params( &(currentAtomType->w0), &(currentAtomType->v0));
789 >      set_sticky_params( &(currentAtomType->w0), &(currentAtomType->v0),
790 >                         &(currentAtomType->v0p),
791 >                         &(currentAtomType->rl), &(currentAtomType->ru),
792 >                         &(currentAtomType->rlp), &(currentAtomType->rup));
793      }
794  
795      if( currentAtomType->name[0] != '\0' ){
# Line 773 | Line 799 | void DUFF::readParams( void ){
799                   &(currentAtomType->isSSD),
800                   &(currentAtomType->isDipole),
801                   &isGB,
802 +                 &isEAM,
803 +                 &isCharge,
804                   &(currentAtomType->epslon),
805                   &(currentAtomType->sigma),
806 +                 &charge,
807                   &(currentAtomType->dipole),
808                   &isError );
809        if( isError ){
# Line 865 | Line 894 | void DUFF::readParams( void ){
894      MPIcheckPoint();
895  
896      headBondType = new LinkedBondType;
897 <    recieveFrcStruct( &bondInfo, mpiBondStructType );
897 >    receiveFrcStruct( &bondInfo, mpiBondStructType );
898      while( !bondInfo.last ){
899  
900        headBondType->add( bondInfo );
901 <      recieveFrcStruct( &bondInfo, mpiBondStructType );
901 >      receiveFrcStruct( &bondInfo, mpiBondStructType );
902      }
903    }
904  
# Line 893 | Line 922 | void DUFF::readParams( void ){
922      fastForward( "BendTypes", "initializeBends" );
923  
924      // we are now at the bendTypes section
925 <
925 >    
926      eof_test =  fgets( readLine, sizeof(readLine), frcFile );
927      lineNum++;
928          
# Line 948 | Line 977 | void DUFF::readParams( void ){
977      MPIcheckPoint();
978  
979      headBendType = new LinkedBendType;
980 <    recieveFrcStruct( &bendInfo, mpiBendStructType );
980 >    receiveFrcStruct( &bendInfo, mpiBendStructType );
981      while( !bendInfo.last ){
982  
983        headBendType->add( bendInfo );
984 <      recieveFrcStruct( &bendInfo, mpiBendStructType );
984 >      receiveFrcStruct( &bendInfo, mpiBendStructType );
985      }
986    }
987  
# Line 1033 | Line 1062 | void DUFF::readParams( void ){
1062      MPIcheckPoint();
1063  
1064      headTorsionType = new LinkedTorsionType;
1065 <    recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1065 >    receiveFrcStruct( &torsionInfo, mpiTorsionStructType );
1066      while( !torsionInfo.last ){
1067  
1068        headTorsionType->add( torsionInfo );
1069 <      recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1069 >      receiveFrcStruct( &torsionInfo, mpiTorsionStructType );
1070      }
1071    }
1072  
# Line 1091 | Line 1120 | void DUFF::initializeAtoms( int nAtoms, Atom** the_ato
1120    // initialize the atoms
1121    
1122    DirectionalAtom* dAtom;
1123 +  double ji[3];
1124  
1125    for(int i=0; i<nAtoms; i++ ){
1126  
# Line 1104 | Line 1134 | void DUFF::initializeAtoms( int nAtoms, Atom** the_ato
1134      }
1135      
1136      the_atoms[i]->setMass( currentAtomType->mass );
1107 –    the_atoms[i]->setEpslon( currentAtomType->epslon );
1108 –    the_atoms[i]->setSigma( currentAtomType->sigma );
1137      the_atoms[i]->setIdent( currentAtomType->ident );
1110 –    the_atoms[i]->setLJ();
1138  
1139      if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
1140  
# Line 1115 | Line 1142 | void DUFF::initializeAtoms( int nAtoms, Atom** the_ato
1142        if( the_atoms[i]->isDirectional() ){
1143          
1144          dAtom = (DirectionalAtom *) the_atoms[i];
1118 –        dAtom->setMu( currentAtomType->dipole );
1145          dAtom->setHasDipole( 1 );
1146 <        dAtom->setJx( 0.0 );
1147 <        dAtom->setJy( 0.0 );
1148 <        dAtom->setJz( 0.0 );
1146 >
1147 >        ji[0] = 0.0;
1148 >        ji[1] = 0.0;
1149 >        ji[2] = 0.0;
1150 >
1151 >        dAtom->setJ( ji );
1152          
1153          if(!strcmp("SSD",the_atoms[i]->getType())){
1154            dAtom->setI( waterI );
1126 –          dAtom->setSSD( 1 );
1155          }
1156          else if(!strcmp("HEAD",the_atoms[i]->getType())){
1157            dAtom->setI( headI );
1130 –          dAtom->setSSD( 0 );
1158          }
1159          else{
1160            sprintf(painCave.errMsg,
# Line 1151 | Line 1178 | void DUFF::initializeAtoms( int nAtoms, Atom** the_ato
1178      else{
1179        if( the_atoms[i]->isDirectional() ){
1180          sprintf( painCave.errMsg,
1181 <                 "DUFF error: Atom \"%s\" was given a standard"
1181 >                 "DUFF error: Atom \"%s\" was given a standard "
1182                   "orientation in the BASS file, yet it is not a dipole.\n",
1183                   currentAtomType->name);
1184          painCave.isFatal = 1;
# Line 1189 | Line 1216 | void DUFF::initializeBonds( int nBonds, Bond** bondArr
1216        simError();
1217      }
1218      
1219 <    if( !strcmp( currentBondType->type, "fixed" ) ){
1220 <      
1219 >    switch( currentBondType->type ){
1220 >
1221 >    case FIXED_BOND:
1222 >            
1223        bondArray[i] = new ConstrainedBond( *the_atoms[a],
1224                                            *the_atoms[b],
1225                                            currentBondType->d0 );
1226        entry_plug->n_constraints++;
1227 +      break;
1228 +
1229 +    case HARMONIC_BOND:
1230 +      
1231 +      bondArray[i] = new HarmonicBond( *the_atoms[a],
1232 +                                       *the_atoms[b],
1233 +                                       currentBondType->d0,
1234 +                                       currentBondType->k0 );
1235 +      break;
1236 +      
1237 +    default:
1238 +
1239 +      break;
1240 +      // do nothing
1241      }
1242    }
1243   }
# Line 1259 | Line 1302 | void DUFF::initializeBends( int nBends, Bend** bendArr
1302          }
1303          
1304          gBend = new GhostBend( *the_atoms[a],
1305 <                               *the_atoms[b] );                        
1305 >                               *the_atoms[b]);
1306 >                                                                      
1307          gBend->setConstants( currentBendType->k1,
1308                               currentBendType->k2,
1309                               currentBendType->k3,
# Line 1275 | Line 1319 | void DUFF::initializeBends( int nBends, Bend** bendArr
1319                               currentBendType->k3,
1320                               currentBendType->t0 );
1321          bendArray[i] = qBend;
1322 <      }
1322 >      }      
1323      }
1324    }
1325   }
# Line 1382 | Line 1426 | int TPE::parseAtom( char *lineBuffer, int lineNum, ato
1426   }
1427  
1428  
1429 < int TPE::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1429 > int DUFF_NS::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1430  
1431    char* the_token;
1432    
# Line 1468 | Line 1512 | int TPE::parseAtom( char *lineBuffer, int lineNum, ato
1512        }
1513        
1514        info.v0 = atof( the_token );
1515 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1516 +        sprintf( painCave.errMsg,
1517 +                 "Error parseing AtomTypes: line %d\n", lineNum );
1518 +        painCave.isFatal = 1;
1519 +        simError();
1520 +      }
1521 +      
1522 +      info.v0p = atof( the_token );
1523 +
1524 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1525 +        sprintf( painCave.errMsg,
1526 +                 "Error parseing AtomTypes: line %d\n", lineNum );
1527 +        painCave.isFatal = 1;
1528 +        simError();
1529 +      }
1530 +      
1531 +      info.rl = atof( the_token );
1532 +
1533 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1534 +        sprintf( painCave.errMsg,
1535 +                 "Error parseing AtomTypes: line %d\n", lineNum );
1536 +        painCave.isFatal = 1;
1537 +        simError();
1538 +      }
1539 +      
1540 +      info.ru = atof( the_token );
1541 +
1542 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1543 +        sprintf( painCave.errMsg,
1544 +                 "Error parseing AtomTypes: line %d\n", lineNum );
1545 +        painCave.isFatal = 1;
1546 +        simError();
1547 +      }
1548 +      
1549 +      info.rlp = atof( the_token );
1550 +
1551 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1552 +        sprintf( painCave.errMsg,
1553 +                 "Error parseing AtomTypes: line %d\n", lineNum );
1554 +        painCave.isFatal = 1;
1555 +        simError();
1556 +      }
1557 +      
1558 +      info.rup = atof( the_token );
1559      }
1560 <    else info.v0 = info.w0 = 0.0;
1560 >    else info.v0 = info.w0 = info.v0p = info.rl = info.ru = info.rlp = info.rup = 0.0;
1561  
1562      return 1;
1563    }
1564    else return 0;
1565   }
1566  
1567 < int TPE::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1567 > int DUFF_NS::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1568  
1569    char* the_token;
1570 +  char bondType[30];
1571    
1572    the_token = strtok( lineBuffer, " \n\t,;" );
1573    if( the_token != NULL ){
# Line 1501 | Line 1590 | int TPE::parseBond( char *lineBuffer, int lineNum, bon
1590        simError();
1591      }
1592      
1593 <    strcpy( info.type, the_token );
1593 >    strcpy( bondType, the_token );
1594      
1595 <    if( !strcmp( info.type, "fixed" ) ){
1595 >    if( !strcmp( bondType, "fixed" ) ){
1596 >      info.type = FIXED_BOND;
1597 >      
1598 >      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1599 >        sprintf( painCave.errMsg,
1600 >                 "Error parseing BondTypes: line %d\n", lineNum );
1601 >        painCave.isFatal = 1;
1602 >        simError();
1603 >      }
1604 >      
1605 >      info.d0 = atof( the_token );
1606 >      
1607 >      info.k0=0.0;
1608 >    }
1609 >    else if( !strcmp( bondType, "harmonic" ) ){
1610 >      info.type = HARMONIC_BOND;
1611 >      
1612        if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1613          sprintf( painCave.errMsg,
1614                   "Error parseing BondTypes: line %d\n", lineNum );
# Line 1512 | Line 1617 | int TPE::parseBond( char *lineBuffer, int lineNum, bon
1617        }
1618        
1619        info.d0 = atof( the_token );
1620 +
1621 +      if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1622 +        sprintf( painCave.errMsg,
1623 +                 "Error parseing BondTypes: line %d\n", lineNum );
1624 +        painCave.isFatal = 1;
1625 +        simError();
1626 +      }
1627 +      
1628 +      info.k0 = atof( the_token );
1629      }
1630 +
1631      else{
1632        sprintf( painCave.errMsg,
1633                 "Unknown DUFF bond type \"%s\" at line %d\n",
1634 <               info.type,
1634 >               bondType,
1635                 lineNum );
1636        painCave.isFatal = 1;
1637        simError();
# Line 1528 | Line 1643 | int TPE::parseBend( char *lineBuffer, int lineNum, ben
1643   }
1644  
1645  
1646 < int TPE::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1646 > int DUFF_NS::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1647  
1648    char* the_token;
1649    
# Line 1616 | Line 1731 | int TPE::parseTorsion( char *lineBuffer, int lineNum,
1731    else return 0;
1732   }
1733  
1734 < int TPE::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1734 > int DUFF_NS::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1735    
1736    char*  the_token;
1737  

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines