ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-3.0/src/UseTheForce/ForceField.cpp
(Generate patch)

Comparing branches/new_design/OOPSE-3.0/src/UseTheForce/ForceField.cpp (file contents):
Revision 1716, Fri Nov 5 21:04:13 2004 UTC vs.
Revision 1783 by tim, Wed Nov 24 21:33:00 2004 UTC

# Line 1 | Line 1
1 + /*
2 + * Copyright (C) 2000-2004  Object Oriented Parallel Simulation Engine (OOPSE) project
3 + *
4 + * Contact: oopse@oopse.org
5 + *
6 + * This program is free software; you can redistribute it and/or
7 + * modify it under the terms of the GNU Lesser General Public License
8 + * as published by the Free Software Foundation; either version 2.1
9 + * of the License, or (at your option) any later version.
10 + * All we ask is that proper credit is given for our work, which includes
11 + * - but is not limited to - adding the above copyright notice to the beginning
12 + * of your source code files, and to any copyright notice that you may distribute
13 + * with programs based on this work.
14 + *
15 + * This program is distributed in the hope that it will be useful,
16 + * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18 + * GNU Lesser General Public License for more details.
19 + *
20 + * You should have received a copy of the GNU Lesser General Public License
21 + * along with this program; if not, write to the Free Software
22 + * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
23 + *
24 + */
25 +
26 + /**
27 +  * @file ForceField.cpp
28 +  * @author tlin
29 +  * @date 11/04/2004
30 +  * @time 22:51am
31 +  * @version 1.0
32 +  */
33 +  
34   #include "UseTheForce/ForceField.hpp"
35 + #include "utils/simError.h"
36 + namespace oopse {
37  
38 < AtomType* ForceField::getMatchingAtomType(const string &at) {
38 > ForceField::ForceField() : hasVariant_(false){
39 >    char* tempPath;
40 >    tempPath = getenv("FORCE_PARAM_PATH");
41  
42 <  map<string, AtomType*>::iterator iter;
43 <  
44 <  iter = atomTypeMap.find(at);
45 <  if (iter != atomTypeMap.end()) {
46 <    return iter->second;
47 <  } else {
11 <    return NULL;
12 <  }
42 >    if (tempPath == NULL) {
43 >        //convert a macro from compiler to a string in c++
44 >        STR_DEFINE(ffPath_, FRC_PATH );
45 >    } else {
46 >        ffPath_ = tempPath;
47 >    }
48   }
49  
50 < BondType* ForceField::getMatchingBondType(const string &at1,
51 <                                          const string &at2) {
50 > AtomType* ForceField::getAtomType(const std::string &at) {
51 >    std::vector<std::string> keys;
52 >    keys.push_back(at);
53 >    return atomTypeCont_.find(keys);
54 > }
55  
56 <  map<pair<string,string>, BondType*>::iterator iter;
57 <  vector<BondType*> foundTypes;
56 > BondType* ForceField::getBondType(const std::string &at1, const std::string &at2) {
57 >    std::vector<std::string> keys;
58 >    keys.push_back(at1);
59 >    keys.push_back(at2);    
60 >    return bondTypeCont_.find(keys, wildCardAtomTypeName_);
61  
62 <  iter = bondTypeMap.find(pair<at1, at2>);
22 <  if (iter != bondTypeMap.end()) {
23 <    // exact match, so just return it
24 <    return iter->second;
25 <  }
62 > }
63  
64 <  iter = bondTypeMap.find(pair<at2, at1>);
65 <  if (iter != bondTypeMap.end()) {
66 <    // exact match in reverse order, so just return it
67 <    return iter->second;
68 <  }
64 > BendType* ForceField::getBendType(const std::string &at1, const std::string &at2,
65 >                        const std::string &at3) {
66 >    std::vector<std::string> keys;
67 >    keys.push_back(at1);
68 >    keys.push_back(at2);    
69 >    keys.push_back(at3);    
70 >    return bendTypeCont_.find(keys, wildCardAtomTypeName_);
71 > }
72  
73 <  iter = bondTypeMap.find(pair<at1, wildCardAtomTypeName>);
74 <  if (iter != bondTypeMap.end()) {
75 <    foundTypes.push_back(iter->second);
76 <  }
73 > TorsionType* ForceField::getTorsionType(const std::string &at1, const std::string &at2,
74 >                              const std::string &at3, const std::string &at4) {
75 >    std::vector<std::string> keys;
76 >    keys.push_back(at1);
77 >    keys.push_back(at2);    
78 >    keys.push_back(at3);    
79 >    keys.push_back(at4);    
80 >    return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
81  
82 <  iter = bondTypeMap.find(pair<at2, wildCardAtomTypeName>);
39 <  if (iter != bondTypeMap.end()) {
40 <    foundTypes.push_back(iter->second);
41 <  }
82 > }
83  
84 <  iter = bondTypeMap.find(pair<wildCardAtomTypeName, at1>);
85 <  if (iter != bondTypeMap.end()) {
86 <    foundTypes.push_back(iter->second);
87 <  }
84 > BondType* ForceField::getExactBondType(const std::string &at1, const std::string &at2){
85 >    std::vector<std::string> keys;
86 >    keys.push_back(at1);
87 >    keys.push_back(at2);    
88 >    return bondTypeCont_.find(keys);
89 > }
90  
91 <  iter = bondTypeMap.find(pair<wildCardAtomTypeName, at2>);
92 <  if (iter != bondTypeMap.end()) {
93 <    foundTypes.push_back(iter->second);
94 <  }
95 <  
96 <  if (foundTypes.empty()) {
97 <    return NULL;
98 <  } else {
56 <    
91 > BendType* ForceField::getExactBendType(const std::string &at1, const std::string &at2,
92 >                            const std::string &at3){
93 >    std::vector<std::string> keys;
94 >    keys.push_back(at1);
95 >    keys.push_back(at2);    
96 >    keys.push_back(at3);    
97 >    return bendTypeCont_.find(keys);
98 > }
99  
100 <
100 > TorsionType* ForceField::getExactTorsionType(const std::string &at1, const std::string &at2,
101 >                                  const std::string &at3, const std::string &at4){
102 >    std::vector<std::string> keys;
103 >    keys.push_back(at1);
104 >    keys.push_back(at2);    
105 >    keys.push_back(at3);    
106 >    keys.push_back(at4);  
107 >    return torsionTypeCont_.find(keys);
108 > }
109 > bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
110 >    std::vector<std::string> keys;
111 >    keys.push_back(at);
112 >    return atomTypeCont_.add(keys, atomType);
113 > }
114  
115 + bool ForceField::addBondType(const std::string &at1, const std::string &at2, BondType* bondType) {
116 +    std::vector<std::string> keys;
117 +    keys.push_back(at1);
118 +    keys.push_back(at2);    
119 +    return bondTypeCont_.add(keys, bondType);
120  
121 + }
122  
123 <  
123 > bool ForceField::addBendType(const std::string &at1, const std::string &at2,
124 >                        const std::string &at3, BendType* bendType) {
125 >    std::vector<std::string> keys;
126 >    keys.push_back(at1);
127 >    keys.push_back(at2);    
128 >    keys.push_back(at3);    
129 >    return bendTypeCont_.add(keys, bendType);
130 > }
131  
132 + bool ForceField::addTorsionType(const std::string &at1, const std::string &at2,
133 +                              const std::string &at3, const std::string &at4, TorsionType* torsionType) {
134 +    std::vector<std::string> keys;
135 +    keys.push_back(at1);
136 +    keys.push_back(at2);    
137 +    keys.push_back(at3);    
138 +    keys.push_back(at4);    
139 +    return torsionTypeCont_.add(keys, torsionType);
140 + }
141  
142 < BendType* ForceField::getMatchingBendType(const string &at1, const string &at2,
143 <                                          const string &at3);
144 < TorsionType* ForceField::getMatchingTorsionType(const string &at1, const string &at2,
145 <                                                const string &at3, const string &at4);
142 > double ForceField::getRcutFromAtomType(AtomType* at) {
143 >    /**@todo */
144 >    GenericData* data;
145 >    double rcut = 0.0;
146  
147 < double ForceField::getRcutForAtomType(AtomType* at);
147 >    if (at->isLennardJones()) {
148 >        data = at->getPropertyByName("LennardJones");
149 >        if (data != NULL) {
150 >            LJParamGenericData* ljData = dynamic_cast<LJParamGenericData*>(data);
151  
152 +            if (ljData != NULL) {
153 +                LJParam ljParam = ljData->getData();
154  
155 < vector<vector<string> > generateWildcardSequence(const vector<string> atomTypes) {
156 <  
157 <   vector<vector<string> > results;
155 >                //by default use 2.5*sigma as cutoff radius
156 >                rcut = 2.5 * ljParam.sigma;
157 >                
158 >            } else {
159 >                    sprintf( painCave.errMsg,
160 >                           "Can not cast GenericData to LJParam\n");
161 >                    painCave.severity = OOPSE_ERROR;
162 >                    painCave.isFatal = 1;
163 >                    simError();          
164 >            }            
165 >        } else {
166 >            sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
167 >            painCave.severity = OOPSE_ERROR;
168 >            painCave.isFatal = 1;
169 >            simError();          
170 >        }
171 >    }
172  
173 <  
173 >    return rcut;    
174 > }
175  
176  
177 <   vector<vector< string> > getAllWildcardPermutations(const vector<string> myAts) {
178 <    
179 <     int nStrings;
180 <     vector<string> oneResult;
181 <     vector<vector<string> > allResults;
177 > ifstrstream* ForceField::openForceFieldFile(const std::string& filename) {
178 >    std::string forceFieldFilename(filename);
179 >    ifstrstream* ffStream = new ifstrstream();
180 >    
181 >    //try to open the force filed file in current directory first    
182 >    ffStream->open(forceFieldFilename.c_str());
183 >    if(!ffStream->is_open()){
184  
185 <     nStrings = myAts.size();
185 >        forceFieldFilename = ffPath_ + "/" + forceFieldFilename;
186 >        ffStream->open( forceFieldFilename.c_str() );
187  
188 <     if (nStrings == 1) {
189 <       oneResult.push_back(wildcardCharacter);
190 <       allResults.push_back(oneResult);
191 <       return allResults;
192 <     } else {
193 <      
194 <       for (i=0; i < nStrings; i++) {
195 <         oneResult = myAts;
196 <         replace(oneResult.begin(), oneResult.end(),
188 >        //if current directory does not contain the force field file,
189 >        //try to open it in the path        
190 >        if(!ffStream->is_open()){
191 >
192 >            sprintf( painCave.errMsg,
193 >               "Error opening the force field parameter file:\n"
194 >               "\t%s\n"
195 >               "\tHave you tried setting the FORCE_PARAM_PATH environment "
196 >               "variable?\n",
197 >               forceFieldFilename.c_str() );
198 >            painCave.severity = OOPSE_ERROR;
199 >            painCave.isFatal = 1;
200 >            simError();
201 >        }
202 >    }  
203 >
204 >    return ffStream;
205 >
206 > }
207 >
208 > } //end namespace oopse

Diff Legend

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