1 |
#!@PERLINTERP@ -w |
2 |
|
3 |
# program that builds water boxes |
4 |
|
5 |
# author = "Chris Fennell |
6 |
# version = "$Revision: 1.3 $" |
7 |
# date = "$Date: 2006-10-02 23:27:40 $" |
8 |
# 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 |
$invalidWater = 0; |
26 |
$waterCase = -1; |
27 |
|
28 |
# get our options |
29 |
getopts('hmrvd:l:n:w:'); |
30 |
|
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 |
print " -m : print out a water.md file (file with all water models)\n"; |
41 |
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 |
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 |
if ($opt_r){ |
68 |
$doRandomize = $opt_r; |
69 |
} |
70 |
|
71 |
if (defined($opt_w)){ |
72 |
$waterName = $opt_w; |
73 |
} else { |
74 |
$waterName = 'SPCE'; |
75 |
} |
76 |
validateWater(); |
77 |
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 |
if ($opt_m || $invalidWater){ |
200 |
printWaterMD(); |
201 |
} |
202 |
|
203 |
# 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 |
sub writeOutFile { |
212 |
# 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 |
print OUTFILE " <MetaData>\n"; |
252 |
|
253 |
# 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 |
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 |
sub findCutoff { |
290 |
$boxLength2 = 0.5*$boxLength; |
291 |
if ($boxLength2 > $cutoff){ |
292 |
# the default is good |
293 |
} else { |
294 |
$cutoff = int($boxLength2); |
295 |
} |
296 |
} |
297 |
|
298 |
sub findAlpha { |
299 |
$alpha = $alphaInt - $cutoff*$alphaSlope; |
300 |
} |
301 |
|
302 |
sub printFrameData { |
303 |
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 |
|
311 |
sub printWaterMD { |
312 |
open(WATERMD, ">./water.md") || die "\tError: can't open file water.md\n"; |
313 |
$waterFileHandle = 'WATERMD'; |
314 |
|
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 |
sub printCl { |
335 |
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 |
sub printNa { |
346 |
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 |
sub printSSD_E { |
357 |
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 |
sub printSSD_RF { |
369 |
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 |
sub printSSD { |
381 |
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 |
sub printSSD1 { |
393 |
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 |
sub printTRED { |
405 |
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 |
sub printTIP3P { |
429 |
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 |
sub printTIP4P { |
456 |
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 |
sub printTIP4PEw { |
487 |
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 |
sub printTIP5P { |
518 |
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 |
sub printTIP5PE { |
553 |
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 |
sub printSPCE { |
588 |
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 |
sub printSPC { |
615 |
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 |
sub printDPD { |
642 |
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 |
sub printFakeWater { |
654 |
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 |
sub validateWater { |
666 |
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 |
sub printWaterModel { |
684 |
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 |
} |