ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/visitors/ZconsVisitor.cpp
(Generate patch)

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines