궤도 최적화를 통한 SQD 추정치 개선

샘플 기반 양자 대각화(SQD)는 전자 배열의 고정된 부분공간에서 해밀토니안을 대각화함으로써 기저 상태 에너지를 근사화한다. 이러한 추정치는 해밀토니안이 표현되는 궤도 기저에 따라 달라지며, 궤도 최적화 (OO)는 이러한 자유도를 활용하여 부분공간을 확대하지 않고도 에너지를 낮추는 방식이다.

이 가이드에서는 $N_2$ 분자에 대해 SQD를 실행한 다음, 궤도 최적화를 통해 결과를 개선합니다. 이때 ffsim 를 사용하여 해밀토니안을 표현하고, 에너지를 최소화하는 궤도 회전 방향을 구합니다.

SQD 실행

우리는 분자 궤도(MO) 기저에서 $N_2$ 에 대한 분자 적분을 구축하고, 균일한 무작위 표본을 생성한 다음, SQD를 실행하여 기저 상태 근사치를 구합니다.

import numpy as np
import pyscf
import pyscf.cc
import pyscf.mcscf
from qiskit_addon_sqd.counts import generate_bit_array_uniform
from qiskit_addon_sqd.fermion import diagonalize_fermionic_hamiltonian

# Specify molecule properties
num_orbitals = 16
num_elec_a = num_elec_b = 5
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

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

# Run SQD
result = diagonalize_fermionic_hamiltonian(
    hcore,
    eri,
    bit_array,
    samples_per_batch=100,
    norb=num_orbitals,
    nelec=(num_elec_a, num_elec_b),
    num_batches=1,
    max_iterations=5,
    symmetrize_spin=True,
    seed=rng,
)

Output:

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

sqd_energy = result.energy + nuclear_repulsion_energy
print(f"Exact energy:  {exact_energy:.8f}")
print(f"SQD energy:    {sqd_energy:.8f}")

Output:

Exact energy:  -109.04667178
SQD energy:    -108.98469255

궤도를 최적화한다

궤도 최적화는 변분 에너지를 낮추는 궤도 회전 방식을 찾는 과정입니다

E = \langle \psi | \mathcal{U}^\dagger\, H\, \mathcal{U} | \psi \rangle

SQD 기저 상태 근사 $|\psi\rangle$ 의 경우, 궤도 회전은 $N \times N$ 의 유니터리 행렬 $\mathbf{U}$ (여기서 $N$ 는 공간 궤도의 개수)로 지정되며, 이는 연산자를 통해 다체 상태에 작용한다.

\mathcal{U} = \exp\left[\sum_{pq, \sigma} \log(\mathbf{U})_{pq}\, a^\dagger_{p\sigma} a_{q\sigma}\right].

ffsim.optimize_orbitals $\mathbf{U}$ 행렬 을 반환하며, 이를 궤도 기저에 (를 통해 hamiltonian.rotated) 적용하는 것은 $\mathcal{U}$ 를 상태에 적용하는 것과 동일하다. 자세한 내용은 ffsim의 궤도-회전 설명을 참고하십시오.

궤도를 회전시키면 부분공간에서 인식되는 해밀토니안이 변하므로, 에너지가 더 이상 개선되지 않을 때까지 다음 두 단계를 번갈아 가며 수행합니다:

고정된 구성 집합에 대해 현재 기저에서 해밀토니안을 대각화하라.
결과 상태의 에너지를 최소화하는 회전 각도를 구하여 궤도를 최적화한 다음, 적분식을 새로운 기저로 변환한다.

우리는 궤도 회전 단계를 다음으로 위임하며 ffsim.optimize_orbitals, 이 단계에서는 상태의 1체 및 2체 환원 밀도 행렬(RDM)을 바탕으로 에너지를 최소화하는 회전 각도를 구합니다. 해당 조항 참조 . 자세한 내용은 II A 4를 참조하십시오.

이 경우 궤도 최적화가 도움이 되는 이유

SCF 분자 궤도(MO) 기저 함수는 전체 CI 문제에 대해 궤도 회전에 대해 정적 입니다. 그러나 SQD는 작은 절단 부분공간(여기서는 약 1,900만 개의 전체 CI 결정자 중 수백 개의 CI 문자열)에서 작동하는데, 이 부분공간에 대해서는 MO 기저가 일반적으로 최적이 아니기 때문에, 궤도를 회전시키면 해당 부분공간이 표현할 수 있는 에너지가 낮아진다.

import ffsim
from pyscf import fci

# ffsim's ``MolecularHamiltonian`` uses the same "chemist" ordering for the two-body
# tensor as PySCF's ``eri``, and stores the nuclear repulsion energy as the constant
# term so that expectation values come out as total energies.

교대 대각화와 궤도 최적화

우리는 대각화 부분공간을 앞서 SQD를 통해 발견된 구성으로 고정 하여, 각 반복 과정에서 궤도 회전 효과만을 분리하여 분석할 수 있도록 합니다. 각 반복마다:

PySCF's 의 선택된 CI(Selected-CI) 솔버를 사용하여, 현재 기저에서 고정된 부분공간에 대해 해밀토니안을 대각화합니다.
결과 상태의 RDM을 생성하며 , 이것이ffsim.optimize_orbitals 필요한 전부입니다.
궤도를 최적화합니다 : 에너지를 최소화하는 회전 각도를 반환하며ffsim.optimize_orbitals, 이를 적분에 적용하여 개선된 기저로 전환합니다.

우리는 각 최적화 단계에 앞서 에너지 값을 기록합니다. 기저가 매 반복마다 개선되기 때문에, 이 수열은 고정된 부분공간에서 달성 가능한 최상의 에너지 값을 향해 단조 감소합니다.

# Fix the diagonalization subspace to the configurations found by SQD.
ci_strings = (result.sci_state.ci_strs_a, result.sci_state.ci_strs_b)
nelec = (num_elec_a, num_elec_b)

# Start from the MO basis in which we ran SQD.
hamiltonian_opt = ffsim.MolecularHamiltonian(
    hcore, eri, constant=nuclear_repulsion_energy
)

num_iters = 10
for i in range(num_iters):
    # Diagonalize over the fixed subspace in the current basis.
    myci = fci.selected_ci.SelectedCI()
    myci = fci.addons.fix_spin_(myci, ss=spin_sq)
    _, amplitudes = fci.selected_ci.kernel_fixed_space(
        myci,
        hamiltonian_opt.one_body_tensor,
        hamiltonian_opt.two_body_tensor,
        num_orbitals,
        nelec,
        ci_strs=ci_strings,
    )

    # Build the RDMs and record the energy before re-optimizing the orbitals.
    dm1, dm2 = myci.make_rdm12(amplitudes, num_orbitals, nelec)
    rdm = ffsim.ReducedDensityMatrix(dm1, dm2)
    energy = rdm.expectation(hamiltonian_opt).real
    print(f"Iteration {i}: energy = {energy:.8f}")

    # Rotate the Hamiltonian into the energy-minimizing basis for the next iteration.
    # optimize_orbitals returns the unitary matrix U minimizing
    # rdm.rotated(U).expectation(hamiltonian), equivalently
    # rdm.expectation(hamiltonian.rotated(U.conj().T)), so we rotate by U^dagger.
    orbital_rotation = ffsim.optimize_orbitals(rdm, hamiltonian_opt)
    hamiltonian_opt = hamiltonian_opt.rotated(orbital_rotation.T.conj())

Output:

Iteration 0: energy = -108.98452447
Iteration 1: energy = -108.99981993
Iteration 2: energy = -109.00585329
Iteration 3: energy = -109.00816569
Iteration 4: energy = -109.00936616
Iteration 5: energy = -109.01014322
Iteration 6: energy = -109.01069439
Iteration 7: energy = -109.01109308
Iteration 8: energy = -109.01138928
Iteration 9: energy = -109.01161411

결과를 비교해 보세요

궤도 최적화는 고정 부분공간 추정치를 개선하여, 정확한 에너지 값과의 격차를 상당 부분 좁히면서도 그 값보다 높은 수준을 유지합니다.

# Diagonalize once more in the final optimized basis to report the improved energy.
myci = fci.selected_ci.SelectedCI()
myci = fci.addons.fix_spin_(myci, ss=spin_sq)
_, amplitudes = fci.selected_ci.kernel_fixed_space(
    myci,
    hamiltonian_opt.one_body_tensor,
    hamiltonian_opt.two_body_tensor,
    num_orbitals,
    nelec,
    ci_strs=ci_strings,
)
dm1, dm2 = myci.make_rdm12(amplitudes, num_orbitals, nelec)
energy_after_oo = (
    ffsim.ReducedDensityMatrix(dm1, dm2).expectation(hamiltonian_opt).real
)

print(f"Exact energy:      {exact_energy:.8f}")
print(f"SQD energy (MO):   {sqd_energy:.8f}")
print(f"Energy after OO:   {energy_after_oo:.8f}")

Output:

Exact energy:      -109.04667178
SQD energy (MO):   -108.98469255
Energy after OO:   -109.01178727

이 페이지가 도움이 되었습니까?

GitHub에서 버그, 오타를 보고하거나 컨텐츠를 요청하십시오.