1 |
chuckv |
678 |
#include "latticeBuilder.hpp" |
2 |
|
|
#include <cmath> |
3 |
|
|
#include <cstdlib> |
4 |
|
|
#include "simError.h" |
5 |
|
|
|
6 |
|
|
|
7 |
chuckv |
700 |
Lattice::Lattice(int latticeType, double latticeParameter){ |
8 |
chuckv |
678 |
int hasError; |
9 |
|
|
|
10 |
|
|
latticePosX = NULL; |
11 |
|
|
latticePosY = NULL; |
12 |
|
|
latticePosZ = NULL; |
13 |
|
|
|
14 |
gezelter |
817 |
startX = 0.0; |
15 |
|
|
startY = 0.0; |
16 |
|
|
startZ = 0.0; |
17 |
chuckv |
678 |
|
18 |
|
|
switch(latticeType){ |
19 |
|
|
|
20 |
|
|
case FCC_LATTICE_TYPE: |
21 |
chuckv |
700 |
hasError = createFccLattice(latticeParameter); |
22 |
chuckv |
678 |
break; |
23 |
|
|
|
24 |
|
|
case BCC_LATTICE_TYPE: |
25 |
chuckv |
700 |
hasError = createBccLattice(latticeParameter); |
26 |
chuckv |
678 |
break; |
27 |
|
|
|
28 |
|
|
case HCP_LATTICE_TYPE: |
29 |
chuckv |
700 |
hasError = createHcpLattice(latticeParameter); |
30 |
chuckv |
678 |
break; |
31 |
|
|
|
32 |
|
|
case HCPWATER_LATTICE_TYPE: |
33 |
chuckv |
700 |
hasError = createHcpWaterLattice(latticeParameter); |
34 |
chuckv |
678 |
break; |
35 |
|
|
|
36 |
gezelter |
817 |
case ORTHORHOMBIC_LATTICE_TYPE: |
37 |
|
|
hasError = createFccLattice(latticeParameter); |
38 |
|
|
break; |
39 |
|
|
|
40 |
|
|
|
41 |
chuckv |
678 |
default: |
42 |
chuckv |
700 |
//simerror here..... |
43 |
|
|
; |
44 |
chuckv |
678 |
} |
45 |
chuckv |
700 |
|
46 |
|
|
latticePosX = new double[nCellSites]; |
47 |
|
|
latticePosY = new double[nCellSites]; |
48 |
|
|
latticePosZ = new double[nCellSites]; |
49 |
|
|
|
50 |
chuckv |
678 |
return; |
51 |
|
|
} |
52 |
|
|
|
53 |
gezelter |
817 |
Lattice::Lattice(int latticeType, double latticeParameters[3]){ |
54 |
|
|
int hasError; |
55 |
chuckv |
678 |
|
56 |
gezelter |
817 |
latticePosX = NULL; |
57 |
|
|
latticePosY = NULL; |
58 |
|
|
latticePosZ = NULL; |
59 |
|
|
|
60 |
|
|
startX = 0.0; |
61 |
|
|
startY = 0.0; |
62 |
|
|
startZ = 0.0; |
63 |
|
|
|
64 |
|
|
switch(latticeType){ |
65 |
|
|
|
66 |
|
|
case ORTHORHOMBIC_LATTICE_TYPE: |
67 |
|
|
hasError = createOrthorhombicLattice(latticeParameters); |
68 |
|
|
break; |
69 |
|
|
|
70 |
|
|
|
71 |
|
|
default: |
72 |
|
|
//simerror here..... |
73 |
|
|
; |
74 |
|
|
} |
75 |
|
|
|
76 |
|
|
latticePosX = new double[nCellSites]; |
77 |
|
|
latticePosY = new double[nCellSites]; |
78 |
|
|
latticePosZ = new double[nCellSites]; |
79 |
|
|
|
80 |
|
|
return; |
81 |
|
|
} |
82 |
|
|
|
83 |
|
|
|
84 |
chuckv |
678 |
Lattice::~Lattice(void){ |
85 |
|
|
if (latticePosX != NULL) delete[] latticePosX; |
86 |
|
|
if (latticePosY != NULL) delete[] latticePosY; |
87 |
|
|
if (latticePosZ != NULL) delete[] latticePosZ; |
88 |
|
|
} |
89 |
|
|
|
90 |
|
|
|
91 |
chuckv |
700 |
int Lattice::createFccLattice(double latticeParameter){ |
92 |
chuckv |
678 |
double cell2; |
93 |
|
|
double rroot3; |
94 |
|
|
|
95 |
gezelter |
817 |
cellLengthX = latticeParameter; |
96 |
|
|
cellLengthY = latticeParameter; |
97 |
|
|
cellLengthZ = latticeParameter; |
98 |
chuckv |
700 |
|
99 |
|
|
cell2 = 0.5 * latticeParameter; |
100 |
chuckv |
678 |
rroot3 = 1.0 / sqrt(3.0); |
101 |
|
|
|
102 |
|
|
nCellSites = 4; |
103 |
|
|
// create new unit cells |
104 |
|
|
thisUnitCell.sx = new double[nCellSites]; |
105 |
|
|
thisUnitCell.sy = new double[nCellSites]; |
106 |
|
|
thisUnitCell.sz = new double[nCellSites]; |
107 |
|
|
thisUnitCell.s_ex = new double[nCellSites]; |
108 |
|
|
thisUnitCell.s_ey = new double[nCellSites]; |
109 |
|
|
thisUnitCell.s_ez = new double[nCellSites]; |
110 |
|
|
|
111 |
|
|
// add members to each unit cell |
112 |
|
|
// Molecule 1 |
113 |
|
|
thisUnitCell.sx[0] = 0.0; |
114 |
|
|
thisUnitCell.sy[0] = 0.0; |
115 |
|
|
thisUnitCell.sz[0] = 0.0; |
116 |
|
|
thisUnitCell.s_ex[0] = rroot3; |
117 |
|
|
thisUnitCell.s_ey[0] = rroot3; |
118 |
|
|
thisUnitCell.s_ez[0] = rroot3; |
119 |
|
|
|
120 |
|
|
// Molecule 2 |
121 |
|
|
thisUnitCell.sx[1] = 0.0; |
122 |
|
|
thisUnitCell.sy[1] = cell2; |
123 |
|
|
thisUnitCell.sz[1] = cell2; |
124 |
|
|
thisUnitCell.s_ex[1] = -rroot3; |
125 |
|
|
thisUnitCell.s_ey[1] = rroot3; |
126 |
|
|
thisUnitCell.s_ez[1] = -rroot3; |
127 |
|
|
|
128 |
|
|
// Molecule 3 |
129 |
|
|
thisUnitCell.sx[2] = cell2; |
130 |
|
|
thisUnitCell.sy[2] = cell2; |
131 |
|
|
thisUnitCell.sz[2] = 0.0; |
132 |
|
|
thisUnitCell.s_ex[2] = rroot3; |
133 |
|
|
thisUnitCell.s_ey[2] = -rroot3; |
134 |
|
|
thisUnitCell.s_ez[2] = -rroot3; |
135 |
|
|
|
136 |
|
|
// Molecule 4 |
137 |
|
|
thisUnitCell.sx[3] = cell2; |
138 |
|
|
thisUnitCell.sy[3] = 0.0; |
139 |
|
|
thisUnitCell.sz[3] = cell2; |
140 |
|
|
thisUnitCell.s_ex[3] = -rroot3; |
141 |
|
|
thisUnitCell.s_ey[3] = -rroot3; |
142 |
|
|
thisUnitCell.s_ez[3] = rroot3; |
143 |
|
|
|
144 |
|
|
return 0; |
145 |
|
|
} |
146 |
|
|
|
147 |
|
|
|
148 |
gezelter |
817 |
int Lattice::createOrthorhombicLattice(double latticeParameters[3]){ |
149 |
mmeineke |
821 |
double cellx2, celly2, cellz2; |
150 |
gezelter |
817 |
double rroot3; |
151 |
|
|
|
152 |
mmeineke |
821 |
cellLengthX = latticeParameters[0]; |
153 |
|
|
cellLengthY = latticeParameters[1]; |
154 |
|
|
cellLengthZ = latticeParameters[2]; |
155 |
gezelter |
817 |
|
156 |
mmeineke |
821 |
cellx2 = 0.5 * latticeParameters[0]; |
157 |
|
|
celly2 = 0.5 * latticeParameters[1]; |
158 |
|
|
cellz2 = 0.5 * latticeParameters[2]; |
159 |
gezelter |
817 |
rroot3 = 1.0 / sqrt(3.0); |
160 |
|
|
|
161 |
|
|
nCellSites = 4; |
162 |
|
|
// create new unit cells |
163 |
|
|
thisUnitCell.sx = new double[nCellSites]; |
164 |
|
|
thisUnitCell.sy = new double[nCellSites]; |
165 |
|
|
thisUnitCell.sz = new double[nCellSites]; |
166 |
|
|
thisUnitCell.s_ex = new double[nCellSites]; |
167 |
|
|
thisUnitCell.s_ey = new double[nCellSites]; |
168 |
|
|
thisUnitCell.s_ez = new double[nCellSites]; |
169 |
|
|
|
170 |
|
|
// add members to each unit cell |
171 |
|
|
// Molecule 1 |
172 |
|
|
thisUnitCell.sx[0] = 0.0; |
173 |
|
|
thisUnitCell.sy[0] = 0.0; |
174 |
|
|
thisUnitCell.sz[0] = 0.0; |
175 |
|
|
thisUnitCell.s_ex[0] = rroot3; |
176 |
|
|
thisUnitCell.s_ey[0] = rroot3; |
177 |
|
|
thisUnitCell.s_ez[0] = rroot3; |
178 |
|
|
|
179 |
|
|
// Molecule 2 |
180 |
|
|
thisUnitCell.sx[1] = 0.0; |
181 |
|
|
thisUnitCell.sy[1] = celly2; |
182 |
|
|
thisUnitCell.sz[1] = cellz2; |
183 |
|
|
thisUnitCell.s_ex[1] = -rroot3; |
184 |
|
|
thisUnitCell.s_ey[1] = rroot3; |
185 |
|
|
thisUnitCell.s_ez[1] = -rroot3; |
186 |
|
|
|
187 |
|
|
// Molecule 3 |
188 |
|
|
thisUnitCell.sx[2] = cellx2; |
189 |
|
|
thisUnitCell.sy[2] = celly2; |
190 |
|
|
thisUnitCell.sz[2] = 0.0; |
191 |
|
|
thisUnitCell.s_ex[2] = rroot3; |
192 |
|
|
thisUnitCell.s_ey[2] = -rroot3; |
193 |
|
|
thisUnitCell.s_ez[2] = -rroot3; |
194 |
|
|
|
195 |
|
|
// Molecule 4 |
196 |
|
|
thisUnitCell.sx[3] = cellx2; |
197 |
|
|
thisUnitCell.sy[3] = 0.0; |
198 |
|
|
thisUnitCell.sz[3] = cellz2; |
199 |
|
|
thisUnitCell.s_ex[3] = -rroot3; |
200 |
|
|
thisUnitCell.s_ey[3] = -rroot3; |
201 |
|
|
thisUnitCell.s_ez[3] = rroot3; |
202 |
|
|
|
203 |
|
|
return 0; |
204 |
|
|
} |
205 |
|
|
|
206 |
|
|
|
207 |
chuckv |
678 |
// Body centered cubic lattice |
208 |
chuckv |
700 |
int Lattice::createBccLattice(double latticeParameter){ |
209 |
chuckv |
678 |
return 0; |
210 |
|
|
} |
211 |
|
|
|
212 |
|
|
|
213 |
|
|
// Standard HCP lattice |
214 |
chuckv |
700 |
int Lattice::createHcpLattice(double latticeParameter){ |
215 |
chuckv |
678 |
return 0; |
216 |
|
|
} |
217 |
|
|
|
218 |
|
|
// HCP contains tetrahedral sites for waters |
219 |
chuckv |
700 |
int Lattice::createHcpWaterLattice(double latticeParameter){ |
220 |
|
|
|
221 |
|
|
return 0; |
222 |
|
|
|
223 |
chuckv |
678 |
double rroot3; |
224 |
|
|
double cell; |
225 |
|
|
double cell2; |
226 |
|
|
double cell4; |
227 |
|
|
double cell16; |
228 |
|
|
double cell23; |
229 |
|
|
double cell34; |
230 |
|
|
double cell09; |
231 |
|
|
double cell59; |
232 |
|
|
|
233 |
|
|
nCellSites = 16; |
234 |
|
|
|
235 |
|
|
// create new unit cells |
236 |
|
|
thisUnitCell.sx = new double[nCellSites]; |
237 |
|
|
thisUnitCell.sy = new double[nCellSites]; |
238 |
|
|
thisUnitCell.sz = new double[nCellSites]; |
239 |
|
|
thisUnitCell.s_ex = new double[nCellSites]; |
240 |
|
|
thisUnitCell.s_ey = new double[nCellSites]; |
241 |
|
|
thisUnitCell.s_ez = new double[nCellSites]; |
242 |
|
|
|
243 |
|
|
rroot3 = 1.0 / sqrt(3.0); |
244 |
chuckv |
700 |
cell = latticeParameter; |
245 |
chuckv |
678 |
cell2 = 0.5 * cell; |
246 |
|
|
cell4 = 0.25 * cell; |
247 |
|
|
cell16 = cell / 6.0; |
248 |
|
|
cell23 = (2.0 * cell) / 3.0; |
249 |
|
|
cell34 = 0.75 * cell; |
250 |
|
|
cell09 = 0.11238 * cell; |
251 |
|
|
cell59 = 0.61238 * cell; |
252 |
|
|
|
253 |
|
|
// ** build the unit cell ** |
254 |
|
|
|
255 |
|
|
// molecule 1 |
256 |
|
|
|
257 |
|
|
thisUnitCell.sx[0] = 0.0; |
258 |
|
|
thisUnitCell.sy[0] = 0.0; |
259 |
|
|
thisUnitCell.sz[0] = 0.0; |
260 |
|
|
|
261 |
|
|
// molecule 2 |
262 |
|
|
|
263 |
|
|
thisUnitCell.sx[1] = cell4; |
264 |
|
|
thisUnitCell.sy[1] = cell2; |
265 |
|
|
thisUnitCell.sz[1] = 0.0; |
266 |
|
|
|
267 |
|
|
// molecule 3 |
268 |
|
|
|
269 |
|
|
thisUnitCell.sx[2] = cell4; |
270 |
|
|
thisUnitCell.sy[2] = cell16; |
271 |
|
|
thisUnitCell.sz[2] = cell09; |
272 |
|
|
|
273 |
|
|
// molecule 4 |
274 |
|
|
|
275 |
|
|
thisUnitCell.sx[3] = 0.0; |
276 |
|
|
thisUnitCell.sy[3] = cell23; |
277 |
|
|
thisUnitCell.sz[3] = cell09; |
278 |
|
|
|
279 |
|
|
// molecule 5 |
280 |
|
|
|
281 |
|
|
thisUnitCell.sx[4] = cell4; |
282 |
|
|
thisUnitCell.sy[4] = cell16; |
283 |
|
|
thisUnitCell.sz[4] = cell2; |
284 |
|
|
|
285 |
|
|
// molecule 6 |
286 |
|
|
|
287 |
|
|
thisUnitCell.sx[5] = 0.0; |
288 |
|
|
thisUnitCell.sy[5] = cell23; |
289 |
|
|
thisUnitCell.sz[5] = cell2; |
290 |
|
|
|
291 |
|
|
// molecule 7 |
292 |
|
|
|
293 |
|
|
thisUnitCell.sx[6] = 0.0; |
294 |
|
|
thisUnitCell.sy[6] = 0.0; |
295 |
|
|
thisUnitCell.sz[6] = cell59; |
296 |
|
|
|
297 |
|
|
|
298 |
|
|
// molecule 8 |
299 |
|
|
|
300 |
|
|
thisUnitCell.sx[7] = cell4; |
301 |
|
|
thisUnitCell.sy[7] = cell2; |
302 |
|
|
thisUnitCell.sz[7] = cell59; |
303 |
|
|
|
304 |
|
|
// molecule 9 |
305 |
|
|
|
306 |
|
|
thisUnitCell.sz[8] = cell2; |
307 |
|
|
thisUnitCell.sy[8] = 0.0; |
308 |
|
|
thisUnitCell.sz[8] = 0.0; |
309 |
|
|
|
310 |
|
|
// molecule 10 |
311 |
|
|
|
312 |
|
|
thisUnitCell.sz[9] = cell34; |
313 |
|
|
thisUnitCell.sy[9] = cell2; |
314 |
|
|
thisUnitCell.sz[9] = 0.0; |
315 |
|
|
|
316 |
|
|
// molecule 11 |
317 |
|
|
|
318 |
|
|
thisUnitCell.sz[10] = cell34; |
319 |
|
|
thisUnitCell.sy[10] = cell16; |
320 |
|
|
thisUnitCell.sz[10] = cell09; |
321 |
|
|
|
322 |
|
|
// molecule 12 |
323 |
|
|
|
324 |
|
|
thisUnitCell.sz[11] = cell2; |
325 |
|
|
thisUnitCell.sy[11] = cell23; |
326 |
|
|
thisUnitCell.sz[11] = cell09; |
327 |
|
|
|
328 |
|
|
// molecule 13 |
329 |
|
|
|
330 |
|
|
thisUnitCell.sz[12] = cell34; |
331 |
|
|
thisUnitCell.sy[12] = cell16; |
332 |
|
|
thisUnitCell.sz[12] = cell2; |
333 |
|
|
|
334 |
|
|
// molecule 14 |
335 |
|
|
|
336 |
|
|
thisUnitCell.sz[13]= cell2; |
337 |
|
|
thisUnitCell.sy[13]= cell23; |
338 |
|
|
thisUnitCell.sz[13]= cell2; |
339 |
|
|
|
340 |
|
|
// molecule 15 |
341 |
|
|
|
342 |
|
|
thisUnitCell.sz[14] = cell2; |
343 |
|
|
thisUnitCell.sy[14] = 0.0; |
344 |
|
|
thisUnitCell.sz[14] = cell59; |
345 |
|
|
|
346 |
|
|
// molecule 16 |
347 |
|
|
|
348 |
|
|
thisUnitCell.sz[15] = cell34; |
349 |
|
|
thisUnitCell.sy[15] = cell2; |
350 |
|
|
thisUnitCell.sz[15] = cell59; |
351 |
|
|
|
352 |
|
|
return 0; |
353 |
|
|
} |
354 |
|
|
|
355 |
|
|
|
356 |
|
|
//Returns lattice points when called repeatedly |
357 |
|
|
int Lattice::getLatticePoints(double** thePosX, double** thePosY, |
358 |
chuckv |
700 |
double** thePosZ, |
359 |
|
|
int ix, int iy, int iz){ |
360 |
chuckv |
678 |
|
361 |
|
|
int iref; |
362 |
|
|
|
363 |
|
|
for( iref=0;iref < nCellSites;iref++){ |
364 |
|
|
|
365 |
gezelter |
817 |
latticePosX[iref] = startX + thisUnitCell.sx[iref] + |
366 |
|
|
cellLengthX * (double( ix ) - 0.5); |
367 |
|
|
latticePosY[iref] = startY + thisUnitCell.sy[iref] + |
368 |
|
|
cellLengthY * (double( iy ) - 0.5); |
369 |
|
|
latticePosZ[iref] = startZ + thisUnitCell.sz[iref] + |
370 |
|
|
cellLengthZ * (double( iz ) - 0.5); |
371 |
chuckv |
678 |
|
372 |
|
|
} |
373 |
|
|
|
374 |
|
|
*thePosX = latticePosX; |
375 |
|
|
*thePosY = latticePosY; |
376 |
|
|
*thePosZ = latticePosZ; |
377 |
|
|
|
378 |
|
|
return nCellSites; |
379 |
|
|
} |
380 |
|
|
|