Quantum Package and Champ

Create a file called h2o.xyz: with the geometry of the water molecule:

3
Water
O       0.                     0.   0.
H      -0.756950272703377558   0.  -0.585882234512562827
H       0.756950272703377558   0.  -0.585882234512562827

Store the pseudopotential parameters in a file named PSEUDO:

H GEN 0 0
3
 1.000000000000 1 25.000000000000
25.000000000000 3 10.821821902641
-8.228005709676 2  9.368618758833

O GEN 2 1
3
6.00000000 1 9.29793903
55.78763416 3 8.86492204
-38.81978498 2 8.62925665
1
38.41914135 2 8.71924452

Store the basis set parameters in a file named BASIS:

HYDROGEN
s 3
1 6.46417546   0.063649375945
2 1.13891461   0.339233210576
3 0.28003249   0.702654522063
s 1
1 0.05908405   1.00000000
p 1
1 0.51368060   1.00000000

OXYGEN
s 9
1 0.125346     0.055741
2 0.268022     0.304848
3 0.573098     0.453752
4 1.225429     0.295926
5 2.620277     0.019567
6 5.602818     -0.128627
7 11.980245     0.012024
8 25.616801     0.000407
9 54.775216     -0.000076
s 1
1 0.258551     1.000000
p 9
1 0.083598     0.044958
2 0.167017     0.150175
3 0.333673     0.255999
4 0.666627     0.281879
5 1.331816     0.242835
6 2.660761     0.161134
7 5.315785     0.082308
8 10.620108     0.039899
9 21.217318     0.004679
p 1
1 0.267865     1.000000
d 1
1 1.232753     1.000000

First, we can compute with QP the energies of the single-determinant wave functions with the 2 different sets of MOs.

qp set_file h2o_hf
qp run print_energy

qp set_file h2o_dft
qp run print_energy

These commands return the energy of the wavefunction contained in the EZFIO database. These values will be useful for checking that the QMC setup is OK. You should obtain the energies:

We will now convert the TREXIO files into input files suitable for CHAMP:

Create a new directory named H2O_HF and copy the TREXIO file h2o_hf.trexio into it. Go inside this directory and run

python3 ~filippi/Tutorial-QMC-School/trex2champ.py --trex "h2o_hf.trexio" \
                       --motype  "Canonical" \
                       --backend "HDF5" \
                       --basis_prefix "BFD-cc-pVDZ" \
                       --lcao \
                       --geom \
                       --basis \
                       --ecp \
                       --det

Many files were created. Now, create a directory named pool, and move some files into the pool:

mkdir pool
mv *.xyz *bfinfo BFD-* ECP* pool

You can now create an input file for CHAMP vmc_h2o_hf.inp :

%module general
    title           'H2O HF calculation'
    pool            './pool/'
    pseudopot       ECP
    basis           BFD-cc-pVDZ
    mode            'vmc_one_mpi1'
%endmodule


load molecule        $pool/champ_v2_h2o_hf_geom.xyz
load basis_num_info  $pool/champ_v2_h2o_hf_with_g.bfinfo

load orbitals        champ_v2_h2o_hf_orbitals.lcao
load determinants    champ_v2_h2o_hf_determinants.det
load jastrow         jastrow.start

%module electrons
    nup           4
    nelec         8
%endmodule


%module blocking_vmc
    vmc_nstep     20
    vmc_nblk      20000
    vmc_nblkeq    1
    vmc_nconf_new 0
%endmodule

Create the file for the Jastrow factor as follows, and save it as jastrow.start:

jastrow_parameter   1
  0  0  0           norda,nordb,nordc
   0.60000000   0.00000000     scalek,a21
   0.00000000   0.00000000   (a(iparmj),iparmj=1,nparma)
   0.00000000   0.00000000   (a(iparmj),iparmj=1,nparma)
   0.00000000   1.00000000   (b(iparmj),iparmj=1,nparmb)
 (c(iparmj),iparmj=1,nparmc)
 (c(iparmj),iparmj=1,nparmc)
end

This files implies that there is no Jastrow factor (\(\exp(J)=1\)).

Create the submission script as presented in section 3, and submit the job. You should obtain the Hartree-Fock energy.

