Psithon: Structuring an Input File¶

To allow arbitrarily complex computations to be performed, PSI4 was built upon the Python interpreter. However, to make the input syntax simpler, some pre-processing of the input file is performed before it is interpreted, resulting in Python syntax that is customized for PSI, termed Psithon. In this section we will describe the essential features of the Psithon language. PSI4 is distributed with an extensive test suite, described in section Test Suite and Sample Inputs; the input files for these test cases can be found in the samples subdirectory of the top-level PSI4 source directory, and should serve as useful examples.

Physical Constants¶

For convenience, the Python interpreter will execute the contents of the ~/.psi4rc file in the current user’s home area (if present) before performing any tasks in the input file. This allows frequently used python variables to be automatically defined in all input files. For example, if we repeatedly make use of the universal gravitational constant, the following line could be placed in the ~/.psi4rc file

UGC = 6.67384E-11 # m^3 / kg^-1 s^-2

which would make the variable UGC available in all PSI4 input files. For convenience, the physical constants used within the PSI4 code (which are obtained from the 3rd edition of the IUPAC Green book [Cohen:GreenBook:2008]) are also automatically loaded as Psithon variables (before ~/.psi4rc is loaded, so that ~/.psi4rc values can be overridden by the user).

The physical constants used within PSI4, which are automatically made available within all PSI4 input files.

psi_h                        = 6.62606896E-34       # The Planck constant (Js)               
psi_c                        = 2.99792458E8         # Speed of light (ms$^{-1}$)             
psi_kb                       = 1.3806504E-23        # The Boltzmann constant (JK$^{-1}$)     
psi_R                        = 8.314472             # Universal gas constant (JK$^{-1}$mol$^{-1}$) 
psi_bohr2angstroms           = 0.52917720859        # Bohr to Angstroms conversion factor    
psi_bohr2m                   = 0.52917720859E-10    # Bohr to meters conversion factor       
psi_bohr2cm                  = 0.52917720859E-8     # Bohr to centimeters conversion factor  
psi_amu2g                    = 1.660538782E-24      # Atomic mass units to grams conversion factor 
psi_amu2kg                   = 1.660538782E-27      # Atomic mass units to kg conversion factor 
psi_au2amu                   = 5.485799097E-4       # Atomic units (m$@@e$) to atomic mass units conversion factor 
psi_hartree2J                = 4.359744E-18         # Hartree to joule conversion factor     
psi_hartree2aJ               = 4.359744             # Hartree to attojoule (10$^{-18}$J) conversion factor 
psi_cal2J                    = 4.184                # Calorie to joule conversion factor     
psi_dipmom_au2si             = 8.47835281E-30       # Atomic units to SI units (Cm) conversion factor for dipoles 
psi_dipmom_au2debye          = 2.54174623           # Atomic units to Debye conversion factor for dipoles 
psi_dipmom_debye2si          = 3.335640952E-30      # Debye to SI units (Cm) conversion factor for dipoles 
psi_c_au                     = 137.035999679        # Speed of light in atomic units         
psi_hartree2ev               = 27.21138             # Hartree to eV conversion factor        
psi_hartree2wavenumbers      = 219474.6             # Hartree to cm$^{-1}$ conversion factor 
psi_hartree2kcalmol          = 627.5095             # Hartree to kcal mol$^{-1}$ conversion factor 
psi_hartree2MHz              = 6.579684E9           # Hartree to MHz conversion factor       
psi_kcalmol2wavenumbers      = 349.7551             # kcal mol$^{-1}$ to cm$^{-1}$ conversion factor 
psi_e0                       = 8.854187817E-12      # Vacuum permittivity (Fm$^{-1}$)        
psi_na                       = 6.02214179E23        # Avagadro's number                      
psi_me                       = 9.10938215E-31       # Electron rest mass (in kg)             

The psi_ prefix is to prevent clashes with user-defined variables in PSI4 input files.

Molecule and Geometry Specification¶

PSI4 has a very flexible input parser that allows the user to provide geometries as Cartesian coordinates, Z-matrix variables, or a combination of both. The use of fixed values and variables are supported for both. For example, the geometry for H₂ can be specified a number of ways, using the molecule keyword:

molecule{
  H
  H 1 0.9
}

 or

molecule{
  H
  H 1 r
  r = 0.9
}

 or

molecule{
  H1
  H2 H1 0.9
}

 or

molecule{
  H 0.0 0.0 0.0
  H 0.0 0.0 0.9
}

 or

molecule{
  H 0.0 0.0 0.0
  H 0.0 0.0 r
  r = 0.9
}

 or

molecule{
  H 0.0 0.0 -r
  H 0.0 0.0 r
  r = 0.45
}

Blank lines are ignored and, unlike regular Python syntax, indentation within the molecule block does not matter, although the molecule keyword itself must be aligned within the input according to standard Python syntax. For more examples of geometry specification, see the mints1 input file in the samples folder. It is also possible to mix Cartesian and Z-matrix geometry specifications, as demonstrated in the mints4 and mints6 sample input files. For example, consider the following geometry specification, taken from the mints6 input:

molecule alanine {
    N           -1.527107413251     0.745960643462     0.766603000356
    C           -0.075844098953     0.811790225041     0.711418672248
    C            0.503195220163    -0.247849447550    -0.215671574613
    O           -0.351261319421    -0.748978309671    -1.089590304723
    O            1.639498336738    -0.571249748886    -0.174705953194
    H           -1.207655674855    -0.365913941094    -0.918035522052
    # First, remove the H from the alpha carbon.  This line could be deleted
    # and is only included for completeness
    #H            0.429560656538     0.717651915252     1.673774709694
    # Now patch it, using a Z Matrix specification.  This patch can be applied
    # anywhere in the coord specification, as long as it appears lower than
    # the atoms referenced, as is usual for Z-Matrices
    C  2  rCC   3  aCCC   1  dCCCN
    H  7  rCH1  2  aHCC1  3  dHCCC1
    H  7  rCH2  2  aHCC2  3  dHCCC2
    H  7  rCH3  2  aHCC3  3  dHCCC3
    H            0.221781602033     1.772570540211     0.286988509018
    H           -1.833601608592     0.108401996052     1.481873213172
    H           -1.925572581453     1.640882152784     0.986471814808

    aCCC = 108.0
    rCC = 1.4
    dCCCN = 120
    rCH1 = 1.08
    rCH2 = 1.08
    rCH3 = 1.08
    aHCC1 = 109.0
    aHCC2 = 109.0
    aHCC3 = 109.0
    dHCCC1 = 0.0
    dHCCC2 = 120.0
    dHCCC3 = 240.0
}

Here, we remove the hydrogen from the alpha carbon of glycine and replace it with a methyl group. Applying this patch using Cartesian coordinates is difficult, because it depends on the orientation of the existing glycine unit. In this example, we use Z-Matrix coordinates to define the methyl group, and define the orientation in terms of the existing glycine Cartesian coordinates which is much easier to visualize than the corresponding Cartesian-only approach.

Multiple Molecules¶

To facilitate more elaborate computations, it is possible to provide a name for each molecule, and tell PSI4 which one should be used in a given calculation. For example, consider the following input file:

molecule h2{
  H
  H 1 0.9
}


set basis cc-pvdz
set reference rhf
energy('scf')

molecule h{
  H
}

set basis cc-pvdz
set reference uhf
energy('scf')

Here, two separate jobs are performed on two different molecules; the first is performed on H₂, while the second is for H atom. The last molecule to be specified is the “active” molecule by default. To explicitly activate a named molecule, the activate keyword is provided. Using this keyword, the above input file can be equivalently written as follows:

molecule h2{
  H
  H 1 0.9
}

molecule h{
  H
}

activate(h2)
set basis cc-pvdz
set reference rhf
energy('scf')

activate(h)
set basis cc-pvdz
set reference uhf
energy('scf')

Note that whenever the molecule is changed, the basis set must be specified again. The following section provides more details about the job control keywords used in the above examples.

Molecule Keywords¶

In addition to specifying the geometry, additional information can be provided in the molecule optional_molecule_name {...} block. If two integers charge multiplicity are encountered on any line of the molecule block, they are interpreted as the molecular charge and multiplicity ( $2 \times M_s + 1$ ), respectively. The symmetry can be specified by a line reading symmetry symbol, where symbol is the Schönflies symbol of the (Abelian) point group to use for the computation. This need not be specified, as the molecular symmetry is automatically detected by PSI4. Certain computations require that the molecule is not reoriented; this can be achieved by adding either no_reorient or noreorient. By default, Ångström units are used; this is changed by adding a line that reads units spec, where spec is one of ang, angstrom, a.u., au, or bohr.

