How to get potential energy surface from GAMESS?

Question

I'm trying to learn GAMESS and performing a simple PES scan of oxygen molecule. After the completion of the calculations I get a .log file with lots of sections like that:

 ---- SURFACE MAPPING GEOMETRY ----
 COORD 1= 1.500 COORD 2= 0.000
 HAS ENERGY VALUE    -149.438497
 O      0.00000   0.00000   1.25000
 O      0.00000   0.00000   0.25000
 ----------------------------------

So, is it possible to get PES scan result as one table with all of coordinates?

Here is my input:

 $BASIS GBASIS=N31 NGAUSS=6 $END
 $CONTRL SCFTYP=UHF RUNTYP=SURFACE MAXIT=50 MULT=3 COORD=ZMT NZVAR=0 $END
 $SYSTEM MWORDS=937 MEMDDI=937 PARALL=.TRUE. $END
 $SURF IVEC1(1)=1,2 IGRP1(1)=1 ORIG1=1 DISP1=0.5 NDISP1=10 $END

 $DATA 
O2
C1
O  
O   1   R12

R12 = 0.5
 $END
 $ZMAT IZMAT(1)=1,2,1,  $END

Answer 1

I know it is easy to get the PES tabulated very neatly in MOLPRO with a command like:

{table,r,scf,ccsd,ccsd_t 
head, R,HF-SCF,CCSD,CCSD(T) 
sort,1,2,3} 

which for a diatomic molecule like yours, gives the following output:

  R      HF-SCF        CCSD       CCSD(T)
 1.5 -108.3566620 -108.6007993 -108.6060512
 1.6 -108.6053845 -108.8602358 -108.8662569
 1.7 -108.7675654 -109.0332410 -109.0401624
 1.8 -108.8668845 -109.1435633 -109.1515330
 1.9 -108.9206732 -109.2086077 -109.2177903
 2.0 -108.9417215 -109.2412265 -109.2518027
 2.1 -108.9395089 -109.2509444 -109.2631091
 2.2 -108.9210615 -109.2448104 -109.2587715 

Someone that uses GAMESS more than me might give an analogous solution for GAMESS.

However you might want a very quick solution so that you can get back to working on your project, so below is a solution you can use immediately.

The command:

grep "COORD 1=" gamess.log 

would give the following type of following output:

COORD 1= 1.5 COORD 2= 0.000
COORD 1= 1.6 COORD 2= 0.000
COORD 1= 1.7 COORD 2= 0.000
COORD 1= 1.8 COORD 2= 0.000
COORD 1= 1.9 COORD 2= 0.000
COORD 1= 2.0 COORD 2= 0.000
COORD 1= 2.1 COORD 2= 0.000
COORD 1= 2.2 COORD 2= 0.000

and the command:

grep "HAS ENERGY VALUE" gamess.log 

would give the following type of output:

HAS ENERGY VALUE   -108.3566620 
HAS ENERGY VALUE   -108.6053845 
HAS ENERGY VALUE   -108.7675654 
HAS ENERGY VALUE   -108.8668845 
HAS ENERGY VALUE   -108.9206732 
HAS ENERGY VALUE   -108.9417215 
HAS ENERGY VALUE   -108.9395089 
HAS ENERGY VALUE   -108.9210615 

You then have the R and V(R) values printed in a more convenient way, and you can then chop the rest off easily in VIM, or if you're unfamiliar with VIM you can import this data into Excel or Google Sheets or MATLAB and have the software automatically separate the columns, from which point you can copy and paste the results into a table.

The grep command can be improved, and combined with clever use of awk to make the table directly from the command line, but it would be more complicated and my above answer is probably the simplest and most immediate solution without figuring out how to get GAMESS to make the table the way MOLPRO does (if GAMESS can even do it).

If anyone wishes to test the above grep commands or to test their own improvements of them, the above grep results were done on this file which is now in the "Modeling Matters" GitHub repository.

Answer 2

It's not the simplest alternative, but sometimes to learn a bit of regular expressions and one scripting language can help a lot in tasks like this. It's a really worth time investment, if you have a little time to spare.

For example, I ran your input in my machine, saving the results to PES_scan_oxygen.log, and in Python 3 wrote a script to read this log and extract the coordinates and energies to a CSV file, that can be opened with a spreadsheet. The lines preceded by # are comments, to help explain what is done at each step of the script:

# Imports regular expression (regex) library
import re

def pes2csv(gamess_output):
    '''
    Function takes the PES scan output as its input, and returns .CSV file
    with matching values of radius and energy.
    '''
    # Creates regex pattern for line with energy value
    energy_pattern = re.compile(r"HAS ENERGY VALUE\s+(\S+)")
    # Creates regex pattern for line with radius value
    coordinate_pattern = re.compile(r"COORD 1= (\S+) ")

    # Opens file with the raw data
    with open(gamess_output, "r") as f:
        # Splits the file in lines, makes a list with them
        line_list = f.readlines()

    # Initiate list where we'll save our R and Energy values with headers
    pair_list = ["R;Energy\n"]

    # For every line i in the file
    for i in range(len(line_list)):
        # If line i is the start of surface mapping geometry section
        if line_list[i] == " ---- SURFACE MAPPING GEOMETRY ----\n":
            # Searches the pattern for radius in the next line
            R_match = coordinate_pattern.search(line_list[i+1])
            # Searches the pattern for energy in the second next line
            Energy_match = energy_pattern.search(line_list[i+2])
            # Assembles one line for our final CSV using ";" as delimiter, puts it on pair_list
            pair_list.append(R_match.group(1) + ";" + Energy_match.group(1) + "\n")

    # Creates a CSV file named PES.csv, in the same folder the script is saved
    with open("PES.csv", "w") as g:
        # Save the contents of pair_list to the CSV file we just created
        g.writelines(pair_list)
    print("Data saved to CSV file. Look your folder")
    
# Now we can call the function on our gamess output file, to extract the data
pes2csv("PES_scan_oxygen.log")

If you have Python 3 installed in your machine, save this script to a .py file in the same folder where the gamess output file is located, and runs it, now there's a PES.csv file in the folder, with just the data you want, as seen in the spreadsheet:

The decimal separator looks wrong, but it's just because I have the language of Libreoffice set to Portuguese, where dots are used as the thousands delimiter and comma to separate the decimal part of a number (the reverse of what is done in English).

For this particular example, to write down the script was probably overkill, as the final table has only 11 lines, and Nike's suggestion of using grep + manual import would be faster. But the problem with approaches involving manual steps is that they don't scale well. In large data sets it becomes tedious and time-consuming. In large spreadsheets you can spend a lot of time just scrolling up and down. So the time investment learning to automate it pays when you have large files or lots of files, what eventually happens with everybody.