Index: ThirdParty/mpqc_open/lib/perl/AtomicBases.pm
===================================================================
--- ThirdParty/mpqc_open/lib/perl/AtomicBases.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/AtomicBases.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,320 @@
+#
+eval 'exec perl $0 $*'
+    if 0;
+
+require AtomicBasis;
+
+package AtomicBases;
+
+$symbol_to_name{"H"}="hydrogen";
+$symbol_to_name{"He"}="helium";
+$symbol_to_name{"Li"}="lithium";
+$symbol_to_name{"Be"}="beryllium";
+$symbol_to_name{"B"}="boron";
+$symbol_to_name{"C"}="carbon";
+$symbol_to_name{"N"}="nitrogen";
+$symbol_to_name{"O"}="oxygen";
+$symbol_to_name{"F"}="fluorine";
+$symbol_to_name{"Ne"}="neon";
+$symbol_to_name{"Na"}="sodium";
+$symbol_to_name{"Mg"}="magnesium";
+$symbol_to_name{"Al"}="aluminum";
+$symbol_to_name{"Si"}="silicon";
+$symbol_to_name{"P"}="phosphorus";
+$symbol_to_name{"S"}="sulfur";
+$symbol_to_name{"Cl"}="chlorine";
+$symbol_to_name{"Ar"}="argon";
+$symbol_to_name{"K"}="potassium";
+$symbol_to_name{"Ca"}="calcium";
+$symbol_to_name{"Sc"}="scandium";
+$symbol_to_name{"Ti"}="titanium";
+$symbol_to_name{"V"}="vanadium";
+$symbol_to_name{"Cr"}="chromium";
+$symbol_to_name{"Mn"}="manganese";
+$symbol_to_name{"Fe"}="iron";
+$symbol_to_name{"Co"}="cobalt";
+$symbol_to_name{"Ni"}="nickel";
+$symbol_to_name{"Cu"}="copper";
+$symbol_to_name{"Zn"}="zinc";
+$symbol_to_name{"Ga"}="gallium";
+$symbol_to_name{"Ge"}="germanium";
+$symbol_to_name{"As"}="arsenic";
+$symbol_to_name{"Se"}="selenium";
+$symbol_to_name{"Br"}="bromine";
+$symbol_to_name{"Kr"}="krypton";
+$symbol_to_name{"Rb"}="rubidium";
+$symbol_to_name{"Sr"}="strontium";
+$symbol_to_name{"Y"}="yttrium";
+$symbol_to_name{"Zr"}="zirconium";
+$symbol_to_name{"Nb"}="niobium";
+$symbol_to_name{"Mo"}="molybdenum";
+$symbol_to_name{"Tc"}="technetium";
+$symbol_to_name{"Ru"}="ruthenium";
+$symbol_to_name{"Rh"}="rhodium";
+$symbol_to_name{"Pd"}="palladium";
+$symbol_to_name{"Ag"}="silver";
+$symbol_to_name{"Cd"}="cadminium";
+$symbol_to_name{"In"}="indium";
+$symbol_to_name{"Sn"}="tin";
+$symbol_to_name{"Sb"}="antimony";
+$symbol_to_name{"Te"}="tellurium";
+$symbol_to_name{"I"}="iodine";
+$symbol_to_name{"Xe"}="xenon";
+$symbol_to_name{"Cs"}="cesium";
+$symbol_to_name{"Ba"}="barium";
+$symbol_to_name{"La"}="lanthanium";
+$symbol_to_name{"Ce"}="cerium";
+$symbol_to_name{"Pr"}="praseodymium";
+$symbol_to_name{"Nd"}="neodymium";
+$symbol_to_name{"Pm"}="promethium";
+$symbol_to_name{"Sm"}="samarium";
+$symbol_to_name{"Eu"}="europium";
+$symbol_to_name{"Gd"}="gadolinium";
+$symbol_to_name{"Tb"}="terbium";
+$symbol_to_name{"Dy"}="dysprosium";
+$symbol_to_name{"Ho"}="holmium";
+$symbol_to_name{"Er"}="erbium";
+$symbol_to_name{"Tm"}="thulium";
+$symbol_to_name{"Yb"}="ytterbium";
+$symbol_to_name{"Lu"}="lutetium";
+$symbol_to_name{"Hf"}="hafnium";
+$symbol_to_name{"Ta"}="tantalum";
+$symbol_to_name{"W"}="tungsten";
+$symbol_to_name{"Re"}="rhenium";
+$symbol_to_name{"Os"}="osmium";
+$symbol_to_name{"Ir"}="iridium";
+$symbol_to_name{"Pt"}="platinum";
+$symbol_to_name{"Au"}="gold";
+$symbol_to_name{"Hg"}="mercury";
+$symbol_to_name{"Tl"}="thallium";
+$symbol_to_name{"Pb"}="lead";
+$symbol_to_name{"Bi"}="bismuth";
+$symbol_to_name{"Po"}="polonium";
+$symbol_to_name{"At"}="astatine";
+$symbol_to_name{"Rn"}="radon";
+$symbol_to_name{"Fr"}="francium";
+$symbol_to_name{"Ra"}="radium";
+$symbol_to_name{"Ac"}="actinium";
+$symbol_to_name{"Th"}="thorium";
+$symbol_to_name{"Pa"}="protactinium";
+$symbol_to_name{"U"}="uranium";
+$symbol_to_name{"Np"}="neptunium";
+$symbol_to_name{"Pu"}="plutonium";
+$symbol_to_name{"Am"}="americium";
+$symbol_to_name{"Cm"}="curium";
+$symbol_to_name{"Bk"}="berkelium";
+$symbol_to_name{"Cf"}="californium";
+$symbol_to_name{"Es"}="einsteinum";
+$symbol_to_name{"Fm"}="fermium";
+$symbol_to_name{"Md"}="mendelevium";
+$symbol_to_name{"No"}="nobelium";
+$symbol_to_name{"Lr"}="lawrencium";
+
+@number_to_symbol = ();
+$number_to_symbol[1] = "H";
+$number_to_symbol[2] = "He";
+$number_to_symbol[3] = "Li";
+$number_to_symbol[4] = "Be";
+$number_to_symbol[5] = "B";
+$number_to_symbol[6] = "C";
+$number_to_symbol[7] = "N";
+$number_to_symbol[8] = "O";
+$number_to_symbol[9] = "F";
+$number_to_symbol[10] = "Ne";
+$number_to_symbol[11] = "Na";
+$number_to_symbol[12] = "Mg";
+$number_to_symbol[13] = "Al";
+$number_to_symbol[14] = "Si";
+$number_to_symbol[15] = "P";
+$number_to_symbol[16] = "S";
+$number_to_symbol[17] = "Cl";
+$number_to_symbol[18] = "Ar";
+$number_to_symbol[19] = "K";
+$number_to_symbol[20] = "Ca";
+$number_to_symbol[21] = "Sc";
+$number_to_symbol[22] = "Ti";
+$number_to_symbol[23] = "V";
+$number_to_symbol[24] = "Cr";
+$number_to_symbol[25] = "Mn";
+$number_to_symbol[26] = "Fe";
+$number_to_symbol[27] = "Co";
+$number_to_symbol[28] = "Ni";
+$number_to_symbol[29] = "Cu";
+$number_to_symbol[30] = "Zn";
+$number_to_symbol[31] = "Ga";
+$number_to_symbol[32] = "Ge";
+$number_to_symbol[33] = "As";
+$number_to_symbol[34] = "Se";
+$number_to_symbol[35] = "Br";
+$number_to_symbol[36] = "Kr";
+$number_to_symbol[37] = "Rb";
+$number_to_symbol[38] = "Sr";
+$number_to_symbol[39] = "Y";
+$number_to_symbol[40] = "Zr";
+$number_to_symbol[41] = "Nb";
+$number_to_symbol[42] = "Mo";
+$number_to_symbol[43] = "Tc";
+$number_to_symbol[44] = "Ru";
+$number_to_symbol[45] = "Rh";
+$number_to_symbol[46] = "Pd";
+$number_to_symbol[47] = "Ag";
+$number_to_symbol[48] = "Cd";
+$number_to_symbol[49] = "In";
+$number_to_symbol[50] = "Sn";
+$number_to_symbol[51] = "Sb";
+$number_to_symbol[52] = "Te";
+$number_to_symbol[53] = "I";
+$number_to_symbol[54] = "Xe";
+$number_to_symbol[55] = "Cs";
+$number_to_symbol[56] = "Ba";
+$number_to_symbol[57] = "La";
+$number_to_symbol[58] = "Ce";
+$number_to_symbol[59] = "Pr";
+$number_to_symbol[60] = "Nd";
+$number_to_symbol[61] = "Pm";
+$number_to_symbol[62] = "Sm";
+$number_to_symbol[63] = "Eu";
+$number_to_symbol[64] = "Gd";
+$number_to_symbol[65] = "Tb";
+$number_to_symbol[66] = "Dy";
+$number_to_symbol[67] = "Ho";
+$number_to_symbol[68] = "Er";
+$number_to_symbol[69] = "Tm";
+$number_to_symbol[70] = "Yb";
+$number_to_symbol[71] = "Lu";
+$number_to_symbol[72] = "Hf";
+$number_to_symbol[73] = "Ta";
+$number_to_symbol[74] = "W";
+$number_to_symbol[75] = "Re";
+$number_to_symbol[76] = "Os";
+$number_to_symbol[77] = "Ir";
+$number_to_symbol[78] = "Pt";
+$number_to_symbol[79] = "Au";
+$number_to_symbol[80] = "Hg";
+$number_to_symbol[81] = "Tl";
+$number_to_symbol[82] = "Pb";
+$number_to_symbol[83] = "Bi";
+$number_to_symbol[84] = "Po";
+$number_to_symbol[85] = "At";
+$number_to_symbol[86] = "Rn";
+$number_to_symbol[87] = "Fr";
+$number_to_symbol[88] = "Ra";
+$number_to_symbol[89] = "Ac";
+$number_to_symbol[90] = "Th";
+$number_to_symbol[91] = "Pa";
+$number_to_symbol[92] = "U";
+$number_to_symbol[93] = "Np";
+$number_to_symbol[94] = "Pu";
+$number_to_symbol[95] = "Am";
+$number_to_symbol[96] = "Cm";
+$number_to_symbol[97] = "Bk";
+$number_to_symbol[98] = "Cf";
+$number_to_symbol[99] = "Es";
+$number_to_symbol[100] = "Fm";
+$number_to_symbol[101] = "Md";
+$number_to_symbol[102] = "No";
+$number_to_symbol[103] = "Lr";
+$number_to_symbol[104] = "Rf";
+$number_to_symbol[105] = "Ha";
+
+%symbol_to_number = {};
+foreach my $i (1..105) {
+    $symbol_to_number{$number_to_symbol[$i]} = $i;
+}
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+
+    $self->{"name"} = "";
+    $self->{"atoms"} = {};
+    $self->{"comment"} = {};
+}
+
+sub name {
+    my $self = shift;
+    return $self->{"name"};
+}
+
+sub set_name {
+    my $self = shift;
+    my $name = shift;
+    $self->{"name"} = $name;
+}
+
+sub read_gaussian_file {
+    my $self = shift;
+    my $filename = shift;
+    open(FILE,"<$filename");
+    $self->read_gaussian(*FILE);
+    close(FILE);
+}
+
+
+sub read_gaussian {
+    my $self = shift;
+    my $file = shift;
+    while (<$file>) {
+        if (/([A-Za-z]+)/) {
+            my $sym = canonical_symbol($1);
+            if (!exists($symbol_to_name{$sym})) {
+                die "bad atomic symbol \"$sym\"";
+            }
+            $self->{"atoms"}->{$sym} = new AtomicBasis();
+            $self->{"atoms"}->{$sym}->read_gaussian($file);
+        }
+    }
+}
+
+sub atoms {
+    my $self = shift;
+    return keys(%{$self->{"atoms"}});
+}
+
+sub comment {
+    my $self = shift;
+    return $self->{"comment"};
+}
+
+sub atom {
+    my $self = shift;
+    my $sym = canonical_symbol(shift);
+    return $self->{"atoms"}->{$sym};
+}
+
+sub canonical_symbol {
+    my $sym = shift;
+    return uc(substr($sym,0,1)) . lc(substr($sym,1,length($sym)));
+}
+
+sub write_keyval_file {
+    my $self = shift;
+    my $filename = shift;
+    open(FILE,">$filename");
+    $self->write_keyval(*FILE);
+    close(FILE);
+};
+
+sub write_keyval {
+    my $self = shift;
+    my $file = shift;
+    my @atoms = $self->atoms();
+    print $file "basis: (\n";
+    foreach my $atom (sort { $symbol_to_number{$a} <=> $symbol_to_number{$b} }
+                      @atoms) {
+        printf $file " %s: \"%s\": [\n", $symbol_to_name{$atom}, $self->name();
+        $self->atom($atom)->write_keyval($file);
+        print $file " ]\n";
+    }
+    print $file ")\n";
+}
+
+1;
Index: ThirdParty/mpqc_open/lib/perl/AtomicBasis.pm
===================================================================
--- ThirdParty/mpqc_open/lib/perl/AtomicBasis.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/AtomicBasis.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,225 @@
+#
+eval 'exec perl $0 $*'
+    if 0;
+
+package AtomicBasis;
+
+$fltrx = "((?:-?\\d+|-?\\d+\\.\\d*|-?\\d*\\.\\d+)(?:[eEdD][+-]?\\d+)?)";
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+
+    $self->{"pure_d"} = 1;
+    $self->{"pure_f_plus"} = 1;
+    $self->{"exponents"} = [];
+    $self->{"coefficients"} = [];
+    $self->{"am"} = [];
+    $self->{"pure"} = [];
+}
+
+sub default_pure {
+    my $self = shift;
+    my $am = shift;
+    if ($am == 2 && $self->{"pure_d"}) { return 1; }
+    elsif ($am > 2 && $self->{"pure_f_plus"}) { return 1; }
+    return 0;
+}
+
+my $amtypes = "SPDFGHIKLMN";
+
+sub am_string_to_number {
+    my $amstr = uc(shift);
+    if (length($amstr) != 1) {
+        die "invalid am string \"$amstr\"";
+    }
+    my $index = index($amtypes,$amstr);
+    if ($index == -1) {
+        die "am string \"$amstr\" not found";
+    }
+    return $index;
+}
+
+sub am_number_to_string {
+    my $am = shift;
+    return substr($amtypes,$am,1);
+}
+
+sub stdflt {
+    my $num = shift;
+    $num =~ s/D/e/;
+    $num =~ s/d/e/;
+    $num =~ s/E/e/;
+    return $num;
+}
+
+sub read_gaussian {
+    my $self = shift;
+    my $file = shift;
+    my $ishell = 0;
+    while (<$file>) {
+        if (/\*\*\*\*/) {
+            last;
+        }
+        elsif (/([A-Za-z]+) +([0-9]+) +([0-9.]+)/) {
+            my $amstr = uc($1);
+            my $nprim = $2;
+            my $scale = $3;
+            if ($scale != 1) {
+                die "cannot handle scale = $scale (must be 1)";
+            }
+            if ($amstr eq "SP") {
+                $self->{"am"}->[$ishell]->[0] = 0;
+                $self->{"am"}->[$ishell]->[1] = 1;
+                $self->{"pure"}->[$ishell]->[0] = 0;
+                $self->{"pure"}->[$ishell]->[1] = 0;
+                foreach my $i (0..$nprim-1) {
+                    while (<$file>) { last; }
+                    if (/$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                        my $exp = $1;
+                        my $coefs = $2;
+                        my $coefp = $3;
+                        $self->{"exponents"}->[$ishell]->[$i]
+                            = stdflt($exp);
+                        $self->{"coefficients"}->[$ishell]->[$i]->[0]
+                            = stdflt($coefs);
+                        $self->{"coefficients"}->[$ishell]->[$i]->[1]
+                            = stdflt($coefp);
+                    }
+                    else {
+                        die "bad exponent coefficient line";
+                    }
+                }
+            }
+            else {
+                my $am = am_string_to_number($amstr);
+                $self->{"am"}->[$ishell]->[0] = $am;
+                $self->{"pure"}->[$ishell]->[0] = $self->default_pure($am);
+                foreach my $i (0..$nprim-1) {
+                    while (<$file>) { last; }
+                    if (/$fltrx\s+$fltrx\s*$/) {
+                        my $exp = $1;
+                        my $coef = $2;
+                        $self->{"exponents"}->[$ishell]->[$i]
+                            = stdflt($exp);
+                        $self->{"coefficients"}->[$ishell]->[$i]->[0]
+                            = stdflt($coef);
+                    }
+                    else {
+                        die "bad exponent coefficient line";
+                    }
+                }
+            }
+            $ishell++;
+        }
+        else {
+            die "could not parse line $_";
+        }
+    }
+}
+
+sub nshell {
+    my $self = shift;
+    my @exp = @{$self->{"exponents"}};
+    return $#exp + 1;
+}
+
+sub nprim {
+    my $self = shift;
+    my $ishell = shift;
+    my @exp = @{$self->{"exponents"}->[$ishell]};
+    return $#exp + 1;
+}
+
+sub ncon {
+    my $self = shift;
+    my $ishell = shift;
+    my @exp = @{$self->{"am"}->[$ishell]};
+    return $#exp + 1;
+}
+
+sub exp {
+    my $self = shift;
+    my $ishell = shift;
+    my $iprim = shift;
+    my $exp = $self->{"exponents"}->[$ishell]->[$iprim];
+    return $exp;
+}
+
+sub coef {
+    my $self = shift;
+    my $ishell = shift;
+    my $iprim = shift;
+    my $icon = shift;
+    my $coef = $self->{"coefficients"}->[$ishell]->[$iprim]->[$icon];
+    return $coef;
+}
+
+sub am {
+    my $self = shift;
+    my $ishell = shift;
+    my $icon = shift;
+    my $am = $self->{"am"}->[$ishell]->[$icon];
+    return $am;
+}
+
+sub pure {
+    my $self = shift;
+    my $ishell = shift;
+    my $icon = shift;
+    my $pure = $self->{"pure"}->[$ishell]->[$icon];
+    return $pure;
+}
+
+sub amstr {
+    my $self = shift;
+    my $ishell = shift;
+    my $icon = shift;
+    my $am = $self->{"am"}->[$ishell]->[$icon];
+    return am_number_to_string($am);
+}
+
+sub write_keyval {
+    my $self = shift;
+    my $file = shift;
+    foreach my $ishell (0..$self->nshell()-1) {
+        print $file "  (type: [";
+        # write out am (and puream)
+        foreach my $icon (0..$self->ncon($ishell)-1) {
+            if ($icon > 0) { print $file " "; }
+            if ($self->pure($ishell,$icon)) {
+                printf $file "(am=%s puream=1)",
+                lc($self->amstr($ishell,$icon));
+            }
+            else {
+                printf $file "am=%s", lc($self->amstr($ishell,$icon));
+            }
+        }
+        print $file "]\n";
+        print $file "   {exp";
+        # write out coef:0...
+        foreach my $icon (0..$self->ncon($ishell)-1) {
+            print $file " coef:$icon";
+        }
+        print "} = {\n";
+        foreach my $iprim (0..$self->nprim($ishell)-1) {
+            printf $file "   %s", $self->exp($ishell,$iprim);
+            # write out coefficients
+            foreach my $icon (0..$self->ncon($ishell)-1) {
+                printf $file " %s", $self->coef($ishell,$iprim,$icon);
+            }
+            printf $file "\n";
+        }
+        print $file "  })\n";
+    }
+}
+
+1;
Index: ThirdParty/mpqc_open/lib/perl/Makefile
===================================================================
--- ThirdParty/mpqc_open/lib/perl/Makefile	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/Makefile	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,12 @@
+TOPDIR=../..
+ifndef SRCDIR
+  SRCDIR=$(shell pwd)
+endif
+
+include $(SRCDIR)/$(TOPDIR)/lib/GlobalMakefile
+
+include $(SRCDIR)/$(TOPDIR)/lib/GlobalRules
+
+install::
+	$(INSTALL) $(INSTALLDIROPT) $(installroot)$(scdatadir)/perl
+	$(INSTALL) $(INSTALLLIBOPT) $(SRCDIR)/*.pm $(installroot)$(scdatadir)/perl
Index: ThirdParty/mpqc_open/lib/perl/Molecule.pm
===================================================================
--- ThirdParty/mpqc_open/lib/perl/Molecule.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/Molecule.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,382 @@
+#
+eval 'exec perl $0 $*'
+    if 0;
+
+use POSIX; # for acos
+require QCParse;
+
+##########################################################################
+
+package Molecule;
+$debug = 0;
+$fltrx = "((?:-?\\d+|-?\\d+\\.\\d*|-?\\d*\\.\\d+)(?:[eE][+-]?\\d+)?)";
+
+%z_to_sym = (
+    "1" => "h",
+    "2" => "he",
+    "3" => "li",
+    "4" => "be",
+    "5" => "b",
+    "6" => "c",
+    "7" => "n",
+    "8" => "o",
+    "9" => "f",
+    "10" => "ne",
+    "11" => "na",
+    "12" => "mg",
+    "13" => "al",
+    "14" => "si",
+    "15" => "p",
+    "16" => "s",
+    "17" => "cl",
+    "18" => "ar",
+    "19" => "k",
+    "20" => "ca",
+    "21" => "sc",
+    "22" => "ti",
+    "23" => "v",
+    "24" => "cr",
+    "25" => "mn",
+    "26" => "fe",
+    "27" => "co",
+    "28" => "ni",
+    "29" => "cu",
+    "30" => "zn",
+    "31" => "ga",
+    "32" => "ge",
+    "33" => "as",
+    "34" => "se",
+    "35" => "br",
+    "36" => "kr",
+    "37" => "rb",
+    "38" => "sr",
+    "39" => "y",
+    "40" => "zr",
+    "41" => "nb",
+    "42" => "mo",
+    "43" => "tc",
+    "44" => "ru",
+    "45" => "rh",
+    "46" => "pd",
+    "47" => "ag",
+    "48" => "cd",
+    "49" => "in",
+    "50" => "sn",
+    "51" => "sb",
+    "52" => "te",
+    "53" => "i",
+    "54" => "xe",
+    "55" => "cs",
+    "56" => "ba",
+    "57" => "la",
+    "58" => "ce",
+    "59" => "pr",
+    "60" => "nd",
+    "61" => "pm",
+    "62" => "sm",
+    "63" => "eu",
+    "64" => "gd",
+    "65" => "tb",
+    "66" => "dy",
+    "67" => "ho",
+    "68" => "er",
+    "69" => "tm",
+    "70" => "yb",
+    "71" => "lu",
+    "72" => "hf",
+    "73" => "ta",
+    "74" => "w",
+    "75" => "re",
+    "76" => "os",
+    "77" => "ir",
+    "78" => "pt",
+    "79" => "au",
+    "80" => "hg",
+    "81" => "tl",
+    "82" => "pb",
+    "83" => "bi",
+    "84" => "po",
+    "85" => "at",
+    "86" => "rn",
+    "87" => "fr",
+    "88" => "ra",
+    "89" => "ac",
+    "90" => "th",
+    "91" => "pa",
+    "92" => "u",
+    "93" => "np",
+    "94" => "pu",
+    "95" => "am",
+    "96" => "cm",
+    "97" => "bk",
+    "98" => "cf",
+    "99" => "es",
+    "100" => "fm",
+    "101" => "md",
+    "102" => "no",
+    "103" => "lr",
+    "104" => "rf",
+    "105" => "ha"
+);
+
+%sym_to_z = (
+    "h" => "1",
+    "he" => "2",
+    "li" => "3",
+    "be" => "4",
+    "b" => "5",
+    "c" => "6",
+    "n" => "7",
+    "o" => "8",
+    "f" => "9",
+    "ne" => "10",
+    "na" => "11",
+    "mg" => "12",
+    "al" => "13",
+    "si" => "14",
+    "p" => "15",
+    "s" => "16",
+    "cl" => "17",
+    "ar" => "18",
+    "k" => "19",
+    "ca" => "20",
+    "sc" => "21",
+    "ti" => "22",
+    "v" => "23",
+    "cr" => "24",
+    "mn" => "25",
+    "fe" => "26",
+    "co" => "27",
+    "ni" => "28",
+    "cu" => "29",
+    "zn" => "30",
+    "ga" => "31",
+    "ge" => "32",
+    "as" => "33",
+    "se" => "34",
+    "br" => "35",
+    "kr" => "36",
+    "rb" => "37",
+    "sr" => "38",
+    "y" => "39",
+    "zr" => "40",
+    "nb" => "41",
+    "mo" => "42",
+    "tc" => "43",
+    "ru" => "44",
+    "rh" => "45",
+    "pd" => "46",
+    "ag" => "47",
+    "cd" => "48",
+    "in" => "49",
+    "sn" => "50",
+    "sb" => "51",
+    "te" => "52",
+    "i" => "53",
+    "xe" => "54",
+    "cs" => "55",
+    "ba" => "56",
+    "la" => "57",
+    "ce" => "58",
+    "pr" => "59",
+    "nd" => "60",
+    "pm" => "61",
+    "sm" => "62",
+    "eu" => "63",
+    "gd" => "64",
+    "tb" => "65",
+    "dy" => "66",
+    "ho" => "67",
+    "er" => "68",
+    "tm" => "69",
+    "yb" => "70",
+    "lu" => "71",
+    "hf" => "72",
+    "ta" => "73",
+    "w" => "74",
+    "re" => "75",
+    "os" => "76",
+    "ir" => "77",
+    "pt" => "78",
+    "au" => "79",
+    "hg" => "80",
+    "tl" => "81",
+    "pb" => "82",
+    "bi" => "83",
+    "po" => "84",
+    "at" => "85",
+    "rn" => "86",
+    "fr" => "87",
+    "ra" => "88",
+    "ac" => "89",
+    "th" => "90",
+    "pa" => "91",
+    "u" => "92",
+    "np" => "93",
+    "pu" => "94",
+    "am" => "95",
+    "cm" => "96",
+    "bk" => "97",
+    "cf" => "98",
+    "es" => "99",
+    "fm" => "100",
+    "md" => "101",
+    "no" => "102",
+    "lr" => "103",
+    "rf" => "104",
+    "ha" => "105"
+);
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+    my $mol = shift;
+
+    $self->{"position"} = [];
+    $self->{"element"} = [];
+    my $i = 0;
+    while ($mol =~ s/^\s*(\w+)\s+$fltrx\s+$fltrx\s+$fltrx\s*\n//) {
+        my $sym = $1;
+        my $x = $2;
+        my $y = $3;
+        my $z = $4;
+        $self->{"element"}[$i] = $sym;
+        $self->{"position"}[$i] = [ $x, $y, $z ];
+        $i++;
+    }
+
+    $self->{"natom"} = $i;
+}
+
+sub string {
+    my $self = shift;
+    my $mol = "";
+    for ($i = 0; $i < $self->n_atom(); $i++) {
+        $mol = sprintf "%s  %s  %14.10f %14.10f %14.10f\n",
+                       $mol, $self->element($i),
+                       $self->position($i,0),
+                       $self->position($i,1),
+                       $self->position($i,2);
+    }
+
+    $mol;
+}
+
+sub n_atom {
+    my $self = shift;
+    printf "Molecule: returning natom = %d\n", $self->{"natom"} if ($debug);
+    $self->{"natom"};
+}
+
+sub element {
+    my $self = shift;
+    my $i = shift;
+    $self->{"element"}[$i];
+}
+
+sub z {
+    my $self = shift;
+    my $i = shift;
+    $sym_to_z{lc($self->{"element"}[$i])};
+}
+
+sub position {
+    my $self = shift;
+    my $i = shift;
+    my $xyz = shift;
+    $self->{"position"}[$i][$xyz];
+}
+
+sub geometry {
+    my $self = shift;
+    my $geom = [];
+    my $ngeom = 0;
+    my $fence = $self->{"natom"}-1;
+    foreach $i (0..$fence) {
+        foreach $xyz (0..2) {
+            $geom->[$ngeom] = $self->{"position"}[$i][$xyz];
+            $ngeom = $ngeom + 1;
+        }
+    }
+    $geom;
+}
+
+sub atom_xyz {
+    my $self = shift;
+    my $i = shift;
+    my $geom = [];
+    my $ngeom = 0;
+    foreach $xyz (0..2) {
+        $geom->[$ngeom] = $self->{"position"}[$i][$xyz];
+        $ngeom = $ngeom + 1;
+    }
+    $geom;
+}
+
+sub dot {
+    my $v1 = shift;
+    my $v2 = shift;
+    my $d = 0.0;
+    foreach $xyz (0..2) {
+        $d = $d + $v1->[$xyz] * $v2->[$xyz];
+    }
+    $d;
+}
+
+sub diff {
+    my $v1 = shift;
+    my $v2 = shift;
+    my $diff = [];
+    foreach $xyz (0..2) {
+        $diff->[$xyz] = $v1->[$xyz] - $v2->[$xyz];
+    }
+    $diff;
+}
+
+sub vecstr {
+    my $v = shift;
+    sprintf "%12.8f %12.8f %12.8f", $v->[0], $v->[1], $v->[2];
+}
+
+# numbering starts at 1 for bond
+sub bond {
+    my $self = shift;
+    my $a1 = shift() - 1;
+    my $a2 = shift() - 1;
+    my $d = diff($self->atom_xyz($a1),$self->atom_xyz($a2));
+    #printf "v1 = %s\n", vecstr($self->atom_xyz($a1));
+    #printf "v2 = %s\n", vecstr($self->atom_xyz($a2));
+    #printf "diff = %s\n", vecstr($d);
+    sqrt(dot($d,$d));
+}
+
+# numbering starts at 1 for bend
+sub bend {
+    my $self = shift;
+    my $a1 = shift() - 1;
+    my $a2 = shift() - 1;
+    my $a3 = shift() - 1;
+    my $diff12 = diff($self->atom_xyz($a1), $self->atom_xyz($a2));
+    my $diff32 = diff($self->atom_xyz($a3), $self->atom_xyz($a2));
+    POSIX::acos(dot($diff12,$diff32)/sqrt(dot($diff12,$diff12)*dot($diff32,$diff32))) * 180.0 / 3.14159265358979323846;
+}
+
+# numbering starts at 1 for tors
+sub tors {
+    my $self = shift;
+    my $a1 = shift() - 1;
+    my $a2 = shift() - 1;
+    my $a3 = shift() - 1;
+    my $a4 = shift() - 1;
+    die "Molecule::tors not available";
+}
+
+1;
Index: ThirdParty/mpqc_open/lib/perl/QCParse.pm
===================================================================
--- ThirdParty/mpqc_open/lib/perl/QCParse.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/QCParse.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,1232 @@
+#
+eval 'exec perl $0 $*'
+    if 0;
+
+require Molecule;
+
+##########################################################################
+
+package QCParse;
+$debug = 0;
+
+sub testparse {
+    my $parse = new QCParse;
+
+    my $string = "x:
+ xval
+test_basis: STO-3G 6-311G**
+charge: 1
+method: scf
+basis: sto-3g
+state: 3b1
+molecule:
+  H 0 0.0000001 1.00000001
+  H 0 0 -1
+gradient: yes
+optimize: no
+frequencies: yes
+properties: NPA
+y:
+yval
+z: zval1 zval2
+zval3
+h:
+0 a
+1
+ 2  c";
+
+    print "string:\n--------------\n$string\n--------------\n";
+
+    $parse->parse_string($string);
+    $parse->doprint();
+
+    my @t = $parse->value_as_array('h');
+    print "-----------------\n";
+    for ($i = 0; $i <= $#t; $i++) {
+        print "$i: $t[$i]\n";
+    }
+    print "-----------------\n";
+
+    @t = $parse->value_as_lines('h');
+    print "-----------------\n";
+    for ($i = 0; $i <= $#t; $i++) {
+        print "$i: $t[$i]\n";
+    }
+    print "-----------------\n";
+
+    my $qcinp = new QCInput($parse);
+    my $test_basis = $parse->value("test_basis");
+    my @test_basis_a = $parse->value_as_array("test_basis");
+    my $state = $qcinp->state();
+    my $mult = $qcinp->mult();
+    my $method = $qcinp->method();
+    my $charge = $qcinp->charge();
+    my $basis = $qcinp->basis();
+    my $gradient = $qcinp->gradient();
+    my $frequencies = $qcinp->frequencies();
+    my $optimize = $qcinp->optimize();
+    my $natom = $qcinp->n_atom();
+    foreach $i (@test_basis_a) {
+        print "test_basis_a: $i\n";
+    }
+    print "test_basis = $test_basis\n";
+    print "state = $state\n";
+    print "mult = $mult\n";
+    print "method = $method\n";
+    print "basis = $basis\n";
+    print "optimize = $optimize\n";
+    print "gradient = $gradient\n";
+    print "frequencies = $frequencies\n";
+    print "natom = $natom\n";
+    for ($i = 0; $i < $natom; $i++) {
+        printf "%s %14.8f %14.8f %14.8f\n", $qcinp->element($i),
+                                $qcinp->position($i,0),
+                                $qcinp->position($i,1),
+                                $qcinp->position($i,2);
+    }
+    printf "qcinp errors: %s\n", $qcinp->error();
+
+    my $inpwr = new MPQCInputWriter($qcinp);
+    printf "MPQC input:\n%s", $inpwr->input_string();
+}
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize();
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+    $self->{'keyval'} = {};
+    $self->{'error'} = "";
+}
+
+sub parse_file {
+    my $self = shift;
+    my $file = shift;
+    if (! -f "$file") {
+        $self->{"ok"} = 0;
+        $self->error("File $file not found.");
+        return;
+    }
+    open(INPUT, "<$file");
+    my $string = "";
+    while (<INPUT>) {
+        $string = "$string$_";
+    }
+    close(INPUT);
+    #print "Got file:\n$string\n";
+    $self->parse_string($string);
+    $self->{"ok"} = 1;
+}
+
+sub write_file {
+    my $self = shift;
+    my $file = shift;
+    my $keyval = $self->{'keyval'};
+    my @keys = keys(%$keyval);
+    open(OUTPUT, ">$file");
+    foreach $key (@keys) {
+        my $value = $keyval->{$key};
+        print OUTPUT "${key}:\n";
+        print OUTPUT "$value\n";
+    }
+    close(OUTPUT);
+}
+
+sub parse_string {
+    my $self = shift;
+    my $string = shift;
+    my $value = "";
+    my $keyword = "";
+    $string = "$string\n";
+    while ($string) {
+        $string =~ s/^[^\n]*\n//;
+        $_ = $&;
+        s/#.*//;
+        if (/^\s*(\w+)\s*:\s*(.*)\s*$/) {
+            $self->add($keyword, $value);
+            $keyword = $1;
+            $value = $2;
+        }
+        elsif (/^\s*$/) {
+            $self->add($keyword, $value);
+            $keyword = "";
+            $value = "";
+        }
+        else {
+            $value = "$value$_";
+        }
+    }
+    $self->add($keyword, $value);
+}
+
+sub add {
+    my $self = shift;
+    my $keyword = shift;
+    my $value = shift;
+    if ($keyword ne "") {
+        $self->{'keyval'}{$keyword} = $value;
+        printf("%s = %s\n", $keyword, $value) if ($debug);
+    }
+}
+
+# returns the value of the keyword
+sub value {
+    my $self = shift;
+    my $keyword = shift;
+    my $keyval = $self->{'keyval'};
+    my $value = $keyval->{$keyword};
+    return $value;
+}
+
+# sets the value of the keyword
+sub set_value {
+    my $self = shift;
+    my $keyword = shift;
+    my $value = shift;
+    my $keyval = $self->{'keyval'};
+    $keyval->{$keyword} = $value;
+    return $value;
+}
+
+# returns the value of the keyword
+sub boolean_value {
+    my $self = shift;
+    my $keyword = shift;
+    my $keyval = $self->{'keyval'};
+    $_ = $keyval->{$keyword};
+    return "1" if (/^\s*(y|yes|1|true|t)\s*$/i);
+    return "0" if (/^\s*(n|no|0|false|f|)\s*$/i);
+    "";
+}
+
+# returns an array of whitespace delimited tokens
+sub value_as_array {
+    my $self = shift;
+    my $keyword = shift;
+    my $keyval = $self->{'keyval'};
+    my $value = $keyval->{$keyword};
+    my @array = ();
+    $i = 0;
+    $value =~ s/^\s+$//;
+    while ($value ne '') {
+        $value =~ s/^\s*(\S+)\s*//s;
+        $array[$i] = $1;
+        $i++;
+    }
+    return @array;
+}
+
+# returns an array reference of whitespace delimited tokens
+sub value_as_arrayref {
+    my $self = shift;
+    my $keyword = shift;
+    my $keyval = $self->{'keyval'};
+    my $value = $keyval->{$keyword};
+    my $array = [];
+    $i = 0;
+    $value =~ s/^\s+$//;
+    while ($value ne '') {
+        $value =~ s/^\s*(\S+)\s*//s;
+        $array->[$i] = $1;
+        $i++;
+    }
+    return $array;
+}
+
+# returns an array of lines
+sub value_as_lines {
+    my $self = shift;
+    my $keyword = shift;
+    my $keyval = $self->{'keyval'};
+    my $value = $keyval->{$keyword};
+    my @array = ();
+    $i = 0;
+    while ($value) {
+        $value =~ s/^\s*(.*)\s*\n//;
+        $array[$i] = $1;
+        $i++;
+    }
+    return @array;
+}
+
+# returns 1 if the input file existed
+sub ok {
+    my $self = shift;
+    $self->{"ok"};
+}
+
+sub display {
+    my $self = shift;
+    my @keys = @_ ? @_ : sort keys %$self;
+    foreach $key (@keys) {
+        print "\t$key => $self->{$key}\n";
+    }
+}
+
+sub doprint {
+    my $self = shift;
+    print "QCParse:\n";
+    my $keyval = $self->{'keyval'};
+    foreach $i (keys %$keyval) {
+        my $val = $keyval->{$i};
+        $val =~ s/\n/\\n/g;
+        print "keyword = $i, value = $val\n";
+    }
+}
+
+sub error {
+    my $self = shift;
+    my $msg = shift;
+    $self->{"error"} = "$self->{'error'}$msg";
+}
+
+##########################################################################
+
+package QCInput;
+$debug = 0;
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+    my $parser = shift;
+    if ($parser eq "") {
+        $parser = new QCParse;
+    }
+    $self->{"parser"} = $parser;
+    $self->{"error"} = $parser->error();
+
+    $self->{"molecule"} = new Molecule($parser->value("molecule"));
+}
+
+sub error {
+    my $self = shift;
+    my $msg = shift;
+    $self->{"error"} = "$self->{'error'}$msg";
+}
+
+sub checkpoint {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("checkpoint");
+    my $val = $self->{"parser"}->value("checkpoint");
+    if ($val ne "" && $bval eq "") {
+        $self->error("Bad value for checkpoint: $val");
+        $bval = "0";
+    }
+    elsif ($val eq "") {
+        $bval = "1";
+    }
+    $bval;
+}
+
+sub restart {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("restart");
+    my $val = $self->{"parser"}->value("restart");
+    if ($val ne "" && $bval eq "") {
+        $self->error("Bad value for restart: $val");
+        $bval = "0";
+    }
+    elsif ($val eq "") {
+        $bval = "1";
+    }
+    $bval;
+}
+
+sub label {
+    my $self = shift;
+    $self->{"parser"}->value("label");
+}
+
+sub charge {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("charge");
+    s/^\s+//;
+    s/\s+$//;
+    s/^\+//;
+    if (/^$/) { $_ = "0"; }
+    if (! /^-?\d+$/) {
+        $self->error("Bad charge: $_ (using 0)\n");
+        $_ = "0";
+    }
+    $_;
+}
+
+sub method {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("method");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $self->error("No method given (using default).\n");
+        $_ = "SCF";
+    }
+    tr/a-z/A-Z/;
+    $_;
+}
+
+sub symmetry {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("symmetry");
+    s/^\s*//;
+    s/\s*$//;
+    uc $_;
+}
+
+sub memory {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("memory");
+    s/^\s*//;
+    s/\s*$//;
+    if ($_ eq "") {
+        $_ = 32000000;
+    }
+    $_;
+}
+
+sub state {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("state");
+    s/^\s*//;
+    s/\s*$//;
+    uc $_;
+}
+
+sub mult {
+    my $self = shift;
+    $_ = $self->state();
+    s/^\s*(\d+)/\1/;
+    if (/^\s*$/) {
+        $_ = 1;
+    }
+    $_;
+}
+
+sub basis {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("basis");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $self->error("No basis given (using default).\n");
+        $_ = "STO-3G";
+    }
+    $_;
+}
+
+sub auxbasis {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("auxbasis");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $self->error("No auxiliary basis given (using default).\n");
+        $_ = "STO-3G";
+    }
+    $_;
+}
+
+sub grid {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("grid");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $_ = "default";
+    }
+    $_;
+}
+
+sub gradient {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("gradient");
+    if ($bval eq "") {
+        my $val = $self->{"parser"}->value("gradient");
+        $self->error("Bad value for gradient: $val");
+    }
+    $bval;
+}
+
+sub fzc {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("fzc");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $_ = 0;
+    }
+    $_;
+}
+
+sub fzv {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("fzv");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "") {
+        $_ = 0;
+    }
+    $_;
+}
+
+sub docc {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("docc");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "" || $_ eq "-") {
+        $_ = "auto";
+    }
+    $_;
+}
+
+sub socc {
+    my $self = shift;
+    $_ = $self->{"parser"}->value("socc");
+    s/^\s+//;
+    s/\s+$//;
+    if ($_ eq "" || $_ eq "-") {
+        $_ = "auto";
+    }
+    $_;
+}
+
+sub optimize {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("optimize");
+    if ($bval eq "") {
+        my $val = $self->{"parser"}->value("optimize");
+        $self->error("Bad value for optimize: $val");
+    }
+    $bval;
+}
+
+# returns "" if orthog_method not set
+sub orthog_method {
+    my $self = shift;
+    my $bval = $self->{"parser"}->value("orthog_method");
+    $bval;
+}
+
+# returns "" if lindep_tol not set
+sub lindep_tol {
+    my $self = shift;
+    my $bval = $self->{"parser"}->value("lindep_tol");
+    $bval;
+}
+
+sub transition_state {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("transition_state");
+    if ($bval eq "") {
+        my $val = $self->{"parser"}->value("transition_state");
+        $self->error("Bad value for transtion_state: $val");
+    }
+    $bval;
+}
+
+sub frequencies {
+    my $self = shift;
+    my $bval = $self->{"parser"}->boolean_value("frequencies");
+    if ($bval eq "") {
+        my $val = $self->{"parser"}->value("frequencies");
+        $self->error("Bad value for frequencies: $val");
+    }
+    $bval;
+}
+
+sub axyz_lines {
+    my $self = shift;
+    $self->molecule()->string();
+}
+
+sub molecule() {
+    my $self = shift;
+    return $self->{"molecule"};
+}
+
+sub n_atom {
+    my $self = shift;
+    printf "QCInput: returning natom = %d\n", $self->{"natom"} if ($debug);
+    $self->molecule()->n_atom();
+}
+
+sub element {
+    my $self = shift;
+    $self->molecule()->element(@_);
+}
+
+sub position {
+    my $self = shift;
+    $self->molecule()->position(@_);
+}
+
+sub write_file {
+    my $self = shift;
+    my $file = shift;
+    my $parser = $self->{'parser'};
+    $parser->write_file($file);
+}
+
+sub mode_following() {
+    my $self = shift;
+    return scalar($self->{"parser"}->value_as_array("followed")) != 0;
+}
+
+# returns 1 if the input file existed
+sub ok {
+    my $self = shift;
+    $self->{"parser"}->{"ok"};
+}
+
+##########################################################################
+
+package InputWriter;
+
+# Input Writer is abstract
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize() {
+    my $self = shift;
+    my $qcinput = shift;
+    $self->{"qcinput"} = $qcinput;
+}
+
+# this should be overridden
+sub input_string() {
+    "";
+}
+
+sub write_input() {
+    my $self = shift;
+    my $file = shift;
+    my $input = $self->input_string();
+    open(OUTPUT,">$file");
+    printf OUTPUT "%s", $input;
+    close(OUTPUT);
+}
+
+sub write_qcinput {
+    my $self = shift;
+    my $file = shift;
+    my $qcinput = $self->{'qcinput'};
+    $qcinput->write_file($file);
+}
+
+##########################################################################
+
+package MPQCInputWriter;
+@ISA = qw( InputWriter );
+%methodmap = ("MP2-R12/A" => "MBPT2_R12",
+              "MP2-R12/A'" => "MBPT2_R12",
+              "MP2" => "MBPT2",
+              "OPT1[2]" => "MBPT2",
+              "OPT2[2]" => "MBPT2",
+              "ZAPT2" => "MBPT2",
+              "MP2V1" => "MBPT2",
+              "OPT1[2]V1" => "MBPT2",
+              "OPT2[2]V1" => "MBPT2",
+              "ZAPT2V1" => "MBPT2",
+              "MP2V2" => "MBPT2",
+              "OPT1[2]V2" => "MBPT2",
+              "OPT2[2]V2" => "MBPT2",
+              "ZAPT2V2" => "MBPT2",
+              "MP2V2LB" => "MBPT2",
+              "OPT1[2]V2LB" => "MBPT2",
+              "OPT2[2]V2LB" => "MBPT2",
+              "ZAPT2V2LB" => "MBPT2",
+              "ROSCF" => "SCF",
+              "SCF" => "SCF",
+              "UHF" => "UHF",
+              "CLHF" => "CLHF",
+              "HSOSHF" => "HSOSHF",
+              "HF" => "SCF",
+              "HFK" => "DFT",
+              "XALPHA" => "DFT",
+              "HFS"    => "DFT",
+              "HFB"    => "DFT",
+              "HFG96"  => "DFT",
+              "BLYP"   => "DFT",
+              "B3LYP"  => "DFT",
+              "KMLYP"  => "DFT",
+              "B3PW91"  => "DFT",
+              "PBE"    => "DFT",
+              "PW91"    => "DFT",
+              "SPZ81"   => "DFT",
+              "B3P86"   => "DFT",
+              "BP86"    => "DFT",
+              "BPW91"   => "DFT",
+              "CLHFK" => "DFT",
+              "CLXALPHA" => "DFT",
+              "CLHFS"    => "DFT",
+              "CLHFB"    => "DFT",
+              "CLHFG96"  => "DFT",
+              "CLBLYP"   => "DFT",
+              "CLB3LYP"  => "DFT",
+              "CLKMLYP"  => "DFT",
+              "CLB3PW91"  => "DFT",
+              "CLPBE"    => "DFT",
+              "CLPW91"    => "DFT",
+              "SPZ81"   => "DFT",
+              "B3P86"   => "DFT",
+              "BP86"    => "DFT",
+              "BPW91"   => "DFT",
+              "HSOSHFK" => "DFT",
+              "HSOSXALPHA" => "DFT",
+              "HSOSHFS"    => "DFT",
+              "HSOSHFB"    => "DFT",
+              "HSOSHFG96"  => "DFT",
+              "HSOSBLYP"   => "DFT",
+              "HSOSB3LYP"  => "DFT",
+              "HSOSKMLYP"  => "DFT",
+              "HSOSB3PW91"  => "DFT",
+              "HSOSPBE"    => "DFT",
+              "HSOSPW91"    => "DFT",
+              "HSOSSPZ81"   => "DFT",
+              "HSOSB3P86"   => "DFT",
+              "HSOSBP86"    => "DFT",
+              "HSOSBPW91"   => "DFT",
+              "UHFK" => "DFT",
+              "UXALPHA" => "DFT",
+              "UHFS"    => "DFT",
+              "UHFB"    => "DFT",
+              "UHFG96"  => "DFT",
+              "UBLYP"   => "DFT",
+              "UB3LYP"  => "DFT",
+              "UKMLYP"  => "DFT",
+              "UB3PW91"  => "DFT",
+              "UPBE"    => "DFT",
+              "UPW91"    => "DFT",
+              "USPZ81"   => "DFT",
+              "UB3P86"   => "DFT",
+              "UBP86"    => "DFT",
+              "UBPW91"   => "DFT",
+             );
+%mbpt2r12stdapproxmap = ("MP2-R12/A" => "A",
+                         "MP2-R12/A'" => "A'",
+                         );
+%mbpt2map = ("MP2" => "mp",
+             "OPT1[2]" => "opt1",
+             "OPT2[2]" => "opt2",
+             "ZAPT2" => "zapt",
+             "MP2V1" => "mp",
+             "OPT1[2]V1" => "opt1",
+             "OPT2[2]V1" => "opt2",
+             "ZAPT2V1" => "zapt",
+             "MP2V2" => "mp",
+             "OPT1[2]V2" => "opt1",
+             "OPT2[2]V2" => "opt2",
+             "ZAPT2V2" => "zapt",
+             "MP2V2LB" => "mp",
+             "OPT1[2]V2LB" => "opt1",
+             "OPT2[2]V2LB" => "opt2",
+             "ZAPT2V2LB" => "zapt");
+%mbpt2algmap = ("MP2" => "",
+                "OPT1[2]" => "",
+                "OPT2[2]" => "",
+                "ZAPT2" => "",
+                "MP2V1" => "v1",
+                "OPT1[2]V1" => "v1",
+                "OPT2[2]V1" => "v1",
+                "ZAPT2V1" => "v1",
+                "MP2V2" => "v2",
+                "OPT1[2]V2" => "v2",
+                "OPT2[2]V2" => "v2",
+                "ZAPT2V2" => "v2",
+                "MP2V2LB" => "v2lb",
+                "OPT1[2]V2LB" => "v2lb",
+                "OPT2[2]V2LB" => "v2lb",
+                "ZAPT2V2LB" => "v2lb");
+$debug = 0;
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize() {
+    my $self = shift;
+    my $qcinput = shift;
+    $self->{"qcinput"} = $qcinput;
+}
+
+sub docc_string() {
+    my $self = shift;
+    my $qcinput = $self->{"qcinput"};
+    my $occs = $qcinput->docc();
+    if ($occs eq "auto") { return ""; }
+    $occs =~ s/,/ /g;
+    "docc = [ $occs ]";
+}
+
+sub socc_string() {
+    my $self = shift;
+    my $qcinput = $self->{"qcinput"};
+    my $occs = $qcinput->socc();
+    if ($occs eq "auto") { return ""; }
+    $occs =~ s/,/ /g;
+    "socc = [ $occs ]";
+}
+
+sub input_string() {
+    my $self = shift;
+    my $qcinput = $self->{"qcinput"};
+    my $qcparse = $qcinput->{"parser"};
+
+    my $use_cints = 0;
+    my $do_cca = $qcparse->value("do_cca");
+
+    printf "molecule = %s\n", $qcparse->value("molecule") if ($debug);
+
+    my $symmetry = $qcinput->symmetry();
+    my $mol = "% molecule specification";
+    $mol = "$mol\nmolecule<Molecule>: (";
+    $symmetry = lc $symmetry if ($symmetry eq "AUTO");
+    if ($qcinput->frequencies()) {
+        $mol = "$mol\n  symmetry = C1";
+    }
+    else {
+        $mol = "$mol\n  symmetry = $symmetry";
+    }
+    $mol = "$mol\n  unit = angstrom";
+    $mol = "$mol\n  { atoms geometry } = {";
+    printf "MPQCInputWriter: natom = %d\n", $qcinput->n_atom() if ($debug);
+    my $i;
+    for ($i = 0; $i < $qcinput->n_atom(); $i++) {
+        $mol = sprintf "%s\n    %2s     [ %18.12f %18.12f %18.12f ]",
+                       $mol, $qcinput->element($i),
+                       $qcinput->position($i,0),
+                       $qcinput->position($i,1),
+                       $qcinput->position($i,2);
+    }
+    $mol = "$mol\n  }";
+    $mol = "$mol\n)\n";
+
+    my $basis = "% basis set specification";
+    $basis = "$basis\nbasis<GaussianBasisSet>: (";
+    $basis = sprintf "%s\n  name = \"%s\"", $basis, $qcinput->basis();
+    $basis = "$basis\n  molecule = \$:molecule";
+    $basis = "$basis\n)\n";
+
+    my $integrals = "";
+    if($do_cca) {
+      $integrals = "% using cca integrals";
+      $integrals = "$integrals\nintegrals<IntegralCCA>: (";
+      my $buffer_type = $qcparse->value("integral_buffer");
+      if( $buffer_type ne "opaque" && $buffer_type ne "array" ) {
+        $buffer_type = "opaque";
+      }     
+      my $int_package = $qcparse->value("integral_package");
+      if( $int_package ne "intv3" && $int_package ne "cints" ) {
+        $int_package = "intv3";
+      }
+      $integrals = "$integrals\n  integral_buffer = $buffer_type";
+      $integrals = "$integrals\n  integral_package = $int_package";
+      $integrals = "$integrals\n  evaluator_factory = MPQC.IntegralEvaluatorFactory";
+      $integrals = "$integrals\n  molecule = \$:molecule";
+      $integrals = "$integrals\n)\n";
+    }
+
+    my $fixed = $qcparse->value_as_arrayref("fixed");
+    my $followed = $qcparse->value_as_arrayref("followed");
+    if (scalar(@{$fixed}) != 0) {
+        $fixed = $self->mpqc_fixed_coor($fixed);
+    }
+    else {
+        $fixed = "";
+    }
+    if (scalar(@{$followed}) != 0) {
+        $followed = $self->mpqc_followed_coor($followed);
+    }
+    else {
+        $followed = "";
+    }
+
+    my $coor = "  % molecular coordinates for optimization";
+    $coor = "$coor\n  coor<SymmMolecularCoor>: (";
+    $coor = "$coor\n    molecule = \$:molecule";
+    $coor = "$coor\n    generator<IntCoorGen>: (";
+    $coor = "$coor\n      molecule = \$:molecule";
+    $coor = "$coor\n    )";
+    $coor = "$coor$followed";
+    $coor = "$coor$fixed";
+    $coor = "$coor\n  )\n";
+
+    my $charge = $qcinput->charge();
+    my $mult = $qcinput->mult();
+    my $docc = $self->docc_string();
+    my $socc = $self->socc_string();
+
+    my $grid = $qcinput->grid();
+
+    my $memory = $qcinput->memory();
+    my $inputmethod = $methodmap{uc($qcinput->method())};
+    my $method = "$inputmethod";
+    $method = "SCF" if ($method eq "");
+    my $openmethod = substr(uc($qcinput->method()),0,4);
+    if (substr($openmethod,0,2) eq "CL") { $openmethod = "CL"; }
+    if (substr($openmethod,0,1) eq "U") { $openmethod = "U"; }
+    if ($method eq "SCF") {
+        if ($openmethod eq "U") {
+            $method = "UHF";
+        }
+        elsif ($openmethod eq "CL") {
+            $method = "CLHF";
+        }
+        elsif ($openmethod eq "HSOS") {
+            $method = "HSOSHF";
+        }
+        elsif ($qcinput->mult() == 1) {
+            $method = "CLHF";
+            $openmethod = "CL";
+        }
+        else {
+            $method = "HSOSHF";
+            $openmethod = "HSOS";
+        }
+    }
+    my $functional;
+    if ($method eq "DFT") {
+        $functional = uc($qcinput->method());
+        if ($openmethod eq "U") {
+            $method = "UKS";
+            $functional = substr($functional,1);
+        }
+        elsif ($openmethod eq "CL") {
+            $method = "CLKS";
+            $functional = substr($functional,2);
+        }
+        elsif ($openmethod eq "HSOS") {
+            $method = "HSOSKS";
+            $functional = substr($functional,4);
+        }
+        elsif ($qcinput->mult() == 1) {
+            $method = "CLKS";
+            $openmethod = "CL";
+        }
+        else {
+            $method = "UKS";
+            $openmethod = "U";
+        }
+    }
+    my $orthog_method = $qcinput->orthog_method();
+    my $lindep_tol = $qcinput->lindep_tol();
+    my $mole = "  do_energy = yes";
+    if ($qcinput->gradient()) {
+        $mole = "$mole\n  do_gradient = yes";
+    }
+    else {
+        $mole = "$mole\n  do_gradient = no";
+    }
+    if($do_cca) {
+      $mole = "$mole\n  do_cca = yes";
+    }
+    $mole = "$mole\n  % method for computing the molecule's energy";
+    $mole = "$mole\n  mole<$method>: (";
+    $mole = "$mole\n    molecule = \$:molecule";
+    $mole = "$mole\n    basis = \$:basis";
+    $mole = "$mole\n    coor = \$..:coor";
+    $mole = "$mole\n    memory = $memory";
+    if($do_cca) {
+      $mole = "$mole\n    integrals = \$:integrals";
+    }
+    if ($inputmethod eq "SCF" || $inputmethod eq "UHF"
+        || $method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS") {
+        $mole = "$mole\n    total_charge = $charge";
+        $mole = "$mole\n    multiplicity = $mult";
+        $mole = "$mole\n    print_npa = yes";
+        if ($docc ne "") {$mole = "$mole\n    $docc";}
+        if ($socc ne "") {$mole = "$mole\n    $socc";}
+        if ($orthog_method ne "" ) {
+            $mole = "$mole\n    orthog_method = $orthog_method";
+        }
+        if ($lindep_tol ne "" ) {
+            $mole = "$mole\n    lindep_tol = $lindep_tol";
+        }
+    }
+    if ($method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS") {
+        $mole = "$mole\n    functional<StdDenFunctional>: name = \"$functional\"";
+    }
+    if (($method eq "CLKS" || $method eq "UKS" || $method eq "HSOSKS")
+        && $grid ne "default") {
+        $mole = "$mole\n    integrator<RadialAngularIntegrator>: (grid = $grid)";
+    }
+    if ($method eq "MBPT2_R12") {
+        my $stdapprox = $mbpt2r12stdapproxmap{uc($qcinput->method())};
+        my $auxbasis = $qcinput->auxbasis();
+        my $fzc = $qcinput->fzc();
+
+        $mole = sprintf "%s\n    stdapprox = \"%s\"", $mole, $stdapprox;
+        $mole = "$mole\n    integrals<IntegralCints>: ()";
+        $mole = "$mole\n    nfzc = $fzc";
+        # don't write an auxbasis if the auxbasis is the same as the basis set.
+        # this will speed up the calculation
+        if ("$auxbasis" ne "" && "$auxbasis" ne $qcinput->basis()) {
+            $mole = "$mole\n    aux_basis<GaussianBasisSet>: (";
+            $mole = sprintf "%s\n      name = \"%s\"", $mole, $auxbasis;
+            $mole = "$mole\n      molecule = \$:molecule";
+            $mole = "$mole\n    )\n";
+        }
+        $mole = append_reference($mole,"CLHF",$charge,$mult,$memory,$orthog_method,
+                                 $lindep_tol,$docc,$socc,"DZ (Dunning)");
+        $use_cints = 1;
+    }
+    elsif ($method eq "MBPT2") {
+        my $fzc = $qcinput->fzc();
+        my $fzv = $qcinput->fzv();
+        my $mbpt2method = $mbpt2map{uc($qcinput->method())};
+        my $mbpt2algorithm = $mbpt2algmap{uc($qcinput->method())};
+        $mole = "$mole\n    method = $mbpt2method";
+        if ($mbpt2algorithm ne "") {
+            $mole = "$mole\n    algorithm = $mbpt2algorithm";
+        }
+        $mole = "$mole\n    nfzc = $fzc";
+        $mole = "$mole\n    nfzv = $fzv";
+        my $refmethod = "";
+        if ($qcinput->mult() == 1) {
+            $refmethod = "CLHF";
+        }
+        else {
+            $refmethod = "HSOSHF";
+        }
+        $mole = append_reference($mole,$refmethod,$charge,$mult,$memory,$orthog_method,
+                                 $lindep_tol,$docc,$socc,"STO-3G");
+    }
+    elsif (! ($basis =~ /^STO/
+              || $basis =~ /^MI/
+              || $basis =~ /^\d-\d1G$/) && ! $do_cca ) {
+        my $guessmethod = "${openmethod}HF";
+        $mole = "$mole\n    guess_wavefunction<$guessmethod>: (";
+        $mole = "$mole\n      molecule = \$:molecule";
+        $mole = "$mole\n      total_charge = $charge";
+        $mole = "$mole\n      multiplicity = $mult";
+        if ($docc ne "") {$mole = "$mole\n      $docc";}
+        if ($socc ne "") {$mole = "$mole\n      $socc";}
+        $mole = "$mole\n      basis<GaussianBasisSet>: (";
+        $mole = "$mole\n        molecule = \$:molecule";
+        $mole = "$mole\n        name = \"STO-3G\"";
+        $mole = "$mole\n      )";
+        $mole = "$mole\n      memory = $memory";
+        if($do_cca) {
+          $mole = "$mole\n      integrals = \$:integrals";
+        }
+        $mole = "$mole\n    )";
+    }
+    if ($qcinput->frequencies()) {
+        $mole = "$mole\n    hessian<FinDispMolecularHessian>: (";
+        if ($symmetry ne "C1") {
+            $mole="$mole\n      point_group<PointGroup>: symmetry = $symmetry";
+        }
+        $mole = "$mole\n      checkpoint = no";
+        $mole = "$mole\n      restart = no";
+        $mole = "$mole\n    )";
+    }
+    $mole = "$mole\n  )\n";
+
+    my $opt;
+    if ($qcinput->optimize()) {
+        $opt = "  optimize = yes";
+    }
+    else {
+        $opt = "  optimize = no";
+    }
+    my $optclass, $updateclass;
+    if ($qcinput->transition_state()) {
+        $optclass = "EFCOpt";
+        $updateclass = "PowellUpdate";
+    }
+    else {
+        $optclass = "QNewtonOpt";
+        $updateclass = "BFGSUpdate";
+    }
+    $opt = "$opt\n  % optimizer object for the molecular geometry";
+    $opt = "$opt\n  opt<$optclass>: (";
+    $opt = "$opt\n    max_iterations = 20";
+    $opt = "$opt\n    function = \$..:mole";
+    if ($qcinput->transition_state()) {
+        $opt = "$opt\n    transition_state = yes";
+        if ($qcinput->mode_following()) {
+            $opt = "$opt\n    hessian = [ [ -0.1 ] ]";
+            $opt = "$opt\n    mode_following = yes";
+        }
+    }
+    $opt = "$opt\n    update<$updateclass>: ()";
+    $opt = "$opt\n    convergence<MolEnergyConvergence>: (";
+    $opt = "$opt\n      cartesian = yes";
+    $opt = "$opt\n      energy = \$..:..:mole";
+    $opt = "$opt\n    )";
+    $opt = "$opt\n  )\n";
+
+    my $freq = "";
+    if ($qcinput->frequencies()) {
+        $freq = "% vibrational frequency input";
+        $freq = "$freq\n  freq<MolecularFrequencies>: (";
+        if ($symmetry ne "C1") {
+            $freq = "$freq\n    point_group<PointGroup>: symmetry = $symmetry";
+        }
+        $freq = "$freq\n    molecule = \$:molecule";
+        $freq = "$freq\n  )\n";
+    }
+
+    my $mpqcstart = sprintf ("mpqc: (\n  checkpoint = %s\n",
+                             bool_to_yesno($qcinput->checkpoint()));
+    $mpqcstart = sprintf ("%s  savestate = %s\n",
+                          $mpqcstart,bool_to_yesno($qcinput->checkpoint()));
+    $mpqcstart = sprintf ("%s  restart = %s\n",
+                          $mpqcstart,bool_to_yesno($qcinput->restart()));
+    if ($use_cints) {
+        $mpqcstart = "$mpqcstart  integrals<IntegralCints>: ()\n";
+    }
+    my $mpqcstop = ")\n";
+    my $emacs = "% Emacs should use -*- KeyVal -*- mode\n";
+    my $warn = "% this file was automatically generated\n";
+    my $lab = $qcinput->label();
+    my $label = "";
+    if (! $lab =~ /^\s*$/) {
+        $label = "% label: $lab";
+        $label =~ s/\n/\n% label: /g;
+        $label = "$label\n";
+    }
+    "$emacs$warn$label$mol$basis$integrals$mpqcstart$coor$mole$opt$freq$mpqcstop";
+}
+
+sub mpqc_fixed_coor {
+    my $self = shift;
+    my $coorref = shift;
+    my $result = "";
+    $result = "\n    fixed<SetIntCoor>: [";
+    while (scalar(@{$coorref}) != 0) {
+        my $nextcoor = $self->mpqc_sum_coor("      ","",$coorref);
+        $result = "$result\n$nextcoor";
+    }
+    $result = "$result\n    ]";
+}
+
+sub mpqc_followed_coor {
+    my $self = shift;
+    my $coorref = shift;
+    sprintf "\n%s", $self->mpqc_sum_coor("    ","followed",$coorref);
+}
+
+sub mpqc_sum_coor {
+    my $self = shift;
+    my $space = shift;
+    my $name = shift;
+    my $coor = shift;
+    my $result = "$space$name<SumIntCoor>:(";
+    $result = "$result\n$space  coor: [";
+    my @coef = ();
+    do {
+        $coef[$ncoor] = shift @{$coor};
+        my $simple = $self->mpqc_coor($coor);
+        $result = "$result\n$space    $simple";
+        $ncoor = $ncoor + 1;
+    } while($coor->[0] eq "+" && shift @{$coor} eq "+");
+    $result = "$result\n$space  ]";
+    $result = "$result\n$space  coef = [";
+    my $i;
+    foreach $i (0..$#coef) {
+        $result = "$result $coef[$i]";
+    }
+    $result = "$result]";
+    $result = "$result\n$space)";
+    $result;
+}
+
+sub mpqc_coor {
+    my $self = shift;
+    my $coor = shift;
+    my $type = shift @{$coor};
+    if ($type eq "TORS") {
+        return sprintf "<TorsSimpleCo>:(atoms = [%d %d %d %d])",
+                       shift @{$coor},shift @{$coor},
+                       shift @{$coor},shift @{$coor};
+    }
+    if ($type eq "BEND") {
+        return sprintf "<BendSimpleCo>:(atoms = [%d %d %d])",
+                       shift @{$coor},shift @{$coor},
+                       shift @{$coor};
+    }
+    if ($type eq "STRE") {
+        return sprintf "<StreSimpleCo>:(atoms = [%d %d])",
+                        shift @{$coor},shift @{$coor};
+    }
+}
+
+sub bool_to_yesno {
+    if (shift) { return "yes"; }
+    else { return "no"; }
+}
+
+sub append_reference {
+    my $mole = shift;
+    my $refmethod = shift;
+    my $charge = shift;
+    my $mult = shift;
+    my $memory = shift;
+    my $orthog_method = shift;
+    my $lindep_tol = shift;
+    my $docc = shift;
+    my $socc = shift;
+    my $guessbasis = shift;
+    $mole = "$mole\n    reference<$refmethod>: (";
+    $mole = "$mole\n      molecule = \$:molecule";
+    $mole = "$mole\n      basis = \$:basis";
+    $mole = "$mole\n      total_charge = $charge";
+    $mole = "$mole\n      multiplicity = $mult";
+    $mole = "$mole\n      memory = $memory";
+    if ($orthog_method ne "" ) {
+        $mole = "$mole\n      orthog_method = $orthog_method";
+    }
+    if ($lindep_tol ne "" ) {
+        $mole = "$mole\n      lindep_tol = $lindep_tol";
+    }
+    if ($docc ne "") {$mole = "$mole\n      $docc";}
+    if ($socc ne "") {$mole = "$mole\n      $socc";}
+    if (! ($basis =~ /^STO/
+           || $basis =~ /^MI/
+           || $basis =~ /^\d-\d1G$/)) {
+        $mole = "$mole\n      guess_wavefunction<$refmethod>: (";
+        $mole = "$mole\n        molecule = \$:molecule";
+        $mole = "$mole\n        total_charge = $charge";
+        $mole = "$mole\n        multiplicity = $mult";
+        if ($docc ne "") {$mole = "$mole\n        $docc";}
+        if ($socc ne "") {$mole = "$mole\n        $socc";}
+        $mole = "$mole\n        basis<GaussianBasisSet>: (";
+        $mole = "$mole\n          molecule = \$:molecule";
+        $mole = "$mole\n          name = \"$guessbasis\"";
+        $mole = "$mole\n        )";
+        $mole = "$mole\n        memory = $memory";
+        $mole = "$mole\n      )";
+    }
+    $mole = "$mole\n    )";
+    return $mole;
+}
+
+1;
+
Index: ThirdParty/mpqc_open/lib/perl/QCResult.pm
===================================================================
--- ThirdParty/mpqc_open/lib/perl/QCResult.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
+++ ThirdParty/mpqc_open/lib/perl/QCResult.pm	(revision 860145ede0311b3539e419277d8ca4a3091f40b2)
@@ -0,0 +1,692 @@
+#
+eval 'exec perl $0 $*'
+    if 0;
+
+require QCParse;
+require Molecule;
+
+##########################################################################
+
+package QCResult;
+
+# this seems to not work as expected
+$fltrx = "([-+]?(?:\\d+\\.\\d*|\\d+|\\.\\d+)(?:[eEdD][+-]?\\d+)?)";
+
+$have_nodenum = 0;
+
+sub test {
+    my $i;
+    @nums = ("-1.0", "-1", "-.1", "-1.",
+             "-1.0E-04", "-1E-04", "-.1E-04", "-1.E-04");
+    foreach $i ( @nums ) {
+        $i =~ /$fltrx/;
+        my $flt = $1;
+        printf "%10s %10s %12.8f\n", $i, $flt, $flt;
+    }
+}
+
+sub new {
+    my $this = shift;
+    my $class = ref($this) || $this;
+    my $self = {};
+    bless $self, $class;
+    $self->initialize(@_);
+    return $self;
+}
+
+sub initialize {
+    my $self = shift;
+    my $infile = shift;
+    my $outfile = shift;
+    my $parse = new QCParse();
+    $parse->parse_file($infile);
+    $self->{"qcinput"} = QCInput::new QCInput($parse);
+
+    $self->{"exists"} = 0;
+    $self->{"ok"} = 0;
+    if (-e $outfile) {
+        $self->{"exists"} = 1;
+        open(OUTFILE,"<$outfile");
+        my $first = 1;
+        while (<OUTFILE>) {
+            if ($first && /^ *[0-9]+:/) {
+                $have_nodenum = 1;
+                $first = 0;
+            }
+            s/^ *[0-9]+:// if ($have_nodenum);
+            if (/^\s*MPQC:/) {
+                $self->parse_mpqc(\*OUTFILE);
+                last;
+            }
+            elsif (/^\s*Entering Gaussian System/) {
+                $self->parse_g94(\*OUTFILE);
+                last;
+            }
+        }
+        close(OUTFILE);
+    }
+}
+
+sub parse_g94 {
+    my $self = shift;
+    my $out = shift;
+    my $scfenergy = "";
+    my $mp2energy = "";
+    my $ccsd_tenergy = "";
+    my $t1norm = "";
+    my $optconverged = 0;
+    my $molecule;
+    my $havefreq = 0;
+    my $freq = [];
+    my $ifreq = 0;
+    my $b3pw91energy = "";
+    while (<$out>) {
+        s/^ *[0-9]+:// if ($have_nodenum);
+        if (/^\s*SCF Done:  E\(RHF\) =\s*$fltrx\s/) {
+            $scfenergy = $1;
+        }
+        elsif (/^\s*E2\s*=\s*$fltrx\s*EUMP2\s*=\s*$fltrx/) {
+            $mp2energy = $2;
+        }
+        elsif (/^\s*CCSD\(T\)\s*=\s*$fltrx/) {
+            $ccsd_tenergy = $1;
+            $ccsd_tenergy =~ s/[DdE]/e/;
+        }
+        elsif (/E\(RB\+HF-PW91\)\s*=\s*$fltrx/) {
+            $b3pw91energy = $1;
+            $b3pw91energy =~ s/[DdE]/e/;
+        }
+        elsif (/^\s*T1 Diagnostic\s*=\s*$fltrx/) {
+            $t1norm = $1;
+        }
+        elsif (/^\s*-- Stationary point found./
+               || /CONVERGENCE CRITERIA APPARENTLY SATISFIED/) {
+            $optconverged = 1;
+        }
+        elsif ($optconverged
+               && (/^\s*Input orientation:/ || /^\s*Standard orientation:/)) {
+            <$out>; <$out>; <$out>; <$out>;
+            my $molstr = "";
+            while (<$out>) {
+                if (! /^\s+\d+\s+(\d+)\s+$fltrx\s+$fltrx\s+$fltrx\s+/) {
+                    last;
+                }
+                $molstr = "${molstr} $1 $2 $3 $4\n";
+            }
+            $molecule = new Molecule($molstr);
+        }
+        elsif (/^\s*Frequencies --\s+$fltrx\s+$fltrx\s+$fltrx/) {
+            $freq->[$ifreq] = $1; $ifreq++;
+            $freq->[$ifreq] = $2; $ifreq++;
+            $freq->[$ifreq] = $3; $ifreq++;
+            $havefreq = 1;
+        }
+        elsif (/^\s*Frequencies --\s+$fltrx\s+$fltrx/) {
+            $freq->[$ifreq] = $1; $ifreq++;
+            $freq->[$ifreq] = $2; $ifreq++;
+            $havefreq = 1;
+        }
+        elsif (/^\s*Frequencies --\s+$fltrx\s/) {
+            $freq->[$ifreq] = $1; $ifreq++;
+            $havefreq = 1;
+        }
+    }
+    $self->{"scfenergy"} = $scfenergy;
+    $self->{"mp2energy"} = $mp2energy;
+    $self->{"b3pw91energy"} = $b3pw91energy;
+    $self->{"optconverged"} = $optconverged;
+    if ($optconverged) {
+        $self->{"optmolecule"} = $molecule;
+    }
+    $self->{"have_frequencies"} = $havefreq;
+    if ($havefreq) {
+        my @tmp = sort(@$freq);
+        $freq = \@tmp;
+        $self->{"freq"} = $freq;
+    }
+    my $qcinput = $self->{"qcinput"};
+    my $method = $qcinput->method();
+    if ($method eq "MP2") {
+        $self->{"energy"} = $mp2energy;
+    }
+    elsif ($method eq "CCSD(T)") {
+        $self->{"energy"} = $ccsd_tenergy;
+    }
+    elsif ($method eq "B3-PW91") {
+        $self->{"energy"} = $b3pw91energy;
+    }
+    elsif ($method eq "SCF"
+           || $method eq "ROSCF") {
+        $self->{"energy"} = $scfenergy;
+    }
+
+    $self->{"t1norm"} = $t1norm;
+
+    $self->{"ok"} = 0;
+    if ($self->{"energy"} ne "") {
+        if ($qcinput->optimize()) {
+            if ($self->{"optconverged"}) {
+                $self->{"ok"} = 1
+                }
+            else {
+                #printf "not ok because not converged\n";
+            }
+        }
+        else {
+            $self->{"ok"} = 1;
+        }
+        if ($qcinput->frequencies() && ! $havefreq) {
+            $self->{"ok"} = 0;
+            #printf "not ok because no freq\n";
+        }
+    }
+    else {
+        #printf "not ok because no energy\n";
+    }
+}
+
+sub parse_mpqc {
+    my $self = shift;
+    my $out = shift;
+    my $scfenergy = "";
+    my $opt1energy = "";
+    my $opt2energy = "";
+    my $zapt2energy = "";
+    my $mp2energy = "";
+    my $optconverged = 0;
+    my $molecule;
+    my $havefreq = 0;
+    my $freq;
+    my $wante = 1;
+    my $ifreq = 0;
+    my $grad;
+    my $ngrad;
+    my $state = "none";
+    my $molecularenergy = "";
+    my $s2norm = "";
+    my $to_angstrom = 0.52917706;
+    my $geometry_conversion = $to_angstrom;
+    my $error = "";
+    my $ccsd_energy = "";
+    my $ccsd_t_energy = "";
+    my $t1norm = "";
+    my $t12norm = "";
+    my $psioutput = "";
+    while (<$out>) {
+        s/^ *[0-9]+:// if ($have_nodenum);
+        if ($state eq "read grad" && $wante) {
+            if (/^\s*([0-9]+\s+[A-Za-z]+\s+)?([\-\.0-9]+)\s+([\-\.0-9]+)\s+([\-\.0-9]+)/) {
+                $grad->[$ngrad + 0] = $2;
+                $grad->[$ngrad + 1] = $3;
+                $grad->[$ngrad + 2] = $4;
+                $ngrad = $ngrad + 3;
+            }
+            elsif (/^\s*([0-9]+\s+)?([\-\.0-9]+)/) {
+                $grad->[$ngrad] = $2;
+                $ngrad = $ngrad + 1;
+            }
+            else {
+                $self->{"grad"} = $grad;
+                $state = "none";
+            }
+        }
+
+        if ($wante && /total scf energy =\s+$fltrx/) {
+            $scfenergy = $1;
+        }
+        elsif ($wante && /OPT1 energy .*:\s+$fltrx/) {
+            $opt1energy = $1;
+        }
+        elsif ($wante && /OPT2 energy .*:\s+$fltrx/) {
+            $opt2energy = $1;
+        }
+        elsif ($wante && /ZAPT2 energy .*:\s+$fltrx/) {
+            $zapt2energy = $1;
+        }
+        elsif ($wante && /MP2 energy .*:\s+$fltrx/) {
+            $mp2energy = $1;
+        }
+        elsif (/CSCF: An SCF program written in C/) {
+            $psioutput = 1;
+        }
+        elsif ($psioutput && $wante
+               && /total energy       =\s*$fltrx/) {
+            $scfenergy = $1;
+        }
+        elsif ($psioutput && $wante
+               && /\s*ITER\s+CORRELATION ENERGY\s+T1 2-NORM\s+T1 DIAG/) {
+            # this is a PSI CC output embedded in an MPQC output
+            # grab iteration 0
+            my $ecorr;
+            <$out>;
+            while (<$out>) {
+                if (/^\s*\d+\s+$fltrx\s+$fltrx\s+$fltrx/) {
+                    $ecorr = $1;
+                    $t12norm = $2;
+                    $t1norm = $3;
+                }
+                else {
+                    last;
+                }
+            }
+            $ccsd_energy = $ecorr + $scfenergy;
+        }
+# IMBN: Added the following elsif to handle new cc output 
+        elsif ($psioutput && $wante
+               && /\s*ITER\s+CORRELATION ENERGY\s+T1 DIAG\s+D1\(CCSD\)/) {
+            # this is a PSI CC output embedded in an MPQC output
+            # grab iteration 0
+            my $ecorr;
+            <$out>;
+            while (<$out>) {
+                if (/^\s*\d+\s+$fltrx\s+$fltrx\s+$fltrx/) {
+                    $ecorr = $1;
+                    $t1norm = $2;
+                    $t12norm = $3;
+                }
+                else {
+                    last;
+                }
+            }
+            $ccsd_energy = $ecorr + $scfenergy;
+        }
+        elsif ($wante && /Value of the MolecularEnergy:\s+$fltrx/) {
+            $molecularenergy = $1;
+        }
+        elsif ($wante && /S2\(ov\) matrix 1-norm\s*=\s*$fltrx/) {
+            $self->{"s2matrix1norm"} = $1;
+        }
+        elsif ($wante && /There are degenerate orbitals within an irrep/) {
+            $self->{"degenerate"} = 1;
+        }
+        elsif ($wante && /D1\(MP2\)\s*=\s*$fltrx/) {
+            $self->{"d1mp2"} = $1;
+        }
+        elsif ($wante && /D2\(MP1\)\s*=\s*$fltrx/) {
+            $self->{"d2mp1"} = $1;
+        }
+        elsif ($wante && /S2\(ov\) matrix inf-norm\s*=\s*$fltrx/) {
+            $self->{"s2matrixinfnorm"} = $1;
+        }
+        elsif ($wante && /S2 norm\s*=\s*$fltrx/) {
+            $s2norm = $1;
+        }
+        elsif ($wante && /Largest S2 values.*:/) {
+            my $s2large_coef = [];
+            my $s2large_i = [];
+            my $s2large_a = [];
+            my $is2 = 0;
+            while (<$out>) {
+                if (/^\s*([0-9]+)\s+$fltrx\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s+->\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s*$/) {
+                    $s2large_coef->[$is2] = $2;
+                    $s2large_i->[$is2] = "$3$4";
+                    $s2large_a->[$is2] = "$5$6";
+                    $is2 = $is2 + 1;
+                }
+                else {
+                    last;
+                }
+            }
+            $self->{"s2large_coef"} = $s2large_coef;
+            $self->{"s2large_i"} = $s2large_i;
+            $self->{"s2large_a"} = $s2large_a;
+        }
+        elsif ($wante && /Largest first order coefficients.*:/) {
+            my $d1large_coef = [];
+            my $d1large_i = [];
+            my $d1large_j = [];
+            my $d1large_a = [];
+            my $d1large_b = [];
+            my $d1large_spin = [];
+            my $id1 = 0;
+            while (<$out>) {
+                if (/^\s*([0-9]+)\s+$fltrx\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s+->\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s+([0-9]+)\s+([A-Za-z0-9\'\"]+)\s+\((....)\)\s*$/) {
+                    $d1large_coef->[$id1] = $2;
+                    $d1large_i->[$id1] = "$3$4";
+                    $d1large_j->[$id1] = "$5$6";
+                    $d1large_a->[$id1] = "$7$8";
+                    $d1large_b->[$id1] = "$9$10";
+                    $d1large_spin->[$id1] = $11;
+                    $id1 = $id1 + 1;
+                }
+                else {
+                    last;
+                }
+            }
+            $self->{"d1large_coef"} = $d1large_coef;
+            $self->{"d1large_i"} = $d1large_i;
+            $self->{"d1large_j"} = $d1large_j;
+            $self->{"d1large_a"} = $d1large_a;
+            $self->{"d1large_b"} = $d1large_b;
+            $self->{"d1large_spin"} = $d1large_spin;
+        }
+        elsif ($wante && /^\s*Natural\s+Population\s+Analysis:\s*$/) {
+            my $npacharge = [];
+            my $npashellpop = [];
+            <$out>;
+            my $iatom = 0;
+            while (<$out>) {
+                if (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2 ];
+                }
+                elsif (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2, $3 ];
+                }
+                elsif (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2, $3, $4 ];
+                }
+                elsif (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2, $3, $4, $5 ];
+                }
+                elsif (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2, $3, $4, $5, $6 ];
+                }
+                elsif (/^\s*\d+\s+[A-Za-z]+\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s+$fltrx\s*$/) {
+                    $npacharge->[$iatom] = $1;
+                    $npashellpop->[$iatom] = [ $2, $3, $4, $5, $6, $7 ];
+                }
+                elsif (/^\s*$/) {
+                    last;
+                }
+                else {
+                    die "AM too high to read NPA shell populations (line=$_)";
+                }
+                $iatom = $iatom + 1;
+            }
+            $self->{"npacharge"} = $npacharge;
+            $self->{"npashellpop"} = $npashellpop;
+        }
+        elsif (/The optimization has converged/) {
+            $optconverged = 1;
+        }
+        elsif (/^\s*Beginning displacement/) {
+            # don't grab energies for displaced geometries
+            $wante = 0;
+        }
+        elsif ($optconverged
+               && /^\s+n\s+atom\s+label\s+x\s+y\s+z\s+mass\s*$/) {
+            # old style geometry
+            my $molstr = "";
+            while (<$out>) {
+                s/^ *[0-9]+:// if ($have_nodenum);
+                if (! /^\s+\d+\s+(\w+)\s+$fltrx\s+$fltrx\s+$fltrx\s+/
+                 && ! /^\s+\d+\s+(\w+)\s+\S+\s+$fltrx\s+$fltrx\s+$fltrx\s+/) {
+                    last;
+                }
+                $molstr = sprintf "%s %s %16.14f %16.14f %16.14f\n",
+                    ${molstr}, $1, $2 * $to_angstrom,
+                    $3 * $to_angstrom, $4 * $to_angstrom;
+            }
+            $molecule = new Molecule($molstr);
+        }
+        elsif (/^\s*molecule<Molecule>:/) {
+            $geometry_conversion = $to_angstrom;
+        }
+        elsif (/^\s*unit = \"angstrom\"\s*/) { # "
+            $geometry_conversion = 1.0;
+        }
+        elsif ($optconverged
+               && /n\s+atoms\s+(atom_labels\s+)?geometry/) {
+            # new style geometry
+            my $molstr = "";
+            while (<$out>) {
+                s/^ *[0-9]+:// if ($have_nodenum);
+                if (! /^\s+\d+\s+(\w+)\s+\[\s*$fltrx\s+$fltrx\s+$fltrx\s*\]/
+                 && ! /^\s+\d+\s+(\w+)\s+\"[^\"]*\"+\s+\[\s*$fltrx\s+$fltrx\s+$fltrx\s*\]/) { # " (unconfuse emacs)
+                    last;
+                }
+                $molstr = sprintf "%s %s %16.14f %16.14f %16.14f\n",
+                    ${molstr}, $1, $2 * $geometry_conversion,
+                    $3 * $geometry_conversion, $4 * $geometry_conversion;
+            }
+            $molecule = new Molecule($molstr);
+        }
+        elsif (/^\s+Total (MP2 )?[Gg]radient/) {
+            $state = "read grad";
+            $grad = [];
+            $ngrad = 0;
+        }
+        elsif (/^\s+Frequencies .*:\s*$/) {
+            # read the frequencies
+            while (<$out>) {
+                s/^ *[0-9]+:// if ($have_nodenum);
+                if (/^\s+\d+\s+$fltrx\s*$/) {
+                    $freq->[$ifreq] = $1;
+                    $ifreq++;
+                }
+                elsif (/THERMODYNAMIC ANALYSIS/) {
+                    last;
+                }
+            }
+            $havefreq = 1;
+        }
+    }
+    $self->{"t1norm"} = $t1norm;
+    $self->{"t1matrix2norm"} = $t12norm;
+    $self->{"s2norm"} = $s2norm;
+    $self->{"ccsdenergy"} = $ccsd_energy;
+    $self->{"ccsd_tenergy"} = $ccsd_t_energy;
+    $self->{"scfenergy"} = $scfenergy;
+    $self->{"opt1energy"} = $opt1energy;
+    $self->{"opt2energy"} = $opt2energy;
+    $self->{"zapt2energy"} = $zapt2energy;
+    $self->{"mp2energy"} = $mp2energy;
+    $self->{"optconverged"} = $optconverged;
+    $self->{"molecularenergy"} = $molecularenergy;
+    if ($optconverged) {
+        $self->{"optmolecule"} = $molecule;
+    }
+    $self->{"have_frequencies"} = $havefreq;
+    if ($havefreq) {
+        my @tmp = sort(@$freq);
+        $freq = \@tmp;
+        $self->{"freq"} = $freq;
+    }
+    my $qcinput = $self->{"qcinput"};
+    my $method = $qcinput->method();
+    #printf "qcinput ok = %d\n", $qcinput->ok();
+    if ($method eq "MP2") {
+        if ($mp2energy ne "") {
+            $self->{"energy"} = $mp2energy;
+        }
+        elsif ($qcinput->mult() == 1) {
+            $self->{"energy"} = $zapt2energy;
+        }
+    }
+    elsif ($method eq "OPT1[2]") {
+        $self->{"energy"} = $opt1energy;
+    }
+    elsif ($method eq "OPT2[2]") {
+        $self->{"energy"} = $opt2energy;
+    }
+    elsif ($method eq "CCSD") {
+        $self->{"energy"} = $ccsd_energy;
+    }
+    elsif ($method eq "CCSD(T)") {
+        $self->{"energy"} = $ccsd_t_energy;
+    }
+    elsif ($qcinput->ok()
+           && ($method eq "SCF" || $method eq "ROSCF")) {
+        $self->{"energy"} = $scfenergy;
+    }
+
+    if ($self->{"energy"} eq "") {
+        $self->{"energy"} = $self->{"molecularenergy"};
+    }
+
+    $self->{"ok"} = 0;
+    if ($self->{"energy"} ne "") {
+        if ($qcinput->optimize()) {
+            if ($self->{"optconverged"}) {
+                $self->{"ok"} = 1;
+            }
+            else {
+                $error = "$error geom not converged\n";
+            }
+        }
+        else {
+            $self->{"ok"} = 1;
+        }
+        if ($qcinput->frequencies() && ! $havefreq) {
+            $self->{"ok"} = 0;
+            $error = "$error no freq\n";
+        }
+    }
+    else {
+        $error = "$error no energy\n";
+    }
+
+    if (! $qcinput->ok()) {
+        $self->{"ok"} = 0;
+        $error = "$error qcinput error: $qcinput->{'error'}";
+    }
+
+    $self->{"error"} = $error;
+}
+
+sub ok {
+    my $self = shift;
+    $self->{"ok"}
+}
+
+sub error {
+    my $self = shift;
+    $self->{"error"};
+}
+
+sub exists {
+    my $self = shift;
+    $self->{"exists"}
+}
+
+sub input {
+    my $self = shift;
+    $self->{"qcinput"}
+}
+
+sub inputok {
+    my $self = shift;
+    $self->{"qcinput"}->ok();
+}
+
+sub energy {
+    my $self = shift;
+    $self->{"energy"}
+}
+
+sub s2norm {
+    my $self = shift;
+    $self->{"s2norm"}
+}
+
+sub s2matrix1norm {
+    my $self = shift;
+    $self->{"s2matrix1norm"}
+}
+
+sub degenerate {
+    my $self = shift;
+    $self->{"degenerate"}
+}
+
+sub d1mp2 {
+    my $self = shift;
+    $self->{"d1mp2"}
+}
+
+sub d2mp1 {
+    my $self = shift;
+    $self->{"d2mp1"}
+}
+
+sub s2matrixinfnorm {
+    my $self = shift;
+    $self->{"s2matrixinfnorm"}
+}
+
+sub t1norm {
+    my $self = shift;
+    $self->{"t1norm"}
+}
+
+sub t1matrix2norm {
+    my $self = shift;
+    $self->{"t1matrix2norm"}
+}
+
+sub optmolecule {
+    my $self = shift;
+    $self->{"optmolecule"};
+}
+
+sub gradient {
+    my $self = shift;
+    $self->{"grad"};
+}
+
+sub frequencies {
+    my $self = shift;
+    $self->{"freq"}
+}
+
+sub d1large_coef {
+    my $self = shift;
+    $self->{"d1large_coef"}
+}
+
+sub d1large_i {
+    my $self = shift;
+    $self->{"d1large_i"}
+}
+
+sub d1large_j {
+    my $self = shift;
+    $self->{"d1large_j"}
+}
+
+sub d1large_a {
+    my $self = shift;
+    $self->{"d1large_a"}
+}
+
+sub d1large_b {
+    my $self = shift;
+    $self->{"d1large_b"}
+}
+
+sub d1large_spin {
+    my $self = shift;
+    $self->{"d1large_spin"}
+}
+
+sub s2large_coef {
+    my $self = shift;
+    $self->{"s2large_coef"}
+}
+
+sub s2large_i {
+    my $self = shift;
+    $self->{"s2large_i"}
+}
+
+sub s2large_a {
+    my $self = shift;
+    $self->{"s2large_a"}
+}
+
+sub npacharge {
+    my $self = shift;
+    $self->{"npacharge"}
+}
+
+sub npashellpop {
+    my $self = shift;
+    $self->{"npashellpop"}
+}
+
+1;