Ghost Atoms¶

While many common computations, such as SAPT and counterpoise corrections, can be greatly simplified using the notation described in Non-Covalently Bonded Molecule Fragments, manual specification of ghost atoms is sometimes required. Either

molecule he2 {
    He
    Gh(He) 1 2.0
}

molecule he2 {
    He
    @He 1 2.0
}

will generate a helium dimer, with the second atom ghosted, i.e., possessing basis functions but no electrons or nuclear charge. See dfmp2_1 and ghosts for a demonstration of both mechanisms for specifying ghost atoms.

Geometries from the PubChem Database¶

Obtaining rough starting guess geometries can be burdensome. The Z-matrix coordinate system was designed to provide chemists with an intuitive method for guessing structures in terms of bond lengths and angles. While Z-matrix input is intuitive for small molecules with few degrees of freedom, it quickly becomes laborious as the system size grows. To obtain a reasonable starting guess geometry, PSI4 can take a chemical name as input; this is then used to attempt to retrieve Cartesian coordinates from the [PubChem] database.

For example, to run a computation on benzene, we can use the following molecule specification:

molecule benzene {
    pubchem:benzene
}

If the computer is connected to the internet, the above code will instruct PSI4 to search PubChem for a starting structure. The search is actually performed for compounds whose name contains “benzene”, so multiple entries will be returned. If the name provided (“benzene” in the above example) exactly matches one of the results, that entry will be used. If no exact match is found the results, along with a unique chemical identifier (CID), are printed to the output file, prompting the user to provide a more specific name. For example, if we know that we want to run a computation on a compound whose name(s) contain “benzene”, but we’re not sure of the exact IUPAC name, the following input can be used:

molecule benzene {
    pubchem:benzene*
}

Appending the “*” prevents an exact match from being found and, at the time of writing, the following results are displayed in the output file:

Chemical ID     IUPAC Name
         241   benzene
        7371   benzenesulfonic acid
       91526   benzenesulfonate
         244   phenylmethanol
         727   1,2,3,4,5,6-hexachlorocyclohexane
         240   benzaldehyde
       65723   benzenesulfonohydrazide
       74296   N-phenylbenzenesulfonamide
         289   benzene-1,2-diol
         243   benzoic acid
        7370   benzenesulfonamide
      636822   1,2,4-trimethoxy-5-[(E)-prop-1-enyl]benzene
        7369   benzenesulfonyl chloride
       12932   N-[2-di(propan-2-yloxy)phosphinothioylsulfanylethyl]benzenesulfonamide
        7505   benzonitrile
       78438   N-[anilino(phenyl)phosphoryl]aniline
       12581   3-phenylpropanenitrile
      517327   sodium benzenesulfonate
      637563   1-methoxy-4-[(E)-prop-1-enyl]benzene
      252325   [(E)-prop-1-enyl]benzene

Note that some of these results do not contain the string “benzene”; these compounds have synonyms containing that text. We can now replace the “benzene*” in the input file with one of the above compounds using either the IUPAC name or the CID provided in the list, viz:

molecule benzene {
    pubchem:637563
}

 or

molecule benzene {
    pubchem:1-methoxy-4-[(E)-prop-1-enyl]benzene
}

Some of the structures in the database are quite loosely optimized and do not have the correct symmetry. Before starting the computation, PSI4 will check to see if the molecule is close to having each of the possible symmetries, and will adjust the structure accordingly so that the maximum symmetry is utilized.

The standard keywords, described in Sec. Molecule Keywords, can be used in conjuction to specify charge, multiplicity, symmetry to use, etc. .

Symmetry¶

For efficiency, PSI4 can utilize the largest Abelian subgroup of the full point group of the molecule. Concomitantly a number of quantities, such as SOCC and DOCC, are arrays whose entries pertain to irreducible representations (irreps) of the molecular point group. Ordering of irreps follows the convention used in Cotton’s Chemical Applications of Group Theory, as detailed in Table Irreps. We refer to this convention as “Cotton Ordering” hereafter.