1 |
#! /usr/bin/env python |
2 |
|
3 |
"""Dipolar Lattice Builder |
4 |
|
5 |
Creates cubic lattices of dipoles to test the dipole-dipole |
6 |
interaction code. |
7 |
|
8 |
Usage: buildDipolarArray |
9 |
|
10 |
Options: |
11 |
-h, --help show this help |
12 |
-a, --array-type-A use array type "A" (default) |
13 |
-b, --array-type-B use array type "B" |
14 |
-l, --lattice=... use the specified lattice ( SC, FCC, or BCC ) |
15 |
-d, --direction=... use dipole orientation (001, 111, or 011) |
16 |
-c, --constant=... use the specified lattice constant |
17 |
-n use the specified number of unit cells |
18 |
-o, --output-file=... use specified output (.xyz) file |
19 |
|
20 |
Type "A" arrays have nearest neighbor strings of antiparallel dipoles. |
21 |
|
22 |
Type "B" arrays have nearest neighbor strings of antiparallel dipoles |
23 |
if the dipoles are contained in a plane perpendicular to the dipole |
24 |
direction that passes through the dipole. |
25 |
|
26 |
Example: |
27 |
buildDipolarArray -a -l fcc -d 001 -c 5 -n 3 -o A_fcc_001.xyz |
28 |
|
29 |
""" |
30 |
|
31 |
__author__ = "Dan Gezelter (gezelter@nd.edu)" |
32 |
__version__ = "$Revision: 1639 $" |
33 |
__date__ = "$Date: 2011-09-24 16:18:07 -0400 (Sat, 24 Sep 2011) $" |
34 |
|
35 |
__copyright__ = "Copyright (c) 2013 by the University of Notre Dame" |
36 |
__license__ = "OpenMD" |
37 |
|
38 |
import sys |
39 |
import getopt |
40 |
import string |
41 |
import math |
42 |
import numpy |
43 |
|
44 |
def usage(): |
45 |
print __doc__ |
46 |
|
47 |
def createLattice(latticeType, latticeNumber, latticeConstant, dipoleDirection, doTypeA, outputFileName): |
48 |
a1 = numpy.array([1.0,0.0,0.0]) |
49 |
a2 = numpy.array([0.0,1.0,0.0]) |
50 |
a3 = numpy.array([0.0,0.0,1.0]) |
51 |
|
52 |
sc = numpy.array([[0.0,0.0,0.0]]) |
53 |
|
54 |
bcc = numpy.array([[0.0,0.0,0.0], |
55 |
[0.5,0.5,0.5]]) |
56 |
|
57 |
fcc = numpy.array([[0.0,0.0,0.0], |
58 |
[0.5,0.5,0.0], |
59 |
[0.0,0.5,0.5], |
60 |
[0.5,0.0,0.5]]) |
61 |
|
62 |
dipole001 = numpy.array([0.0,0.0,1.0]) |
63 |
dipole011 = numpy.array([0.0,1.0/math.sqrt(2.0),1.0/math.sqrt(2.0)]) |
64 |
dipole111 = numpy.array([1.0/math.sqrt(3.0),1.0/math.sqrt(3.0),1.0/math.sqrt(3.0)]) |
65 |
|
66 |
if (latticeType.lower() == 'sc'): |
67 |
basis_atoms = sc |
68 |
elif(latticeType.lower() == 'fcc'): |
69 |
basis_atoms = fcc |
70 |
elif(latticeType.lower() == 'bcc'): |
71 |
basis_atoms = bcc |
72 |
else: |
73 |
print "unknown lattice type" |
74 |
print __doc__ |
75 |
|
76 |
ijklvector = numpy.array([1,1,1,1]) |
77 |
|
78 |
if (dipoleDirection.lower() == '001'): |
79 |
basis_dipoles = numpy.array([dipole001, -dipole001]) |
80 |
elif(dipoleDirection.lower() == '011'): |
81 |
basis_dipoles = numpy.array([dipole011, -dipole011]) |
82 |
elif(dipoleDirection.lower() == '111'): |
83 |
basis_dipoles = numpy.array([dipole111, -dipole111]) |
84 |
else: |
85 |
print "unknown lattice type" |
86 |
print __doc__ |
87 |
|
88 |
if (not doTypeA): |
89 |
if (latticeType.lower() == 'sc'): |
90 |
if (dipoleDirection.lower() == '001'): |
91 |
ijkl = numpy.array([1,1,0,0]) |
92 |
elif(dipoleDirection.lower() == '011'): |
93 |
ijkl = numpy.array([1,1,-1,0]) |
94 |
elif(dipoleDirection.lower() == '111'): |
95 |
ijkl = numpy.array([1,1,-2,0]) |
96 |
elif(latticeType.lower() == 'fcc'): |
97 |
if (dipoleDirection.lower() == '001'): |
98 |
ijkl = numpy.array([1,1,0,1]) |
99 |
elif(dipoleDirection.lower() == '011'): |
100 |
ijkl = numpy.array([1,1,-1,2]) |
101 |
elif(dipoleDirection.lower() == '111'): |
102 |
ijkl = numpy.array([1,1,-2,1]) |
103 |
elif(latticeType.lower() == 'bcc'): |
104 |
if (dipoleDirection.lower() == '001'): |
105 |
ijkl = numpy.array([1,1,0,1]) |
106 |
elif(dipoleDirection.lower() == '011'): |
107 |
ijkl = numpy.array([1,1,-1,1]) |
108 |
elif(dipoleDirection.lower() == '111'): |
109 |
ijkl = numpy.array([1,1,-2,0]) |
110 |
else: |
111 |
if (latticeType.lower() == 'sc'): |
112 |
ijkl = numpy.array([1,1,1,0]) |
113 |
elif(latticeType.lower() == 'fcc'): |
114 |
ijkl = numpy.array([1,1,1,1]) |
115 |
elif(latticeType.lower() == 'bcc'): |
116 |
ijkl = numpy.array([0,0,0,1]) |
117 |
|
118 |
bravais_lattice = [] |
119 |
bravais_dipoles = [] |
120 |
|
121 |
for i in range(latticeNumber): |
122 |
for j in range(latticeNumber): |
123 |
for k in range(latticeNumber): |
124 |
for l in range(len(basis_atoms)): |
125 |
bravais_lattice.append(i*a1 + j*a2 + k*a3 + basis_atoms[l]) |
126 |
tester = numpy.array([i,j,k,l]) |
127 |
flip = numpy.dot(tester, ijkl) % 2 |
128 |
print flip |
129 |
bravais_dipoles.append(basis_dipoles[flip]) |
130 |
|
131 |
outputFile = open(outputFileName,'w') |
132 |
|
133 |
outputFile.write(repr(len(bravais_lattice)) + '\n') |
134 |
Hxx = latticeConstant * latticeNumber |
135 |
Hyy = latticeConstant * latticeNumber |
136 |
Hzz = latticeConstant * latticeNumber |
137 |
|
138 |
outputFile.write('Hmat: {{%d, 0, 0}, {0, %d, 0}, {0, 0, %d}}\n' % (Hxx, Hyy, Hzz)) |
139 |
|
140 |
for i in range(len(bravais_lattice)): |
141 |
xcart = latticeConstant*(bravais_lattice[i][0]) |
142 |
ycart = latticeConstant*(bravais_lattice[i][1]) |
143 |
zcart = latticeConstant*(bravais_lattice[i][2]) |
144 |
dx = bravais_dipoles[i][0] |
145 |
dy = bravais_dipoles[i][1] |
146 |
dz = bravais_dipoles[i][2] |
147 |
outputFile.write('Ar ' + repr(xcart) + ' ' + repr(ycart) + ' ' + repr(zcart) + ' ' + repr(dx) + ' ' + repr(dy) + ' ' + repr(dz) + '\n') |
148 |
|
149 |
outputFile.close() |
150 |
|
151 |
def main(argv): |
152 |
|
153 |
doTypeA = True |
154 |
haveOutputFileName = False |
155 |
latticeType = "fcc" |
156 |
dipoleDirection = "001" |
157 |
latticeNumber = 3 |
158 |
latticeConstant = 4 |
159 |
try: |
160 |
opts, args = getopt.getopt(argv, "habl:d:c:n:o:", ["help","array-type-A", "array-type-B", "lattice=" "direction=", "constant=", "output-file="]) |
161 |
except getopt.GetoptError: |
162 |
usage() |
163 |
sys.exit(2) |
164 |
for opt, arg in opts: |
165 |
if opt in ("-h", "--help"): |
166 |
usage() |
167 |
sys.exit() |
168 |
elif opt in ("-a", "--array-type-A"): |
169 |
doTypeA = True |
170 |
elif opt in ("-b", "--array-type-B"): |
171 |
doTypeA = False |
172 |
elif opt in ("-l", "--lattice"): |
173 |
latticeType = arg |
174 |
elif opt in ("-d", "--direction"): |
175 |
dipoleDirection = arg |
176 |
elif opt in ("-c", "--constant"): |
177 |
latticeConstant = float(arg) |
178 |
elif opt in ("-n"): |
179 |
latticeNumber = int(arg) |
180 |
elif opt in ("-o", "--output-file"): |
181 |
outputFileName = arg |
182 |
haveOutputFileName = True |
183 |
if (not haveOutputFileName): |
184 |
usage() |
185 |
print "No output file was specified" |
186 |
sys.exit() |
187 |
|
188 |
createLattice(latticeType, latticeNumber, latticeConstant, dipoleDirection, doTypeA, outputFileName); |
189 |
|
190 |
if __name__ == "__main__": |
191 |
if len(sys.argv) == 1: |
192 |
usage() |
193 |
sys.exit() |
194 |
main(sys.argv[1:]) |