OSDN Git Service

Start implementing 'additional exponent' for JANPA output
[molby/Molby.git] / Scripts / gamess.rb
index 065df14..94c74fe 100755 (executable)
@@ -1,3 +1,4 @@
+# coding: utf-8
 #
 #  gamess.rb
 #
@@ -30,6 +31,132 @@ end
 
 class Molecule
 
+  #  Import nbo log.  The argument is either an array of Strings
+  #  or an IO object.
+  def import_nbo_log(lines)
+    if lines.is_a?(IO)
+         getline = lambda { lines.gets }
+       else
+         getline = lambda { lines.shift }
+    end
+    nbo = Hash.new
+       while (ln = getline.call) != nil
+         if ln =~ /done with NBO analysis/
+           break
+         end
+         if ln =~ /NATURAL POPULATIONS:  Natural atomic orbital occupancies/
+           #  List of Natural AOs
+               getline.call
+               getline.call
+               getline.call
+               nbo["nao"] = []
+               while (ln = getline.call) != nil
+                 if ln =~ /^\s*$/
+                   #  Skip a blank line
+                   ln = getline.call
+                       if ln =~ /^\s*$/
+                         #  Double blank lines indicates the end of the table
+                         break
+                       end
+                 end
+                 ln.chomp!
+                 next unless ln =~ /^ *(\d+) *([A-Za-z]+) *(\d+) *([a-z]+) +([A-Za-z]+)\( *(\d+)([a-z]+)\)/
+                 i = $1.to_i - 1
+                 an = $3.to_i - 1
+                 atom_label = $2 + $3.to_s
+                 label = $4
+                 type = $5
+                 pn = $6
+                 vals = $~.post_match.split.map { |e| e.to_f }
+                 #  atom_index, type, pn+label, occupancy[, energy]
+                 nbo["nao"].push([atom_label, type, pn + label] + vals)
+               end
+         elsif ln =~ / NATURAL BOND ORBITALS \(Summary\):/
+           nbo["nbo"] = []
+               getline.call
+           while (ln = getline.call) != nil
+                 break if ln =~ /Total Lewis/
+                 if ln =~ /^\s*$/
+                       #  Skip a blank line
+                       ln = getline.call
+                       if ln =~ /^\s*$/
+                         #  Double blank lines indicates the end of the table
+                         break
+                       end
+                 end
+                 if ln =~ /^\s*(\d+)\. *([A-Za-z]+\*?) *\( *(\d+)\) *([A-Za-z]+) *(\d+)(- *([A-Za-z]+) *(\d+))? *(\d\.\d+)/
+                       idx = $1.to_i
+                       orb_kind = $2
+                       orb_idx = $3
+                       atom_label = $4 + ($5.to_i - 1).to_s
+                       if $6 != nil
+                         atom_label += "-" + $7 + ($8.to_i - 1).to_s
+                       end
+                       occ = $9.to_f
+                       nbo["nbo"].push([orb_kind + "_" + orb_idx, atom_label, occ])
+                 end
+               end
+         elsif ln =~ /([A-Z]+)s in the ([A-Z]+) basis:/ || ln =~ /([A-Z]+) Fock matrix:/
+           #  Read matrix
+           dst = $1
+               src = $2
+               if src == nil
+                 src = "F"   #  Fock
+               end
+               key = src + "/" + dst
+               getline.call
+               getline.call
+               getline.call
+               elements = []
+               labels = []
+               idx = 0
+               block = 0
+               while (ln = getline.call) != nil
+                 if ln =~ /^\s*$/
+                   #  Blank line: end of one block
+                   ln = getline.call
+                       if ln =~ /^\s*$/
+                         #  Double blank lines indicates the end of the table
+                         break
+                       else
+                         #  Begin next section
+                         ln = getline.call
+                         idx = 0
+                         block += 1
+                         next
+                       end
+                 end
+                 ln.chomp!
+                 ln =~ /-?\d+\.\d+/
+                 lab = $~.pre_match.strip
+                 a = ([$~[0]] + $~.post_match.split).map { |e| e.to_f }
+                 if block == 0
+                   lab =~ /^[0-9]+\. *(.*)$/
+                   labels.push($1)
+                 end
+                 (elements[idx] ||= []).concat(a)
+                 idx += 1
+               end
+               nbo[key] = LAMatrix.new(elements).transpose!
+               key2 = (src == "F" ? dst : src) + "_L"
+               if nbo[key2] == nil
+                 nbo[key2] = labels
+               end
+         end
+       end
+       @nbo = nbo
+       if @nbo["AO/NAO"]
+         #  Convert NAO based matrix into AO based matrix
+         @nbo.keys.each { |key|
+           if key[0..3] == "NAO/"
+                 key2 = "AO/" + key[4..-1]
+                 @nbo[key2] = @nbo["AO/NAO"] * @nbo[key]
+               end
+         }
+       end
+       # puts @nbo.inspect
+  end
+  
   def Molecule.read_gamess_basis_sets(fname)
     # $gamess_basis = Hash.new unless $gamess_basis
        $gamess_ecp = Hash.new unless $gamess_ecp
@@ -122,6 +249,9 @@ class Molecule
     gmsname = get_global_settings("gamess.executable_path")
     gmsdir = nil
     gmsvers = nil
+
+       cygwin_version = false  #  For windows: old cygwin version
+       
     while 1
       if gmsname == nil || !File.exist?(gmsname)
         gmsname = Dialog.open_panel("Please locate the GAMESS executable")
@@ -189,6 +319,13 @@ class Molecule
       sep = "/"
     end
 
+       #  Old (obsolete) cygwin version, using ddikick
+    if $platform == "win"
+         if gmsvers == "11"
+           cygwin_version = true
+         end
+       end
+
     #  Get the host name etc.
     hostname = backquote("hostname").chomp
     if $platform == "win"
@@ -208,6 +345,10 @@ class Molecule
     #  Redirect standard output to the log file
     logname = scrdir + sep + logbase
     fpout = File.open(logname, "w")
+       if cygwin_version
+         #  The cygwin version uses LF as the eol character
+         fpout.binmode
+       end
     fpout.print "----- GAMESS execution script -----\n"
     fpout.print "This job is running on host #{hostname}\n"
     fpout.print "under operating system #{uname} at #{Time.now.to_s}\n"
@@ -378,7 +519,7 @@ class Molecule
        #  if ncpus < 2
        #    ncpus = 2
        #  end
-         if gmsvers == "11"
+         if cygwin_version
            #  Old (obsolete) cygwin version, using ddikick
         fpout.print "ddikick will run #{ncpus} compute process\n"
                ENV["CYGWIN"] = "nodosfilewarning"
@@ -403,9 +544,15 @@ class Molecule
     fplog = File.open(logname, "r")
     size = 0
     lines = []
+       lineno = 0
     last_line = ""
     search_mode = 0
