1 |
chrisfen |
2996 |
#!@PERLINTERP@ -w |
2 |
|
|
|
3 |
|
|
# program that builds water boxes |
4 |
|
|
|
5 |
|
|
# author = "Chris Fennell |
6 |
chrisfen |
3031 |
# version = "$Revision: 1.4 $" |
7 |
|
|
# date = "$Date: 2006-10-02 23:54:30 $" |
8 |
chrisfen |
2996 |
# copyright = "Copyright (c) 2006 by the University of Notre Dame" |
9 |
|
|
# license = "OOPSE" |
10 |
|
|
|
11 |
|
|
use Getopt::Std; |
12 |
|
|
|
13 |
|
|
$tolerance = 1.0E-8; |
14 |
|
|
$mass = 2.99151E-23; # mass of H2O in grams |
15 |
|
|
$cm3ToAng3 = 1E24; # convert cm^3 to angstroms^3 |
16 |
|
|
$densityConvert = $mass*$cm3ToAng3; |
17 |
|
|
$lattice = 0; |
18 |
|
|
$nMol = 500; |
19 |
|
|
$density = 1.0; |
20 |
|
|
$doRandomize = 0; |
21 |
|
|
$cutoff = 12; |
22 |
|
|
$alpha = 0.2125; |
23 |
|
|
$alphaInt = 0.5125; |
24 |
|
|
$alphaSlope = 0.025; |
25 |
chrisfen |
2997 |
$invalidWater = 0; |
26 |
chrisfen |
3030 |
$waterCase = -1; |
27 |
chrisfen |
2996 |
|
28 |
|
|
# get our options |
29 |
chrisfen |
3031 |
getopts('hmrvd:l:n:o:w:'); |
30 |
chrisfen |
2996 |
|
31 |
|
|
# if we don't have a filename, drop to -h |
32 |
chrisfen |
3031 |
$opt_h = 'true' if $#ARGV != -1; |
33 |
chrisfen |
2996 |
|
34 |
|
|
# our option output |
35 |
|
|
if ($opt_h){ |
36 |
|
|
print "waterBoxer: builds water boxes\n\n"; |
37 |
chrisfen |
3031 |
print "usage: waterBoxer [-hmrv] [-d density] [-l lattice] [-n # waters]\n"; |
38 |
|
|
print "\t[-o file name] [-w water name] \n\n"; |
39 |
chrisfen |
2996 |
print " -h : show this message\n"; |
40 |
chrisfen |
2997 |
print " -m : print out a water.md file (file with all water models)\n"; |
41 |
chrisfen |
2996 |
print " -r : randomize orientations\n"; |
42 |
|
|
print " -v : verbose output\n\n"; |
43 |
|
|
print " -d real : density in g/cm^3\n"; |
44 |
|
|
print " (default: 1)\n"; |
45 |
|
|
print " -l integer : 0 - face centered cubic, 1 - simple cubic\n"; |
46 |
|
|
print " (default: 0)\n"; |
47 |
|
|
print " -n integer : # of water molecules\n"; |
48 |
|
|
print " (default: 500)\n"; |
49 |
chrisfen |
3031 |
print " -o char : output file name\n"; |
50 |
|
|
print " (default: freshWater.md)\n"; |
51 |
chrisfen |
2996 |
print " -w char : name of the water stunt double\n"; |
52 |
|
|
print " (default: SPCE)\n"; |
53 |
|
|
print "Example:\n"; |
54 |
chrisfen |
3031 |
die " waterBoxer -d 0.997 -n 864 -w SSD_RF -o ssdrfWater.md\n"; |
55 |
chrisfen |
2996 |
} |
56 |
|
|
|
57 |
|
|
# set some variables to be used in the code |
58 |
chrisfen |
3031 |
if (defined($opt_o)){ |
59 |
|
|
$fileName = $opt_o; |
60 |
chrisfen |
2996 |
} else { |
61 |
chrisfen |
3031 |
$fileName = 'freshWater.md'; |
62 |
chrisfen |
2996 |
} |
63 |
chrisfen |
2997 |
if ($opt_m){ |
64 |
|
|
die "$fileName cannot be \"water.md\"\n\tPlease choose a different name\n" if $fileName eq 'water.md'; |
65 |
|
|
$waterFileHandle = 'WATERMD'; |
66 |
|
|
} else { |
67 |
|
|
$waterFileHandle = 'OUTFILE'; |
68 |
|
|
} |
69 |
chrisfen |
2996 |
if ($opt_r){ |
70 |
|
|
$doRandomize = $opt_r; |
71 |
|
|
} |
72 |
chrisfen |
2997 |
|
73 |
chrisfen |
2996 |
if (defined($opt_w)){ |
74 |
|
|
$waterName = $opt_w; |
75 |
|
|
} else { |
76 |
|
|
$waterName = 'SPCE'; |
77 |
|
|
} |
78 |
chrisfen |
2997 |
validateWater(); |
79 |
chrisfen |
3031 |
if ($invalidWater == 1){ |
80 |
|
|
print "Warning: \'$waterName\' is not a recognized water model name.\n"; |
81 |
|
|
print " Use the \'-m\' option to generate a \'water.md\' with the\n"; |
82 |
|
|
print " recognized water model geometries.\n"; |
83 |
|
|
} |
84 |
chrisfen |
2996 |
if (defined($opt_d)){ |
85 |
|
|
if ($opt_d =~ /^[0-9]/) { |
86 |
|
|
$density = $opt_d; |
87 |
|
|
} else { |
88 |
|
|
die "\t-d value ($opt_d) is not a valid number\n\tPlease choose a positive real # value\n"; |
89 |
|
|
} |
90 |
|
|
} |
91 |
|
|
if (defined($opt_l)){ |
92 |
|
|
if ($opt_l =~ /^[0-9]/) { |
93 |
|
|
$lattice = $opt_l; |
94 |
|
|
if ($lattice != 0 && $lattice != 1){ |
95 |
|
|
die "\t-l value ($opt_l) is not a valid number\n\tPlease choose 0 or 1\n"; |
96 |
|
|
} |
97 |
|
|
} else { |
98 |
|
|
die "\t-l value ($opt_l) is not a valid number\n\tPlease choose 0 or 1\n"; |
99 |
|
|
} |
100 |
|
|
} |
101 |
|
|
if (defined($opt_n)){ |
102 |
|
|
if ($opt_n =~ /^[0-9]/) { |
103 |
|
|
$nMol = $opt_n; |
104 |
|
|
} else { |
105 |
|
|
die "\t-n value ($opt_n) is not a valid number\n\tPlease choose a non-negative integer\n"; |
106 |
|
|
} |
107 |
|
|
} |
108 |
|
|
|
109 |
|
|
# open the file writer |
110 |
|
|
open(OUTFILE, ">./$fileName") || die "\tError: can't open file $fileName\n"; |
111 |
|
|
|
112 |
|
|
# check to set magic lattice numbers |
113 |
|
|
if ($lattice == 0){ |
114 |
|
|
$crystalNumReal = ($nMol/4.0)**(1.0/3.0); |
115 |
|
|
$crystalNum = int($crystalNumReal + $tolerance); |
116 |
|
|
$remainder = $crystalNumReal - $crystalNum; |
117 |
|
|
|
118 |
|
|
# if crystalNumReal wasn't an integer, we bump the crystal to the next |
119 |
|
|
# magic number |
120 |
|
|
if ($remainder > $tolerance){ |
121 |
|
|
$crystalNum = $crystalNum + 1; |
122 |
|
|
$newMol = 4 * $crystalNum**3; |
123 |
|
|
print "WARNING: The number chosen ($nMol) failed to build a clean\n"; |
124 |
|
|
print "fcc lattice. The number of molecules has been increased to\n"; |
125 |
|
|
print "the next magic number ($newMol).\n"; |
126 |
|
|
$nMol = $newMol; |
127 |
|
|
} |
128 |
|
|
} elsif ($lattice == 1){ |
129 |
|
|
$crystalNumReal = ($nMol/1.0)**(1.0/3.0); |
130 |
|
|
$crystalNum = int($crystalNumReal); |
131 |
|
|
$remainder = $crystalNumReal - $crystalNum; |
132 |
|
|
|
133 |
|
|
# again, if crystalNumReal wasn't an integer, we bump the crystal to the next |
134 |
|
|
# magic number |
135 |
|
|
if ($remainder > $tolerance){ |
136 |
|
|
$crystalNum = $crystalNum + 1; |
137 |
|
|
$newMol = $crystalNum**3; |
138 |
|
|
print "WARNING: The number chosen ($nMol) failed to build a clean\n"; |
139 |
|
|
print "simple cubic lattice. The number of molecules has been\n"; |
140 |
|
|
print "increased to the next magic number ($newMol).\n"; |
141 |
|
|
$nMol = $newMol; |
142 |
|
|
} |
143 |
|
|
} |
144 |
|
|
|
145 |
|
|
# now we can start building the crystals |
146 |
|
|
$boxLength = ($nMol*$densityConvert/$density)**(1.0/3.0); |
147 |
|
|
$cellLength = $boxLength / $crystalNum; |
148 |
|
|
$cell2 = $cellLength*0.5; |
149 |
|
|
|
150 |
|
|
if ($lattice == 0) { |
151 |
|
|
# build the unit cell |
152 |
|
|
# molecule 0 |
153 |
|
|
$xCorr[0] = 0.0; |
154 |
|
|
$yCorr[0] = 0.0; |
155 |
|
|
$zCorr[0] = 0.0; |
156 |
|
|
# molecule 1 |
157 |
|
|
$xCorr[1] = 0.0; |
158 |
|
|
$yCorr[1] = $cell2; |
159 |
|
|
$zCorr[1] = $cell2; |
160 |
|
|
# molecule 2 |
161 |
|
|
$xCorr[2] = $cell2; |
162 |
|
|
$yCorr[2] = $cell2; |
163 |
|
|
$zCorr[2] = 0.0; |
164 |
|
|
# molecule 3 |
165 |
|
|
$xCorr[3] = $cell2; |
166 |
|
|
$yCorr[3] = 0.0; |
167 |
|
|
$zCorr[3] = $cell2; |
168 |
|
|
# assemble the lattice |
169 |
|
|
$counter = 0; |
170 |
|
|
for ($z = 0; $z < $crystalNum; $z++) { |
171 |
|
|
for ($y = 0; $y < $crystalNum; $y++) { |
172 |
|
|
for ($x = 0; $x < $crystalNum; $x++) { |
173 |
|
|
for ($uc = 0; $uc < 4; $uc++) { |
174 |
|
|
$xCorr[$uc+$counter] = $xCorr[$uc] + $cellLength*$x; |
175 |
|
|
$yCorr[$uc+$counter] = $yCorr[$uc] + $cellLength*$y; |
176 |
|
|
$zCorr[$uc+$counter] = $zCorr[$uc] + $cellLength*$z; |
177 |
|
|
} |
178 |
|
|
$counter = $counter + 4; |
179 |
|
|
} |
180 |
|
|
} |
181 |
|
|
} |
182 |
|
|
|
183 |
|
|
} elsif ($lattice == 1) { |
184 |
|
|
# build the unit cell |
185 |
|
|
# molecule 0 |
186 |
|
|
$xCorr[0] = $cell2; |
187 |
|
|
$yCorr[0] = $cell2; |
188 |
|
|
$zCorr[0] = $cell2; |
189 |
|
|
#assemble the lattice |
190 |
|
|
$counter = 0; |
191 |
|
|
for ($z = 0; $z < $crystalNum; $z++) { |
192 |
|
|
for ($y = 0; $y < $crystalNum; $y++) { |
193 |
|
|
for ($x = 0; $x < $crystalNum; $x++) { |
194 |
|
|
$xCorr[$counter] = $xCorr[0] + $cellLength*$x; |
195 |
|
|
$yCorr[$counter] = $yCorr[0] + $cellLength*$y; |
196 |
|
|
$zCorr[$counter] = $zCorr[0] + $cellLength*$z; |
197 |
|
|
|
198 |
|
|
$counter++; |
199 |
|
|
} |
200 |
|
|
} |
201 |
|
|
} |
202 |
|
|
} |
203 |
|
|
|
204 |
|
|
writeOutFile(); |
205 |
chrisfen |
3031 |
print "Water box \"$fileName\" generated.\n"; |
206 |
chrisfen |
2996 |
|
207 |
chrisfen |
3031 |
if ($opt_m){ |
208 |
chrisfen |
2997 |
printWaterMD(); |
209 |
chrisfen |
3031 |
print "File \"water.md\" generated for inclusion in \"$fileName\"\n"; |
210 |
chrisfen |
2997 |
} |
211 |
|
|
|
212 |
chrisfen |
2996 |
# this marks the end of the main program, below is subroutines |
213 |
|
|
|
214 |
|
|
sub acos { |
215 |
|
|
my ($rad) = @_; |
216 |
|
|
my $ret = atan2(sqrt(1 - $rad*$rad), $rad); |
217 |
|
|
return $ret; |
218 |
|
|
} |
219 |
|
|
|
220 |
chrisfen |
3030 |
sub writeOutFile { |
221 |
chrisfen |
2996 |
# write out the header |
222 |
|
|
print OUTFILE "<OOPSE version=4>\n"; |
223 |
|
|
findCutoff(); |
224 |
|
|
findAlpha(); |
225 |
|
|
printMetaData(); |
226 |
|
|
printFrameData(); |
227 |
|
|
print OUTFILE " <StuntDoubles>\n"; |
228 |
|
|
|
229 |
|
|
# shift the box center to the origin and write out the coordinates |
230 |
|
|
for ($i = 0; $i < $nMol; $i++) { |
231 |
|
|
$xCorr[$i] -= 0.5*$boxLength; |
232 |
|
|
$yCorr[$i] -= 0.5*$boxLength; |
233 |
|
|
$zCorr[$i] -= 0.5*$boxLength; |
234 |
|
|
|
235 |
|
|
$q0 = 1.0; |
236 |
|
|
$q1 = 0.0; |
237 |
|
|
$q2 = 0.0; |
238 |
|
|
$q3 = 0.0; |
239 |
|
|
|
240 |
|
|
if ($doRandomize == 1){ |
241 |
|
|
$cosTheta = 2.0*rand() - 1.0; |
242 |
|
|
$theta = acos($cosTheta); |
243 |
|
|
$phi = 2.0*3.14159265359*rand(); |
244 |
|
|
$psi = 2.0*3.14159265359*rand(); |
245 |
|
|
|
246 |
|
|
$q0 = cos(0.5*$theta)*cos(0.5*($phi + $psi)); |
247 |
|
|
$q1 = sin(0.5*$theta)*cos(0.5*($phi - $psi)); |
248 |
|
|
$q2 = sin(0.5*$theta)*sin(0.5*($phi - $psi)); |
249 |
|
|
$q3 = cos(0.5*$theta)*sin(0.5*($phi + $psi)); |
250 |
|
|
} |
251 |
|
|
|
252 |
|
|
print OUTFILE "$i\tpq\t$xCorr[$i] $yCorr[$i] $zCorr[$i] "; |
253 |
|
|
print OUTFILE "$q0 $q1 $q2 $q3\n"; |
254 |
|
|
} |
255 |
|
|
|
256 |
|
|
print OUTFILE " </StuntDoubles>\n </Snapshot>\n</OOPSE>\n"; |
257 |
|
|
} |
258 |
|
|
|
259 |
|
|
sub printMetaData { |
260 |
chrisfen |
2997 |
print OUTFILE " <MetaData>\n"; |
261 |
chrisfen |
2996 |
|
262 |
chrisfen |
2997 |
# print the water model or includes |
263 |
|
|
if ($opt_m){ |
264 |
|
|
print OUTFILE "#include \"water.md\""; |
265 |
|
|
} else { |
266 |
|
|
printWaterModel(); |
267 |
|
|
} |
268 |
|
|
printFakeWater() if $invalidWater == 1; |
269 |
|
|
|
270 |
|
|
# now back to the metaData output |
271 |
|
|
print OUTFILE "\n\ncomponent{ |
272 |
chrisfen |
2996 |
type = \"$waterName\"; |
273 |
|
|
nMol = $nMol; |
274 |
|
|
} |
275 |
|
|
|
276 |
|
|
ensemble = NVE; |
277 |
|
|
forceField = \"DUFF\"; |
278 |
|
|
electrostaticSummationMethod = \"shifted_force\"; |
279 |
|
|
electrostaticScreeningMethod = \"damped\"; |
280 |
|
|
cutoffRadius = $cutoff; |
281 |
|
|
|
282 |
|
|
targetTemp = 300; |
283 |
|
|
targetPressure = 1.0; |
284 |
|
|
|
285 |
|
|
tauThermostat = 1e3; |
286 |
|
|
tauBarostat = 1e4; |
287 |
|
|
|
288 |
|
|
dt = 2.0; |
289 |
|
|
runTime = 1e3; |
290 |
|
|
|
291 |
|
|
tempSet = \"true\"; |
292 |
|
|
thermalTime = 10; |
293 |
|
|
sampleTime = 100; |
294 |
|
|
statusTime = 2; |
295 |
|
|
</MetaData>\n"; |
296 |
|
|
} |
297 |
|
|
|
298 |
chrisfen |
3030 |
sub findCutoff { |
299 |
chrisfen |
2996 |
$boxLength2 = 0.5*$boxLength; |
300 |
|
|
if ($boxLength2 > $cutoff){ |
301 |
|
|
# the default is good |
302 |
|
|
} else { |
303 |
|
|
$cutoff = int($boxLength2); |
304 |
|
|
} |
305 |
|
|
} |
306 |
|
|
|
307 |
chrisfen |
3030 |
sub findAlpha { |
308 |
chrisfen |
2996 |
$alpha = $alphaInt - $cutoff*$alphaSlope; |
309 |
|
|
} |
310 |
|
|
|
311 |
chrisfen |
3030 |
sub printFrameData { |
312 |
chrisfen |
2996 |
print OUTFILE |
313 |
|
|
" <Snapshot> |
314 |
|
|
<FrameData> |
315 |
|
|
Time: 0 |
316 |
|
|
Hmat: {{ $boxLength, 0, 0 }, { 0, $boxLength, 0 }, { 0, 0, $boxLength }} |
317 |
|
|
</FrameData>\n"; |
318 |
|
|
} |
319 |
chrisfen |
2997 |
|
320 |
chrisfen |
3030 |
sub printWaterMD { |
321 |
chrisfen |
2997 |
open(WATERMD, ">./water.md") || die "\tError: can't open file water.md\n"; |
322 |
chrisfen |
3030 |
$waterFileHandle = 'WATERMD'; |
323 |
chrisfen |
2997 |
|
324 |
|
|
print WATERMD "#ifndef _WATER_MD_\n#define _WATER_MD_\n"; |
325 |
|
|
printCl(); |
326 |
|
|
printNa(); |
327 |
|
|
printSSD_E(); |
328 |
|
|
printSSD_RF(); |
329 |
|
|
printSSD(); |
330 |
|
|
printSSD1(); |
331 |
|
|
printTRED(); |
332 |
|
|
printTIP3P(); |
333 |
|
|
printTIP4P(); |
334 |
|
|
printTIP4PEw(); |
335 |
|
|
printTIP5P(); |
336 |
|
|
printTIP5PE(); |
337 |
|
|
printSPCE(); |
338 |
|
|
printSPC(); |
339 |
|
|
printDPD(); |
340 |
|
|
print WATERMD "\n\n#endif"; |
341 |
|
|
} |
342 |
|
|
|
343 |
chrisfen |
3030 |
sub printCl { |
344 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
345 |
|
|
name = \"Cl-\"; |
346 |
|
|
|
347 |
|
|
atom[0]{ |
348 |
|
|
type = \"Cl-\"; |
349 |
|
|
position(0.0, 0.0, 0.0); |
350 |
|
|
} |
351 |
|
|
}" |
352 |
|
|
} |
353 |
|
|
|
354 |
chrisfen |
3030 |
sub printNa { |
355 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
356 |
|
|
name = \"Na+\"; |
357 |
|
|
|
358 |
|
|
atom[0]{ |
359 |
|
|
type = \"Na+\"; |
360 |
|
|
position(0.0, 0.0, 0.0); |
361 |
|
|
} |
362 |
|
|
}" |
363 |
|
|
} |
364 |
|
|
|
365 |
chrisfen |
3030 |
sub printSSD_E { |
366 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
367 |
|
|
name = \"SSD_E\"; |
368 |
|
|
|
369 |
|
|
atom[0]{ |
370 |
|
|
type = \"SSD_E\"; |
371 |
|
|
position( 0.0, 0.0, 0.0 ); |
372 |
|
|
orientation( 0.0, 0.0, 0.0 ); |
373 |
|
|
} |
374 |
|
|
}" |
375 |
|
|
} |
376 |
|
|
|
377 |
chrisfen |
3030 |
sub printSSD_RF { |
378 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
379 |
|
|
name = \"SSD_RF\"; |
380 |
|
|
|
381 |
|
|
atom[0]{ |
382 |
|
|
type = \"SSD_RF\"; |
383 |
|
|
position( 0.0, 0.0, 0.0 ); |
384 |
|
|
orientation( 0.0, 0.0, 0.0 ); |
385 |
|
|
} |
386 |
|
|
}" |
387 |
|
|
} |
388 |
|
|
|
389 |
chrisfen |
3030 |
sub printSSD { |
390 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
391 |
|
|
name = \"SSD\"; |
392 |
|
|
|
393 |
|
|
atom[0]{ |
394 |
|
|
type = \"SSD\"; |
395 |
|
|
position( 0.0, 0.0, 0.0 ); |
396 |
|
|
orientation( 0.0, 0.0, 0.0 ); |
397 |
|
|
} |
398 |
|
|
}" |
399 |
|
|
} |
400 |
|
|
|
401 |
chrisfen |
3030 |
sub printSSD1 { |
402 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
403 |
|
|
name = \"SSD1\"; |
404 |
|
|
|
405 |
|
|
atom[0]{ |
406 |
|
|
type = \"SSD1\"; |
407 |
|
|
position( 0.0, 0.0, 0.0 ); |
408 |
|
|
orientation( 0.0, 0.0, 0.0 ); |
409 |
|
|
} |
410 |
|
|
}" |
411 |
|
|
} |
412 |
|
|
|
413 |
chrisfen |
3030 |
sub printTRED { |
414 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
415 |
|
|
name = \"TRED\"; |
416 |
|
|
|
417 |
|
|
atom[0]{ |
418 |
|
|
type = \"TRED\"; |
419 |
|
|
position( 0.0, 0.0, 0.0 ); |
420 |
|
|
orientation( 0.0, 0.0, 0.0 ); |
421 |
|
|
} |
422 |
|
|
atom[1]{ |
423 |
|
|
type = \"EP_TRED\"; |
424 |
|
|
position( 0.0, 0.0, 0.5 ); |
425 |
|
|
} |
426 |
|
|
|
427 |
|
|
rigidBody[0]{ |
428 |
|
|
members(0, 1); |
429 |
|
|
} |
430 |
|
|
|
431 |
|
|
cutoffGroup{ |
432 |
|
|
members(0, 1); |
433 |
|
|
} |
434 |
|
|
}" |
435 |
|
|
} |
436 |
|
|
|
437 |
chrisfen |
3030 |
sub printTIP3P { |
438 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
439 |
|
|
name = \"TIP3P\"; |
440 |
|
|
|
441 |
|
|
atom[0]{ |
442 |
|
|
type = \"O_TIP3P\"; |
443 |
|
|
position( 0.0, 0.0, -0.06556 ); |
444 |
|
|
} |
445 |
|
|
atom[1]{ |
446 |
|
|
type = \"H_TIP3P\"; |
447 |
|
|
position( 0.0, 0.75695, 0.52032 ); |
448 |
|
|
} |
449 |
|
|
atom[2]{ |
450 |
|
|
type = \"H_TIP3P\"; |
451 |
|
|
position( 0.0, -0.75695, 0.52032 ); |
452 |
|
|
} |
453 |
|
|
|
454 |
|
|
rigidBody[0]{ |
455 |
|
|
members(0, 1, 2); |
456 |
|
|
} |
457 |
|
|
|
458 |
|
|
cutoffGroup{ |
459 |
|
|
members(0, 1, 2); |
460 |
|
|
} |
461 |
|
|
}" |
462 |
|
|
} |
463 |
|
|
|
464 |
chrisfen |
3030 |
sub printTIP4P { |
465 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
466 |
|
|
name = \"TIP4P\"; |
467 |
|
|
|
468 |
|
|
atom[0]{ |
469 |
|
|
type = \"O_TIP4P\"; |
470 |
|
|
position( 0.0, 0.0, -0.06556 ); |
471 |
|
|
} |
472 |
|
|
atom[1]{ |
473 |
|
|
type = \"H_TIP4P\"; |
474 |
|
|
position( 0.0, 0.75695, 0.52032 ); |
475 |
|
|
} |
476 |
|
|
atom[2]{ |
477 |
|
|
type = \"H_TIP4P\"; |
478 |
|
|
position( 0.0, -0.75695, 0.52032 ); |
479 |
|
|
} |
480 |
|
|
atom[3]{ |
481 |
|
|
type = \"EP_TIP4P\"; |
482 |
|
|
position( 0.0, 0.0, 0.08444 ); |
483 |
|
|
} |
484 |
|
|
|
485 |
|
|
rigidBody[0]{ |
486 |
|
|
members(0, 1, 2, 3); |
487 |
|
|
} |
488 |
|
|
|
489 |
|
|
cutoffGroup{ |
490 |
|
|
members(0, 1, 2, 3); |
491 |
|
|
} |
492 |
|
|
}" |
493 |
|
|
} |
494 |
|
|
|
495 |
chrisfen |
3030 |
sub printTIP4PEw { |
496 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
497 |
|
|
name = \"TIP4P-Ew\"; |
498 |
|
|
|
499 |
|
|
atom[0]{ |
500 |
|
|
type = \"O_TIP4P-Ew\"; |
501 |
|
|
position( 0.0, 0.0, -0.06556 ); |
502 |
|
|
} |
503 |
|
|
atom[1]{ |
504 |
|
|
type = \"H_TIP4P-Ew\"; |
505 |
|
|
position( 0.0, 0.75695, 0.52032 ); |
506 |
|
|
} |
507 |
|
|
atom[2]{ |
508 |
|
|
type = \"H_TIP4P-Ew\"; |
509 |
|
|
position( 0.0, -0.75695, 0.52032 ); |
510 |
|
|
} |
511 |
|
|
atom[3]{ |
512 |
|
|
type = \"EP_TIP4P-Ew\"; |
513 |
|
|
position( 0.0, 0.0, 0.05944 ); |
514 |
|
|
} |
515 |
|
|
|
516 |
|
|
rigidBody[0]{ |
517 |
|
|
members(0, 1, 2, 3); |
518 |
|
|
} |
519 |
|
|
|
520 |
|
|
cutoffGroup{ |
521 |
|
|
members(0, 1, 2, 3); |
522 |
|
|
} |
523 |
|
|
}" |
524 |
|
|
} |
525 |
|
|
|
526 |
chrisfen |
3030 |
sub printTIP5P { |
527 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
528 |
|
|
name = \"TIP5P\"; |
529 |
|
|
|
530 |
|
|
atom[0]{ |
531 |
|
|
type = \"O_TIP5P\"; |
532 |
|
|
position( 0.0, 0.0, -0.06556 ); |
533 |
|
|
} |
534 |
|
|
atom[1]{ |
535 |
|
|
type = \"H_TIP5P\"; |
536 |
|
|
position( 0.0, 0.75695, 0.52032 ); |
537 |
|
|
} |
538 |
|
|
atom[2]{ |
539 |
|
|
type = \"H_TIP5P\"; |
540 |
|
|
position( 0.0, -0.75695, 0.52032 ); |
541 |
|
|
} |
542 |
|
|
atom[3]{ |
543 |
|
|
type = \"EP_TIP5P\"; |
544 |
|
|
position( 0.57154, 0.0, -0.46971 ); |
545 |
|
|
} |
546 |
|
|
atom[4]{ |
547 |
|
|
type = \"EP_TIP5P\"; |
548 |
|
|
position( -0.57154, 0.0, -0.46971 ); |
549 |
|
|
} |
550 |
|
|
|
551 |
|
|
rigidBody[0]{ |
552 |
|
|
members(0, 1, 2, 3, 4); |
553 |
|
|
} |
554 |
|
|
|
555 |
|
|
cutoffGroup{ |
556 |
|
|
members(0, 1, 2, 3, 4); |
557 |
|
|
} |
558 |
|
|
}" |
559 |
|
|
} |
560 |
|
|
|
561 |
chrisfen |
3030 |
sub printTIP5PE { |
562 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
563 |
|
|
name = \"TIP5P-E\"; |
564 |
|
|
|
565 |
|
|
atom[0]{ |
566 |
|
|
type = \"O_TIP5P-E\"; |
567 |
|
|
position( 0.0, 0.0, -0.06556 ); |
568 |
|
|
} |
569 |
|
|
atom[1]{ |
570 |
|
|
type = \"H_TIP5P\"; |
571 |
|
|
position( 0.0, 0.75695, 0.52032 ); |
572 |
|
|
} |
573 |
|
|
atom[2]{ |
574 |
|
|
type = \"H_TIP5P\"; |
575 |
|
|
position( 0.0, -0.75695, 0.52032 ); |
576 |
|
|
} |
577 |
|
|
atom[3]{ |
578 |
|
|
type = \"EP_TIP5P\"; |
579 |
|
|
position( 0.57154, 0.0, -0.46971 ); |
580 |
|
|
} |
581 |
|
|
atom[4]{ |
582 |
|
|
type = \"EP_TIP5P\"; |
583 |
|
|
position( -0.57154, 0.0, -0.46971 ); |
584 |
|
|
} |
585 |
|
|
|
586 |
|
|
rigidBody[0]{ |
587 |
|
|
members(0, 1, 2, 3, 4); |
588 |
|
|
} |
589 |
|
|
|
590 |
|
|
cutoffGroup{ |
591 |
|
|
members(0, 1, 2, 3, 4); |
592 |
|
|
} |
593 |
|
|
}" |
594 |
|
|
} |
595 |
|
|
|
596 |
chrisfen |
3030 |
sub printSPCE { |
597 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
598 |
|
|
name = \"SPCE\"; |
599 |
|
|
|
600 |
|
|
atom[0]{ |
601 |
|
|
type = \"O_SPCE\"; |
602 |
|
|
position( 0.0, 0.0, -0.06461 ); |
603 |
|
|
} |
604 |
|
|
atom[1]{ |
605 |
|
|
type = \"H_SPCE\"; |
606 |
|
|
position( 0.0, 0.81649, 0.51275 ); |
607 |
|
|
} |
608 |
|
|
atom[2]{ |
609 |
|
|
type = \"H_SPCE\"; |
610 |
|
|
position( 0.0, -0.81649, 0.51275 ); |
611 |
|
|
} |
612 |
|
|
|
613 |
|
|
rigidBody[0]{ |
614 |
|
|
members(0, 1, 2); |
615 |
|
|
} |
616 |
|
|
|
617 |
|
|
cutoffGroup{ |
618 |
|
|
members(0, 1, 2); |
619 |
|
|
} |
620 |
|
|
}" |
621 |
|
|
} |
622 |
|
|
|
623 |
chrisfen |
3030 |
sub printSPC { |
624 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
625 |
|
|
name = \"SPC\"; |
626 |
|
|
|
627 |
|
|
atom[0]{ |
628 |
|
|
type = \"O_SPC\"; |
629 |
|
|
position( 0.0, 0.0, -0.06461 ); |
630 |
|
|
} |
631 |
|
|
atom[1]{ |
632 |
|
|
type = \"H_SPC\"; |
633 |
|
|
position( 0.0, 0.81649, 0.51275 ); |
634 |
|
|
} |
635 |
|
|
atom[2]{ |
636 |
|
|
type = \"H_SPC\"; |
637 |
|
|
position( 0.0, -0.81649, 0.51275 ); |
638 |
|
|
} |
639 |
|
|
|
640 |
|
|
rigidBody[0]{ |
641 |
|
|
members(0, 1, 2); |
642 |
|
|
} |
643 |
|
|
|
644 |
|
|
cutoffGroup{ |
645 |
|
|
members(0, 1, 2); |
646 |
|
|
} |
647 |
|
|
}" |
648 |
|
|
} |
649 |
|
|
|
650 |
chrisfen |
3030 |
sub printDPD { |
651 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
652 |
|
|
name = \"DPD\"; |
653 |
|
|
|
654 |
|
|
atom[0]{ |
655 |
|
|
type = \"DPD\"; |
656 |
|
|
position(0.0, 0.0, 0.0); |
657 |
|
|
} |
658 |
|
|
}" |
659 |
|
|
} |
660 |
|
|
|
661 |
|
|
|
662 |
chrisfen |
3030 |
sub printFakeWater { |
663 |
chrisfen |
2997 |
print $waterFileHandle "\n\nmolecule{ |
664 |
|
|
name = \"$waterName\"; |
665 |
|
|
|
666 |
|
|
atom[0]{ |
667 |
|
|
type = \"$waterName\"; |
668 |
|
|
position(0.0, 0.0, 0.0); |
669 |
|
|
} |
670 |
|
|
}" |
671 |
|
|
} |
672 |
|
|
|
673 |
|
|
|
674 |
chrisfen |
3030 |
sub validateWater { |
675 |
chrisfen |
2997 |
if ($waterName eq 'Cl-') { $waterCase = 0; } |
676 |
|
|
elsif ($waterName eq 'Na+') { $waterCase = 1; } |
677 |
|
|
elsif ($waterName eq 'SSD_E') { $waterCase = 2; } |
678 |
|
|
elsif ($waterName eq 'SSD_RF') { $waterCase = 3; } |
679 |
|
|
elsif ($waterName eq 'SSD') { $waterCase = 4; } |
680 |
|
|
elsif ($waterName eq 'SSD1') { $waterCase = 5; } |
681 |
|
|
elsif ($waterName eq 'TIP3P') { $waterCase = 6; } |
682 |
|
|
elsif ($waterName eq 'TIP4P') { $waterCase = 7; } |
683 |
|
|
elsif ($waterName eq 'TIP4P-Ew') { $waterCase = 8; } |
684 |
|
|
elsif ($waterName eq 'TIP5P') { $waterCase = 9; } |
685 |
|
|
elsif ($waterName eq 'TIP5P-E') { $waterCase = 10; } |
686 |
|
|
elsif ($waterName eq 'SPCE') { $waterCase = 11; } |
687 |
|
|
elsif ($waterName eq 'SPC') { $waterCase = 12; } |
688 |
|
|
elsif ($waterName eq 'DPD') { $waterCase = 13; } |
689 |
|
|
else { $invalidWater = 1; } |
690 |
|
|
} |
691 |
|
|
|
692 |
chrisfen |
3030 |
sub printWaterModel { |
693 |
chrisfen |
2997 |
if ($waterCase == 0) { printCl(); } |
694 |
|
|
elsif ($waterCase == 1) { printNa(); } |
695 |
|
|
elsif ($waterCase == 2) { printSSD_E(); } |
696 |
|
|
elsif ($waterCase == 3) { printSSD_RF(); } |
697 |
|
|
elsif ($waterCase == 4) { printSSD(); } |
698 |
|
|
elsif ($waterCase == 5) { printSSD1(); } |
699 |
|
|
elsif ($waterCase == 6) { printTIP3P(); } |
700 |
|
|
elsif ($waterCase == 7) { printTIP4P(); } |
701 |
|
|
elsif ($waterCase == 8) { printTIP4PEw(); } |
702 |
|
|
elsif ($waterCase == 9) { printTIP5P(); } |
703 |
|
|
elsif ($waterCase == 10) { printTIP5PE(); } |
704 |
|
|
elsif ($waterCase == 11) { printSPCE(); } |
705 |
|
|
elsif ($waterCase == 12) { printSPC(); } |
706 |
|
|
elsif ($waterCase == 13) { printDPD(); } |
707 |
|
|
} |