ใช้ trim SQD เพื่อเลือกคอนฟิกูเรชันได้ดีขึ้น

พฤติกรรมเริ่มต้นของ diagonalize_fermionic_hamiltonian คือทำ diagonalization ในแต่ละแบตช์ย่อยที่สุ่มมาและเก็บผลลัพธ์ที่ดีที่สุด ส่วน Trim SQD จะมองแบตช์เป็นตัวเลือกที่ต้องคัดกรอง: แต่ละแบตช์ถูกทำ diagonalization แล้ว ตัดแต่ง เหลือเฉพาะคอนฟิกูเรชันที่มีน้ำหนักมากที่สุด จากนั้นนำส่วนที่เหลือของทุกแบตช์มารวมกันแล้วทำ diagonalization ร่วมกัน diagonalization แบบรวมนี้คือตัวที่พลังงานของมันถูกรายงานในรอบนั้น และผลที่ตัดแต่งแล้วของมันจะเป็นจุดเริ่มต้นของรอบถัดไป

การทำงานของ trim SQD อิงจากงานวิจัยร่วมของ Cleveland Clinic Foundation, RIKEN และ IBM® หากคุณใช้ฟีเจอร์นี้ใน งานวิจัยของคุณ โปรดอ้างอิงบทความดังกล่าวด้วย คุณสามารถหาการอ้างอิงได้ในไฟล์ CITATION.bib โดยมีคีย์การอ้างอิงคือ merz2026crossing

สำหรับรายละเอียดเพิ่มเติมเกี่ยวกับโค้ด SQD ที่ใช้ในตัวอย่างนี้ ดูบทเรียนแฮมิลโทเนียนทางเคมี

ตั้งค่าโมเลกุล

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

