ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-3.0/src/types/MoleculeStamp.cpp
(Generate patch)

Comparing trunk/OOPSE-3.0/src/types/MoleculeStamp.cpp (file contents):
Revision 1490 by gezelter, Fri Sep 24 04:16:43 2004 UTC vs.
Revision 2469 by tim, Fri Dec 2 15:38:03 2005 UTC

# Line 1 | Line 1
1 + /*
2 + * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3 + *
4 + * The University of Notre Dame grants you ("Licensee") a
5 + * non-exclusive, royalty free, license to use, modify and
6 + * redistribute this software in source and binary code form, provided
7 + * that the following conditions are met:
8 + *
9 + * 1. Acknowledgement of the program authors must be made in any
10 + *    publication of scientific results based in part on use of the
11 + *    program.  An acceptable form of acknowledgement is citation of
12 + *    the article in which the program was described (Matthew
13 + *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14 + *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15 + *    Parallel Simulation Engine for Molecular Dynamics,"
16 + *    J. Comput. Chem. 26, pp. 252-271 (2005))
17 + *
18 + * 2. Redistributions of source code must retain the above copyright
19 + *    notice, this list of conditions and the following disclaimer.
20 + *
21 + * 3. Redistributions in binary form must reproduce the above copyright
22 + *    notice, this list of conditions and the following disclaimer in the
23 + *    documentation and/or other materials provided with the
24 + *    distribution.
25 + *
26 + * This software is provided "AS IS," without a warranty of any
27 + * kind. All express or implied conditions, representations and
28 + * warranties, including any implied warranty of merchantability,
29 + * fitness for a particular purpose or non-infringement, are hereby
30 + * excluded.  The University of Notre Dame and its licensors shall not
31 + * be liable for any damages suffered by licensee as a result of
32 + * using, modifying or distributing the software or its
33 + * derivatives. In no event will the University of Notre Dame or its
34 + * licensors be liable for any lost revenue, profit or data, or for
35 + * direct, indirect, special, consequential, incidental or punitive
36 + * damages, however caused and regardless of the theory of liability,
37 + * arising out of the use of or inability to use software, even if the
38 + * University of Notre Dame has been advised of the possibility of
39 + * such damages.
40 + */
41 +
42   #include <stdlib.h>
43   #include <stdio.h>
44   #include <string.h>
45   #include <iostream>
46  
47 < #include "MoleculeStamp.hpp"
47 > #include "types/MoleculeStamp.hpp"
48  
49 < char MoleculeStamp::errMsg[500];
49 > namespace oopse {
50 > MoleculeStamp::MoleculeStamp() {
51 >    DefineParameter(Name, "name");
52 >    
53 >    deprecatedKeywords_.insert("nAtoms");
54 >    deprecatedKeywords_.insert("nBonds");
55 >    deprecatedKeywords_.insert("nBends");
56 >    deprecatedKeywords_.insert("nTorsions");
57 >    deprecatedKeywords_.insert("nRigidBodies");
58 >    deprecatedKeywords_.insert("nCutoffGroups");
59 >    
60 > }
61  
62 < MoleculeStamp::MoleculeStamp(){
11 <  
12 <  n_atoms = 0;
13 <  n_bonds = 0;
14 <  n_bends = 0;
15 <  n_torsions = 0;
16 <  n_rigidbodies = 0;
17 <  n_cutoffgroups = 0;
18 <  n_integrable = 0;
62 > MoleculeStamp::~MoleculeStamp() {
63  
64 <  unhandled = NULL;
21 <  atoms = NULL;
22 <  bonds = NULL;
23 <  bends = NULL;
24 <  torsions = NULL;
25 <  rigidBodies = NULL;
26 <  cutoffGroups = NULL;
64 > }
65  
66 <  have_name = 0;
67 <  have_atoms = 0;
68 <  have_bonds = 0;
69 <  have_bends = 0;
70 <  have_torsions = 0;
71 <  have_rigidbodies = 0;
72 <  have_cutoffgroups = 0;
66 > bool MoleculeStamp::addAtomStamp( AtomStamp* atom) {
67 >    bool ret = addIndexSensitiveStamp(atomStamps_, atom);
68 >    if (!ret) {
69 >        std::cout << "multiple atoms have the same index: " << atom->getIndex() <<" in " << getName()  << " Molecule\n";
70 >    }
71 >    return ret;
72 >    
73 > }
74  
75 + bool MoleculeStamp::addBondStamp( BondStamp* bond) {
76 +    bondStamps_.push_back(bond);
77 +    return true;
78   }
79  
80 < MoleculeStamp::~MoleculeStamp(){
81 <  int i;
82 <  
83 <  if( unhandled != NULL) delete unhandled;
42 <  
43 <  if( rigidBodies != NULL ) {
44 <    for( i=0; i<n_rigidbodies; i++ ) delete rigidBodies[i];
45 <  }
46 <  delete[] rigidBodies;
80 > bool MoleculeStamp::addBendStamp( BendStamp* bend) {
81 >    bendStamps_.push_back(bend);
82 >    return true;
83 > }
84  
85 <  if( cutoffGroups != NULL ) {
86 <    for( i=0; i<n_cutoffgroups; i++ ) delete cutoffGroups[i];
87 <  }
88 <  delete[] cutoffGroups;
52 <  
53 <  if( atoms != NULL ){
54 <    for( i=0; i<n_atoms; i++ ) delete atoms[i];
55 <  }
56 <  delete[] atoms;
57 <  
58 <  if( bonds != NULL ){
59 <    for( i=0; i<n_bonds; i++ ) delete bonds[i];
60 <  }
61 <  delete[] bonds;
62 <  
63 <  if( bends != NULL ){
64 <    for( i=0; i<n_bends; i++ ) delete bends[i];
65 <  }
66 <  delete[] bends;
67 <  
68 <  if( torsions != NULL ){
69 <    for( i=0; i<n_torsions; i++ ) delete torsions[i];
70 <  }
71 <  delete[] torsions;
72 <  
85 > bool MoleculeStamp::addTorsionStamp( TorsionStamp* torsion) {
86 >    torsionStamps_.push_back(torsion);
87 >    return true;
88 > }
89  
90 + bool MoleculeStamp::addRigidBodyStamp( RigidBodyStamp* rigidbody) {
91 +    bool ret = addIndexSensitiveStamp(rigidBodyStamps_, rigidbody);
92 +    if (!ret) {
93 +        std::cout << "multiple rigidbodies have the same index: " << rigidbody->getIndex() <<" in " << getName()  << " Molecule\n";
94 +    }
95 +    return ret;
96 + }
97  
98 <
98 > bool MoleculeStamp::addCutoffGroupStamp( CutoffGroupStamp* cutoffgroup) {
99 >    cutoffGroupStamps_.push_back(cutoffgroup);
100 >    return true;
101   }
102  
103 < char* MoleculeStamp::assignString( char* lhs, char* rhs ){
104 <  
80 <  if( !strcmp( lhs, "name" ) ){
81 <    strcpy( name, rhs );
82 <    have_name = 1;
83 <  }
84 <  else{
85 <    if( unhandled == NULL ) unhandled = new LinkedAssign( lhs, rhs );
86 <    else unhandled->add( lhs, rhs );
87 <    have_extras = 1;
88 <  }
89 <  return NULL;
103 > bool MoleculeStamp::addFragmentStamp( FragmentStamp* fragment) {
104 >    return addIndexSensitiveStamp(fragmentStamps_, fragment);
105   }
106 +    
107 + void MoleculeStamp::validate() {
108 +    DataHolder::validate();
109  
110 < char* MoleculeStamp::assignDouble( char* lhs, double rhs ){
111 <  int i;
110 >    std::vector<AtomStamp*>::iterator ai = std::find(atomStamps_.begin(), atomStamps_.end(), static_cast<AtomStamp*>(NULL));
111 >    if (ai != atomStamps_.end()) {
112 >        std::cout << "Error in Molecule " << getName() << ": atom[" << ai - atomStamps_.begin()<< "] is missing\n";
113 >    }
114  
115 <  if( !strcmp( lhs, "nAtoms" ) ){
116 <    n_atoms = (int)rhs;
115 >     std::vector<RigidBodyStamp*>::iterator ri = std::find(rigidBodyStamps_.begin(), rigidBodyStamps_.end(), static_cast<RigidBodyStamp*>(NULL));
116 >     if (ri != rigidBodyStamps_.end()) {
117 >         std::cout << "Error in Molecule " << getName() << ":rigidBody[" <<  ri - rigidBodyStamps_.begin()<< "] is missing\n";
118 >     }
119      
120 <    if( have_atoms ){
121 <      sprintf( errMsg,
122 <               "MoleculeStamp error, n_atoms already declared"
101 <               " for molecule: %s\n",
102 <               name);
103 <      return strdup( errMsg );
120 >    std::vector<FragmentStamp*>::iterator fi = std::find(fragmentStamps_.begin(), fragmentStamps_.end(), static_cast<FragmentStamp*>(NULL));
121 >    if (fi != fragmentStamps_.end()) {
122 >        std::cout << "Error in Molecule " << getName() << ":fragment[" <<  fi - fragmentStamps_.begin()<< "] is missing\n";
123      }
105    have_atoms = 1;
106    atoms = new AtomStamp*[n_atoms];
107    for( i=0; i<n_atoms; i++ ) atoms[i] = NULL;
108  }
109  
110  else if( !strcmp( lhs, "nBonds" ) ){
111    n_bonds = (int)rhs;
124  
125 <    if( have_bonds ){
126 <      sprintf( errMsg,  
127 <               "MoleculeStamp error, n_bonds already declared for"
128 <               " molecule: %s\n",
129 <               name);
130 <      return strdup( errMsg );
125 >    //make sure index is not out of range
126 >    int natoms = getNAtoms();
127 >    for(int i = 0; i < getNBonds(); ++i) {
128 >        BondStamp* bondStamp = getBondStamp(i);
129 >        if (bondStamp->getA() >=  natoms && bondStamp->getB() >= natoms) {
130 >            std::cout << "Error in Molecule " << getName() <<  ": bond between " << bondStamp->getA() << " and " << bondStamp->getB() << " is invalid\n";
131 >        }
132      }
133 <    have_bonds = 1;
134 <    bonds = new BondStamp*[n_bonds];
135 <    for( i=0; i<n_bonds; i++ ) bonds[i] = NULL;
136 <  }
137 <  
138 <  else if( !strcmp( lhs, "nBends" ) ){
139 <    n_bends = (int)rhs;
133 >    for(int i = 0; i < getNBends(); ++i) {
134 >        BendStamp* bendStamp = getBendStamp(i);
135 >        std::vector<int> bendAtoms =  bendStamp->getMembers();
136 >        std::vector<int>::iterator j =std::find_if(bendAtoms.begin(), bendAtoms.end(), std::bind2nd(std::greater<int>(), natoms-1));
137 >        if (j != bendAtoms.end()) {
138 >            std::cout << "Error in Molecule " << getName();
139 >        }
140  
141 <    if( have_bends ){
142 <      sprintf( errMsg,
143 <               "MoleculeStamp error, n_bends already declared for"
144 <               " molecule: %s\n",
145 <               name);
146 <      return strdup( errMsg );
141 >        if (bendAtoms.size() == 2 && !bendStamp->haveGhostVectorSource()) {
142 >            std::cout << "Error in Molecule " << getName() << ": ghostVectorSouce is missing";
143 >        }
144 >    }    
145 >    for(int i = 0; i < getNBends(); ++i) {
146 >        TorsionStamp* torsionStamp = getTorsionStamp(i);
147 >        std::vector<int> torsionAtoms =  torsionStamp ->getMembers();
148 >        std::vector<int>::iterator j =std::find_if(torsionAtoms.begin(), torsionAtoms.end(), std::bind2nd(std::greater<int>(), natoms-1));
149 >        if (j != torsionAtoms.end()) {
150 >            std::cout << "Error in Molecule " << getName();
151 >        }
152      }
153 <    have_bends = 1;
154 <    bends = new BendStamp*[n_bends];
155 <    for( i=0; i<n_bends; i++ ) bends[i] = NULL;
156 <  }
157 <  
158 <  else if( !strcmp( lhs, "nTorsions" ) ){
159 <    n_torsions = (int)rhs;
142 <
143 <    if( have_torsions ){
144 <      sprintf( errMsg,
145 <               "MoleculeStamp error, n_torsions already declared for"
146 <               " molecule: %s\n",
147 <               name );
148 <      return strdup( errMsg );
153 >    for(int i = 0; i < getNCutoffGroups(); ++i) {
154 >        CutoffGroupStamp* cutoffGroupStamp = getCutoffGroupStamp(i);
155 >        std::vector<int> cutoffGroupAtoms =  cutoffGroupStamp ->getMembers();
156 >        std::vector<int>::iterator j =std::find_if(cutoffGroupAtoms.begin(), cutoffGroupAtoms.end(), std::bind2nd(std::greater<int>(), natoms-1));
157 >        if (j != cutoffGroupAtoms.end()) {
158 >            std::cout << "Error in Molecule " << getName();
159 >        }
160      }
161 <    have_torsions = 1;
162 <    torsions = new TorsionStamp*[n_torsions];
163 <    for( i=0; i<n_torsions; i++ ) torsions[i] = NULL;
164 <  }
165 <
166 <  else if( !strcmp( lhs, "nRigidBodies" ) ){
156 <    n_rigidbodies = (int)rhs;
157 <
158 <    if( have_rigidbodies ){
159 <      sprintf( errMsg,
160 <               "MoleculeStamp error, n_rigidbodies already declared for"
161 <               " molecule: %s\n",
162 <               name );
163 <      return strdup( errMsg );
161 >        
162 >    atom2Rigidbody.resize(natoms);
163 >    // negative number means atom is a free atom, does not belong to rigidbody
164 >    //every element in atom2Rigidbody has unique negative number at the very beginning
165 >    for(int i = 0; i < atom2Rigidbody.size(); ++i) {
166 >        atom2Rigidbody[i] = -1 - i;
167      }
165    have_rigidbodies = 1;
166    rigidBodies = new RigidBodyStamp*[n_rigidbodies];
167    for( i=0; i<n_rigidbodies; i++ ) rigidBodies[i] = NULL;
168  }
168  
169 <  else if( !strcmp( lhs, "nCutoffGroups" ) ){
170 <    n_cutoffgroups = (int)rhs;
171 <
172 <    if( have_cutoffgroups ){
173 <      sprintf( errMsg,
174 <               "MoleculeStamp error, n_cutoffgroups already declared for"
176 <               " molecule: %s\n",
177 <               name );
178 <      return strdup( errMsg );
169 >    for (int i = 0; i < getNRigidBodies(); ++i) {
170 >        RigidBodyStamp* rbStamp = getRigidBodyStamp(i);
171 >        std::vector<int> members = rbStamp->getMembers();
172 >        for(std::vector<int>::iterator j = members.begin(); j != members.end(); ++j) {
173 >            atom2Rigidbody[*j] = i;                
174 >        }
175      }
176 <    have_cutoffgroups = 1;
177 <    cutoffGroups = new CutoffGroupStamp*[n_cutoffgroups];
178 <    for( i=0; i<n_cutoffgroups; i++ ) cutoffGroups[i] = NULL;
179 <  }
180 <  
181 <  else{
182 <    if( unhandled == NULL ) unhandled = new LinkedAssign( lhs, rhs );
183 <    else unhandled->add( lhs, rhs );
184 <    have_extras = 1;
185 <  }
186 <  return NULL;
187 < }
188 <
189 < char*  MoleculeStamp::assignInt( char* lhs, int rhs ){
190 <  int i;
191 <  
192 <  if( !strcmp( lhs, "nAtoms" ) ){
193 <    n_atoms = rhs;
194 <    
195 <    if( have_atoms ){
196 <      sprintf( errMsg,
197 <               "MoleculeStamp error, n_atoms already declared for"
198 <               " molecule: %s\n",
199 <               name);
200 <      return strdup( errMsg );
201 <    }
202 <    have_atoms = 1;
203 <    atoms = new AtomStamp*[n_atoms];
204 <    for( i=0; i<n_atoms; i++ ) atoms[i] = NULL;
205 <  }
206 <  
207 <  else if( !strcmp( lhs, "nBonds" ) ){
208 <    n_bonds = rhs;
176 >    //make sure atoms belong to same rigidbody do not bond to each other
177 >    for(int i = 0; i < getNBonds(); ++i) {
178 >        BondStamp* bondStamp = getBondStamp(i);
179 >        if (atom2Rigidbody[bondStamp->getA()] == atom2Rigidbody[bondStamp->getB()])
180 >            std::cout << "Error in Molecule " << getName() << ": "<<"bond between " << bondStamp->getA() << " and " << bondStamp->getB() << "belong to same rigidbody " << atom2Rigidbody[bondStamp->getA()] << "\n";
181 >        }
182 >        
183 >    for(int i = 0; i < getNBends(); ++i) {
184 >        BendStamp* bendStamp = getBendStamp(i);
185 >        std::vector<int> bendAtoms =  bendStamp->getMembers();
186 >        std::vector<int> rigidSet(getNRigidBodies(), 0);
187 >        std::vector<int>::iterator j;
188 >        for( j = bendAtoms.begin(); j != bendAtoms.end(); ++j) {
189 >            int rigidbodyIndex = atom2Rigidbody[*j];
190 >            if (rigidbodyIndex >= 0) {
191 >                ++rigidSet[rigidbodyIndex];
192 >                if (rigidSet[rigidbodyIndex] > 1) {
193 >                    std::cout << "Error in Molecule " << getName() << ": ";
194 >                    //std::cout << "atoms of bend " <<  << "belong to same rigidbody " << rigidbodyIndex << "\n";                    
195 >                }
196 >            }
197 >        }
198 >    }      
199 >    for(int i = 0; i < getNTorsions(); ++i) {
200 >        TorsionStamp* torsionStamp = getTorsionStamp(i);
201 >        std::vector<int> torsionAtoms =  torsionStamp->getMembers();
202 >        std::vector<int> rigidSet(getNRigidBodies(), 0);
203 >        std::vector<int>::iterator j;
204 >        for( j = torsionAtoms.begin(); j != torsionAtoms.end(); ++j) {
205 >            int rigidbodyIndex = atom2Rigidbody[*j];
206 >            if (rigidbodyIndex >= 0) {
207 >                ++rigidSet[rigidbodyIndex];
208 >                if (rigidSet[rigidbodyIndex] > 1) {
209 >                    std::cout << "Error in Molecule " << getName() << ": ";
210 >                    //std::cout << "atoms of torsion " <<  << "belong to same rigidbody " << rigidbodyIndex << "\n";                    
211 >                }
212 >            }
213 >        }
214 >    }
215  
214    if( have_bonds ){
215      sprintf( errMsg,
216               "MoleculeStamp error, n_bonds already declared for"
217               " molecule: %s\n",
218               name);
219      return strdup( errMsg );
220    }
221    have_bonds = 1;
222    bonds = new BondStamp*[n_bonds];
223    for( i=0; i<n_bonds; i++ ) bonds[i] = NULL;
224  }
225  
226  else if( !strcmp( lhs, "nBends" ) ){
227    n_bends = rhs;
216  
217 <    if( have_bends ){
218 <      sprintf( errMsg,
219 <               "MoleculeStamp error, n_bends already declared for"
220 <               " molecule: %s\n",
233 <               name );
234 <      return strdup( errMsg );
235 <    }
236 <    have_bends = 1;
237 <    bends = new BendStamp*[n_bends];
238 <    for( i=0; i<n_bends; i++ ) bends[i] = NULL;
239 <  }
240 <  
241 <  else if( !strcmp( lhs, "nTorsions" ) ){
242 <    n_torsions = rhs;
217 >    //fill in bond information into atom
218 >    fillBondInfo();
219 >    findBends();
220 >    findTorsions();
221  
222 <    if( have_torsions ){
223 <      sprintf( errMsg,
224 <               "MoleculeStamp error, n_torsions already declared for"
225 <               " molecule: %s\n",
248 <               name);
249 <      return strdup( errMsg );
222 >    int nrigidAtoms = 0;
223 >    for (int i = 0; i < getNRigidBodies(); ++i) {
224 >        RigidBodyStamp* rbStamp = getRigidBodyStamp(i);
225 >        nrigidAtoms += rbStamp->getNMembers();
226      }
227 <    have_torsions = 1;
252 <    torsions = new TorsionStamp*[n_torsions];
253 <    for( i=0; i<n_torsions; i++ ) torsions[i] = NULL;
254 <  }
227 >    nintegrable_ = getNAtoms()+ getNRigidBodies() - nrigidAtoms;
228  
229 <  else if( !strcmp( lhs, "nRigidBodies" ) ){
257 <    n_rigidbodies = rhs;
229 > }
230  
231 <    if( have_rigidbodies ){
260 <      sprintf( errMsg,
261 <               "MoleculeStamp error, n_rigidbodies already declared for"
262 <               " molecule: %s\n",
263 <               name);
264 <      return strdup( errMsg );
265 <    }
266 <    have_rigidbodies = 1;
267 <    rigidBodies = new RigidBodyStamp*[n_rigidbodies];
268 <    for( i=0; i<n_rigidbodies; i++ ) rigidBodies[i] = NULL;
269 <  }
270 <  else if( !strcmp( lhs, "nCutoffGroups" ) ){
271 <    n_cutoffgroups = rhs;
272 <
273 <    if( have_cutoffgroups ){
274 <      sprintf( errMsg,
275 <               "MoleculeStamp error, n_cutoffgroups already declared for"
276 <               " molecule: %s\n",
277 <               name);
278 <      return strdup( errMsg );
279 <    }
280 <    have_cutoffgroups = 1;
281 <    cutoffGroups = new CutoffGroupStamp*[n_cutoffgroups];
282 <    for( i=0; i<n_cutoffgroups; i++ ) cutoffGroups[i] = NULL;
283 <  }
284 <  else{
285 <    if( unhandled == NULL ) unhandled = new LinkedAssign( lhs, rhs );
286 <    else unhandled->add( lhs, rhs );
287 <    have_extras = 1;
288 <  }
289 <  return NULL;
231 > void MoleculeStamp::fillBondInfo() {
232   }
233  
234 < char* MoleculeStamp::addAtom( AtomStamp* the_atom, int atomIndex ){
293 <  
294 <  if( have_atoms && atomIndex < n_atoms ) atoms[atomIndex] = the_atom;
295 <  else {
296 <    if( have_atoms ){
297 <      sprintf( errMsg, "MoleculeStamp error, %d out of nAtoms range",
298 <               atomIndex );
299 <      return strdup( errMsg );
300 <    }
301 <    else return strdup("MoleculeStamp error, nAtoms not given before"
302 <                       " first atom declaration." );
303 <  }
234 > void MoleculeStamp::findBends() {
235  
305  return NULL;
236   }
237  
238 < char* MoleculeStamp::addRigidBody( RigidBodyStamp* the_rigidbody,
309 <                                   int rigidBodyIndex ){
310 <  
311 <  if( have_rigidbodies && rigidBodyIndex < n_rigidbodies )
312 <    rigidBodies[rigidBodyIndex] = the_rigidbody;
313 <  else {
314 <    if( have_rigidbodies ){
315 <      sprintf( errMsg, "MoleculeStamp error, %d out of nRigidBodies range",
316 <               rigidBodyIndex );
317 <      return strdup( errMsg );
318 <    }
319 <    else return strdup("MoleculeStamp error, nRigidBodies not given before"
320 <                       " first rigidBody declaration." );
321 <  }
322 <  
323 <  return NULL;
324 < }
238 > void MoleculeStamp::findTorsions() {
239  
326 char* MoleculeStamp::addCutoffGroup( CutoffGroupStamp* the_cutoffgroup,
327                                     int cutoffGroupIndex ){
328  
329  if( have_cutoffgroups && cutoffGroupIndex < n_cutoffgroups )
330    cutoffGroups[cutoffGroupIndex] = the_cutoffgroup;
331  else {
332    if( have_cutoffgroups ){
333      sprintf( errMsg, "MoleculeStamp error, %d out of nCutoffGroups range",
334               cutoffGroupIndex );
335      return strdup( errMsg );
336    }
337    else return strdup("MoleculeStamp error, nCutoffGroups not given before"
338                       " first CutoffGroup declaration." );
339  }
340  
341  return NULL;
240   }
241  
344 char* MoleculeStamp::addBond( BondStamp* the_bond, int bondIndex ){
345  
346  
347  if( have_bonds && bondIndex < n_bonds ) bonds[bondIndex] = the_bond;
348  else{
349    if( have_bonds ){
350      sprintf( errMsg, "MoleculeStamp error, %d out of nBonds range",
351               bondIndex );
352      return strdup( errMsg );
353    }
354    else return strdup("MoleculeStamp error, nBonds not given before"
355                       "first bond declaration." );
356  }
357
358  return NULL;
359 }
360
361 char* MoleculeStamp::addBend( BendStamp* the_bend, int bendIndex ){
362  
363  
364  if( have_bends && bendIndex < n_bends ) bends[bendIndex] = the_bend;
365  else{
366    if( have_bends ){
367      sprintf( errMsg, "MoleculeStamp error, %d out of nBends range",
368               bendIndex );
369      return strdup( errMsg );
370    }
371    else return strdup("MoleculeStamp error, nBends not given before"
372                       "first bend declaration." );
373  }
374
375  return NULL;
376 }
377
378 char* MoleculeStamp::addTorsion( TorsionStamp* the_torsion, int torsionIndex ){
379  
380  
381  if( have_torsions && torsionIndex < n_torsions )
382    torsions[torsionIndex] = the_torsion;
383  else{
384    if( have_torsions ){
385      sprintf( errMsg, "MoleculeStamp error, %d out of nTorsions range",
386               torsionIndex );
387      return strdup( errMsg );
388    }
389    else return strdup("MoleculeStamp error, nTorsions not given before"
390                       "first torsion declaration." );
391  }
392
393  return NULL;
394 }
395
396 char* MoleculeStamp::checkMe( void ){
397  
398  int i;
399  short int no_atom, no_rigidbody, no_cutoffgroup;
400
401  if( !have_name ) return strdup( "MoleculeStamp error. Molecule's name"
402                                  " was not given.\n" );
403  
404  if( !have_atoms ){
405    return strdup( "MoleculeStamp error. Molecule contains no atoms." );
406  }
407  
408  no_rigidbody = 0;
409  for( i=0; i<n_rigidbodies; i++ ){
410    if( rigidBodies[i] == NULL ) no_rigidbody = 1;
411  }
412
413  if( no_rigidbody ){
414    sprintf( errMsg,
415             "MoleculeStamp error. Not all of the RigidBodies were"
416             " declared in molecule \"%s\".\n", name );
417    return strdup( errMsg );
418  }
419
420  no_cutoffgroup = 0;
421  for( i=0; i<n_cutoffgroups; i++ ){
422    if( cutoffGroups[i] == NULL ) no_cutoffgroup = 1;
423  }
424
425  if( no_cutoffgroup ){
426    sprintf( errMsg,
427             "MoleculeStamp error. Not all of the CutoffGroups were"
428             " declared in molecule \"%s\".\n", name );
429    return strdup( errMsg );
430  }
431  
432  no_atom = 0;
433  for( i=0; i<n_atoms; i++ ){
434    if( atoms[i] == NULL ) no_atom = 1;
435  }
436
437  if( no_atom ){
438    sprintf( errMsg,
439             "MoleculeStamp error. Not all of the atoms were"
440             " declared in molecule \"%s\".\n", name );
441    return strdup( errMsg );
442  }
443
444  n_integrable = n_atoms;
445  for (i = 0; i < n_rigidbodies; i++)
446    n_integrable = n_integrable - rigidBodies[i]->getNMembers() + 1; //rigidbody is an integrable object
447  
448  if (n_integrable <= 0 || n_integrable > n_atoms) {
449    sprintf( errMsg,
450             "MoleculeStamp error. n_integrable is either <= 0 or"
451             " greater than n_atoms in molecule \"%s\".\n", name );
452    return strdup( errMsg );
453  }
454
455  return NULL;
456 }  
457
458
242   //Function Name: isBondInSameRigidBody
243   //Return true is both atoms of the bond belong to the same rigid body, otherwise return false
244   bool MoleculeStamp::isBondInSameRigidBody(BondStamp* bond){
# Line 503 | Line 286 | bool MoleculeStamp::isAtomInRigidBody(int atomIndex, i
286    numRb = this->getNRigidBodies();
287    
288    for(int i = 0 ; i < numRb; i++){
289 <    rbStamp = this->getRigidBody(i);
289 >    rbStamp = this->getRigidBodyStamp(i);
290      numAtom = rbStamp->getNMembers();
291      for(int j = 0; j < numAtom; j++)
292 <      if (rbStamp->getMember(j) == atomIndex){
292 >      if (rbStamp->getMemberAt(j) == atomIndex){
293          whichRigidBody = i;
294          consAtomIndex = j;
295          return true;
# Line 520 | Line 303 | vector<pair<int, int> > MoleculeStamp::getJointAtoms(i
303   //return the position of joint atom apears in  rigidbody's definition
304   //for the time being, we will use the most inefficient algorithm, the complexity is O(N2)
305   //actually we could improve the complexity to O(NlgN) by sorting the atom index in rigid body first
306 < vector<pair<int, int> > MoleculeStamp::getJointAtoms(int rb1, int rb2){
306 > std::vector<std::pair<int, int> > MoleculeStamp::getJointAtoms(int rb1, int rb2){
307    RigidBodyStamp* rbStamp1;
308    RigidBodyStamp* rbStamp2;
309    int natomInRb1;
310    int natomInRb2;
311    int atomIndex1;
312    int atomIndex2;
313 <  vector<pair<int, int> > jointAtomIndexPair;
313 >  std::vector<std::pair<int, int> > jointAtomIndexPair;
314    
315 <  rbStamp1 = this->getRigidBody(rb1);
315 >  rbStamp1 = this->getRigidBodyStamp(rb1);
316    natomInRb1 =rbStamp1->getNMembers();
317  
318 <  rbStamp2 = this->getRigidBody(rb2);
318 >  rbStamp2 = this->getRigidBodyStamp(rb2);
319    natomInRb2 =rbStamp2->getNMembers();
320  
321    for(int i = 0; i < natomInRb1; i++){
322 <    atomIndex1 = rbStamp1->getMember(i);
322 >    atomIndex1 = rbStamp1->getMemberAt(i);
323        
324      for(int j= 0; j < natomInRb1; j++){
325 <      atomIndex2 = rbStamp2->getMember(j);
325 >      atomIndex2 = rbStamp2->getMemberAt(j);
326  
327        if(atomIndex1 == atomIndex2){
328 <        jointAtomIndexPair.push_back(make_pair(i, j));
328 >        jointAtomIndexPair.push_back(std::make_pair(i, j));
329          break;
330        }
331        
332 <    }//end for(j =0)
332 >    }
333  
334 <  }//end for (i = 0)
334 >  }
335  
336    return jointAtomIndexPair;
337   }
338 +
339 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines