تحديد بُعد الفضاء الجزئي

في هذا الشرح العملي، سنوضح تأثير بُعد الفضاء الجزئي في تقنية استرجاع التكوين المتسقة ذاتيًا.

مسبقًا، لا نعرف ما هو بُعد الفضاء الجزئي الصحيح للحصول على مستوى دقة مستهدف. ومع ذلك، نعلم أن زيادة بُعد الفضاء الجزئي تزيد من دقة الطريقة. لذلك يمكننا دراسة دقة التنبؤات كدالة في بُعد الفضاء الجزئي. حدّد الجزيء وخصائصه.

# 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

ولّد بعض سلاسل البتات العشوائية لتحاكي عينات 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
)

استدعِ SQD مع أحجام دفعات متزايدة.

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))

يُظهر هذا الرسم أن زيادة بُعد الفضاء الجزئي تؤدي إلى نتائج أكثر دقة.

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()

تحديد حد صريح باستخدام max_dim

في المسح أعلاه، كان بُعد الفضاء الجزئي يُتحكَّم به بشكل غير مباشر فقط. يحدد samples_per_batch عدد سلاسل البتات المسحوبة في كل دفعة، لكن بُعد الفضاء الجزئي الناتج يتحدد بعدد سلاسل السبين ألفا والسبين بيتا المتمايزة التي تحتويها سلاسل البتات تلك، كما أن مجموعة التكوينات المتاحة لأخذ العينات الجزئية تنمو من تكرار لاسترجاع التكوين إلى التالي. ولهذا السبب فإن أبعاد الفضاء الجزئي على المحور الأفقي أعلاه أعداد غير منتظمة لم نخترها نحن.

يحدد الوسيط max_dim في diagonalize_fermionic_hamiltonian البُعد مباشرة. فهو يحدّ من عدد سلاسل CI المحتفَظ بها في كل قطاع سبين، وبما أن الفضاء الجزئي يمتد بالجداء الديكارتي للقطاعين، فإن بُعده هو حاصل ضرب بُعدي القطاعين:

• يحدّ max_dim=d كلا القطاعين بـ d من السلاسل، وبالتالي يكون بُعد الفضاء الجزئي على الأكثر d ** 2.

• يحدّ max_dim=(d_a, d_b) قطاعي السبين ألفا والسبين بيتا كلًا على حدة، وبالتالي يكون بُعد الفضاء الجزئي على الأكثر d_a * d_b. وعندما يكون symmetrize_spin=True، يجب أن يتساوى الحدّان.

• max_dim=None، وهو الافتراضي، لا يفرض أي حد.

عندما تحتوي الدفعة على سلاسل أكثر مما يسمح به الحد، تُستبعد الأقل أهمية: تُحفظ أولًا السلاسل الممرَّرة عبر include_configurations، تليها السلاسل المنقولة من التكرار السابق (انظر carryover_threshold)، وتُملأ الخانات المتبقية بالسلاسل المأخوذة عينةً بترتيب تنازلي حسب عدد مرات أخذها.

ولأن max_dim يحدّ من حجم أكبر عملية تقطير قطري سيُطلب من حلّال الحالات الذاتية تنفيذها، فهو وسيلة مريحة لإبقاء الحلّال ضمن ميزانية الوقت والذاكرة: اختر samples_per_batch بسخاء، ودع max_dim يقرر الحجم الذي يُسمح للفضاء الجزئي بالنمو إليه.

أدناه، نكرر الدراسة، لكن هذه المرة نثبّت samples_per_batch ونزيد max_dim.

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]

أصبح بُعد الفضاء الجزئي الآن بالضبط max_dim ** 2 في كل تشغيل: فقد تم اختياره بدلًا من اكتشافه. ويُظهر رسم المسحين معًا أن الدقة تحكمها أبعاد الفضاء الجزئي نفسه، أيًّا كان الوسيط المستخدم للوصول إليه. ومع ذلك، ليس بالضرورة أن ينطبق المنحنيان تمامًا. فتشغيلات max_dim تسحب 600 سلسلة بتات لكل دفعة ثم تحتفظ بالسلاسل الأكثر أخذًا للعينات، لذلك يمكنها عند بُعد معين أن تحتفظ بمجموعة تكوينات أفضل من تشغيل وصل إلى ذلك البُعد بأخذ عدد أقل من سلاسل البتات أصلًا.

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()