OOPSE/libmdtools/AtomVisitor.cpp

#include <cstring>
#include "AtomVisitor.hpp"
#include "DirectionalAtom.hpp"
#include "MatVec3.h"
#include "RigidBody.hpp"

void BaseAtomVisitor::visit(RigidBody* rb){
  //vector<Atom*> myAtoms;
  //vector<Atom*>::iterator atomIter;

  //myAtoms = rb->getAtoms();
  
  //for(atomIter = myAtoms.begin(); atomIter != myAtoms.end(); ++atomIter)
  //  (*atomIter)->accept(this);
}

void BaseAtomVisitor::setVisited(Atom* atom){
  GenericData* data;
  data = atom->getProperty("VISITED");

  //if visited property is not existed, add it as new property
  if(data == NULL){
    data = new GenericData();
    data->setID("VISITED");
    atom->addProperty(data);  
  }
}

bool BaseAtomVisitor::isVisited(Atom* atom){
  GenericData* data;
  data = atom->getProperty("VISITED");
  return data == NULL ?  false : true;
}

void SSDAtomVisitor::visit(DirectionalAtom* datom){

  vector<AtomInfo*> atoms;

  //we need to convert SSD into 4 differnet atoms
  //one oxygen atom, two hydrogen atoms and one pseudo atom which is the center of the mass
  //of the water with a dipole moment
  double h1[3] = {0.0, -0.75695, 0.5206};
  double h2[3] = {0.0, 0.75695, 0.5206};
  double ox[3] = {0.0, 0.0, -0.0654};
  double u[3] = {0, 0, 1};
  double rotMatrix[3][3];
  double rotTrans[3][3];
  AtomInfo* atomInfo;
  double pos[3];
  double vel[3];
  double newVec[3];
  double q[4];
  AtomData* atomData;
  GenericData* data;
  bool haveAtomData;
  
  //if atom is not SSD atom, just skip it
  if(strcmp(datom->getType(), "SSD"))
    return;

  data = datom->getProperty("ATOMDATA");
  if(data != NULL){

    atomData = dynamic_cast<AtomData*>(data);  
    if(atomData == NULL){
      cerr << "can not get Atom Data from " << datom->getType() << endl;
      atomData = new AtomData; 
      haveAtomData = false;      
    }
    else
      haveAtomData = true;
  }
  else{
    atomData = new AtomData;
    haveAtomData = false;
  }
   
  
  datom->getPos(pos);
  datom->getQ(q);
  datom->getA(rotMatrix);

  // We need A^T to convert from body-fixed to space-fixed:
  transposeMat3(rotMatrix, rotTrans);
  
  //center of mass of the water molecule
  matVecMul3(rotTrans, u, newVec);
  atomInfo = new AtomInfo;
  atomInfo->AtomType = "X";
  atomInfo->pos[0] = pos[0];
  atomInfo->pos[1] = pos[1];
  atomInfo->pos[2] = pos[2];
  atomInfo->dipole[0] = newVec[0];
  atomInfo->dipole[1] = newVec[1];
  atomInfo->dipole[2] = newVec[2];

  atomData->addAtomInfo(atomInfo);

  //oxygen
  matVecMul3(rotTrans, ox, newVec);
  atomInfo = new AtomInfo;
  atomInfo->AtomType = "O";
  atomInfo->pos[0] = pos[0] + newVec[0];
  atomInfo->pos[1] = pos[1] + newVec[1];
  atomInfo->pos[2] = pos[2] + newVec[2];
  atomInfo->dipole[0] = 0.0;
  atomInfo->dipole[1] = 0.0;
  atomInfo->dipole[2] = 0.0;
  atomData->addAtomInfo(atomInfo);


  //hydrogen1
    matVecMul3(rotTrans, h1, newVec);
  atomInfo = new AtomInfo;
  atomInfo->AtomType = "H";
  atomInfo->pos[0] = pos[0] + newVec[0];
  atomInfo->pos[1] = pos[1] + newVec[1];
  atomInfo->pos[2] = pos[2] + newVec[2];
  atomInfo->dipole[0] = 0.0;
  atomInfo->dipole[1] = 0.0;
  atomInfo->dipole[2] = 0.0;
  atomData->addAtomInfo(atomInfo);

  //hydrogen2
  matVecMul3(rotTrans, h2, newVec);
  atomInfo = new AtomInfo;
  atomInfo->AtomType = "H";
  atomInfo->pos[0] = pos[0] + newVec[0];
  atomInfo->pos[1] = pos[1] + newVec[1];
  atomInfo->pos[2] = pos[2] + newVec[2];
  atomInfo->dipole[0] = 0.0;
  atomInfo->dipole[1] = 0.0;
  atomInfo->dipole[2] = 0.0;
  atomData->addAtomInfo(atomInfo);

  //add atom data into atom's property

  if(!haveAtomData){
    atomData->setID("ATOMDATA");
    datom->addProperty(atomData);
  }

  setVisited(datom);

}

