source: ThirdParty/mpqc_open/lib/perl/QCParse.pm@ ef7b4b

Action_Thermostats Add_AtomRandomPerturbation Add_RotateAroundBondAction Add_SelectAtomByNameAction Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_StructOpt_integration_tests AutomationFragmentation_failures Candidate_v1.6.0 Candidate_v1.6.1 Candidate_v1.7.0 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph_documentation Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_ChronosMutex Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion GeometryObjects Gui_displays_atomic_force_velocity IndependentFragmentGrids_IntegrationTest JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool PartialCharges_OrthogonalSummation PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks RotateToPrincipalAxisSystem_UndoRedo StoppableMakroAction TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps Ubuntu_1604_changes stable
Last change on this file since ef7b4b was 860145, checked in by Frederik Heber <heber@…>, 9 years ago

Merge commit '0b990dfaa8c6007a996d030163a25f7f5fc8a7e7' as 'ThirdParty/mpqc_open'
  • Property mode set to 100755
File size: 32.5 KB
Line 
1#
2eval 'exec perl $0 $*'
3 if 0;
4
5require Molecule;
6
7##########################################################################
8
9package QCParse;
10$debug = 0;
11
12sub testparse {
13 my $parse = new QCParse;
14
15 my $string = "x:
16 xval
17test_basis: STO-3G 6-311G**
18charge: 1
19method: scf
20basis: sto-3g
21state: 3b1
22molecule:
23 H 0 0.0000001 1.00000001
24 H 0 0 -1
25gradient: yes
26optimize: no
27frequencies: yes
28properties: NPA
29y:
30yval
31z: zval1 zval2
32zval3
33h:
340 a
351
36 2 c";
37
38 print "string:\n--------------\n$string\n--------------\n";
39
40 $parse->parse_string($string);
41 $parse->doprint();
42
43 my @t = $parse->value_as_array('h');
44 print "-----------------\n";
45 for ($i = 0; $i <= $#t; $i++) {
46 print "$i: $t[$i]\n";
47 }
48 print "-----------------\n";
49
50 @t = $parse->value_as_lines('h');
51 print "-----------------\n";
52 for ($i = 0; $i <= $#t; $i++) {
53 print "$i: $t[$i]\n";
54 }
55 print "-----------------\n";
56
57 my $qcinp = new QCInput($parse);
58 my $test_basis = $parse->value("test_basis");
59 my @test_basis_a = $parse->value_as_array("test_basis");
60 my $state = $qcinp->state();
61 my $mult = $qcinp->mult();
62 my $method = $qcinp->method();
63 my $charge = $qcinp->charge();
64 my $basis = $qcinp->basis();
65 my $gradient = $qcinp->gradient();
66 my $frequencies = $qcinp->frequencies();
67 my $optimize = $qcinp->optimize();
68 my $natom = $qcinp->n_atom();
69 foreach $i (@test_basis_a) {
70 print "test_basis_a: $i\n";
71 }
72 print "test_basis = $test_basis\n";
73 print "state = $state\n";
74 print "mult = $mult\n";
75 print "method = $method\n";
76 print "basis = $basis\n";
77 print "optimize = $optimize\n";
78 print "gradient = $gradient\n";
79 print "frequencies = $frequencies\n";
80 print "natom = $natom\n";
81 for ($i = 0; $i < $natom; $i++) {
82 printf "%s %14.8f %14.8f %14.8f\n", $qcinp->element($i),
83 $qcinp->position($i,0),
84 $qcinp->position($i,1),
85 $qcinp->position($i,2);
86 }
87 printf "qcinp errors: %s\n", $qcinp->error();
88
89 my $inpwr = new MPQCInputWriter($qcinp);
90 printf "MPQC input:\n%s", $inpwr->input_string();
91}
92
93sub new {
94 my $this = shift;
95 my $class = ref($this) || $this;
96 my $self = {};
97 bless $self, $class;
98 $self->initialize();
99 return $self;
100}
101
102sub initialize {
103 my $self = shift;
104 $self->{'keyval'} = {};
105 $self->{'error'} = "";
106}
107
108sub parse_file {
109 my $self = shift;
110 my $file = shift;
111 if (! -f "$file") {
112 $self->{"ok"} = 0;
113 $self->error("File $file not found.");
114 return;
115 }
116 open(INPUT, "<$file");
117 my $string = "";
118 while (<INPUT>) {
119 $string = "$string$_";
120 }
121 close(INPUT);
122 #print "Got file:\n$string\n";
123 $self->parse_string($string);
124 $self->{"ok"} = 1;
125}
126
127sub write_file {
128 my $self = shift;
129 my $file = shift;
130 my $keyval = $self->{'keyval'};
131 my @keys = keys(%$keyval);
132 open(OUTPUT, ">$file");
133 foreach $key (@keys) {
134 my $value = $keyval->{$key};
135 print OUTPUT "${key}:\n";
136 print OUTPUT "$value\n";
137 }
138 close(OUTPUT);
139}
140
141sub parse_string {
142 my $self = shift;
143 my $string = shift;
144 my $value = "";
145 my $keyword = "";
146 $string = "$string\n";
147 while ($string) {
148 $string =~ s/^[^\n]*\n//;
149 $_ = $&;
150 s/#.*//;
151 if (/^\s*(\w+)\s*:\s*(.*)\s*$/) {
152 $self->add($keyword, $value);
153 $keyword = $1;
154 $value = $2;
155 }
156 elsif (/^\s*$/) {
157 $self->add($keyword, $value);
158 $keyword = "";
159 $value = "";
160 }
161 else {
162 $value = "$value$_";
163 }
164 }
165 $self->add($keyword, $value);
166}
167
168sub add {
169 my $self = shift;
170 my $keyword = shift;
171 my $value = shift;
172 if ($keyword ne "") {
173 $self->{'keyval'}{$keyword} = $value;
174 printf("%s = %s\n", $keyword, $value) if ($debug);
175 }
176}
177
178# returns the value of the keyword
179sub value {
180 my $self = shift;
181 my $keyword = shift;
182 my $keyval = $self->{'keyval'};
183 my $value = $keyval->{$keyword};
184 return $value;
185}
186
187# sets the value of the keyword
188sub set_value {
189 my $self = shift;
190 my $keyword = shift;
191 my $value = shift;
192 my $keyval = $self->{'keyval'};
193 $keyval->{$keyword} = $value;
194 return $value;
195}
196
197# returns the value of the keyword
198sub boolean_value {
199 my $self = shift;
200 my $keyword = shift;
201 my $keyval = $self->{'keyval'};
202 $_ = $keyval->{$keyword};
203 return "1" if (/^\s*(y|yes|1|true|t)\s*$/i);
204 return "0" if (/^\s*(n|no|0|false|f|)\s*$/i);
205 "";
206}
207
208# returns an array of whitespace delimited tokens
209sub value_as_array {
210 my $self = shift;
211 my $keyword = shift;
212 my $keyval = $self->{'keyval'};
213 my $value = $keyval->{$keyword};
214 my @array = ();
215 $i = 0;
216 $value =~ s/^\s+$//;
217 while ($value ne '') {
218 $value =~ s/^\s*(\S+)\s*//s;
219 $array[$i] = $1;
220 $i++;
221 }
222 return @array;
223}
224
225# returns an array reference of whitespace delimited tokens
226sub value_as_arrayref {
227 my $self = shift;
228 my $keyword = shift;
229 my $keyval = $self->{'keyval'};
230 my $value = $keyval->{$keyword};
231 my $array = [];
232 $i = 0;
233 $value =~ s/^\s+$//;
234 while ($value ne '') {
235 $value =~ s/^\s*(\S+)\s*//s;
236 $array->[$i] = $1;
237 $i++;
238 }
239 return $array;
240}
241
242# returns an array of lines
243sub value_as_lines {
244 my $self = shift;
245 my $keyword = shift;
246 my $keyval = $self->{'keyval'};
247 my $value = $keyval->{$keyword};
248 my @array = ();
249 $i = 0;
250 while ($value) {
251 $value =~ s/^\s*(.*)\s*\n//;
252 $array[$i] = $1;
253 $i++;
254 }
255 return @array;
256}
257
258# returns 1 if the input file existed
259sub ok {
260 my $self = shift;
261 $self->{"ok"};
262}
263
264sub display {
265 my $self = shift;
266 my @keys = @_ ? @_ : sort keys %$self;
267 foreach $key (@keys) {
268 print "\t$key => $self->{$key}\n";
269 }
270}
271
272sub doprint {
273 my $self = shift;
274 print "QCParse:\n";
275 my $keyval = $self->{'keyval'};
276 foreach $i (keys %$keyval) {
277 my $val = $keyval->{$i};
278 $val =~ s/\n/\\n/g;
279 print "keyword = $i, value = $val\n";
280 }
281}
282
283sub error {
284 my $self = shift;
285 my $msg = shift;
286 $self->{"error"} = "$self->{'error'}$msg";
287}
288
289##########################################################################
290
291package QCInput;
292$debug = 0;
293
294sub new {
295 my $this = shift;
296 my $class = ref($this) || $this;
297 my $self = {};
298 bless $self, $class;
299 $self->initialize(@_);
300 return $self;
301}
302
303sub initialize {
304 my $self = shift;
305 my $parser = shift;
306 if ($parser eq "") {
307 $parser = new QCParse;
308 }
309 $self->{"parser"} = $parser;
310 $self->{"error"} = $parser->error();
311
312 $self->{"molecule"} = new Molecule($parser->value("molecule"));
313}
314
315sub error {
316 my $self = shift;
317 my $msg = shift;
318 $self->{"error"} = "$self->{'error'}$msg";
319}
320
321sub checkpoint {
322 my $self = shift;
323 my $bval = $self->{"parser"}->boolean_value("checkpoint");
324 my $val = $self->{"parser"}->value("checkpoint");
325 if ($val ne "" && $bval eq "") {
326 $self->error("Bad value for checkpoint: $val");
327 $bval = "0";
328 }
329 elsif ($val eq "") {
330 $bval = "1";
331 }
332 $bval;
333}
334
335sub restart {
336 my $self = shift;
337 my $bval = $self->{"parser"}->boolean_value("restart");
338 my $val = $self->{"parser"}->value("restart");
339 if ($val ne "" && $bval eq "") {
340 $self->error("Bad value for restart: $val");
341 $bval = "0";
342 }
343 elsif ($val eq "") {
344 $bval = "1";
345 }
346 $bval;
347}
348
349sub label {
350 my $self = shift;
351 $self->{"parser"}->value("label");
352}
353
354sub charge {
355 my $self = shift;
356 $_ = $self->{"parser"}->value("charge");
357 s/^\s+//;
358 s/\s+$//;
359 s/^\+//;
360 if (/^$/) { $_ = "0"; }
361 if (! /^-?\d+$/) {
362 $self->error("Bad charge: $_ (using 0)\n");
363 $_ = "0";
364 }
365 $_;
366}
367
368sub method {
369 my $self = shift;
370 $_ = $self->{"parser"}->value("method");
371 s/^\s+//;
372 s/\s+$//;
373 if ($_ eq "") {
374 $self->error("No method given (using default).\n");
375 $_ = "SCF";
376 }
377 tr/a-z/A-Z/;
378 $_;
379}
380
381sub symmetry {
382 my $self = shift;
383 $_ = $self->{"parser"}->value("symmetry");
384 s/^\s*//;
385 s/\s*$//;
386 uc $_;
387}
388
389sub memory {
390 my $self = shift;
391 $_ = $self->{"parser"}->value("memory");
392 s/^\s*//;
393 s/\s*$//;
394 if ($_ eq "") {
395 $_ = 32000000;
396 }
397 $_;
398}
399
400sub state {
401 my $self = shift;
402 $_ = $self->{"parser"}->value("state");
403 s/^\s*//;
404 s/\s*$//;
405 uc $_;
406}
407
408sub mult {
409 my $self = shift;
410 $_ = $self->state();
411 s/^\s*(\d+)/\1/;
412 if (/^\s*$/) {
413 $_ = 1;
414 }
415 $_;
416}
417
418sub basis {
419 my $self = shift;
420 $_ = $self->{"parser"}->value("basis");
421 s/^\s+//;
422 s/\s+$//;
423 if ($_ eq "") {
424 $self->error("No basis given (using default).\n");
425 $_ = "STO-3G";
426 }
427 $_;
428}
429
430sub auxbasis {
431 my $self = shift;
432 $_ = $self->{"parser"}->value("auxbasis");
433 s/^\s+//;
434 s/\s+$//;
435 if ($_ eq "") {
436 $self->error("No auxiliary basis given (using default).\n");
437 $_ = "STO-3G";
438 }
439 $_;
440}
441
442sub grid {
443 my $self = shift;
444 $_ = $self->{"parser"}->value("grid");
445 s/^\s+//;
446 s/\s+$//;
447 if ($_ eq "") {
448 $_ = "default";
449 }
450 $_;
451}
452
453sub gradient {
454 my $self = shift;
455 my $bval = $self->{"parser"}->boolean_value("gradient");
456 if ($bval eq "") {
457 my $val = $self->{"parser"}->value("gradient");
458 $self->error("Bad value for gradient: $val");
459 }
460 $bval;
461}
462
463sub fzc {
464 my $self = shift;
465 $_ = $self->{"parser"}->value("fzc");
466 s/^\s+//;
467 s/\s+$//;
468 if ($_ eq "") {
469 $_ = 0;
470 }
471 $_;
472}
473
474sub fzv {
475 my $self = shift;
476 $_ = $self->{"parser"}->value("fzv");
477 s/^\s+//;
478 s/\s+$//;
479 if ($_ eq "") {
480 $_ = 0;
481 }
482 $_;
483}
484
485sub docc {
486 my $self = shift;
487 $_ = $self->{"parser"}->value("docc");
488 s/^\s+//;
489 s/\s+$//;
490 if ($_ eq "" || $_ eq "-") {
491 $_ = "auto";
492 }
493 $_;
494}
495
496sub socc {
497 my $self = shift;
498 $_ = $self->{"parser"}->value("socc");
499 s/^\s+//;
500 s/\s+$//;
501 if ($_ eq "" || $_ eq "-") {
502 $_ = "auto";
503 }
504 $_;
505}
506
507sub optimize {
508 my $self = shift;
509 my $bval = $self->{"parser"}->boolean_value("optimize");
510 if ($bval eq "") {
511 my $val = $self->{"parser"}->value("optimize");
512 $self->error("Bad value for optimize: $val");
513 }
514 $bval;
515}
516
517# returns "" if orthog_method not set
518sub orthog_method {
519 my $self = shift;
520 my $bval = $self->{"parser"}->value("orthog_method");
521 $bval;
522}
523
524# returns "" if lindep_tol not set
525sub lindep_tol {
526 my $self = shift;
527 my $bval = $self->{"parser"}->value("lindep_tol");
528 $bval;
529}
530
531sub transition_state {
532 my $self = shift;
533 my $bval = $self->{"parser"}->boolean_value("transition_state");
534 if ($bval eq "") {
535 my $val = $self->{"parser"}->value("transition_state");
536 $self->error("Bad value for transtion_state: $val");
537 }
538 $bval;
539}
540
541sub frequencies {
542 my $self = shift;
543 my $bval = $self->{"parser"}->boolean_value("frequencies");
544 if ($bval eq "") {
545 my $val = $self->{"parser"}->value("frequencies");
546 $self->error("Bad value for frequencies: $val");
547 }
548 $bval;
549}
550
551sub axyz_lines {
552 my $self = shift;
553 $self->molecule()->string();
554}
555
556sub molecule() {
557 my $self = shift;
558 return $self->{"molecule"};
559}
560
561sub n_atom {
562 my $self = shift;
563 printf "QCInput: returning natom = %d\n", $self->{"natom"} if ($debug);
564 $self->molecule()->n_atom();
565}
566
567sub element {
568 my $self = shift;
569 $self->molecule()->element(@_);
570}
571
572sub position {
573 my $self = shift;
574 $self->molecule()->position(@_);
575}
576
577sub write_file {
578 my $self = shift;
579 my $file = shift;
580 my $parser = $self->{'parser'};
581 $parser->write_file($file);
582}
583
584sub mode_following() {
585 my $self = shift;
586 return scalar($self->{"parser"}->value_as_array("followed")) != 0;
587}
588
589# returns 1 if the input file existed
590sub ok {
591 my $self = shift;
592 $self->{"parser"}->{"ok"};
593}
594
595##########################################################################
596
597package InputWriter;
598
599# Input Writer is abstract
600sub new {
601 my $this = shift;
602 my $class = ref($this) || $this;
603 my $self = {};
604 bless $self, $class;
605 $self->initialize(@_);
606 return $self;
607}
608
609sub initialize() {
610 my $self = shift;
611 my $qcinput = shift;
612 $self->{"qcinput"} = $qcinput;
613}
614
615# this should be overridden
616sub input_string() {
617 "";
618}
619
620sub write_input() {
621 my $self = shift;
622 my $file = shift;
623 my $input = $self->input_string();
624 open(OUTPUT,">$file");
625 printf OUTPUT "%s", $input;
626 close(OUTPUT);
627}
628
629sub write_qcinput {
630 my $self = shift;
631 my $file = shift;
632 my $qcinput = $self->{'qcinput'};
633 $qcinput->write_file($file);
634}
635
636##########################################################################
637
638package MPQCInputWriter;
639@ISA = qw( InputWriter );
640%methodmap = ("MP2-R12/A" => "MBPT2_R12",
641 "MP2-R12/A'" => "MBPT2_R12",
642 "MP2" => "MBPT2",
643 "OPT1[2]" => "MBPT2",
644 "OPT2[2]" => "MBPT2",
645 "ZAPT2" => "MBPT2",
646 "MP2V1" => "MBPT2",
647 "OPT1[2]V1" => "MBPT2",
648 "OPT2[2]V1" => "MBPT2",
649 "ZAPT2V1" => "MBPT2",
650 "MP2V2" => "MBPT2",
651 "OPT1[2]V2" => "MBPT2",
652 "OPT2[2]V2" => "MBPT2",
653 "ZAPT2V2" => "MBPT2",
654 "MP2V2LB" => "MBPT2",
655 "OPT1[2]V2LB" => "MBPT2",
656 "OPT2[2]V2LB" => "MBPT2",
657 "ZAPT2V2LB" => "MBPT2",
658 "ROSCF" => "SCF",
659 "SCF" => "SCF",
660 "UHF" => "UHF",
661 "CLHF" => "CLHF",
662 "HSOSHF" => "HSOSHF",
663 "HF" => "SCF",
664 "HFK" => "DFT",
665 "XALPHA" => "DFT",
666 "HFS" => "DFT",
667 "HFB" => "DFT",
668 "HFG96" => "DFT",
669 "BLYP" => "DFT",
670 "B3LYP" => "DFT",
671 "KMLYP" => "DFT",
672 "B3PW91" => "DFT",
673 "PBE" => "DFT",
674 "PW91" => "DFT",
675 "SPZ81" => "DFT",
676 "B3P86" => "DFT",
677 "BP86" => "DFT",
678 "BPW91" => "DFT",
679 "CLHFK" => "DFT",
680 "CLXALPHA" => "DFT",
681 "CLHFS" => "DFT",
682 "CLHFB" => "DFT",
683 "CLHFG96" => "DFT",
684 "CLBLYP" => "DFT",
685 "CLB3LYP" => "DFT",
686 "CLKMLYP" => "DFT",
687 "CLB3PW91" => "DFT",
688 "CLPBE" => "DFT",
689 "CLPW91" => "DFT",
690 "SPZ81" => "DFT",
691 "B3P86" => "DFT",
692 "BP86" => "DFT",
693 "BPW91" => "DFT",
694 "HSOSHFK" => "DFT",
695 "HSOSXALPHA" => "DFT",
696 "HSOSHFS" => "DFT",
697 "HSOSHFB" => "DFT",
698 "HSOSHFG96" => "DFT",
699 "HSOSBLYP" => "DFT",
700 "HSOSB3LYP" => "DFT",
701 "HSOSKMLYP" => "DFT",
702 "HSOSB3PW91" => "DFT",
703 "HSOSPBE" => "DFT",
704 "HSOSPW91" => "DFT",
705 "HSOSSPZ81" => "DFT",
706 "HSOSB3P86" => "DFT",
707 "HSOSBP86" => "DFT",
708 "HSOSBPW91" => "DFT",
709 "UHFK" => "DFT",
710 "UXALPHA" => "DFT",
711 "UHFS" => "DFT",
712 "UHFB" => "DFT",
713 "UHFG96" => "DFT",
714 "UBLYP" => "DFT",
715 "UB3LYP" => "DFT",
716 "UKMLYP" => "DFT",
717 "UB3PW91" => "DFT",
718 "UPBE" => "DFT",
719 "UPW91" => "DFT",
720 "USPZ81" => "DFT",
721 "UB3P86" => "DFT",
722 "UBP86" => "DFT",
723 "UBPW91" => "DFT",
724 );
725%mbpt2r12stdapproxmap = ("MP2-R12/A" => "A",
726 "MP2-R12/A'" => "A'",
727 );
728%mbpt2map = ("MP2" => "mp",
729 "OPT1[2]" => "opt1",
730 "OPT2[2]" => "opt2",
731 "ZAPT2" => "zapt",
732 "MP2V1" => "mp",
733 "OPT1[2]V1" => "opt1",
734 "OPT2[2]V1" => "opt2",
735 "ZAPT2V1" => "zapt",
736 "MP2V2" => "mp",
737 "OPT1[2]V2" => "opt1",
738 "OPT2[2]V2" => "opt2",
739 "ZAPT2V2" => "zapt",
740 "MP2V2LB" => "mp",
741 "OPT1[2]V2LB" => "opt1",
742 "OPT2[2]V2LB" => "opt2",
743 "ZAPT2V2LB" => "zapt");
744%mbpt2algmap = ("MP2" => "",
745 "OPT1[2]" => "",
746 "OPT2[2]" => "",
747 "ZAPT2" => "",
748 "MP2V1" => "v1",
749 "OPT1[2]V1" => "v1",
750 "OPT2[2]V1" => "v1",
751 "ZAPT2V1" => "v1",
752 "MP2V2" => "v2",
753 "OPT1[2]V2" => "v2",
754 "OPT2[2]V2" => "v2",
755 "ZAPT2V2" => "v2",
756 "MP2V2LB" => "v2lb",
757 "OPT1[2]V2LB" => "v2lb",
758 "OPT2[2]V2LB" => "v2lb",
759 "ZAPT2V2LB" => "v2lb");
760$debug = 0;
761
762sub new {
763 my $this = shift;
764 my $class = ref($this) || $this;
765 my $self = {};
766 bless $self, $class;
767
768 $self->initialize(@_);
769 return $self;
770}
771
772sub initialize() {
773 my $self = shift;
774 my $qcinput = shift;
775 $self->{"qcinput"} = $qcinput;
776}
777
778sub docc_string() {
779 my $self = shift;
780 my $qcinput = $self->{"qcinput"};
781 my $occs = $qcinput->docc();
782 if ($occs eq "auto") { return ""; }
783 $occs =~ s/,/ /g;
784 "docc = [ $occs ]";
785}
786
787sub socc_string() {
788 my $self = shift;
789 my $qcinput = $self->{"qcinput"};
790 my $occs = $qcinput->socc();
791 if ($occs eq "auto") { return ""; }
792 $occs =~ s/,/ /g;
793 "socc = [ $occs ]";
794}
795
796sub input_string() {
797 my $self = shift;
798 my $qcinput = $self->{"qcinput"};
799 my $qcparse = $qcinput->{"parser"};
800
801 my $use_cints = 0;
802 my $do_cca = $qcparse->value("do_cca");
803
804 printf "molecule = %s\n", $qcparse->value("molecule") if ($debug);
805
806 my $symmetry = $qcinput->symmetry();
807 my $mol = "% molecule specification";
808 $mol = "$mol\nmolecule<Molecule>: (";
809 $symmetry = lc $symmetry if ($symmetry eq "AUTO");
810 if ($qcinput->frequencies()) {
811 $mol = "$mol\n symmetry = C1";
812 }
813 else {
814 $mol = "$mol\n symmetry = $symmetry";
815 }
816 $mol = "$mol\n unit = angstrom";
817 $mol = "$mol\n { atoms geometry } = {";
818 printf "MPQCInputWriter: natom = %d\n", $qcinput->n_atom() if ($debug);
819 my $i;
820 for ($i = 0; $i < $qcinput->n_atom(); $i++) {
821 $mol = sprintf "%s\n %2s [ %18.12f %18.12f %18.12f ]",
822 $mol, $qcinput->element($i),
823 $qcinput->position($i,0),
824 $qcinput->position($i,1),
825 $qcinput->position($i,2);
826 }
827 $mol = "$mol\n }";
828 $mol = "$mol\n)\n";
829
830 my $basis = "% basis set specification";
831 $basis = "$basis\nbasis<GaussianBasisSet>: (";
832 $basis = sprintf "%s\n name = \"%s\"", $basis, $qcinput->basis();
833 $basis = "$basis\n molecule = \$:molecule";
834 $basis = "$basis\n)\n";
835
836 my $integrals = "";
837 if($do_cca) {
838 $integrals = "% using cca integrals";
839 $integrals = "$integrals\nintegrals<IntegralCCA>: (";
840 my $buffer_type = $qcparse->value("integral_buffer");
841 if( $buffer_type ne "opaque" && $buffer_type ne "array" ) {
842 $buffer_type = "opaque";
843 }
844 my $int_package = $qcparse->value("integral_package");
845 if( $int_package ne "intv3" && $int_package ne "cints" ) {
846 $int_package = "intv3";
847 }
848 $integrals = "$integrals\n integral_buffer = $buffer_type";
849 $integrals = "$integrals\n integral_package = $int_package";
850 $integrals = "$integrals\n evaluator_factory = MPQC.IntegralEvaluatorFactory";
851 $integrals = "$integrals\n molecule = \$:molecule";
852 $integrals = "$integrals\n)\n";
853 }
854
855 my $fixed = $qcparse->value_as_arrayref("fixed");
856 my $followed = $qcparse->value_as_arrayref("followed");
857 if (scalar(@{$fixed}) != 0) {
858 $fixed = $self->mpqc_fixed_coor($fixed);
859 }
860 else {
861 $fixed = "";
862 }
863 if (scalar(@{$followed}) != 0) {
864 $followed = $self->mpqc_followed_coor($followed);
865 }
866 else {
867 $followed = "";
868 }
869
870 my $coor = " % molecular coordinates for optimization";
871 $coor = "$coor\n coor<SymmMolecularCoor>: (";
872 $coor = "$coor\n molecule = \$:molecule";
873 $coor = "$coor\n generator<IntCoorGen>: (";
874 $coor = "$coor\n molecule = \$:molecule";
875 $coor = "$coor\n )";
876 $coor = "$coor$followed";
877 $coor = "$coor$fixed";
878 $coor = "$coor\n )\n";
879
880 my $charge = $qcinput->charge();
881 my $mult = $qcinput->mult();
882 my $docc = $self->docc_string();
883 my $socc = $self->socc_string();
884
885 my $grid = $qcinput->grid();
886
887 my $memory = $qcinput->memory();
888 my $inputmethod = $methodmap{uc($qcinput->method())};
889 my $method = "$inputmethod";
890 $method = "SCF" if ($method eq "");
891 my $openmethod = substr(uc($qcinput->method()),0,4);
892 if (substr($openmethod,0,2) eq "CL") { $openmethod = "CL"; }
893 if (substr($openmethod,0,1) eq "U") { $openmethod = "U"; }
894 if ($method eq "SCF") {
895 if ($openmethod eq "U") {
896 $method = "UHF";
897 }
898 elsif ($openmethod eq "CL") {
899 $method = "CLHF";
900 }
901 elsif ($openmethod eq "HSOS") {
902 $method = "HSOSHF";
903 }
904 elsif ($qcinput->mult() == 1) {
905 $method = "CLHF";
906 $openmethod = "CL";
907 }
908 else {
909 $method = "HSOSHF";
910 $openmethod = "HSOS";
911 }
912 }
913 my $functional;
914 if ($method eq "DFT") {
915 $functional = uc($qcinput->method());
916 if ($openmethod eq "U") {
917 $method = "UKS";
918 $functional = substr($functional,1);
919 }
920 elsif ($openmethod eq "CL") {
921 $method = "CLKS";
922 $functional = substr($functional,2);
923 }
924 elsif ($openmethod eq "HSOS") {
925 $method = "HSOSKS";
926 $functional = substr($functional,4);
927 }
928 elsif ($qcinput->mult() == 1) {
929 $method = "CLKS";
930 $openmethod = "CL";
931 }
932 else {
933 $method = "UKS";
934 $openmethod = "U";
935 }
936 }
937 my $orthog_method = $qcinput->orthog_method();
938 my $lindep_tol = $qcinput->lindep_tol();
939 my $mole = " do_energy = yes";
940 if ($qcinput->gradient()) {
941 $mole = "$mole\n do_gradient = yes";
942 }
943 else {
944 $mole = "$mole\n do_gradient = no";
945 }
946 if($do_cca) {
947 $mole = "$mole\n do_cca = yes";
948 }
949 $mole = "$mole\n % method for computing the molecule's energy";
950 $mole = "$mole\n mole<$method>: (";
951 $mole = "$mole\n molecule = \$:molecule";
952 $mole = "$mole\n basis = \$:basis";
953 $mole = "$mole\n coor = \$..:coor";
954 $mole = "$mole\n memory = $memory";
955 if($do_cca) {
956 $mole = "$mole\n integrals = \$:integrals";
957 }
958 if ($inputmethod eq "SCF" || $inputmethod eq "UHF"
959 || $method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS") {
960 $mole = "$mole\n total_charge = $charge";
961 $mole = "$mole\n multiplicity = $mult";
962 $mole = "$mole\n print_npa = yes";
963 if ($docc ne "") {$mole = "$mole\n $docc";}
964 if ($socc ne "") {$mole = "$mole\n $socc";}
965 if ($orthog_method ne "" ) {
966 $mole = "$mole\n orthog_method = $orthog_method";
967 }
968 if ($lindep_tol ne "" ) {
969 $mole = "$mole\n lindep_tol = $lindep_tol";
970 }
971 }
972 if ($method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS") {
973 $mole = "$mole\n functional<StdDenFunctional>: name = \"$functional\"";
974 }
975 if (($method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS")
976 && $grid ne "default") {
977 $mole = "$mole\n integrator<RadialAngularIntegrator>: (grid = $grid)";
978 }
979 if ($method eq "MBPT2_R12") {
980 my $stdapprox = $mbpt2r12stdapproxmap{uc($qcinput->method())};
981 my $auxbasis = $qcinput->auxbasis();
982 my $fzc = $qcinput->fzc();
983
984 $mole = sprintf "%s\n stdapprox = \"%s\"", $mole, $stdapprox;
985 $mole = "$mole\n integrals<IntegralCints>: ()";
986 $mole = "$mole\n nfzc = $fzc";
987 # don't write an auxbasis if the auxbasis is the same as the basis set.
988 # this will speed up the calculation
989 if ("$auxbasis" ne "" && "$auxbasis" ne $qcinput->basis()) {
990 $mole = "$mole\n aux_basis<GaussianBasisSet>: (";
991 $mole = sprintf "%s\n name = \"%s\"", $mole, $auxbasis;
992 $mole = "$mole\n molecule = \$:molecule";
993 $mole = "$mole\n )\n";
994 }
995 $mole = append_reference($mole,"CLHF",$charge,$mult,$memory,$orthog_method,
996 $lindep_tol,$docc,$socc,"DZ (Dunning)");
997 $use_cints = 1;
998 }
999 elsif ($method eq "MBPT2") {
1000 my $fzc = $qcinput->fzc();
1001 my $fzv = $qcinput->fzv();
1002 my $mbpt2method = $mbpt2map{uc($qcinput->method())};
1003 my $mbpt2algorithm = $mbpt2algmap{uc($qcinput->method())};
1004 $mole = "$mole\n method = $mbpt2method";
1005 if ($mbpt2algorithm ne "") {
1006 $mole = "$mole\n algorithm = $mbpt2algorithm";
1007 }
1008 $mole = "$mole\n nfzc = $fzc";
1009 $mole = "$mole\n nfzv = $fzv";
1010 my $refmethod = "";
1011 if ($qcinput->mult() == 1) {
1012 $refmethod = "CLHF";
1013 }
1014 else {
1015 $refmethod = "HSOSHF";
1016 }
1017 $mole = append_reference($mole,$refmethod,$charge,$mult,$memory,$orthog_method,
1018 $lindep_tol,$docc,$socc,"STO-3G");
1019 }
1020 elsif (! ($basis =~ /^STO/
1021 || $basis =~ /^MI/
1022 || $basis =~ /^\d-\d1G$/) && ! $do_cca ) {
1023 my $guessmethod = "${openmethod}HF";
1024 $mole = "$mole\n guess_wavefunction<$guessmethod>: (";
1025 $mole = "$mole\n molecule = \$:molecule";
1026 $mole = "$mole\n total_charge = $charge";
1027 $mole = "$mole\n multiplicity = $mult";
1028 if ($docc ne "") {$mole = "$mole\n $docc";}
1029 if ($socc ne "") {$mole = "$mole\n $socc";}
1030 $mole = "$mole\n basis<GaussianBasisSet>: (";
1031 $mole = "$mole\n molecule = \$:molecule";
1032 $mole = "$mole\n name = \"STO-3G\"";
1033 $mole = "$mole\n )";
1034 $mole = "$mole\n memory = $memory";
1035 if($do_cca) {
1036 $mole = "$mole\n integrals = \$:integrals";
1037 }
1038 $mole = "$mole\n )";
1039 }
1040 if ($qcinput->frequencies()) {
1041 $mole = "$mole\n hessian<FinDispMolecularHessian>: (";
1042 if ($symmetry ne "C1") {
1043 $mole="$mole\n point_group<PointGroup>: symmetry = $symmetry";
1044 }
1045 $mole = "$mole\n checkpoint = no";
1046 $mole = "$mole\n restart = no";
1047 $mole = "$mole\n )";
1048 }
1049 $mole = "$mole\n )\n";
1050
1051 my $opt;
1052 if ($qcinput->optimize()) {
1053 $opt = " optimize = yes";
1054 }
1055 else {
1056 $opt = " optimize = no";
1057 }
1058 my $optclass, $updateclass;
1059 if ($qcinput->transition_state()) {
1060 $optclass = "EFCOpt";
1061 $updateclass = "PowellUpdate";
1062 }
1063 else {
1064 $optclass = "QNewtonOpt";
1065 $updateclass = "BFGSUpdate";
1066 }
1067 $opt = "$opt\n % optimizer object for the molecular geometry";
1068 $opt = "$opt\n opt<$optclass>: (";
1069 $opt = "$opt\n max_iterations = 20";
1070 $opt = "$opt\n function = \$..:mole";
1071 if ($qcinput->transition_state()) {
1072 $opt = "$opt\n transition_state = yes";
1073 if ($qcinput->mode_following()) {
1074 $opt = "$opt\n hessian = [ [ -0.1 ] ]";
1075 $opt = "$opt\n mode_following = yes";
1076 }
1077 }
1078 $opt = "$opt\n update<$updateclass>: ()";
1079 $opt = "$opt\n convergence<MolEnergyConvergence>: (";
1080 $opt = "$opt\n cartesian = yes";
1081 $opt = "$opt\n energy = \$..:..:mole";
1082 $opt = "$opt\n )";
1083 $opt = "$opt\n )\n";
1084
1085 my $freq = "";
1086 if ($qcinput->frequencies()) {
1087 $freq = "% vibrational frequency input";
1088 $freq = "$freq\n freq<MolecularFrequencies>: (";
1089 if ($symmetry ne "C1") {
1090 $freq = "$freq\n point_group<PointGroup>: symmetry = $symmetry";
1091 }
1092 $freq = "$freq\n molecule = \$:molecule";
1093 $freq = "$freq\n )\n";
1094 }
1095
1096 my $mpqcstart = sprintf ("mpqc: (\n checkpoint = %s\n",
1097 bool_to_yesno($qcinput->checkpoint()));
1098 $mpqcstart = sprintf ("%s savestate = %s\n",
1099 $mpqcstart,bool_to_yesno($qcinput->checkpoint()));
1100 $mpqcstart = sprintf ("%s restart = %s\n",
1101 $mpqcstart,bool_to_yesno($qcinput->restart()));
1102 if ($use_cints) {
1103 $mpqcstart = "$mpqcstart integrals<IntegralCints>: ()\n";
1104 }
1105 my $mpqcstop = ")\n";
1106 my $emacs = "% Emacs should use -*- KeyVal -*- mode\n";
1107 my $warn = "% this file was automatically generated\n";
1108 my $lab = $qcinput->label();
1109 my $label = "";
1110 if (! $lab =~ /^\s*$/) {
1111 $label = "% label: $lab";
1112 $label =~ s/\n/\n% label: /g;
1113 $label = "$label\n";
1114 }
1115 "$emacs$warn$label$mol$basis$integrals$mpqcstart$coor$mole$opt$freq$mpqcstop";
1116}
1117
1118sub mpqc_fixed_coor {
1119 my $self = shift;
1120 my $coorref = shift;
1121 my $result = "";
1122 $result = "\n fixed<SetIntCoor>: [";
1123 while (scalar(@{$coorref}) != 0) {
1124 my $nextcoor = $self->mpqc_sum_coor(" ","",$coorref);
1125 $result = "$result\n$nextcoor";
1126 }
1127 $result = "$result\n ]";
1128}
1129
1130sub mpqc_followed_coor {
1131 my $self = shift;
1132 my $coorref = shift;
1133 sprintf "\n%s", $self->mpqc_sum_coor(" ","followed",$coorref);
1134}
1135
1136sub mpqc_sum_coor {
1137 my $self = shift;
1138 my $space = shift;
1139 my $name = shift;
1140 my $coor = shift;
1141 my $result = "$space$name<SumIntCoor>:(";
1142 $result = "$result\n$space coor: [";
1143 my @coef = ();
1144 do {
1145 $coef[$ncoor] = shift @{$coor};
1146 my $simple = $self->mpqc_coor($coor);
1147 $result = "$result\n$space $simple";
1148 $ncoor = $ncoor + 1;
1149 } while($coor->[0] eq "+" && shift @{$coor} eq "+");
1150 $result = "$result\n$space ]";
1151 $result = "$result\n$space coef = [";
1152 my $i;
1153 foreach $i (0..$#coef) {
1154 $result = "$result $coef[$i]";
1155 }
1156 $result = "$result]";
1157 $result = "$result\n$space)";
1158 $result;
1159}
1160
1161sub mpqc_coor {
1162 my $self = shift;
1163 my $coor = shift;
1164 my $type = shift @{$coor};
1165 if ($type eq "TORS") {
1166 return sprintf "<TorsSimpleCo>:(atoms = [%d %d %d %d])",
1167 shift @{$coor},shift @{$coor},
1168 shift @{$coor},shift @{$coor};
1169 }
1170 if ($type eq "BEND") {
1171 return sprintf "<BendSimpleCo>:(atoms = [%d %d %d])",
1172 shift @{$coor},shift @{$coor},
1173 shift @{$coor};
1174 }
1175 if ($type eq "STRE") {
1176 return sprintf "<StreSimpleCo>:(atoms = [%d %d])",
1177 shift @{$coor},shift @{$coor};
1178 }
1179}
1180
1181sub bool_to_yesno {
1182 if (shift) { return "yes"; }
1183 else { return "no"; }
1184}
1185
1186sub append_reference {
1187 my $mole = shift;
1188 my $refmethod = shift;
1189 my $charge = shift;
1190 my $mult = shift;
1191 my $memory = shift;
1192 my $orthog_method = shift;
1193 my $lindep_tol = shift;
1194 my $docc = shift;
1195 my $socc = shift;
1196 my $guessbasis = shift;
1197 $mole = "$mole\n reference<$refmethod>: (";
1198 $mole = "$mole\n molecule = \$:molecule";
1199 $mole = "$mole\n basis = \$:basis";
1200 $mole = "$mole\n total_charge = $charge";
1201 $mole = "$mole\n multiplicity = $mult";
1202 $mole = "$mole\n memory = $memory";
1203 if ($orthog_method ne "" ) {
1204 $mole = "$mole\n orthog_method = $orthog_method";
1205 }
1206 if ($lindep_tol ne "" ) {
1207 $mole = "$mole\n lindep_tol = $lindep_tol";
1208 }
1209 if ($docc ne "") {$mole = "$mole\n $docc";}
1210 if ($socc ne "") {$mole = "$mole\n $socc";}
1211 if (! ($basis =~ /^STO/
1212 || $basis =~ /^MI/
1213 || $basis =~ /^\d-\d1G$/)) {
1214 $mole = "$mole\n guess_wavefunction<$refmethod>: (";
1215 $mole = "$mole\n molecule = \$:molecule";
1216 $mole = "$mole\n total_charge = $charge";
1217 $mole = "$mole\n multiplicity = $mult";
1218 if ($docc ne "") {$mole = "$mole\n $docc";}
1219 if ($socc ne "") {$mole = "$mole\n $socc";}
1220 $mole = "$mole\n basis<GaussianBasisSet>: (";
1221 $mole = "$mole\n molecule = \$:molecule";
1222 $mole = "$mole\n name = \"$guessbasis\"";
1223 $mole = "$mole\n )";
1224 $mole = "$mole\n memory = $memory";
1225 $mole = "$mole\n )";
1226 }
1227 $mole = "$mole\n )";
1228 return $mole;
1229}
1230
12311;
1232
Note: See TracBrowser for help on using the repository browser.