{ "cells": [ { "cell_type": "markdown", "metadata": { "id": "wD4caskLRUBB" }, "source": [ "# MPS Import (Renormalizer)\n", "\n", "[![Open in Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/block-hczhai/block2-preview/blob/master/docs/source/tutorial/mps-import-renormalizer.ipynb)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "PB7VbXEYPFfb", "outputId": "fa3ca563-23af-4ddf-9025-c3587f12e325" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " Preparing metadata (setup.py) ... \u001b[?25l\u001b[?25hdone\n", " Building wheel for renormalizer (setup.py) ... \u001b[?25l\u001b[?25hdone\n" ] } ], "source": [ "!pip install block2==0.5.4rc6 -qq --progress-bar off --extra-index-url=https://block-hczhai.github.io/block2-preview/pypi/\n", "!pip install pyscf==2.8.0 -qq --progress-bar off\n", "!pip install renormalizer==0.0.11 -qq --progress-bar off" ] }, { "cell_type": "markdown", "metadata": { "id": "3i9j5QaF8vSW" }, "source": [ "## Introduction\n", "\n", "In this tutorial we explain how to import an MPS from ``Renormalizer`` to ``block2``. ``Renormalizer`` (https://github.com/shuaigroup/Renormalizer) is a Python tensor network package focused on electron-phonon quantum dynamics. This can be used for computing the overlap between a ``block2`` MPS (in the ``SZ`` or ``SU2`` mode) and an external MPS, which may be needed in quantum information or in Quantum Monte Carlo calculations (https://arxiv.org/abs/2405.05440).\n", "\n", "## The external MPS\n", "\n", "We use the hydrogen chain as a test system. The required integrals for quantum chemistry Hamiltonian can be generated using ``pyscf``." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "id": "hOEjXXRl94Tk" }, "outputs": [], "source": [ "from pyscf import gto, scf, lo, tools\n", "import numpy as np\n", "from pyblock2._pyscf.ao2mo import integrals as itg\n", "\n", "N, R = 6, 3.0\n", "mf = scf.RHF(gto.M(atom=[['H', (0, 0, i * R)] for i in range(N)], basis='sto3g', verbose=3, symmetry='c1', unit='bohr'))\n", "mf.mo_coeff = lo.orth.lowdin(mf.mol.intor('cint1e_ovlp_sph'))\n", "ncas, n_elec, spin, ecore, h1e, h2e, orb_sym = itg.get_rhf_integrals(mf, g2e_symm=8)" ] }, { "cell_type": "markdown", "metadata": { "id": "xQA7GKYq_Cxx" }, "source": [ "We can find the ground state using the ``Renormalizer`` package with ``SZ`` symmetry, which is used as an external MPS." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "b2DTsUKaAAKg", "outputId": "cc72483c-ba00-41e0-f3fd-20b18c1cebe2" }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "2025-06-14 20:16:30,306[INFO] nuclear repulsion: 2.899999999999999\n", "INFO:renormalizer.model.h_qc:nuclear repulsion: 2.899999999999999\n", "2025-06-14 20:16:30,309[INFO] spin norbs: 12\n", "INFO:renormalizer.model.h_qc:spin norbs: 12\n", "2025-06-14 20:16:32,492[INFO] optimization method: 2site\n", "INFO:renormalizer.mps.gs:optimization method: 2site\n", "2025-06-14 20:16:32,495[INFO] e_rtol: 1e-06\n", "INFO:renormalizer.mps.gs:e_rtol: 1e-06\n", "2025-06-14 20:16:32,497[INFO] e_atol: 1e-08\n", "INFO:renormalizer.mps.gs:e_atol: 1e-08\n", "2025-06-14 20:16:32,499[INFO] procedure: [[30, 0.4], [30, 0.2], [30, 0.1], [30, 0], [30, 0], [30, 0], [30, 0]]\n", "INFO:renormalizer.mps.gs:procedure: [[30, 0.4], [30, 0.2], [30, 0.1], [30, 0], [30, 0], [30, 0], [30, 0]]\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "mpo_bond_dims:[1, 4, 16, 37, 54, 68, 81, 68, 54, 37, 16, 4, 1]\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "2025-06-14 20:16:35,500[INFO] DMRG has converged!\n", "INFO:renormalizer.mps.gs:DMRG has converged!\n", "2025-06-14 20:16:35,575[INFO] mps current size: 48.4KiB, Matrix product bond dim:[1, 2, 4, 8, 16, 30, 30, 30, 16, 8, 4, 2, 1]\n", "INFO:renormalizer.mps.gs:mps current size: 48.4KiB, Matrix product bond dim:[1, 2, 4, 8, 16, 30, 30, 30, 16, 8, 4, 2, 1]\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "lowest energy from Renormalizer MPS: -2.957646039005119\n" ] } ], "source": [ "import tempfile\n", "from renormalizer import Model, Mps, Mpo, optimize_mps\n", "from renormalizer.model import h_qc\n", "from renormalizer.utils import log\n", "from logging import INFO\n", "\n", "np.product = np.prod\n", "logger = log.logging.getLogger(\"renormalizer\")\n", "logger.setLevel(INFO)\n", "\n", "ftmp = tempfile.NamedTemporaryFile()\n", "tools.fcidump.from_mo(mf.mol, ftmp.name, mf.mo_coeff)\n", "h1e_spin, h2e_spin, nuc = h_qc.read_fcidump(ftmp.name, ncas)\n", "basis, ham_terms = h_qc.qc_model(h1e_spin, h2e_spin)\n", "model = Model(basis, ham_terms)\n", "mpo = Mpo(model)\n", "print(f\"mpo_bond_dims:{mpo.bond_dims}\")\n", "\n", "M = 30\n", "procedure = [[M, 0.4], [M, 0.2], [M, 0.1], [M, 0], [M, 0], [M, 0], [M, 0]]\n", "mps = Mps.random(model, mf.mol.nelec, M, percent=1.0)\n", "mps.optimize_config.procedure = procedure\n", "mps.optimize_config.method = \"2site\"\n", "energies, mps = optimize_mps(mps.copy(), mpo)\n", "gs_e = min(energies)+nuc\n", "print(f'lowest energy from Renormalizer MPS: {gs_e}')" ] }, { "cell_type": "markdown", "metadata": { "id": "bV0aslrSAhfY" }, "source": [ "The following code is used to convert the external MPS from a spin-orbital basis to fermionic spatial orbtial basis with the left-canonicalized format." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "gS_QmkjBBFP2", "outputId": "24722d90-8b5c-4502-8238-a707abdd0ec0" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "bond dimensions: [1, 4, 16, 30, 16, 4]\n", "quantum numbers: [[[0, 0]], [[1, 0], [0, 1], [1, 1], [0, 0]], [[0, 1], [0, 1], [1, 2], [1, 2], [2, 1], [2, 1], [0, 0], [1, 1], [1, 1], [1, 1], [1, 1], [2, 0], [0, 2], [2, 2], [1, 0], [1, 0]], [[0, 1], [1, 2], [1, 2], [1, 2], [1, 2], [2, 1], [2, 1], [2, 1], [2, 1], [3, 1], [3, 1], [1, 1], [1, 1], [1, 1], [1, 1], [0, 3], [2, 0], [2, 0], [3, 0], [2, 3], [0, 2], [0, 2], [2, 2], [2, 2], [2, 2], [2, 2], [1, 0], [3, 2], [1, 3], [1, 3]], [[1, 2], [1, 2], [2, 1], [2, 1], [3, 1], [1, 1], [2, 3], [2, 3], [3, 3], [2, 2], [2, 2], [2, 2], [2, 2], [3, 2], [3, 2], [1, 3]], [[2, 3], [3, 3], [2, 2], [3, 2]], [[3, 3]]]\n" ] } ], "source": [ "mps.ensure_left_canonical()\n", "\n", "def spin_to_fermion_mps(mps_spin, qn):\n", " nsites_spin = len(mps_spin)\n", " assert nsites_spin % 2 == 0\n", " nsites_fermion = nsites_spin // 2\n", " qnl_fermion = []\n", " mps_fermion = []\n", "\n", " for i in range(nsites_fermion):\n", " ml = mps_spin[2*i].shape[0]\n", " mr = mps_spin[2*i+1].shape[-1]\n", " merged_array = np.einsum(\"ipj, jql->ipql\", mps_spin[2*i], mps_spin[2*i+1])\n", " merged_array = np.transpose(merged_array, (0, 2, 1, 3))\n", " merged_array = merged_array.reshape(ml, 4, mr)\n", " mps_fermion.append(merged_array)\n", " qnl_fermion.append(qn[2*i])\n", " qnl_fermion.append(qn[-1])\n", " return mps_fermion, qnl_fermion\n", "\n", "qn = []\n", "for i_mps in range(len(mps)):\n", " iqn = [x.tolist() for x in mps.qn[i_mps]]\n", " qn.append(iqn)\n", "qn[0] = [[0, 0]]\n", "qn.append([[mf.mol.nelec[0], mf.mol.nelec[1]]])\n", "\n", "mps_fermion, qn_fermion = spin_to_fermion_mps(mps, qn)\n", "mps_tensors = [mps_fermion[i_mps] for i_mps in range(len(mps_fermion))]\n", "print(\"bond dimensions:\", [ts.shape[0] for ts in mps_tensors])\n", "print(\"quantum numbers:\", qn_fermion)" ] }, { "cell_type": "markdown", "metadata": { "id": "A9LKUpdABqQJ" }, "source": [ "## MPS overlap\n", "\n", "Here we show how to transform this external MPS to a ``block2`` MPS and compute the overlap between this MPS and the ground state MPS optimized in ``block2`` for the same system." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "id": "Lt_Z8MILCbzh" }, "outputs": [], "source": [ "from pyblock2.driver.core import DMRGDriver, SymmetryTypes\n", "from pyblock2.algebra.io import MPSTools\n", "\n", "bond_dims = [20] * 4 + [30] * 4\n", "noises = [1e-4] * 4 + [1e-5] * 4 + [0]\n", "thrds = [1e-10] * 8\n", "\n", "driver = DMRGDriver(scratch=\"./tmp\", symm_type=SymmetryTypes.SZ)\n", "driver.initialize_system(n_sites=ncas, n_elec=n_elec, spin=spin, orb_sym=orb_sym)\n", "\n", "pybra_sz = MPSTools.from_renormalizer_dense_sz(mps_tensors, qn_fermion)\n", "pybra_su2 = MPSTools.trans_sz_to_su2(pybra_sz, driver.basis, driver.target, target_twos=0)" ] }, { "cell_type": "markdown", "metadata": { "id": "FJyM9ze-CYq0" }, "source": [ "The overlap computation can be done in either the ``SZ`` or ``SU2`` mode in ``block2``.\n", "\n", "In the ``SZ`` mode:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "3W7vnmBVDGnn", "outputId": "8c463e4e-bea5-4d30-be17-0c96433933e5" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "integral symmetrize error = 0.0\n", "integral cutoff error = 0.0\n", "mpo terms = 2286\n", "\n", "Build MPO | Nsites = 6 | Nterms = 2286 | Algorithm = FastBIP | Cutoff = 1.00e-20\n", " Site = 0 / 6 .. Mmpo = 26 DW = 0.00e+00 NNZ = 26 SPT = 0.0000 Tmvc = 0.001 T = 0.009\n", " Site = 1 / 6 .. Mmpo = 66 DW = 0.00e+00 NNZ = 243 SPT = 0.8584 Tmvc = 0.001 T = 0.010\n", " Site = 2 / 6 .. Mmpo = 110 DW = 0.00e+00 NNZ = 459 SPT = 0.9368 Tmvc = 0.001 T = 0.011\n", " Site = 3 / 6 .. Mmpo = 66 DW = 0.00e+00 NNZ = 1147 SPT = 0.8420 Tmvc = 0.001 T = 0.012\n", " Site = 4 / 6 .. Mmpo = 26 DW = 0.00e+00 NNZ = 243 SPT = 0.8584 Tmvc = 0.000 T = 0.008\n", " Site = 5 / 6 .. Mmpo = 1 DW = 0.00e+00 NNZ = 26 SPT = 0.0000 Tmvc = 0.000 T = 0.006\n", "Ttotal = 0.057 Tmvc-total = 0.004 MPO bond dimension = 110 MaxDW = 0.00e+00\n", "NNZ = 2144 SIZE = 18004 SPT = 0.8809\n", "\n", "Rank = 0 Ttotal = 0.120 MPO method = FastBipartite bond dimension = 110 NNZ = 2144 SIZE = 18004 SPT = 0.8809\n", "\n", "Sweep = 0 | Direction = forward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.663 | E = -2.9576460853 | DW = 6.66567e-05\n", "\n", "Sweep = 1 | Direction = backward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.724 | E = -2.9576460853 | DE = 1.04e-11 | DW = 5.92236e-05\n", "\n", "Sweep = 2 | Direction = forward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.779 | E = -2.9576460853 | DE = -3.64e-11 | DW = 6.66526e-05\n", "\n", "Sweep = 3 | Direction = backward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.837 | E = -2.9576460853 | DE = -2.38e-12 | DW = 5.92239e-05\n", "\n", "Sweep = 4 | Direction = forward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.896 | E = -2.9576460853 | DE = 2.36e-12 | DW = 5.02379e-06\n", "\n", "Sweep = 5 | Direction = backward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.959 | E = -2.9576460853 | DE = -3.69e-11 | DW = 2.73774e-06\n", "\n", "Sweep = 6 | Direction = forward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 1.024 | E = -2.9576460853 | DE = -1.95e-11 | DW = 5.02112e-06\n", "\n", "Sweep = 7 | Direction = backward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 1.084 | E = -2.9576460853 | DE = 1.93e-11 | DW = 2.73775e-06\n", "\n", "Sweep = 8 | Direction = forward | Bond dimension = 30 | Noise = 0.00e+00 | Dav threshold = 1.00e-09\n", "Time elapsed = 1.144 | E = -2.9576460852 | DE = 1.66e-10 | DW = 1.91390e-06\n", "\n", "Overlap (SZ): 0.9999999659573686\n" ] } ], "source": [ "driver.symm_type = SymmetryTypes.SZ\n", "driver.initialize_system(n_sites=ncas, n_elec=n_elec, spin=spin, orb_sym=orb_sym)\n", "zbra = MPSTools.to_block2(pybra_sz, driver.basis, tag='ZBRA')\n", "mpo = driver.get_qc_mpo(h1e=h1e, g2e=h2e, ecore=ecore, iprint=1)\n", "\n", "zket = driver.get_random_mps(tag=\"ZKET\", bond_dim=20, nroots=1)\n", "energies = driver.dmrg(mpo, zket, n_sweeps=20, bond_dims=bond_dims, noises=noises, thrds=thrds, iprint=1)\n", "zket = driver.adjust_mps(zket, dot=1)[0]\n", "\n", "impo = driver.get_identity_mpo()\n", "print('Overlap (SZ):', driver.expectation(zbra, impo, zket))" ] }, { "cell_type": "markdown", "metadata": { "id": "7SRKmVeNDgdL" }, "source": [ "In the ``SU2`` mode:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "colab": { "base_uri": "https://localhost:8080/" }, "id": "g189S736DiNx", "outputId": "df485345-7011-4b1d-d4cd-589f24478ba6" }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "integral symmetrize error = 0.0\n", "integral cutoff error = 0.0\n", "mpo terms = 863\n", "\n", "Build MPO | Nsites = 6 | Nterms = 863 | Algorithm = FastBIP | Cutoff = 1.00e-20\n", " Site = 0 / 6 .. Mmpo = 13 DW = 0.00e+00 NNZ = 13 SPT = 0.0000 Tmvc = 0.000 T = 0.003\n", " Site = 1 / 6 .. Mmpo = 34 DW = 0.00e+00 NNZ = 97 SPT = 0.7805 Tmvc = 0.000 T = 0.003\n", " Site = 2 / 6 .. Mmpo = 56 DW = 0.00e+00 NNZ = 180 SPT = 0.9055 Tmvc = 0.001 T = 0.005\n", " Site = 3 / 6 .. Mmpo = 34 DW = 0.00e+00 NNZ = 405 SPT = 0.7873 Tmvc = 0.000 T = 0.008\n", " Site = 4 / 6 .. Mmpo = 14 DW = 0.00e+00 NNZ = 100 SPT = 0.7899 Tmvc = 0.000 T = 0.003\n", " Site = 5 / 6 .. Mmpo = 1 DW = 0.00e+00 NNZ = 14 SPT = 0.0000 Tmvc = 0.000 T = 0.002\n", "Ttotal = 0.024 Tmvc-total = 0.002 MPO bond dimension = 56 MaxDW = 0.00e+00\n", "NNZ = 809 SIZE = 4753 SPT = 0.8298\n", "\n", "Rank = 0 Ttotal = 0.038 MPO method = FastBipartite bond dimension = 56 NNZ = 809 SIZE = 4753 SPT = 0.8298\n", "\n", "Sweep = 0 | Direction = forward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.046 | E = -2.9576460853 | DW = 1.81931e-06\n", "\n", "Sweep = 1 | Direction = backward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.081 | E = -2.9576460854 | DE = -5.92e-11 | DW = 9.59991e-07\n", "\n", "Sweep = 2 | Direction = forward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.116 | E = -2.9576460854 | DE = -3.03e-12 | DW = 1.81905e-06\n", "\n", "Sweep = 3 | Direction = backward | Bond dimension = 20 | Noise = 1.00e-04 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.150 | E = -2.9576460854 | DE = 4.31e-12 | DW = 9.60002e-07\n", "\n", "Sweep = 4 | Direction = forward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.597 | E = -2.9576460854 | DE = -4.28e-12 | DW = 4.15116e-09\n", "\n", "Sweep = 5 | Direction = backward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.718 | E = -2.9576460854 | DE = -3.61e-12 | DW = 7.70782e-10\n", "\n", "Sweep = 6 | Direction = forward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 0.865 | E = -2.9576460854 | DE = -4.10e-12 | DW = 4.13558e-09\n", "\n", "Sweep = 7 | Direction = backward | Bond dimension = 30 | Noise = 1.00e-05 | Dav threshold = 1.00e-10\n", "Time elapsed = 1.129 | E = -2.9576460854 | DE = 1.87e-12 | DW = 7.70814e-10\n", "\n", "Sweep = 8 | Direction = forward | Bond dimension = 30 | Noise = 0.00e+00 | Dav threshold = 1.00e-09\n", "Time elapsed = 1.398 | E = -2.9576460851 | DE = 2.61e-10 | DW = 6.52891e-18\n", "\n", "Overlap (SU2): 0.9999990300129581\n" ] } ], "source": [ "driver.symm_type = SymmetryTypes.SU2\n", "driver.initialize_system(n_sites=ncas, n_elec=n_elec, spin=spin, orb_sym=orb_sym)\n", "bra = MPSTools.to_block2(pybra_su2, driver.basis, tag='BRA')\n", "mpo = driver.get_qc_mpo(h1e=h1e, g2e=h2e, ecore=ecore, iprint=1)\n", "\n", "ket = driver.get_random_mps(tag=\"KET\", bond_dim=20, nroots=1)\n", "energies = driver.dmrg(mpo, ket, n_sweeps=20, bond_dims=bond_dims, noises=noises, thrds=thrds, iprint=1)\n", "ket = driver.adjust_mps(ket, dot=1)[0]\n", "\n", "impo = driver.get_identity_mpo()\n", "print('Overlap (SU2):', driver.expectation(bra, impo, ket))" ] } ], "metadata": { "colab": { "provenance": [] }, "kernelspec": { "display_name": "Python 3", "name": "python3" }, "language_info": { "name": "python" } }, "nbformat": 4, "nbformat_minor": 0 }