const string SSDAtomVisitor::toString(){
  char buffer[65535];
  string result;
  
  sprintf(buffer ,"------------------------------------------------------------------\n");
  result += buffer;

  sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
  result += buffer;

  sprintf(buffer , "Visitor Description: Convert SSD into 4 different atoms\n");
  result += buffer;

  sprintf(buffer ,"------------------------------------------------------------------\n");
  result += buffer;

  return result;
}

//----------------------------------------------------------------------------//

void DefaultAtomVisitor::visit(Atom* atom){
  AtomData* atomData;
  AtomInfo* atomInfo;
  double pos[3];

  if(isVisited(atom))
    return;

 atomInfo =new AtomInfo;

  atomData = new AtomData; 
  atomData->setID("ATOMDATA");
 
  atom->getPos(pos);
  atomInfo->AtomType = atom->getType();
  atomInfo->pos[0] = pos[0];
  atomInfo->pos[1] = pos[1];
  atomInfo->pos[2] = pos[2];
  atomInfo->dipole[0] = 0.0;
  atomInfo->dipole[1] = 0.0;
  atomInfo->dipole[2] = 0.0;


  atomData->addAtomInfo(atomInfo);
  
  atom->addProperty(atomData);

  setVisited(atom);
}
void DefaultAtomVisitor::visit(DirectionalAtom* datom){
  AtomData* atomData;
  AtomInfo* atomInfo;
  double pos[3];
  double u[3];

  if(isVisited(datom))
    return;
  
  datom->getPos(pos);
  datom->getU(u);

  atomData = new AtomData; 
  atomData->setID("ATOMDATA");
  atomInfo =new AtomInfo;
  
  atomInfo->AtomType = datom->getType();
  atomInfo->pos[0] = pos[0];
  atomInfo->pos[1] = pos[1];
  atomInfo->pos[2] = pos[2];
  atomInfo->dipole[0] = u[0];
  atomInfo->dipole[1] = u[1];
  atomInfo->dipole[2] = u[2];  

  atomData->addAtomInfo(atomInfo);

  datom->addProperty(atomData);

  setVisited(datom);
}


