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 (21 years ago) by mmeineke
File size:	4777 byte(s)
Log Message:	Most of the integrator and NVT seem to be working now.
#	Content
1	#include <iostream>
2
3	using namespace std;
4
5	#include "Atom.hpp"
6
7	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	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	void Atom::createArrays (int the_nElements) {
42	int i;
43
44	nElements = the_nElements;
45
46	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	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	}
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	new_Amat[j] = AAmat[i];
149	}
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	new_mu[j] = Amu[i];
158	}
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	}