_needed = {"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 numpy as np
import pyscf
import pyscf.mcscf
from qiskit_addon_sqd.counts import generate_bit_array_uniform
from qiskit_addon_sqd.fermion import diagonalize_fermionic_hamiltonian
from qiskit_addon_sqd.trim import TrimPolicy

# 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

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

WARN: Unable to to identify input symmetry using original axes.
Different symmetry axes will be used.

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

รันลูป

trim SQD หนึ่งรอบทำสิ่งต่อไปนี้:

1. รอบคัดกรอง แบ่งบิตสตริงเป็นแบตช์ ทำ diagonalization แต่ละแบตช์ และเก็บสัดส่วน trim_ratio ของสตริง CI ในแต่ละภาคสปิน โดยจัดอันดับตามน้ำหนัก

2. รอบรวม ทำ diagonalization ส่วนที่เหลือที่รวมกันแล้วหนึ่งครั้ง นี่คือผลลัพธ์ที่รอบนั้นรายงาน และสตริงที่ตัดแต่งแล้วของมันจะเป็นจุดเริ่มต้นของรอบถัดไป

sci_solver ของลูปทำทั้งสองรอบ และค่าเริ่มต้นคือตัวแก้ PySCF ของแพ็กเกจนี้

เปิดใช้ trim SQD โดยส่ง TrimPolicy ดังนี้:

result = diagonalize_fermionic_hamiltonian(
    hcore,
    eri,
    bit_array,
    samples_per_batch=300,
    norb=num_orbitals,
    nelec=(num_elec_a, num_elec_b),
    num_batches=5,
    max_iterations=5,
    policy=TrimPolicy(trim_ratio=0.5, max_strings_per_trim=250),
    symmetrize_spin=True,
    max_dim=250,
    seed=rng,
)

sci_state = result.sci_state
print(
    f"Merged subspace: {len(sci_state.ci_strs_a)} x {len(sci_state.ci_strs_b)} CI strings"
)
print(f"Exact energy:     {exact_energy}")
print(f"Estimated energy: {result.energy + nuclear_repulsion_energy}")

Merged subspace: 148 x 148 CI strings
Exact energy:     -109.04667177808024
Estimated energy: -109.034369147063

ตรวจสอบว่าการคัดกรองช่วยหรือไม่

การคัดกรองให้ประโยชน์มากที่สุดเมื่อคอนฟิกูเรชันที่สุ่มมามีจำนวนมากกว่าที่ทำ diagonalization ได้ในคราวเดียวมาก ทำให้การจัดอันดับพูลที่ใหญ่กว่าพบคอนฟิกูเรชันที่แบตช์เดียวอาจพลาดไปได้จริง และให้ประโยชน์น้อยที่สุดเมื่อปริภูมิย่อยขนาดที่ทำได้ครอบคลุมฟังก์ชันคลื่นส่วนใหญ่อยู่แล้ว

การคัดกรองขึ้นอยู่กับว่าแบตช์ต่างกันออกไป TrimPolicy จัดการเรื่องนี้โดยปรับพฤติกรรมของ diagonalize_fermionic_hamiltonian ให้สุ่มพูลเดียวของบิตสตริง samples_per_batch * num_batches ตัว แล้วแบ่งให้แต่ละแบตช์ ดังนั้นจะไม่มีบิตสตริงใดปรากฏมากกว่าหนึ่งแบตช์

จำกัดขอบเขตปริภูมิย่อย

ปริภูมิย่อยในที่นี้ถูกขึงด้วย ผลคูณคาร์ทีเซียน ของอาร์เรย์สตริง CI สปินแอลฟาและสปินเบตา

การตัดแต่งจัดอันดับคอนฟิกูเรชัน แต่สิ่งที่ส่งต่อไปคือรายการสตริงต่อสปินสองรายการ การเก็บคอนฟิกูเรชันน้ำหนักสูงสุด 150 ตัวของปริภูมิย่อยมักเกี่ยวข้องกับสตริงแอลฟาที่ต่างกันราว 150 ตัวและสตริงเบตาที่ต่างกันราว 150 ตัว ซึ่งผลคูณของมันครอบคลุมคอนฟิกูเรชันประมาณ 22,500 ตัว ไม่ใช่ 150 ปริภูมิย่อยที่ส่งให้รอบถัดไปจึงมักใหญ่กว่าชุดที่จัดอันดับไว้มาก และโตแบบกำลังสองตามจำนวนที่คุณเก็บไว้

การตัดแต่งและการส่งต่อถูกจำกัดแยกกัน เพราะเกิดขึ้นคนละจุด: การตัดแต่งลดขนาดแต่ละแบตช์ของรอบคัดกรอง ส่วนการส่งต่อคือสิ่งที่ข้ามไปยังรอบถัดไป

มีอาร์กิวเมนต์สี่ตัวที่จำกัดมัน และควรตั้งอย่างน้อยหนึ่งตัว:

• trim_ratio คือสัดส่วนของสตริงในแต่ละแบตช์ ต่อภาคสปิน ที่ผ่านการคัดกรอง ค่าที่ต่ำลงทำให้ปริภูมิย่อยรวมเล็กลง แต่ทิ้งคอนฟิกูเรชันที่รอบถัดไปต้องค้นพบใหม่ด้วยการสุ่มตัวอย่าง

• max_strings_per_trim จำกัดจำนวนสตริงที่แต่ละแบตช์ให้ ต่อภาค โดยไม่คำนึงถึงสัดส่วน มีประโยชน์เมื่อใช้กับ symmetrize_spin ซึ่งรวมทั้งสองภาคสปินเป็นรายการเดียวจึงทำให้สิ่งที่แต่ละแบตช์ให้เพิ่มเป็นประมาณสองเท่า ขอบเขตนี้เป็นต่อแบตช์: ปริภูมิย่อยรวมอาจมีสตริงได้มากถึง num_batches เท่า หักส่วนที่แบตช์ต่าง ๆ มีร่วมกัน ซึ่งในทางปฏิบัติมีไม่มาก

• carryover_threshold, carryover_ratio และ max_carryover จำกัดสตริงที่เป็นจุดเริ่มต้นของรอบถัดไป โดยใช้กับปริภูมิย่อยรวมทั้งหมด ไม่ใช่ต่อแบตช์ carryover_threshold คือเกณฑ์แอมพลิจูดที่ StandardPolicy ใช้ โดยมีค่าเริ่มต้นเดียวกัน อีกสองตัวจำกัดจำนวน ซึ่งเกณฑ์เพียงอย่างเดียวทำไม่ได้

• max_dim บน diagonalize_fermionic_hamiltonian เอง จำกัดแต่ละภาคสปินของทุกปริภูมิย่อยที่ลูปทำ diagonalization รวมถึงปริภูมิย่อยรวมด้วย ทำหน้าที่เป็นตัวสำรอง และใช้กับสิ่งที่ policy คืนมา หลังจากรวมการส่งต่อเข้ากับตัวอย่างใหม่แล้ว

การรันข้างต้นตั้ง max_strings_per_trim=250 และ max_dim=250 ร่วมกัน ดังนั้นทั้งสิ่งที่แต่ละแบตช์ให้และปริภูมิย่อยที่สร้างจากมันจะไม่บานปลาย

ใช้ตัวแก้ค่าเฉพาะภายนอก

คุณสามารถใช้ trim SQD กับตัวแก้ค่าเฉพาะภายนอกได้โดยส่งอาร์กิวเมนต์ sci_solver:

from sbd.sbd_solver import create_sbd_solver

result = diagonalize_fermionic_hamiltonian(
    hcore,
    eri,
    bit_array,
    samples_per_batch=300,
    norb=num_orbitals,
    nelec=(num_elec_a, num_elec_b),
    num_batches=5,
    policy=TrimPolicy(trim_ratio=0.1),
    sci_solver=create_sbd_solver(sbd_config={"method": 0, "eps": 1e-10}),
    seed=rng,
)

สนิปเพ็ตก่อนหน้าต้องติดตั้งแพ็กเกจ sbd-eigensolver จาก PyPI

กระจายงาน

เมื่อใช้ TrimPolicy ลูปจะเรียก sci_solver สองครั้งต่อรอบ: ครั้งหนึ่งกับทุกแบตช์ และอีกครั้งกับปริภูมิย่อยรวม ทั้งสองการเรียกผ่านอินเทอร์เฟซ sci_solver ปกติ ซึ่งการรองรับหลายโพรเซสของแพ็กเกจนี้จะเรียกร่วมกันบนทุก MPI rank

มีการประมวลผลแบบขนานที่เป็นไปได้สองระดับ และเป็นอิสระต่อกัน:

• ภายใน diagonalization หนึ่งครั้ง ตัวแก้เช่น SBD สามารถกระจาย diagonalization หนึ่งครั้งไปยังหลาย rank

• ข้ามแบตช์ ตัวแก้อาจทำ diagonalization หลายแบตช์พร้อมกันได้ ในหลักการ

เมื่อได้เซตของแบตช์ ตัวแก้จะเป็นผู้กำหนดวิธีมอบหมายงานให้โพรเซสทำงานของมัน