Vibrational Hamiltonians

Open in Colab

[1]:
!pip install block2==0.5.3rc16 -qq --progress-bar off --extra-index-url=https://block-hczhai.github.io/block2-preview/pypi/

In this tutorial, we explain how to perform DMRG for vibrational Hamiltonians to find vibrational energy levels and vibrational mode assignments. We use the Watson Hamiltonian and expand the Potential Energy Surface (PES) as a sextic force field in terms of normal mode coordinates about the equilibrium geometry:

\[\begin{split}\hat{H} = -\frac{1}{2}\sum_{i=1} \frac{\partial^2}{\partial Q_i^2} + \frac{1}{2!} \sum_{ij} F_{ij} Q_i Q_j + \frac{1}{3!} \sum_{ijk} F_{ijk} Q_i Q_j Q_k + \frac{1}{4!} \sum_{ijkl} F_{ijkl} Q_i Q_j Q_k Q_l \\ + \frac{1}{5!} \sum_{ijklm} F_{ijklm} Q_i Q_j Q_k Q_l Q_m + \frac{1}{6!} \sum_{ijklmn} F_{ijklmn} Q_i Q_j Q_k Q_l Q_m Q_n\end{split}\]

For details of the definition of the Hamiltonian and the wavefunction ansatz, see M. Sibaev, D. L. Crittenden. PyVCI: A flexible open-source code for calculating accurate molecular infrared spectra. Computer Physics Communications, 203, 290-297 (2016).

Force Field Parameters

