CASINO workflow

This workflow performs a quantum monte carlo (QMC) calculation using "trial" wave function (WFN) generated by QCHEM and ORCA.

Authors

• Vladimir Konjkov

Usage

Step 1: Install workflow and additional utilities

Download CASINO QMC program and install it.

Download CASINO workflow from repository

Download molden2qmc.py utility required to convert "trial" WFN from MOLDEN output format to CASINO input format.

Download multideterminant.py utility required to convert multi determinant information available in CASSCF and post-HF methods to CASINO format.

Choose program to generate "trial" WFN (QCHEM > 4.0 or ORCA > 4.0.1) and install it.

Make sure all these programs are accessible through the PATH variable.

Step 2: Configure workflow

Configure the workflow according to your needs via editing the file ORCA/Snakefile or QCHEM/Snakefile depending on what program you will use.

The basic information necessary for calculations is contained in config.yaml file:

• STD_ERR: desireable accuracy of DMC energy

• VMC_NCONFIG is a number of VMC config in VMC optimization step

• MOLECULES is a list of molecular geometry file names in xyz-format (without extension) located in the chem_database directory for which the calculation will be performed.

• METHODS is a list of quantum chemistry methods to obtain a "trial" WFN.

• BASES is a list of bases in which the "trial" function will be represented.

• JASTROWS is a filename of template (without extension) defining JASTROW factor

Depending on the program used to generate "trial" WFN, the following rules are available:

ORCA rules:

• generation ORCA output file.

rule ALL_ORCA:
    input: '{method}/{basis}/{molecule}/mol.out'

where:

• method - method available in ORCA to calculate "trial" WFN like HF, any DFT methods (i.e. B3LYP, CAM-BLYP, PBE0), OO-RI-MP2, OO-RI-SCS-MP2, CASSCF(N.M) for multideterminant extension.
• basis - any basis available in ORCA (i.e. cc-pVDZ, aug-cc-pVQZ, def2-SVP).
• molecule - molecular geometry file name in xyz-format (without extension) located in the chem_database directory.

QCHEM rules

• generation QCHEM output file.

rule ALL_QCHEM:
    input: '{method}/{basis}/{molecule}/mol.out'

where:

• method - method available in QCHEM to calculate "trial" WFN including HF, any DFT methods (i.e. B3LYP, CAM-BLYP, PBE0), any orbital optimized methods from the list (OD, OD(2), VOD, VOD(2), QCCD, QCCD(2), VQCCD). in case of OO-method T2-amplitudes subsequently used as determinant's weights but some type of active space truncation should be specified.
  It should be done with two ways:
  • OD all available T2-amplitudes were taken.
  • OD_10 first 10 active orbitals were taken.
  • OD_0.01 all orbitals with T2-amplitudes greater than 0.01 were taken.
• basis - any basis available in QCHEM (i.e. cc-pVDZ, aug-cc-pVQZ, pc-1).
• molecule - molecular geometry file names in xyz-format (without extension) located in the chem_database directory.

QMC rules

• pure VMC calculation (without JASTROW)

  rule ALL_VMC:
      input: '{method}/{basis}/{molecule}/VMC/{nstep}/out'
  

  where:

  • nstep - number of VMC statistic accumulation steps.

  This rule intended to check whether the conversion of the "trial" WFN to the CASINO format is correct and HF energy is equal to pure VMC one. It is not required to calculate the DMC energy.

• JASTROW coefficients optimization using some optimization plan

  rule ALL_VMC_OPT:
      input: '{method}/{basis}/{molecule}/VMC_OPT/{opt_plan}/{jastrow}/out'
  

  where:

  • opt_plan - name (without extension) of CASINO input file in opt_plan directory which define JASTROW optimization plan.
  • jastrow - name (without extension) of template for parameters.casl file in casl directory.

• VMC energy calculation with optimized JASTROW coefficients.

  rule ALL_VMC_OPT_ENERGY:
      input: '{method}/{basis}/{molecule}/VMC_OPT_ENERGY/{opt_plan}/{jastrow}/{nstep}/out'
  

  where:

  • nstep - number of VMC statistic accumulation steps.

  This rule is not required to calculate the DMC energy, but this may be necessary if you need the VMC energy with high accuracy.

• FN-DMC energy calculation to achieve desired accuracy (specified in the variable STD_ERR).

  rule ALL_VMC_DMC:
      input: '{method}/{basis}/{molecule}/VMC_DMC/{opt_plan}/{jastrow}/tmax_{dt}_{nconfig}_{i}/out'
  

  where:

  • dt - part of denominator (integer) to calculate the DMC time step using the formula 1.0/(max_Z² * 3.0 * dt).
  • nconfig - number of configuration in DMC calculation (1024 is recommended)
  • i - stage of DMC calculation (1 - only DMC equilibration and fixed step (50000) of DMC accumulation run, 2 - additional DMC accumulation run to achieve desired accuracy)

• JASTROW coefficients optimization using some optimization plan with BACKFLOW transformed WFN.

  rule ALL_VMC_OPT_BF:
      input: '{method}/{basis}/{molecule}/VMC_OPT_BF/{opt_plan}/{jastrow}__{backflow}/out'
  

  where:

  • backflow - combination of {ETA-term}_{MU-term}_{PHI-term-electron-nucleus}{PHI-term-electron-electron} expansion orders:
    • 3_3_00 - ETA-term expansion of order 3 and MU-term of order 3, no PHI-term
    • 9_9_33 - ETA-term expansion of order 9 and MU-term of order 9, PHI-term with electron-nucleus and electron-electron expansion of order 3 both.

  It should be noted that the BACKFLOW transformation in CASINO is implemented only up to f-orbitals.

• VMC energy calculation with optimized JASTROW coefficients with BACKFLOW transformation.

  rule ALL_VMC_OPT_ENERGY_BF:
      input: '{method}/{basis}/{molecule}/VMC_OPT_ENERGY_BF/{opt_plan}/{jastrow}__{backflow}/{nstep}/out'
  

  where:

  • nstep - number of VMC statistic accumulation steps.

  This rule is not required to calculate the DMC energy.

• FN-DMC energy calculation to achieve desired accuracy (specified in the variable STD_ERR) with BACKFLOW transformed WFN.

  rule ALL_VMC_DMC_BF:
      input: '{method}/{basis}/{molecule}/VMC_DMC_BF/{opt_plan}/{jastrow}__{backflow}/tmax_{dt}_{nconfig}_{i}/out'
  

Step 3: Execute workflow

All you need to execute this workflow is to install Snakemake via the Conda package manager. Software needed by this workflow is automatically deployed into isolated environments by Snakemake.

Test your configuration by performing a dry-run via

snakemake <rule> -n

If you took a files from previous calculations on another server or from this repository and want to continue the calculations,
you need to touch all files to prevent recalculations with the command

snakemake <rule> -t

Execute the workflow locally via

snakemake <rule>

Examples

To demonstrate the possibilities of this workflow, examples of calculations are given.

ORCA

(U)HF/def2-QZVP "trial" WFN for H-Ar atoms

Since the DMC calculation time for achieving a given accuracy increases ~ Z^5.5, I used several configuration files to perform calculations: config_00001.yaml for He-Li, config_0001.yaml for Be-N, config_003.yaml for O-Si, config_001.yaml for P-Ar.

The VMC optimization plan was used in the first step performs varmin optimization and then emin optimization for the remaining eight steps. This is a fairly effective optimization plan, as shown in the graphs, in all cases, it would be enough only 4-5 steps.

Other useful graphs summarizing the results of this calculations are given.

• H (²S_1/2)

  The Hartree–Fock (HF) "trial" WFN for the ground state of H-atom is nodeless, so DMC energy is exact, also as no electron correlations are present VMC energy is also exact.

  Ψ_HF(R) = ψ_1s(r₁)

• He (¹S₀)

  Spacial part of The Hartree–Fock (HF) "trial" WFN for the ground state of He-atom is symmetric so has no nodal surface, thus DMC energy is exact.

  Ψ_HF(R) = ψ_1s(r₁) • ψ_1s(r₂)

• Li (²S_1/2)

  For The Hartree–Fock (HF) "trial" WFN electrons nodal surface is determined by equation r₁ = r₂ when 1 and 2 label the electrons in the same spin channel [1].

  Ψ_HF(R) = det[ψ_1s(r₁), ψ_2s(r₂)] • ψ_1s(r₃)

  The electron 1 therefore “sees” the node as a sphere which passes through the position of electron 2 and is centered around the nucleus. The wave function will be equal to zero if electron 1 occupies any point on the spherical nodal surface. For the correlated electrons this is not strictly exact, as the correlation with the electron in the spin-down channel will cause deformations away from a perfect sphere. For example, the excitation 1s²2s¹2p² will have a contribution to the exact ground state and would in principle lead to a departure from the single particle node (i.e., the sphere will slightly deform to ellipsoid or perhaps a more complicated surface that would depend on the position of the minority spin electron). It is therefore quite remarkable that the HF nodal surface seems to be so accurate: the total energy with the HF nodes, is accurate to ~ 0.05 mHa.

• Be (¹S₀)

  The Hartree–Fock (HF) "trial" WFN for Be is given by a Slater determinant which is block diagonal in spin so it can be broken into a product of the spin channels [2].

  Ψ_HF(R) = det[ψ_1s(r₁), ψ_2s(r₂)] • det[ψ_1s(r₃), ψ_2s(r₄)]

  nodal surface is determined by equation (r₁ - r₂)(r₃ - r₄) = 0 which clearly shows that there are 2 * 2 = 4 nodal pockets. However, it has been found that for ground state the correct number of nodal domains is two. The accurate nodal surface for this system is actually remarkably well described by a two configuration wave function where HF WFN is augmented by adding 1s²2p² double excitation which corresponds to a near-degeneracy effect.

• B (²P_1/2)

  The Hartree–Fock (HF) "trial" WFN for B also should be improved by adding 1s²2p³ configuration, however, even after this improvement fixed node error is ~ 5 mHa and additional configurations is required to achieve chemical accuracy. The article [3] chapter 5. Quantum Monte Carlo Calculations with Multi-Reference Trial Wave Functions (H-J Flad et al.) points out that next important configuration is 1s²2s¹2p¹3d¹ which reduces the error to ~ 1.7 mHa.

• C (³P₀)

  The Hartree–Fock (HF) "trial" WFN for C also should be improved by adding 1s²2p⁴ configuration, but this again does not give us an exact nodal surface.

• N (⁴S_3/2)

  Only single determinant approximation was used.

• O (³P₂)

  Only single determinant approximation was used.

• F (²P_3/2)

  Only single determinant approximation was used.

• Ne (¹S₀)

  Only single determinant approximation was used.

• Na (¹S₀)

  Only single determinant approximation was used.

QCHEM