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