Skip to content
Yuan Yao edited this page Dec 25, 2021 · 9 revisions

Welcome to Arrow's wiki!

Installation and Compilation

Make sure you have installed MPI.

git clone --recursive https://github.com/QMC-Cornell/shci
cd shci
make -j

If using intel compiler, compile with make -f Makefile.intel -j.

Example Run

All you need as input are config.json and FCIDUMP files.

Example inputs and outputs for carbon and chromium atoms and nitrogen molecule are in the examples/ directory.

To run other systems, you will have to obtain an integrals file, FCIDUMP, and modify the values in config.json accordingly. Many software packages can generate FCIDUMP, such as PySCF and Molpro.

If running a system with > 128 orbitals, open shci/src/det/half_det.h and change the line #define N_CHUNKS 2 to a number so that N_CHUNKS * 64 >= n_orb before compiling.

Arrow uses hybrid MPI+OpenMP parallelization. Run 1 MPI process per node and number of OpenMP threads per node = # cores/ # nodes.

This command runs with one node:

mpirun -n 1 <path_to_executable>

Note: When using one node, make sure there is only one MPI process. When using multiple nodes, make sure all threads are running on each node.

config.json arguments

A minimal config.json file looks like this (this one is for the carbon atom provided in the examples/ directory):

{
  "system": "chem",
  "n_up": 4,
  "n_dn": 2,
  "eps_vars": [
    5e-5,
    2e-5,
    1e-5
  ],
  "eps_vars_schedule": [
    2e-3,
    1e-3,
    5e-4,
    2e-4,
    1e-4
  ],
  "chem": {
    "point_group": "d2h"
  }
}

Here, we have specified the number of electrons of both spins, the point group symmetry of the orbitals in FCIDUMP, and the eps_vars threshold values for the variational wave function. The eps_vars_schedule argument specified how the threshold should be gradually reduced.

The lists below contain some more input variables that can be added.

General

  • n_up, n_dn (required): number of up / down electrons.
  • system (required): only support chem and heg for now.
  • occs_up: occupied orbitals for the starting up determinant, default: lowest ones.
  • occs_dn: occupied orbitals for the starting dn determinant, default: lowest ones.
  • eps_vars (required): an array of variational epsilons from big to small.
  • eps_vars_schedule: an array of additional variational epsilons run beforehand for a better selection of determinants, default: an empty array.
  • eps_pt_dtm: 🌴 deterministic perturbation epsilon, default: max(eps_var*eps_pt_dtm_ratio,2e-6).
  • eps_pt_psto: 🌴 pseudo stochastic perturbation epsilon, default: max(eps_var*eps_pt_psto_ratio,1e-7)
  • eps_pt: 🌴 perturbation epsilon, default: eps_var*eps_pt_ratio.
  • eps_pt_dtm_ratio: 🌴 ratio eps_pt_dtm/eps_var, default: eps_var / 10.
  • eps_pt_psto_ratio: 🌴 ratio eps_pt_psto/eps_var, default: eps_var / 100.
  • eps_pt_ratio: 🌴 ratio eps_pt/eps_var, default: eps_var / 1000.
  • max_pt_iterations: 🌴 maximum stochastic perturbation iterations, default: 100.
  • n_batches_pt_sto: 🌴 number of batches for stochastic perturbation, default: 16.
  • n_samples_pt_sto: 🌴 number of samples for stochastic perturbation, default: choose based on available system memory.
  • random_seed: for stochastic perturbation, default: 347634253.
  • target_error: target error for stochastic perturbation, default: 1.0e-5.
  • var_only: run variation only, useful e.g. when optimizing orbs, default: false.
  • force_var: run variation even if valid wavefunction files already exist, useful e.g. when optimizing orbs, default: false.
  • skip_var: skip the extra read of the wavefunction when wavefunction files already exist, default: false.
  • var_sd: 🌴 include all singles and doubles excitation, i.e. at least CISD, default: false.
  • get_pair_contrib: 🌴 calculate occupied pair contribution, default: false.
  • time_sym: 🌴 uses time-reversal symmetry, recommended for singlets, default: false.
  • s2: 🌴 calculates s squared, default: false.
  • optorb: generates optimized orbitals FCIDUMP, default: false.
  • second_rejection: it uses 2nd criterion for choosing dets, useful when core excit allowed, default: false.
  • second_rejection_factor: default: false.
  • load_integrals_cache: loads FCIDUMP information from integrals_cache, default: false.
  • get_1rdm_csv, get_2rdm_csv: 🌱 calculates the density matrices, default: false. If it is true, the program outputs density matrices for the smallest eps_var in the csv format. For the two body density matrix, the p q r s columns represent a+_p a+_q a_r a_s.
  • get_green: 🌱 calculates the Green's function, default: false. If it is true, w_green is the real part of the frequency and n_green is the imaginary part. advanced_green selects G- (true) or G+ (false), default: false. The Green's function matrix is returned in csv format.
  • hc_server_mode: 🌱 operates as an H * c server, default: false. If it is true, the program serves as an RPC server for performing H * c matrix-vector multiplication after finishing the matrix reconstruction. It can work with any language that supports direct socket IO. A python client interface / demo is provided via hc_client.py. hc_client.py exposes a class called HcClient, which has three public methods: getN for getting the number of determinants, getCoefs for getting the coefficients array as a numpy array, and Hc(arr) which performs the matrix-vector multiplication on a numpy array of dtype either np.float64 or np.complex64 and returns the resulting numpy array of the same type. The HcClient accepts several optional construction options:
    • nProcs: default: 1, where the program uses all the cores on the master node. To run the hc_server across nodes, set nProcs to the number of nodes allocated to the job when creating the HcClient instance.
    • runtimePath: default: the current working directory. The config.json and FCIDUMP shall exist in the runtime path.
    • shciPath: default: the shci under the current working directory, probably needs to be changed to the actual path of the program.
    • port: default: 2018. Change it together with the value in src/solver/hc_server.h in case of a port conflict.
    • verbose: default: true.

