OOPSE/libmdtools/Atom.cpp

#include <iostream>

using namespace std;

#include "Atom.hpp"

double* Atom::pos; // the position array
double* Atom::vel; // the velocity array
double* Atom::frc; // the forc array
double* Atom::trq; // the torque vector  ( space fixed )
double* Atom::Amat; // the rotation matrix
double* Atom::mu;   // the array of dipole moments
double* Atom::ul;   // the lab frame unit directional vector
int Atom::nElements;

Atom::Atom(int theIndex) {
  c_n_hyd = 0; 
  has_dipole = 0;
  is_VDW = 0;
  is_LJ = 0;
  
  index = theIndex;
  offset = 3 * index;
  offsetX = offset;
  offsetY = offset+1;
  offsetZ = offset+2;
  
  Axx = index*9;
  Axy = Axx+1;
  Axz = Axx+2;
  
  Ayx = Axx+3;
  Ayy = Axx+4;
  Ayz = Axx+5;
  
  Azx = Axx+6;
  Azy = Axx+7;
  Azz = Axx+8;
}

void Atom::createArrays (int the_nElements) {
  int i;
  
  nElements = the_nElements;

  pos = new double[nElements*3];
  vel = new double[nElements*3];
  frc = new double[nElements*3];
  trq = new double[nElements*3];
  Amat = new double[nElements*9];
  mu = new double[nElements];
  ul = new double[nElements*3];
  
  // init directional values to zero
  
  for( i=0; i<nElements; i++){
    trq[i] = 0.0;
    trq[i+1] = 0.0;
    trq[i+2] = 0.0;
    
    Amat[i] = 1.0;
    Amat[i+1] = 0.0;
    Amat[i+2] = 0.0;
    
    Amat[i+3] = 0.0;
    Amat[i+4] = 1.0;
    Amat[i+5] = 0.0;
    
    Amat[i+6] = 0.0;
    Amat[i+7] = 0.0;
    Amat[i+8] = 1.0;
    
    mu[i] = 0.0;    
    
    ul[i] = 1.0;
    ul[i+1] = 0.0;
    ul[i+2] = 0.0;
  }
}

void Atom::destroyArrays(void) {
  delete[] pos;
  delete[] vel;
  delete[] frc;
  delete[] trq;
  delete[] Amat;
  delete[] mu;
}

void Atom::setIndex(int theIndex) {
  index = theIndex; 
  offset = index*3;
  offsetX = offset;
  offsetY = offset+1;
  offsetZ = offset+2;
  
  Axx = index*9;
  Axy = Axx+1;
  Axz = Axx+2;
  
  Ayx = Axx+3;
  Ayy = Axx+4;
  Ayz = Axx+5;
  
  Azx = Axx+6;
  Azy = Axx+7;
  Azz = Axx+8;
}

void Atom::addAtoms(int nAdded, double* Apos, double* Avel, double* Afrc, 
                   double* Atrq, double* AAmat, double* Amu, 
                   double* Aul) {

  int nNew = nElements+nAdded;

  double* new_pos = new double[nNew*3];
  double* new_vel = new double[nNew*3];
  double* new_frc = new double[nNew*3];
  double* new_trq = new double[nNew*3];
  double* new_Amat = new double[nNew*9];
  double* new_mu = new double[nNew];
  double* new_ul = new double[nNew*3];
  int i, j;
  
  for (i = 0; i < 3*nElements; i++) {
    new_pos[i] = pos[i];
    new_vel[i] = vel[i];
    new_frc[i] = frc[i];
    new_trq[i] = trq[i];
    new_ul[i] = ul[i];
  }

  for(i = 0; i < 3*nAdded; i++) {
    j = i + 3*nElements;
    new_pos[j] = Apos[i];
    new_vel[j] = Avel[i];
    new_frc[j] = Afrc[i];
    new_trq[j] = Atrq[i];
    new_ul[j] = Aul[i]; 
  }

  for (i = 0; i < 9*nElements; i++) {
    new_Amat[i] = Amat[i];
  }

  for(i = 0; i < 9*nAdded; i++) {
    j = i + 9*nElements;
    new_Amat[j] = AAmat[i];
  }

  for (i = 0; i < nElements; i++) {
    new_mu[i] = mu[i];
  }

  for(i = 0; i < nAdded; i++) {
    j = i + nElements;
    new_mu[j] = Amu[i];
  }

  delete[] pos;
  delete[] vel;
  delete[] frc;
  delete[] trq;
  delete[] Amat;
  delete[] mu;

  pos = new_pos;
  vel = new_vel;
  frc = new_frc;
  trq = new_trq;
  ul = new_ul;
  Amat = new_Amat;
  mu = new_mu;

  nElements = nNew;
}

void Atom::deleteAtom(int theIndex) {
  deleteRange(theIndex, theIndex);
}

void Atom::deleteRange(int startIndex, int stopIndex) {
 
  int nNew = nElements-(stopIndex-startIndex+1);

  double* new_pos = new double[nNew*3];
  double* new_vel = new double[nNew*3];
  double* new_frc = new double[nNew*3];
  double* new_trq = new double[nNew*3];
  double* new_Amat = new double[nNew*9];
  double* new_mu = new double[nNew];
  double* new_ul = new double[nNew*3];
  int i, j;
  
  for (i = 0; i < 3*startIndex; i++) {
    new_pos[i] = pos[i];
    new_vel[i] = vel[i];
    new_frc[i] = frc[i];
    new_trq[i] = trq[i];
    new_ul[i] = ul[i];
  }

  for(i = 3*(stopIndex + 1); i < 3*nElements; i++) {
    j = i - 3*startIndex + 1;
    new_pos[j] = pos[i];
    new_vel[j] = vel[i];
    new_frc[j] = frc[i];
    new_trq[j] = trq[i];
    new_ul[j] = ul[i]; 
  }

  for (i = 0; i < 9*startIndex; i++) {
    new_Amat[i] = Amat[i];
  }

  for(i = 9*(stopIndex + 1); i < 9*nElements; i++) {
    j = i - 9*startIndex + 1;
    new_Amat[j] = Amat[i];
  }

  for (i = 0; i < startIndex; i++) {
    new_mu[i] = mu[i];
  }

  for(i = (stopIndex+1); i < nElements; i++) {
    j = i - startIndex + 1;
    new_mu[j] = mu[i];
  }

  delete[] pos;
  delete[] vel;
  delete[] frc;
  delete[] trq;
  delete[] Amat;
  delete[] mu;

  pos = new_pos;
  vel = new_vel;
  frc = new_frc;
  trq = new_trq;
  ul = new_ul;
  Amat = new_Amat;
  mu = new_mu;

  nElements = nNew;
}
Revision:	561
Committed:	Fri Jun 20 20:29:36 2003 UTC (22 years, 4 months ago) by mmeineke
Original Path:	*trunk/OOPSE/libmdtools/Atom.cpp*
File size:	4777 byte(s)
Log Message:	Most of the integrator and NVT seem to be working now.
#	User	Rev	Content
1	mmeineke	561	#include <iostream>
2
3			using namespace std;
4
5	gezelter	409	#include "Atom.hpp"
6
7	mmeineke	413	double* Atom::pos; // the position array
8			double* Atom::vel; // the velocity array
9			double* Atom::frc; // the forc array
10			double* Atom::trq; // the torque vector ( space fixed )
11			double* Atom::Amat; // the rotation matrix
12			double* Atom::mu; // the array of dipole moments
13			double* Atom::ul; // the lab frame unit directional vector
14			int Atom::nElements;
15
16	gezelter	409	Atom::Atom(int theIndex) {
17			c_n_hyd = 0;
18			has_dipole = 0;
19			is_VDW = 0;
20			is_LJ = 0;
21
22			index = theIndex;
23			offset = 3 * index;
24			offsetX = offset;
25			offsetY = offset+1;
26			offsetZ = offset+2;
27
28			Axx = index*9;
29			Axy = Axx+1;
30			Axz = Axx+2;
31
32			Ayx = Axx+3;
33			Ayy = Axx+4;
34			Ayz = Axx+5;
35
36			Azx = Axx+6;
37			Azy = Axx+7;
38			Azz = Axx+8;
39			}
40
41	mmeineke	413	void Atom::createArrays (int the_nElements) {
42	gezelter	409	int i;
43
44	mmeineke	413	nElements = the_nElements;
45
46	gezelter	409	pos = new double[nElements*3];
47			vel = new double[nElements*3];
48			frc = new double[nElements*3];
49			trq = new double[nElements*3];
50			Amat = new double[nElements*9];
51			mu = new double[nElements];
52			ul = new double[nElements*3];
53
54			// init directional values to zero
55
56			for( i=0; i<nElements; i++){
57			trq[i] = 0.0;
58			trq[i+1] = 0.0;
59			trq[i+2] = 0.0;
60
61			Amat[i] = 1.0;
62			Amat[i+1] = 0.0;
63			Amat[i+2] = 0.0;
64
65			Amat[i+3] = 0.0;
66			Amat[i+4] = 1.0;
67			Amat[i+5] = 0.0;
68
69			Amat[i+6] = 0.0;
70			Amat[i+7] = 0.0;
71			Amat[i+8] = 1.0;
72
73			mu[i] = 0.0;
74
75			ul[i] = 1.0;
76			ul[i+1] = 0.0;
77			ul[i+2] = 0.0;
78			}
79			}
80
81			void Atom::destroyArrays(void) {
82			delete[] pos;
83			delete[] vel;
84			delete[] frc;
85			delete[] trq;
86			delete[] Amat;
87			delete[] mu;
88			}
89
90			void Atom::setIndex(int theIndex) {
91			index = theIndex;
92			offset = index*3;
93			offsetX = offset;
94			offsetY = offset+1;
95			offsetZ = offset+2;
96
97			Axx = index*9;
98			Axy = Axx+1;
99			Axz = Axx+2;
100
101			Ayx = Axx+3;
102			Ayy = Axx+4;
103			Ayz = Axx+5;
104
105			Azx = Axx+6;
106			Azy = Axx+7;
107			Azz = Axx+8;
108			}
109
110			void Atom::addAtoms(int nAdded, double* Apos, double* Avel, double* Afrc,
111			double* Atrq, double* AAmat, double* Amu,
112			double* Aul) {
113
114			int nNew = nElements+nAdded;
115
116			double* new_pos = new double[nNew*3];
117			double* new_vel = new double[nNew*3];
118			double* new_frc = new double[nNew*3];
119			double* new_trq = new double[nNew*3];
120			double* new_Amat = new double[nNew*9];
121			double* new_mu = new double[nNew];
122			double* new_ul = new double[nNew*3];
123			int i, j;
124
125			for (i = 0; i < 3*nElements; i++) {
126			new_pos[i] = pos[i];
127			new_vel[i] = vel[i];
128			new_frc[i] = frc[i];
129			new_trq[i] = trq[i];
130			new_ul[i] = ul[i];
131			}
132
133			for(i = 0; i < 3*nAdded; i++) {
134			j = i + 3*nElements;
135	gezelter	410	new_pos[j] = Apos[i];
136			new_vel[j] = Avel[i];
137			new_frc[j] = Afrc[i];
138			new_trq[j] = Atrq[i];
139			new_ul[j] = Aul[i];
140	gezelter	409	}
141
142			for (i = 0; i < 9*nElements; i++) {
143			new_Amat[i] = Amat[i];
144			}
145
146			for(i = 0; i < 9*nAdded; i++) {
147			j = i + 9*nElements;
148	gezelter	410	new_Amat[j] = AAmat[i];
149	gezelter	409	}
150
151			for (i = 0; i < nElements; i++) {
152			new_mu[i] = mu[i];
153			}
154
155			for(i = 0; i < nAdded; i++) {
156			j = i + nElements;
157	gezelter	410	new_mu[j] = Amu[i];
158	gezelter	409	}
159
160			delete[] pos;
161			delete[] vel;
162			delete[] frc;
163			delete[] trq;
164			delete[] Amat;
165			delete[] mu;
166
167			pos = new_pos;
168			vel = new_vel;
169			frc = new_frc;
170			trq = new_trq;
171			ul = new_ul;
172			Amat = new_Amat;
173			mu = new_mu;
174
175			nElements = nNew;
176			}
177
178			void Atom::deleteAtom(int theIndex) {
179			deleteRange(theIndex, theIndex);
180			}
181
182			void Atom::deleteRange(int startIndex, int stopIndex) {
183
184			int nNew = nElements-(stopIndex-startIndex+1);
185
186			double* new_pos = new double[nNew*3];
187			double* new_vel = new double[nNew*3];
188			double* new_frc = new double[nNew*3];
189			double* new_trq = new double[nNew*3];
190			double* new_Amat = new double[nNew*9];
191			double* new_mu = new double[nNew];
192			double* new_ul = new double[nNew*3];
193			int i, j;
194
195			for (i = 0; i < 3*startIndex; i++) {
196			new_pos[i] = pos[i];
197			new_vel[i] = vel[i];
198			new_frc[i] = frc[i];
199			new_trq[i] = trq[i];
200			new_ul[i] = ul[i];
201			}
202
203			for(i = 3(stopIndex + 1); i < 3nElements; i++) {
204			j = i - 3*startIndex + 1;
205			new_pos[j] = pos[i];
206			new_vel[j] = vel[i];
207			new_frc[j] = frc[i];
208			new_trq[j] = trq[i];
209			new_ul[j] = ul[i];
210			}
211
212			for (i = 0; i < 9*startIndex; i++) {
213			new_Amat[i] = Amat[i];
214			}
215
216			for(i = 9(stopIndex + 1); i < 9nElements; i++) {
217			j = i - 9*startIndex + 1;
218			new_Amat[j] = Amat[i];
219			}
220
221			for (i = 0; i < startIndex; i++) {
222			new_mu[i] = mu[i];
223			}
224
225			for(i = (stopIndex+1); i < nElements; i++) {
226			j = i - startIndex + 1;
227			new_mu[j] = mu[i];
228			}
229
230			delete[] pos;
231			delete[] vel;
232			delete[] frc;
233			delete[] trq;
234			delete[] Amat;
235			delete[] mu;
236
237			pos = new_pos;
238			vel = new_vel;
239			frc = new_frc;
240			trq = new_trq;
241			ul = new_ul;
242			Amat = new_Amat;
243			mu = new_mu;
244
245			nElements = nNew;
246			}