First we need to load the force field parameters in the Hamiltonian. We consider \(\mathrm{N_2H_2}\) as an example, where the PES can be found in the PyPES package (https://github.com/dlc62/pypes-lib) and the reference energy levels computed using Vibrational Configuration Interaction (VCI) are given in M. Sibaev, D. L. Crittenden. The PyPES library of high quality semi‐global potential energy surfaces.” Journal of Computational Chemistry 29, 2200-2207 (2015).

We first download the force field file force_field.tar (without decompressing it). It will be automaticaly decompressed in the following Python script.

[2]:
!wget https://github.com/dlc62/pypes-lib/raw/master/force_field.tar
--2024-08-01 02:09:02--  https://github.com/dlc62/pypes-lib/raw/master/force_field.tar
Resolving github.com (github.com)... 140.82.113.4
Connecting to github.com (github.com)|140.82.113.4|:443... connected.
HTTP request sent, awaiting response... 302 Found
Location: https://raw.githubusercontent.com/dlc62/pypes-lib/master/force_field.tar [following]
--2024-08-01 02:09:03--  https://raw.githubusercontent.com/dlc62/pypes-lib/master/force_field.tar
Resolving raw.githubusercontent.com (raw.githubusercontent.com)... 185.199.108.133, 185.199.109.133, 185.199.110.133, ...
Connecting to raw.githubusercontent.com (raw.githubusercontent.com)|185.199.108.133|:443... connected.
HTTP request sent, awaiting response... 200 OK
Length: 21156352 (20M) [application/octet-stream]
Saving to: ‘force_field.tar’

force_field.tar     100%[===================>]  20.18M  --.-KB/s    in 0.1s

2024-08-01 02:09:03 (173 MB/s) - ‘force_field.tar’ saved [21156352/21156352]

Perform Vibrational DMRG

The following script will solve the 10 low-energy roots for the vibrational Hamiltonian of the \(\mathrm{N_2H_2}\) molecule (rotational effects are not included).

[3]:
import numpy as np
import math
import tarfile, bz2

pes_name = 'N2H2'
pes_n = None
nroots = 10
nexcit = 4 # max excitation number per mode

# Read force field parameters from PyPES
with tarfile.open(name='force_field.tar') as tf:
    fn = [fn for fn in tf.getnames() if fn.endswith('/molecule_%s.py.bz2' % pes_name)][0]
    pess, ipes, ff = {}, -1, None
    for l in bz2.decompress(tf.extractfile(fn).read()).decode('utf-8').split('\n'):
        if 'if state.pesn ==' in l:
            if ff is not None:
                pess[ipes]['ff'] = ff
            ipes = int(l.split('\"')[-2])
            pess[ipes], ff = {'nm': False, 'cmt': ''}, None
        elif 'state.normal_mode_ff' in l:
            pess[ipes]['nm'] = True
        elif 'state.ff_comment.append' in l:
            pess[ipes]['cmt'] += l.split('\'')[1] + '\n'
        elif 'ff.ff =' in l:
            ff = []
        elif 'else:' in l:
            if ff is not None:
                pess[ipes]['ff'] = ff
            ff = None
        elif ff is not None:
            xs = l.split('[')[1].split(']')[0].split(',')
            ff.append([int(x) for x in xs[:-1]] + [float(xs[-1])])
    pes = pess[pes_n if pes_n is not None else ipes]
    assert pes['nm']

derivs = [{tuple(l[:-1]): l[-1] for l in pes['ff'] if len(l) == i + 1} for i in range(max([len(l) for l in pes['ff']]))]
n_sites = max(derivs[2])[0] + 1
omega = np.array([derivs[2][(i, i)] for i in range(n_sites)]) ** 0.5
au2cm = 219474.6435

print('PES:', pes_name)
print(pes['cmt'])
print('Harmonic Frequencies (cm-1): %s' % ''.join(['%12.4f' % x for x in (omega * au2cm)]))
print('Harmonic ZPE         (cm-1): %12.4f' % np.sum(omega * au2cm / 2))

from pyblock2.driver.core import DMRGDriver, SymmetryTypes
driver = DMRGDriver(scratch="./tmp", symm_type=SymmetryTypes.SAny, n_threads=4)
driver.set_symmetry_groups()
q = driver.bw.SX()

# Table 1 in Computer Physics Communications 203, 290-297 (2016).
site_basis, site_ops = [], []
for k in range(n_sites):
    basis = [(q, nexcit)]
    pmat, qmat = np.zeros((nexcit, nexcit)), np.zeros((nexcit, nexcit))
    for i, j, n, dn in [(i, j, max(i, j), abs(i - j)) for i in range(nexcit) for j in range(nexcit)]:
        if dn == 0:
            pmat[i, j] = -omega[k] * (i + 0.5)
        elif dn == 2:
            pmat[i, j] = 0.5 * omega[k] * (n * (n - 1)) ** 0.5
        elif dn == 1:
            qmat[i, j] = (n  / (2 * omega[k])) ** 0.5
    ops = {"": np.identity(nexcit), "P": pmat, "Q": qmat}
    site_basis.append(basis)
    site_ops.append(ops)

# Set up vibrational Hamiltonian
driver.initialize_system(n_sites=n_sites, vacuum=q, target=q, hamil_init=False)
driver.ghamil = driver.get_custom_hamiltonian(site_basis, site_ops)
b = driver.expr_builder()

b.add_term("P", list(range(n_sites)), -0.5)
for p, idx, f in [(p, idx, f) for p, d in enumerate(derivs) for idx, f in d.items()]:
    fac = np.prod([math.factorial(g) for g in np.unique(idx, return_counts=True)[1]])
    b.add_term("Q" * p, idx, f / fac)

mpo = driver.get_mpo(b.finalize(adjust_order=True, fermionic_ops=''), cutoff=1E-75, iprint=1)

# Run DMRG
kets = [None for _ in range(nroots)]
energies = []
for ir in range(nroots):
    print('\n==== ROOT %4d ====' % ir)
    kets[ir] = driver.get_random_mps(tag="KET-%d" % ir, bond_dim=50, nroots=1)
    energies.append(driver.dmrg(mpo, kets[ir], n_sweeps=10, bond_dims=[50] * 4 + [100] * 4,
        noises=[1e-4] * 4 + [1e-5] * 2 + [0], thrds=[1e-10] * 8, dav_max_iter=300,
        proj_weights=[1.0] * ir, proj_mpss=kets[:ir], iprint=1, tol=1E-11))

# Energy levels and assignment
for ir in range(nroots):
    if ir == 0:
        print('ZPE = %12.4f cm-1' % (energies[ir] * au2cm), end=' ', flush=True)
    else:
        print('VIB = %12.4f cm-1' % ((energies[ir] - energies[0]) * au2cm), end=' ', flush=True)
    kets[ir] = driver.adjust_mps(kets[ir], dot=1)[0]
    csfs, coeffs = driver.get_csf_coefficients(kets[ir], cutoff=0.1, iprint=0)
    csf, weight = csfs[np.argmax(coeffs ** 2)], coeffs[np.argmax(coeffs ** 2)] ** 2
    print(' Max weight = %10.4f Assignment = %s' % (weight, csf))
PES: N2H2
  trans-N2H2 molecule. Saved in normal mode coordinates from PESN 2.
  Reference:
  Martin, J. M. L.; Taylor, P. R. Spectrochim. Acta Part A: Mol. Biomol. Spectrosc. 1997, 53(8), 1039-1050.

Harmonic Frequencies (cm-1):    1328.4271   1350.2811   1558.4338   1621.8428   3269.8469   3301.5852
Harmonic ZPE         (cm-1):    6215.2084

Build MPO | Nsites =     6 | Nterms =        348 | Algorithm = FastBIP | Cutoff = 1.00e-75
 Site =     0 /     6 .. Mmpo =     4 DW = 0.00e+00 NNZ =        4 SPT = 0.0000 Tmvc = 0.000 T = 0.014
 Site =     1 /     6 .. Mmpo =    11 DW = 0.00e+00 NNZ =       16 SPT = 0.6364 Tmvc = 0.000 T = 0.027
 Site =     2 /     6 .. Mmpo =    23 DW = 0.00e+00 NNZ =       66 SPT = 0.7391 Tmvc = 0.000 T = 0.034
 Site =     3 /     6 .. Mmpo =    19 DW = 0.00e+00 NNZ =      104 SPT = 0.7620 Tmvc = 0.000 T = 0.013
 Site =     4 /     6 .. Mmpo =     7 DW = 0.00e+00 NNZ =       36 SPT = 0.7293 Tmvc = 0.000 T = 0.016
 Site =     5 /     6 .. Mmpo =     1 DW = 0.00e+00 NNZ =        7 SPT = 0.0000 Tmvc = 0.000 T = 0.007
Ttotal =      0.112 Tmvc-total = 0.001 MPO bond dimension =    23 MaxDW = 0.00e+00
NNZ =          233 SIZE =          878 SPT = 0.7346

Rank =     0 Ttotal =      0.196 MPO method = FastBipartite bond dimension =      23 NNZ =          233 SIZE =          878 SPT = 0.7346

==== ROOT    0 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.162 | E =       0.0278175746 | DW = 1.50629e-20

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.230 | E =       0.0278175746 | DE = -4.86e-17 | DW = 2.94238e-20

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.303 | E =       0.0278175746 | DE = -4.54e-15 | DW = 9.66411e-21

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.406 | E =       0.0278175746 | DE = 4.55e-15 | DW = 1.94996e-20

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.591 | E =       0.0278175746 | DE = -4.53e-15 | DW = 1.11346e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.681 | E =       0.0278175746 | DE = 4.52e-15 | DW = 5.49483e-21

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      0.789 | E =       0.0278175746 | DE = -4.55e-15 | DW = 1.43430e-20


==== ROOT    1 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.510 | E =       0.0337897813 | DW = 8.91046e-16

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.609 | E =       0.0337897813 | DE = -5.29e-12 | DW = 7.20577e-18

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.717 | E =       0.0337897813 | DE = 0.00e+00 | DW = 1.43698e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.819 | E =       0.0337897813 | DE = -2.08e-17 | DW = 3.82610e-18

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.923 | E =       0.0337897813 | DE = 1.39e-17 | DW = 1.50304e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.030 | E =       0.0337897813 | DE = 9.71e-17 | DW = 2.46714e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.124 | E =       0.0337897813 | DE = 1.39e-17 | DW = 1.36100e-20


==== ROOT    2 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.372 | E =       0.0336294440 | DW = 6.32649e-12

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.482 | E =       0.0336294440 | DE = -3.26e-12 | DW = 1.32102e-17

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.612 | E =       0.0336294440 | DE = 0.00e+00 | DW = 8.63491e-21

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.718 | E =       0.0336294440 | DE = -6.94e-18 | DW = 1.09199e-17

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.840 | E =       0.0336294440 | DE = -6.94e-18 | DW = 1.24511e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.955 | E =       0.0336294440 | DE = 6.94e-18 | DW = 1.33919e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.074 | E =       0.0336294440 | DE = 1.39e-17 | DW = 9.75384e-21


==== ROOT    3 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.336 | E =       0.0347691843 | DW = 2.71296e-15

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.453 | E =       0.0347691843 | DE = -9.19e-13 | DW = 4.24397e-20

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.546 | E =       0.0347691843 | DE = -6.94e-18 | DW = 1.57276e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.633 | E =       0.0347691843 | DE = -3.47e-17 | DW = 2.62031e-20

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.740 | E =       0.0347691843 | DE = 0.00e+00 | DW = 1.64701e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.824 | E =       0.0347691843 | DE = 3.47e-17 | DW = 1.24025e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      0.913 | E =       0.0347691843 | DE = -1.39e-17 | DW = 1.95745e-20


==== ROOT    4 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.221 | E =       0.0350118953 | DW = 2.73902e-17

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.320 | E =       0.0350118953 | DE = -1.27e-12 | DW = 2.21093e-19

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.426 | E =       0.0350118953 | DE = 0.00e+00 | DW = 8.90919e-21

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.520 | E =       0.0350118953 | DE = -6.94e-18 | DW = 2.29959e-18

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.624 | E =       0.0350118953 | DE = 0.00e+00 | DW = 1.89815e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.738 | E =       0.0350118953 | DE = 0.00e+00 | DW = 1.07261e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      0.846 | E =       0.0350118953 | DE = 6.94e-18 | DW = 9.05114e-21


==== ROOT    5 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.313 | E =       0.0393996116 | DW = 2.56369e-15

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.424 | E =       0.0393996116 | DE = -5.93e-12 | DW = 4.63254e-18

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.538 | E =       0.0393996116 | DE = 0.00e+00 | DW = 1.62961e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.647 | E =       0.0393996116 | DE = 4.16e-17 | DW = 1.40207e-17

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.773 | E =       0.0393996116 | DE = -6.94e-18 | DW = 1.72877e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.890 | E =       0.0393996116 | DE = -2.78e-17 | DW = 1.85333e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      0.998 | E =       0.0393996116 | DE = 0.00e+00 | DW = 1.03031e-20


==== ROOT    6 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.511 | E =       0.0397070198 | DW = 1.26937e-13

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.628 | E =       0.0397070198 | DE = -9.73e-12 | DW = 1.34647e-17

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.770 | E =       0.0397070198 | DE = 1.39e-17 | DW = 1.20041e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.884 | E =       0.0397070198 | DE = -6.94e-18 | DW = 1.12444e-17

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.003 | E =       0.0397070198 | DE = 6.94e-18 | DW = 1.36156e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.124 | E =       0.0397070198 | DE = 0.00e+00 | DW = 1.99115e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.240 | E =       0.0397070198 | DE = -6.94e-18 | DW = 1.60524e-20


==== ROOT    7 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.430 | E =       0.0398511574 | DW = 3.70941e-15

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.561 | E =       0.0398511574 | DE = -4.99e-12 | DW = 2.45878e-20

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.691 | E =       0.0398511574 | DE = -1.39e-17 | DW = 2.44316e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.823 | E =       0.0398511574 | DE = -1.39e-17 | DW = 3.23054e-19

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      0.957 | E =       0.0398511574 | DE = -6.94e-18 | DW = 2.15729e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.083 | E =       0.0398511574 | DE = 3.47e-17 | DW = 1.74024e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.218 | E =       0.0398511574 | DE = 0.00e+00 | DW = 1.91998e-20


==== ROOT    8 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.476 | E =       0.0406183665 | DW = 2.45735e-15

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.607 | E =       0.0406183665 | DE = -5.63e-12 | DW = 2.86690e-17

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.754 | E =       0.0406183665 | DE = 6.94e-18 | DW = 1.81005e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.888 | E =       0.0406183665 | DE = -6.94e-18 | DW = 2.33180e-17

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.028 | E =       0.0406183665 | DE = 2.08e-17 | DW = 2.37587e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.194 | E =       0.0406183665 | DE = -9.02e-17 | DW = 2.19444e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.340 | E =       0.0406183665 | DE = 1.39e-17 | DW = 2.59923e-20


==== ROOT    9 ====

Sweep =    0 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.588 | E =       0.0395978303 | DW = 5.53041e-11

Sweep =    1 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.748 | E =       0.0395978303 | DE = -8.83e-12 | DW = 4.70079e-20

Sweep =    2 | Direction =  forward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      0.903 | E =       0.0395978303 | DE = 6.94e-18 | DW = 1.44256e-20

Sweep =    3 | Direction = backward | Bond dimension =   50 | Noise =  1.00e-04 | Dav threshold =  1.00e-10
Time elapsed =      1.048 | E =       0.0395978303 | DE = 6.94e-18 | DW = 3.29833e-20

Sweep =    4 | Direction =  forward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.203 | E =       0.0395978303 | DE = -1.39e-17 | DW = 1.81447e-20

Sweep =    5 | Direction = backward | Bond dimension =  100 | Noise =  1.00e-05 | Dav threshold =  1.00e-10
Time elapsed =      1.353 | E =       0.0395978303 | DE = 0.00e+00 | DW = 1.83984e-20

Sweep =    6 | Direction =  forward | Bond dimension =  100 | Noise =  0.00e+00 | Dav threshold =  1.00e-10
Time elapsed =      1.505 | E =       0.0395978303 | DE = 0.00e+00 | DW = 1.61162e-20

ZPE =    6105.2523 cm-1  Max weight =     0.9705 Assignment = [0 0 0 0 0 0]
VIB =    1310.7479 cm-1  Max weight =     0.9699 Assignment = [0 1 0 0 0 0]
VIB =    1275.5580 cm-1  Max weight =     0.9732 Assignment = [1 0 0 0 0 0]
VIB =    1525.7021 cm-1  Max weight =     0.9334 Assignment = [0 0 1 0 0 0]
VIB =    1578.9710 cm-1  Max weight =     0.9524 Assignment = [0 0 0 1 0 0]
VIB =    2541.9634 cm-1  Max weight =     0.9322 Assignment = [2 0 0 0 0 0]
VIB =    2609.4317 cm-1  Max weight =     0.9226 Assignment = [0 2 0 0 0 0]
VIB =    2641.0663 cm-1  Max weight =     0.8964 Assignment = [0 3 0 0 0 0]
VIB =    2809.4492 cm-1  Max weight =     0.7351 Assignment = [0 1 1 0 0 0]
VIB =    2585.4674 cm-1  Max weight =     0.9687 Assignment = [1 1 0 0 0 0]