Chemistry Block chem

The chem block contains input arguments specific to quantum chemistry systems.

  • point_group (required): supports C1, C2, Cs, Ci, C2v, C2h, Coov, D2, D2h, and Dooh.
  • irreps: an array of irreducible representations. If occupations are also given, they together determine the starting determinant, otherwise, the lowest orbitals are filled. occs_up and occs_dn when specified (outside chem block) have priority over irreps.
  • irrep_occs_up and irrep_occs_dn: occupation of each irreducible representation for up and down electrons respectively in the starting determinant. The lowest orbitals of each irrep. are filled. Ignored if irreps is not given.

Orbital Optimization

The optimization routine performs orthogonal transformations among the molecular orbitals to lower the variational energy. The arguments for orbital optimization should be specified in the optimization block in config.json.

  • natorb_iter: number of natural orbital iterations. Default: 1.
  • optorb_iter: number of iterations of full orbital optimization. Default: 20.
  • micro_iter: number of micro-iterations where the determinants are not re-selected. Default: 0.
  • rotation_matrix: whether to write out the orbital rotation matrix. Default: false.
  • method: what optimization method to use. See comparison of the various methods in this paper. Implementation details are in src/chem/optimization.cc.
  • amsgrad: AMSGrad method.
  • newton: Newton's method.
  • appnewton: Newton's method with the orbital Hessian matrix approximated by its diagonal.
  • bfgs: block matrix inversion with the BFGS method.
  • accelerate: use overshooting to accelerate convergence in optimization, default: false. Specific to method: app_newton and newton.
  • parameters block: optimization parameters specific to each method that can be fine tuned.

An example input that uses the accelerated diagonal Newton's method is:

  "optimization":{
          "method":"appnewton",
          "accelerate": true,
          "natorb_iter": 1,
          "optorb_iter":20
  }

An example input that uses the BFGS method with alternating macro- and micro-iterations is:

  "optimization":{
          "method":"bfgs",
          "natorb_iter": 1,
          "optorb_iter":10,
          "micro_iter": 1
  }

The FCIDUMP file is written out for natural orbitals as FCIDUMP_natorb and for optimized orbitals as FCIDUMP_optorb. By default, Arrow prints out FCIDUMP_optorb every 4 iterations if the energy goes down monotonically. FCIDUMP_optorb will not be written out if the current energy is higher than the lowest energy so far and will instead keep on iterating until it reaches some energy lower than the lowest so far. See the implementation logic in optimization_run() in src/solver/solver.h.

To see the energy convergence, grep for 'dE:' in the output:

Iter 1: natural orbitals E: -114.31880655 ndet: 256,228 dE: -0.01465866
Iter 2: optimized orbitals (bfgs) E: -114.31890043 ndet: 256,533 dE: -0.00009387
Iter 3: optimized orbitals (bfgs) E: -114.32011843 ndet: 267,198 dE: -0.00121801
Iter 4: optimized orbitals (bfgs) E: -114.32041848 ndet: 263,495 dE: -0.00030005
Iter 5: optimized orbitals (bfgs) E: -114.32052269 ndet: 263,439 dE: -0.00010421
Iter 6: optimized orbitals (bfgs) E: -114.32137590 ndet: 267,548 dE: -0.00085321
Iter 7: optimized orbitals (bfgs) E: -114.32164856 ndet: 265,856 dE: -0.00027266
Iter 8: optimized orbitals (bfgs) E: -114.32174660 ndet: 265,582 dE: -0.00009804
Iter 9: optimized orbitals (bfgs) E: -114.32219101 ndet: 264,690 dE: -0.00044441
Iter 10: optimized orbitals (bfgs) E: -114.32238104 ndet: 264,452 dE: -0.00019003

To see which iteration the printed FCIDUMP file corresponds to, grep for 'FCIDUMP' in addition to 'dE:'.

Active Space Calculations

One can partition the set of orbitals into active and virtual orbitals and restrict all excitations to be within the active space. To do this, set chem/active_space to true and specify the number of active orbitals in each irrep in chem/irrep_active_space and the program will automatically set the lowest energy orbitals in each irrep to be the active orbitals.

  "chem": {
    "point_group": "D2h",
    "active_space": true,
    "irrep_active_space": [6,3,3,1,6,3,3,1]
  },
  "reorder_orbs": false,

In the output, the active orbitals in each irrep will be listed.

Note that currently reorder_orbs needs to be set to false manually to ensure that the orbital order (based on one-determinant orbital energies) stays the same for the optimization run and the final run.

Since active-virtual rotations are included here, accelerated Newton's methods do not work well, but BFGS should still work. In addition, it is usually beneficial to calculate natural orbitals for all orbitals and not just for the active space before doing active space orbital optimization, since the starting HF orbitals can be far from optimal for the virtual space.