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 2 by gezelter, Fri Sep 24 04:16:43 2004 UTC vs.
branches/development/src/visitors/ZconsVisitor.cpp (file contents), Revision 1665 by gezelter, Tue Nov 22 20:38:56 2011 UTC

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

Comparing:
trunk/src/visitors/ZconsVisitor.cpp (property svn:keywords), Revision 2 by gezelter, Fri Sep 24 04:16:43 2004 UTC vs.
branches/development/src/visitors/ZconsVisitor.cpp (property svn:keywords), Revision 1665 by gezelter, Tue Nov 22 20:38:56 2011 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines