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

Comparing trunk/OOPSE-2.0/src/visitors/ZconsVisitor.cpp (file contents):
Revision 1625 by tim, Thu Oct 21 16:22:01 2004 UTC vs.
Revision 1930 by gezelter, Wed Jan 12 22:41:40 2005 UTC

# Line 1 | Line 1
1 < #include "visitors/ZconsVisitor.hpp"
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 <cmath>
43 <
43 > #include "visitors/ZconsVisitor.hpp"
44 > #include "primitives/Molecule.hpp"
45 > #include "utils/StringUtils.hpp"
46   namespace oopse {
47  
48 < ZConsVisitor::ZConsVisitor(SimInfo* info) : BaseVisitor(), zconsReader(NULL){
7 <  GenericData* data;
8 <  DoubleGenericData* tolerance;  
9 <  ZConsParaData* zConsParaData;
10 <  StringGenericData* filename;
11 <  DoubleGenericData* sampleTime;
12 <  vector<ZConsParaItem>* parameters;
48 > ZConsVisitor::ZConsVisitor(SimInfo* info) : BaseVisitor(), info_(info), zconsReader_(NULL){
49  
50 <  this->info = info;
51 <  visitorName = "ZConsVisitor";
52 <  
17 <  //retrieve tolerance for z-constraint molecuels
18 <  data = info->getProperty(ZCONSTOL_ID);
50 >    visitorName = "ZConsVisitor";
51 >    currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
52 >    Globals* simParam = info_->getSimParams();
53  
54 <  if (!data){
55 <    cerr << "Can not get zconstraint tolerance from SimInfo" << endl;
22 <    haveZcons = false;  
23 <    return;
24 <  }
25 <  else{
26 <    tolerance = dynamic_cast<DoubleGenericData*>(data);
27 <
28 <    if (!tolerance){
29 <      cerr << "Can not get zconstraint tolerance  from SimInfo" << endl;
30 <      haveZcons = false;    
31 <      return;
54 >    if (simParam->haveZconstraintTime()){
55 >        zconsTime_ = simParam->getZconsTime();
56      }
57      else{
58 <      zconsTol = tolerance->getData();
58 >        sprintf(painCave.errMsg,
59 >        "ZConstraint error: If you use a ZConstraint,\n"
60 >        "\tyou must set zconsTime.\n");
61 >        painCave.isFatal = 1;
62 >        simError();
63      }
36  }
64  
65 <  //retrieve sample time of z-contraint
66 <  data = info->getProperty(ZCONSTIME_ID);
40 <
41 <  if (!data){
42 <    cerr << "Can not get zcons time  from SimInfo" << endl;
43 <    haveZcons = false;    
44 <    return;
45 <  }
46 <  else{
47 <    sampleTime = dynamic_cast<DoubleGenericData*>(data);
48 <
49 <    if (!sampleTime){
50 <      cerr << "Can not get zcons time  from SimInfo" << endl;
51 <      haveZcons = false;    
52 <      return;
65 >    if (simParam->haveZconsTol()){
66 >        zconsTol_ = simParam->getZconsTol();
67      }
68      else{
69 <     zconsTime = sampleTime->getData();
69 >        zconsTol_ = 0.01;
70 >        sprintf(painCave.errMsg,
71 >            "ZConstraint Warning: Tolerance for z-constraint method is not specified.\n"
72 >            "\tOOPSE will use a default value of %f.\n"
73 >            "\tTo set the tolerance, use the zconsTol variable.\n",
74 >            zconsTol_);
75 >        painCave.isFatal = 0;
76 >        simError();      
77 >    }    
78 >        
79 >    int nZconstraints = simParam->getNzConstraints();
80 >    ZconStamp** stamp = simParam->getZconStamp();
81 >    for (int i = 0; i < nZconstraints; i++){
82 >        int zmolIndex = stamp[i]->getMolIndex();
83 >        zmolStates_.insert(std::make_pair(zmolIndex, zsMoving));
84      }
57  }
85  
59  //retrieve index of z-constraint molecules
60  data = info->getProperty(ZCONSPARADATA_ID);
61  if (!data){
62    cerr << "Can not get index of z-constraint molecules from SimInfo" << endl;
63    haveZcons = false;
64    return;
65  }
66  else{
67    zConsParaData = dynamic_cast<ZConsParaData*>(data);
86  
87 <    if (!zConsParaData){
88 <      cerr << "Can not get index of z-constraint molecules from SimInfo" << endl;
89 <      haveZcons = false;
90 <      return;
87 >    //fill zatomToZmol_ array
88 >    /** @todo only works for single version now*/
89 >    std::map<int, ZConsState>::iterator j;
90 >    for (j = zmolStates_.begin(); j != zmolStates_.end(); ++j) {
91 >        Molecule* mol = info_->getMoleculeByGlobalIndex(j->first);
92 >        assert(mol != NULL);
93 >        Molecule::AtomIterator ai;
94 >        Atom* at;
95 >        for (at = mol->beginAtom(ai); at != NULL; at = mol->nextAtom(ai)) {
96 >            zatomToZmol_.insert(std::make_pair(at->getGlobalIndex(), mol->getGlobalIndex()));
97 >        }
98      }
74    else{
75      vector<ZConsParaItem>::iterator i;
76      Molecule* mol;
99  
100 <      parameters = zConsParaData->getData();
101 <      for(i = parameters->begin(); i != parameters->end(); ++i){
100 >    zconsFilename_ = getPrefix(info_->getFinalConfigFileName()) + ".fz";
101 >    
102 >    zconsReader_ = new ZConsReader(info);
103  
104 <        mol = findZconsMol(i->zconsIndex);
105 <        if(mol != NULL)
83 <          zconsMol.push_back(mol);
84 <
85 <      }
86 <
87 <      if(zconsMol.size() < 1){
88 <        cerr << "number of zconstraint molecules is less than one" << endl;
89 <        haveZcons = false;
90 <        return;
91 <      }
92 <      
93 <    }//end if (!zConsParaData)
94 <
95 <  }//end  if (!data  
104 >    if (zconsReader_->hasNextFrame())
105 >        zconsReader_->readNextFrame();
106    
97  //retrieve output filename of z force
98  data = info->getProperty(ZCONSFILENAME_ID);
99  if (!data){
100    cerr << "Can not get filename of z-constraint output from SimInfo" << endl;
101    haveZcons = false;
102    return;
103  }
104  else{
105    filename = dynamic_cast<StringGenericData*>(data);
106
107    if (!filename){
108      cerr << "Can not get filename of z-constraint output from SimInfo" << endl;
109      haveZcons = false;
110      return;
111    }
112    else{
113      zconsFilename = filename->getData();
114    }
115  }
116
117  zconsReader = new ZConsReader(info);
118
119  if (zconsReader->hasNextFrame())
120  zconsReader->readNextFrame();
121  
122  haveZcons = true;
107   }
108  
109   ZConsVisitor::~ZConsVisitor(){
110 <  if(!zconsReader)
111 <    delete zconsReader;
110 >  if(!zconsReader_)
111 >    delete zconsReader_;
112    
113   }
114  
115   void ZConsVisitor::visit(Atom* atom){
116 <  string prefix;
117 <  if(isZconstraint(atom->getIndex(), prefix))
118 <    internalVisit(atom, prefix);
116 >    std::string prefix;
117 >    if(isZconstraint(atom->getGlobalIndex(), prefix))
118 >        internalVisit(atom, prefix);
119   }
120  
121   void ZConsVisitor::visit(DirectionalAtom* datom){
122 <  string prefix;
123 <  
124 <  if(isZconstraint(datom->getIndex(), prefix))
125 <    internalVisit(datom, prefix);
122 >    std::string prefix;
123 >
124 >    if(isZconstraint(datom->getGlobalIndex(), prefix))
125 >        internalVisit(datom, prefix);
126   }
127  
128   void ZConsVisitor::visit(RigidBody* rb){
129 <  string prefix;
130 <  vector<Atom*> atoms;
131 <  
132 <  atoms = rb->getAtoms();
133 <    
134 <  if(isZconstraint(atoms[0]->getIndex(), prefix))
135 <    internalVisit(rb, prefix);
129 >    std::string prefix;
130 >    std::vector<Atom*> atoms;
131 >
132 >    atoms = rb->getAtoms();
133 >
134 >    if(isZconstraint(atoms[0]->getGlobalIndex(), prefix))
135 >        internalVisit(rb, prefix);
136   }
137  
138 < Molecule* ZConsVisitor::findZconsMol(int index){
139 <  Molecule* mols;
140 <  mols = info->molecules;
141 <  for(int i =0; i < info->n_mol; i++)
142 <    if(index == mols[i].getMyIndex())
143 <      return &mols[i];
138 > void ZConsVisitor::update(){
139 >    Molecule* mol;
140 >    Vector3d com;
141 >    ZConsState state;
142 >    std::map<int, ZConsState>::iterator i;
143 >    for ( i = zmolStates_.begin(); i != zmolStates_.end(); ++i) {
144 >        i->second = zsMoving;
145 >    }
146 >    
147 >    readZconsFile(currSnapshot_->getTime());
148  
149 <  return NULL;
149 >    const std::vector<ZconsData>& fixedZmolData = zconsReader_->getFixedZMolData();
150 >    std::vector<ZconsData>::const_iterator j;
151 >    for (j = fixedZmolData.begin(); j != fixedZmolData.end(); ++j) {
152 >        std::map<int, ZConsState>::iterator k = zmolStates_.find(j->zmolIndex);
153 >        assert(k != zmolStates_.end());
154 >        k->second = zsFixed;
155 >    }
156 >    
157   }
158  
159 < void ZConsVisitor::update(){
160 <  Molecule* mol;
161 <  double com[3];
162 <  ZConsState state;
163 <  Atom** atoms;
164 <  int natoms;
159 > void ZConsVisitor::readZconsFile(double time) {
160 >    double tempTime;
161 >    while(zconsReader_->hasNextFrame()){
162 >        tempTime = zconsReader_->getCurTime();
163 >        if(tempTime >= time) {
164 >            return;
165 >        }
166 >        
167 >        zconsReader_->readNextFrame();
168 >    }
169 > }
170  
171 <  getZconsPos(info->currentTime);
172 <  
173 <  for(size_t i = 0; i < zconsMol.size(); i++){
174 <    zconsMol[i]->getCOM(com);
171 > void ZConsVisitor::internalVisit(StuntDouble* sd, const std::string& prefix){
172 >    GenericData* data;
173 >    AtomData* atomData;
174 >    AtomInfo* atomInfo;
175 >    std::vector<AtomInfo*>::iterator iter;
176  
177 <    if(fabs(com[2] - zconsPos[i]) < zconsTol)
178 <      state = zsFixed;
177 >    //if there is not atom data, just skip it
178 >    data = sd->getPropertyByName("ATOMDATA");
179 >    if(data != NULL){
180 >        atomData = dynamic_cast<AtomData*>(data);  
181 >        if(atomData == NULL)
182 >            return;
183 >        }
184      else
185 <      state = zsMoving;
185 >        return;
186  
187 <    //update the state of zconstraint atom
188 <    natoms = zconsMol[i]->getNAtoms();
183 <    atoms = zconsMol[i]->getMyAtoms();    
184 <    for(int j =0; j < natoms; j++)
185 <      zconsState[atoms[j]->getIndex()] = state;
186 <  }
187 <  
187 >    for(atomInfo  = atomData->beginAtomInfo(iter); atomInfo; atomInfo = atomData->nextAtomInfo(iter))
188 >        (atomInfo->AtomType).insert(0, prefix);
189   }
190  
190 void ZConsVisitor::getZconsPos(double time){
191  
192 <  vector<double> tempPos;
193 <  vector<double> prevPos;
194 <  double tempTime;
192 > bool ZConsVisitor::isZconstraint(int atomIndex, std::string& prefix){
193 >    std::string prefixString[] = {"ZF", "ZM"};
194 >    std::map<int, int>::iterator i = zatomToZmol_.find(atomIndex);
195 >    if (i ==  zatomToZmol_.end() ){
196 >        prefix = "";
197 >        return false;
198 >    } else {
199  
200 <  while(true){
201 <    tempTime = zconsReader->getCurTime();
202 <
203 <    zconsPos = zconsReader->getZConsPos();
204 <    
201 <    if(tempTime >= time)
202 <      return;
203 <    
204 <    zconsReader->readNextFrame();
205 <
206 <  }
207 <
200 >        std::map<int, ZConsState>::iterator j = zmolStates_.find(i->second);
201 >        assert(j !=zmolStates_.end());
202 >        prefix = prefixString[j->second];
203 >        return true;
204 >    }
205   }
206  
207 < void ZConsVisitor::internalVisit(StuntDouble* sd, const string& prefix){
208 <  GenericData* data;
209 <  AtomData* atomData;
213 <  AtomInfo* atomInfo;
214 <  vector<AtomInfo*>::iterator iter;
207 > const std::string ZConsVisitor::toString(){
208 >    char buffer[65535];
209 >    std::string result;
210  
211 <  
212 <  //if there is not atom data, just skip it
218 <  data = sd->getProperty("ATOMDATA");
219 <  if(data != NULL){
220 <    atomData = dynamic_cast<AtomData*>(data);  
221 <    if(atomData == NULL)
222 <      return;
223 <  }
224 <  else
225 <    return;
211 >    sprintf(buffer ,"------------------------------------------------------------------\n");
212 >    result += buffer;
213  
214 <  for(atomInfo  = atomData->beginAtomInfo(iter); atomInfo; atomInfo = atomData->nextAtomInfo(iter))
215 <    (atomInfo->AtomType).insert(0, prefix);
229 < }
214 >    sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
215 >    result += buffer;
216  
217 +    sprintf(buffer , "number of zconstraint molecule: %d\n", zmolStates_.size());
218 +    result += buffer;
219  
220 < bool ZConsVisitor::isZconstraint(int index, string& prefix){
221 <  map<int, ZConsState>::iterator i;
234 <  string prefixString[] = {"ZF", "ZM"};
235 <  
236 <  i = zconsState.find(index);
237 <  if(i !=zconsState.end()){
238 <    prefix = prefixString[(*i).second];
239 <    return true;
240 <  }
241 <  else{
242 <    prefix = "";
243 <    return false;
244 <  }
245 < }
220 >    sprintf(buffer , "zconstraint tolerance = %lf\n", zconsTol_);
221 >    result += buffer;
222  
223 < const string ZConsVisitor::toString(){
224 <  char buffer[65535];
249 <  string result;
250 <  
251 <  sprintf(buffer ,"------------------------------------------------------------------\n");
252 <  result += buffer;
223 >    sprintf(buffer , "zconstraint sample time = %lf\n", zconsTime_);
224 >    result += buffer;
225  
226 <  sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
227 <  result += buffer;
256 <  
257 <  sprintf(buffer , "number of zconstraint molecule: %d\n", zconsMol.size());
258 <  result += buffer;
226 >    sprintf(buffer , "zconstraint output filename = %s\n", zconsFilename_.c_str());
227 >    result += buffer;
228  
229 <  sprintf(buffer , "zconstraint tolerance = %lf\n", zconsTol);
230 <  result += buffer;
231 <
232 <  sprintf(buffer , "zconstraint sample time = %lf\n", zconsTime);
233 <  result += buffer;
234 <
266 <  sprintf(buffer , "zconstraint output filename = %s\n", zconsFilename.c_str());
267 <  result += buffer;
229 >    std::map<int, ZConsState>::iterator i;
230 >    int j = 0;
231 >    for ( i = zmolStates_.begin(); i != zmolStates_.end(); ++i) {
232 >        sprintf(buffer ,"zconstraint molecule[%d] = %d\n", j++, i->first);
233 >        result += buffer;
234 >    }
235    
236 <  for(size_t i = 0; i < zconsMol.size(); ++i){
270 <    sprintf(buffer ,"zconstraint molecule[%d] = %d\n", i, zconsMol[i]->getMyIndex());
236 >    sprintf(buffer ,"------------------------------------------------------------------\n");
237      result += buffer;
238  
239 <  }
274 <
275 <  
276 <  sprintf(buffer ,"------------------------------------------------------------------\n");
277 <  result += buffer;
278 <
279 <  return result;
239 >    return result;
240   }
241  
242  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines