Determinant Tools
N2(10o, 7e) (STO3G)
Ground-state DMRG (N2 STO3G) with C++ optimized core functions:
import numpy as np
from pyblock3.algebra.mpe import MPE
from pyblock3.hamiltonian import Hamiltonian
from pyblock3.fcidump import FCIDUMP
fd = 'data/N2.STO3G.FCIDUMP'
bond_dim = 250
hamil = Hamiltonian(FCIDUMP(pg='d2h').read(fd), flat=True)
mpo = hamil.build_qc_mpo()
mpo, _ = mpo.compress(cutoff=1E-9, norm_cutoff=1E-9)
mps = hamil.build_mps(bond_dim)
dmrg = MPE(mps, mpo, mps).dmrg(bdims=[bond_dim], noises=[1E-6, 0], dav_thrds=[1E-4], iprint=2, n_sweeps=10)
ener = dmrg.energies[-1]
print("Energy = %20.12f" % ener)
Check MPO:
print('MPO = ', mpo.show_bond_dims())
mpo, error = mpo.compress(cutoff=1E-12)
print('MPO = ', mpo.show_bond_dims())
Check MPS:
print('MPS = ', mps.show_bond_dims())
print(mps.norm())
Check ground-state energy:
mps @ (mpo @ mps)
MPO-MPS Contraction
mps.opts = {}
hmps = mpo @ mps
print(hmps.show_bond_dims())
print(hmps.norm())
The result MPS can be compressed:
hmps, _ = hmps.compress(cutoff=1E-12)
print(hmps.show_bond_dims())
print(hmps.norm())
MPO truncation to bond dimension 50:
cmpo, cps_error = mpo.compress(max_bond_dim=50, cutoff=1E-12)
print('error = ', cps_error)
print(cmpo.show_bond_dims())
Apply contracted MPO to MPS:
hmps = cmpo @ mps
print(hmps.show_bond_dims())
print(hmps.norm())
Determinants
Using SliceableTensor:
smps = mps.to_non_flat().to_sliceable()
print(smps[0])
print('-'*20)
print(smps[0][:, 2:, 2])
print('-'*20)
print(smps[0][:, :2, 2].infos)
print('-'*20)
print(smps[0])
print('-'*20)
print(smps.amplitude([3, 3, 0, 3, 0, 3, 3, 3, 3, 0]))
If the determinant belongs to another symmetry sector, the overlap should be zero:
print(smps.amplitude([3, 3, 0, 0, 0, 3, 3, 3, 3, 0]))
Check the overlap for all doubly occupied determinants:
import itertools
coeffs = []
for ocp in itertools.combinations(range(10), 7):
det = [0] * mps.n_sites
for t in ocp:
det[t] = 3
tt = time.perf_counter()
coeffs.append(smps.amplitude(det))
print(np.array(det), "%10.5f" % coeffs[-1])
Check the sum of probabilities:
print((np.array(coeffs) ** 2).sum())