ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/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 2071 by gezelter, Sat Mar 7 21:41:51 2015 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. Redistributions of source code must retain the above copyright
10 > *    notice, this list of conditions and the following disclaimer.
11 > *
12 > * 2. Redistributions in binary form must reproduce the above copyright
13 > *    notice, this list of conditions and the following disclaimer in the
14 > *    documentation and/or other materials provided with the
15 > *    distribution.
16 > *
17 > * This software is provided "AS IS," without a warranty of any
18 > * kind. All express or implied conditions, representations and
19 > * warranties, including any implied warranty of merchantability,
20 > * fitness for a particular purpose or non-infringement, are hereby
21 > * excluded.  The University of Notre Dame and its licensors shall not
22 > * be liable for any damages suffered by licensee as a result of
23 > * using, modifying or distributing the software or its
24 > * derivatives. In no event will the University of Notre Dame or its
25 > * licensors be liable for any lost revenue, profit or data, or for
26 > * direct, indirect, special, consequential, incidental or punitive
27 > * damages, however caused and regardless of the theory of liability,
28 > * arising out of the use of or inability to use software, even if the
29 > * University of Notre Dame has been advised of the possibility of
30 > * such damages.
31 > *
32 > * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
33 > * research, please cite the appropriate papers when you publish your
34 > * work.  Good starting points are:
35 > *                                                                      
36 > * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).            
37 > * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
38 > * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
39 > * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 > * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41 > */
42 >
43   #include <cmath>
44 + #include "visitors/ZconsVisitor.hpp"
45 + #include "primitives/Molecule.hpp"
46 + #include "utils/StringUtils.hpp"
47 + #include "types/ZconsStamp.hpp"
48 + namespace OpenMD {
49  
50 < namespace oopse {
50 >  ZConsVisitor::ZConsVisitor(SimInfo* info) : BaseVisitor(),
51 >                                              zconsReader_(NULL), info_(info){
52 >    
53 >    visitorName = "ZConsVisitor";
54 >    currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
55 >    Globals* simParam = info_->getSimParams();
56  
57 < ZConsVisitor::ZConsVisitor(SimInfo* info) : BaseVisitor(), zconsReader(NULL){
58 <  GenericData* data;
8 <  DoubleGenericData* tolerance;  
9 <  ZConsParaData* zConsParaData;
10 <  StringGenericData* filename;
11 <  DoubleGenericData* sampleTime;
12 <  vector<ZConsParaItem>* parameters;
13 <
14 <  this->info = info;
15 <  visitorName = "ZConsVisitor";
16 <  
17 <  //retrieve tolerance for z-constraint molecuels
18 <  data = info->getProperty(ZCONSTOL_ID);
19 <
20 <  if (!data){
21 <    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;
57 >    if (simParam->haveZconsTime()){
58 >      zconsTime_ = simParam->getZconsTime();
59      }
60      else{
61 <      zconsTol = tolerance->getData();
61 >      sprintf(painCave.errMsg,
62 >              "ZConstraint error: If you use a ZConstraint,\n"
63 >              "\tyou must set zconsTime.\n");
64 >      painCave.isFatal = 1;
65 >      simError();
66      }
36  }
67  
68 <  //retrieve sample time of z-contraint
69 <  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;
68 >    if (simParam->haveZconsTol()){
69 >      zconsTol_ = simParam->getZconsTol();
70      }
71      else{
72 <     zconsTime = sampleTime->getData();
72 >      zconsTol_ = 0.01;
73 >      sprintf(painCave.errMsg,
74 >              "ZConstraint Warning: Tolerance for z-constraint method is not specified.\n"
75 >              "\tOpenMD will use a default value of %f.\n"
76 >              "\tTo set the tolerance, use the zconsTol variable.\n",
77 >              zconsTol_);
78 >      painCave.isFatal = 0;
79 >      simError();      
80 >    }    
81 >        
82 >    int nZconstraints = simParam->getNZconsStamps();
83 >    std::vector<ZConsStamp*> stamp = simParam->getZconsStamps();
84 >    for (int i = 0; i < nZconstraints; i++){
85 >      int zmolIndex = stamp[i]->getMolIndex();
86 >      zmolStates_.insert(std::make_pair(zmolIndex, zsMoving));
87      }
57  }
88  
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);
89  
90 <    if (!zConsParaData){
91 <      cerr << "Can not get index of z-constraint molecules from SimInfo" << endl;
92 <      haveZcons = false;
93 <      return;
90 >    //fill zatomToZmol_ array
91 >    /** @todo only works for single version now*/
92 >    std::map<int, ZConsState>::iterator j;
93 >    for (j = zmolStates_.begin(); j != zmolStates_.end(); ++j) {
94 >      Molecule* mol = info_->getMoleculeByGlobalIndex(j->first);
95 >      assert(mol != NULL);
96 >      Molecule::AtomIterator ai;
97 >      Atom* at;
98 >      for (at = mol->beginAtom(ai); at != NULL; at = mol->nextAtom(ai)) {
99 >        zatomToZmol_.insert(std::make_pair(at->getGlobalIndex(), mol->getGlobalIndex()));
100 >      }
101      }
74    else{
75      vector<ZConsParaItem>::iterator i;
76      Molecule* mol;
102  
103 <      parameters = zConsParaData->getData();
104 <      for(i = parameters->begin(); i != parameters->end(); ++i){
103 >    zconsFilename_ = getPrefix(info_->getFinalConfigFileName()) + ".fz";
104 >    
105 >    zconsReader_ = new ZConsReader(info);
106  
107 <        mol = findZconsMol(i->zconsIndex);
108 <        if(mol != NULL)
109 <          zconsMol.push_back(mol);
107 >    if (zconsReader_->hasNextFrame())
108 >      zconsReader_->readNextFrame();
109 >  
110 >  }
111  
112 <      }
113 <
114 <      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  
112 >  ZConsVisitor::~ZConsVisitor(){
113 >    if(!zconsReader_)
114 >      delete zconsReader_;
115    
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;
116    }
104  else{
105    filename = dynamic_cast<StringGenericData*>(data);
117  
118 <    if (!filename){
119 <      cerr << "Can not get filename of z-constraint output from SimInfo" << endl;
120 <      haveZcons = false;
121 <      return;
111 <    }
112 <    else{
113 <      zconsFilename = filename->getData();
114 <    }
118 >  void ZConsVisitor::visit(Atom* atom){
119 >    std::string prefix;
120 >    if(isZconstraint(atom->getGlobalIndex(), prefix))
121 >      internalVisit(atom, prefix);
122    }
123  
124 <  zconsReader = new ZConsReader(info);
124 >  void ZConsVisitor::visit(DirectionalAtom* datom){
125 >    std::string prefix;
126  
127 <  if (zconsReader->hasNextFrame())
128 <  zconsReader->readNextFrame();
129 <  
122 <  haveZcons = true;
123 < }
127 >    if(isZconstraint(datom->getGlobalIndex(), prefix))
128 >      internalVisit(datom, prefix);
129 >  }
130  
131 < ZConsVisitor::~ZConsVisitor(){
132 <  if(!zconsReader)
133 <    delete zconsReader;
128 <  
129 < }
131 >  void ZConsVisitor::visit(RigidBody* rb){
132 >    std::string prefix;
133 >    std::vector<Atom*> atoms;
134  
135 < void ZConsVisitor::visit(Atom* atom){
132 <  string prefix;
133 <  if(isZconstraint(atom->getIndex(), prefix))
134 <    internalVisit(atom, prefix);
135 < }
135 >    atoms = rb->getAtoms();
136  
137 < void ZConsVisitor::visit(DirectionalAtom* datom){
138 <  string prefix;
139 <  
140 <  if(isZconstraint(datom->getIndex(), prefix))
141 <    internalVisit(datom, prefix);
142 < }
137 >    if(isZconstraint(atoms[0]->getGlobalIndex(), prefix))
138 >      internalVisit(rb, prefix);
139 >  }
140  
141 < void ZConsVisitor::visit(RigidBody* rb){
142 <  string prefix;
143 <  vector<Atom*> atoms;
144 <  
145 <  atoms = rb->getAtoms();
141 >  void ZConsVisitor::update(){
142 >    Vector3d com;
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", (int) 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  
244 <  return result;
280 < }
281 <
282 <
283 < }//namespace oopse
244 > }//namespace OpenMD

Comparing trunk/src/visitors/ZconsVisitor.cpp (property svn:keywords):
Revision 132 by tim, Thu Oct 21 16:22:01 2004 UTC vs.
Revision 2071 by gezelter, Sat Mar 7 21:41:51 2015 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines