MPQC(1)

πŸ“¦
version 2.3.1 Β· 2 versions
, in debian bookworm, ubuntu jammy, ubuntu noble

πŸ’– MPQC β€” the Massively Parallel Quantum Chemistry program πŸ’–

okay bestie, sit DOWN, because we're about to talk about mpqc, aka the Massively Parallel Quantum Chemistry program, and she is that girl. this is the General Commands Manual entry, section 1, and yes we're doing the whole thing. grab a boba. πŸ§‹


✨ NAME ✨

mpqc β€” The Massively Parallel Quantum Chemistry program (MPQC) πŸ¦‹


✨ SYNOPSIS ✨

mpqc [options] [filename]

✨ DESCRIPTION ✨

so like, MPQC computes the properties of molecules, ab initio, on a wide variety of computer architectures. like, she gets around. no cap. πŸ‘‘

she can compute closed shell and general restricted openshell HartreeFock energies and gradients, second order openshell perturbation theory (OPT2[2]) and Zaveraged perturbation theory (ZAPT2) energies, and second order closed shell MoellerPlesset perturbation theory energies and gradients. and she also includes methods for optimizing molecules in either Cartesian or internal coordinates. literally so versatile.

MPQC is designed using objectoriented programming techniques and implemented in the C++ programming language. she's giving very much structured, very much class hierarchy, very much slay. πŸ’…


✨ OPTIONS ✨

MPQC can be given options followed by an optional input file name. if you don't give her an input file name, it just defaults to 'mpqc.in' like "fine, i'll pick for us." πŸ™ˆ

these are the command line options she recognizes:

  • -o β€” gives the name of the output file. the default is the console. (yes she will just scream into the terminal, and honestly? iconic.)
  • -i β€” converts a simple input file to an object oriented input file and writes the result to the output. no calculations are done β€” it's a glow-up, not a science experiment. ✨
  • -messagegrp β€” a ParsedKeyVal specification of a MessageGrp object. the default depends on how MPQC was compiled. (iykyk.)
  • -memorygrp β€” a ParsedKeyVal specification of a MemoryGrp object. the default depends on how MPQC was compiled.
  • -threadgrp β€” a ParsedKeyVal specification of a ThreadGrp object. the default depends on how MPQC was compiled.
  • -l β€” sets a limit on the number of basis functions. the default is zero, which means an unlimited number of basis functions. go off, honestly.
  • -W β€” sets the working directory. the default is the current directory.
  • -c β€” check the input and exit. just a little vibe check. πŸ”
  • -v β€” print the version number. πŸ“Œ
  • -w β€” print the warranty information. and bestie… there is no warranty. 😭
  • -d β€” if a debugger object was given in the input, start the debugger running as soon as MPQC is started. straight into it. no hesitation. πŸ’€
  • -h β€” print a list of options. your little menu, your little helper. 🫢
  • -f β€” the name of an object-oriented input file. the default is mpqc.in. this cannot be used if another input file is specified.
⚠️ heads up: this option is deprecated, because both input file formats can be read by giving the input file name on the command line without any option flags. so like… just do that. periodt.

also, real talk: some MPI environments do not pass the command line to slave programs, but supply it when MPI_Init is called. to make MPQC call MPI_Init on start-up, instead of when an MPIMessageGrp is created, name the executable mpqc-mpi. cute little trick, very main character. 🌈


✨ ENVIRONMENTAL VARIABLES ✨

MPQC looks at four environmental variables to set up communication and find library files. machine specific libraries and utilities to run programs in parallel might look at other environment variables as well. but the four that apply on all platforms are:

  • SCLIBDIR β€” the name of the library directory.
  • MESSAGEGRP β€” a ParsedKeyVal specification of a MessageGrp object. the default depends on how MPQC was compiled. see the MessageGrp class documentation for more information.
  • MEMORYGRP β€” a ParsedKeyVal specification of a MemoryGrp object. the default depends on how MPQC was compiled and the MessageGrp in use.
  • THREADGRP β€” a ParsedKeyVal specification of a ThreadGrp object. the default depends on how MPQC was compiled.

by default, MPQC tries to find library files first in the lib sub-directory of the installation directory and then the source code directory. if the library files cannot be found, MPQC must be notified of the new location with the environmental variable SCLIBDIR. πŸ—‚οΈ

the other three keywords specify objects. this is done by giving a mini ParsedKeyVal input in a string. the object is anonymous, that is, no keyword is associated with it. here's an example, bestie:

setenv MESSAGEGRP '<ShmMessageGrp>:(n = 4)'

✨ SHARED MEMORY MULTIPROCESSOR WITH SYSV IPC ✨

by default, MPQC will run on only one CPU. to specify more, you can give a ShmMessageGrp object on the command line. the following would run mpqc in four processes:

mpqc -messagegrp '<ShmMessageGrp>:(n = 4)' input_file

alternately, the ShmMessageGrp object can be given as an environmental variable:

setenv MESSAGEGRP '<ShmMessageGrp>:(n = 4)'
mpqc input_file
😭 DANGER ZONE: if MPQC should unexpectedly die, shared memory segments and semaphores will be left on the machine. these should be promptly cleaned up or other jobs may be prevented from running successfully. like, she leaves her stuff everywhere and the next girl can't get ready. rude.

to see if you have any of these resources allocated, use the ipcs command. the output will look something like:

IPC status from /dev/kmem as of Wed Mar 13 14:42:18 1996
T     ID     KEY        MODE       OWNER    GROUP
Message Queues:
Shared Memory:
m 288800 0x00000000 --rw-------  cljanss     user
Semaphores:
s    390 0x00000000 --ra-------  cljanss     user
s    391 0x00000000 --ra-------  cljanss     user

to remove the IPC resources used by cljanss in the above example on IRIX, type:

ipcrm -m 288800
ipcrm -s 390
ipcrm -s 391

and on Linux, type:

ipcrm shm 288800
ipcrm sem 390
ipcrm sem 391

(notice the syntax is different on the two, bestie β€” that's not a typo, that's just how it is. 🫠)


✨ SHARED MEMORY MULTIPROCESSOR WITH POSIX THREADS ✨

by default, MPQC will run with only one thread. to specify more, you can give a PthreadThreadGrp object on the command line.

real talk: MPQC is not parallelized to as large an extent with threads as it is with the more conventional distributed memory model, so you might not get the best performance using this technique. on the other hand, the memory overhead is lower and no interprocess communication is needed. trade-offs, bestie, trade-offs. βš–οΈ

the following would run MPQC in four threads:

mpqc -threadgrp '<PthreadThreadGrp>:(num_threads = 4)' input_file

alternately, the PthreadThreadGrp object can be given as an environmental variable:

setenv THREADGRP '<PthreadThreadGrp>:(n = 4)'
mpqc input_file

✨ SHARED OR DISTRIBUTED MEMORY MULTIPROCESSOR WITH MPI ✨

a MPIMessageGrp object is used to run using MPI. the number of nodes used is determined by the MPI run-time and is not specified as input data to MPIMessageGrp. so don't even try to hand her that number, she doesn't want it. πŸ™…β€β™€οΈ

mpqc -messagegrp '<MPIMessageGrp>:()' input_file

alternately, the MPIMessageGrp object can be given as an environmental variable:

setenv MESSAGEGRP '<MPIMessageGrp>:()'
mpqc input_file

usually, a special command is needed to start MPI jobs; typically it is named mpirun. πŸš€


✨ INPUT ✨

MPQC supports two input formats.

  1. the primary input is an object oriented format which gives users access to all of MPQCs options. the full menu. the deluxe. πŸ‘‘
  2. the second format allows access to a subset of MPQCs capabilities, but is more intuitive and easier to learn. honestly? a great starter pack.

new users are advised to start with the simplified format. πŸ’– and here's the slay part β€” MPQC can be used to convert the simplified format to the full object-oriented format with the -i option. glow-up unlocked. ✨

✨ Simple Input ✨

Okay bestie, so the simple input format is literally just keywords followed by a : and then a value. The keywords are case sensitive β€” no cap. And the values might get modified by options that live in parentheses. For example, this iconic input does an optimization of water using density functional theory with the B3LYP exchange-correlation functional:

% B3LYP optimization of water
optimize: yes
method: KS (xc = B3LYP)
basis: 3-21G*
molecule:
    O    0.172   0.000   0.000
    H    0.745   0.000   0.754
    H    0.745   0.000  -0.754

Comments start with a % and go all the way to the end of the line. Basis set names containing special characters, such as a space or parentheses, must be quoted inside a pair of double quotes. The accepted keywords are:

  • molecule β€” gives the atom types and coordinates. The following options can be used:
    • bohr β€” the coordinates are given in Bohr.
    • angstrom β€” the coordinates are given in Angstroms.
    • charge β€” this option can be given after an element x y z quadruple. This will override the charge on the atom. For example, (charge = 0) can be given for the ghost atoms in a counterpoise correction calculation.
    • multiplicity β€” gives the multiplicity of the molecule. The default is 1.
  • optimize β€” if yes, then an optimization will be performed. The default is no. The following options can be given:
    • cartesian β€” use Cartesian coordinates.
    • internal β€” use internal coordinates.
    • redundant β€” use redundant internal coordinates.
  • gradient β€” if yes, then a gradient calculation will be performed. The default is no.
  • frequencies β€” if yes, then the frequencies will be obtained. The default is no.
  • charge β€” specifies the charge on the molecule. The default is 0.
  • method β€” specifies the method. There is no default and the possible values are:
    • HF β€” Hartree-Fock. Unrestricted HF is used if multiplicity > 1.
    • RHF β€” Restricted Hartree-Fock.
    • UHF β€” Unrestricted Hartree-Fock.
    • KS β€” Kohn-Sham. Unrestricted KS is used if multiplicity > 1.
    • RKS β€” Restricted Kohn-Sham.
    • UKS β€” Unrestricted Kohn-Sham.
    • MP2 β€” Second order Moeller-Plesset perturbation theory. Only available for multiplicity = 1.
    • ZAPT2 β€” Z-averaged perturbation theory. Only available for multiplicity > 1. No gradient, optimization, or frequencies are possible.
  • The following options are valid with the KS, RKS, and UKS methods:
    • grid β€” specifies the grid to be used for numerical integrations. The following values can be given:
      • xcoarse
      • coarse
      • medium
      • fine
      • xfine
      • ultrafine
    • xc β€” specifies the exchange-correlation functional. There is no default. See the table in the StdDenFunctional class documentation for the possible values.
    • basis β€” specifies the basis set. There is no default. See the table in the GaussianBasisSet class documentation for the available basis sets.
  • restart β€” set to yes to restart an optimization. The default is no.
  • checkpoint β€” set to no to not save checkpoint files during an optimization. The default is yes.
  • symmetry β€” specifies the Schoenflies symbol of the point group of the molecule. The default is auto, which will cause the program to find the highest order Abelian subgroup of the molecule.
  • docc β€” gives the number of doubly occupied orbitals in each irreducible representation in a parenthesized list. The symmetry must be specified and not be auto. The method must be restricted.
  • socc β€” gives the number of single occupied orbitals in each irreducible representation in a parenthesized list. The symmetry must be specified and not be auto. The method must be restricted.
  • alpha β€” gives the number of alpha occupied orbitals in each irreducible representation in a parenthesized list. The symmetry must be specified and not be auto. The method must be unrestricted.
  • beta β€” gives the number of beta occupied orbitals in each irreducible representation in a parenthesized list. The symmetry must be specified and not be auto. The method must be unrestricted.
  • frozen_docc β€” gives the number of frozen core orbitals. Can be either a single integer or a parenthesized list giving the frozen core orbitals in each irreducible representation. In the latter case the symmetry must be given and not be auto.
  • frozen_uocc β€” gives the number of frozen virtual orbitals. Can be either a single integer or a parenthesized list giving the frozen virtual orbitals in each irreducible representation. In the latter case the symmetry must be given and not be auto.

πŸ¦‹ Object-Oriented Input πŸ¦‹

MPQC is an object-oriented program that directly allows the user to specify objects that MPQC then manipulates to obtain energies, properties, etc. This makes the input very flexible, but very complex. However, most calculations should be quite similar to the one of the examples given later in this chapter. The best way to get started is to use one of the example input files and modify it to meet your needs.

MPQC starts off by creating a ParsedKeyVal object that parses the input file specified on the command line. The format of the input file is documented in . It is basically a free format input that associates keywords and logical groupings of keywords with values. The values can be scalars, arrays, or objects.

The keywords recognized by MPQC begin with the mpqc prefix. That is, they must be nested between an mpqc:( and a ). Alternately, each keyword can be individually prefixed by mpqc:. The primary keywords are given below. Some of the keywords specify objects, in which case the object will require more ParsedKeyVal input. These objects are created from the input by using their ParsedKeyVal constructors. These constructors are documented with the source code documentation for the class.

  • mole β€” This is the most important keyword for MPQC. It specifies the MolecularEnergy object. This is an object that knows how to compute the energy of a molecule. The specializations of MolecularEnergy that are most commonly used are CLKS, HSOSKS, UKS, CLHF, HSOSHF, UHF, and MBPT2.
  • opt β€” This keyword must be specified for optimizations. It specifies an Optimize object. Usually, QNewtonOpt is best for finding minima and EFCOpt is best for transition states.
  • freq β€” This keyword must be specified to compute frequencies. It specifies a MolecularFrequencies object.
  • thread β€” This specifies an object of type ThreadGrp that can be used to advantage on shared-memory multiprocessor machines for certain types of calculations. This keyword can be overridden by giving the ThreadGrp in the environment or command line. See the section on running MPQC for more information.
  • checkpoint β€” The value of this keyword is Boolean. If true, then optimizations will be checkpointed after each iteration. The checkpoint file suffix is .ckpt. The default is to checkpoint.
  • savestate β€” The value of this keyword is Boolean. If true, then the states of the optimizer and wavefunction objects will be saved after the calculation completes. The output file suffix is .wfn. The default is to save state.
  • restart β€” The value of this keyword is Boolean. If true, mpqc will attempt to restart the calculation. If the checkpoint file is not found, the calculation will continue as if the value were false. The default is true.
  • restart_file β€” This gives the name of a file from which restart information is read. If the file name ends in .wfn the MolecularEnergy object will be restored. Otherwise, the Optimize object will be restored. The default file name is formed by appending .ckpt to the input file name with the extension removed.
  • do_energy β€” The value of this keyword is Boolean. If true a single point energy calculation will be done for the MolecularEnergy object given with the mole keyword. The default is true.
  • do_gradient β€” The value of this keyword is Boolean. If true a single point gradient calculation will be done for the MolecularEnergy object given with the mole keyword. The default is false.
  • optimize β€” The value of this keyword is Boolean. If true and the opt keyword was set to a valid value, then an optimization will be performed. The default is true.
  • write_pdb β€” The value of this keyword is Boolean. If true a PDB file with the molecular coordinates will be written.
  • filename β€” The value of this keyword is a string that gives a name from which checkpoint and other filenames are constructed. The default is the basename of the input file.
  • print_timings β€” If this is true, timing information is printed at the end of the run. The default is true.

There are also some utility keywords that tell mpqc some technical details about how to do the calculation:

  • debug β€” This optional keyword gives a Debugger object which can be used to help find the problem if MPQC encounters a catastrophic error.
  • matrixkit β€” This optional keyword gives a SCMatrixKit specialization which is used to produce matrices of the desired type. The default is a ReplSCMatrixKit which replicates matrices on all of the nodes. Other choices are not thoroughly tested.

✨ EXAMPLES ✨

Okay bestie, real talk β€” this is the part where it all comes together 🫢. Let's do a Hartree-Fock calculation on water, and I'm gonna walk you through every single line because you deserve to actually get it, not just copy-paste and pray πŸ™. First, the whole input, then the breakdown with descriptions, mwah πŸ’‹.

% This input does a Hartree-Fock calculation on water.
molecule<Molecule>: (
  symmetry = C2V
  unit = angstrom
  { atoms geometry } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
  )
)

We start with a descriptive comment. Comments begin with a %. Everything from the % to the end of the line is ignored. Literally, the computer does not care. It's just for us πŸ’….

% This input does a Hartree-Fock calculation on water.

Now lets set up a Molecule object. The name of the object comes first β€” it is molecule. Then, in angle brackets, comes the type of the molecule, which is the class Molecule. The keyword and class name are followed by a : and then several pieces of input grouped between a pair of matching parentheses. These parentheses contain the information that will be given to Molecule KeyVal constructor. So it's like, "hey bestie, here's everything you need to know about me" πŸ’–.

