| 1 |
gezelter |
1889 |
#! /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:]) |
