Menghadkan dimensi subruang

Dalam panduan ini, kami akan menunjukkan kesan dimensi subruang dalam teknik pemulihan konfigurasi swakonsisten.

A priori, kami tidak tahu apakah dimensi subruang yang betul untuk mencapai tahap ketepatan sasaran. Walau bagaimanapun, kami tahu bahawa meningkatkan dimensi subruang meningkatkan ketepatan kaedah ini. Oleh itu, kami boleh mengkaji ketepatan ramalan sebagai fungsi dimensi subruang. Tentukan molekul dan sifat-sifatnya.

# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "pyscf": "pyscf", "qiskit_addon_sqd": "qiskit-addon-sqd"}
_missing = [pip for module, pip in _needed.items()
            if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
    %pip install -q {_pip}
if not _missing:
    print("\u2713 All required packages are installed")

import warnings

import pyscf
import pyscf.cc
import pyscf.mcscf

warnings.filterwarnings("ignore")

# Specify molecule properties
open_shell = False
spin_sq = 0

# Build N2 molecule
mol = pyscf.gto.Mole()
mol.build(
    atom=[["N", (0, 0, 0)], ["N", (1.0, 0, 0)]],
    basis="6-31g",
    symmetry="Dooh",
)

# Define active space
n_frozen = 2
active_space = range(n_frozen, mol.nao_nr())

# Get molecular integrals
scf = pyscf.scf.RHF(mol).run()
num_orbitals = len(active_space)
n_electrons = int(sum(scf.mo_occ[active_space]))
num_elec_a = (n_electrons + mol.spin) // 2
num_elec_b = (n_electrons - mol.spin) // 2
cas = pyscf.mcscf.CASCI(scf, num_orbitals, (num_elec_a, num_elec_b))
mo = cas.sort_mo(active_space, base=0)
hcore, nuclear_repulsion_energy = cas.get_h1cas(mo)
eri = pyscf.ao2mo.restore(1, cas.get_h2cas(mo), num_orbitals)

# Compute exact energy
exact_energy = cas.run().e_tot

converged SCF energy = -108.835236570775
CASCI E = -109.046671778080  E(CI) = -32.8155692383188  S^2 = 0.0000000

Jana beberapa rentetan bit rawak sebagai proksi sampel QPU.

import numpy as np
from qiskit_addon_sqd.counts import generate_bit_array_uniform

# Create a seed to control randomness throughout this workflow
rng = np.random.default_rng(24)

# Generate random samples
bit_array = generate_bit_array_uniform(
    10_000, num_orbitals * 2, rand_seed=rng
)

Panggil SQD dengan saiz kelompok yang semakin meningkat.

from qiskit_addon_sqd.fermion import diagonalize_fermionic_hamiltonian

list_samples_per_batch = [50, 200, 400, 600]

# SQD options
max_iterations = 5

# Eigenstate solver options
num_batches = 10
max_davidson_cycles = 200

energies = []
subspace_dimensions = []

for samples_per_batch in list_samples_per_batch:
    result = diagonalize_fermionic_hamiltonian(
        hcore,
        eri,
        bit_array,
        samples_per_batch=samples_per_batch,
        norb=num_orbitals,
        nelec=(num_elec_a, num_elec_b),
        num_batches=num_batches,
        max_iterations=max_iterations,
        symmetrize_spin=True,
        seed=rng,
    )
    energies.append(result.energy)
    subspace_dimensions.append(np.prod(result.sci_state.amplitudes.shape))

Plot ini menunjukkan bahawa peningkatan dimensi subruang menghasilkan keputusan yang lebih tepat.

import matplotlib.pyplot as plt

# Data for energies plot
x1 = subspace_dimensions
y1 = np.array(energies) + nuclear_repulsion_energy

fig, axs = plt.subplots(1, 1, figsize=(12, 6))

# Plot energies
axs.plot(x1, y1, marker=".", markersize=20, label="Estimated")
axs.set_xticks(x1)
axs.set_xticklabels(x1)
axs.axhline(y=exact_energy, color="red", linestyle="--", label="Exact")
axs.set_title("Approximated Ground State Energy vs subspace dimension")
axs.set_xlabel("Subspace dimension")
axs.set_ylabel("Energy (Ha)")
axs.legend()

plt.tight_layout()
plt.show()

Menetapkan had yang jelas dengan max_dim

Dalam sapuan di atas, dimensi subruang dikawal secara tidak langsung sahaja. samples_per_batch menetapkan berapa banyak rentetan bit diambil ke dalam setiap kelompok, tetapi dimensi subruang yang terhasil ditentukan oleh bilangan rentetan spin-alfa dan spin-beta yang berbeza yang terkandung dalam rentetan bit tersebut, dan kumpulan konfigurasi yang tersedia untuk subsampel bertambah daripada satu lelaran pemulihan konfigurasi kepada yang seterusnya. Itulah sebabnya dimensi subruang pada paksi mengufuk di atas ialah nombor yang tidak sekata yang tidak kami pilih.

Argumen max_dim bagi diagonalize_fermionic_hamiltonian menghadkan dimensi secara langsung. Ia mengehadkan bilangan rentetan CI yang dikekalkan dalam setiap sektor spin, dan memandangkan subruang direntangkan oleh hasil darab Cartesian bagi kedua-dua sektor, dimensinya ialah hasil darab kedua-dua dimensi sektor:

• max_dim=d mengehadkan kedua-dua sektor kepada d rentetan, jadi dimensi subruang paling banyak d ** 2.

• max_dim=(d_a, d_b) mengehadkan sektor spin-alfa dan spin-beta secara berasingan, jadi dimensi subruang paling banyak d_a * d_b. Apabila symmetrize_spin=True, kedua-dua had mestilah sama.

• max_dim=None, iaitu lalai, tidak mengenakan sebarang had.

Apabila sesuatu kelompok mengandungi lebih banyak rentetan daripada yang dibenarkan oleh had, rentetan yang paling kurang penting dibuang: rentetan yang diberikan melalui include_configurations dikekalkan dahulu, diikuti rentetan yang dibawa dari lelaran sebelumnya (lihat carryover_threshold), dan slot yang selebihnya diisi dengan rentetan yang disampel mengikut tertib menurun kekerapan ia disampel.

Oleh sebab max_dim menghadkan saiz pepenjuruan terbesar yang akan diminta daripada penyelesai keadaan eigen, ia merupakan cara yang mudah untuk mengekalkan penyelesai dalam bajet masa jalan dan memori: pilih samples_per_batch dengan murah hati, dan biarkan max_dim menentukan sebesar mana subruang dibenarkan membesar.

Di bawah, kami mengulangi kajian ini, kali ini dengan samples_per_batch ditetapkan dan max_dim ditingkatkan.

list_max_dim = [100, 200, 400, 800]

# Use a generous number of samples per batch, so that max_dim is what limits
# the dimension of the subspace
samples_per_batch = 600

energies_max_dim = []
subspace_dimensions_max_dim = []

for max_dim in list_max_dim:
    result = diagonalize_fermionic_hamiltonian(
        hcore,
        eri,
        bit_array,
        samples_per_batch=samples_per_batch,
        norb=num_orbitals,
        nelec=(num_elec_a, num_elec_b),
        num_batches=num_batches,
        max_iterations=max_iterations,
        symmetrize_spin=True,
        max_dim=max_dim,
        seed=rng,
    )
    energies_max_dim.append(result.energy)
    subspace_dimensions_max_dim.append(
        int(np.prod(result.sci_state.amplitudes.shape))
    )

print(f"max_dim:            {list_max_dim}")
print(f"subspace dimension: {subspace_dimensions_max_dim}")

max_dim:            [100, 200, 400, 800]
subspace dimension: [10000, 40000, 160000, 640000]

Dimensi subruang kini tepat max_dim ** 2 dalam setiap larian: ia dipilih dan bukannya ditemui. Memplot kedua-dua sapuan bersama menunjukkan bahawa ketepatan ditentukan oleh dimensi subruang itu sendiri, apa pun argumen yang digunakan untuk mencapainya. Walau bagaimanapun, kedua-dua lengkung itu tidak semestinya bertindih sepenuhnya. Larian max_dim mengambil 600 rentetan bit bagi setiap kelompok dan kemudian mengekalkan rentetan yang paling kerap disampel, jadi pada dimensi tertentu ia boleh mengekalkan set konfigurasi yang lebih baik berbanding larian yang mencapai dimensi itu dengan mensampel lebih sedikit rentetan bit sejak awal.

fig, axs = plt.subplots(1, 1, figsize=(12, 6))

axs.plot(
    subspace_dimensions,
    np.array(energies) + nuclear_repulsion_energy,
    marker=".",
    markersize=20,
    label="Varying samples_per_batch",
)
axs.plot(
    subspace_dimensions_max_dim,
    np.array(energies_max_dim) + nuclear_repulsion_energy,
    marker="^",
    markersize=12,
    label="Varying max_dim",
)
axs.axhline(y=exact_energy, color="red", linestyle="--", label="Exact")
axs.set_xscale("log")
axs.ticklabel_format(axis="y", useOffset=False)
axs.set_title("Approximated Ground State Energy vs subspace dimension")
axs.set_xlabel("Subspace dimension")
axs.set_ylabel("Energy (Ha)")
axs.legend()

plt.tight_layout()
plt.show()