Scripts/gamess.rb

   1 #
   2 #  gamess.rb
   3 #
   4 #  Created by Toshi Nagata on 2009/11/22.
   5 #  Copyright 2009 Toshi Nagata. All rights reserved.
   6 #
   7 # This program is free software; you can redistribute it and/or modify
   8 # it under the terms of the GNU General Public License as published by
   9 # the Free Software Foundation version 2 of the License.
  10 #
  11 # This program is distributed in the hope that it will be useful,
  12 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14 # GNU General Public License for more details.
  15
  16 class MemoryIO < IO
  17   def initialize
  18     @buffer = ""
  19   end
  20   def write(s)
  21         @buffer << s
  22   end
  23   def buffer
  24     @buffer
  25   end
  26   def empty
  27     @buffer = ""
  28   end
  29 end
  30
  31 class Molecule
  32
  33   #  Import nbo log.  The argument is either an array of Strings
  34   #  or an IO object.
  35   def import_nbo_log(lines)
  36     if lines.is_a?(IO)
  37           getline = lambda { lines.gets }
  38         else
  39           getline = lambda { lines.shift }
  40     end
  41     nbo = Hash.new
  42         while (ln = getline.call) != nil
  43           if ln =~ /done with NBO analysis/
  44             break
  45           end
  46           if ln =~ /NATURAL POPULATIONS:  Natural atomic orbital occupancies/
  47             #  List of Natural AOs
  48                 getline.call
  49                 getline.call
  50                 getline.call
  51                 nbo["nao"] = []
  52                 while (ln = getline.call) != nil
  53                   if ln =~ /^\s*$/
  54                     #  Skip a blank line
  55                     ln = getline.call
  56                         if ln =~ /^\s*$/
  57                           #  Double blank lines indicates the end of the table
  58                           break
  59                         end
  60                   end
  61                   ln.chomp!
  62                   next unless ln =~ /^ *(\d+) *([A-Za-z]+) *(\d+) *([a-z]+) +([A-Za-z]+)\( *(\d+)([a-z]+)\)/
  63                   i = $1.to_i - 1
  64                   an = $3.to_i - 1
  65                   atom_label = $2 + $3.to_s
  66                   label = $4
  67                   type = $5
  68                   pn = $6
  69                   vals = $~.post_match.split.map { |e| e.to_f }
  70                   #  atom_index, type, pn+label, occupancy[, energy]
  71                   nbo["nao"].push([atom_label, type, pn + label] + vals)
  72                 end
  73           elsif ln =~ / NATURAL BOND ORBITALS \(Summary\):/
  74             nbo["nbo"] = []
  75                 getline.call
  76             while (ln = getline.call) != nil
  77                   break if ln =~ /Total Lewis/
  78                   if ln =~ /^\s*$/
  79                         #  Skip a blank line
  80                         ln = getline.call
  81                         if ln =~ /^\s*$/
  82                           #  Double blank lines indicates the end of the table
  83                           break
  84                         end
  85                   end
  86                   if ln =~ /^\s*(\d+)\. *([A-Za-z]+\*?) *\( *(\d+)\) *([A-Za-z]+) *(\d+)(- *([A-Za-z]+) *(\d+))? *(\d\.\d+)/
  87                         idx = $1.to_i
  88                         orb_kind = $2
  89                         orb_idx = $3
  90                         atom_label = $4 + ($5.to_i - 1).to_s
  91                         if $6 != nil
  92                           atom_label += "-" + $7 + ($8.to_i - 1).to_s
  93                         end
  94                         occ = $9.to_f
  95                         nbo["nbo"].push([orb_kind + "_" + orb_idx, atom_label, occ])
  96                   end
  97                 end
  98           elsif ln =~ /([A-Z]+)s in the ([A-Z]+) basis:/ || ln =~ /([A-Z]+) Fock matrix:/
  99             #  Read matrix
 100             dst = $1
 101                 src = $2
 102                 if src == nil
 103                   src = "F"   #  Fock
 104                 end
 105                 key = src + "/" + dst
 106                 getline.call
 107                 getline.call
 108                 getline.call
 109                 elements = []
 110                 labels = []
 111                 idx = 0
 112                 block = 0
 113                 while (ln = getline.call) != nil
 114                   if ln =~ /^\s*$/
 115                     #  Blank line: end of one block
 116                     ln = getline.call
 117                         if ln =~ /^\s*$/
 118                           #  Double blank lines indicates the end of the table
 119                           break
 120                         else
 121                           #  Begin next section
 122                           ln = getline.call
 123                           idx = 0
 124                           block += 1
 125                           next
 126                         end
 127                   end
 128                   ln.chomp!
 129                   ln =~ /-?\d+\.\d+/
 130                   lab = $~.pre_match.strip
 131                   a = ([$~[0]] + $~.post_match.split).map { |e| e.to_f }
 132                   if block == 0
 133                     lab =~ /^[0-9]+\. *(.*)$/
 134                     labels.push($1)
 135                   end
 136                   (elements[idx] ||= []).concat(a)
 137                   idx += 1
 138                 end
 139                 nbo[key] = LAMatrix.new(elements).transpose!
 140                 key2 = (src == "F" ? dst : src) + "_L"
 141                 if nbo[key2] == nil
 142                   nbo[key2] = labels
 143                 end
 144           end
 145         end
 146         @nbo = nbo
 147         if @nbo["AO/NAO"]
 148           #  Convert NAO based matrix into AO based matrix
 149           @nbo.keys.each { |key|
 150             if key[0..3] == "NAO/"
 151                   key2 = "AO/" + key[4..-1]
 152                   @nbo[key2] = @nbo["AO/NAO"] * @nbo[key]
 153                 end
 154           }
 155         end
 156         # puts @nbo.inspect
 157   end
 158
 159   def Molecule.read_gamess_basis_sets(fname)
 160     # $gamess_basis = Hash.new unless $gamess_basis
 161         $gamess_ecp = Hash.new unless $gamess_ecp
 162         # $gamess_basis_desc = Hash.new unless $gamess_basis_desc
 163         # $gamess_basis_keys = [] unless $gamess_basis_keys
 164         basename = File.basename(fname, ".*")
 165         keys = []
 166         descname = nil
 167     File.open(fname, "r") { |fp|
 168           while (s = fp.gets)
 169             if s =~ /^\s*!\s*/ && descname == nil
 170               #  Get the descriptive name from the first comment line
 171                   s = Regexp.last_match.post_match
 172                   if s =~ /  EMSL/
 173                     descname = Regexp.last_match.pre_match
 174                   else
 175                     descname = (s.split)[0]
 176                   end
 177                   $gamess_basis_desc[basename] = descname
 178                   next
 179                 end
 180                 ss, bas = (s.split)[0..1]     #  Tokens delimited by whitespaces
 181                 next if ss == nil || ss == ""
 182                 if ss == "$ECP"
 183                 #  Read the ECP section
 184                   ecp = $gamess_ecp[basename]
 185                   if !ecp
 186                     ecp = []
 187                         $gamess_ecp[basename] = ecp
 188                   end
 189                   keys << basename unless keys.include?(basename)
 190                   while (s = fp.gets)
 191                     break if s=~ /\$END/
 192                         ary = s.split  #  (PNAME, PTYPE, IZCORE, LMAX+1)
 193                         elem = ary[0].match(/\w{1,2}/).to_a[0]  #  ary[0] should begin with an element name
 194                         ap = Parameter.builtin.elements.find { |p| p.name.casecmp(elem) == 0 }
 195                         raise MolbyError, "the ECP definition does not begin with an element name: #{s}" if !ap
 196                         ecpdef = s
 197                         ln = 1
 198                         (0..Integer(ary[3])).each {
 199                           s = fp.gets
 200                           raise MolbyError, "the ECP definition ends unexpectedly at line #{ln} for: #{s}" if !s
 201                           ln += 1
 202                           ary = s.split  #  (NGPOT, comment)
 203                           ecpdef += s
 204                           (1..Integer(ary[0])).each {
 205                             s = fp.gets
 206                                 raise MolbyError, "the ECP definition ends unexpectedly at line #{ln} for: #{s}" if !s
 207                                 ln += 1
 208                                 ecpdef += s
 209                           }
 210                         }
 211                     ecp[ap.index] = ecpdef
 212                   end
 213                 elsif ss =~ /\W/
 214                   #  Comments or other unrecognizable lines
 215                   next
 216                 elsif (ap = Parameter.builtin.elements.find { |p| p.name.casecmp(ss) == 0 || p.fullname.casecmp(ss) == 0 })
 217                   #  Valid basis definition
 218                   if bas == nil || bas =~ /\W/
 219                     bas = basename
 220                   end
 221                   basis = $gamess_basis[bas]
 222                   if !basis
 223                     basis = []
 224                         $gamess_basis[bas] = basis
 225                   end
 226                   keys << bas unless keys.include?(bas)
 227                   basdef = ""
 228                   while (s = fp.gets) && s =~ /\S/
 229                     basdef += s
 230                   end
 231                   basis[ap.index] = basdef
 232                 else
 233                   raise MolbyError, "ss is not a valid element symbol or name: #{s}"
 234                 end
 235           end
 236     }
 237         unless $gamess_basis_keys.include?(basename)
 238           $gamess_basis_keys.push(basename)
 239         end
 240   end
 241
 242   #  Execute GAMESS (inferior copy of rungms script)
 243   #  inpname is the input file
 244   #  mol (optional) is the molecule from which the GAMESS input was built.
 245   #  If mol is specified and RUNTYP=OPTIMIZE, then the intermediate structures are
 246   #  displayed real-time.
 247   def Molecule.execute_gamess(inpname, mol = nil)
 248     gmsname = get_global_settings("gamess.executable_path")
 249     gmsdir = nil
 250     gmsvers = nil
 251
 252         cygwin_version = false  #  For windows: old cygwin version
 253
 254     while 1
 255       if gmsname == nil || !File.exist?(gmsname)
 256         gmsname = Dialog.open_panel("Please locate the GAMESS executable")
 257         exit if gmsname == nil
 258       end
 259       gmsbase = File.basename(gmsname)
 260       gmsdir = File.dirname(gmsname)
 261       if gmsbase =~ /gamess\.(.*)\.(exe|x)$/i
 262         gmsvers = $1
 263         break
 264       else
 265         gmsname = nil
 266         error_message_box(gmsbase + " does not look like a GAMESS executable!")
 267       end
 268     end
 269
 270 #       if mol == nil
 271 #               mol = Molecule.open(inpname)
 272 #               if mol == nil
 273 #                       error_message_box("Cannot open #{inpname} as GAMESS input")
 274 #                       return
 275 #               end
 276 #       end
 277
 278     inpbase = File.basename(inpname)
 279     inpdir = File.dirname(inpname)
 280     inpbody = inpbase.sub(/\.inp$/, "")
 281     logbase = inpbody + ".log"
 282
 283     set_global_settings("gamess.executable_path", gmsname)
 284
 285     ncpus = get_global_settings("gamess.ncpus").to_i
 286         if ncpus == 0
 287           ncpus = 1
 288         end
 289
 290     #  Prepare the scratch directory
 291         begin
 292           scrdir = create_temp_dir("gamess", inpbody)
 293         rescue
 294           error_message_box($!.to_s)
 295           return
 296         end
 297 #    scrdir = $home_directory + "/Molby/gamess"
 298 #       begin
 299 #         mkdir_recursive(scrdir)
 300 #       rescue
 301 #         error_message_box("Cannot create directory #{scrdir}: " + $!.to_s)
 302 #         return
 303 #       end
 304
 305 #    scrdir = scrdir + "/" + inpbody + "." + $$.to_s + ".0"
 306 #    n = 0
 307 #    while File.exist?(scrdir)
 308 #      scrdir = scrdir.sub(".#{n}", ".#{n + 1}")
 309 #      n += 1
 310 #    end
 311 #    Dir.mkdir(scrdir)
 312
 313     if $platform == "win"
 314       sep = "\\"
 315       scrdir.gsub!("/", sep)
 316       gmsdir.gsub!("/", sep)
 317     else
 318       sep = "/"
 319     end
 320
 321         #  Old (obsolete) cygwin version, using ddikick
 322     if $platform == "win"
 323           if gmsvers == "11"
 324             cygwin_version = true
 325           end
 326         end
 327
 328     #  Get the host name etc.
 329     hostname = backquote("hostname").chomp
 330     if $platform == "win"
 331           s = backquote("cmd.exe /c dir \"#{scrdir}\"")
 332       freebytes = s.split("\n").pop.match(/([0-9,]+)[^0-9]*$/).to_a[1]
 333       if freebytes
 334         freebytes = (freebytes.gsub(",","").to_i / 1024).to_s + " Kbytes"
 335       else
 336         freebytes = "(unknown)"
 337       end
 338       uname = backquote("cmd.exe /c ver").to_s.gsub("\n", "")
 339     else
 340       freebytes = `df -k "#{scrdir}"`
 341       uname = `uname`.chomp
 342     end
 343
 344     #  Redirect standard output to the log file
 345     logname = scrdir + sep + logbase
 346     fpout = File.open(logname, "w")
 347         if cygwin_version
 348           #  The cygwin version uses LF as the eol character
 349           fpout.binmode
 350         end
 351     fpout.print "----- GAMESS execution script -----\n"
 352     fpout.print "This job is running on host #{hostname}\n"
 353     fpout.print "under operating system #{uname} at #{Time.now.to_s}\n"
 354     fpout.print "Available scratch disk space (Kbyte units) at beginning of the job is\n"
 355     fpout.print freebytes + "\n"
 356
 357     #  Copy the input file
 358         scrprefix = scrdir + sep + inpbody
 359         if $platform == "win"
 360           scrbody = inpbody
 361         else
 362       scrbody = scrprefix
 363         end
 364     filecopy(inpname, scrprefix + ".F05")
 365     File.open("#{scrdir}/.in_use", "w") { |fp| }
 366
 367     #  Prepare environmental variables
 368     auxdir = "#{gmsdir}#{sep}auxdata"
 369     ENV["ERICFMT"] = "#{auxdir}#{sep}ericfmt.dat"
 370     ENV["MCPPATH"] = "#{auxdir}#{sep}MCP"
 371     ENV["BASPATH"] = "#{auxdir}#{sep}BASES"
 372     ENV["QUANPOL"] = "#{auxdir}#{sep}QUANPOL"
 373     ENV["EXTBAS"] = "/dev/null"
 374     ENV["IRCDATA"] = "#{scrbody}.irc"
 375         ENV["RESTART"] = "#{scrbody}.rst"
 376         ENV["TRAJECT"] = "#{scrbody}.trj"
 377     ENV["PUNCH"] = "#{scrbody}.dat"
 378     ENV["INPUT"] = "#{scrbody}.F05"
 379     ENV["AOINTS"] = "#{scrbody}.F08"
 380     ENV["MOINTS"] = "#{scrbody}.F09"
 381     ENV["DICTNRY"] = "#{scrbody}.F10"
 382     ENV["DRTFILE"] = "#{scrbody}.F11"
 383     ENV["CIVECTR"] = "#{scrbody}.F12"
 384     ENV["CASINTS"] = "#{scrbody}.F13"
 385     ENV["CIINTS"] = "#{scrbody}.F14"
 386     ENV["WORK15"] = "#{scrbody}.F15"
 387     ENV["WORK16"] = "#{scrbody}.F16"
 388     ENV["CSFSAVE"] = "#{scrbody}.F17"
 389     ENV["FOCKDER"] = "#{scrbody}.F18"
 390     ENV["WORK19"] = "#{scrbody}.F19"
 391     ENV["DASORT"] = "#{scrbody}.F20"
 392     ENV["DFTINTS"] = "#{scrbody}.F21"
 393     ENV["DFTGRID"] = "#{scrbody}.F22"
 394     ENV["JKFILE"] = "#{scrbody}.F23"
 395     ENV["ORDINT"] = "#{scrbody}.F24"
 396     ENV["EFPIND"] = "#{scrbody}.F25"
 397     ENV["PCMDATA"] = "#{scrbody}.F26"
 398     ENV["PCMINTS"] = "#{scrbody}.F27"
 399     ENV["MLTPL"] = "#{scrbody}.F28"
 400     ENV["MLTPLT"] = "#{scrbody}.F29"
 401     ENV["DAFL30"] = "#{scrbody}.F30"
 402     ENV["SOINTX"] = "#{scrbody}.F31"
 403     ENV["SOINTY"] = "#{scrbody}.F32"
 404     ENV["SOINTZ"] = "#{scrbody}.F33"
 405     ENV["SORESC"] = "#{scrbody}.F34"
 406     ENV["SIMEN"] = "#{scrbody}.simen"
 407     ENV["SIMCOR"] = "#{scrbody}.simcor"
 408     ENV["GCILIST"] = "#{scrbody}.F37"
 409     ENV["HESSIAN"] = "#{scrbody}.F38"
 410     ENV["SOCCDAT"] = "#{scrbody}.F40"
 411     ENV["AABB41"] = "#{scrbody}.F41"
 412     ENV["BBAA42"] = "#{scrbody}.F42"
 413     ENV["BBBB43"] = "#{scrbody}.F43"
 414     ENV["MCQD50"] = "#{scrbody}.F50"
 415     ENV["MCQD51"] = "#{scrbody}.F51"
 416     ENV["MCQD52"] = "#{scrbody}.F52"
 417     ENV["MCQD53"] = "#{scrbody}.F53"
 418     ENV["MCQD54"] = "#{scrbody}.F54"
 419     ENV["MCQD55"] = "#{scrbody}.F55"
 420     ENV["MCQD56"] = "#{scrbody}.F56"
 421     ENV["MCQD57"] = "#{scrbody}.F57"
 422     ENV["MCQD58"] = "#{scrbody}.F58"
 423     ENV["MCQD59"] = "#{scrbody}.F59"
 424     ENV["MCQD60"] = "#{scrbody}.F60"
 425     ENV["MCQD61"] = "#{scrbody}.F61"
 426     ENV["MCQD62"] = "#{scrbody}.F62"
 427     ENV["MCQD63"] = "#{scrbody}.F63"
 428     ENV["MCQD64"] = "#{scrbody}.F64"
 429     ENV["NMRINT1"] = "#{scrbody}.F61"
 430     ENV["NMRINT2"] = "#{scrbody}.F62"
 431     ENV["NMRINT3"] = "#{scrbody}.F63"
 432     ENV["NMRINT4"] = "#{scrbody}.F64"
 433     ENV["NMRINT5"] = "#{scrbody}.F65"
 434     ENV["NMRINT6"] = "#{scrbody}.F66"
 435     ENV["DCPHFH2"] = "#{scrbody}.F67"
 436     ENV["DCPHF21"] = "#{scrbody}.F68"
 437     ENV["GVVPT"] = "#{scrbody}.F69"
 438
 439     #    next files are used only during coupled cluster runs, so let's
 440     #    display the numerous definitions only if they are to be used.
 441         #    Don't set these variables on windows, where the memory for environmental variables
 442         #    is limited.
 443         if $platform != "win"
 444     ENV["CCREST"] = "#{scrbody}.F70"
 445     ENV["CCDIIS"] = "#{scrbody}.F71"
 446     ENV["CCINTS"] = "#{scrbody}.F72"
 447     ENV["CCT1AMP"] = "#{scrbody}.F73"
 448     ENV["CCT2AMP"] = "#{scrbody}.F74"
 449     ENV["CCT3AMP"] = "#{scrbody}.F75"
 450     ENV["CCVM"] = "#{scrbody}.F76"
 451     ENV["CCVE"] = "#{scrbody}.F77"
 452     ENV["EOMSTAR"] = "#{scrbody}.F80"
 453     ENV["EOMVEC1"] = "#{scrbody}.F81"
 454     ENV["EOMVEC2"] = "#{scrbody}.F82"
 455     ENV["EOMHC1"] = "#{scrbody}.F83"
 456     ENV["EOMHC2"] = "#{scrbody}.F84"
 457     ENV["EOMHHHH"] = "#{scrbody}.F85"
 458     ENV["EOMPPPP"] = "#{scrbody}.F86"
 459     ENV["EOMRAMP"] = "#{scrbody}.F87"
 460     ENV["EOMRTMP"] = "#{scrbody}.F88"
 461     ENV["EOMDG12"] = "#{scrbody}.F89"
 462     ENV["MMPP"] = "#{scrbody}.F90"
 463     ENV["MMHPP"] = "#{scrbody}.F91"
 464     ENV["MMCIVEC"] = "#{scrbody}.F92"
 465     ENV["MMCIVC1"] = "#{scrbody}.F93"
 466     ENV["MMCIITR"] = "#{scrbody}.F94"
 467     ENV["MMNEXM"] = "#{scrbody}.F95"
 468     ENV["MMNEXE"] = "#{scrbody}.F96"
 469     ENV["MMNREXM"] = "#{scrbody}.F97"
 470     ENV["MMNREXE"] = "#{scrbody}.F98"
 471         end
 472     #
 473     #     next are for TDHFX code, not used by current GAMESS
 474     #
 475         if false
 476     ENV["OLI201"] = "#{scrbody}.F201"
 477     ENV["OLI202"] = "#{scrbody}.F202"
 478     ENV["OLI203"] = "#{scrbody}.F203"
 479     ENV["OLI204"] = "#{scrbody}.F204"
 480     ENV["OLI205"] = "#{scrbody}.F205"
 481     ENV["OLI206"] = "#{scrbody}.F206"
 482     ENV["OLI207"] = "#{scrbody}.F207"
 483     ENV["OLI208"] = "#{scrbody}.F208"
 484     ENV["OLI209"] = "#{scrbody}.F209"
 485     ENV["OLI210"] = "#{scrbody}.F210"
 486     ENV["OLI211"] = "#{scrbody}.F211"
 487     ENV["OLI212"] = "#{scrbody}.F212"
 488     ENV["OLI213"] = "#{scrbody}.F213"
 489     ENV["OLI214"] = "#{scrbody}.F214"
 490     ENV["OLI215"] = "#{scrbody}.F215"
 491     ENV["OLI216"] = "#{scrbody}.F216"
 492     ENV["OLI217"] = "#{scrbody}.F217"
 493     ENV["OLI218"] = "#{scrbody}.F218"
 494     ENV["OLI219"] = "#{scrbody}.F219"
 495     ENV["OLI220"] = "#{scrbody}.F220"
 496     ENV["OLI221"] = "#{scrbody}.F221"
 497     ENV["OLI222"] = "#{scrbody}.F222"
 498     ENV["OLI223"] = "#{scrbody}.F223"
 499     ENV["OLI224"] = "#{scrbody}.F224"
 500     ENV["OLI225"] = "#{scrbody}.F225"
 501     ENV["OLI226"] = "#{scrbody}.F226"
 502     ENV["OLI227"] = "#{scrbody}.F227"
 503     ENV["OLI228"] = "#{scrbody}.F228"
 504     ENV["OLI229"] = "#{scrbody}.F229"
 505     ENV["OLI230"] = "#{scrbody}.F230"
 506     ENV["OLI231"] = "#{scrbody}.F231"
 507     ENV["OLI232"] = "#{scrbody}.F232"
 508     ENV["OLI233"] = "#{scrbody}.F233"
 509     ENV["OLI234"] = "#{scrbody}.F234"
 510     ENV["OLI235"] = "#{scrbody}.F235"
 511     ENV["OLI236"] = "#{scrbody}.F236"
 512     ENV["OLI237"] = "#{scrbody}.F237"
 513     ENV["OLI238"] = "#{scrbody}.F238"
 514     ENV["OLI239"] = "#{scrbody}.F239"
 515     end
 516
 517     if $platform == "win"
 518         #  if ncpus < 2
 519         #    ncpus = 2
 520         #  end
 521           if cygwin_version
 522             #  Old (obsolete) cygwin version, using ddikick
 523         fpout.print "ddikick will run #{ncpus} compute process\n"
 524                 ENV["CYGWIN"] = "nodosfilewarning"
 525           else
 526         fpout.print "Microsoft MPI will be running GAMESS on 1 node.\n"
 527         fpout.print "The binary kicked off by 'mpiexec' is gamess.#{gmsvers}.exe\n"
 528         fpout.print "MS-MPI will run #{ncpus*2} compute process\n"
 529         #  File containing environmental variables
 530         envfil = "#{scrprefix}.GMS.ENV"
 531         fp = File.open(envfil, "w")
 532         ENV.each { |k, v| fp.print "#{k}=#{v}\n" }
 533         fp.close
 534         #  File containing arguments to mpiexec
 535         procfil = "#{scrprefix}.processes.mpd"
 536         fp = File.open(procfil, "w")
 537         fp.print "-env ENVFIL \"#{envfil}\" -wdir \"#{scrdir}\" -n #{ncpus*2} \"#{gmsdir}#{sep}gamess.#{gmsvers}.exe\"\n"
 538         fp.close
 539           end
 540     end
 541
 542     fpout.close
 543     fplog = File.open(logname, "r")
 544     size = 0
 545     lines = []
 546         lineno = 0
 547     last_line = ""
 548     search_mode = 0
 549         nserch = -1
 550         rflag = 0
 551         ne_alpha = ne_beta = 0
 552         mo_count = 0
 553         ncomps = 0
 554
 555     #  Callback procs
 556         term_callback = lambda { |m, n|
 557           msg = "GAMESS execution on #{inpbase} "
 558           hmsg = "GAMESS "
 559           if n == 0
 560             msg += "succeeded."
 561                 hmsg += "Completed"
 562                 icon = :info
 563           else
 564             msg += "failed with status #{n}."
 565                 hmsg += "Failed"
 566                 icon = :error
 567           end
 568           msg += "\n(In directory #{inpdir})"
 569
 570           File.delete("#{scrdir}/.in_use")
 571
 572           ext_to_keep = [".dat", ".rst", ".trj", ".efp", ".gamma", ".log"]
 573           ext_to_keep.each { |ex|
 574             if File.exists?("#{scrprefix}#{ex}")
 575                   filecopy("#{scrprefix}#{ex}", "#{inpdir}#{sep}#{inpbody}#{ex}")
 576             end
 577           }
 578           Dir.foreach(scrdir) { |file|
 579             if file != "." && file != ".." && !ext_to_keep.include?(File.extname(file))
 580                   File.delete("#{scrdir}#{sep}#{file}")
 581             end
 582           }
 583
 584           begin
 585             erase_old_logs(scrdir, "latest", 5)
 586           rescue
 587             error_message_box($!.to_s)
 588                 return
 589           end
 590
 591           if (script = get_global_settings("gamess.postfix_script")) != nil && script != ""
 592             eval(script)
 593           end
 594
 595           if mol != nil
 596             message_box(msg, hmsg, :ok, icon)
 597       end
 598     }
 599
 600     timer_callback = lambda { |m, n|
 601       fplog.seek(0, IO::SEEK_END)
 602       sizec = fplog.tell
 603       if sizec > size
 604         #  Read new lines
 605         fplog.seek(size, IO::SEEK_SET)
 606         fplog.each_line { |line|
 607           if line[-1, 1] == "\n"
 608             lines.push(last_line + line)
 609             last_line = ""
 610           else
 611             last_line += line
 612             break
 613           end
 614         }
 615         size = fplog.tell
 616         last_i = nil
 617         i = 0
 618         while i < lines.count
 619           line = lines[i]
 620           if line =~ /GEOMETRY SEARCH POINT NSERCH= *(\d+)/
 621             nserch = $1.to_i
 622             last_i = i
 623                         search_mode = 1
 624                         energy = nil
 625                   elsif line =~ /START OF DRC CALCULATION/
 626                     search_mode = 3
 627                         nserch = 1
 628                         energy = nil
 629                   elsif line =~ /CONSTRAINED OPTIMIZATION POINT/
 630                     search_mode = 2
 631                         energy = nil
 632                   elsif line =~ /POINT *([0-9]+) *ON THE REACTION PATH/
 633                     nserch = $1.to_i
 634                         last_i = i
 635                   elsif line =~ /ATOMIC BASIS SET/
 636                         j = i
 637                         while j < lines.count
 638                           line = lines[j]
 639                           break if line =~ /TOTAL NUMBER OF BASIS SET/
 640                           j += 1
 641                         end
 642                         if j < lines.count
 643                           #  Found
 644                           bs_lines = []
 645                           ii = i
 646                           while ii <= j
 647                             bs_lines.push(lines[ii].chomp)
 648                                 ii += 1
 649                           end
 650                           begin
 651                             if mol
 652                                   ncomps = mol.sub_load_gamess_log_basis_set(bs_lines, lineno + i)
 653                                 end
 654                           rescue
 655                             puts $!.to_s
 656                             puts $!.backtrace.inspect
 657                           end
 658                           last_i = j
 659                         else
 660                           break  #  Wait until all basis set lines are read
 661                         end
 662                   elsif line =~ /NUMBER OF OCCUPIED ORBITALS/
 663                         line =~ /=\s*(\d+)/
 664                         n = Integer($1)
 665                         if line =~ /ALPHA/
 666                           ne_alpha = n
 667                         else
 668                           ne_beta = n
 669                         end
 670                   elsif line =~ /SCFTYP=(\w+)/
 671                         scftyp = $1
 672                         if ne_alpha > 0 || ne_beta > 0
 673                           rflag = 0
 674                           case scftyp
 675                           when "RHF"
 676                                 rflag = 1
 677                           when "ROHF"
 678                                 rflag = 2
 679                           end
 680                         end
 681                   elsif line =~ /^\s*(EIGENVECTORS|MOLECULAR ORBITALS)\s*$/
 682                         if mo_count == 0 && mol
 683                           mol.clear_mo_coefficients
 684                           mol.set_mo_info(:type=>["UHF", "RHF", "ROHF"][rflag], :alpha=>ne_alpha, :beta=>ne_beta)
 685                         end
 686                         i += 2
 687                         j = i
 688                         mo_count += 1
 689                         while j < lines.count
 690                           line = lines[j]
 691                           break if line =~ /\.\.\.\.\.\./ || line =~ /----------------/
 692                           j += 1
 693                         end
 694                         if j == lines.count
 695                           break  #  Wait until complete MO info are read
 696                         end
 697                         ii = i
 698                         mo_lines = []
 699                         while ii < j
 700                           mo_lines.push(lines[ii].chomp)
 701                           ii += 1
 702                         end
 703                         begin
 704                           if mol
 705                             mol.sub_load_gamess_log_mo_coefficients(mo_lines, lineno + i, ncomps)
 706                           end
 707                         rescue
 708                           puts $!.to_s
 709                           puts $!.backtrace.inspect
 710                         end
 711                         last_i = j
 712           elsif line =~ /NSERCH: *([0-9]+)/
 713           #  print line
 714                         if mol
 715                           n = $1
 716                           line =~ /E= +([-.0-9]+)/
 717                           energy = $1.to_f
 718                           line =~ /GRAD\. MAX= +([-.0-9]+)/
 719                           grad = $1
 720                           mol.show_text("Search: #{n}\nGradient: #{grad}")
 721                           # mol.set_property("energy", energy)
 722                         end
 723                         last_i = i
 724                   elsif line =~ /TOTAL ENERGY += *([-.0-9]+)/
 725                     energy = $1
 726                         if mol && search_mode == 2
 727                           mol.show_text("Point: #{nserch}\nEnergy: #{energy}")
 728                         end
 729                         energy = energy.to_f
 730                   elsif line =~ /FINAL .* ENERGY IS +([-.0-9]+) *AFTER/
 731                     energy = $1.to_f
 732                         if mol && search_mode == 0
 733                           mol.set_property("energy", energy)
 734                         end
 735           elsif nserch > 0 && line =~ /ddikick\.x/
 736             last_i = -1
 737             break
 738           elsif mol && nserch > 0 && (line =~ /COORDINATES OF ALL ATOMS/ || line =~ /COORDINATES \(IN ANGSTROM\) FOR \$DATA GROUP ARE/)
 739             #  There should be (natoms) lines
 740                         if line =~ /COORDINATES OF ALL ATOMS/
 741                           #  Skip header lines
 742                           i += 2
 743                     end
 744             if i + mol.natoms + 1 <= lines.count
 745               coords = []
 746               (i + 1...i + 1 + mol.natoms).each { |j|
 747                 name, charge, x, y, z = lines[j].split
 748                 coords.push(Vector3D[x.to_f, y.to_f, z.to_f])
 749               }
 750               mol.create_frame([coords])
 751                           if search_mode == 1 && energy
 752                             mol.set_property("energy", energy)
 753                           end
 754               mol.display
 755               last_i = i + mol.natoms
 756               i = last_i   #  Skip the processed lines
 757             end
 758                   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/
 759                     #  NBO output
 760                     j = i + 1
 761                         while j < lines.count
 762                           break if lines[j] =~ /done with NBO analysis/
 763                           j += 1
 764                         end
 765                         if j < lines.count
 766                           nbo_lines = lines[i..j]
 767                           mol.import_nbo_log(nbo_lines) rescue puts "Error: #{$!}, #{$!.backtrace.inspect}"
 768                           last_i = j
 769                           i = last_i  #  Skip the processed lines
 770                         else
 771                           break  #  Wait until NBO is done
 772                         end
 773           end
 774           i += 1
 775         end
 776         if last_i == -1
 777           lines.clear
 778           break
 779         elsif last_i
 780           lines[0..last_i] = nil
 781                   lineno += last_i + 1
 782         end
 783       end
 784       true
 785     }
 786
 787     if (script = get_global_settings("gamess.prefix_script")) != nil && script != ""
 788           eval(script)
 789         end
 790
 791     if $platform == "win"
 792           if gmsvers == "11"
 793             hosts = "localhost " * ncpus
 794             cmdline = "cmd.exe /c \"#{gmsdir}/ddikick.exe #{gmsdir}/gamess.#{gmsvers}.exe #{inpbody} -ddi #{ncpus} #{ncpus} #{hosts} -scr #{scrdir} < NUL >>#{logname}"
 795           else
 796             cmdline = "cmd.exe /c \"mpiexec -configfile #{procfil} >>#{logname}"
 797           end
 798         else
 799           hosts = "localhost " * ncpus
 800           cmdline = "/bin/sh -c '#{gmsdir}/ddikick.x #{gmsdir}/gamess.#{gmsvers}.x #{inpbody} -ddi #{ncpus} #{ncpus} #{hosts} -scr #{scrdir} < /dev/null >>#{logname}'"
 801         end
 802
 803         if mol
 804           pid = mol.call_subprocess_async(cmdline, term_callback, timer_callback)
 805           if pid < 0
 806             error_message_box("GAMESS failed to start. Please examine GAMESS installation.")
 807             return
 808           end
 809         else
 810           pid = call_subprocess(cmdline, "Running GAMESS")
 811           term_callback.call(nil, pid)
 812           return pid
 813         end
 814
 815   end
 816
 817   def export_gamess(fname, hash)
 818
 819     def reorder_array(ary, ordered_sub_ary)
 820           return (ordered_sub_ary & ary) + (ary - ordered_sub_ary)
 821         end
 822
 823         now = Time.now.to_s
 824         if fname
 825           basename = File.basename(fname, ".*")
 826         else
 827           basename = File.basename(self.name, ".*")
 828         end
 829
 830         #  Various settings
 831         icharg = hash["charge"]
 832         mult = hash["mult"]
 833         runtyp = ["ENERGY", "PROP", "OPTIMIZE"][hash["runtype"].to_i]
 834         scftyp = ["RHF", "ROHF", "UHF"][hash["scftype"].to_i]
 835         bssname = hash["basis"]
 836         bssname2 = hash["secondary_basis"]
 837         if hash["use_secondary_basis"].to_i != 0 && bssname2 != bssname
 838           use_2nd = true
 839           element2 = hash["secondary_elements"].split(/[\s,]+/).map { |name| name.capitalize }
 840         else
 841           use_2nd = false
 842           bssname2 = bssname
 843         end
 844         basis = $gamess_basis[bssname]
 845         basis2 = $gamess_basis[bssname2]
 846
 847         #  Use effective core potentials?
 848         ecp = $gamess_ecp[bssname]
 849         ecp2 = $gamess_ecp[bssname2]
 850         ecp_read = (ecp || ecp2 ? "ECP=READ " : "")
 851
 852         #  Use only one built-in basis set?
 853         gbasis = nil
 854         if !use_2nd
 855           case bssname
 856           when "PM3"
 857                 gbasis = "PM3 NGAUSS=3"
 858           when "STO3G"
 859                 gbasis = "STO NGAUSS=3"
 860           when "321G"
 861                 gbasis = "N21 NGAUSS=3"
 862           when "631G"
 863                 gbasis = "N31 NGAUSS=6"
 864           when "631Gd"
 865                 gbasis = "N31 NGAUSS=6 NDFUNC=1"
 866           when "631Gdp"
 867                 gbasis = "N31 NGAUSS=6 NDFUNC=1 NPFUNC=1"
 868           when "6311G"
 869                 gbasis = "N311 NGAUSS=6"
 870           when "6311Gdp"
 871                 gbasis = "N311 NGAUSS=6 NDFUNC=1 NPFUNC=1"
 872           end
 873         end
 874
 875         #  Count non-dummy atoms
 876         natoms = 0
 877         each_atom { |ap|
 878           natoms += 1 if ap.atomic_number != 0
 879         }
 880
 881         #  Fill hash with default values
 882         h = (hash["CONTRL"] ||= Hash.new)
 883         h["COORD"] ||= "UNIQUE"
 884         h["EXETYP"] ||= "RUN"
 885         h["ICHARG"] ||= (icharg || 0).to_s
 886         h["ICUT"] ||= "20"
 887         h["INTTYP"] ||= "HONDO"
 888         h["ITOL"] ||= "30"
 889         h["MAXIT"] ||= "200"
 890         h["MOLPLT"] ||= ".T."
 891         h["MPLEVL"] ||= "0"
 892         h["MULT"] ||= mult.to_s
 893         h["QMTTOL"] ||= "1e-08"
 894         h["RUNTYP"] ||= runtyp
 895         if (hash["dft"] || 0) != 0 && hash["dfttype"]
 896           h["DFTTYP"] ||= hash["dfttype"]
 897         end
 898         if (hash["use_internal"] || 0) != 0 && (hash["runtype"] == 2 || h["RUNTYP"] == "OPTIMIZE")
 899           nzvar = natoms * 3 - 6  #  TODO: 3N-5 for linear molecules
 900           h["NZVAR"] ||= nzvar.to_s
 901         else
 902           nzvar = 0
 903         end
 904         h["SCFTYP"] ||= scftyp
 905         if ecp_read != ""
 906           h["ECP"] ||= "READ"
 907         end
 908         h["UNITS"] ||= "ANGS"
 909
 910         h = (hash["SCF"] ||= Hash.new)
 911         h["CONV"] ||= "1.0E-06"
 912         h["DIRSCF"] ||= ".T."
 913         h["FDIFF"] ||= ".T."
 914         h["DAMP"] ||= ".T."
 915
 916         h = (hash["STATPT"] ||= Hash.new)
 917         h["NSTEP"] ||= "400"
 918         h["OPTTOL"] ||= "1.0E-06"
 919
 920         h = (hash["SYSTEM"] ||= Hash.new)
 921         h["MEMDDI"] ||= "0"
 922         h["MWORDS"] ||= "16"
 923         h["TIMLIM"] ||= "50000"
 924
 925         h = (hash["GUESS"] ||= Hash.new)
 926         h["GUESS"] ||= "HUCKEL"
 927
 928         if gbasis
 929           h = (hash["BASIS"] ||= Hash.new)
 930           h["GBASIS"] ||= gbasis
 931         end
 932
 933         if nzvar > 0
 934           h = (hash["ZMAT"] ||= Hash.new)
 935           h["DLC"] ||= ".T."
 936           h["AUTO"] ||= ".T."
 937         end
 938
 939         if (hash["esp"] || 0) != 0
 940           h = (hash["ELPOT"] ||= Hash.new)
 941           h["IEPOT"] ||= "1"
 942           h["OUTPUT"] ||= "PUNCH"
 943           h["WHERE"] ||= "PDC"
 944           h = (hash["PDC"] ||= Hash.new)
 945           h["CONSTR"] ||= "NONE"
 946           h["PTSEL"] ||= "CONNOLLY"
 947         end
 948
 949         if fname
 950           fp = File.open(fname, "wb")
 951         else
 952           fp = MemoryIO.new
 953         end
 954     if fp
 955           fp.print "!  GAMESS input\n"
 956           fp.print "!  Generated by Molby at #{now}\n"
 957           fp.print "!  Basis set: " + ($gamess_basis_desc[bssname] || "(not specified)") + "\n"
 958           if use_2nd
 959             fp.print "!  [" + element2.join(", ") + "]: " + ($gamess_basis_desc[bssname2] || "(not specified)") + "\n"
 960           end
 961           ordered = ["CONTRL", "SCF", "STATPT", "SYSTEM", "GUESS", "BASIS", "ZMAT", "ELPOT", "PDC"]
 962           controls = reorder_array(hash.keys.select { |k| k != "key_order" && hash[k].is_a?(Hash) },
 963                 hash["key_order"] || ordered)
 964           controls.each { |k|
 965             h = hash[k]
 966                 next if h == nil || h.size == 0 || (h["key_order"] && h.size == 1)
 967                 if k == "CONTRL"
 968                   ordered = ["COORD", "EXETYP", "ICHARG", "ICUT", "INTTYP", "ITOL", "MAXIT", "MOLPLT", "MPLEVL",
 969                              "MULT", "QMTTOL", "RUNTYP", "SCFTYP", "ECP", "UNITS", "DFTTYP", "NZVAR"]
 970                 elsif k == "SCF"
 971                   ordered = ["CONV", "DIRSCF", "FDIFF", "DAMP"]
 972                 elsif k == "STATPT"
 973                   ordered = ["NSTEP", "OPTTOL"]
 974                 elsif k == "SYSTEM"
 975                   ordered = ["MEMDDI", "MWORDS", "TIMLIM"]
 976                 elsif k == "GUESS"
 977                   ordered = ["GUESS"]
 978                 elsif k == "BASIS"
 979                   ordered = ["GBASIS"]
 980                 elsif k == "ZMAT"
 981                   ordered = ["DLC", "AUTO"]
 982                 elsif k == "ELPOT"
 983                   ordered = ["IEPOT", "OUTPUT", "WHERE"]
 984                 elsif k == "PDC"
 985                   ordered = ["CONSTR", "PTSEL"]
 986                 else
 987                   ordered = []
 988                 end
 989             keys = reorder_array(h.keys, h["key_order"] || ordered)
 990                 n = 0
 991                 keys.each_with_index { |kk, i|
 992                   v = h[kk]
 993                   if n == 0
 994                     fp.printf " $%-6s", k
 995                   elsif n % 3 == 0
 996                     fp.print "\n        "
 997                   end
 998                   if v == nil || v == ""
 999                     #  No value
1000                         fp.print " " + kk
1001                   else
1002                     fp.printf " %s=%s", kk, h[kk].to_s
1003                   end
1004                   n += 1
1005                 }
1006                 if n > 0
1007                   fp.print " $END\n"
1008                 end
1009           }
1010           fp.print " $DATA\n#{basename}\nC1 0\n"
1011           secondary = []
1012           each_atom { |ap|
1013             next if ap.atomic_number == 0
1014                 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
1015                 if use_2nd && element2.include?(ap.element)
1016                   secondary[ap.index] = true
1017                 end
1018             if !gbasis
1019                   #  Basis specification followed by a blank line
1020                   bas = (secondary[ap.index] ? basis2 : basis)
1021                   if bas.is_a?(Array)
1022                     bas = bas[ap.atomic_number]
1023                   end
1024                   if bas
1025                     fp.print bas
1026                     fp.print "\n"
1027                   else
1028                     puts "Warning: basis set is not defined for atom #{ap.index}, element #{ap.element}"
1029                   end
1030                 end
1031           }
1032           fp.print " $END\n"
1033           ecp_ary = []
1034           if ecp || ecp2
1035             fp.print " $ECP\n"
1036                 each_atom { |ap|
1037                   an = ap.atomic_number
1038                   next if an == 0
1039                   ecpp = (secondary[ap.index] ? ecp2 : ecp)
1040                   e = ecp_ary[an] || (ecpp && ecpp[an])
1041                   if e
1042                     #  Cache the PNAME of the $ECP entry and re-use it
1043                     ecp_ary[an] ||= (e.split)[0] + "\n"
1044                   else
1045                     e = ap.element.upcase + "-ECP NONE\n"
1046                   end
1047                   fp.print e
1048                 }
1049             fp.print " $END\n"
1050           end
1051         end
1052         if fname == nil
1053           s = fp.buffer
1054       fp.empty
1055           return s
1056         else
1057           fp.close
1058           return fname
1059         end
1060   end
1061
1062   def cmd_edit_gamess_input(s)
1063     h = Dialog.run("Edit GAMESS Input", "OK", "Cancel", :resizable=>true) {
1064           layout(1,
1065             item(:textview, :value=>s, :tag=>"edit", :width=>400, :height=>400, :flex=>[0,0,0,0,1,1]),
1066             :flex=>[0,0,0,0,1,1]
1067           )
1068           set_min_size(300, 300)
1069         }
1070         if h[:status] == 0
1071           return h["edit"]
1072         else
1073           return nil
1074         end
1075   end
1076
1077   def cmd_create_gamess_input
1078
1079     mol = self
1080
1081     if natoms == 0
1082       raise MolbyError, "cannot create GAMESS input; the molecule is empty"
1083     end
1084
1085         #  Descriptive text and internal string for popup menus
1086 #    bset_desc = ["PM3", "STO-3G", "3-21G", "6-31G", "6-31G(d)", "6-31G(d,p)", "6-311G", "6-311G(d,p)", "LanL2DZ"]
1087 #       bset_internal = ["PM3", "STO3G", "321G", "631G", "631Gd", "631Gdp", "6311G", "6311Gdp", "LanL2DZ"]
1088     bset_desc = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
1089         dft_desc = ["B3LYP"]
1090         dft_internal = ["B3LYP"]
1091
1092         defaults = {"scftype"=>0, "runtype"=>0, "charge"=>"0", "mult"=>"1",
1093           "basis"=>4, "use_secondary_basis"=>0, "secondary_elements"=>"",
1094           "secondary_basis"=>8, "esp"=>0, "ncpus"=>"1"}
1095
1096         gamess_input_direct = nil
1097
1098     user_input = Hash.new
1099         ["CONTRL", "SCF", "STATPT", "SYSTEM", "GUESS", "BASIS"].each { |k|
1100           user_input[k] = Hash.new
1101         }
1102     if ".inp".casecmp(self.name[-4, 4]) == 0
1103           #  Read input and analyze commands
1104           fp = open(self.path, "r")
1105           key = hh = nil
1106           if fp
1107             while ln = fp.gets
1108                   ln.strip!
1109                   break if "$DATA".casecmp(ln[0, 5]) == 0
1110                   ln.split.each { |s|
1111                     if s[0] == "$"
1112                           #  New key
1113                           key = s[1..-1].upcase
1114                           hh = user_input[key] = Hash.new
1115                           (user_input["key_order"] ||= []).push(key)
1116                         else
1117                           k, v = s.split("=")
1118                           if key && hh
1119                             k.upcase!
1120                             hh[k] = v
1121                                 (hh["key_order"] ||= []).push(k)
1122                           end
1123                         end
1124                   }
1125                 end  #  end while
1126                 fp.close
1127                 user_input["charge"] = user_input["CONTRL"]["ICHARG"]
1128                 user_input["mult"] = user_input["CONTRL"]["MULT"]
1129                 user_input["runtype"] = ((i = ["ENERGY", "PROP", "OPTIMIZE"].find_index(user_input["CONTRL"]["RUNTYP"])) ? i.to_s : nil)
1130                 user_input["scftype"] = ((i = ["RHF", "ROHF", "UHF"].find_index(user_input["CONTRL"]["SCFTYP"])) ? i.to_s : nil)
1131                 dft_type = dft_internal.find_index(user_input["CONTRL"]["DFTTYP"])
1132                 if dft_type
1133                   user_input["dfttype"] = dft_type.to_s
1134                   user_input["dft"] = 1
1135                 end
1136                 bssname = nil
1137                 user_input["basis"] = "-1"
1138                 case user_input["BASIS"]["GBASIS"]
1139                 when "PM3"
1140                   bssname = "PM3"
1141                 when "STO"
1142                   bssname = "STO3G"
1143                 when "N21"
1144                   bssname = "321G"
1145                 when "N31"
1146                   if user_input["NDFUNC"] == "1"
1147                     if user_input["NPFUNC"] == "1"
1148                           bssname = "631Gdp"
1149                         else
1150                           bssname = "631Gd"
1151                         end
1152                   else
1153                     bssname = "631G"
1154                   end
1155                 when "N311"
1156                   if user_input["NDFUNC"] == "1" && user_input["NPFUNC"] == "1"
1157                     bssname = "6311Gdp"
1158                   else
1159                     bssname = "6311G"
1160                   end
1161                 end
1162                 if bssname
1163                   user_input["basis"] = $gamess_basis_keys.find_index(bssname).to_s
1164                 end
1165           #  puts user_input.inspect
1166           end
1167         end
1168
1169     hash = Dialog.run("GAMESS Export") {
1170       def load_basis_set_sub(item)
1171             fname = Dialog.open_panel("Select a file containing GAMESS basis set:")
1172                 if fname
1173                   Molecule.read_gamess_basis_sets(fname)
1174                   bset_desc_new = $gamess_basis_keys.map { |key| $gamess_basis_desc[key] }
1175                   sel1 = attr("basis", :value)
1176                   sel2 = attr("secondary_basis", :value)
1177                   set_attr("basis", :subitems=>bset_desc_new)
1178                   set_attr("basis", :value=>sel1)
1179                   set_attr("secondary_basis", :subitems=>bset_desc_new)
1180                   set_attr("secondary_basis", :value=>sel2)
1181                 end
1182       end
1183           def select_gamess_path(item)
1184             while 1
1185               fname = Dialog.open_panel("Locate GAMESS executable:")
1186                   return if fname == nil
1187                   bname = File.basename(fname)
1188                   if bname =~ /gamess\.(.*)\.(exe|x)$/i
1189                     set_attr("executable_path", :value=>fname)
1190                     return
1191                   else
1192                     error_message_box("\"#{bname}\" does not look like a GAMESS executable!  Please try again.")
1193                   end
1194                 end
1195           end
1196           def set_optional_scripts(item)
1197             h = Dialog.run("GAMESS Optional Scripts") {
1198                   s_pre = get_global_settings("gamess.prefix_script")
1199                   s_post = get_global_settings("gamess.postfix_script")
1200                   layout(1,
1201                     item(:text, :title=>"Script to run before GAMESS execution:"),
1202                         item(:textview, :width=>400, :height=>200, :value=>s_pre, :tag=>"prefix"),
1203                     item(:text, :title=>"Script to run after GAMESS execution:"),
1204                         item(:textview, :width=>400, :height=>200, :value=>s_pre, :tag=>"postfix"))
1205                 }
1206                 if h[:status] == 0
1207                   set_global_settings("gamess.prefix_script", h["prefix"])
1208                   set_global_settings("gamess.postfix_script", h["postfix"])
1209                 end
1210           end
1211           layout(4,
1212                 item(:text, :title=>"SCF type"),
1213                 item(:popup, :subitems=>["RHF", "ROHF", "UHF"], :tag=>"scftype"),
1214             item(:text, :title=>"Run type"),
1215                 item(:popup, :subitems=>["Energy", "Property", "Optimize"], :tag=>"runtype",
1216                   :action=>lambda { |it| set_attr("use_internal", :enabled=>(it[:value] == 2)) } ),
1217
1218                 item(:checkbox, :title=>"Use internal coordinates for structure optimization", :tag=>"use_internal"),
1219                 -1, -1, -1,
1220
1221                 item(:text, :title=>"Charge"),
1222                 item(:textfield, :width=>80, :tag=>"charge"),
1223                 item(:text, :title=>"Multiplicity"),
1224                 item(:textfield, :width=>80, :tag=>"mult"),
1225
1226                 item(:checkbox, :title=>"Use DFT", :tag=>"dft",
1227                   :action=>lambda { |it| set_attr("dfttype", :enabled=>(it[:value] != 0)) } ),
1228                 -1,
1229                 item(:text, :title=>"DFT type"),
1230                 item(:popup, :subitems=>dft_desc, :tag=>"dfttype"),
1231
1232                 item(:line),
1233                 -1, -1, -1,
1234
1235                 item(:text, :title=>"Basis set"),
1236                 item(:popup, :subitems=>bset_desc + ["(no select)"], :tag=>"basis"),
1237                 -1,
1238                 -1,
1239
1240                 item(:button, :title=>"Load Basis Set...", :action=>:load_basis_set_sub),
1241                 -1, -1, -1,
1242
1243                 item(:checkbox, :title=>"Use secondary basis set", :tag=>"use_secondary_basis",
1244                   :action=>lambda { |it|
1245                     flag = (it[:value] != 0)
1246                         set_attr("secondary_elements", :enabled=>flag)
1247                         set_attr("secondary_basis", :enabled=>flag)
1248                   }),
1249                 -1, -1, -1,
1250
1251                 item(:text, :title=>"   Elements"),
1252                 item(:textfield, :width=>80, :tag=>"secondary_elements"),
1253                 item(:text, :title=>"Basis set"),
1254                 item(:popup, :subitems=>bset_desc, :tag=>"secondary_basis"),
1255
1256                 item(:line),
1257                 -1, -1, -1,
1258
1259                 item(:checkbox, :title=>"Calculate electrostatic potential (ESP)", :tag=>"esp"),
1260                 -1, -1, -1,
1261
1262                 item(:line),
1263                 -1, -1, -1,
1264
1265                 item(:checkbox, :title=>"Execute GAMESS on this machine", :tag=>"execute_local",
1266                   :action=>lambda { |it|
1267                     flag = (it[:value] != 0)
1268                         set_attr("executable_path", :enabled=>flag)
1269                         set_attr("select_path", :enabled=>flag)
1270                         set_attr("ncpus", :enabled=>flag)
1271                   }),
1272                 -1, -1, -1,
1273
1274                 item(:text, :title=>"   Path"),
1275                 item(:textfield, :width=>300, :tag=>"executable_path"),
1276                 -1, -1,
1277
1278                 -1,
1279                 item(:button, :title=>"Select Path...", :tag=>"select_path", :action=>:select_gamess_path),
1280                 -1,
1281                 item(:button, :title=>"Optional Scripts...", :action=>:set_optional_scripts),
1282
1283                 item(:text, :title=>"   N of CPUs"),
1284                 item(:textfield, :width=>80, :tag=>"ncpus"),
1285                 -1, -1,
1286
1287                 item(:line),
1288                 -1, -1, -1,
1289
1290                 item(:button, :title=>"Edit GAMESS Input and Go", :action=>lambda { |it|
1291                   h = Hash.new
1292                   each_item { |it2|
1293                     if (tag = it2[:tag]) != nil
1294                           h[tag] = it2[:value]
1295                         end
1296                   }
1297                   h["basis"] = $gamess_basis_keys[h["basis"]]
1298                   h["secondary_basis"] = $gamess_basis_keys[h["secondary_basis"]]
1299                   h["dfttype"] = dft_internal[h["dfttype"]]
1300                   gamess_input_direct = mol.cmd_edit_gamess_input(mol.export_gamess(nil, h))
1301                   if gamess_input_direct
1302                         end_modal(0)
1303                   end
1304                 }),
1305                 -1, -1, -1
1306           )
1307           values = Hash.new
1308           each_item { |it|
1309             tag = it[:tag]
1310                 if tag
1311                   values[tag] = (user_input[tag] || get_global_settings("gamess.#{tag}") || defaults[tag])
1312                   if tag == "basis" && values[tag] == "-1"
1313                     values[tag] = (bset_desc.count).to_s
1314                   end
1315                   it[:value] = values[tag]
1316                 end
1317           }
1318           set_attr("secondary_elements", :enabled=>(values["use_secondary_basis"] == 1))
1319           set_attr("secondary_basis", :enabled=>(values["use_secondary_basis"] == 1))
1320           set_attr("dfttype", :enabled=>(values["dft"] == 1))
1321           set_attr("use_internal", :enabled=>(values["runtype"] == 2))
1322           set_attr("executable_path", :enabled=>(values["execute_local"] == 1))
1323           set_attr("select_path", :enabled=>(values["execute_local"] == 1))
1324           set_attr("ncpus", :enabled=>(values["execute_local"] == 1))
1325         }
1326         hash.each_pair { |key, value|
1327           next if key == :status
1328           set_global_settings("gamess.#{key}", value)
1329         }
1330         if hash[:status] == 0
1331           #  Specify basis by internal keys
1332           hash["basis"] = $gamess_basis_keys[hash["basis"]]
1333           hash["secondary_basis"] = $gamess_basis_keys[hash["secondary_basis"]]
1334           hash["dfttype"] = dft_internal[hash["dfttype"]]
1335           basename = (self.path ? File.basename(self.path, ".*") : self.name)
1336           fname = Dialog.save_panel("Export GAMESS input file:", self.dir, basename + ".inp", "GAMESS input file (*.inp)|*.inp|All files|*.*")
1337           return nil if !fname
1338           if gamess_input_direct
1339           #  puts "gamess_input_direct = \"#{gamess_input_direct}\""
1340             File.open(fname, "w") { |fp| fp.print(gamess_input_direct) }
1341           else
1342             export_gamess(fname, hash)
1343           end
1344           if hash["execute_local"] == 1
1345             Molecule.execute_gamess(fname, self)
1346           end
1347         else
1348           nil
1349         end
1350
1351   end
1352
1353   def cmd_create_cube
1354     grid = default_MO_grid
1355         if grid == nil
1356           Dialog.run {
1357             layout(1,
1358                   item(:text, :title=>"This molecule does not contain MO information."))
1359           }
1360           return
1361         end
1362
1363     mos = selected_MO
1364         if mos == nil || mos.length == 0
1365       Dialog.run {
1366             layout(1,
1367                   item(:text, :title=>"Please select MO(s) in the MO Info table."))
1368       }
1369           return
1370         end
1371         hash = Dialog.run {
1372           layout(1,
1373             item(:text, :title=>"Please specify cube dimensions (in angstrom units):"),
1374             layout(4,
1375                   item(:text, :title=>"Origin"),
1376                   item(:textfield, :width=>100, :height=>20, :tag=>"originx", :value=>sprintf("%.6f", grid[0].x)),
1377                   item(:textfield, :width=>100, :height=>20, :tag=>"originy", :value=>sprintf("%.6f", grid[0].y)),
1378                   item(:textfield, :width=>100, :height=>20, :tag=>"originz", :value=>sprintf("%.6f", grid[0].z)),
1379                   item(:text, :title=>"Delta"),
1380                   item(:textfield, :width=>100, :height=>20, :tag=>"deltax", :value=>sprintf("%.6f", grid[1])),
1381                   item(:textfield, :width=>100, :height=>20, :tag=>"deltay", :value=>sprintf("%.6f", grid[2])),
1382                   item(:textfield, :width=>100, :height=>20, :tag=>"deltaz", :value=>sprintf("%.6f", grid[3])),
1383                   item(:text, :title=>"Step"),
1384                   item(:textfield, :width=>100, :height=>20, :tag=>"stepx", :value=>grid[4].to_s),
1385                   item(:textfield, :width=>100, :height=>20, :tag=>"stepy", :value=>grid[5].to_s),
1386                   item(:textfield, :width=>100, :height=>20, :tag=>"stepz", :value=>grid[6].to_s)))
1387         }
1388         if hash[:status] == 0
1389           path = self.path || self.name
1390           dir = self.dir || Dir.pwd
1391           origin = Vector3D[hash["originx"], hash["originy"], hash["originz"]]
1392           dx = hash["deltax"]
1393           dy = hash["deltay"]
1394           dz = hash["deltaz"]
1395           nx = hash["stepx"]
1396           ny = hash["stepy"]
1397           nz = hash["stepz"]
1398           basename = File.basename(path, ".*")
1399           filenames = []
1400           mo_type = self.mo_type
1401           mos.each { |n|
1402             fname1 = fname2 = nil
1403             alpha = (mo_type != "UHF" ? "" : "alpha ")
1404                 a = (mo_type != "UHF" ? "" : "a")
1405             fname1 = Dialog.save_panel("Cube file name for #{alpha}MO #{n}", dir, basename + "_#{n}#{a}.cube", "Gaussian cube file (*.cube)|*.cube")
1406                 if (mo_type == "UHF")
1407                   fname2 = Dialog.save_panel("Cube file name for beta MO #{n}", dir, basename + "_#{n}b.cube", "Gaussian cube file (*.cube)|*.cube")
1408                 end
1409                 filenames.push([n, fname1, fname2])
1410           }
1411           filenames.each { |pair|
1412             n = pair[0]
1413                 alpha = (mo_type != "UHF" ? "" : "alpha ")
1414             show_progress_panel("Creating cube file for #{alpha}MO #{n}...")
1415                 if pair[1]
1416                   cubegen(pair[1], n, origin, dx, dy, dz, nx, ny, nz, true)
1417                 end
1418                 if pair[2] && mo_type == "UHF"
1419                   set_progress_message("Creating cube file for beta MO #{n}...")
1420                   cubegen(pair[2], n, origin, dx, dy, dz, nx, ny, nz, true, true)
1421                 end
1422             hide_progress_panel
1423           }
1424         end
1425   end
1426
1427 end
1428
1429 $gamess_basis = {
1430   "PM3"   => " PM3 0\n",
1431   "STO3G" => " STO 3\n",
1432   "321G"  => " N21 3\n",
1433   "631G"  => " N31 6\n" }
1434 $gamess_basis_desc = {
1435   "PM3"   => "PM3",
1436   "STO3G" => "STO-3G",
1437   "321G"  => "3-21G",
1438   "631G"  => "6-31G" }
1439 $gamess_basis_keys = ["PM3", "STO3G", "321G", "631G"]
1440
1441 ["631Gd", "631Gdp", "631+Gd", "631++Gdp", "6311Gdp", "6311+Gd", "6311++Gdp", "6311++G2d2p", "6311++G3df3pd", "LanL2DZ"].each { |n|
1442   Molecule.read_gamess_basis_sets("basis_sets/#{n}.txt")
1443 }
1444
1445 register_menu("QChem\tCreate GAMESS Input...",
1446   :cmd_create_gamess_input, :non_empty)
1447 register_menu("QChem\tCreate MOPAC6 Input...",
1448   :cmd_create_mopac_input, :non_empty)    # mopac6.rb
1449 register_menu("QChem\tCreate MO Cube...",
1450   :cmd_create_cube, lambda { |m| m && m.mo_type } )