Improving convergence in CheMPS2

[quanp]

Hi Sebastian,

I'm trying to do a DMRG-CI calculation with 22 electrons in 22 orbitals. I used both Block and CheMPS2. In Block, I have an energy of -4686.37279 Ha, while in CheMPS2 I only got -4686.114229 Ha, indicating that the calculation with CheMPS2 was trapped at a local minimum. How can I improve the result of CheMPS2?

You can find the FCIDUMP, input and output files at: https://github.com/quanp/NiFe

The input of CheMPS2 is FCIDUMP = /scratch/leuven/305/vsc30513/sing_Ni_tpss.work/FCIDUMP GROUP = 0

MULTIPLICITY = 1 NELECTRONS = 22 IRREP = 0 EXCITATION = 0

SWEEP_STATES = 250, 500, 500 SWEEP_ENERGY_CONV = 1e-8, 1e-8, 1e-8 SWEEP_MAX_SWEEPS = 10, 10, 10 SWEEP_NOISE_PREFAC = 0.05, 0.05, 0.00 SWEEP_DVDSON_RTOL = 1e-5, 1e-10, 1e-10

NOCC = 0 NACT = 22 NVIR = 0

SCF_STATE_AVG = FALSE SCF_DIIS_THR = 1e-2 SCF_GRAD_THR = 1e-7 SCF_MAX_ITER = 1 SCF_ACTIVE_SPACE = L

CASPT2_CALC = FALSE CASPT2_ORBS = A CASPT2_IPEA = 0.25 CASPT2_IMAG = 0.1 CASPT2_CHECKPT = TRUE

PRINT_CORR = TRUE

TMP_FOLDER = /scratch/leuven/305/vsc30513/

The input of Block is nelec 22 orbitals FCIDUMP-nife spin 0 irrep 1 schedule default maxM 500 maxiter 10 point_group c1 hf_occ integral

[SebWouters]

Hi Quan,

I don't think you localize the orbitals in block, or perhaps with a split localization, although I'm not sure. With the 'L' option for SCF_ACTIVE_SPACE, you allow all orbitals to be localized. So I changed SCF_ACTIVE_SPACE to

SCF_ACTIVE_SPACE = I

so that the input orbitals are used. A second thing is that I think that block doesn't use the D=250, although I might be wrong. So I left that sweep instruction out as well:

SWEEP_STATES = 500, 500
SWEEP_ENERGY_CONV = 1e-8, 1e-8
SWEEP_MAX_SWEEPS = 10, 10
SWEEP_NOISE_PREFAC = 0.05, 0.00
SWEEP_DVDSON_RTOL = 1e-5, 1e-10

and then after the first couple of sweeps I already get an energy of -4686.33233792072 and it's still decreasing.

Best wishes, Sebastian

[SebWouters]

Hi Quan,

So I looked at your block output, and it indeed uses D=250 as well. I don't think they localize, but they reorder the orbitals always according to the Fiedler vector of the exchange matrix. At first sight, the remaining difference -4686.33 (chemps2) vs. -4686.37 (block) is due to the ordering, although I have to check.

Best wishes, Sebastian

[SebWouters]

Hi Quan,

In commit 4389aa3017afdb9d5ed6082f7e4df6cb729d5af0 (and onwards) you can now specify

SCF_ACTIVE_SPACE = F

This option does not rotate the orbitals, but only reorders them. For your test case, it gives (initially) the output output-quick.txt which yields the same energy. I didn't let it run fully, but it shows that the difference between chemps2 and block came from the orbital ordering.

Best wishes, Sebastian

[quanp]

Hi Sebastian,

Thank you for your support. This shows that using localized orbitals, which was used in many studies, is not always the best option.

Best regards, Quan

[quanp]

Hi Sebastian,