molecule<Molecule>: (

The point group of the molecule is needed. This is done by assigning symmetry to a case insensitive Schoenflies symbol that is used to initialize a PointGroup object. An Abelian point group should be used. Noted πŸŽ€.

  symmetry = C2V

The default unit for the Cartesian coordinates is Bohr. You can specify other units by assigned unit to a string that will be used to initialize a Units object. Iconic of us to give options ✨.

  unit = angstrom

Finally, the atoms and coordinates are given. This can be given in the shorthand table syntax shown below. The headings of the table are the keywords between the first pair of brackets. These are followed by an = and another pair of brackets that contain the data. The first datum is assigned to the first element of the array that corresponds to the first heading, atom. The second datum is assigned to the first element of the array associated with the second heading, geometry, and so on. Here the second datum is actually a vector: the x, y and z coordinates of the first atom. Okay this part is giving spreadsheet energy but make it quantum 🧚.

  { atoms                       geometry                   } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)

Next, a basis set object is given.

basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)

Now we will give the main body of input. All the subsequent keywords will be grouped in the mpqc section of the input (that is, each keyword will be prefixed with mpqc:). Lowkey the mpqc: block is the main character of this whole file πŸ‘‘.

mpqc: (

Next we give the mole keyword which provides a specialization of the MolecularEnergy class. In this case we will do a closed-shell Hartree-Fock calculation. That is done with an object of type CLHF. The keywords that CLHF accepts are given with the documentation for the CLHF class, usually in the description of the const RefKeyVal& constructor for the class. Also with the CLHF documentation is a list of parent classes. Each of the parent classes may also have input. This input is included with the rest of the input for the child class. It's giving inheritance, and honestly, obsessed πŸ’œ.

  mole<CLHF>: (

The next line specifies the molecule to be used. There are two things to note, first that this is actually a reference to complete molecule specification elsewhere in the input file. The $ indicates that this is a reference and the keyword following the $ is the actual location of the molecule. The : in front of the keyword means that the keyword is not relative to the current location in the input, but rather relative to the root of the tree of keywords. Thus, this line grabs the molecule that was specified above. The molecule object could have been placed here, but frequently it is necessary that several objects refer to the exact same object and this can only be done using references. It's like tagging someone in a group chat instead of retyping their whole bio 🫢.

The second point is that if you look at the documentation for CLHF, you will see that it doesn't read molecule keyword. However, if you follow its parent classes up to MolecularEnergy, you'll find that molecule is indeed read. No cap β€” always check the parents 😭.

    molecule = $:molecule

Just as we gave molecule, specify the basis set with the basis keyword as follows:

    basis = $:basis

Now we close off the parentheses we opened above and we are finished. Periodt πŸ’….

  )
)

🌸 Sample Object-Oriented Input Files 🌸

The easiest way to get started with mpqc is to start with one of sample inputs that most nearly matches your problem. All of the samples inputs shown here can be found in the directory src/bin/mpqc/samples. Go forth and steal shamelessly, bestie πŸ“.

πŸ’– Hartree-Fock Energy πŸ’–

The following input will compute the Hartree-Fock energy of water.

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C2V
  unit = angstrom
  { atoms geometry } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    memory = 16000000
  )
)

🌈 MP2 Energy 🌈

The following input will compute the MP2 energy of water.

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C2V
  unit = angstrom
  { atoms geometry } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % method for computing the molecule's energy
  mole<MBPT2>: (
    molecule = $:molecule
    basis = $:basis
    memory = 16000000
    % reference wavefunction
    reference<CLHF>: (
      molecule = $:molecule
      basis = $:basis
      memory = 16000000
    )
  )
)

πŸ¦‹ Hartree-Fock Optimization πŸ¦‹

