5 # Created by Toshi Nagata on 2009/11/22.
6 # Copyright 2009 Toshi Nagata. All rights reserved.
8 # This program is free software; you can redistribute it and/or modify
9 # it under the terms of the GNU General Public License as published by
10 # the Free Software Foundation version 2 of the License.
12 # This program is distributed in the hope that it will be useful,
13 # but WITHOUT ANY WARRANTY; without even the implied warranty of
14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 # GNU General Public License for more details.
34 # Import nbo log. The argument is either an array of Strings
36 def import_nbo_log(lines)
38 getline = lambda { lines.gets }
40 getline = lambda { lines.shift }
43 while (ln = getline.call) != nil
44 if ln =~ /done with NBO analysis/
47 if ln =~ /NATURAL POPULATIONS: Natural atomic orbital occupancies/
53 while (ln = getline.call) != nil
58 # Double blank lines indicates the end of the table
63 next unless ln =~ /^ *(\d+) *([A-Za-z]+) *(\d+) *([a-z]+) +([A-Za-z]+)\( *(\d+)([a-z]+)\)/
66 atom_label = $2 + $3.to_s
70 vals = $~.post_match.split.map { |e| e.to_f }
71 # atom_index, type, pn+label, occupancy[, energy]
72 nbo["nao"].push([atom_label, type, pn + label] + vals)
74 elsif ln =~ / NATURAL BOND ORBITALS \(Summary\):/
77 while (ln = getline.call) != nil
78 break if ln =~ /Total Lewis/
83 # Double blank lines indicates the end of the table
87 if ln =~ /^\s*(\d+)\. *([A-Za-z]+\*?) *\( *(\d+)\) *([A-Za-z]+) *(\d+)(- *([A-Za-z]+) *(\d+))? *(\d\.\d+)/
91 atom_label = $4 + ($5.to_i - 1).to_s
93 atom_label += "-" + $7 + ($8.to_i - 1).to_s
96 nbo["nbo"].push([orb_kind + "_" + orb_idx, atom_label, occ])
99 elsif ln =~ /([A-Z]+)s in the ([A-Z]+) basis:/ || ln =~ /([A-Z]+) Fock matrix:/
106 key = src + "/" + dst
114 while (ln = getline.call) != nil
116 # Blank line: end of one block
119 # Double blank lines indicates the end of the table
131 lab = $~.pre_match.strip
132 a = ([$~[0]] + $~.post_match.split).map { |e| e.to_f }
134 lab =~ /^[0-9]+\. *(.*)$/
137 (elements[idx] ||= []).concat(a)
140 nbo[key] = LAMatrix.new(elements).transpose!
141 key2 = (src == "F" ? dst : src) + "_L"
149 # Convert NAO based matrix into AO based matrix
150 @nbo.keys.each { |key|
151 if key[0..3] == "NAO/"
152 key2 = "AO/" + key[4..-1]
153 @nbo[key2] = @nbo["AO/NAO"] * @nbo[key]
160 def Molecule.read_gamess_basis_sets(fname)
161 # $gamess_basis = Hash.new unless $gamess_basis
162 $gamess_ecp = Hash.new unless $gamess_ecp
163 # $gamess_basis_desc = Hash.new unless $gamess_basis_desc
164 # $gamess_basis_keys = [] unless $gamess_basis_keys
165 basename = File.basename(fname, ".*")
168 File.open(fname, "r") { |fp|
170 if s =~ /^\s*!\s*/ && descname == nil
171 # Get the descriptive name from the first comment line
172 s = Regexp.last_match.post_match
174 descname = Regexp.last_match.pre_match
176 descname = (s.split)[0]
178 $gamess_basis_desc[basename] = descname
181 ss, bas = (s.split)[0..1] # Tokens delimited by whitespaces
182 next if ss == nil || ss == ""
184 # Read the ECP section
185 ecp = $gamess_ecp[basename]
188 $gamess_ecp[basename] = ecp
190 keys << basename unless keys.include?(basename)
193 ary = s.split # (PNAME, PTYPE, IZCORE, LMAX+1)
194 elem = ary[0].match(/\w{1,2}/).to_a[0] # ary[0] should begin with an element name
195 ap = Parameter.builtin.elements.find { |p| p.name.casecmp(elem) == 0 }
196 raise MolbyError, "the ECP definition does not begin with an element name: #{s}" if !ap
199 (0..Integer(ary[3])).each {
201 raise MolbyError, "the ECP definition ends unexpectedly at line #{ln} for: #{s}" if !s
203 ary = s.split # (NGPOT, comment)
205 (1..Integer(ary[0])).each {
207 raise MolbyError, "the ECP definition ends unexpectedly at line #{ln} for: #{s}" if !s
212 ecp[ap.index] = ecpdef
215 # Comments or other unrecognizable lines
217 elsif (ap = Parameter.builtin.elements.find { |p| p.name.casecmp(ss) == 0 || p.fullname.casecmp(ss) == 0 })
218 # Valid basis definition
219 if bas == nil || bas =~ /\W/
222 basis = $gamess_basis[bas]
225 $gamess_basis[bas] = basis
227 keys << bas unless keys.include?(bas)
229 while (s = fp.gets) && s =~ /\S/
232 basis[ap.index] = basdef
234 raise MolbyError, "ss is not a valid element symbol or name: #{s}"
238 unless $gamess_basis_keys.include?(basename)
239 $gamess_basis_keys.push(basename)
243 # Execute GAMESS (inferior copy of rungms script)
244 # inpname is the input file
245 # mol (optional) is the molecule from which the GAMESS input was built.
246 # If mol is specified and RUNTYP=OPTIMIZE, then the intermediate structures are
247 # displayed real-time.
248 def Molecule.execute_gamess(inpname, mol = nil)
249 gmsname = get_global_settings("gamess.executable_path")
253 cygwin_version = false # For windows: old cygwin version
256 if gmsname == nil || !File.exist?(gmsname)
257 gmsname = Dialog.open_panel("Please locate the GAMESS executable")
258 exit if gmsname == nil
260 gmsbase = File.basename(gmsname)
261 gmsdir = File.dirname(gmsname)
262 if gmsbase =~ /gamess\.(.*)\.(exe|x)$/i
267 error_message_box(gmsbase + " does not look like a GAMESS executable!")
272 # mol = Molecule.open(inpname)
274 # error_message_box("Cannot open #{inpname} as GAMESS input")
279 inpbase = File.basename(inpname)
280 inpdir = File.dirname(inpname)
281 inpbody = inpbase.sub(/\.inp$/, "")
282 logbase = inpbody + ".log"
284 set_global_settings("gamess.executable_path", gmsname)
286 ncpus = get_global_settings("gamess.ncpus").to_i
291 # Prepare the scratch directory
293 scrdir = create_temp_dir("gamess", inpbody)
295 error_message_box($!.to_s)
298 # scrdir = $home_directory + "/Molby/gamess"
300 # mkdir_recursive(scrdir)
302 # error_message_box("Cannot create directory #{scrdir}: " + $!.to_s)
306 # scrdir = scrdir + "/" + inpbody + "." + $$.to_s + ".0"
308 # while File.exist?(scrdir)
309 # scrdir = scrdir.sub(".#{n}", ".#{n + 1}")
314 if $platform == "win"
316 scrdir.gsub!("/", sep)
317 gmsdir.gsub!("/", sep)
322 # Old (obsolete) cygwin version, using ddikick
323 if $platform == "win"
325 cygwin_version = true
329 # Get the host name etc.
330 hostname = backquote("hostname").chomp
331 if $platform == "win"
332 s = backquote("cmd.exe /c dir \"#{scrdir}\"")
333 freebytes = s.split("\n").pop.match(/([0-9,]+)[^0-9]*$/).to_a[1]
335 freebytes = (freebytes.gsub(",","").to_i / 1024).to_s + " Kbytes"
337 freebytes = "(unknown)"
339 uname = backquote("cmd.exe /c ver").to_s.gsub("\n", "")
341 freebytes = `df -k "#{scrdir}"`
342 uname = `uname`.chomp
345 # Redirect standard output to the log file
346 logname = scrdir + sep + logbase
347 fpout = File.open(logname, "w")
349 # The cygwin version uses LF as the eol character
352 fpout.print "----- GAMESS execution script -----\n"
353 fpout.print "This job is running on host #{hostname}\n"
354 fpout.print "under operating system #{uname} at #{Time.now.to_s}\n"
355 fpout.print "Available scratch disk space (Kbyte units) at beginning of the job is\n"
356 fpout.print freebytes + "\n"
358 # Copy the input file
359 scrprefix = scrdir + sep + inpbody
360 if $platform == "win"
365 filecopy(inpname, scrprefix + ".F05")
366 File.open("#{scrdir}/.in_use", "w") { |fp| }
368 # Prepare environmental variables
369 auxdir = "#{gmsdir}#{sep}auxdata"
370 ENV["ERICFMT"] = "#{auxdir}#{sep}ericfmt.dat"
371 ENV["MCPPATH"] = "#{auxdir}#{sep}MCP"
372 ENV["BASPATH"] = "#{auxdir}#{sep}BASES"
373 ENV["QUANPOL"] = "#{auxdir}#{sep}QUANPOL"
374 ENV["EXTBAS"] = "/dev/null"
375 ENV["IRCDATA"] = "#{scrbody}.irc"
376 ENV["RESTART"] = "#{scrbody}.rst"
377 ENV["TRAJECT"] = "#{scrbody}.trj"
378 ENV["PUNCH"] = "#{scrbody}.dat"
379 ENV["INPUT"] = "#{scrbody}.F05"
380 ENV["AOINTS"] = "#{scrbody}.F08"
381 ENV["MOINTS"] = "#{scrbody}.F09"
382 ENV["DICTNRY"] = "#{scrbody}.F10"
383 ENV["DRTFILE"] = "#{scrbody}.F11"
384 ENV["CIVECTR"] = "#{scrbody}.F12"
385 ENV["CASINTS"] = "#{scrbody}.F13"
386 ENV["CIINTS"] = "#{scrbody}.F14"
387 ENV["WORK15"] = "#{scrbody}.F15"
388 ENV["WORK16"] = "#{scrbody}.F16"
389 ENV["CSFSAVE"] = "#{scrbody}.F17"
390 ENV["FOCKDER"] = "#{scrbody}.F18"
391 ENV["WORK19"] = "#{scrbody}.F19"
392 ENV["DASORT"] = "#{scrbody}.F20"
393 ENV["DFTINTS"] = "#{scrbody}.F21"
394 ENV["DFTGRID"] = "#{scrbody}.F22"
395 ENV["JKFILE"] = "#{scrbody}.F23"
396 ENV["ORDINT"] = "#{scrbody}.F24"
397 ENV["EFPIND"] = "#{scrbody}.F25"
398 ENV["PCMDATA"] = "#{scrbody}.F26"
399 ENV["PCMINTS"] = "#{scrbody}.F27"
400 ENV["MLTPL"] = "#{scrbody}.F28"
401 ENV["MLTPLT"] = "#{scrbody}.F29"
402 ENV["DAFL30"] = "#{scrbody}.F30"
403 ENV["SOINTX"] = "#{scrbody}.F31"
404 ENV["SOINTY"] = "#{scrbody}.F32"
405 ENV["SOINTZ"] = "#{scrbody}.F33"
406 ENV["SORESC"] = "#{scrbody}.F34"
407 ENV["SIMEN"] = "#{scrbody}.simen"
408 ENV["SIMCOR"] = "#{scrbody}.simcor"
409 ENV["GCILIST"] = "#{scrbody}.F37"
410 ENV["HESSIAN"] = "#{scrbody}.F38"
411 ENV["SOCCDAT"] = "#{scrbody}.F40"
412 ENV["AABB41"] = "#{scrbody}.F41"
413 ENV["BBAA42"] = "#{scrbody}.F42"
414 ENV["BBBB43"] = "#{scrbody}.F43"
415 ENV["MCQD50"] = "#{scrbody}.F50"
416 ENV["MCQD51"] = "#{scrbody}.F51"
417 ENV["MCQD52"] = "#{scrbody}.F52"
418 ENV["MCQD53"] = "#{scrbody}.F53"
419 ENV["MCQD54"] = "#{scrbody}.F54"
420 ENV["MCQD55"] = "#{scrbody}.F55"
421 ENV["MCQD56"] = "#{scrbody}.F56"
422 ENV["MCQD57"] = "#{scrbody}.F57"
423 ENV["MCQD58"] = "#{scrbody}.F58"
424 ENV["MCQD59"] = "#{scrbody}.F59"
425 ENV["MCQD60"] = "#{scrbody}.F60"
426 ENV["MCQD61"] = "#{scrbody}.F61"
427 ENV["MCQD62"] = "#{scrbody}.F62"
428 ENV["MCQD63"] = "#{scrbody}.F63"
429 ENV["MCQD64"] = "#{scrbody}.F64"
430 ENV["NMRINT1"] = "#{scrbody}.F61"
431 ENV["NMRINT2"] = "#{scrbody}.F62"
432 ENV["NMRINT3"] = "#{scrbody}.F63"
433 ENV["NMRINT4"] = "#{scrbody}.F64"
434 ENV["NMRINT5"] = "#{scrbody}.F65"
435 ENV["NMRINT6"] = "#{scrbody}.F66"
436 ENV["DCPHFH2"] = "#{scrbody}.F67"
437 ENV["DCPHF21"] = "#{scrbody}.F68"
438 ENV["GVVPT"] = "#{scrbody}.F69"
440 # next files are used only during coupled cluster runs, so let's
441 # display the numerous definitions only if they are to be used.
442 # Don't set these variables on windows, where the memory for environmental variables
444 if $platform != "win"
445 ENV["CCREST"] = "#{scrbody}.F70"
446 ENV["CCDIIS"] = "#{scrbody}.F71"
447 ENV["CCINTS"] = "#{scrbody}.F72"
448 ENV["CCT1AMP"] = "#{scrbody}.F73"
449 ENV["CCT2AMP"] = "#{scrbody}.F74"
450 ENV["CCT3AMP"] = "#{scrbody}.F75"
451 ENV["CCVM"] = "#{scrbody}.F76"
452 ENV["CCVE"] = "#{scrbody}.F77"
453 ENV["EOMSTAR"] = "#{scrbody}.F80"
454 ENV["EOMVEC1"] = "#{scrbody}.F81"
455 ENV["EOMVEC2"] = "#{scrbody}.F82"
456 ENV["EOMHC1"] = "#{scrbody}.F83"
457 ENV["EOMHC2"] = "#{scrbody}.F84"
458 ENV["EOMHHHH"] = "#{scrbody}.F85"
459 ENV["EOMPPPP"] = "#{scrbody}.F86"
460 ENV["EOMRAMP"] = "#{scrbody}.F87"
461 ENV["EOMRTMP"] = "#{scrbody}.F88"
462 ENV["EOMDG12"] = "#{scrbody}.F89"
463 ENV["MMPP"] = "#{scrbody}.F90"
464 ENV["MMHPP"] = "#{scrbody}.F91"
465 ENV["MMCIVEC"] = "#{scrbody}.F92"
466 ENV["MMCIVC1"] = "#{scrbody}.F93"
467 ENV["MMCIITR"] = "#{scrbody}.F94"
468 ENV["MMNEXM"] = "#{scrbody}.F95"
469 ENV["MMNEXE"] = "#{scrbody}.F96"
470 ENV["MMNREXM"] = "#{scrbody}.F97"
471 ENV["MMNREXE"] = "#{scrbody}.F98"
474 # next are for TDHFX code, not used by current GAMESS
477 ENV["OLI201"] = "#{scrbody}.F201"
478 ENV["OLI202"] = "#{scrbody}.F202"
479 ENV["OLI203"] = "#{scrbody}.F203"
480 ENV["OLI204"] = "#{scrbody}.F204"
481 ENV["OLI205"] = "#{scrbody}.F205"
482 ENV["OLI206"] = "#{scrbody}.F206"
483 ENV["OLI207"] = "#{scrbody}.F207"
484 ENV["OLI208"] = "#{scrbody}.F208"
485 ENV["OLI209"] = "#{scrbody}.F209"
486 ENV["OLI210"] = "#{scrbody}.F210"
487 ENV["OLI211"] = "#{scrbody}.F211"
488 ENV["OLI212"] = "#{scrbody}.F212"
489 ENV["OLI213"] = "#{scrbody}.F213"
490 ENV["OLI214"] = "#{scrbody}.F214"
491 ENV["OLI215"] = "#{scrbody}.F215"
492 ENV["OLI216"] = "#{scrbody}.F216"
493 ENV["OLI217"] = "#{scrbody}.F217"
494 ENV["OLI218"] = "#{scrbody}.F218"
495 ENV["OLI219"] = "#{scrbody}.F219"
496 ENV["OLI220"] = "#{scrbody}.F220"
497 ENV["OLI221"] = "#{scrbody}.F221"
498 ENV["OLI222"] = "#{scrbody}.F222"
499 ENV["OLI223"] = "#{scrbody}.F223"
500 ENV["OLI224"] = "#{scrbody}.F224"
501 ENV["OLI225"] = "#{scrbody}.F225"
502 ENV["OLI226"] = "#{scrbody}.F226"
503 ENV["OLI227"] = "#{scrbody}.F227"
504 ENV["OLI228"] = "#{scrbody}.F228"
505 ENV["OLI229"] = "#{scrbody}.F229"
506 ENV["OLI230"] = "#{scrbody}.F230"
507 ENV["OLI231"] = "#{scrbody}.F231"
508 ENV["OLI232"] = "#{scrbody}.F232"
509 ENV["OLI233"] = "#{scrbody}.F233"
510 ENV["OLI234"] = "#{scrbody}.F234"
511 ENV["OLI235"] = "#{scrbody}.F235"
512 ENV["OLI236"] = "#{scrbody}.F236"
513 ENV["OLI237"] = "#{scrbody}.F237"
514 ENV["OLI238"] = "#{scrbody}.F238"
515 ENV["OLI239"] = "#{scrbody}.F239"
518 if $platform == "win"
523 # Old (obsolete) cygwin version, using ddikick
524 fpout.print "ddikick will run #{ncpus} compute process\n"
525 ENV["CYGWIN"] = "nodosfilewarning"
527 fpout.print "Microsoft MPI will be running GAMESS on 1 node.\n"
528 fpout.print "The binary kicked off by 'mpiexec' is gamess.#{gmsvers}.exe\n"
529 fpout.print "MS-MPI will run #{ncpus*2} compute process\n"
530 # File containing environmental variables
531 envfil = "#{scrprefix}.GMS.ENV"
532 fp = File.open(envfil, "w")
533 ENV.each { |k, v| fp.print "#{k}=#{v}\n" }
535 # File containing arguments to mpiexec
536 procfil = "#{scrprefix}.processes.mpd"
537 fp = File.open(procfil, "w")
538 fp.print "-env ENVFIL \"#{envfil}\" -wdir \"#{scrdir}\" -n #{ncpus*2} \"#{gmsdir}#{sep}gamess.#{gmsvers}.exe\"\n"
544 fplog = File.open(logname, "r")
552 ne_alpha = ne_beta = 0
557 term_callback = lambda { |m, n|
558 msg = "GAMESS execution on #{inpbase} "
565 msg += "failed with status #{n}."
569 msg += "\n(In directory #{inpdir})"
571 File.delete("#{scrdir}/.in_use")
573 ext_to_keep = [".dat", ".rst", ".trj", ".efp", ".gamma", ".log"]
574 ext_to_keep.each { |ex|
575 if File.exists?("#{scrprefix}#{ex}")
576 filecopy("#{scrprefix}#{ex}", "#{inpdir}#{sep}#{inpbody}#{ex}")
579 Dir.foreach(scrdir) { |file|
580 if file != "." && file != ".." && !ext_to_keep.include?(File.extname(file))
581 File.delete("#{scrdir}#{sep}#{file}")
586 erase_old_logs(scrdir, "latest", 5)
588 error_message_box($!.to_s)
592 if (script = get_global_settings("gamess.postfix_script")) != nil && script != ""
597 message_box(msg, hmsg, :ok, icon)
601 timer_callback = lambda { |m, n|
602 fplog.seek(0, IO::SEEK_END)
606 fplog.seek(size, IO::SEEK_SET)
607 fplog.each_line { |line|
608 if line[-1, 1] == "\n"
609 lines.push(last_line + line)
619 while i < lines.count
621 if line =~ /GEOMETRY SEARCH POINT NSERCH= *(\d+)/
626 elsif line =~ /START OF DRC CALCULATION/
630 elsif line =~ /CONSTRAINED OPTIMIZATION POINT/
633 elsif line =~ /POINT *([0-9]+) *ON THE REACTION PATH/
636 elsif line =~ /ATOMIC BASIS SET/
638 while j < lines.count
640 break if line =~ /TOTAL NUMBER OF BASIS SET/
648 bs_lines.push(lines[ii].chomp)
653 ncomps = mol.sub_load_gamess_log_basis_set(bs_lines, lineno + i)
656 $stderr.write($!.to_s + "\n")
657 $stderr.write($!.backtrace.inspect + "\n")
661 break # Wait until all basis set lines are read
663 elsif line =~ /NUMBER OF OCCUPIED ORBITALS/
671 elsif line =~ /SCFTYP=(\w+)/
673 if ne_alpha > 0 || ne_beta > 0
682 elsif line =~ /^\s*(EIGENVECTORS|MOLECULAR ORBITALS)\s*$/
683 if mo_count == 0 && mol
684 mol.clear_mo_coefficients
685 mol.set_mo_info(:type=>["UHF", "RHF", "ROHF"][rflag], :alpha=>ne_alpha, :beta=>ne_beta)
690 while j < lines.count
692 break if line =~ /\.\.\.\.\.\./ || line =~ /----------------/
696 break # Wait until complete MO info are read
701 mo_lines.push(lines[ii].chomp)
706 mol.sub_load_gamess_log_mo_coefficients(mo_lines, lineno + i, ncomps)
709 $stderr.write($!.to_s + "\n")
710 $stderr.write($!.backtrace.inspect + "\n")
713 elsif line =~ /NSERCH: *([0-9]+)/
717 line =~ /E= +([-.0-9]+)/
719 line =~ /GRAD\. MAX= +([-.0-9]+)/
721 mol.show_text("Search: #{n}\nGradient: #{grad}")
722 # mol.set_property("energy", energy)
725 elsif line =~ /TOTAL ENERGY += *([-.0-9]+)/
727 if mol && search_mode == 2
728 mol.show_text("Point: #{nserch}\nEnergy: #{energy}")
731 elsif line =~ /FINAL .* ENERGY IS +([-.0-9]+) *AFTER/
733 if mol && search_mode == 0
734 mol.set_property("energy", energy)
736 elsif nserch > 0 && line =~ /ddikick\.x/
739 elsif mol && nserch > 0 && (line =~ /COORDINATES OF ALL ATOMS/ || line =~ /COORDINATES \(IN ANGSTROM\) FOR \$DATA GROUP ARE/)
740 # There should be (natoms) lines
741 if line =~ /COORDINATES OF ALL ATOMS/
745 if i + mol.natoms + 1 <= lines.count
747 (i + 1...i + 1 + mol.natoms).each { |j|
748 name, charge, x, y, z = lines[j].split
749 coords.push(Vector3D[x.to_f, y.to_f, z.to_f])
751 mol.create_frame([coords])
752 if search_mode == 1 && energy
753 mol.set_property("energy", energy)
756 last_i = i + mol.natoms
757 i = last_i # Skip the processed lines
759 elsif line =~ /N A T U R A L B O N D O R B I T A L A N A L Y S I S/
762 while j < lines.count
763 break if lines[j] =~ /done with NBO analysis/
767 nbo_lines = lines[i..j]
768 mol.import_nbo_log(nbo_lines) rescue puts "Error: #{$!}, #{$!.backtrace.inspect}"
770 i = last_i # Skip the processed lines
772 break # Wait until NBO is done
781 lines[0..last_i] = nil
788 if (script = get_global_settings("gamess.prefix_script")) != nil && script != ""
792 if $platform == "win"
794 hosts = "localhost " * ncpus
795 cmdline = "cmd.exe /c \"#{gmsdir}/ddikick.exe #{gmsdir}/gamess.#{gmsvers}.exe #{inpbody} -ddi #{ncpus} #{ncpus} #{hosts} -scr #{scrdir} < NUL >>#{logname}"
797 cmdline = "cmd.exe /c \"mpiexec -configfile #{procfil} >>#{logname}"
800 hosts = "localhost " * ncpus
801 cmdline = "/bin/sh -c '#{gmsdir}/ddikick.x #{gmsdir}/gamess.#{gmsvers}.x #{inpbody} -ddi #{ncpus} #{ncpus} #{hosts} -scr #{scrdir} < /dev/null >>#{logname}'"
805 pid = mol.call_subprocess_async(cmdline, term_callback, timer_callback)
807 error_message_box("GAMESS failed to start. Please examine GAMESS installation.")
811 pid = call_subprocess(cmdline, "Running GAMESS")
812 term_callback.call(nil, pid)
818 def export_gamess(fname, hash)
820 def reorder_array(ary, ordered_sub_ary)
821 return (ordered_sub_ary & ary) + (ary - ordered_sub_ary)
826 basename = File.basename(fname, ".*")
828 basename = File.basename(self.name, ".*")
832 icharg = hash["charge"]
834 runtyp = ["ENERGY", "PROP", "OPTIMIZE", "SADPOINT", "IRC"][hash["runtype"].to_i]
835 scftyp = ["RHF", "ROHF", "UHF"][hash["scftype"].to_i]
836 bssname = hash["basis"]
837 bssname2 = hash["secondary_basis"]
838 if hash["use_secondary_basis"].to_i != 0 && bssname2 != bssname
840 element2 = hash["secondary_elements"].split(/[\s,]+/).map { |name| name.capitalize }
845 basis = $gamess_basis[bssname]
846 basis2 = $gamess_basis[bssname2]
848 # Use effective core potentials?
849 ecp = $gamess_ecp[bssname]
850 ecp2 = $gamess_ecp[bssname2]
851 ecp_read = (ecp || ecp2 ? "ECP=READ " : "")
853 # Use only one built-in basis set?
858 gbasis = "PM3 NGAUSS=3"
860 gbasis = "STO NGAUSS=3"
862 gbasis = "N21 NGAUSS=3"
864 gbasis = "N31 NGAUSS=6"
866 gbasis = "N31 NGAUSS=6 NDFUNC=1"
868 gbasis = "N31 NGAUSS=6 NDFUNC=1 NPFUNC=1"
870 gbasis = "N311 NGAUSS=6"
872 gbasis = "N311 NGAUSS=6 NDFUNC=1 NPFUNC=1"
876 # Count non-dummy atoms
879 natoms += 1 if ap.atomic_number != 0
882 # Fill hash with default values
883 h = (hash["CONTRL"] ||= Hash.new)
884 h["COORD"] ||= "UNIQUE"
885 h["EXETYP"] ||= "RUN"
886 h["ICHARG"] ||= (icharg || 0).to_s
888 h["INTTYP"] ||= "HONDO"
891 h["MOLPLT"] ||= ".T."
893 h["MULT"] ||= mult.to_s
894 h["QMTTOL"] ||= "1e-08"
895 h["RUNTYP"] ||= runtyp
896 if (hash["dft"] || 0) != 0 && hash["dfttype"]
897 h["DFTTYP"] ||= hash["dfttype"]
899 if (hash["use_internal"] || 0) != 0 && (hash["runtype"] >= 2 || h["RUNTYP"] == "OPTIMIZE" || h["RUNTYP"] == "SADPOINT" || h["RUNTYP"] == "IRC")
900 nzvar = natoms * 3 - 6 # TODO: 3N-5 for linear molecules
901 h["NZVAR"] ||= nzvar.to_s
905 h["SCFTYP"] ||= scftyp
909 h["UNITS"] ||= "ANGS"
911 h = (hash["SCF"] ||= Hash.new)
912 h["CONV"] ||= "1.0E-06"
913 h["DIRSCF"] ||= ".T."
917 h = (hash["STATPT"] ||= Hash.new)
919 if runtyp == "SADPOINT"
920 h["OPTTOL"] ||= "1.0E-05"
922 h["HSSEND"] ||= ".T."
923 elsif runtyp == "IRC"
924 h["OPTTOL"] ||= "1.0E-05"
926 h["HSSEND"] ||= ".T."
928 h["OPTTOL"] ||= "1.0E-06"
930 if hash["eliminate_freedom"] == 0
931 h["PROJCT"] ||= ".F."
934 h = (hash["SYSTEM"] ||= Hash.new)
936 if runtyp == "SADPOINT" || runtyp == "IRC"
941 h["TIMLIM"] ||= "50000"
943 h = (hash["GUESS"] ||= Hash.new)
944 h["GUESS"] ||= "HUCKEL"
947 h = (hash["BASIS"] ||= Hash.new)
948 h["GBASIS"] ||= gbasis
952 h = (hash["ZMAT"] ||= Hash.new)
958 h = (hash["IRC"] ||= Hash.new)
964 h["OPTTOL"] = "1.0E-05"
967 if (hash["esp"] || 0) != 0
968 h = (hash["ELPOT"] ||= Hash.new)
970 h["OUTPUT"] ||= "PUNCH"
972 h = (hash["PDC"] ||= Hash.new)
973 h["CONSTR"] ||= "NONE"
974 h["PTSEL"] ||= "CONNOLLY"
977 if (hash["include_nbo"] || 0) != 0
978 h = (hash["NBO"] ||= Hash.new)
980 ["nao", "nbo", "nho", "nlmo", "pnao", "pnbo", "pnho", "pnlmo"].each { |nao|
981 if (hash[nao] || 0) != 0
982 s += " f" + nao + " ao" + nao
989 fp = File.open(fname, "wb")
994 fp.print "! GAMESS input\n"
995 fp.print "! Generated by Molby at #{now}\n"
996 fp.print "! Basis set: " + ($gamess_basis_desc[bssname] || "(not specified)") + "\n"
998 fp.print "! [" + element2.join(", ") + "]: " + ($gamess_basis_desc[bssname2] || "(not specified)") + "\n"
1000 ordered = ["CONTRL", "SCF", "STATPT", "SYSTEM", "GUESS", "BASIS", "ZMAT", "ELPOT", "PDC"]
1001 controls = reorder_array(hash.keys.select { |k| k != "key_order" && hash[k].is_a?(Hash) },
1002 hash["key_order"] || ordered)
1005 next if h == nil || h.size == 0 || (h["key_order"] && h.size == 1)
1007 ordered = ["COORD", "EXETYP", "ICHARG", "ICUT", "INTTYP", "ITOL", "MAXIT", "MOLPLT", "MPLEVL",
1008 "MULT", "QMTTOL", "RUNTYP", "SCFTYP", "ECP", "UNITS", "DFTTYP", "NZVAR"]
1010 ordered = ["CONV", "DIRSCF", "FDIFF", "DAMP"]
1012 ordered = ["NSTEP", "OPTTOL", "HESS", "HSSEND"]
1014 ordered = ["MEMDDI", "MWORDS", "TIMLIM"]
1018 ordered = ["GBASIS"]
1020 ordered = ["DLC", "AUTO"]
1022 ordered = ["IEPOT", "OUTPUT", "WHERE"]
1024 ordered = ["CONSTR", "PTSEL"]
1026 ordered = ["SADDLE", "TSENGY", "FORWRD", "NPOINT", "STRIDE", "OPTTOL"]
1030 keys = reorder_array(h.keys, h["key_order"] || ordered)
1032 keys.each_with_index { |kk, i|
1035 fp.printf " $%-6s", k
1039 if v == nil || v == ""
1043 fp.printf " %s=%s", kk, h[kk].to_s
1051 fp.print " $DATA\n#{basename}\nC1 0\n"
1054 next if ap.atomic_number == 0
1055 fp.printf "%-6s %4d %10.6f %10.6f %10.6f\n", ap.name, ap.atomic_number, ap.r.x, ap.r.y, ap.r.z
1056 if use_2nd && element2.include?(ap.element)
1057 secondary[ap.index] = true
1060 # Basis specification followed by a blank line
1061 bas = (secondary[ap.index] ? basis2 : basis)
1063 bas = bas[ap.atomic_number]
1069 puts "Warning: basis set is not defined for atom #{ap.index}, element #{ap.element}"
1078 an = ap.atomic_number
1080 ecpp = (secondary[ap.index] ? ecp2 : ecp)
1081 e = ecp_ary[an] || (ecpp && ecpp[an])
1083 # Cache the PNAME of the $ECP entry and re-use it
1084 ecp_ary[an] ||= (e.split)[0] + "\n"
1086 e = ap.element.upcase + "-ECP NONE\n"
1103 def copy_section_from_gamess_output
1104 if @gamess_output_fname
1105 dir = File.dirname(@gamess_output_fname)
1109 fname = Dialog.open_panel("Select GAMESS Output:", dir, "*.log;*.dat")
1111 fp = open(fname, "rb")
1114 if @gamess_output_fname == fname && @gamess_output_mtime == mtime
1115 # Use the cached value
1116 k = @gamess_output_sections
1121 if ln.start_with?(" $")
1123 if keyword !~ /\$END/
1129 # Cache the information
1130 @gamess_output_fname = fname
1131 @gamess_output_mtime = mtime
1132 @gamess_output_sections = k
1134 keywords = k.keys.sort { |a, b| k[a] <=> k[b] }
1135 h = Dialog.run("Select GAMESS section to copy", "Copy", "Cancel") {
1137 item(:popup, :subitems=>keywords, :tag=>"section"))
1140 fp.seek(k[keywords[h["section"]]])
1144 break if ln =~ /\$END/
1146 export_to_clipboard(s)
1153 def cmd_edit_gamess_input(s)
1155 h = Dialog.run("Edit GAMESS Input", "OK", "Cancel", :resizable=>true) {
1157 item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
1158 item(:button, :title=>"Copy Section from GAMESS Output...", :action=>lambda { |it| mol.copy_section_from_gamess_output } ),
1159 :flex=>[0,0,0,0,1,1]
1161 set_min_size(300, 300)
1170 def cmd_create_gamess_input
1175 raise MolbyError, "cannot create GAMESS input; the molecule is empty"
1178 # Descriptive text and internal string for popup menus
1179 # bset_desc = ["PM3", "STO-3G", "3-21G", "6-31G", "6-31G(d)", "6-31G(d,p)", "6-311G", "6-311G(d,p)", "LanL2DZ"]
1180 # bset_internal = ["PM3", "STO3G", "321G", "631G", "631Gd", "631Gdp", "6311G", "6311Gdp", "LanL2DZ"]
1181 bset_desc = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
1182 dft_desc = ["B3LYP"]
1183 dft_internal = ["B3LYP"]
1185 defaults = {"scftype"=>0, "runtype"=>0, "use_internal"=>1, "eliminate_freedom"=>1, "charge"=>"0", "mult"=>"1",
1186 "basis"=>4, "use_secondary_basis"=>0, "secondary_elements"=>"",
1187 "secondary_basis"=>8, "esp"=>0, "ncpus"=>"1"}
1189 gamess_input_direct = nil
1191 user_input = Hash.new
1192 ["CONTRL", "SCF", "STATPT", "SYSTEM", "GUESS", "BASIS"].each { |k|
1193 user_input[k] = Hash.new
1195 if ".inp".casecmp(self.name[-4, 4]) == 0
1196 # Read input and analyze commands
1197 fp = open(self.path, "r")
1202 break if "$DATA".casecmp(ln[0, 5]) == 0
1206 key = s[1..-1].upcase
1207 hh = user_input[key] = Hash.new
1208 (user_input["key_order"] ||= []).push(key)
1214 (hh["key_order"] ||= []).push(k)
1220 user_input["charge"] = user_input["CONTRL"]["ICHARG"]
1221 user_input["mult"] = user_input["CONTRL"]["MULT"]
1222 user_input["runtype"] = ["ENERGY", "PROP", "OPTIMIZE"].find_index(user_input["CONTRL"]["RUNTYP"])
1223 user_input["scftype"] = ["RHF", "ROHF", "UHF"].find_index(user_input["CONTRL"]["SCFTYP"])
1224 dft_type = dft_internal.find_index(user_input["CONTRL"]["DFTTYP"])
1226 user_input["dfttype"] = dft_type
1227 user_input["dft"] = 1
1230 user_input["basis"] = "-1"
1231 case user_input["BASIS"]["GBASIS"]
1239 if user_input["NDFUNC"] == "1"
1240 if user_input["NPFUNC"] == "1"
1249 if user_input["NDFUNC"] == "1" && user_input["NPFUNC"] == "1"
1256 user_input["basis"] = $gamess_basis_keys.find_index(bssname).to_s
1258 # puts user_input.inspect
1262 hash = Dialog.run("GAMESS Export") {
1263 def load_basis_set_sub(item)
1264 fname = Dialog.open_panel("Select a file containing GAMESS basis set:")
1266 Molecule.read_gamess_basis_sets(fname)
1267 bset_desc_new = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
1268 sel1 = attr("basis", :value)
1269 sel2 = attr("secondary_basis", :value)
1270 set_attr("basis", :subitems=>bset_desc_new)
1271 set_attr("basis", :value=>sel1)
1272 set_attr("secondary_basis", :subitems=>bset_desc_new)
1273 set_attr("secondary_basis", :value=>sel2)
1276 def select_gamess_path(item)
1278 fname = Dialog.open_panel("Locate GAMESS executable:")
1279 return if fname == nil
1280 bname = File.basename(fname)
1281 if bname =~ /gamess\.(.*)\.(exe|x)$/i
1282 set_attr("executable_path", :value=>fname)
1285 error_message_box("\"#{bname}\" does not look like a GAMESS executable! Please try again.")
1289 def set_optional_scripts(item)
1290 h = Dialog.run("GAMESS Optional Scripts") {
1291 s_pre = get_global_settings("gamess.prefix_script")
1292 s_post = get_global_settings("gamess.postfix_script")
1294 item(:text, :title=>"Script to run before GAMESS execution:"),
1295 item(:textview, :width=>400, :height=>200, :value=>s_pre, :tag=>"prefix"),
1296 item(:text, :title=>"Script to run after GAMESS execution:"),
1297 item(:textview, :width=>400, :height=>200, :value=>s_post, :tag=>"postfix"))
1300 set_global_settings("gamess.prefix_script", h["prefix"])
1301 set_global_settings("gamess.postfix_script", h["postfix"])
1304 nbos = ["nao", "nbo", "nho", "nlmo", "pnao", "pnbo", "pnho", "pnlmo"]
1306 item(:text, :title=>"SCF type"),
1307 item(:popup, :subitems=>["RHF", "ROHF", "UHF"], :tag=>"scftype"),
1308 item(:text, :title=>"Run type"),
1309 item(:popup, :subitems=>["Energy", "Property", "Optimize", "Sadpoint", "IRC"], :tag=>"runtype",
1310 :action=>lambda { |it| set_attr("use_internal", :enabled=>(it[:value] >= 2)) } ),
1311 item(:checkbox, :title=>"Use internal coordinates for structure optimization", :tag=>"use_internal"),
1313 item(:checkbox, :title=>"Eliminate translation and rotational degrees of freedom", :tag=>"eliminate_freedom"),
1316 item(:text, :title=>"Charge"),
1317 item(:textfield, :width=>80, :tag=>"charge"),
1318 item(:text, :title=>"Multiplicity"),
1319 item(:textfield, :width=>80, :tag=>"mult"),
1321 item(:checkbox, :title=>"Use DFT", :tag=>"dft",
1322 :action=>lambda { |it| set_attr("dfttype", :enabled=>(it[:value] != 0)) } ),
1324 item(:text, :title=>"DFT type"),
1325 item(:popup, :subitems=>dft_desc, :tag=>"dfttype"),
1330 item(:text, :title=>"Basis set"),
1331 item(:popup, :subitems=>bset_desc + ["(no select)"], :tag=>"basis"),
1335 item(:button, :title=>"Load Basis Set...", :action=>:load_basis_set_sub),
1338 item(:checkbox, :title=>"Use secondary basis set", :tag=>"use_secondary_basis",
1339 :action=>lambda { |it|
1340 flag = (it[:value] != 0)
1341 set_attr("secondary_elements", :enabled=>flag)
1342 set_attr("secondary_basis", :enabled=>flag)
1346 item(:text, :title=>" Elements"),
1347 item(:textfield, :width=>80, :tag=>"secondary_elements"),
1348 item(:text, :title=>"Basis set"),
1349 item(:popup, :subitems=>bset_desc, :tag=>"secondary_basis"),
1354 item(:checkbox, :title=>"Calculate electrostatic potential (ESP)", :tag=>"esp"),
1360 item(:checkbox, :title=>"Include NBO instructions", :tag=>"include_nbo",
1361 :action=>lambda { |it|
1362 flag = (it[:value] != 0)
1363 nbos.each { |nbo| set_attr(nbo, :enabled=>flag) }
1366 item(:checkbox, :title=>"NAO", :tag=>"nao"),
1367 item(:checkbox, :title=>"NBO", :tag=>"nbo"),
1368 item(:checkbox, :title=>"NHO", :tag=>"nho"),
1369 item(:checkbox, :title=>"NLMO", :tag=>"nlmo"),
1370 item(:checkbox, :title=>"PNAO", :tag=>"pnao"),
1371 item(:checkbox, :title=>"PNBO", :tag=>"pnbo"),
1372 item(:checkbox, :title=>"PNHO", :tag=>"pnho"),
1373 item(:checkbox, :title=>"PNLMO", :tag=>"pnlmo"),
1378 item(:checkbox, :title=>"Execute GAMESS on this machine", :tag=>"execute_local",
1379 :action=>lambda { |it|
1380 flag = (it[:value] != 0)
1381 set_attr("executable_path", :enabled=>flag)
1382 set_attr("select_path", :enabled=>flag)
1383 set_attr("ncpus", :enabled=>flag)
1387 item(:text, :title=>" Path"),
1388 item(:textfield, :width=>300, :tag=>"executable_path"),
1392 item(:button, :title=>"Select Path...", :tag=>"select_path", :action=>:select_gamess_path),
1394 item(:button, :title=>"Optional Scripts...", :action=>:set_optional_scripts),
1396 item(:text, :title=>" N of CPUs"),
1397 item(:textfield, :width=>80, :tag=>"ncpus"),
1403 item(:button, :title=>"Edit GAMESS Input and Go", :action=>lambda { |it|
1406 if (tag = it2[:tag]) != nil
1407 h[tag] = it2[:value]
1410 h["basis"] = $gamess_basis_keys[h["basis"]]
1411 h["secondary_basis"] = $gamess_basis_keys[h["secondary_basis"]]
1412 h["dfttype"] = dft_internal[h["dfttype"]]
1413 gamess_input_direct = mol.cmd_edit_gamess_input(mol.export_gamess(nil, h))
1414 if gamess_input_direct
1424 values[tag] = (user_input[tag] || get_global_settings("gamess.#{tag}") || defaults[tag])
1425 if tag == "basis" && values[tag] == "-1"
1426 values[tag] = (bset_desc.count).to_s
1428 it[:value] = values[tag]
1431 set_attr("secondary_elements", :enabled=>(values["use_secondary_basis"] == 1))
1432 set_attr("secondary_basis", :enabled=>(values["use_secondary_basis"] == 1))
1433 set_attr("dfttype", :enabled=>(values["dft"] == 1))
1434 set_attr("use_internal", :enabled=>(values["runtype"] >= 2))
1435 set_attr("executable_path", :enabled=>(values["execute_local"] == 1))
1436 set_attr("select_path", :enabled=>(values["execute_local"] == 1))
1437 set_attr("ncpus", :enabled=>(values["execute_local"] == 1))
1439 set_attr(nao, :enabled=>(values["include_nbo"] == 1))
1442 hash.each_pair { |key, value|
1443 next if key == :status
1444 set_global_settings("gamess.#{key}", value)
1446 if hash[:status] == 0
1447 # Specify basis by internal keys
1448 hash["basis"] = $gamess_basis_keys[hash["basis"]]
1449 hash["secondary_basis"] = $gamess_basis_keys[hash["secondary_basis"]]
1450 hash["dfttype"] = dft_internal[hash["dfttype"]]
1451 basename = (self.path ? File.basename(self.path, ".*") : self.name)
1452 fname = Dialog.save_panel("Export GAMESS input file:", self.dir, basename + ".inp", "GAMESS input file (*.inp)|*.inp|All files|*.*")
1453 return nil if !fname
1454 if gamess_input_direct
1455 # puts "gamess_input_direct = \"#{gamess_input_direct}\""
1456 File.open(fname, "w") { |fp| fp.print(gamess_input_direct) }
1458 export_gamess(fname, hash)
1460 if hash["execute_local"] == 1
1461 Molecule.execute_gamess(fname, self)
1470 grid = default_MO_grid
1474 item(:text, :title=>"This molecule does not contain MO information."))
1480 if mos == nil || mos.length == 0
1483 item(:text, :title=>"Please select MO(s) in the MO Info table."))
1489 item(:text, :title=>"Please specify cube dimensions (in angstrom units):"),
1491 item(:text, :title=>"Origin"),
1492 item(:textfield, :width=>100, :height=>20, :tag=>"originx", :value=>sprintf("%.6f", grid[0].x)),
1493 item(:textfield, :width=>100, :height=>20, :tag=>"originy", :value=>sprintf("%.6f", grid[0].y)),
1494 item(:textfield, :width=>100, :height=>20, :tag=>"originz", :value=>sprintf("%.6f", grid[0].z)),
1495 item(:text, :title=>"Delta"),
1496 item(:textfield, :width=>100, :height=>20, :tag=>"deltax", :value=>sprintf("%.6f", grid[1])),
1497 item(:textfield, :width=>100, :height=>20, :tag=>"deltay", :value=>sprintf("%.6f", grid[2])),
1498 item(:textfield, :width=>100, :height=>20, :tag=>"deltaz", :value=>sprintf("%.6f", grid[3])),
1499 item(:text, :title=>"Step"),
1500 item(:textfield, :width=>100, :height=>20, :tag=>"stepx", :value=>grid[4].to_s),
1501 item(:textfield, :width=>100, :height=>20, :tag=>"stepy", :value=>grid[5].to_s),
1502 item(:textfield, :width=>100, :height=>20, :tag=>"stepz", :value=>grid[6].to_s)))
1504 if hash[:status] == 0
1505 path = self.path || self.name
1506 dir = self.dir || Dir.pwd
1507 origin = Vector3D[hash["originx"], hash["originy"], hash["originz"]]
1514 basename = File.basename(path, ".*")
1516 mo_type = self.mo_type
1518 fname1 = fname2 = nil
1519 alpha = (mo_type != "UHF" ? "" : "alpha ")
1520 a = (mo_type != "UHF" ? "" : "a")
1521 fname1 = Dialog.save_panel("Cube file name for #{alpha}MO #{n}", dir, basename + "_#{n}#{a}.cube", "Gaussian cube file (*.cube)|*.cube")
1522 if (mo_type == "UHF")
1523 fname2 = Dialog.save_panel("Cube file name for beta MO #{n}", dir, basename + "_#{n}b.cube", "Gaussian cube file (*.cube)|*.cube")
1525 filenames.push([n, fname1, fname2])
1527 filenames.each { |pair|
1529 alpha = (mo_type != "UHF" ? "" : "alpha ")
1530 show_progress_panel("Creating cube file for #{alpha}MO #{n}...")
1532 cubegen(pair[1], n, origin, dx, dy, dz, nx, ny, nz, true)
1534 if pair[2] && mo_type == "UHF"
1535 set_progress_message("Creating cube file for beta MO #{n}...")
1536 cubegen(pair[2], n, origin, dx, dy, dz, nx, ny, nz, true, true)
1543 def export_psi4(fname, hash)
1546 fp = File.open(fname, "wb")
1551 fp.print "# Psi4 input\n"
1552 fp.print "# Generated by Molby at #{now}\n"
1553 fp.print "import sys\n"
1554 fp.print "import re\n"
1555 fp.print "base = re.sub('\\.\\w*$', '', sys.argv[1]) # Basename of the input file\n"
1556 fp.print "molecule mol {\n"
1557 charge = Integer(hash["charge"])
1558 mult = Integer(hash["mult"])
1559 if charge != 0 || mult != 1
1560 fp.printf " %d %d\n", charge, mult
1563 fp.printf "%-8s %16.12f %16.12f %16.12f\n", ap.element, ap.x, ap.y, ap.z
1565 use_symmetry = Integer(hash["use_symmetry"])
1566 move_to_com = Integer(hash["move_to_com"])
1567 do_reorient = Integer(hash["do_reorient"])
1572 fp.print "noreorient\n"
1574 if use_symmetry == 0
1575 fp.print "symmetry c1\n"
1577 fp.print "units angstrom\n"
1579 fp.print "set basis #{hash['basis']}\n"
1580 fp.print "set reference #{hash['scftype']}\n"
1581 options = "return_wfn=True"
1583 options += ", dft_functional='#{hash['dfttype']}'"
1585 runtype = hash['runtype'].downcase
1586 fp.print "energy, wfn = #{runtype}('scf', #{options})\n"
1587 fp.print "wfn.write_molden(f'{base}.molden')\n"
1588 if hash['run_junpa'] != 0
1590 fp.print "# Interface for JANPA\n"
1591 fp.print "# cf. https://sourceforge.net/p/janpa/wiki/psi4Examples/\n"
1592 fp.print "d = wfn.Da().to_array()\n"
1593 fp.print "s = wfn.S().to_array()\n"
1594 fp.print "c = wfn.Ca().to_array()\n"
1595 fp.print "occs = c.T.dot(s.dot(d).dot(s).dot(c))\n"
1596 fp.print "molden(wfn, f'{base}.da.molden', density_a = psi4.core.Matrix.from_array(occs))\n"
1609 def Molecule.is_java_available
1610 if $platform == "win"
1611 f = get_global_settings("java_home")
1613 ENV["JAVA_HOME"] = f
1614 if !ENV["PATH"].split(";").find { |e| e == "#{f}\\bin" }
1615 ENV["PATH"] = "#{f}\\bin;" + ENV["PATH"]
1619 return (call_subprocess("java -version", nil) == 0)
1622 def Molecule.make_java_available
1623 if $platform == "win"
1624 fname = Dialog.open_panel("Locate JDK Folder (if you have one):", "c:\\", nil, true)
1625 return false if fname == nil
1626 fname.sub!(/\//, "\\")
1627 if File.exists?("#{fname}\\bin\\java.exe")
1628 set_global_settings("java_home", fname)
1629 if Molecule.is_java_available()
1633 error_message_box("Cannot run Java. Please examine your installation again.")
1635 elsif $platform == "mac"
1636 message_box("Please download OpenJDK, and move it into /Library/Java/JavaVirtualMachines folder.", "Install Java", :ok)
1639 message_box("Please install Java virtual machine.", "Install Java", :ok)
1645 # inppath is the input file minus extention
1646 # mol is the molecule (may be nil)
1647 def Molecule.execute_janpa(inppath, mol, spherical)
1648 # nbo_desc = # JANPA
1649 janpa_dir = "#{ResourcePath}/JANPA"
1651 outfile = "#{inppath}.janpa.log"
1653 cmd1 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-NormalizeBF", "-i", "#{inppath}.da.molden", "-o", "#{inppath}.in.molden"]
1656 cmd1 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-frompsi4v1mo", "-NormalizeBF", "-cart2pure", "-i", "#{inppath}.da.molden", "-o", "#{inppath}.in.molden"]
1657 cmd2 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-frompsi4v1mo", "-NormalizeBF", "-cart2pure", "-i", "#{inppath}.molden", "-o", "#{inppath}.spherical.molden"]
1659 cmd3 = ["java", "-jar", "#{janpa_dir}/janpa.jar", "-i", "#{inppath}.in.molden"]
1660 ["nao", "pnao", "aho", "lho", "lpo", "clpo"].each { |type|
1661 generate = (get_global_settings("psi4.#{type}").to_i != 0)
1662 if type == "pnao" || type == "nao" || type == "lho"
1663 # PLHO is generated within Molby from JANPA NAO/PNAO/LHO
1664 generate ||= (get_global_settings("psi4.plho").to_i != 0)
1667 cmd3.push("-#{type.upcase}_Molden_File", "#{inppath}.#{type.upcase}.molden")
1670 # show_progress_panel("Executing JANPA...")
1671 procname = "molden2molden"
1672 flag = call_subprocess(cmd1, procname, nil, outfile)
1673 if flag && cmd2 != nil
1674 procname = "molden2molden"
1675 flag = call_subprocess(cmd2, procname, nil, ">>#{outfile}")
1679 flag = call_subprocess(cmd3, procname, nil, ">>#{outfile}")
1683 # import JANPA log and molden output
1684 # Files: inppath.janpa.log, inppath.{NAO,PNAO,AHO,LHO,LPO,CLPO,spherical}.molden
1685 mol.sub_load_janpa_log(inppath)
1689 status = $?.exitstatus
1691 message_box("Execution of #{procname} failed with status #{status}.", "JANPA Failed")
1697 # inpname is the input file
1698 # mol (optional) is the molecule from which the Psi4 input was built.
1699 def Molecule.execute_psi4(inpname, mol = nil)
1701 inpbase = File.basename(inpname)
1702 inpbody = File.basename(inpname, ".*")
1703 inpdir = File.dirname(inpname)
1705 # Set PATH for psi4conda
1706 psi4folder = get_global_settings("psi4.psi4conda_folder")
1707 ncpus = get_global_settings("psi4.ncpus").to_i
1708 orgpath = ENV["PATH"]
1710 if $platform == "win"
1711 ENV["PATH"] = "#{psi4folder};#{psi4folder}\\Library\\mingw-w64\\bin;#{psi4folder}\\Library\\usr\\bin;#{psi4folder}\\Library\\bin;#{psi4folder}\\Scripts;#{psi4folder}\\bin;#{psi4folder}\\condabin;" + ENV["PATH"]
1713 ENV["PATH"] = "#{psi4folder}/bin:#{psi4folder}/condabin:" + ENV["PATH"]
1716 cmdargv = ["psi4", "#{inpbase}"]
1718 cmdargv.push("-n", "#{ncpus}")
1725 outfile = inpdir + "/" + inpbody + ".out"
1726 if File.exists?(outfile)
1729 outbackfile = inpdir + "/" + inpbody + "~" + (n == 1 ? "" : "#{n}") + ".out"
1730 break if !File.exists?(outbackfile)
1733 File.rename(outfile, outbackfile)
1745 timer_callback = lambda { |m, n|
1748 if timer_count < 10 # Only handle every 1 seconds
1753 if File.exists?(outfile)
1754 fplog = Kernel.open(outfile, "rt")
1760 return true # Skip until outfile is available
1763 fplog.seek(0, IO::SEEK_END)
1764 current_size = fplog.tell
1765 if current_size > last_size
1767 fplog.seek(last_size, IO::SEEK_SET)
1768 fplog.each_line { |line|
1769 if line[-1, 1] == "\n"
1770 lines.push(last_line + line)
1777 last_size = fplog.tell
1780 getline = lambda { if li < lines.length; li += 1; return lines[li - 1]; else return nil; end }
1782 while (line = getline.call)
1783 if line =~ /==> Geometry <==/
1784 # Skip until line containing "------"
1785 while line = getline.call
1786 break if line =~ /------/
1790 # Read atom positions
1791 while line = getline.call
1793 break if line =~ /^\s*$/
1794 tokens = line.split(' ')
1795 vecs.push(Vector3D[Float(tokens[1]), Float(tokens[2]), Float(tokens[3])])
1798 if vecs.length < mol.natoms
1799 break # Log file is incomplete
1801 # Does this geometry differ from the last one?
1802 vecs.length.times { |i|
1803 if (mol.atoms[i].r - vecs[i]).length2 > 1.0e-20
1804 # Create a new frame and break
1806 vecs.length.times { |j|
1807 mol.atoms[j].r = vecs[j]
1814 elsif line =~ /Final Energy: +([-.0-9]+)/
1815 # Energy for this geometry
1817 mol.set_property("energy", energy)
1818 mol.show_text("Frame: #{mol.nframes - 1}\nEnergy: #{energy}")
1819 print("Frame: #{mol.nframes - 1} Energy: #{energy}\n")
1824 elsif line =~ /ROHF/
1828 elsif line =~ /^ *Nalpha *= *(\d+)/
1829 nalpha = Integer($1)
1830 elsif line =~ /^ *Nbeta *= *(\d+)/
1832 elsif line =~ /^ *Spherical Harmonics\?: *(\w+)/
1833 spherical = ($1 == "true")
1837 lines.slice!(0, next_index)
1842 # Some error occurred during callback. Show the message in the console and abort.
1843 $stderr.write("#{e.message}\n")
1844 $stderr.write("#{e.backtrace.inspect}\n")
1849 # Terminate callback
1850 term_callback = lambda { |m, n|
1853 msg = "Psi4 execution of #{inpbase} "
1856 msg += "cannot be started. Please examine Psi4 installation."
1865 msg += "failed with status #{n}."
1869 msg += "\n(In directory #{inpdir})"
1872 # Try to load final lines of the logfile
1874 timer_callback.call(m, n)
1875 # Try to load molden file if available
1877 mol.clear_mo_coefficients
1878 mol.set_mo_info(:type => hf_type, :alpha => nalpha, :beta => nbeta)
1879 molden = inpdir + "/" + inpbody +".molden"
1880 if File.exists?(molden)
1881 fp = Kernel.open(molden, "rt")
1882 mol.instance_eval { @lineno = 0 }
1884 # mol.@hf_type should be set before calling sub_load_molden
1885 mol.sub_load_molden(fp)
1888 $stderr.write("#{e.message}\n")
1889 $stderr.write("#{e.backtrace.inspect}\n")
1892 if (get_global_settings("psi4.run_janpa").to_i == 1)
1894 Molecule.execute_janpa(inpdir + "/" + inpbody, mol, spherical)
1897 # The child process actually did not start
1898 # Restore the old file if outbackfile is not nil
1899 if outbackfile && !File.exists?(outfile)
1900 File.rename(outbackfile, outfile)
1903 if outbackfile && File.exists?(outbackfile)
1904 File.delete(outbackfile)
1906 ENV["PATH"] = orgpath
1909 message_box(msg, hmsg, :ok, icon)
1912 $stderr.write("#{e.message}\n")
1913 $stderr.write("#{e.backtrace.inspect}\n")
1920 pid = mol.call_subprocess_async(cmdargv, term_callback, timer_callback)
1922 # This may not happen on OSX or Linux (don't know for MSW)
1923 error_message_box("Psi4 failed to start. Please examine Psi4 installation.")
1927 status = call_subprocess(cmdargv, "Running Psi4")
1928 term_callback.call(nil, status)
1933 def cmd_edit_psi4_input(s)
1935 h = Dialog.run("Edit Psi4 Input", "OK", "Cancel", :resizable=>true) {
1937 item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
1939 set_min_size(300, 300)
1948 def cmd_create_psi4_input
1953 raise MolbyError, "cannot create Psi4 input; the molecule is empty"
1956 # Basis sets Cf. https://psicode.org/psi4manual/master/basissets_tables.html#apdx-basistables
1957 bset_desc = ["STO-3G", "3-21G", "6-31G", "6-31G(d)", "6-31G(d,p)", "6-311G", "6-311G(d)", "6-311G(d,p)", "LanL2DZ"]
1958 dfttype_desc = ["B3LYP"]
1959 runtype_desc = ["Energy", "Optimize"]
1960 scftype_desc = ["RHF", "ROHF", "UHF"]
1961 nbo_desc = ["nao", "pnao", "aho", "lho", "lpo", "clpo"] # JANPA
1962 user_input = Hash.new
1963 defaults = {"scftype"=>0, "runtype"=>0, "move_to_com"=>0, "do_reorient"=>0,
1964 "use_symmetry"=>0, "charge"=>"0", "mult"=>"1",
1965 "basis"=>1, "use_secondary_basis"=>0, "secondary_elements"=>"",
1966 "secondary_basis"=>8, "esp"=>0, "ncpus"=>"1"}
1967 psi4_input_direct = nil
1969 hash = Dialog.run("Psi4 Export") {
1970 def load_basis_set_sub(item)
1971 fname = Dialog.open_panel("Select a file containing GAMESS basis set:")
1973 Molecule.read_gamess_basis_sets(fname)
1974 bset_desc_new = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
1975 sel1 = attr("basis", :value)
1976 sel2 = attr("secondary_basis", :value)
1977 set_attr("basis", :subitems=>bset_desc_new)
1978 set_attr("basis", :value=>sel1)
1979 set_attr("secondary_basis", :subitems=>bset_desc_new)
1980 set_attr("secondary_basis", :value=>sel2)
1983 def select_psi4_folder(item)
1984 # By default, psi4conda is installed in the home directory
1986 fname = Dialog.open_panel("Locate 'psi4conda' Folder:", ENV["HOME"] || ENV["USERPROFILE"], nil, true)
1987 return if fname == nil
1988 bname = File.basename(fname)
1989 if bname == "psi4conda"
1992 h = Dialog.run("", "OK", "Choose Again") {
1994 item(:text, :title=>"The folder does not have the name 'psi4conda'."),
1995 item(:text, :title=>"Do you want to use it anyway?"))
2002 set_attr("psi4conda_folder", :value=>fname)
2006 item(:text, :title=>"SCF type"),
2007 item(:popup, :subitems=>scftype_desc, :tag=>"scftype"),
2008 item(:text, :title=>"Run type"),
2009 item(:popup, :subitems=>runtype_desc, :tag=>"runtype"),
2011 item(:checkbox, :title=>"Detect symmetry", :tag=>"use_symmetry"),
2014 item(:checkbox, :title=>"Move the molecule to the center of mass", :tag=>"move_to_com"),
2017 item(:checkbox, :title=>"Rotate the molecule to the symmetry principle axes", :tag=>"do_reorient"),
2020 item(:text, :title=>"Charge"),
2021 item(:textfield, :width=>80, :tag=>"charge"),
2022 item(:text, :title=>"Multiplicity"),
2023 item(:textfield, :width=>80, :tag=>"mult"),
2025 item(:checkbox, :title=>"Use DFT", :tag=>"dft",
2026 :action=>lambda { |it| set_attr("dfttype", :enabled=>(it[:value] != 0)) } ),
2028 item(:text, :title=>"DFT type"),
2029 item(:popup, :subitems=>dfttype_desc, :tag=>"dfttype"),
2034 item(:text, :title=>"Basis set"),
2035 item(:popup, :subitems=>bset_desc, :tag=>"basis"),
2039 #item(:button, :title=>"Load Basis Set...", :action=>:load_basis_set_sub),
2042 #item(:checkbox, :title=>"Use secondary basis set", :tag=>"use_secondary_basis",
2043 # :action=>lambda { |it|
2044 # flag = (it[:value] != 0)
2045 # set_attr("secondary_elements", :enabled=>flag)
2046 # set_attr("secondary_basis", :enabled=>flag)
2050 #item(:text, :title=>" Elements"),
2051 #item(:textfield, :width=>80, :tag=>"secondary_elements"),
2052 #item(:text, :title=>"Basis set"),
2053 #item(:popup, :subitems=>bset_desc, :tag=>"secondary_basis"),
2058 #item(:checkbox, :title=>"Calculate electrostatic potential (ESP)", :tag=>"esp"),
2064 item(:checkbox, :title=>"Run JANPA after Psi4", :tag=>"run_janpa",
2065 :action=>lambda { |it|
2066 flag = (it[:value] != 0)
2067 nbo_desc.each { |nbo| set_attr(nbo, :enabled=>flag) }
2071 item(:checkbox, :title=>"NAO", :tag=>"nao"),
2072 item(:checkbox, :title=>"PNAO", :tag=>"pnao"),
2073 item(:checkbox, :title=>"AHO", :tag=>"aho"),
2074 item(:checkbox, :title=>"LHO", :tag=>"lho"),
2076 item(:checkbox, :title=>"PLHO*", :tag=>"plho"),
2077 item(:checkbox, :title=>"LPO", :tag=>"lpo"),
2078 item(:checkbox, :title=>"CLPO", :tag=>"clpo"),
2081 item(:text, :title=>"* Not JANPA original; Molby extension", :font=>[9]),
2087 item(:checkbox, :title=>"Execute Psi4 on this machine", :tag=>"execute_local",
2088 :action=>lambda { |it|
2089 flag = (it[:value] != 0)
2090 set_attr("psi4conda_folder", :enabled=>flag)
2091 set_attr("select_folder", :enabled=>flag)
2092 set_attr("ncpus", :enabled=>flag)
2096 item(:text, :title=>" Psi4 Folder"),
2097 item(:textfield, :width=>300, :tag=>"psi4conda_folder"),
2101 item(:button, :title=>"Select Folder...", :tag=>"select_folder", :action=>:select_psi4_folder),
2103 # item(:button, :title=>"Optional Scripts...", :action=>:set_optional_scripts),
2105 item(:text, :title=>" N of CPUs"),
2106 item(:textfield, :width=>80, :tag=>"ncpus"),
2112 item(:button, :title=>"Edit Psi4 Input and Go", :action=>lambda { |it|
2115 if (tag = it2[:tag]) != nil
2116 h[tag] = it2[:value]
2119 h["basis"] = bset_desc[h["basis"] || 0]
2120 h["secondary_basis"] = bset_desc[h["secondary_basis"] || 0]
2121 h["dfttype"] = dfttype_desc[h["dfttype"] || 0]
2122 h["runtype"] = runtype_desc[h["runtype"] || 0]
2123 h["scftype"] = scftype_desc[h["scftype"] || 0]
2124 psi4_input_direct = mol.cmd_edit_psi4_input(mol.export_psi4(nil, h))
2125 if psi4_input_direct
2135 values[tag] = (user_input[tag] || get_global_settings("psi4.#{tag}") || defaults[tag])
2137 it[:value] = values[tag]
2141 #set_attr("secondary_elements", :enabled=>(values["use_secondary_basis"] == 1))
2142 #set_attr("secondary_basis", :enabled=>(values["use_secondary_basis"] == 1))
2143 set_attr("dfttype", :enabled=>(values["dft"] == 1))
2144 set_attr("psi4conda_folder", :enabled=>(values["execute_local"] == 1))
2145 set_attr("select_folder", :enabled=>(values["execute_local"] == 1))
2146 set_attr("ncpus", :enabled=>(values["execute_local"] == 1))
2148 # set_attr(nao, :enabled=>(values["include_nbo"] == 1))
2152 hash.each_pair { |key, value|
2153 next if key == :status
2154 set_global_settings("psi4.#{key}", value)
2156 if hash[:status] == 0
2157 # Specify basis by internal keys
2158 hash["basis"] = bset_desc[hash["basis"] || 0]
2159 hash["secondary_basis"] = bset_desc[hash["secondary_basis"] || 0]
2160 hash["dfttype"] = dfttype_desc[hash["dfttype"] || 0]
2161 hash["runtype"] = runtype_desc[hash["runtype"] || 0]
2162 hash["scftype"] = scftype_desc[hash["scftype"] || 0]
2163 basename = (self.path ? File.basename(self.path, ".*") : self.name)
2164 fname = Dialog.save_panel("Export Psi4 input file:", self.dir, basename + ".in", "Psi4 input file (*.in)|*.in|All files|*.*")
2165 return nil if !fname
2166 if psi4_input_direct
2167 File.open(fname, "w") { |fp| fp.print(psi4_input_direct) }
2169 export_psi4(fname, hash)
2171 if hash["execute_local"] == 1
2172 if hash["run_janpa"] == 1
2173 # Check if Java is available
2174 if !Molecule.is_java_available()
2175 if !Molecule.make_java_available()
2180 @hf_type = hash["scftype"]
2181 Molecule.execute_psi4(fname, self)
2187 end # End def create_psi4_input
2193 "PM3" => " PM3 0\n",
2194 "STO3G" => " STO 3\n",
2195 "321G" => " N21 3\n",
2196 "631G" => " N31 6\n" }
2197 $gamess_basis_desc = {
2199 "STO3G" => "STO-3G",
2202 $gamess_basis_keys = ["PM3", "STO3G", "321G", "631G"]
2204 ["631Gd", "631Gdp", "631+Gd", "631++Gdp", "6311Gdp", "6311+Gd", "6311++Gdp", "6311++G2d2p", "6311++G3df3pd", "LanL2DZ"].each { |n|
2205 Molecule.read_gamess_basis_sets("basis_sets/#{n}.txt")
2208 register_menu("QChem\tCreate GAMESS Input...",
2209 :cmd_create_gamess_input, :non_empty)
2210 register_menu("QChem\tCreate Psi4 Input...",
2211 :cmd_create_psi4_input, :non_empty)
2212 register_menu("QChem\tCreate MOPAC6 Input...",
2213 :cmd_create_mopac_input, :non_empty) # mopac6.rb
2214 register_menu("QChem\tCreate MO Cube...",
2215 :cmd_create_cube, lambda { |m| m && m.mo_type } )