Sorry for reopening this question as I feel that there is still room for improvement for the convergence of CheMPS2, especially when using large active space. For example, I did a calculation with 24 electrons in 24 orbitals and I always obtain much higher total energy with CheMPS2 than with Block. Increasing D to 1000 (and even 2000) doesn't help. You can see the input and output files at https://github.com/quanp/NiFe

Best regards, Quan

[SebWouters]

Hi Quan,

I think I know what the matter is. The Fiedler vector is in CheMPS2 currently only for orbitals within the same irrep. My guess is that it's globally in block. I implemented it this way because my localization scheme doesn't mix orbitals between different irreps, but that doesn't matter for the orbital reordering.

I'll notify you once it's adjusted. I'll also add a print statement of the new ordering so we can check the actual orbital reordering.

Best wishes, Sebastian

[SebWouters]

Hi Quan,

It should be fixed with the latest commit 79fb7c6f2b4cc2d628e5593fc8aa548c572a9e83.

I get the orbital reordering [ 6, 13, 17, 22, 9, 2, 11, 5, 14, 19, 8, 1, 12, 7, 0, 4, 23, 18, 16, 20, 15, 21, 10, 3 ], which is exactly the same as block, if you take point group order in chemps2 (Ag Bg Au Bu) and block (Ag Au Bu Bg) into account, and notice that my orbital reordering is then the reverse of block's ordering.

Also add a little bit more noise please to avoid getting stuck in local minima. I used the input:

FCIDUMP = FCIDUMP
GROUP   = 6

MULTIPLICITY = 1
NELECTRONS   = 24
IRREP        = 0
EXCITATION   = 0

SWEEP_STATES       = 250,  500,  500  
SWEEP_ENERGY_CONV  = 1e-9, 1e-9, 1e-9 
SWEEP_MAX_SWEEPS   = 10,   10,   10
SWEEP_NOISE_PREFAC = 0.2,  0.2,  0.0
SWEEP_DVDSON_RTOL  = 1e-3, 1e-3, 1e-10

NOCC =  0,  0,  0,  0
NACT =  7,  7,  5,  5
NVIR =  0,  0,  0,  0

SCF_STATE_AVG    = FALSE
SCF_DIIS_THR     = 1e-2
SCF_GRAD_THR     = 1e-7
SCF_MAX_ITER     = 1
SCF_ACTIVE_SPACE = F

CASPT2_CALC    = FALSE
CASPT2_ORBS    = A
CASPT2_IPEA    = 0.25
CASPT2_IMAG    = 0.2
CASPT2_CHECKPT = TRUE

PRINT_CORR = TRUE

TMP_FOLDER = /tmp

Which gives the output energies:

******************************************************************
***  Information on completed instruction 0:
***     The reduced virtual dimension DSU(2)               = 250
***     Minimum energy encountered during all instructions = -3568.52268030624
***     Minimum energy encountered during the last sweep   = -3568.52268030624
***     Maximum discarded weight during the last sweep     = 1.3985333969913e-05
******************************************************************
***  Information on completed instruction 1:
***     The reduced virtual dimension DSU(2)               = 500
***     Minimum energy encountered during all instructions = -3568.52361986675
***     Minimum energy encountered during the last sweep   = -3568.52361986675
***     Maximum discarded weight during the last sweep     = 5.60788687687778e-06
******************************************************************
***  Information on completed instruction 2:
***     The reduced virtual dimension DSU(2)               = 500
***     Minimum energy encountered during all instructions = -3568.52363556337
***     Minimum energy encountered during the last sweep   = -3568.52363554665
***     Maximum discarded weight during the last sweep     = 6.18133592786532e-06
******************************************************************

Please also note that block's davidson tolerance parameter is the square of mine. While mine is the norm, in block it's the norm squared. So if you say SWEEP_DVDSON_RTOL = 1e-10, it corresponds to 1e-20 in block!

Best wishes, Sebastian

[quanp]

Hi Sebastian,

Thank you for fixing the problem.

Best regards, Quan