const string DefaultAtomVisitor::toString(){
  char buffer[65535];
  string result;
  
  sprintf(buffer ,"------------------------------------------------------------------\n");
  result += buffer;

  sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
  result += buffer;

  sprintf(buffer , "Visitor Description: copy atom infomation into atom data\n");
  result += buffer;

  sprintf(buffer ,"------------------------------------------------------------------\n");
  result += buffer;

  return result;
}    
Revision:	1132
Committed:	Sat Apr 24 04:31:36 2004 UTC (20 years, 2 months ago) by tim
File size:	5942 byte(s)
Log Message:	add reaction field correction to charge-charge interaction
#	User	Rev	Content
1	tim	1118	#include <cstring>
2			#include "AtomVisitor.hpp"
3			#include "DirectionalAtom.hpp"
4			#include "MatVec3.h"
5	tim	1119	#include "RigidBody.hpp"
6	tim	1118
7	tim	1119	void BaseAtomVisitor::visit(RigidBody* rb){
8			//vector<Atom*> myAtoms;
9			//vector<Atom*>::iterator atomIter;
10
11			//myAtoms = rb->getAtoms();
12
13			//for(atomIter = myAtoms.begin(); atomIter != myAtoms.end(); ++atomIter)
14			// (*atomIter)->accept(this);
15			}
16
17	tim	1118	void BaseAtomVisitor::setVisited(Atom* atom){
18			GenericData* data;
19			data = atom->getProperty("VISITED");
20
21			//if visited property is not existed, add it as new property
22			if(data == NULL){
23			data = new GenericData();
24			data->setID("VISITED");
25			atom->addProperty(data);
26			}
27			}
28
29			bool BaseAtomVisitor::isVisited(Atom* atom){
30			GenericData* data;
31			data = atom->getProperty("VISITED");
32			return data == NULL ? false : true;
33			}
34
35			void SSDAtomVisitor::visit(DirectionalAtom* datom){
36
37			vector<AtomInfo*> atoms;
38
39			//we need to convert SSD into 4 differnet atoms
40			//one oxygen atom, two hydrogen atoms and one pseudo atom which is the center of the mass
41			//of the water with a dipole moment
42			double h1[3] = {0.0, -0.75695, 0.5206};
43			double h2[3] = {0.0, 0.75695, 0.5206};
44			double ox[3] = {0.0, 0.0, -0.0654};
45			double u[3] = {0, 0, 1};
46			double rotMatrix[3][3];
47	tim	1129	double rotTrans[3][3];
48	tim	1118	AtomInfo* atomInfo;
49			double pos[3];
50			double vel[3];
51			double newVec[3];
52			double q[4];
53			AtomData* atomData;
54			GenericData* data;
55			bool haveAtomData;
56
57			//if atom is not SSD atom, just skip it
58	tim	1126	if(strcmp(datom->getType(), "SSD"))
59	tim	1118	return;
60
61			data = datom->getProperty("ATOMDATA");
62			if(data != NULL){
63
64			atomData = dynamic_cast<AtomData*>(data);
65			if(atomData == NULL){
66			cerr << "can not get Atom Data from " << datom->getType() << endl;
67			atomData = new AtomData;
68			haveAtomData = false;
69			}
70			else
71			haveAtomData = true;
72			}
73			else{
74			atomData = new AtomData;
75			haveAtomData = false;
76			}
77
78
79			datom->getPos(pos);
80			datom->getQ(q);
81			datom->getA(rotMatrix);
82	tim	1129
83			// We need A^T to convert from body-fixed to space-fixed:
84			transposeMat3(rotMatrix, rotTrans);
85	tim	1118
86			//center of mass of the water molecule
87	tim	1129	matVecMul3(rotTrans, u, newVec);
88	tim	1118	atomInfo = new AtomInfo;
89			atomInfo->AtomType = "X";
90			atomInfo->pos[0] = pos[0];
91			atomInfo->pos[1] = pos[1];
92			atomInfo->pos[2] = pos[2];
93			atomInfo->dipole[0] = newVec[0];
94			atomInfo->dipole[1] = newVec[1];
95			atomInfo->dipole[2] = newVec[2];
96
97			atomData->addAtomInfo(atomInfo);
98
99			//oxygen
100	tim	1129	matVecMul3(rotTrans, ox, newVec);
101	tim	1118	atomInfo = new AtomInfo;
102			atomInfo->AtomType = "O";
103			atomInfo->pos[0] = pos[0] + newVec[0];
104			atomInfo->pos[1] = pos[1] + newVec[1];
105			atomInfo->pos[2] = pos[2] + newVec[2];
106			atomInfo->dipole[0] = 0.0;
107			atomInfo->dipole[1] = 0.0;
108			atomInfo->dipole[2] = 0.0;
109			atomData->addAtomInfo(atomInfo);
110
111
112			//hydrogen1
113	tim	1129	matVecMul3(rotTrans, h1, newVec);
114	tim	1118	atomInfo = new AtomInfo;
115			atomInfo->AtomType = "H";
116			atomInfo->pos[0] = pos[0] + newVec[0];
117			atomInfo->pos[1] = pos[1] + newVec[1];
118			atomInfo->pos[2] = pos[2] + newVec[2];
119			atomInfo->dipole[0] = 0.0;
120			atomInfo->dipole[1] = 0.0;
121			atomInfo->dipole[2] = 0.0;
122			atomData->addAtomInfo(atomInfo);
123
124			//hydrogen2
125	tim	1129	matVecMul3(rotTrans, h2, newVec);
126	tim	1118	atomInfo = new AtomInfo;
127			atomInfo->AtomType = "H";
128			atomInfo->pos[0] = pos[0] + newVec[0];
129			atomInfo->pos[1] = pos[1] + newVec[1];
130			atomInfo->pos[2] = pos[2] + newVec[2];
131			atomInfo->dipole[0] = 0.0;
132			atomInfo->dipole[1] = 0.0;
133			atomInfo->dipole[2] = 0.0;
134			atomData->addAtomInfo(atomInfo);
135
136			//add atom data into atom's property
137
138			if(!haveAtomData){
139			atomData->setID("ATOMDATA");
140			datom->addProperty(atomData);
141			}
142
143			setVisited(datom);
144
145			}
146
147	tim	1126	const string SSDAtomVisitor::toString(){
148			char buffer[65535];
149			string result;
150
151			sprintf(buffer ,"------------------------------------------------------------------\n");
152			result += buffer;
153
154			sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
155			result += buffer;
156
157	tim	1132	sprintf(buffer , "Visitor Description: Convert SSD into 4 different atoms\n");
158	tim	1126	result += buffer;
159
160			sprintf(buffer ,"------------------------------------------------------------------\n");
161			result += buffer;
162
163			return result;
164			}
165
166			//----------------------------------------------------------------------------//
167
168	tim	1118	void DefaultAtomVisitor::visit(Atom* atom){
169			AtomData* atomData;
170			AtomInfo* atomInfo;
171			double pos[3];
172
173			if(isVisited(atom))
174			return;
175
176			atomInfo =new AtomInfo;
177	tim	1119
178			atomData = new AtomData;
179			atomData->setID("ATOMDATA");
180	tim	1118
181			atom->getPos(pos);
182			atomInfo->AtomType = atom->getType();
183			atomInfo->pos[0] = pos[0];
184			atomInfo->pos[1] = pos[1];
185			atomInfo->pos[2] = pos[2];
186			atomInfo->dipole[0] = 0.0;
187			atomInfo->dipole[1] = 0.0;
188			atomInfo->dipole[2] = 0.0;
189
190	tim	1119
191			atomData->addAtomInfo(atomInfo);
192
193	tim	1118	atom->addProperty(atomData);
194
195			setVisited(atom);
196			}
197			void DefaultAtomVisitor::visit(DirectionalAtom* datom){
198			AtomData* atomData;
199			AtomInfo* atomInfo;
200			double pos[3];
201			double u[3];
202
203			if(isVisited(datom))
204			return;
205
206			datom->getPos(pos);
207			datom->getU(u);
208
209	tim	1119	atomData = new AtomData;
210			atomData->setID("ATOMDATA");
211	tim	1118	atomInfo =new AtomInfo;
212
213			atomInfo->AtomType = datom->getType();
214			atomInfo->pos[0] = pos[0];
215			atomInfo->pos[1] = pos[1];
216			atomInfo->pos[2] = pos[2];
217			atomInfo->dipole[0] = u[0];
218			atomInfo->dipole[1] = u[1];
219			atomInfo->dipole[2] = u[2];
220
221	tim	1119	atomData->addAtomInfo(atomInfo);
222
223	tim	1118	datom->addProperty(atomData);
224
225			setVisited(datom);
226			}
227	tim	1126
228
229			const string DefaultAtomVisitor::toString(){
230			char buffer[65535];
231			string result;
232
233			sprintf(buffer ,"------------------------------------------------------------------\n");
234			result += buffer;
235
236			sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
237			result += buffer;
238
239			sprintf(buffer , "Visitor Description: copy atom infomation into atom data\n");
240			result += buffer;
241
242			sprintf(buffer ,"------------------------------------------------------------------\n");
243			result += buffer;
244
245			return result;
246			}