The following input will optimize the geometry of water using the quasi-Newton method.

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C2V
  unit = angstrom
  { atoms geometry } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '6-31G*'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
  )
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
  )
  % optimizer object for the molecular geometry
  opt<QNewtonOpt>: (
    function = $..:mole
    update<BFGSUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

✨ Optimization with a Computed Guess Hessian ✨

The following input will optimize the geometry of water using the quasi-Newton method. The guess Hessian will be computed at a lower level of theory. It's giving "do the cheap version first, then glow up" πŸ’….

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C2V
  unit = angstrom
  { atoms geometry } = {
    O     [     0.00000000     0.00000000     0.37000000 ]
    H     [     0.78000000     0.00000000    -0.18000000 ]
    H     [    -0.78000000     0.00000000    -0.18000000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '6-31G*'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
  )
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
    guess_hessian<FinDispMolecularHessian>: (
      molecule = $:molecule
      only_totally_symmetric = yes
      eliminate_cubic_terms = no
      checkpoint = no
      energy<CLHF>: (
        molecule = $:molecule
        memory = 16000000
        basis<GaussianBasisSet>: (
          name = '3-21G'
          molecule = $:molecule
        )
      )
    )
  )
  % optimizer object for the molecular geometry
  opt<QNewtonOpt>: (
    function = $..:mole
    update<BFGSUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

πŸ”₯ Optimization Using Newton's Method πŸ”₯

The following input will optimize the geometry of water using the Newton's method. The Hessian will be computed at each step in the optimization. However, Hessian recomputation is usually not worth the cost; try using the computed Hessian as a guess Hessian for a quasi-Newton method before resorting to a Newton optimization.

Okay so this is the "expensive but thorough" option β€” the manpage itself is telling you it's usually not worth the cost, so listen to it before you tank your runtime bestie 😭.

okay so we're mid-flow through the example files, no big deal, just sliding right in like we never left πŸ’…

% Emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = c2v
  unit = angstrom
  { atoms geometry } = {
     O     [     0.00000000     0.00000000     0.36937294 ]
     H     [     0.78397590     0.00000000    -0.18468647 ]
     H     [    -0.78397590     0.00000000    -0.18468647 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '3-21G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  restart = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
  )
  do_energy = no
  do_gradient = no
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    memory = 16000000
    coor = $..:coor
    guess_wavefunction<CLHF>: (
      molecule = $:molecule
      total_charge = 0
      basis<GaussianBasisSet>: (
        molecule = $:molecule
        name = 'STO-3G'
      )
      memory = 16000000
    )
    hessian<FinDispMolecularHessian>: (
      only_totally_symmetric = yes
      eliminate_cubic_terms = no
      checkpoint = no
    )
  )
  optimize = yes
  % optimizer object for the molecular geometry
  opt<NewtonOpt>: (
    print_hessian = yes
    max_iterations = 20
    function = $..:mole
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

πŸ¦‹ Hartree-Fock Frequencies πŸ¦‹

okay so this next input file? it computes Hartree-Fock frequencies via finite displacements. bestie it also runs a thermodynamic analysis for you, because why not ✨ And here's the tea: if you also hand it optimization input, it'll run the optimization first, and then do the frequencies. Order matters, no cap.

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C1
  { atoms geometry } = {
    O     [  0.0000000000    0.0000000000    0.8072934188 ]
    H     [  1.4325589285    0.0000000000   -0.3941980761 ]
    H     [ -1.4325589285    0.0000000000   -0.3941980761 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    memory = 16000000
  )
% vibrational frequency input
  freq<MolecularFrequencies>: (
    molecule = $:molecule
  )
)

πŸŽ€ Giving Coordinates and a Guess Hessian πŸŽ€

so this example is literally serving several features that are really independent all at once. iconic. let's break it down:

  • the variable coordinates are explicitly given instead of being generated automatically for you. that's highkey useful when you're handing over a guess Hessian, which is exactly what's happening here πŸ‘€
  • that Hessian, as given by the user, is not complete β€” and QNewtonOpt will fill in the missing values using a guess the Hessian provided by the MolecularEnergy object. so don't stress about blanks, she's got it 🫢
  • also, fixed coordinates are given in this sample input. yes bestie, both at once. multitasking queen.
% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C1
  { atoms geometry } = {
      H    [ 0.088    2.006    1.438 ]
      O    [ 0.123    3.193    0.000 ]
      H    [ 0.088    2.006   -1.438 ]
      O    [ 4.502    5.955   -0.000 ]
      H    [ 2.917    4.963   -0.000 ]
      H    [ 3.812    7.691   -0.000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
  )
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
      extra_bonds = [ 2 5 ]
    )
    % use these instead of generated coordinates
    variable<SetIntCoor>: [
      <StreSimpleCo>:( atoms = [ 2 5 ] )
      <BendSimpleCo>:( atoms = [ 2 5 4 ] )
      <OutSimpleCo>: ( atoms = [ 5 2 1 3 ] )
      <SumIntCoor>: (
        coor: [
          <StreSimpleCo>:( atoms = [ 1 2 ] )
          <StreSimpleCo>:( atoms = [ 2 3 ] )
          ]
        coef = [ 1.0 1.0 ]
        )
      <SumIntCoor>: (
        coor: [
          <StreSimpleCo>:( atoms = [ 4 5 ] )
          <StreSimpleCo>:( atoms = [ 4 6 ] )
          ]
        coef = [ 1.0 1.0 ]
        )
      <BendSimpleCo>:( atoms = [ 1 2 3 ] )
      <BendSimpleCo>:( atoms = [ 5 4 6 ] )
    ]
    % these are fixed by symmetry anyway,
    fixed<SetIntCoor>: [
      <SumIntCoor>: (
        coor: [
          <StreSimpleCo>:( atoms = [ 1 2 ] )
          <StreSimpleCo>:( atoms = [ 2 3 ] )
          ]
        coef = [ 1.0 -1.0 ]
        )
      <SumIntCoor>: (
        coor: [
          <StreSimpleCo>:( atoms = [ 4 5 ] )
          <StreSimpleCo>:( atoms = [ 4 6 ] )
          ]
        coef = [ 1.0 -1.0 ]
        )
      <TorsSimpleCo>:( atoms = [ 2 5 4 6] )
      <OutSimpleCo>:( atoms = [ 3 2 6 4 ] )
      <OutSimpleCo>:( atoms = [ 1 2 6 4 ] )
    ]
  )
  % optimizer object for the molecular geometry
  opt<QNewtonOpt>: (
    function = $..:mole
    update<BFGSUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
    % give a partial guess hessian in internal coordinates
    % the missing elements will be filled in automatically
    hessian = [
        [  0.0109261670 ]
        [ -0.0004214845    0.0102746106  ]
        [ -0.0008600592    0.0030051330    0.0043149957 ]
        [  0.0             0.0             0.0          ]
        [  0.0             0.0             0.0          ]
        [  0.0             0.0             0.0          ]
        [  0.0             0.0             0.0          ]
     ]
  )
)

🌸 Optimization with a Hydrogen Bond 🌸

okay important heads-up bestie: the automatic internal coordinate generator will literally fail if it can't find enough redundant internal coordinates. it's giving picky. 😭

so what do you do? you explicitly create the internal coordinate generator in the input and hand it the extra connectivity info it's missing, like so πŸ‘‡

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = C1
  { atoms geometry } = {
      H    [ 0.088    2.006    1.438 ]
      O    [ 0.123    3.193    0.000 ]
      H    [ 0.088    2.006   -1.438 ]
      O    [ 4.502    5.955   -0.000 ]
      H    [ 2.917    4.963   -0.000 ]
      H    [ 3.812    7.691   -0.000 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = 'STO-3G'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
  )
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    % give an internal coordinate generator that knows about the
    % hydrogen bond between atoms 2 and 5
    generator<IntCoorGen>: (
      molecule = $:molecule
      extra_bonds = [ 2 5 ]
    )
  )
  % optimizer object for the molecular geometry
  opt<QNewtonOpt>: (
    function = $..:mole
    update<BFGSUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

πŸ’œ Fixed Coordinate Optimization πŸ’œ

this one shows you how to selectively fix internal coordinates in an optimization. you can give any number of linearly independent coordinates β€” but plot twist, they have to stay linearly independent the whole way through the optimization. and that condition might not hold, because the coordinates can be nonlinear. so like. be aware 🫠

by default, the initial fixed coordinates' values are taken from the cartesian geometry given by the Molecule object. however, the molecule will be displaced to the internal coordinate values given with the fixed internal coordinates if the have_fixed_values keyword is set to true β€” which is exactly what this example does ✨

in that case, the initial Cartesian geometry should be reasonably close to the desired initial geometry, and all of the variable coordinates will be frozen to their original values during the initial displacement. noted, bestie πŸ“

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = CS
  { atoms geometry } = {
    H [  3.04 -0.69 -1.59 ]
    H [  3.04 -0.69  1.59 ]
    N [  2.09 -0.48 -0.00 ]
    C [ -0.58 -0.15  0.00 ]
    H [ -1.17  1.82  0.00 ]
    H [ -1.41 -1.04 -1.64 ]
    H [ -1.41 -1.04  1.64 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '4-31G*'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
    have_fixed_values = yes
    fixed<SetIntCoor>: [
      <OutSimpleCo>: ( value = -0.1
                       label = 'N-inversion'
                       atoms = [4 3 2 1] )
      ]
  )
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
  )
  % optimizer object for the molecular geometry
  opt<QNewtonOpt>: (
    max_iterations = 20
    function = $..:mole
    update<BFGSUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

πŸ‘‘ Transition State Optimization πŸ‘‘

this example shows a transition state optimization of the N-inversion in using mode following. the initial geometry was obtained by doing a few fixed coordinate optimizations along the inversion coordinate.

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = CS
  { atoms geometry } = {
    H [  3.045436 -0.697438 -1.596748 ]
    H [  3.045436 -0.697438  1.596748 ]
    N [  2.098157 -0.482779 -0.000000 ]
    C [ -0.582616 -0.151798  0.000000 ]
    H [ -1.171620  1.822306  0.000000 ]
    H [ -1.417337 -1.042238 -1.647529 ]
    H [ -1.417337 -1.042238  1.647529 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '4-31G*'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
    followed<OutSimpleCo> = [ 'N-inversion' 4 3 2 1 ]
  )
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
  )
  % optimizer object for the molecular geometry
  opt<EFCOpt>: (
    transition_state = yes
    mode_following = yes
    max_iterations = 20
    function = $..:mole
    update<PowellUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

🌸 Transition State Optimization with a Computed Guess Hessian 🌸

okay so this next one is the glow-up version of that same vibe. πŸ’…

It's a transition state optimization of the N-inversion, and it's still doing the mode following thing. The starting geometry? Got it by running a few fixed coordinate optimizations along the inversion coordinate β€” like, they walked it in first, no cap. πŸšΆβ€β™€οΈβœ¨

Then the magic: an approximate guess Hessian gets computed, and that is what makes the optimization converge so much faster in this case. it's giving "do a little prep work, save a whole lot of time," fr fr. πŸ¦‹

% emacs should use -*- KeyVal -*- mode
% molecule specification
molecule<Molecule>: (
  symmetry = CS
  { atoms geometry } = {
    H [  3.045436 -0.697438 -1.596748 ]
    H [  3.045436 -0.697438  1.596748 ]
    N [  2.098157 -0.482779 -0.000000 ]
    C [ -0.582616 -0.151798  0.000000 ]
    H [ -1.171620  1.822306  0.000000 ]
    H [ -1.417337 -1.042238 -1.647529 ]
    H [ -1.417337 -1.042238  1.647529 ]
  }
)
% basis set specification
basis<GaussianBasisSet>: (
  name = '4-31G*'
  molecule = $:molecule
)
mpqc: (
  checkpoint = no
  savestate = no
  % molecular coordinates for optimization
  coor<SymmMolecularCoor>: (
    molecule = $:molecule
    generator<IntCoorGen>: (
      molecule = $:molecule
    )
    followed<OutSimpleCo> = [ 'N-inversion' 4 3 2 1 ]
  )
  % method for computing the molecule's energy
  mole<CLHF>: (
    molecule = $:molecule
    basis = $:basis
    coor = $..:coor
    memory = 16000000
    guess_hessian<FinDispMolecularHessian>: (
      molecule = $:molecule
      only_totally_symmetric = yes
      eliminate_cubic_terms = no
      checkpoint = no
      energy<CLHF>: (
        molecule = $:molecule
        memory = 16000000
        basis<GaussianBasisSet>: (
          name = '3-21G'
          molecule = $:molecule
        )
      )
    )
  )
  % optimizer object for the molecular geometry
  opt<EFCOpt>: (
    transition_state = yes
    mode_following = yes
    max_iterations = 20
    function = $..:mole
    update<PowellUpdate>: ()
    convergence<MolEnergyConvergence>: (
      cartesian = yes
      energy = $..:..:mole
    )
  )
)

πŸŽ€ note bestie: that inner energy<CLHF> block is running on a whole different basis ('3-21G', not the '4-31G*' from the outer one) β€” that's intentional, it's the cheap little scout run for the guess Hessian. copy it exactly, don't get cute and "unify" them. πŸ™ˆ

and that's the whole page, bestie!! πŸ’–βœ¨ you just read every last option, every environment variable, every shared-memory and MPI setup, every input flavor, and all the example files β€” like, actually all of them. you ate that. πŸ‘‘

This one is Version 2, dated 15 Feb 2008, and it lives under the name MPQC(1) β€” so if you're citing it or looking it up, that's your tag. 🫢 go forth and compute, you gorgeous genius. πŸŒΈπŸ¦‹πŸ’œ

πŸ“¦ where this page comes from