Now reproduce the same steps for the TREXIO file containing the DFT orbitals in directory H2O_DFT.

The energies obtained with VMC without the Jastrow factor should be the same as those computed by QP at the beginning of this section.

The Jastrow factor depends on the electronic (\(\mathbf{r}\)) and nuclear (\(\mathbf{R}\)) coordinates. Its defined as \(\exp(J(\mathbf{r},\mathbf{R}))\), where

\[ J = f_{en} + f_{ee} + f_{een} \]

Electron-nucleus and electron-electron: \(R={1-e^{-\kappa r} \over \kappa}\)

\[ f_{en} = \sum_{i=1}^{N_{\rm elec}} \sum_{\alpha=1}^{N_{\rm nuc}} \left( {a_1 R_{i\alpha} \over 1+a_2R_{i\alpha}} + \sum_{p=2}^{N^a_{\rm ord}} a_{p+1} R_{i\alpha}^p \right) \]

\[ f_{ee} = \sum_{i=2}^{N_{\rm elec}} \sum_{j=1}^{i-1} \left( {b_1 R_{ij} \over 1+b_2R_{ij}} + \sum_{p=2}^{N^b_{\rm ord}} b_{p+1} R_{ij}^p \right) \]

Electron-electron-nucleus: \(R=\exp\left(-\kappa r \right)\)

\[ f_{een} = \sum_{i=2}^{N_{\rm elec}} \sum_{j=1}^{i-1} \sum_{\alpha=1}^{N_{\rm nuc}} \sum_{p=2}^{N^c_{\rm ord}} \sum_{k=p-1}^0 \sum_{l=l_{\rm max}}^0 c_n R_{ij}^k (R_{i\alpha}^l+R_{j\alpha}^l) (R_{i\alpha}R_{j\alpha})^m \]

where \(m={p-k-l \over 2}\)

• Typically \(N^a_{\rm ord}=N^b_{\rm ord}=5\). If \(f_{een}\) is included, \(N^c_{\rm ord}=5\).
• Dependence among \(\{c_n\}\) \(\rightarrow\) \(f_{een}\) does not contribute to cusp-conditions
• \(f_{en}\) and \(f_{een}\): different \(\{a_n\}\) and \(\{c_n\}\) for different atom types

Let us start to run a DMC simulation with the HF orbitals and the optimal Jastrow factor you have just generated.

Create a new directory and copy the wave function TREXIO info and the optimal Jastrow factor (for simplicity, pick the last one).

First, generate an input file as before where you read the wave function files (careful to load the new Jastrow factor) and perform a short VMC calculation to generate the walkers for DMC.

To shorten the VMC run, you can choose a small vmc_nblk in the main input file and modify vmc_nconf_new to be the number of walkers per core you wish. Here, we use the same values as for the starting iterations of the Jastrow factor optimization:

%module blocking_vmc
    vmc_nstep     20
    vmc_nblk      200
    vmc_nblkeq    1
    vmc_nconf_new 100
%endmodule

This will generate 100 walkers per core (vmc_nconf_new) by writing the coordinates of a walker every \(20 \times 200 / 100\) steps. Since the correlation time is less than 2 step in VMC, your walkers will be decorrelated.

A bunch of mc_configs_newX files will appear in your directory, each containing 100 walkers.

cat mc_configs_new* >> mc_configs
rm mc_configs_new*

mc_configs contains now all walkers.

Generate a DMC input

%module blocking_dmc
    dmc_nstep     60
    dmc_nblk      40
    dmc_nblkeq    1
    dmc_nconf     100
%endmodule

%module dmc
    tau           0.05
    etrial      -17.240
    icasula      -1
%endmodule

You also need to change the mode keyword in the input file:

mode            'dmc_one_mpi1'

within the general module.

Some debug files are being created, that you can just erase.

rm problem*
rm mc_configs_new*

To look at the energy, you can do

grep '( 100) =' dmc*out

In the last column, you have the correlation time.

Also perform another calculation with a smaller time step.

Repeat the optimization and DMC calculation with the DFT orbitals. Compare the VMC and DMC energies.

Finally, starting from the DFT orbitals and the optimal two-body Jastrow optimize the full wave function (Jastrow and orbitals).

To this aim, set ioptorb to 1 in the optwf module.

ioptorb     1