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