-
+       nserch = -1
+       rflag = 0
+       ne_alpha = ne_beta = 0
+       mo_count = 0
+       ncomps = 0
+       
     #  Callback procs
        term_callback = lambda { |m, n|
          msg = "GAMESS execution on #{inpbase} "
@@ -441,6 +588,11 @@ class Molecule
            error_message_box($!.to_s)
                return
          end
+         
+         if (script = get_global_settings("gamess.postfix_script")) != nil && script != ""
+           eval(script)
+         end
+
          if mol != nil
            message_box(msg, hmsg, :ok, icon)
       end
@@ -464,7 +616,6 @@ class Molecule
         size = fplog.tell
         last_i = nil
         i = 0
-        nserch = -1
         while i < lines.count
           line = lines[i]
           if line =~ /GEOMETRY SEARCH POINT NSERCH= *(\d+)/
@@ -472,12 +623,93 @@ class Molecule
             last_i = i
                        search_mode = 1
                        energy = nil
+                 elsif line =~ /START OF DRC CALCULATION/
+                   search_mode = 3
+                       nserch = 1
+                       energy = nil
                  elsif line =~ /CONSTRAINED OPTIMIZATION POINT/
                    search_mode = 2
                        energy = nil
                  elsif line =~ /POINT *([0-9]+) *ON THE REACTION PATH/
                    nserch = $1.to_i
                        last_i = i
+                 elsif line =~ /ATOMIC BASIS SET/
+                       j = i
+                       while j < lines.count
+                         line = lines[j]
+                         break if line =~ /TOTAL NUMBER OF BASIS SET/
+                         j += 1
+                       end
+                       if j < lines.count
+                         #  Found
+                         bs_lines = []
+                         ii = i
+                         while ii <= j
+                           bs_lines.push(lines[ii].chomp)
+                               ii += 1
+                         end
+                         begin
+                           if mol
+                                 ncomps = mol.sub_load_gamess_log_basis_set(bs_lines, lineno + i)
+                               end
+                         rescue
+          $stderr.write($!.to_s + "\n")
+          $stderr.write($!.backtrace.inspect + "\n")
+                         end
+                         last_i = j
+                       else
+                         break  #  Wait until all basis set lines are read
+                       end
+                 elsif line =~ /NUMBER OF OCCUPIED ORBITALS/
+                       line =~ /=\s*(\d+)/
+                       n = Integer($1)
+                       if line =~ /ALPHA/
+                         ne_alpha = n
+                       else
+                         ne_beta = n
+                       end
+                 elsif line =~ /SCFTYP=(\w+)/
+                       scftyp = $1
+                       if ne_alpha > 0 || ne_beta > 0
+                         rflag = 0
+                         case scftyp
+                         when "RHF"
+                               rflag = 1
+                         when "ROHF"
+                               rflag = 2
+                         end
+                       end
+                 elsif line =~ /^\s*(EIGENVECTORS|MOLECULAR ORBITALS)\s*$/
+                       if mo_count == 0 && mol
+                         mol.clear_mo_coefficients
+                         mol.set_mo_info(:type=>["UHF", "RHF", "ROHF"][rflag], :alpha=>ne_alpha, :beta=>ne_beta)
+                       end
+                       i += 2
+                       j = i
+                       mo_count += 1
+                       while j < lines.count
+                         line = lines[j]
+                         break if line =~ /\.\.\.\.\.\./ || line =~ /----------------/
+                         j += 1
+                       end
+                       if j == lines.count
+                         break  #  Wait until complete MO info are read
+                       end
+                       ii = i
+                       mo_lines = []
+                       while ii < j
+                         mo_lines.push(lines[ii].chomp)
+                         ii += 1
+                       end
+                       begin
+                         if mol
+                           mol.sub_load_gamess_log_mo_coefficients(mo_lines, lineno + i, ncomps)
+                         end
+                       rescue
+                         $stderr.write($!.to_s + "\n")
+                         $stderr.write($!.backtrace.inspect + "\n")
+                       end
+                       last_i = j
           elsif line =~ /NSERCH: *([0-9]+)/
           #  print line
                        if mol
@@ -487,7 +719,7 @@ class Molecule
                          line =~ /GRAD\. MAX= +([-.0-9]+)/
                          grad = $1
                          mol.show_text("Search: #{n}\nGradient: #{grad}")
-                         mol.set_property("energy", energy)
+                         mol.set_property("energy", energy)
                        end
                        last_i = i
                  elsif line =~ /TOTAL ENERGY += *([-.0-9]+)/
@@ -504,11 +736,15 @@ class Molecule
           elsif nserch > 0 && line =~ /ddikick\.x/
             last_i = -1
             break
-          elsif mol && nserch > 0 && line =~ /COORDINATES OF ALL ATOMS/
-            #  There should be (natoms + 2) lines
-            if i + mol.natoms + 3 <= lines.count
+          elsif mol && nserch > 0 && (line =~ /COORDINATES OF ALL ATOMS/ || line =~ /COORDINATES \(IN ANGSTROM\) FOR \$DATA GROUP ARE/)
+            #  There should be (natoms) lines
+                       if line =~ /COORDINATES OF ALL ATOMS/
+                         #  Skip header lines
+                         i += 2
+                   end
+            if i + mol.natoms + 1 <= lines.count
               coords = []
-              (i + 3...i + 3 + mol.natoms).each { |j|
+              (i + 1...i + 1 + mol.natoms).each { |j|
                 name, charge, x, y, z = lines[j].split
                 coords.push(Vector3D[x.to_f, y.to_f, z.to_f])
               }
@@ -517,9 +753,24 @@ class Molecule
                            mol.set_property("energy", energy)
                          end
               mol.display
-              last_i = i + mol.natoms + 2
+              last_i = i + mol.natoms
               i = last_i   #  Skip the processed lines
             end
+                 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/
+                   #  NBO output
+                   j = i + 1
+                       while j < lines.count
+                         break if lines[j] =~ /done with NBO analysis/
+                         j += 1
+                       end
+                       if j < lines.count
+                         nbo_lines = lines[i..j]
+                         mol.import_nbo_log(nbo_lines) rescue puts "Error: #{$!}, #{$!.backtrace.inspect}"
+                         last_i = j
+                         i = last_i  #  Skip the processed lines
+                       else
+                         break  #  Wait until NBO is done
+                       end
           end
           i += 1
         end
@@ -528,11 +779,16 @@ class Molecule
           break
         elsif last_i
           lines[0..last_i] = nil
+                 lineno += last_i + 1
         end
       end
       true
     }
 
+    if (script = get_global_settings("gamess.prefix_script")) != nil && script != ""
+         eval(script)
+       end
+
     if $platform == "win"
          if gmsvers == "11"
            hosts = "localhost " * ncpus
@@ -562,7 +818,7 @@ class Molecule
   def export_gamess(fname, hash)
 
     def reorder_array(ary, ordered_sub_ary)
-         return ordered_sub_ary + (ary - ordered_sub_ary)
+         return (ordered_sub_ary & ary) + (ary - ordered_sub_ary)
        end
        
        now = Time.now.to_s
@@ -575,7 +831,7 @@ class Molecule
        #  Various settings
        icharg = hash["charge"]
        mult = hash["mult"]
-       runtyp = ["ENERGY", "PROP", "OPTIMIZE"][hash["runtype"].to_i]
+       runtyp = ["ENERGY", "PROP", "OPTIMIZE", "SADPOINT", "IRC"][hash["runtype"].to_i]
        scftyp = ["RHF", "ROHF", "UHF"][hash["scftype"].to_i]
        bssname = hash["basis"]
        bssname2 = hash["secondary_basis"]
@@ -627,7 +883,7 @@ class Molecule
        h = (hash["CONTRL"] ||= Hash.new)
        h["COORD"] ||= "UNIQUE"
        h["EXETYP"] ||= "RUN"
-       h["ICHARG"] ||= icharg.to_s
+       h["ICHARG"] ||= (icharg || 0).to_s
        h["ICUT"] ||= "20"
        h["INTTYP"] ||= "HONDO"
        h["ITOL"] ||= "30"
@@ -637,10 +893,10 @@ class Molecule
        h["MULT"] ||= mult.to_s
        h["QMTTOL"] ||= "1e-08"
        h["RUNTYP"] ||= runtyp
-       if hash["dft"] != 0
+       if (hash["dft"] || 0) != 0 && hash["dfttype"]
          h["DFTTYP"] ||= hash["dfttype"]
        end
-       if hash["use_internal"] != 0 && (hash["runtype"] == 2 || h["RUNTYP"] == "OPTIMIZE")
+       if (hash["use_internal"] || 0) != 0 && (hash["runtype"] >= 2 || h["RUNTYP"] == "OPTIMIZE" || h["RUNTYP"] == "SADPOINT" || h["RUNTYP"] == "IRC")
          nzvar = natoms * 3 - 6  #  TODO: 3N-5 for linear molecules
          h["NZVAR"] ||= nzvar.to_s
        else
@@ -660,11 +916,28 @@ class Molecule
        
        h = (hash["STATPT"] ||= Hash.new)
        h["NSTEP"] ||= "400"
-       h["OPTTOL"] ||= "1.0E-06"
+    if runtyp == "SADPOINT"
+        h["OPTTOL"] ||= "1.0E-05"
+        h["HESS"] ||= "CALC"
+        h["HSSEND"] ||= ".T."
+    elsif runtyp == "IRC"
+        h["OPTTOL"] ||= "1.0E-05"
+        h["HESS"] ||= "READ"
+        h["HSSEND"] ||= ".T."
+    else
+        h["OPTTOL"] ||= "1.0E-06"
+    end
+    if hash["eliminate_freedom"] == 0
+      h["PROJCT"] ||= ".F."
+    end
 
        h = (hash["SYSTEM"] ||= Hash.new)
        h["MEMDDI"] ||= "0"
-       h["MWORDS"] ||= "16"
+    if runtyp == "SADPOINT" || runtyp == "IRC"
+        h["MWORDS"] ||= "32"
+    else
+        h["MWORDS"] ||= "16"
+    end
        h["TIMLIM"] ||= "50000"
        
        h = (hash["GUESS"] ||= Hash.new)
@@ -681,7 +954,17 @@ class Molecule
          h["AUTO"] ||= ".T."
        end
        
-       if hash["esp"] != 0
+    if runtyp == "IRC"
+        h = (hash["IRC"] ||= Hash.new)
+        h["SADDLE"] = ".T."
+        h["TSENGY"] = ".T."
+        h["FORWRD"] = ".T."
+        h["NPOINT"] = "100"
+        h["STRIDE"] = "0.2"
+        h["OPTTOL"] = "1.0E-05"
+    end
+    
+       if (hash["esp"] || 0) != 0
          h = (hash["ELPOT"] ||= Hash.new)
          h["IEPOT"] ||= "1"
          h["OUTPUT"] ||= "PUNCH"
@@ -691,6 +974,17 @@ class Molecule
          h["PTSEL"] ||= "CONNOLLY"
        end
        
+    if (hash["include_nbo"] || 0) != 0
+      h = (hash["NBO"] ||= Hash.new)
+      s = ""
+      ["nao", "nbo", "nho", "nlmo", "pnao", "pnbo", "pnho", "pnlmo"].each { |nao|
+        if (hash[nao] || 0) != 0
+          s += " f" + nao + " ao" + nao
+        end
+      }
+      h[s] = ""
+    end
+
        if fname
          fp = File.open(fname, "wb")
        else
@@ -715,7 +1009,7 @@ class Molecule
                elsif k == "SCF"
                  ordered = ["CONV", "DIRSCF", "FDIFF", "DAMP"]
                elsif k == "STATPT"
-                 ordered = ["NSTEP", "OPTTOL"]
+                 ordered = ["NSTEP", "OPTTOL", "HESS", "HSSEND"]
                elsif k == "SYSTEM"
                  ordered = ["MEMDDI", "MWORDS", "TIMLIM"]
                elsif k == "GUESS"
@@ -728,6 +1022,8 @@ class Molecule
                  ordered = ["IEPOT", "OUTPUT", "WHERE"]
                elsif k == "PDC"
                  ordered = ["CONSTR", "PTSEL"]
+        elsif k == "IRC"
+          ordered = ["SADDLE", "TSENGY", "FORWRD", "NPOINT", "STRIDE", "OPTTOL"]
                else
                  ordered = []
                end
@@ -735,13 +1031,17 @@ class Molecule
                n = 0
                keys.each_with_index { |kk, i|
                  v = h[kk]
-                 next if v == nil || v == ""
                  if n == 0
                    fp.printf " $%-6s", k
                  elsif n % 3 == 0
                    fp.print "\n        "
                  end
-                 fp.printf " %s=%s", kk, h[kk].to_s
+                 if v == nil || v == ""
+                   #  No value
+                       fp.print " " + kk
+                 else
+                   fp.printf " %s=%s", kk, h[kk].to_s
+                 end
                  n += 1
                }
                if n > 0
@@ -800,21 +1100,73 @@ class Molecule
        end
   end
   
-  def cmd_edit_gamess_input(s)
-    h = Dialog.run("Edit GAMESS Input", "OK", "Cancel", :resizable=>true) {
-         layout(1,
-           item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
-           :flex=>[0,0,0,0,1,1]
-         )
-         set_min_size(300, 300)
-       }
-       if h[:status] == 0
-         return h["edit"]
-       else
-         return nil
-       end
+  def copy_section_from_gamess_output
+      if @gamess_output_fname
+          dir = File.dirname(@gamess_output_fname)
+      else
+          dir = self.dir
+      end
+      fname = Dialog.open_panel("Select GAMESS Output:", dir, "*.log;*.dat")
+      if fname
+          fp = open(fname, "rb")
+          if fp
+              mtime = fp.mtime
+              if @gamess_output_fname == fname && @gamess_output_mtime == mtime
+                  #  Use the cached value
+                  k = @gamess_output_sections
+              else
+                  pos = 0
+                  k = Hash.new
+                  fp.each_line { |ln|
+                      if ln.start_with?(" $")
+                          keyword = ln.strip
+                          if keyword !~ /\$END/
+                              k[keyword] = pos
+                          end
+                      end
+                      pos += ln.length
+                  }
+                  #  Cache the information
+                  @gamess_output_fname = fname
+                  @gamess_output_mtime = mtime
+                  @gamess_output_sections = k
+              end
+              keywords = k.keys.sort { |a, b| k[a] <=> k[b] }
+              h = Dialog.run("Select GAMESS section to copy", "Copy", "Cancel") {
+                  layout(1,
+                         item(:popup, :subitems=>keywords, :tag=>"section"))
+              }
+              if h[:status] == 0
+                  fp.seek(k[keywords[h["section"]]])
+                  s = ""
+                  fp.each_line { |ln|
+                      s += ln
+                      break if ln =~ /\$END/
+                  }
+                  export_to_clipboard(s)
+              end
+              fp.close
+          end
+      end
   end
   
+  def cmd_edit_gamess_input(s)
+      mol = self
+      h = Dialog.run("Edit GAMESS Input", "OK", "Cancel", :resizable=>true) {
+          layout(1,
+                 item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
+                 item(:button, :title=>"Copy Section from GAMESS Output...", :action=>lambda { |it| mol.copy_section_from_gamess_output } ),
+                 :flex=>[0,0,0,0,1,1]
+                 )
+                 set_min_size(300, 300)
+      }
+      if h[:status] == 0
+          return h["edit"]
+      else
+          return nil
+      end
+  end
+
   def cmd_create_gamess_input
 
     mol = self
@@ -830,7 +1182,7 @@ class Molecule
        dft_desc = ["B3LYP"]
        dft_internal = ["B3LYP"]
 
-       defaults = {"scftype"=>0, "runtype"=>0, "charge"=>"0", "mult"=>"1",
+       defaults = {"scftype"=>0, "runtype"=>0, "use_internal"=>1, "eliminate_freedom"=>1, "charge"=>"0", "mult"=>"1",
          "basis"=>4, "use_secondary_basis"=>0, "secondary_elements"=>"",
          "secondary_basis"=>8, "esp"=>0, "ncpus"=>"1"}
 
@@ -867,11 +1219,11 @@ class Molecule
                fp.close
                user_input["charge"] = user_input["CONTRL"]["ICHARG"]
                user_input["mult"] = user_input["CONTRL"]["MULT"]
-               user_input["runtype"] = ((i = ["ENERGY", "PROP", "OPTIMIZE"].find_index(user_input["CONTRL"]["RUNTYP"])) ? i.to_s : nil)
-               user_input["scftype"] = ((i = ["RHF", "ROHF", "UHF"].find_index(user_input["CONTRL"]["SCFTYP"])) ? i.to_s : nil)
+               user_input["runtype"] = ["ENERGY", "PROP", "OPTIMIZE"].find_index(user_input["CONTRL"]["RUNTYP"])
+               user_input["scftype"] = ["RHF", "ROHF", "UHF"].find_index(user_input["CONTRL"]["SCFTYP"])
                dft_type = dft_internal.find_index(user_input["CONTRL"]["DFTTYP"])
                if dft_type
-                 user_input["dfttype"] = dft_type.to_s
+                 user_input["dfttype"] = dft_type
                  user_input["dft"] = 1
                end
                bssname = nil
@@ -903,7 +1255,7 @@ class Molecule
                if bssname
                  user_input["basis"] = $gamess_basis_keys.find_index(bssname).to_s
                end
-           puts user_input.inspect
+         #  puts user_input.inspect
          end
        end
        
@@ -934,14 +1286,31 @@ class Molecule
                  end
                end
          end
+         def set_optional_scripts(item)
+           h = Dialog.run("GAMESS Optional Scripts") {
+                 s_pre = get_global_settings("gamess.prefix_script")
+                 s_post = get_global_settings("gamess.postfix_script")
+                 layout(1,
+                   item(:text, :title=>"Script to run before GAMESS execution:"),
+                       item(:textview, :width=>400, :height=>200, :value=>s_pre, :tag=>"prefix"),
+                   item(:text, :title=>"Script to run after GAMESS execution:"),
+                       item(:textview, :width=>400, :height=>200, :value=>s_post, :tag=>"postfix"))
+               }
+               if h[:status] == 0
+                 set_global_settings("gamess.prefix_script", h["prefix"])
+                 set_global_settings("gamess.postfix_script", h["postfix"])
+               end
+         end
+      nbos = ["nao", "nbo", "nho", "nlmo", "pnao", "pnbo", "pnho", "pnlmo"]
          layout(4,
                item(:text, :title=>"SCF type"),
                item(:popup, :subitems=>["RHF", "ROHF", "UHF"], :tag=>"scftype"),
            item(:text, :title=>"Run type"),
-               item(:popup, :subitems=>["Energy", "Property", "Optimize"], :tag=>"runtype",
-                 :action=>lambda { |it| set_attr("use_internal", :enabled=>(it[:value] == 2)) } ),
-
-               item(:checkbox, :title=>"Use internal coordinates for structure optimization", :tag=>"use_internal"),
+               item(:popup, :subitems=>["Energy", "Property", "Optimize", "Sadpoint", "IRC"], :tag=>"runtype",
+                 :action=>lambda { |it| set_attr("use_internal", :enabled=>(it[:value] >= 2)) } ),
+        item(:checkbox, :title=>"Use internal coordinates for structure optimization", :tag=>"use_internal"),
+        -1, -1, -1,
+               item(:checkbox, :title=>"Eliminate translation and rotational degrees of freedom", :tag=>"eliminate_freedom"),
                -1, -1, -1,
 
                item(:text, :title=>"Charge"),
@@ -985,7 +1354,25 @@ class Molecule
                item(:checkbox, :title=>"Calculate electrostatic potential (ESP)", :tag=>"esp"),
                -1, -1, -1,
        
-               item(:line),
+        item(:line),
+        -1, -1, -1,
+
+        item(:checkbox, :title=>"Include NBO instructions", :tag=>"include_nbo",
+             :action=>lambda { |it|
+               flag = (it[:value] != 0)
+               nbos.each { |nbo| set_attr(nbo, :enabled=>flag) }
+             }),
+        -1, -1, -1,
+        item(:checkbox, :title=>"NAO", :tag=>"nao"),
+        item(:checkbox, :title=>"NBO", :tag=>"nbo"),
+        item(:checkbox, :title=>"NHO", :tag=>"nho"),
+        item(:checkbox, :title=>"NLMO", :tag=>"nlmo"),
+        item(:checkbox, :title=>"PNAO", :tag=>"pnao"),
+        item(:checkbox, :title=>"PNBO", :tag=>"pnbo"),
+        item(:checkbox, :title=>"PNHO", :tag=>"pnho"),
+        item(:checkbox, :title=>"PNLMO", :tag=>"pnlmo"),
+             
+        item(:line),
                -1, -1, -1,
                
                item(:checkbox, :title=>"Execute GAMESS on this machine", :tag=>"execute_local",
@@ -1003,7 +1390,8 @@ class Molecule
                
                -1,
                item(:button, :title=>"Select Path...", :tag=>"select_path", :action=>:select_gamess_path),
-               -1, -1,
+               -1,
+               item(:button, :title=>"Optional Scripts...", :action=>:set_optional_scripts),
                
                item(:text, :title=>"   N of CPUs"),
                item(:textfield, :width=>80, :tag=>"ncpus"),
@@ -1043,10 +1431,13 @@ class Molecule
          set_attr("secondary_elements", :enabled=>(values["use_secondary_basis"] == 1))
          set_attr("secondary_basis", :enabled=>(values["use_secondary_basis"] == 1))
          set_attr("dfttype", :enabled=>(values["dft"] == 1))
-         set_attr("use_internal", :enabled=>(values["runtype"] == 2))
+         set_attr("use_internal", :enabled=>(values["runtype"] >= 2))
          set_attr("executable_path", :enabled=>(values["execute_local"] == 1))
          set_attr("select_path", :enabled=>(values["execute_local"] == 1))
          set_attr("ncpus", :enabled=>(values["execute_local"] == 1))
+      nbos.each { |nao|
+        set_attr(nao, :enabled=>(values["include_nbo"] == 1))
+      }
        }
        hash.each_pair { |key, value|
          next if key == :status
@@ -1061,7 +1452,7 @@ class Molecule
          fname = Dialog.save_panel("Export GAMESS input file:", self.dir, basename + ".inp", "GAMESS input file (*.inp)|*.inp|All files|*.*")
          return nil if !fname
          if gamess_input_direct
-           puts "gamess_input_direct = \"#{gamess_input_direct}\""
+         #  puts "gamess_input_direct = \"#{gamess_input_direct}\""
            File.open(fname, "w") { |fp| fp.print(gamess_input_direct) }
          else
            export_gamess(fname, hash)
@@ -1149,6 +1540,653 @@ class Molecule
        end
   end
 
+  def export_psi4(fname, hash)
+    now = Time.now.to_s
+    if fname
+      fp = File.open(fname, "wb")
+    else
+      fp = MemoryIO.new
+    end
+    if fp
+      fp.print "#  Psi4 input\n"
+      fp.print "#  Generated by Molby at #{now}\n"
+      fp.print "import sys\n"
+      fp.print "import re\n"
+      fp.print "base = re.sub('\\.\\w*$', '', sys.argv[1])  #  Basename of the input file\n"
+      fp.print "molecule mol {\n"
+      charge = Integer(hash["charge"])
+      mult = Integer(hash["mult"])
+      if charge != 0 || mult != 1
+        fp.printf " %d %d\n", charge, mult
+      end
+      atoms.each { |ap|
+        fp.printf "%-8s %16.12f %16.12f %16.12f\n", ap.element, ap.x, ap.y, ap.z
+      }
+      use_symmetry = Integer(hash["use_symmetry"])
+      move_to_com = Integer(hash["move_to_com"])
+      do_reorient = Integer(hash["do_reorient"])
+      if move_to_com == 0
+        fp.print "nocom\n"
+      end
+      if do_reorient == 0
+        fp.print "noreorient\n"
+      end
+      if use_symmetry == 0
+        fp.print "symmetry c1\n"
+      end
+      fp.print "units angstrom\n"
+      fp.print "}\n"
+      fp.print "set basis #{hash['basis']}\n"
+      fp.print "set reference #{hash['scftype']}\n"
+      options = "return_wfn=True"
+      if hash['dft'] != 0
+        options += ", dft_functional='#{hash['dfttype']}'"
+      end
+      runtype = hash['runtype'].downcase
+      fp.print "energy, wfn = #{runtype}('scf', #{options})\n"
+      fp.print "wfn.write_molden(f'{base}.molden')\n"
+      if hash['run_junpa'] != 0
+        fp.print "\n"
+        fp.print "#  Interface for JANPA\n"
+        fp.print "#  cf. https://sourceforge.net/p/janpa/wiki/psi4Examples/\n"
+        fp.print "d = wfn.Da().to_array()\n"
+        fp.print "s = wfn.S().to_array()\n"
+        fp.print "c = wfn.Ca().to_array()\n"
+        fp.print "occs = c.T.dot(s.dot(d).dot(s).dot(c))\n"
+        fp.print "molden(wfn, f'{base}.da.molden', density_a = psi4.core.Matrix.from_array(occs))\n"
+      end
+    end
+    if fname == nil
+      s = fp.buffer
+      fp.empty
+      return s
+    else
+      fp.close
+      return fname
+    end
+  end
+
+  def Molecule.is_java_available
+    if $platform == "win"
+      f = get_global_settings("java_home")
+      if f
+        ENV["JAVA_HOME"] = f
+        if !ENV["PATH"].split(";").find { |e| e == "#{f}\\bin" }
+          ENV["PATH"] = "#{f}\\bin;" + ENV["PATH"]
+        end
+      end
+    end
+    return (call_subprocess("java -version", nil) == 0)
+  end
+  
+  def Molecule.make_java_available
+    if $platform == "win"
+      fname = Dialog.open_panel("Locate JDK Folder (if you have one):", "c:\\", nil, true)
+      return false if fname == nil
+      fname.sub!(/\//, "\\")
+      if File.exists?("#{fname}\\bin\\java.exe")
+        set_global_settings("java_home", fname)
+        if Molecule.is_java_available()
+          return true
+        end
+      end
+      error_message_box("Cannot run Java. Please examine your installation again.")
+      return false
+    elsif $platform == "mac"
+      message_box("Please download OpenJDK, and move it into /Library/Java/JavaVirtualMachines folder.", "Install Java", :ok)
+      return false
+    else
+      message_box("Please install Java virtual machine.", "Install Java", :ok)
+      return false
+    end
+  end
+  
+  #  Execute JANPA
+  #  inppath is the input file minus extention
+  #  mol is the molecule (may be nil)
+  def Molecule.execute_janpa(inppath, mol, spherical)
+    #    nbo_desc =   #  JANPA
+    janpa_dir = "#{ResourcePath}/JANPA"
+    status = 0
+    outfile = "#{inppath}.janpa.log"
+    if spherical
+      cmd1 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-NormalizeBF", "-i", "#{inppath}.da.molden", "-o", "#{inppath}.in.molden"]
+      cmd2 = nil
+    else
+      cmd1 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-frompsi4v1mo", "-NormalizeBF", "-cart2pure", "-i", "#{inppath}.da.molden", "-o", "#{inppath}.in.molden"]
+      cmd2 = ["java", "-jar", "#{janpa_dir}/molden2molden.jar", "-frompsi4v1mo", "-NormalizeBF", "-cart2pure", "-i", "#{inppath}.molden", "-o", "#{inppath}.spherical.molden"]
+    end
+    cmd3 = ["java", "-jar", "#{janpa_dir}/janpa.jar", "-i", "#{inppath}.in.molden"]
+    ["nao", "pnao", "aho", "lho", "lpo", "clpo"].each { |type|
+      generate = (get_global_settings("psi4.#{type}").to_i != 0)
+      if type == "pnao" || type == "nao" || type == "lho"
+        #  PLHO is generated within Molby from JANPA NAO/PNAO/LHO
+        generate ||= (get_global_settings("psi4.plho").to_i != 0)
+      end
+      if generate
+        cmd3.push("-#{type.upcase}_Molden_File", "#{inppath}.#{type.upcase}.molden")
+      end
+    }
+    # show_progress_panel("Executing JANPA...")
+    procname = "molden2molden"
+    flag = call_subprocess(cmd1, procname, nil, outfile)
+    if flag && cmd2 != nil
+      procname = "molden2molden"
+      flag = call_subprocess(cmd2, procname, nil, ">>#{outfile}")
+    end
+    if flag
+      procname = "janpa"
+      flag = call_subprocess(cmd3, procname, nil, ">>#{outfile}")
+    end
+    if flag
+      if mol
+        #  import JANPA log and molden output
+        #  Files: inppath.janpa.log, inppath.{NAO,PNAO,AHO,LHO,LPO,CLPO,spherical}.molden
+        mol.sub_load_janpa_log(inppath)
+      end
+      hide_progress_panel
+    else
+      status = $?.exitstatus
+      hide_progress_panel
+      message_box("Execution of #{procname} failed with status #{status}.", "JANPA Failed")
+    end
+    return status
+  end
+  
+  #  Execute Psi4
+  #  inpname is the input file
+  #  mol (optional) is the molecule from which the Psi4 input was built.
+  def Molecule.execute_psi4(inpname, mol = nil)
+
+    inpbase = File.basename(inpname)
+    inpbody = File.basename(inpname, ".*")
+    inpdir = File.dirname(inpname)
+
+    #  Set PATH for psi4conda
+    psi4folder = get_global_settings("psi4.psi4conda_folder")
+    ncpus = get_global_settings("psi4.ncpus").to_i
+    orgpath = ENV["PATH"]
+    orgdir = Dir.pwd
+    if $platform == "win"
+      ENV["PATH"] = "#{psi4folder};#{psi4folder}\\Library\\mingw-w64\\bin;#{psi4folder}\\Library\\usr\\bin;#{psi4folder}\\Library\\bin;#{psi4folder}\\Scripts;#{psi4folder}\\bin;#{psi4folder}\\condabin;" + ENV["PATH"]
+    else
+      ENV["PATH"] = "#{psi4folder}/bin:#{psi4folder}/condabin:" + ENV["PATH"]
+    end
+    Dir.chdir(inpdir)
+    cmdargv = ["psi4", "#{inpbase}"]
+    if ncpus > 0
+      cmdargv.push("-n", "#{ncpus}")
+    end
+    hf_type = nil
+    nalpha = nil
+    nbeta = nil
+    spherical = false
+
+    outfile = inpdir + "/" + inpbody + ".out"
+    if File.exists?(outfile)
+      n = 1
+      while true
+        outbackfile = inpdir + "/" + inpbody + "~" + (n == 1 ? "" : "#{n}") + ".out"
+        break if !File.exists?(outbackfile)
+        n += 1
+      end
+      File.rename(outfile, outbackfile)
+    else
+      outbackfile = nil
+    end
+
+    #  Timer callback
+    timer_count = 0
+    fplog = nil
+    last_size = 0
+    lines = []
+    last_line = ""
+    next_index = 0
+    timer_callback = lambda { |m, n|
+      begin
+        timer_count += 1
+        if timer_count < 10   #  Only handle every 1 seconds
+          return true
+        end
+        timer_count = 0
+        if fplog == nil
+          if File.exists?(outfile)
+            fplog = Kernel.open(outfile, "rt")
+            if fplog == nil
+              return true
+            end
+            last_size = 0
+          else
+            return true  #  Skip until outfile is available
+          end
+        end
+        fplog.seek(0, IO::SEEK_END)
+        current_size = fplog.tell
+        if current_size > last_size
+          #  Read new lines
+          fplog.seek(last_size, IO::SEEK_SET)
+          fplog.each_line { |line|
+            if line[-1, 1] == "\n"
+              lines.push(last_line + line)
+              last_line = ""
+            else
+              last_line += line
+              break
+            end
+          }
+          last_size = fplog.tell
+        end
+        li = next_index
+        getline = lambda { if li < lines.length; li += 1; return lines[li - 1]; else return nil; end }
+        vecs = []
+        while (line = getline.call)
+          if line =~ /==> Geometry <==/
+            #  Skip until line containing "------"
+            while line = getline.call
+              break if line =~ /------/
+            end
+            vecs.clear
+            index = 0
+            #  Read atom positions
+            while line = getline.call
+              line.chomp!
+              break if line =~ /^\s*$/
+              tokens = line.split(' ')
+              vecs.push(Vector3D[Float(tokens[1]), Float(tokens[2]), Float(tokens[3])])
+              index += 1
+            end
+            if vecs.length < mol.natoms
+              break  #  Log file is incomplete
+            end
+            #  Does this geometry differ from the last one?
+            vecs.length.times { |i|
+              if (mol.atoms[i].r - vecs[i]).length2 > 1.0e-20
+                #  Create a new frame and break
+                mol.create_frame
+                vecs.length.times { |j|
+                  mol.atoms[j].r = vecs[j]
+                }
+                break
+              end
+            }
+            next_index = li
+            #  end geometry
+          elsif line =~ /Final Energy: +([-.0-9]+)/
+            #  Energy for this geometry
+            energy = Float($1)
+            mol.set_property("energy", energy)
+            mol.show_text("Frame: #{mol.nframes - 1}\nEnergy: #{energy}")
+            print("Frame: #{mol.nframes - 1} Energy: #{energy}\n")
+            if line =~ /RHF/
+              hf_type = "RHF"
+            elsif line =~ /UHF/
+              hf_type = "UHF"
+            elsif line =~ /ROHF/
+              hf_type = "ROHF"
+            end
+            next_index = li
+          elsif line =~ /^ *Nalpha *= *(\d+)/
+            nalpha = Integer($1)
+          elsif line =~ /^ *Nbeta *= *(\d+)/
+            nbeta = Integer($1)
+          elsif line =~ /^ *Spherical Harmonics\?: *(\w+)/
+            spherical = ($1 == "true")
+          end
+        end
+        if next_index > 0
+          lines.slice!(0, next_index)
+          next_index = 0
+        end
+        return true
+      rescue => e
+        #  Some error occurred during callback. Show the message in the console and abort.
+        $stderr.write("#{e.message}\n")
+        $stderr.write("#{e.backtrace.inspect}\n")
+        return false
+      end
+    }
+
+    #  Terminate callback
+    term_callback = lambda { |m, n|
+      do_janpa = false
+      begin
+        msg = "Psi4 execution of #{inpbase} "
+        hmsg = "Psi4 "
+        if n == -1
+          msg += "cannot be started. Please examine Psi4 installation."
+          hmsg += "Error"
+          icon = :error
+        else
+          if n == 0
+            msg += "succeeded."
+            hmsg += "Completed"
+            icon = :info
+          else
+            msg += "failed with status #{n}."
+            hmsg += "Failed"
+            icon = :error
+          end
+          msg += "\n(In directory #{inpdir})"
+        end
+        if n == 0
+          #  Try to load final lines of the logfile
+          timer_count = 100
+          timer_callback.call(m, n)
+          #  Try to load molden file if available
+          mol.clear_basis_set
+          mol.clear_mo_coefficients
+          mol.set_mo_info(:type => hf_type, :alpha => nalpha, :beta => nbeta)
+          molden = inpdir + "/" + inpbody +".molden"
+          if File.exists?(molden)
+            fp = Kernel.open(molden, "rt")
+            mol.instance_eval { @lineno = 0 }
+            begin
+              #  mol.@hf_type should be set before calling sub_load_molden
+              mol.sub_load_molden(fp)
+              fp.close
+            rescue => e
+              $stderr.write("#{e.message}\n")
+              $stderr.write("#{e.backtrace.inspect}\n")
+            end
+          end
+          if (get_global_settings("psi4.run_janpa").to_i == 1)
+            do_janpa = true
+            Molecule.execute_janpa(inpdir + "/" + inpbody, mol, spherical)
+          end
+        elsif n == -1
+          #  The child process actually did not start
+          #  Restore the old file if outbackfile is not nil
+          if outbackfile && !File.exists?(outfile)
+            File.rename(outbackfile, outfile)
+          end
+        end
+        if outbackfile && File.exists?(outbackfile)
+          File.delete(outbackfile)
+        end
+        ENV["PATH"] = orgpath
+        Dir.chdir(orgdir)
+        if mol != nil
+          message_box(msg, hmsg, :ok, icon)
+        end
+      rescue => e
+        $stderr.write("#{e.message}\n")
+        $stderr.write("#{e.backtrace.inspect}\n")
+      end
+      print("% ")
+      true
+    }
+    
+    if mol
+      pid = mol.call_subprocess_async(cmdargv, term_callback, timer_callback)
+      if pid < 0
+        #  This may not happen on OSX or Linux (don't know for MSW)
+        error_message_box("Psi4 failed to start. Please examine Psi4 installation.")
+        return -1
+      end
+    else
+      status = call_subprocess(cmdargv, "Running Psi4")
+      term_callback.call(nil, status)
+      return status
+    end
+  end
+  
+  def cmd_edit_psi4_input(s)
+    mol = self
+    h = Dialog.run("Edit Psi4 Input", "OK", "Cancel", :resizable=>true) {
+      layout(1,
+        item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
+      )
+      set_min_size(300, 300)
+    }
+    if h[:status] == 0
+        return h["edit"]
+    else
+        return nil
+    end
+  end
+
+  def cmd_create_psi4_input
+  
+    mol = self
+  
+    if natoms == 0
+      raise MolbyError, "cannot create Psi4 input; the molecule is empty"
+    end
+    
+    #  Basis sets Cf. https://psicode.org/psi4manual/master/basissets_tables.html#apdx-basistables
+    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"]
+    dfttype_desc = ["B3LYP"]
+    runtype_desc = ["Energy", "Optimize"]
+    scftype_desc = ["RHF", "ROHF", "UHF"]
+    nbo_desc = ["nao", "pnao", "aho", "lho", "lpo", "clpo"]  #  JANPA
+    user_input = Hash.new
+    defaults = {"scftype"=>0, "runtype"=>0, "move_to_com"=>0, "do_reorient"=>0,
+      "use_symmetry"=>0,  "charge"=>"0", "mult"=>"1",
+      "basis"=>1, "use_secondary_basis"=>0, "secondary_elements"=>"",
+      "secondary_basis"=>8, "esp"=>0, "ncpus"=>"1"}
+    psi4_input_direct = nil
+    
+    hash = Dialog.run("Psi4 Export") {
+      def load_basis_set_sub(item)
+        fname = Dialog.open_panel("Select a file containing GAMESS basis set:")
+        if fname
+          Molecule.read_gamess_basis_sets(fname)
+          bset_desc_new = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
+          sel1 = attr("basis", :value)
+          sel2 = attr("secondary_basis", :value)
+          set_attr("basis", :subitems=>bset_desc_new)
+          set_attr("basis", :value=>sel1)
+          set_attr("secondary_basis", :subitems=>bset_desc_new)
+          set_attr("secondary_basis", :value=>sel2)
+        end
+      end
+      def select_psi4_folder(item)
+        #  By default, psi4conda is installed in the home directory
+        while 1
+          fname = Dialog.open_panel("Locate 'psi4conda' Folder:", ENV["HOME"] || ENV["USERPROFILE"], nil, true)
+          return if fname == nil
+          bname = File.basename(fname)
+          if bname == "psi4conda"
+            break
+          else
+            h = Dialog.run("", "OK", "Choose Again") {
+              layout(1,
+                item(:text, :title=>"The folder does not have the name 'psi4conda'."),
+                item(:text, :title=>"Do you want to use it anyway?"))
+            }
+            if h[:status] == 0
+              break
+            end
+          end
+        end
+        set_attr("psi4conda_folder", :value=>fname)
+      end
+      layout(4,
+        #  ------
+        item(:text, :title=>"SCF type"),
+        item(:popup, :subitems=>scftype_desc, :tag=>"scftype"),
+        item(:text, :title=>"Run type"),
+        item(:popup, :subitems=>runtype_desc, :tag=>"runtype"),
+        #  ------
+        item(:checkbox, :title=>"Detect symmetry", :tag=>"use_symmetry"),
+        -1, -1, -1,
+        #  ------
+        item(:checkbox, :title=>"Move the molecule to the center of mass", :tag=>"move_to_com"),
+        -1, -1, -1,
+        #  ------
+        item(:checkbox, :title=>"Rotate the molecule to the symmetry principle axes", :tag=>"do_reorient"),
+        -1, -1, -1,
+        #  ------
+        item(:text, :title=>"Charge"),
+        item(:textfield, :width=>80, :tag=>"charge"),
+        item(:text, :title=>"Multiplicity"),
+        item(:textfield, :width=>80, :tag=>"mult"),
+        #  ------
+        item(:checkbox, :title=>"Use DFT", :tag=>"dft",
+          :action=>lambda { |it| set_attr("dfttype", :enabled=>(it[:value] != 0)) } ),
+        -1,
+        item(:text, :title=>"DFT type"),
+        item(:popup, :subitems=>dfttype_desc, :tag=>"dfttype"),
+        #  ------
+        item(:line),
+        -1, -1, -1,
+        #  ------
+        item(:text, :title=>"Basis set"),
+        item(:popup, :subitems=>bset_desc, :tag=>"basis"),
+        -1,
+        -1,
+        #  ------
+        #item(:button, :title=>"Load Basis Set...", :action=>:load_basis_set_sub),
+        #-1, -1, -1,
+        #  ------
+        #item(:checkbox, :title=>"Use secondary basis set", :tag=>"use_secondary_basis",
+        #  :action=>lambda { |it|
+        #    flag = (it[:value] != 0)
+        #    set_attr("secondary_elements", :enabled=>flag)
+        #    set_attr("secondary_basis", :enabled=>flag)
+        #  }),
+        #-1, -1, -1,
+        #  ------
+        #item(:text, :title=>"   Elements"),
+        #item(:textfield, :width=>80, :tag=>"secondary_elements"),
+        #item(:text, :title=>"Basis set"),
+        #item(:popup, :subitems=>bset_desc, :tag=>"secondary_basis"),
+        #  ------
+        item(:line),
+        -1, -1, -1,
+        #  ------
+        #item(:checkbox, :title=>"Calculate electrostatic potential (ESP)", :tag=>"esp"),
+        #-1, -1, -1,
+        #  ------
+        #item(:line),
+        #-1, -1, -1,
+        #  ------
+        item(:checkbox, :title=>"Run JANPA after Psi4", :tag=>"run_janpa",
+            :action=>lambda { |it|
+              flag = (it[:value] != 0)
+              nbo_desc.each { |nbo| set_attr(nbo, :enabled=>flag) }
+            }),
+        -1, -1, -1,
+        #  ------
+        item(:checkbox, :title=>"NAO", :tag=>"nao"),
+        item(:checkbox, :title=>"PNAO", :tag=>"pnao"),
+        item(:checkbox, :title=>"AHO", :tag=>"aho"),
+        item(:checkbox, :title=>"LHO", :tag=>"lho"),
+        #  ------
+        item(:checkbox, :title=>"PLHO*", :tag=>"plho"),
+        item(:checkbox, :title=>"LPO", :tag=>"lpo"),
+        item(:checkbox, :title=>"CLPO", :tag=>"clpo"),
+        -1,
+        #  ------
+        item(:text, :title=>"* Not JANPA original; Molby extension", :font=>[9]),
+        -1, -1, -1,
+        #  ------
+        item(:line),
+        -1, -1, -1,
+        #  ------
+        item(:checkbox, :title=>"Execute Psi4 on this machine", :tag=>"execute_local",
+          :action=>lambda { |it|
+            flag = (it[:value] != 0)
+            set_attr("psi4conda_folder", :enabled=>flag)
+            set_attr("select_folder", :enabled=>flag)
+            set_attr("ncpus", :enabled=>flag)
+          }),
+        -1, -1, -1,
+        #  ------
+        item(:text, :title=>" Psi4 Folder"),
+        item(:textfield, :width=>300, :tag=>"psi4conda_folder"),
+        -1, -1,
+        #  ------
+        -1,
+        item(:button, :title=>"Select Folder...", :tag=>"select_folder", :action=>:select_psi4_folder),
+        -1, -1,
+        # item(:button, :title=>"Optional Scripts...", :action=>:set_optional_scripts),
+        #  ------
+        item(:text, :title=>"   N of CPUs"),
+        item(:textfield, :width=>80, :tag=>"ncpus"),
+        -1, -1,
+        #  ------
+        item(:line),
+        -1, -1, -1,
+        #  ------
+        item(:button, :title=>"Edit Psi4 Input and Go", :action=>lambda { |it|
+          h = Hash.new
+          each_item { |it2|
+            if (tag = it2[:tag]) != nil
+              h[tag] = it2[:value]
+            end
+          }
+          h["basis"] = bset_desc[h["basis"] || 0]
+          h["secondary_basis"] = bset_desc[h["secondary_basis"] || 0]
+          h["dfttype"] = dfttype_desc[h["dfttype"] || 0]
+          h["runtype"] = runtype_desc[h["runtype"] || 0]
+          h["scftype"] = scftype_desc[h["scftype"] || 0]
+          psi4_input_direct = mol.cmd_edit_psi4_input(mol.export_psi4(nil, h))
+          if psi4_input_direct
+            end_modal(0)
+          end
+        }),
+        -1, -1, -1
+      )
+      values = Hash.new
+      each_item { |it|
+        tag = it[:tag]
+        if tag
+          values[tag] = (user_input[tag] || get_global_settings("psi4.#{tag}") || defaults[tag])
+          if values[tag]
+            it[:value] = values[tag]
+          end
+        end
+      }
+      #set_attr("secondary_elements", :enabled=>(values["use_secondary_basis"] == 1))
+      #set_attr("secondary_basis", :enabled=>(values["use_secondary_basis"] == 1))
+      set_attr("dfttype", :enabled=>(values["dft"] == 1))
+      set_attr("psi4conda_folder", :enabled=>(values["execute_local"] == 1))
+      set_attr("select_folder", :enabled=>(values["execute_local"] == 1))
+      set_attr("ncpus", :enabled=>(values["execute_local"] == 1))
+      #nbos.each { |nao|
+      #  set_attr(nao, :enabled=>(values["include_nbo"] == 1))
+      #}
+    }  #  end Dialog.run
+
+    hash.each_pair { |key, value|
+      next if key == :status
+      set_global_settings("psi4.#{key}", value)
+    }
+    if hash[:status] == 0
+      #  Specify basis by internal keys
+      hash["basis"] = bset_desc[hash["basis"] || 0]
+      hash["secondary_basis"] = bset_desc[hash["secondary_basis"] || 0]
+      hash["dfttype"] = dfttype_desc[hash["dfttype"] || 0]
+      hash["runtype"] = runtype_desc[hash["runtype"] || 0]
+      hash["scftype"] = scftype_desc[hash["scftype"] || 0]
+      basename = (self.path ? File.basename(self.path, ".*") : self.name)
+      fname = Dialog.save_panel("Export Psi4 input file:", self.dir, basename + ".in", "Psi4 input file (*.in)|*.in|All files|*.*")
+      return nil if !fname
+      if psi4_input_direct
+        File.open(fname, "w") { |fp| fp.print(psi4_input_direct) }
+      else
+        export_psi4(fname, hash)
+      end
+      if hash["execute_local"] == 1
+        if hash["run_janpa"] == 1
+          #  Check if Java is available
+          if !Molecule.is_java_available()
+            if !Molecule.make_java_available()
+              return nil
+            end
+          end
+        end
+        @hf_type = hash["scftype"]
+        Molecule.execute_psi4(fname, self)
+      end
+    else
+      nil
+    end
+  
+  end  #  End def create_psi4_input
+  
+
 end
 
 $gamess_basis = {
@@ -1169,6 +2207,8 @@ $gamess_basis_keys = ["PM3", "STO3G", "321G", "631G"]
 
 register_menu("QChem\tCreate GAMESS Input...",
   :cmd_create_gamess_input, :non_empty)
+register_menu("QChem\tCreate Psi4 Input...",
+  :cmd_create_psi4_input, :non_empty)
 register_menu("QChem\tCreate MOPAC6 Input...",
   :cmd_create_mopac_input, :non_empty)    # mopac6.rb
 register_menu("QChem\tCreate MO Cube...",