From 5d4d8ba58b3ea248084aed6a744b2b32523909a9 Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 18:53:10 -0400 Subject: [PATCH 01/41] add excitation function --- pennylane/qchem/convert.py | 33 +++++++++++++++++++++++++++++++++ 1 file changed, 33 insertions(+) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index e75e804b4d0..48fac426f18 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -354,3 +354,36 @@ def import_operator(qubit_observable, format="openfermion", wires=None, tol=1e01 return qml.dot(*_openfermion_to_pennylane(qubit_observable, wires=wires)) return qml.Hamiltonian(*_openfermion_to_pennylane(qubit_observable, wires=wires)) + + +def _excitations(electrons, orbitals): + r"""Generate all possible single and double excitations from a Hartree-Fock reference state. + + Args: + electrons (int): Number of electrons + orbitals (int) – Number of spin orbitals + + Returns: + tuple(list, list): lists with the indices of the spin orbitals involved in the excitations + + **Example** + + >>> electrons = 2 + >>> orbitals = 4 + >>> _excitations(electrons, orbitals) + ([[[0, 2]], [[0, 3]], [[1, 2]], [[1, 3]]], [[0, 1, 2, 3]]) + """ + singles_p, singles_q = [], [] + doubles_pq, doubles_rs = [], [] + + for i in range(electrons): + singles_p += [i] + doubles_pq += [[k, i] for k in range(i)] + for j in range(electrons, orbitals): + singles_q += [j] + doubles_rs += [[k, j] for k in range(electrons, j)] + + singles = [[[p] + [q]] for p in singles_p for q in singles_q] + doubles = [pq + rs for pq in doubles_pq for rs in doubles_rs] + + return singles, doubles From d75a41530ddc007041317e7ef9e821b72acbb41f Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 19:16:05 -0400 Subject: [PATCH 02/41] add excited states --- pennylane/qchem/convert.py | 54 ++++++++++++++++++++++++++++++++++++-- 1 file changed, 52 insertions(+), 2 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 48fac426f18..a0a4de3be7e 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -360,8 +360,8 @@ def _excitations(electrons, orbitals): r"""Generate all possible single and double excitations from a Hartree-Fock reference state. Args: - electrons (int): Number of electrons - orbitals (int) – Number of spin orbitals + electrons (int): number of electrons + orbitals (int): number of spin orbitals Returns: tuple(list, list): lists with the indices of the spin orbitals involved in the excitations @@ -387,3 +387,53 @@ def _excitations(electrons, orbitals): doubles = [pq + rs for pq in doubles_pq for rs in doubles_rs] return singles, doubles + + +def _excitated_states(electrons, orbitals, excitation): + r"""Generate excited states from a Hartree-Fock reference state. + + Args: + electrons (int): number of electrons + orbitals (int): number of spin orbitals + excitation (int): number of excited electrons + + Returns: + tuple(list, list): lists of excited states and signs obtained by reordering of the creation operators + + **Example** + + >>> electrons = 2 + >>> orbitals = 5 + >>> excitation = 2 + >>> _excitated_states(electrons, orbitals, excitation) + ([28, 26, 25], [ 1, -1, 1]) + """ + hf_state = [1] * electrons + [0] * (orbitals - electrons) + + singles, doubles = _excitations(electrons, orbitals) + + states, signs = [], [] + + if excitation == 1: + for s in singles: + state = hf_state.copy() + state[s[0]], state[s[1]] = state[s[1]], state[s[0]] + states += [state] + order = len(range(s[0], electrons - 1)) + signs.append((-1) ** order) + + if excitation == 2: + for d in doubles: + state = hf_state.copy() + state[d[0]], state[d[2]] = state[d[2]], state[d[0]] + state[d[1]], state[d[3]] = state[d[3]], state[d[1]] + states += [state] + + order_pq = len(range(d[0], electrons - 1)) + order_rs = len(range(d[1], electrons - 1)) + signs.append((-1) ** (order_pq + order_rs + 1)) + + states_str = ["".join([str(i) for i in item]) for item in states] + state_int = [int(bb[::-1], 2) for bb in states_str] + + return np.array(state_int), np.array(signs) From a33369a53d37acf2153ecd03e94cb3fe158e7761 Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 19:25:03 -0400 Subject: [PATCH 03/41] modify var names --- pennylane/qchem/convert.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index a0a4de3be7e..ddab25041cb 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -433,7 +433,7 @@ def _excitated_states(electrons, orbitals, excitation): order_rs = len(range(d[1], electrons - 1)) signs.append((-1) ** (order_pq + order_rs + 1)) - states_str = ["".join([str(i) for i in item]) for item in states] - state_int = [int(bb[::-1], 2) for bb in states_str] + states_str = ["".join([str(i) for i in state]) for state in states] + state_int = [int(state[::-1], 2) for state in states_str] return np.array(state_int), np.array(signs) From 3d2aa331c528a5b3dd1b7898f0a9dad627c5114e Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 19:25:21 -0400 Subject: [PATCH 04/41] modify var names --- pennylane/qchem/convert.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index ddab25041cb..62a40bc2e03 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -434,6 +434,6 @@ def _excitated_states(electrons, orbitals, excitation): signs.append((-1) ** (order_pq + order_rs + 1)) states_str = ["".join([str(i) for i in state]) for state in states] - state_int = [int(state[::-1], 2) for state in states_str] + states_int = [int(state[::-1], 2) for state in states_str] - return np.array(state_int), np.array(signs) + return np.array(states_int), np.array(signs) From d41071bc7598eb94f883676a00129a879d9fada9 Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 22:36:43 -0400 Subject: [PATCH 05/41] add tests --- pennylane/qchem/convert.py | 2 +- tests/qchem/of_tests/test_convert.py | 48 ++++++++++++++++++++++++++-- 2 files changed, 47 insertions(+), 3 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 62a40bc2e03..1753f5a6188 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -383,7 +383,7 @@ def _excitations(electrons, orbitals): singles_q += [j] doubles_rs += [[k, j] for k in range(electrons, j)] - singles = [[[p] + [q]] for p in singles_p for q in singles_q] + singles = [[p] + [q] for p in singles_p for q in singles_q] doubles = [pq + rs for pq in doubles_pq for rs in doubles_rs] return singles, doubles diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 6c2f0f008ab..6b3ce829065 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -12,8 +12,8 @@ # See the License for the specific language governing permissions and # limitations under the License. """ -Unit tests for functions needed for converting observables obtained from external libraries to a -PennyLane observable. +Unit tests for functions needed for converting objects obtained from external libraries to a +PennyLane object. """ # pylint: disable=too-many-arguments,protected-access import os @@ -813,3 +813,47 @@ def dummy_ansatz(phis, wires): res = dummy_cost(phis) assert np.abs(res - expected_cost) < tol["atol"] + + +@pytest.mark.parametrize( + ("electrons", "orbitals", "singles_ref", "doubles_ref"), + [ + # trivial case + (2, 4, [[0, 2], [0, 3], [1, 2], [1, 3]], [[0, 1, 2, 3]]), + ], +) +def test_excitations(electrons, orbitals, singles_ref, doubles_ref): + r"""Test if the _excitations function returns correct single and double excitations.""" + singles, doubles = qchem.convert._excitations(electrons, orbitals) + assert singles == singles_ref + assert doubles == doubles_ref + + +@pytest.mark.parametrize( + ("electrons", "orbitals", "excitation", "states_ref", "signs_ref"), + [ + # reference data computed with pyscf: + # pyscf_addrs, pyscf_signs = tn_addrs_signs(orbitals, electrons, excitation) + # pyscf_state = pyscf.fci.cistring.addrs2str(orbitals, electrons, pyscf_addrs) + # pyscf_state, pyscf_signs + ( + 3, + 8, + 1, + np.array([14, 22, 38, 70, 134, 13, 21, 37, 69, 133, 11, 19, 35, 67, 131]), + np.array([1, 1, 1, 1, 1, -1, -1, -1, -1, -1, 1, 1, 1, 1, 1]), + ), + ( + 3, + 6, + 2, + np.array([28, 44, 52, 26, 42, 50, 25, 41, 49]), + np.array([1, 1, 1, -1, -1, -1, 1, 1, 1]), + ), + ], +) +def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref): + r"""Test if the _excitated_states function returns correct states and signs.""" + states, signs = qchem.convert._excitated_states(electrons, orbitals, excitation) + assert np.allclose(states, states_ref) + assert np.allclose(signs, signs_ref) From 0a8c78d5e93b1579ea1c896f828657e8ea833b30 Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 22:54:52 -0400 Subject: [PATCH 06/41] make code compact --- pennylane/qchem/convert.py | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 1753f5a6188..5c464b2f1cd 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -408,7 +408,7 @@ def _excitated_states(electrons, orbitals, excitation): >>> _excitated_states(electrons, orbitals, excitation) ([28, 26, 25], [ 1, -1, 1]) """ - hf_state = [1] * electrons + [0] * (orbitals - electrons) + hf_state = qml.qchem.hf_state(electrons, orbitals) singles, doubles = _excitations(electrons, orbitals) @@ -417,18 +417,15 @@ def _excitated_states(electrons, orbitals, excitation): if excitation == 1: for s in singles: state = hf_state.copy() - state[s[0]], state[s[1]] = state[s[1]], state[s[0]] + state[s] = state[s[::-1]] states += [state] - order = len(range(s[0], electrons - 1)) - signs.append((-1) ** order) + signs.append((-1) ** len(range(s[0], electrons - 1))) if excitation == 2: for d in doubles: state = hf_state.copy() - state[d[0]], state[d[2]] = state[d[2]], state[d[0]] - state[d[1]], state[d[3]] = state[d[3]], state[d[1]] + state[d] = state[[d[2], d[3], d[0], d[1]]] states += [state] - order_pq = len(range(d[0], electrons - 1)) order_rs = len(range(d[1], electrons - 1)) signs.append((-1) ** (order_pq + order_rs + 1)) @@ -436,4 +433,4 @@ def _excitated_states(electrons, orbitals, excitation): states_str = ["".join([str(i) for i in state]) for state in states] states_int = [int(state[::-1], 2) for state in states_str] - return np.array(states_int), np.array(signs) + return states_int, signs From 52d8ed4cc8b7b712a7902de18371f580f1ac7f5a Mon Sep 17 00:00:00 2001 From: soranjh Date: Mon, 31 Jul 2023 22:57:49 -0400 Subject: [PATCH 07/41] make code compact --- pennylane/qchem/convert.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 5c464b2f1cd..dbc384a624e 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -398,7 +398,7 @@ def _excitated_states(electrons, orbitals, excitation): excitation (int): number of excited electrons Returns: - tuple(list, list): lists of excited states and signs obtained by reordering of the creation operators + tuple(array, array): arrays of excited states and signs obtained by reordering of the creation operators **Example** @@ -406,12 +406,10 @@ def _excitated_states(electrons, orbitals, excitation): >>> orbitals = 5 >>> excitation = 2 >>> _excitated_states(electrons, orbitals, excitation) - ([28, 26, 25], [ 1, -1, 1]) + (array([28, 26, 25]), array([ 1, -1, 1])) """ hf_state = qml.qchem.hf_state(electrons, orbitals) - singles, doubles = _excitations(electrons, orbitals) - states, signs = [], [] if excitation == 1: From e06c7a535cd51cd6f2e8ace74e79154b183aa9b8 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 14:20:49 -0400 Subject: [PATCH 08/41] add ucisd wavefunction constructor --- pennylane/qchem/convert.py | 208 +++++++++++++++++++++++++++++++++++++ 1 file changed, 208 insertions(+) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index dbc384a624e..55acad25d1b 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -432,3 +432,211 @@ def _excitated_states(electrons, orbitals, excitation): states_int = [int(state[::-1], 2) for state in states_str] return states_int, signs + + +def _ucisd_state(cisd_solver, state=0, tol=1e-15): + r"""Construct a wavefunction in the wf_dict format `{ (int_a, int_b): coeff}` from + PySCF's `UCISD` Solver object. + + The format is described in detail in the docstring of :func:`~.cisd_state`. In short, + the wavefunction is returned as a Python `dict` object with the format + `{(int_a, int_b) : coeff}`, where the binary representation of integers `int_a, int_b` + gives the Fock occupation vector over the molecular orbitals of the alpha and beta + electrons, respectively and `coeff` is the value of the CI coefficient. + + The wavefunction is constructed from the configuration interaction coefficients that are + computed by PySCF's `UCISD` object upon calling its `.run()` or `'kernel()` method. In that + computation, the wavefunction is assumed to have the form + + .. math:: + + \ket{\Psi_{CISD}} = c^{(0)} \ket{\text{HF}} + \sum_{i} c^{(1)}_{i} \ket{S_i} + \sum_{i} c^{(2)}_{i} \ket_{D_{i}} + + where :math:`c^{(0,1,2)}` are the + + Parameters + ---------- + PySCF CISD solver + + state + which state to do the conversion for, if multiple were solved for + + Returns + ------- + dict + Dictionary with keys (stra,strb) being binary strings representing + states occupied by alpha and beta electrons, and values being the CI + coefficients. + """ + + norb = mol.nao + nelec = mol.nelectron + nelec_a = int((nelec + mol.spin) / 2) + nelec_b = int((nelec - mol.spin) / 2) + nocc_a, nocc_b = nelec_a, nelec_b + nvir_a, nvir_b = norb - nelec_a, norb - nelec_b + + ## extract the CI coeffs for the chosen state, if multiple were solved for by CISD + assert state in range(cisd_solver.nroots), ( + f"State requested has not " f"been solved for. Re-run with larger nroots." + ) + if cisd_solver.nroots > 1: + cisdvec = cisd_solver.ci[state] + else: + cisdvec = cisd_solver.ci + + # get number of single/double excitations, to partition the cisdvec accordingly + # simple combinatorics (i.e. N orbitals choose n electrons ) gives the answer + size_a = nocc_a * nvir_a + size_b = nocc_b * nvir_b + size_aa = int(nocc_a * (nocc_a - 1) / 2) * int(nvir_a * (nvir_a - 1) / 2) + size_bb = int(nocc_b * (nocc_b - 1) / 2) * int(nvir_b * (nvir_b - 1) / 2) + size_ab = nocc_a * nocc_b * nvir_a * nvir_b + + # cisdvec in PySCF stores excitation coefficients as a flat array [c0, c1a, c1b, c2aa, c2ab, c2bb] + sizes = [1, size_a, size_b, size_aa, size_ab, size_bb] + cumul = np.cumsum(sizes) + idxs = [0] + [slice(cumul[ii], cumul[ii + 1]) for ii in range(len(cumul) - 1)] + + # quick check that the combinatorics is correct + assert np.sum(sizes) == len(cisdvec), ( + f"Number of expected excitations {np.sum(sizes)} does not " + f"correspond to cisdvec size {len(cisdvec)}.\nCheck you combinatorics!" + ) + + # np.cumsum() up to an appropriate size_xx will give the appropriate cisdvec slice + c0, c1a, c1b, c2aa, c2ab, c2bb = [cisdvec[idx] for idx in idxs] + + # with the coefficients c0/c1/c2, can build the wavefunction through helper functions + # use row to store integers representing the Fock vector of alpha electrons, and col for beta + row = [] + col = [] + coeff = [] + + ### Hartree-Fock coefficient ### + # HF vector given by bin(ref_a)[::-1] = 1111...10...0 + ref_a = int(2**nelec_a - 1) + ref_b = int(2**nelec_b - 1) + row.extend([ref_a]) + col.extend([ref_b]) + coeff.extend([c0]) + + ### Single excitations ### + # generate arrays of integers c1x_configs, whose binary form corresponds to + # allowed single excitations, and signs from acting with the excitation operator c^+ c on the HF reference + + # a -> a + # for alpha excitations, beta configs (cols) are unchanged + c1a_configs, c1a_signs = _excitated_states(norb, nelec_a, 1) + row.extend(c1a_configs) + col.extend([ref_b] * size_a) + coeff.extend(c1a * c1a_signs) + + # b -> b + # for beta excitations, alpha configs (rows) are unchanged + c1b_configs, c1b_signs = _excitated_states(norb, nelec_b, 1) + row.extend(c1b_configs) + col.extend([ref_a] * size_b) + coeff.extend(c1b * c1b_signs) + + ### Double excitations ### + # generate arrays of integers c2xx_configs, whose binary form corresponds to + # allowed doubl excitations, and signs from acting with the excitation operator c^+ c c^+ c on the HF reference + + ## aa -> aa + # for alpha excitations, beta configs (cols) are unchanged + c2aa_configs, c2aa_signs = _excitated_states(norb, nelec_a, 2) + row.extend(c2aa_configs) + col.extend([ref_b] * size_aa) + coeff.extend(c2aa * c2aa_signs) + + ## ab -> ab + # generate all possible pairwise combinations of _single_ excitations of alpha and beta sectors + rowvals, colvals = np.array(list(product(c1a_configs, c1b_configs))).T + row.extend(rowvals) + col.extend(colvals) + coeff.extend(c2ab * np.kron(c1a_signs, c1b_signs)) + + ## bb -> bb + # for beta excitations, alpha configs (rows) are unchanged + c2bb_configs, c2bb_signs = _excitated_states(norb, nelec_b, 2) + row.extend([ref_a] * size_bb) + row.extend(c2bb_configs) + coeff.extend(c2bb * c2bb_signs) + + # zip into a Slater det dictionary of the form { (int_a, int_b): coeff } + # where binary form of int_a / int_b gives the Fock occupation vector for alpha/beta sectors + dict_fcimatr = dict(zip(list(zip(row, col)), coeff)) + + # filter based on tolerance cutoff + dict_fcimatr = {key: value for key, value in dict_fcimatr.items() if abs(value) > tol} + + return dict_fcimatr + + +def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): + r"""Wrapper for constructing a wavefunction in the ``wf_dict`` format + ``{ (int_a, int_b): coeff}`` from PySCF's CISD class instance. + + The ``wf_dict`` format is a condensed, unified format for representing + wavefunctions written as a sum of Slater determinants. Each entry in + the dictionary corresponds to a single Slater determinant with coefficient + ``coeff``: the determinant is represented using the integers ``int_a, int_b``, + chosen such that such that their binary representation (i.e. ``bin(int_a)``) + corresponds to strings representing the Fock molecular orbital occupation + vector for the alpha (spin-up) sector and beta (spin-down) sector. + + For example, the Hartree-Fock state of the :math:`\text{H}_2` molecule in + the minimal basis with two electrons, split between the alpha and beta sector, + is written :math:`[1,0] x [1,0]` (high-energy orbitals on the right). In the + `wf_dict` format, it would correspond to integers (1, 1), since + `bin(1) = '0b1' -> [1,0]`. Note that the representation is reversed + (high-energy orbitals are on the left). The doubly excited state + :math:`[0,1] x [0,1]` would correspond to (2, 2), since `bin(2) = '0b10' -> [0,1]`. + + - The wrapper re-directs wavefunction construction to the appropriate method + depending on whether the restricted CISD or unrestricted CISD was used. + + - + + Args: + cisd_solver (PySCF CISD or UCISD Class instance): the class object representing + the CISD / UCISD calculation in PySCF. Must have already carried out the + calculation, e.g. by calling .kernel() or .run(). + + hftype (str): User specifies whether restricted or unrestricted CISD was performed. + The options are 'rhf' for restricted, and 'uhf' for unrestricted. + + state (int): which state to do the conversion for, if several were solved for + when running CISD / UCISD. Default is 0 (the ground state). + + tol (float): the tolerance to which the wavefunction is being built -- Slater + determinants with coefficients smaller than this are discarded. Default + is 1e-15 (all coefficients are included). + + Returns: + dict: Dictionary of the form `{(int_a, int_b) :coeff}`, with integers `int_a, int_b` + having binary represention corresponding to the Fock occupation vector in alpha and beta + spin sectors, respectively, and coeff being the CI coefficients of those configurations. + + **Example** + + >>> from pyscf import gto, scf, ci + >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,1.4)]]) + >>> myhf = scf.UHF(mol).run() + >>> myci = ci.UCISD(myhf).run() + >>> wf_cisd = cisd_state(myci, tol=1e-2) + >>> print(wf_cisd) + {(1, 1): -0.9486830343747924, (2, 2): -0.31622855704290276} + + """ + + if hftype == "uhf": + wf = _ucisd_state(cisd_solver, state=state, tol=tol) + # elif hftype == 'rhf': + # wf = _rcisd_state(cisd_solver, state=state, tol=tol) + else: + print(f"Unknown HF reference character. Exiting") + exit() + + return wf From cebd3488a31d07191726ee631f276ec18505656d Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 14:49:33 -0400 Subject: [PATCH 09/41] added the converter from wf_dict format to pennylane statevector --- pennylane/qchem/convert.py | 47 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 47 insertions(+) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 55acad25d1b..ada0ce6dbca 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -469,10 +469,13 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): coefficients. """ + mol = cisd_solver.mol + norb = mol.nao nelec = mol.nelectron nelec_a = int((nelec + mol.spin) / 2) nelec_b = int((nelec - mol.spin) / 2) + nocc_a, nocc_b = nelec_a, nelec_b nvir_a, nvir_b = norb - nelec_a, norb - nelec_b @@ -640,3 +643,47 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): exit() return wf + +def wfdict_to_statevector(wf_dict, norbs): + ''' + Intermediate converter between Overlapper's wavefunction_dict format + and Pennylane's 2^n statevector. + + Parameters + ---------- + + wf_dict : Overlapper wf_dict format + Dictionary with keys (int_a,int_b) being integers whose binary representation + shows states occupied by alpha and beta electrons; and values being the CI + coefficients. + + For example, int_a = 3 --> bin(3) = "0b11" -> [1100] means the + alpha electrons occuppy the lowest two orbitals. If int_b = 3 also, + then the total state can be written [2200], or [11110000] if showing + spin-orbitals explicitly. + + Returns + ------- + list + normalized statevector of length 2^(number_of_spinorbitals), just as in Pennylane + + ''' + + statevector = np.zeros(2**(2*norbs), dtype=complex) + + for (int_a, int_b), coeff in wf_dict.items(): + # convert to binary + bin_a, bin_b = bin(int_a)[2:], bin(int_b)[2:] + # and re-order (different PySCF vs Pennylane conventions) + bin_a, bin_b = bin_a[::-1], bin_b[::-1] + # interleave to get a spin-orbital form + bin_tot = "".join(i + j for i, j in zip(bin_a, bin_b)) + # pad to get the total number of orbitals + bin_tot_full = bin_tot + "0" * (2*norbs - len(bin_tot)) + idx = int(bin_tot_full, 2) + statevector[idx] += coeff + + # normalize + statevector = statevector / np.sqrt( np.sum( statevector**2 ) ) + + return statevector \ No newline at end of file From 1d6bca9065715d0d5caa145eaddea2a27db67128 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 17:10:11 -0400 Subject: [PATCH 10/41] update cisd_state docstring --- pennylane/qchem/convert.py | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index ada0ce6dbca..068fca4a1cb 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -22,6 +22,7 @@ from pennylane.wires import Wires from pennylane.operation import Tensor, active_new_opmath from pennylane.pauli import pauli_sentence +from itertools import product def _process_wires(wires, n_wires=None): @@ -530,14 +531,14 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): # a -> a # for alpha excitations, beta configs (cols) are unchanged - c1a_configs, c1a_signs = _excitated_states(norb, nelec_a, 1) + c1a_configs, c1a_signs = _excitated_states(nelec_a, norb, 1) row.extend(c1a_configs) col.extend([ref_b] * size_a) coeff.extend(c1a * c1a_signs) # b -> b # for beta excitations, alpha configs (rows) are unchanged - c1b_configs, c1b_signs = _excitated_states(norb, nelec_b, 1) + c1b_configs, c1b_signs = _excitated_states(nelec_b, norb, 1) row.extend(c1b_configs) col.extend([ref_a] * size_b) coeff.extend(c1b * c1b_signs) @@ -548,21 +549,21 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): ## aa -> aa # for alpha excitations, beta configs (cols) are unchanged - c2aa_configs, c2aa_signs = _excitated_states(norb, nelec_a, 2) + c2aa_configs, c2aa_signs = _excitated_states(nelec_a, norb, 2) row.extend(c2aa_configs) col.extend([ref_b] * size_aa) coeff.extend(c2aa * c2aa_signs) ## ab -> ab # generate all possible pairwise combinations of _single_ excitations of alpha and beta sectors - rowvals, colvals = np.array(list(product(c1a_configs, c1b_configs))).T + rowvals, colvals = np.array( list(product(c1a_configs, c1b_configs)), dtype=int ).T.numpy() row.extend(rowvals) col.extend(colvals) - coeff.extend(c2ab * np.kron(c1a_signs, c1b_signs)) + coeff.extend( c2ab * np.kron(c1a_signs, c1b_signs).numpy() ) ## bb -> bb # for beta excitations, alpha configs (rows) are unchanged - c2bb_configs, c2bb_signs = _excitated_states(norb, nelec_b, 2) + c2bb_configs, c2bb_signs = _excitated_states(nelec_b, norb, 2) row.extend([ref_a] * size_bb) row.extend(c2bb_configs) coeff.extend(c2bb * c2bb_signs) @@ -600,7 +601,11 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): - The wrapper re-directs wavefunction construction to the appropriate method depending on whether the restricted CISD or unrestricted CISD was used. - - + - The (private) method converts the CISD object's attribute `.ci`, which stores the + Slater determinant coefficients, into the `wf_dict` described above. Optionally, + if multiple states were computed with PySCF CISD, an excited state's wavefunction + can be obtained by passing `state = K`, where :math:`K` points to the :math:`K`'th + excited state. Args: cisd_solver (PySCF CISD or UCISD Class instance): the class object representing @@ -625,7 +630,7 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): **Example** >>> from pyscf import gto, scf, ci - >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,1.4)]]) + >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]]) >>> myhf = scf.UHF(mol).run() >>> myci = ci.UCISD(myhf).run() >>> wf_cisd = cisd_state(myci, tol=1e-2) From 16de4e4e45e706e5ebcd42809a2a888d945fc248 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 17:26:04 -0400 Subject: [PATCH 11/41] update _ucisd_state docstring --- pennylane/qchem/convert.py | 54 +++++++++++++++++++++++++++++--------- 1 file changed, 42 insertions(+), 12 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 068fca4a1cb..ac05cb83caf 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -453,23 +453,53 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): \ket{\Psi_{CISD}} = c^{(0)} \ket{\text{HF}} + \sum_{i} c^{(1)}_{i} \ket{S_i} + \sum_{i} c^{(2)}_{i} \ket_{D_{i}} - where :math:`c^{(0,1,2)}` are the + where :math:`c^{(0,1,2)}` are the configuration interaction coefficients stored + as the `.ci` attribute of PySCF's `UCISD` object, and :math:`\ket{S_i}, \ket{D_i}` + are the singly and doubly excited Slater determinants, respectively. - Parameters - ---------- - PySCF CISD solver + - The first step is to partition the `.ci` attribute storing the CI coefficients into the :math:`c^{(0,1,2)}` + coefficients. To do this, the number of allowed single and double excitations is computed through combinatorics. + For example, the number of allowed single excitations of :math:`n_e` alpha electrons in :math:`N` orbitals is + :math:`\binom{N}{n_e}`. - state - which state to do the conversion for, if multiple were solved for + - The :func:`~._excited_states` function is then called to generate the configurations corresponding to allowed + single and double excitations, as well as their corresponding signs from permuting the creation operators. These + configurations are again represented by integers whose binary form gives the Fock occupation vector. + + - The configurations for alpha electrons (`row`) and beta electrons (`col`) are collected for all excitations, + and the corresponding coefficients are collected (`coeff`). These three arrays are then used to build the final + `wf_dict`, which is filtered to drop the smallest coefficients according to the `tol` argument. + + Args: + cisd_solver (PySCF UCISD Class instance): the class object representing + the UCISD calculation in PySCF. Must have already carried out the + calculation, e.g. by calling .kernel() or .run(). + + state (int): which state to do the conversion for, if several were solved for + when running UCISD. Default is 0 (the ground state). + + tol (float): the tolerance to which the wavefunction is being built -- Slater + determinants with coefficients smaller than this are discarded. Default + is 1e-15 (all coefficients are included). + + Returns: + dict: Dictionary of the form `{(int_a, int_b) :coeff}`, with integers `int_a, int_b` + having binary represention corresponding to the Fock occupation vector in alpha and beta + spin sectors, respectively, and coeff being the CI coefficients of those configurations. + + **Example** + + >>> from pyscf import gto, scf, ci + >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]]) + >>> myhf = scf.UHF(mol).run() + >>> myci = ci.UCISD(myhf).run() + >>> wf_cisd = cisd_state(myci, tol=1e-2) + >>> print(wf_cisd) + {(1, 1): -0.9486830343747924, (2, 2): -0.31622855704290276} - Returns - ------- - dict - Dictionary with keys (stra,strb) being binary strings representing - states occupied by alpha and beta electrons, and values being the CI - coefficients. """ + mol = cisd_solver.mol norb = mol.nao From 1d29b1f73f86b931449df817d1705bb2b5708933 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 17:30:55 -0400 Subject: [PATCH 12/41] update docstring of final converter to statevector --- pennylane/qchem/convert.py | 32 +++++++++++--------------------- 1 file changed, 11 insertions(+), 21 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index ac05cb83caf..7d02221a89f 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -680,29 +680,19 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): return wf def wfdict_to_statevector(wf_dict, norbs): - ''' - Intermediate converter between Overlapper's wavefunction_dict format - and Pennylane's 2^n statevector. + r"""Final converter between wf_dict format and Pennylane's statevector. - Parameters - ---------- - - wf_dict : Overlapper wf_dict format - Dictionary with keys (int_a,int_b) being integers whose binary representation - shows states occupied by alpha and beta electrons; and values being the CI - coefficients. + Args: + wf_dict (wf_dict format): dictionary with keys (int_a,int_b) being integers whose binary representation + shows the Fock occupation vector for alpha and beta electrons; and values being the CI + coefficients. - For example, int_a = 3 --> bin(3) = "0b11" -> [1100] means the - alpha electrons occuppy the lowest two orbitals. If int_b = 3 also, - then the total state can be written [2200], or [11110000] if showing - spin-orbitals explicitly. - - Returns - ------- - list - normalized statevector of length 2^(number_of_spinorbitals), just as in Pennylane - - ''' + norbs (int): number of molecular orbitals in the system + + Returns: + statevector (list): normalized statevector of length 2^(number_of_spinorbitals), standard in Pennylane + + """ statevector = np.zeros(2**(2*norbs), dtype=complex) From f76de224c604a1278f7e80cec6940f3eb6eb9f6d Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 18:08:30 -0400 Subject: [PATCH 13/41] fixed the example output for cisd_state --- pennylane/qchem/convert.py | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 7d02221a89f..0320b96b8f7 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -490,12 +490,12 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): **Example** >>> from pyscf import gto, scf, ci - >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]]) + >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g') >>> myhf = scf.UHF(mol).run() >>> myci = ci.UCISD(myhf).run() - >>> wf_cisd = cisd_state(myci, tol=1e-2) + >>> wf_cisd = cisd_state(myci, tol=1e-1) >>> print(wf_cisd) - {(1, 1): -0.9486830343747924, (2, 2): -0.31622855704290276} + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} """ @@ -660,13 +660,12 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): **Example** >>> from pyscf import gto, scf, ci - >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]]) + >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g') >>> myhf = scf.UHF(mol).run() >>> myci = ci.UCISD(myhf).run() - >>> wf_cisd = cisd_state(myci, tol=1e-2) + >>> wf_cisd = cisd_state(myci, tol=1e-1) >>> print(wf_cisd) - {(1, 1): -0.9486830343747924, (2, 2): -0.31622855704290276} - + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} """ if hftype == "uhf": @@ -691,7 +690,7 @@ def wfdict_to_statevector(wf_dict, norbs): Returns: statevector (list): normalized statevector of length 2^(number_of_spinorbitals), standard in Pennylane - + """ statevector = np.zeros(2**(2*norbs), dtype=complex) From 4db87786fb6b57ca186b6fbc67430448e1ae5938 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Tue, 1 Aug 2023 21:35:31 -0400 Subject: [PATCH 14/41] added test for private _ucisd func based on h2 molecule --- tests/qchem/of_tests/test_convert.py | 38 ++++++++++++++++++++++++++++ 1 file changed, 38 insertions(+) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 6b3ce829065..89a262e147c 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -857,3 +857,41 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref states, signs = qchem.convert._excitated_states(electrons, orbitals, excitation) assert np.allclose(states, states_ref) assert np.allclose(signs, signs_ref) + +from pyscf import gto +@pytest.mark.parametrize( + ("mol", "hftype", "state", "tol", "wf_ref"), + [ + # reference data computed with pyscf -- identify the FCI matrix entries + # with the correct Fock occupation vector entries + # mycasci = mcscf.CASCI(hf, ncas = orbitals, nelecas = electrons).run() + # print(mycasci.ci) + ( + gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g'), + 'uhf', + 0, + 1e-1, + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} + ), + # ( + # gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='cc-pvdz'), + # 'uhf', + # 0, + # 4e-2, + # {(1, 1): 0.9919704779918753, (2, 2): 0.048530387857156514, (2, 8): 0.04452334895452825,\ + # (4, 4): -0.050035993275715764, (8, 2): -0.04452335710141664, (8, 8): -0.052262296228815036,\ + # (16, 16): -0.04045330124484208, (32, 32): -0.04045330124484209} + # ), + ], +) +def test_private_ucisd_state(mol, hftype, state, tol, wf_ref): + r"""Test the UCISD wavefunction constructor.""" + + from pyscf import scf, ci + + myhf = scf.UHF(mol).run() + myci = ci.UCISD(myhf).run(state=state) + + wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) + + assert { key: round(val, 4) for (key, val) in wf_cisd.items() } == { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file From 611715853c2e7826982f14a725ed96a3babbc3f9 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Wed, 2 Aug 2023 09:05:54 -0400 Subject: [PATCH 15/41] added cc-pvdz test --- tests/qchem/of_tests/test_convert.py | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 89a262e147c..c30ec1717c8 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -873,15 +873,15 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref 1e-1, {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} ), - # ( - # gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='cc-pvdz'), - # 'uhf', - # 0, - # 4e-2, - # {(1, 1): 0.9919704779918753, (2, 2): 0.048530387857156514, (2, 8): 0.04452334895452825,\ - # (4, 4): -0.050035993275715764, (8, 2): -0.04452335710141664, (8, 8): -0.052262296228815036,\ - # (16, 16): -0.04045330124484208, (32, 32): -0.04045330124484209} - # ), + ( + gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='cc-pvdz'), + 'uhf', + 0, + 4e-2, + {(1, 1): 0.9919704779918753, (2, 2): 0.048530387857156514, (2, 8): 0.04452334895452825,\ + (4, 4): -0.050035993275715764, (8, 2): -0.04452335710141664, (8, 8): -0.052262296228815036,\ + (16, 16): -0.04045330124484208, (32, 32): -0.04045330124484209} + ), ], ) def test_private_ucisd_state(mol, hftype, state, tol, wf_ref): From 9f26454bf994dde6a484f8a018ab04c386cd4e0b Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Wed, 2 Aug 2023 09:55:16 -0400 Subject: [PATCH 16/41] added how to generate ref data for test_prive_ucisd --- pennylane/qchem/convert.py | 6 +++--- tests/qchem/of_tests/test_convert.py | 30 ++++++++++++++++------------ 2 files changed, 20 insertions(+), 16 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 0320b96b8f7..f63a390b25e 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -569,8 +569,8 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): # b -> b # for beta excitations, alpha configs (rows) are unchanged c1b_configs, c1b_signs = _excitated_states(nelec_b, norb, 1) - row.extend(c1b_configs) - col.extend([ref_a] * size_b) + row.extend([ref_a] * size_b) + col.extend(c1b_configs) coeff.extend(c1b * c1b_signs) ### Double excitations ### @@ -595,7 +595,7 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): # for beta excitations, alpha configs (rows) are unchanged c2bb_configs, c2bb_signs = _excitated_states(nelec_b, norb, 2) row.extend([ref_a] * size_bb) - row.extend(c2bb_configs) + col.extend(c2bb_configs) coeff.extend(c2bb * c2bb_signs) # zip into a Slater det dictionary of the form { (int_a, int_b): coeff } diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index c30ec1717c8..b85d89f4691 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -858,40 +858,44 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref assert np.allclose(states, states_ref) assert np.allclose(signs, signs_ref) -from pyscf import gto @pytest.mark.parametrize( - ("mol", "hftype", "state", "tol", "wf_ref"), + ("atom", "basis", "hftype", "state", "tol", "wf_ref"), [ # reference data computed with pyscf -- identify the FCI matrix entries # with the correct Fock occupation vector entries - # mycasci = mcscf.CASCI(hf, ncas = orbitals, nelecas = electrons).run() - # print(mycasci.ci) + # mycasci = mcscf.UCASCI(myhf, norbs, (nelec_a, nelec_b)).run() + # sparse_cascimatr = coo_matrix(mycasci.ci, shape=np.shape(cascivec), dtype=float) + # row, col, dat = sparse_cascimatr.row, sparse_cascimatr.col, sparse_cascimatr.data + # strs_row, strs_col = addrs2str(norbs, nelec_a, row), addrs2str(norbs, nelec_b, col) + # wf_ref = dict(zip(list(zip(strs_row, strs_col)), dat)) ( - gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g'), + [['H', (0, 0, 0)], ['H', (0,0,0.71)]], + 'sto6g', 'uhf', 0, 1e-1, {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} ), ( - gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='cc-pvdz'), + [['H', (0, 0, 0)], ['H', (0, 0, 0.71)]], + 'cc-pvdz', 'uhf', 0, 4e-2, - {(1, 1): 0.9919704779918753, (2, 2): 0.048530387857156514, (2, 8): 0.04452334895452825,\ - (4, 4): -0.050035993275715764, (8, 2): -0.04452335710141664, (8, 8): -0.052262296228815036,\ - (16, 16): -0.04045330124484208, (32, 32): -0.04045330124484209} - ), + {(1, 1): 0.9919704801055201, (2, 2): 0.04853035090478132, (2, 8): -0.04452332640264183, \ + (4, 4): -0.05003594588609278, (8, 2): 0.044523334555907366, (8, 8): -0.05226230845395026}), ], ) -def test_private_ucisd_state(mol, hftype, state, tol, wf_ref): +def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): r"""Test the UCISD wavefunction constructor.""" - from pyscf import scf, ci + from pyscf import gto, scf, ci + mol = gto.M(atom=atom, basis=basis) myhf = scf.UHF(mol).run() myci = ci.UCISD(myhf).run(state=state) wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) - assert { key: round(val, 4) for (key, val) in wf_cisd.items() } == { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file + assert { key: round(val, 4) for (key, val) in wf_cisd.items() } ==\ + { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file From d7168fd8e0bf61420a5b2cff79bfe57fcf7e3827 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Wed, 2 Aug 2023 09:56:02 -0400 Subject: [PATCH 17/41] [skip ci] fixed typos --- tests/qchem/of_tests/test_convert.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index b85d89f4691..c7c8c1735ff 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -897,5 +897,5 @@ def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) - assert { key: round(val, 4) for (key, val) in wf_cisd.items() } ==\ - { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file + assert { key: round(val, 4) for (key, val) in wf_cisd.items() } == \ + { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file From 414965a98a8093671e2886480d7cd1d9830d2ae1 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 14:17:51 -0400 Subject: [PATCH 18/41] modify test to account for sign change --- tests/qchem/of_tests/test_convert.py | 30 ++++++++++++++++++---------- 1 file changed, 19 insertions(+), 11 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index c7c8c1735ff..0f100d346ce 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -858,6 +858,7 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref assert np.allclose(states, states_ref) assert np.allclose(signs, signs_ref) + @pytest.mark.parametrize( ("atom", "basis", "hftype", "state", "tol", "wf_ref"), [ @@ -869,21 +870,28 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref # strs_row, strs_col = addrs2str(norbs, nelec_a, row), addrs2str(norbs, nelec_b, col) # wf_ref = dict(zip(list(zip(strs_row, strs_col)), dat)) ( - [['H', (0, 0, 0)], ['H', (0,0,0.71)]], - 'sto6g', - 'uhf', + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "sto6g", + "uhf", 0, 1e-1, - {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159}, ), ( - [['H', (0, 0, 0)], ['H', (0, 0, 0.71)]], - 'cc-pvdz', - 'uhf', + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "cc-pvdz", + "uhf", 0, 4e-2, - {(1, 1): 0.9919704801055201, (2, 2): 0.04853035090478132, (2, 8): -0.04452332640264183, \ - (4, 4): -0.05003594588609278, (8, 2): 0.044523334555907366, (8, 8): -0.05226230845395026}), + { + (1, 1): 0.9919704801055201, + (2, 2): 0.04853035090478132, + (2, 8): -0.04452332640264183, + (4, 4): -0.05003594588609278, + (8, 2): 0.044523334555907366, + (8, 8): -0.05226230845395026, + }, + ), ], ) def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): @@ -897,5 +905,5 @@ def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) - assert { key: round(val, 4) for (key, val) in wf_cisd.items() } == \ - { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file + assert wf_cisd.keys() == wf_ref.keys() + assert np.allclose(abs(np.array(list(wf_cisd.values()))), abs(np.array(list(wf_ref.values())))) From 782768c4212fa45ce2cd7a4a17b78b3bb31e5665 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 15:18:00 -0400 Subject: [PATCH 19/41] remove comments and modify docstring --- pennylane/qchem/convert.py | 230 +++++++++---------------------------- 1 file changed, 55 insertions(+), 175 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index f63a390b25e..9d6c1564fe5 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -22,7 +22,7 @@ from pennylane.wires import Wires from pennylane.operation import Tensor, active_new_opmath from pennylane.pauli import pauli_sentence -from itertools import product +from itertools import product def _process_wires(wires, n_wires=None): @@ -435,52 +435,25 @@ def _excitated_states(electrons, orbitals, excitation): return states_int, signs -def _ucisd_state(cisd_solver, state=0, tol=1e-15): - r"""Construct a wavefunction in the wf_dict format `{ (int_a, int_b): coeff}` from - PySCF's `UCISD` Solver object. +def _ucisd_state(cisd_solver, tol=1e-15): + r"""Construct a wavefunction from PySCF's `UCISD` Solver object. - The format is described in detail in the docstring of :func:`~.cisd_state`. In short, - the wavefunction is returned as a Python `dict` object with the format - `{(int_a, int_b) : coeff}`, where the binary representation of integers `int_a, int_b` - gives the Fock occupation vector over the molecular orbitals of the alpha and beta - electrons, respectively and `coeff` is the value of the CI coefficient. - - The wavefunction is constructed from the configuration interaction coefficients that are - computed by PySCF's `UCISD` object upon calling its `.run()` or `'kernel()` method. In that - computation, the wavefunction is assumed to have the form - - .. math:: - - \ket{\Psi_{CISD}} = c^{(0)} \ket{\text{HF}} + \sum_{i} c^{(1)}_{i} \ket{S_i} + \sum_{i} c^{(2)}_{i} \ket_{D_{i}} - - where :math:`c^{(0,1,2)}` are the configuration interaction coefficients stored - as the `.ci` attribute of PySCF's `UCISD` object, and :math:`\ket{S_i}, \ket{D_i}` - are the singly and doubly excited Slater determinants, respectively. - - - The first step is to partition the `.ci` attribute storing the CI coefficients into the :math:`c^{(0,1,2)}` - coefficients. To do this, the number of allowed single and double excitations is computed through combinatorics. - For example, the number of allowed single excitations of :math:`n_e` alpha electrons in :math:`N` orbitals is - :math:`\binom{N}{n_e}`. - - - The :func:`~._excited_states` function is then called to generate the configurations corresponding to allowed - single and double excitations, as well as their corresponding signs from permuting the creation operators. These - configurations are again represented by integers whose binary form gives the Fock occupation vector. - - - The configurations for alpha electrons (`row`) and beta electrons (`col`) are collected for all excitations, - and the corresponding coefficients are collected (`coeff`). These three arrays are then used to build the final - `wf_dict`, which is filtered to drop the smallest coefficients according to the `tol` argument. + The wavefunction is represented as a dictionary where the keys are tuples representing a + configuration, which corresponds to a Slater determinant, and the values are the CI coefficient + corresponding to that Slater determinant. Each dictionary key is a tuple of two integers. The + binary representation of these integers correspond to a specific configuration, or Slater + determinant. The first number represents the configuration of alpha electrons and the second + number represents the beta electrons. For instance, the Hartree-Fock state + :math:`|1 1 0 0 \rangle` will be represented by the binary string `0011`. This string can be + splited to `01` and `01` for the alpha and beta electrons. The integer corresponding to `01` is + `1`. Then the dictionary representation of a state with `0.99` contribution of the Hartree-Fock + state and `0.01` contribution from the doubly-excited state will be + `{(1, 1): 0.99, (2, 2): 0.01}`. Args: - cisd_solver (PySCF UCISD Class instance): the class object representing - the UCISD calculation in PySCF. Must have already carried out the - calculation, e.g. by calling .kernel() or .run(). - - state (int): which state to do the conversion for, if several were solved for - when running UCISD. Default is 0 (the ground state). - - tol (float): the tolerance to which the wavefunction is being built -- Slater - determinants with coefficients smaller than this are discarded. Default - is 1e-15 (all coefficients are included). + cisd_solver (PySCF UCISD Class instance): the class object representing the UCISD calculation in PySCF + tol (float): the tolerance to which the wavefunction is being built -- Slater determinants + with coefficients smaller than this are discarded. Default is 1e-15 (all coefficients are included). Returns: dict: Dictionary of the form `{(int_a, int_b) :coeff}`, with integers `int_a, int_b` @@ -496,111 +469,54 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): >>> wf_cisd = cisd_state(myci, tol=1e-1) >>> print(wf_cisd) {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} - """ - - mol = cisd_solver.mol + cisdvec = cisd_solver.ci norb = mol.nao nelec = mol.nelectron nelec_a = int((nelec + mol.spin) / 2) nelec_b = int((nelec - mol.spin) / 2) - nocc_a, nocc_b = nelec_a, nelec_b nvir_a, nvir_b = norb - nelec_a, norb - nelec_b - ## extract the CI coeffs for the chosen state, if multiple were solved for by CISD - assert state in range(cisd_solver.nroots), ( - f"State requested has not " f"been solved for. Re-run with larger nroots." - ) - if cisd_solver.nroots > 1: - cisdvec = cisd_solver.ci[state] - else: - cisdvec = cisd_solver.ci - - # get number of single/double excitations, to partition the cisdvec accordingly - # simple combinatorics (i.e. N orbitals choose n electrons ) gives the answer - size_a = nocc_a * nvir_a - size_b = nocc_b * nvir_b - size_aa = int(nocc_a * (nocc_a - 1) / 2) * int(nvir_a * (nvir_a - 1) / 2) - size_bb = int(nocc_b * (nocc_b - 1) / 2) * int(nvir_b * (nvir_b - 1) / 2) - size_ab = nocc_a * nocc_b * nvir_a * nvir_b + size_a, size_b = nelec_a * nvir_a, nelec_b * nvir_b + size_aa = int(nelec_a * (nelec_a - 1) / 2) * int(nvir_a * (nvir_a - 1) / 2) + size_bb = int(nelec_b * (nelec_b - 1) / 2) * int(nvir_b * (nvir_b - 1) / 2) + size_ab = nelec_a * nelec_b * nvir_a * nvir_b - # cisdvec in PySCF stores excitation coefficients as a flat array [c0, c1a, c1b, c2aa, c2ab, c2bb] sizes = [1, size_a, size_b, size_aa, size_ab, size_bb] cumul = np.cumsum(sizes) idxs = [0] + [slice(cumul[ii], cumul[ii + 1]) for ii in range(len(cumul) - 1)] - - # quick check that the combinatorics is correct - assert np.sum(sizes) == len(cisdvec), ( - f"Number of expected excitations {np.sum(sizes)} does not " - f"correspond to cisdvec size {len(cisdvec)}.\nCheck you combinatorics!" - ) - - # np.cumsum() up to an appropriate size_xx will give the appropriate cisdvec slice c0, c1a, c1b, c2aa, c2ab, c2bb = [cisdvec[idx] for idx in idxs] - # with the coefficients c0/c1/c2, can build the wavefunction through helper functions - # use row to store integers representing the Fock vector of alpha electrons, and col for beta - row = [] - col = [] - coeff = [] - - ### Hartree-Fock coefficient ### - # HF vector given by bin(ref_a)[::-1] = 1111...10...0 + # numbers representing the Hartree-Fock vector, e.g., bin(ref_a)[::-1] = 1111...10...0 ref_a = int(2**nelec_a - 1) ref_b = int(2**nelec_b - 1) - row.extend([ref_a]) - col.extend([ref_b]) - coeff.extend([c0]) - ### Single excitations ### - # generate arrays of integers c1x_configs, whose binary form corresponds to - # allowed single excitations, and signs from acting with the excitation operator c^+ c on the HF reference + dict_fcimatr = dict(zip(list(zip([ref_a], [ref_b])), [c0])) - # a -> a - # for alpha excitations, beta configs (cols) are unchanged + # alpha -> alpha excitations c1a_configs, c1a_signs = _excitated_states(nelec_a, norb, 1) - row.extend(c1a_configs) - col.extend([ref_b] * size_a) - coeff.extend(c1a * c1a_signs) + dict_fcimatr.update(dict(zip(list(zip(c1a_configs, [ref_b] * size_a)), c1a * c1a_signs))) - # b -> b - # for beta excitations, alpha configs (rows) are unchanged + # beta -> beta excitations c1b_configs, c1b_signs = _excitated_states(nelec_b, norb, 1) - row.extend([ref_a] * size_b) - col.extend(c1b_configs) - coeff.extend(c1b * c1b_signs) + dict_fcimatr.update(dict(zip(list(zip(c1b_configs, [ref_a] * size_b)), c1b * c1b_signs))) - ### Double excitations ### - # generate arrays of integers c2xx_configs, whose binary form corresponds to - # allowed doubl excitations, and signs from acting with the excitation operator c^+ c c^+ c on the HF reference - - ## aa -> aa - # for alpha excitations, beta configs (cols) are unchanged + # alpha, alpha -> alpha, alpha excitations c2aa_configs, c2aa_signs = _excitated_states(nelec_a, norb, 2) - row.extend(c2aa_configs) - col.extend([ref_b] * size_aa) - coeff.extend(c2aa * c2aa_signs) - - ## ab -> ab - # generate all possible pairwise combinations of _single_ excitations of alpha and beta sectors - rowvals, colvals = np.array( list(product(c1a_configs, c1b_configs)), dtype=int ).T.numpy() - row.extend(rowvals) - col.extend(colvals) - coeff.extend( c2ab * np.kron(c1a_signs, c1b_signs).numpy() ) - - ## bb -> bb - # for beta excitations, alpha configs (rows) are unchanged - c2bb_configs, c2bb_signs = _excitated_states(nelec_b, norb, 2) - row.extend([ref_a] * size_bb) - col.extend(c2bb_configs) - coeff.extend(c2bb * c2bb_signs) + dict_fcimatr.update(dict(zip(list(zip(c2aa_configs, [ref_b] * size_aa)), c2aa * c2aa_signs))) - # zip into a Slater det dictionary of the form { (int_a, int_b): coeff } - # where binary form of int_a / int_b gives the Fock occupation vector for alpha/beta sectors - dict_fcimatr = dict(zip(list(zip(row, col)), coeff)) + # alpha, beta -> alpha, beta excitations + rowvals, colvals = np.array(list(product(c1a_configs, c1b_configs)), dtype=int).T.numpy() + dict_fcimatr.update( + dict(zip(list(zip(rowvals, colvals)), c2ab * np.kron(c1a_signs, c1b_signs).numpy())) + ) + + # beta, beta -> beta, beta excitations + c2bb_configs, c2bb_signs = _excitated_states(nelec_b, norb, 2) + dict_fcimatr.update(dict(zip(list(zip([ref_a] * size_bb, c2bb_configs)), c2bb * c2bb_signs))) # filter based on tolerance cutoff dict_fcimatr = {key: value for key, value in dict_fcimatr.items() if abs(value) > tol} @@ -609,45 +525,15 @@ def _ucisd_state(cisd_solver, state=0, tol=1e-15): def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): - r"""Wrapper for constructing a wavefunction in the ``wf_dict`` format - ``{ (int_a, int_b): coeff}`` from PySCF's CISD class instance. - - The ``wf_dict`` format is a condensed, unified format for representing - wavefunctions written as a sum of Slater determinants. Each entry in - the dictionary corresponds to a single Slater determinant with coefficient - ``coeff``: the determinant is represented using the integers ``int_a, int_b``, - chosen such that such that their binary representation (i.e. ``bin(int_a)``) - corresponds to strings representing the Fock molecular orbital occupation - vector for the alpha (spin-up) sector and beta (spin-down) sector. - - For example, the Hartree-Fock state of the :math:`\text{H}_2` molecule in - the minimal basis with two electrons, split between the alpha and beta sector, - is written :math:`[1,0] x [1,0]` (high-energy orbitals on the right). In the - `wf_dict` format, it would correspond to integers (1, 1), since - `bin(1) = '0b1' -> [1,0]`. Note that the representation is reversed - (high-energy orbitals are on the left). The doubly excited state - :math:`[0,1] x [0,1]` would correspond to (2, 2), since `bin(2) = '0b10' -> [0,1]`. - - - The wrapper re-directs wavefunction construction to the appropriate method - depending on whether the restricted CISD or unrestricted CISD was used. - - - The (private) method converts the CISD object's attribute `.ci`, which stores the - Slater determinant coefficients, into the `wf_dict` described above. Optionally, - if multiple states were computed with PySCF CISD, an excited state's wavefunction - can be obtained by passing `state = K`, where :math:`K` points to the :math:`K`'th - excited state. + r"""Construct a wavefunction from PySCF output. Args: - cisd_solver (PySCF CISD or UCISD Class instance): the class object representing - the CISD / UCISD calculation in PySCF. Must have already carried out the - calculation, e.g. by calling .kernel() or .run(). + cisd_solver (PySCF CISD or UCISD Class instance): the class object representing the + CISD / UCISD calculation in PySCF - hftype (str): User specifies whether restricted or unrestricted CISD was performed. + hftype (str): keyword specifying whether restricted or unrestricted CISD was performed. The options are 'rhf' for restricted, and 'uhf' for unrestricted. - state (int): which state to do the conversion for, if several were solved for - when running CISD / UCISD. Default is 0 (the ground state). - tol (float): the tolerance to which the wavefunction is being built -- Slater determinants with coefficients smaller than this are discarded. Default is 1e-15 (all coefficients are included). @@ -669,45 +555,39 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): """ if hftype == "uhf": - wf = _ucisd_state(cisd_solver, state=state, tol=tol) - # elif hftype == 'rhf': - # wf = _rcisd_state(cisd_solver, state=state, tol=tol) + wf_dict = _ucisd_state(cisd_solver, state=state, tol=tol) else: - print(f"Unknown HF reference character. Exiting") - exit() + raise ValueError("Only restricted, 'rhf', and unrestricted, 'uhf', are supported.") + wf = wfdict_to_statevector(wf_dict, cisd_solver.mol.nao) return wf + def wfdict_to_statevector(wf_dict, norbs): - r"""Final converter between wf_dict format and Pennylane's statevector. + r"""Convert a wavefunction in sparce dictionary format to a Pennylane's statevector. Args: wf_dict (wf_dict format): dictionary with keys (int_a,int_b) being integers whose binary representation shows the Fock occupation vector for alpha and beta electrons; and values being the CI - coefficients. - + coefficients. + norbs (int): number of molecular orbitals in the system - + Returns: statevector (list): normalized statevector of length 2^(number_of_spinorbitals), standard in Pennylane """ - statevector = np.zeros(2**(2*norbs), dtype=complex) + statevector = np.zeros(2 ** (2 * norbs), dtype=complex) for (int_a, int_b), coeff in wf_dict.items(): - # convert to binary bin_a, bin_b = bin(int_a)[2:], bin(int_b)[2:] - # and re-order (different PySCF vs Pennylane conventions) bin_a, bin_b = bin_a[::-1], bin_b[::-1] - # interleave to get a spin-orbital form bin_tot = "".join(i + j for i, j in zip(bin_a, bin_b)) - # pad to get the total number of orbitals - bin_tot_full = bin_tot + "0" * (2*norbs - len(bin_tot)) + bin_tot_full = bin_tot + "0" * (2 * norbs - len(bin_tot)) idx = int(bin_tot_full, 2) statevector[idx] += coeff - # normalize - statevector = statevector / np.sqrt( np.sum( statevector**2 ) ) + statevector = statevector / np.sqrt(np.sum(statevector**2)) - return statevector \ No newline at end of file + return statevector From 4acf9c8e7741ba61c1508c0dd8e0870458aa8464 Mon Sep 17 00:00:00 2001 From: soranjh <40344468+soranjh@users.noreply.github.com> Date: Wed, 2 Aug 2023 15:28:15 -0400 Subject: [PATCH 20/41] Add functions to return excited basis states (#4417) --- pennylane/qchem/convert.py | 30 +++++++++++++++++++++++++++- tests/qchem/of_tests/test_convert.py | 10 +++++----- 2 files changed, 34 insertions(+), 6 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 9d6c1564fe5..91b4827f85c 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -360,6 +360,25 @@ def import_operator(qubit_observable, format="openfermion", wires=None, tol=1e01 def _excitations(electrons, orbitals): r"""Generate all possible single and double excitations from a Hartree-Fock reference state. + This function is a more performant version of `qchem.excitations`, where the order of the + generated excitations is consistent with PySCF. + + Single and double excitations can be generated by acting with the operators + :math:`\hat T_1` and :math:`\hat T_2` on the Hartree-Fock reference state: + + .. math:: + + && \hat{T}_1 = \sum_{p \in \mathrm{occ} \\ q \in \mathrm{unocc}} + \hat{c}_q^\dagger \hat{c}_p \\ + && \hat{T}_2 = \sum_{p>q \in \mathrm{occ} \\ r>s \in + \mathrm{unocc}} \hat{c}_r^\dagger \hat{c}_s^\dagger \hat{c}_p \hat{c}_q. + + + In the equations above the indices :math:`p, q, r, s` run over the + occupied (occ) and unoccupied (unocc) *spin* orbitals and :math:`\hat c` and + :math:`\hat c^\dagger` are the electron annihilation and creation operators, + respectively. + Args: electrons (int): number of electrons orbitals (int): number of spin orbitals @@ -390,9 +409,18 @@ def _excitations(electrons, orbitals): return singles, doubles -def _excitated_states(electrons, orbitals, excitation): +def _excited_configurations(electrons, orbitals, excitation): r"""Generate excited states from a Hartree-Fock reference state. + This function generates excited states in the form of binary strings or integers representing + them, from a Hartree-Fock (HF) reference state. The HF state is assumed to be + :math:`|1 1 ...1 0 ... 0 0 \rangle` where the number of :math:`1` and :math:`0` elements + represents the number of occupied and unoccupied spin orbitals, respectively. The string + representation of the state is obtained by converting the occupation-number vector to a string, + e.g., ``111000`` represents :math:`|1 1 1 0 0 0 \rangle. The number representation of the state + is the integer corresonding to the binary number, e.g., ``111000`` is represented by + :math:`int('111000', 2) = 56`. + Args: electrons (int): number of electrons orbitals (int): number of spin orbitals diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 0f100d346ce..1e27f09c56c 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -818,7 +818,7 @@ def dummy_ansatz(phis, wires): @pytest.mark.parametrize( ("electrons", "orbitals", "singles_ref", "doubles_ref"), [ - # trivial case + # trivial case, e.g., H2/STO-3G (2, 4, [[0, 2], [0, 3], [1, 2], [1, 3]], [[0, 1, 2, 3]]), ], ) @@ -833,7 +833,7 @@ def test_excitations(electrons, orbitals, singles_ref, doubles_ref): ("electrons", "orbitals", "excitation", "states_ref", "signs_ref"), [ # reference data computed with pyscf: - # pyscf_addrs, pyscf_signs = tn_addrs_signs(orbitals, electrons, excitation) + # pyscf_addrs, pyscf_signs = pyscf.ci.cisd.tn_addrs_signs(orbitals, electrons, excitation) # pyscf_state = pyscf.fci.cistring.addrs2str(orbitals, electrons, pyscf_addrs) # pyscf_state, pyscf_signs ( @@ -852,9 +852,9 @@ def test_excitations(electrons, orbitals, singles_ref, doubles_ref): ), ], ) -def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref): - r"""Test if the _excitated_states function returns correct states and signs.""" - states, signs = qchem.convert._excitated_states(electrons, orbitals, excitation) +def test_excited_configurations(electrons, orbitals, excitation, states_ref, signs_ref): + r"""Test if the _excited_configurations function returns correct states and signs.""" + states, signs = qchem.convert._excited_configurations(electrons, orbitals, excitation) assert np.allclose(states, states_ref) assert np.allclose(signs, signs_ref) From 9be65fa9305bfbce0df666d3f80ec1968fcc6a4b Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 15:37:50 -0400 Subject: [PATCH 21/41] correct test --- pennylane/qchem/convert.py | 16 ++++++++-------- tests/qchem/of_tests/test_convert.py | 16 ++++++---------- 2 files changed, 14 insertions(+), 18 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 91b4827f85c..654b3516b68 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -15,14 +15,14 @@ This module contains the functions for converting an external operator to a Pennylane operator. """ import warnings +from itertools import product # pylint: disable=import-outside-toplevel import pennylane as qml from pennylane import numpy as np -from pennylane.wires import Wires from pennylane.operation import Tensor, active_new_opmath from pennylane.pauli import pauli_sentence -from itertools import product +from pennylane.wires import Wires def _process_wires(wires, n_wires=None): @@ -525,15 +525,15 @@ def _ucisd_state(cisd_solver, tol=1e-15): dict_fcimatr = dict(zip(list(zip([ref_a], [ref_b])), [c0])) # alpha -> alpha excitations - c1a_configs, c1a_signs = _excitated_states(nelec_a, norb, 1) + c1a_configs, c1a_signs = _excited_configurations(nelec_a, norb, 1) dict_fcimatr.update(dict(zip(list(zip(c1a_configs, [ref_b] * size_a)), c1a * c1a_signs))) # beta -> beta excitations - c1b_configs, c1b_signs = _excitated_states(nelec_b, norb, 1) + c1b_configs, c1b_signs = _excited_configurations(nelec_b, norb, 1) dict_fcimatr.update(dict(zip(list(zip(c1b_configs, [ref_a] * size_b)), c1b * c1b_signs))) # alpha, alpha -> alpha, alpha excitations - c2aa_configs, c2aa_signs = _excitated_states(nelec_a, norb, 2) + c2aa_configs, c2aa_signs = _excited_configurations(nelec_a, norb, 2) dict_fcimatr.update(dict(zip(list(zip(c2aa_configs, [ref_b] * size_aa)), c2aa * c2aa_signs))) # alpha, beta -> alpha, beta excitations @@ -543,7 +543,7 @@ def _ucisd_state(cisd_solver, tol=1e-15): ) # beta, beta -> beta, beta excitations - c2bb_configs, c2bb_signs = _excitated_states(nelec_b, norb, 2) + c2bb_configs, c2bb_signs = _excited_configurations(nelec_b, norb, 2) dict_fcimatr.update(dict(zip(list(zip([ref_a] * size_bb, c2bb_configs)), c2bb * c2bb_signs))) # filter based on tolerance cutoff @@ -552,7 +552,7 @@ def _ucisd_state(cisd_solver, tol=1e-15): return dict_fcimatr -def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): +def cisd_state(cisd_solver, hftype, tol=1e-15): r"""Construct a wavefunction from PySCF output. Args: @@ -583,7 +583,7 @@ def cisd_state(cisd_solver, hftype, state=0, tol=1e-15): """ if hftype == "uhf": - wf_dict = _ucisd_state(cisd_solver, state=state, tol=tol) + wf_dict = _ucisd_state(cisd_solver, tol=tol) else: raise ValueError("Only restricted, 'rhf', and unrestricted, 'uhf', are supported.") wf = wfdict_to_statevector(wf_dict, cisd_solver.mol.nao) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 1e27f09c56c..0cb4ca0cc9a 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -24,7 +24,7 @@ import pennylane as qml from pennylane import numpy as np from pennylane import qchem -from pennylane.operation import enable_new_opmath, disable_new_opmath +from pennylane.operation import disable_new_opmath, enable_new_opmath openfermion = pytest.importorskip("openfermion") openfermionpyscf = pytest.importorskip("openfermionpyscf") @@ -860,7 +860,7 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig @pytest.mark.parametrize( - ("atom", "basis", "hftype", "state", "tol", "wf_ref"), + ("molecule", "basis", "tol", "wf_ref"), [ # reference data computed with pyscf -- identify the FCI matrix entries # with the correct Fock occupation vector entries @@ -872,16 +872,12 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig ( [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], "sto6g", - "uhf", - 0, 1e-1, {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159}, ), ( [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], "cc-pvdz", - "uhf", - 0, 4e-2, { (1, 1): 0.9919704801055201, @@ -894,16 +890,16 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig ), ], ) -def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): +def test_ucisd_state(molecule, basis, tol, wf_ref): r"""Test the UCISD wavefunction constructor.""" from pyscf import gto, scf, ci - mol = gto.M(atom=atom, basis=basis) + mol = gto.M(atom=molecule, basis=basis) myhf = scf.UHF(mol).run() - myci = ci.UCISD(myhf).run(state=state) + myci = ci.UCISD(myhf).run() - wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) + wf_cisd = qchem.convert._ucisd_state(myci, tol=tol) assert wf_cisd.keys() == wf_ref.keys() assert np.allclose(abs(np.array(list(wf_cisd.values()))), abs(np.array(list(wf_ref.values())))) From a4eceafbfb5f4eb1f54076e699b69642b3b4e355 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 15:39:48 -0400 Subject: [PATCH 22/41] correct import --- tests/qchem/of_tests/test_convert.py | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 0cb4ca0cc9a..d33f9a9bb8a 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -28,7 +28,7 @@ openfermion = pytest.importorskip("openfermion") openfermionpyscf = pytest.importorskip("openfermionpyscf") - +pyscf = pytest.importorskip("pyscf") pauli_ops_and_prod = (qml.PauliX, qml.PauliY, qml.PauliZ, qml.Identity, qml.ops.Prod) pauli_ops_and_tensor = (qml.PauliX, qml.PauliY, qml.PauliZ, qml.Identity, qml.operation.Tensor) @@ -893,11 +893,9 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig def test_ucisd_state(molecule, basis, tol, wf_ref): r"""Test the UCISD wavefunction constructor.""" - from pyscf import gto, scf, ci - - mol = gto.M(atom=molecule, basis=basis) - myhf = scf.UHF(mol).run() - myci = ci.UCISD(myhf).run() + mol = pyscf.gto.M(atom=molecule, basis=basis) + myhf = pyscf.scf.UHF(mol).run() + myci = pyscf.ci.UCISD(myhf).run() wf_cisd = qchem.convert._ucisd_state(myci, tol=tol) From 99b0fee2f2e4a7376efd0ce1240d6ae3e295785d Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 16:02:50 -0400 Subject: [PATCH 23/41] add test for statevector --- tests/qchem/of_tests/test_convert.py | 39 +++++++++++++++++++++++++++- 1 file changed, 38 insertions(+), 1 deletion(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index d33f9a9bb8a..39a131214ac 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -891,7 +891,7 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig ], ) def test_ucisd_state(molecule, basis, tol, wf_ref): - r"""Test the UCISD wavefunction constructor.""" + r"""Test that _ucisd_state returns the correct wavefunction.""" mol = pyscf.gto.M(atom=molecule, basis=basis) myhf = pyscf.scf.UHF(mol).run() @@ -901,3 +901,40 @@ def test_ucisd_state(molecule, basis, tol, wf_ref): assert wf_cisd.keys() == wf_ref.keys() assert np.allclose(abs(np.array(list(wf_cisd.values()))), abs(np.array(list(wf_ref.values())))) + + +@pytest.mark.parametrize( + ("wf_dict", "n_orbitals", "statevector"), + [ + ( # -0.99 |1100> + 0.11 |0011> + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159}, + 2, + np.array( + [ + 0.0, + 0.0, + 0.0, + 0.1066467, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + -0.99429698, + 0.0, + 0.0, + 0.0, + ] + ), + ), + ], +) +def test_wfdict_to_statevector(wf_dict, n_orbitals, statevector): + r"""Test that wfdict_to_statevector returns the correct statevector.""" + + wf_comp = qchem.convert.wfdict_to_statevector(wf_dict, n_orbitals) + + assert np.allclose(wf_comp, statevector) From 41238e618fe1730ca542e51968f2fde835bf4e10 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 16:14:58 -0400 Subject: [PATCH 24/41] add test for cisd --- tests/qchem/of_tests/test_convert.py | 59 ++++++++++++++++++++++------ 1 file changed, 46 insertions(+), 13 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 39a131214ac..6d79b9b1916 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -862,13 +862,6 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig @pytest.mark.parametrize( ("molecule", "basis", "tol", "wf_ref"), [ - # reference data computed with pyscf -- identify the FCI matrix entries - # with the correct Fock occupation vector entries - # mycasci = mcscf.UCASCI(myhf, norbs, (nelec_a, nelec_b)).run() - # sparse_cascimatr = coo_matrix(mycasci.ci, shape=np.shape(cascivec), dtype=float) - # row, col, dat = sparse_cascimatr.row, sparse_cascimatr.col, sparse_cascimatr.data - # strs_row, strs_col = addrs2str(norbs, nelec_a, row), addrs2str(norbs, nelec_b, col) - # wf_ref = dict(zip(list(zip(strs_row, strs_col)), dat)) ( [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], "sto6g", @@ -904,10 +897,10 @@ def test_ucisd_state(molecule, basis, tol, wf_ref): @pytest.mark.parametrize( - ("wf_dict", "n_orbitals", "statevector"), + ("wf_dict", "n_orbitals", "wf_ref"), [ ( # -0.99 |1100> + 0.11 |0011> - {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159}, + {(1, 1): 0.9942969785398778, (2, 2): 0.10664669927602159}, 2, np.array( [ @@ -923,7 +916,7 @@ def test_ucisd_state(molecule, basis, tol, wf_ref): 0.0, 0.0, 0.0, - -0.99429698, + 0.99429698, 0.0, 0.0, 0.0, @@ -932,9 +925,49 @@ def test_ucisd_state(molecule, basis, tol, wf_ref): ), ], ) -def test_wfdict_to_statevector(wf_dict, n_orbitals, statevector): +def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): r"""Test that wfdict_to_statevector returns the correct statevector.""" - wf_comp = qchem.convert.wfdict_to_statevector(wf_dict, n_orbitals) + assert np.allclose(wf_comp, wf_ref) + + +@pytest.mark.parametrize( + ("molecule", "basis", "hftype", "wf_ref"), + [ + ( + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "sto6g", + "uhf", + np.array( + [ + 0.0, + 0.0, + 0.0, + 0.1066467, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.99429698, + 0.0, + 0.0, + 0.0, + ] + ), + ), + ], +) +def test_cisd_state(molecule, basis, hftype, wf_ref): + r"""Test that cisd_state returns the correct wavefunction.""" + + mol = pyscf.gto.M(atom=molecule, basis=basis) + myhf = pyscf.scf.UHF(mol).run() + myci = pyscf.ci.UCISD(myhf).run() + + wf_comp = qchem.convert.cisd_state(myci, hftype) - assert np.allclose(wf_comp, statevector) + assert np.allclose(abs(wf_comp), wf_ref) From 0b8566bfb7a5fceecae221c74a8f906da2fc64ca Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 18:06:30 -0400 Subject: [PATCH 25/41] add error test and fix bug --- pennylane/qchem/convert.py | 6 ++++-- tests/qchem/of_tests/test_convert.py | 6 ++++++ 2 files changed, 10 insertions(+), 2 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 654b3516b68..c676f399999 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -530,7 +530,7 @@ def _ucisd_state(cisd_solver, tol=1e-15): # beta -> beta excitations c1b_configs, c1b_signs = _excited_configurations(nelec_b, norb, 1) - dict_fcimatr.update(dict(zip(list(zip(c1b_configs, [ref_a] * size_b)), c1b * c1b_signs))) + dict_fcimatr.update(dict(zip(list(zip([ref_a] * size_b, c1b_configs)), c1b * c1b_signs))) # alpha, alpha -> alpha, alpha excitations c2aa_configs, c2aa_signs = _excited_configurations(nelec_a, norb, 2) @@ -585,7 +585,9 @@ def cisd_state(cisd_solver, hftype, tol=1e-15): if hftype == "uhf": wf_dict = _ucisd_state(cisd_solver, tol=tol) else: - raise ValueError("Only restricted, 'rhf', and unrestricted, 'uhf', are supported.") + raise ValueError( + "The supported hftype options are 'rhf' for restricted," " and 'uhf' for unrestricted." + ) wf = wfdict_to_statevector(wf_dict, cisd_solver.mol.nao) return wf diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 6d79b9b1916..6b128561b86 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -971,3 +971,9 @@ def test_cisd_state(molecule, basis, hftype, wf_ref): wf_comp = qchem.convert.cisd_state(myci, hftype) assert np.allclose(abs(wf_comp), wf_ref) + + +def test_cisd_state_error(): + r"""Test that an error is raised if a wrong/not-supported hftype symbol is entered.""" + with pytest.raises(ValueError, match="The supported hftype options are"): + _ = qchem.convert.cisd_state("myci", "hf") From ae363b3c9b40750b85cebda4456962806e207cf2 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 18:12:29 -0400 Subject: [PATCH 26/41] add correct input args to test --- tests/qchem/of_tests/test_convert.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 6b128561b86..cdc0eaec40d 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -975,5 +975,9 @@ def test_cisd_state(molecule, basis, hftype, wf_ref): def test_cisd_state_error(): r"""Test that an error is raised if a wrong/not-supported hftype symbol is entered.""" + + myci = pyscf.ci.UCISD + hftype = "wrongtype" + with pytest.raises(ValueError, match="The supported hftype options are"): - _ = qchem.convert.cisd_state("myci", "hf") + _ = qchem.convert.cisd_state(myci, hftype) From 0ddab7bd5281526ae1315dcbee5eaee99ff3c836 Mon Sep 17 00:00:00 2001 From: soranjh Date: Wed, 2 Aug 2023 18:16:22 -0400 Subject: [PATCH 27/41] fix typo --- pennylane/qchem/convert.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index c676f399999..3e17fb710ea 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -586,7 +586,7 @@ def cisd_state(cisd_solver, hftype, tol=1e-15): wf_dict = _ucisd_state(cisd_solver, tol=tol) else: raise ValueError( - "The supported hftype options are 'rhf' for restricted," " and 'uhf' for unrestricted." + "The supported hftype options are 'rhf' for restricted and 'uhf' for unrestricted." ) wf = wfdict_to_statevector(wf_dict, cisd_solver.mol.nao) From 0ceffc85a1d56a3b96ddfd518f21c8e285886015 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Wed, 2 Aug 2023 19:04:14 -0400 Subject: [PATCH 28/41] updating tests with symm argument --- tests/qchem/of_tests/test_convert.py | 135 +++++++++++++++++++++------ 1 file changed, 106 insertions(+), 29 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index c7c8c1735ff..c30276950b9 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -28,7 +28,7 @@ openfermion = pytest.importorskip("openfermion") openfermionpyscf = pytest.importorskip("openfermionpyscf") - +pyscf = pytest.importorskip("pyscf") pauli_ops_and_prod = (qml.PauliX, qml.PauliY, qml.PauliZ, qml.Identity, qml.ops.Prod) pauli_ops_and_tensor = (qml.PauliX, qml.PauliY, qml.PauliZ, qml.Identity, qml.operation.Tensor) @@ -859,43 +859,120 @@ def test_excitated_states(electrons, orbitals, excitation, states_ref, signs_ref assert np.allclose(signs, signs_ref) @pytest.mark.parametrize( - ("atom", "basis", "hftype", "state", "tol", "wf_ref"), + ("molecule", "basis", "symm", "tol", "wf_ref"), [ - # reference data computed with pyscf -- identify the FCI matrix entries - # with the correct Fock occupation vector entries - # mycasci = mcscf.UCASCI(myhf, norbs, (nelec_a, nelec_b)).run() - # sparse_cascimatr = coo_matrix(mycasci.ci, shape=np.shape(cascivec), dtype=float) - # row, col, dat = sparse_cascimatr.row, sparse_cascimatr.col, sparse_cascimatr.data - # strs_row, strs_col = addrs2str(norbs, nelec_a, row), addrs2str(norbs, nelec_b, col) - # wf_ref = dict(zip(list(zip(strs_row, strs_col)), dat)) ( - [['H', (0, 0, 0)], ['H', (0,0,0.71)]], - 'sto6g', - 'uhf', - 0, + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "sto6g", + "d2h", 1e-1, - {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} + {(1, 1): 0.9942969785398776, (2, 2): -0.10664669927602176}, ), ( - [['H', (0, 0, 0)], ['H', (0, 0, 0.71)]], - 'cc-pvdz', - 'uhf', - 0, + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "cc-pvdz", + "d2h", 4e-2, - {(1, 1): 0.9919704801055201, (2, 2): 0.04853035090478132, (2, 8): -0.04452332640264183, \ - (4, 4): -0.05003594588609278, (8, 2): 0.044523334555907366, (8, 8): -0.05226230845395026}), + { + (1, 1): 0.9919704795977625, + (2, 2): -0.048530356564387034, + (2, 8): 0.0445233308500785, + (4, 4): -0.05003594568491194, + (8, 2): 0.04452333085007853, + (8, 8): -0.05226230322043741, + (16, 16): -0.0404759737476627, + (32, 32): -0.0404759737476627 + }, + ), + ], +) +def test_ucisd_state(molecule, basis, symm, tol, wf_ref): + r"""Test that _ucisd_state returns the correct wavefunction.""" + + mol = pyscf.gto.M(atom=molecule, basis=basis, symmetry=symm) + myhf = pyscf.scf.UHF(mol).run() + myci = pyscf.ci.UCISD(myhf).run() + + wf_cisd = qchem.convert._ucisd_state(myci, tol=tol) + + assert wf_cisd.keys() == wf_ref.keys() + assert np.allclose(abs(np.array(list(wf_cisd.values()))), abs(np.array(list(wf_ref.values())))) + + +@pytest.mark.parametrize( + ("wf_dict", "n_orbitals", "wf_ref"), + [ + ( # -0.99 |1100> + 0.11 |0011> + {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602179}, + 2, + np.array( + [ + 0.0, + 0.0, + 0.0, + 0.1066467, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + -0.99429698, + 0.0, + 0.0, + 0.0, + ] + ), + ), ], ) -def test_private_ucisd_state(atom, basis, hftype, state, tol, wf_ref): - r"""Test the UCISD wavefunction constructor.""" +def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): + r"""Test that wfdict_to_statevector returns the correct statevector.""" + wf_comp = qchem.convert.wfdict_to_statevector(wf_dict, n_orbitals) + assert np.allclose(wf_comp, wf_ref) + - from pyscf import gto, scf, ci +@pytest.mark.parametrize( + ("molecule", "basis", "symm", "hftype", "wf_ref"), + [ + ( + [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], + "sto6g", + "d2h", + "uhf", + np.array( + [ + 0.0, + 0.0, + 0.0, + -0.1066467, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.99429698, + 0.0, + 0.0, + 0.0, + ] + ), + ), + ], +) +def test_cisd_state(molecule, basis, symm, hftype, wf_ref): + r"""Test that cisd_state returns the correct wavefunction.""" - mol = gto.M(atom=atom, basis=basis) - myhf = scf.UHF(mol).run() - myci = ci.UCISD(myhf).run(state=state) + mol = pyscf.gto.M(atom=molecule, basis=basis, symmetry=symm) + myhf = pyscf.scf.UHF(mol).run() + myci = pyscf.ci.UCISD(myhf).run() - wf_cisd = qchem.convert.cisd_state(myci, hftype=hftype, state=state, tol=tol) + wf_comp = qchem.convert.cisd_state(myci, hftype) - assert { key: round(val, 4) for (key, val) in wf_cisd.items() } == \ - { key: round(val, 4) for (key, val) in wf_ref.items() } \ No newline at end of file + # overall sign could be +/-1 different + assert np.allclose(wf_comp, wf_ref) or np.allclose(wf_comp, -wf_ref) \ No newline at end of file From eee04106820200757e0fff80611c87f587ec50b1 Mon Sep 17 00:00:00 2001 From: Stepan Fomichev Date: Wed, 2 Aug 2023 19:08:37 -0400 Subject: [PATCH 29/41] reformatted with black --- tests/qchem/of_tests/test_convert.py | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index a72a6a9f239..844a75b33a7 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -875,14 +875,14 @@ def test_excited_configurations(electrons, orbitals, excitation, states_ref, sig "d2h", 4e-2, { - (1, 1): 0.9919704795977625, - (2, 2): -0.048530356564387034, - (2, 8): 0.0445233308500785, - (4, 4): -0.05003594568491194, - (8, 2): 0.04452333085007853, - (8, 8): -0.05226230322043741, - (16, 16): -0.0404759737476627, - (32, 32): -0.0404759737476627 + (1, 1): 0.9919704795977625, + (2, 2): -0.048530356564387034, + (2, 8): 0.0445233308500785, + (4, 4): -0.05003594568491194, + (8, 2): 0.04452333085007853, + (8, 8): -0.05226230322043741, + (16, 16): -0.0404759737476627, + (32, 32): -0.0404759737476627, }, ), ], @@ -978,6 +978,7 @@ def test_cisd_state(molecule, basis, symm, hftype, wf_ref): # overall sign could be +/-1 different assert np.allclose(wf_comp, wf_ref) or np.allclose(wf_comp, -wf_ref) + def test_cisd_state_error(): r"""Test that an error is raised if a wrong/not-supported hftype symbol is entered.""" From 38011470ca1090490b837d9ac8d3ff0461f45040 Mon Sep 17 00:00:00 2001 From: soranjh Date: Thu, 3 Aug 2023 09:43:10 -0400 Subject: [PATCH 30/41] make code compact --- pennylane/qchem/convert.py | 12 +++++------- 1 file changed, 5 insertions(+), 7 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 3e17fb710ea..67cb6aad944 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -594,7 +594,7 @@ def cisd_state(cisd_solver, hftype, tol=1e-15): def wfdict_to_statevector(wf_dict, norbs): - r"""Convert a wavefunction in sparce dictionary format to a Pennylane's statevector. + r"""Convert a wavefunction in sparce dictionary format to a PennyLane statevector. Args: wf_dict (wf_dict format): dictionary with keys (int_a,int_b) being integers whose binary representation @@ -611,12 +611,10 @@ def wfdict_to_statevector(wf_dict, norbs): statevector = np.zeros(2 ** (2 * norbs), dtype=complex) for (int_a, int_b), coeff in wf_dict.items(): - bin_a, bin_b = bin(int_a)[2:], bin(int_b)[2:] - bin_a, bin_b = bin_a[::-1], bin_b[::-1] - bin_tot = "".join(i + j for i, j in zip(bin_a, bin_b)) - bin_tot_full = bin_tot + "0" * (2 * norbs - len(bin_tot)) - idx = int(bin_tot_full, 2) - statevector[idx] += coeff + bin_a, bin_b = bin(int_a)[2:][::-1], bin(int_b)[2:][::-1] + bin_ab = "".join(i + j for i, j in zip(bin_a, bin_b)) + bin_ab_full = bin_ab + "0" * (2 * norbs - len(bin_ab)) + statevector[int(bin_ab_full, 2)] += coeff statevector = statevector / np.sqrt(np.sum(statevector**2)) From c8e69837f4ba4caa7b405f8dbd327d96bdb38c37 Mon Sep 17 00:00:00 2001 From: soranjh Date: Thu, 3 Aug 2023 10:33:55 -0400 Subject: [PATCH 31/41] rename cisd_state to import_state --- pennylane/qchem/convert.py | 20 +++++++++----------- tests/qchem/of_tests/test_convert.py | 22 +++++++++++----------- 2 files changed, 20 insertions(+), 22 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 67cb6aad944..208590c29e7 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -552,15 +552,12 @@ def _ucisd_state(cisd_solver, tol=1e-15): return dict_fcimatr -def cisd_state(cisd_solver, hftype, tol=1e-15): +def import_state(solver, method, tol=1e-15): r"""Construct a wavefunction from PySCF output. Args: - cisd_solver (PySCF CISD or UCISD Class instance): the class object representing the - CISD / UCISD calculation in PySCF - - hftype (str): keyword specifying whether restricted or unrestricted CISD was performed. - The options are 'rhf' for restricted, and 'uhf' for unrestricted. + solver: external wavefunction object that will be converted + method (str): keyword specifying the method of calculation tol (float): the tolerance to which the wavefunction is being built -- Slater determinants with coefficients smaller than this are discarded. Default @@ -582,18 +579,19 @@ def cisd_state(cisd_solver, hftype, tol=1e-15): {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} """ - if hftype == "uhf": - wf_dict = _ucisd_state(cisd_solver, tol=tol) + if method == "ucisd": + wf_dict = _ucisd_state(solver, tol=tol) else: raise ValueError( - "The supported hftype options are 'rhf' for restricted and 'uhf' for unrestricted." + "The supported method options are 'rcisd' for restricted and 'ucisd' for unrestricted" + " CISD calculations." ) - wf = wfdict_to_statevector(wf_dict, cisd_solver.mol.nao) + wf = _wfdict_to_statevector(wf_dict, solver.mol.nao) return wf -def wfdict_to_statevector(wf_dict, norbs): +def _wfdict_to_statevector(wf_dict, norbs): r"""Convert a wavefunction in sparce dictionary format to a PennyLane statevector. Args: diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 844a75b33a7..8b678c5f7f6 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -930,19 +930,19 @@ def test_ucisd_state(molecule, basis, symm, tol, wf_ref): ], ) def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): - r"""Test that wfdict_to_statevector returns the correct statevector.""" - wf_comp = qchem.convert.wfdict_to_statevector(wf_dict, n_orbitals) + r"""Test that _wfdict_to_statevector returns the correct statevector.""" + wf_comp = qchem.convert._wfdict_to_statevector(wf_dict, n_orbitals) assert np.allclose(wf_comp, wf_ref) @pytest.mark.parametrize( - ("molecule", "basis", "symm", "hftype", "wf_ref"), + ("molecule", "basis", "symm", "method", "wf_ref"), [ ( [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], "sto6g", "d2h", - "uhf", + "ucisd", np.array( [ 0.0, @@ -966,24 +966,24 @@ def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): ), ], ) -def test_cisd_state(molecule, basis, symm, hftype, wf_ref): +def test_import_state(molecule, basis, symm, method, wf_ref): r"""Test that cisd_state returns the correct wavefunction.""" mol = pyscf.gto.M(atom=molecule, basis=basis, symmetry=symm) myhf = pyscf.scf.UHF(mol).run() myci = pyscf.ci.UCISD(myhf).run() - wf_comp = qchem.convert.cisd_state(myci, hftype) + wf_comp = qchem.convert.import_state(myci, method) # overall sign could be +/-1 different assert np.allclose(wf_comp, wf_ref) or np.allclose(wf_comp, -wf_ref) -def test_cisd_state_error(): - r"""Test that an error is raised if a wrong/not-supported hftype symbol is entered.""" +def test_import_state_error(): + r"""Test that an error is raised if a wrong/not-supported method symbol is entered.""" myci = pyscf.ci.UCISD - hftype = "wrongtype" + method = "wrongmethod" - with pytest.raises(ValueError, match="The supported hftype options are"): - _ = qchem.convert.cisd_state(myci, hftype) + with pytest.raises(ValueError, match="The supported method options are"): + _ = qchem.convert.import_state(myci, method) From 699373d0477b321f26043448cf81e781811f450e Mon Sep 17 00:00:00 2001 From: soranjh Date: Thu, 3 Aug 2023 13:03:03 -0400 Subject: [PATCH 32/41] modify docs --- pennylane/qchem/convert.py | 65 ++++++++++++++-------------- tests/qchem/of_tests/test_convert.py | 2 +- 2 files changed, 33 insertions(+), 34 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 208590c29e7..9fa3c779086 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -410,16 +410,19 @@ def _excitations(electrons, orbitals): def _excited_configurations(electrons, orbitals, excitation): - r"""Generate excited states from a Hartree-Fock reference state. + r"""Generate excited configurations from a Hartree-Fock reference state. - This function generates excited states in the form of binary strings or integers representing - them, from a Hartree-Fock (HF) reference state. The HF state is assumed to be - :math:`|1 1 ...1 0 ... 0 0 \rangle` where the number of :math:`1` and :math:`0` elements - represents the number of occupied and unoccupied spin orbitals, respectively. The string - representation of the state is obtained by converting the occupation-number vector to a string, - e.g., ``111000`` represents :math:`|1 1 1 0 0 0 \rangle. The number representation of the state - is the integer corresonding to the binary number, e.g., ``111000`` is represented by - :math:`int('111000', 2) = 56`. + This function generates excited configurations in the form of integers representing a binary + string, e.g., :math:`|1 1 0 1 0 0 \rangle is represented by :math:`int('110100', 2) = 52`. + + The excited configurations are generated from a Hartree-Fock (HF) reference state. The HF state + is assumed to have the form :math:`|1 1 ...1 0 ... 0 0 \rangle` where the number of :math:`1` + and :math:`0` elements are the number of occupied and unoccupied spin orbitals, respectively. + The string representation of the state is obtained by converting the occupation-number vector to + a string, e.g., ``111000`` to represent :math:`|1 1 1 0 0 0 \rangle. + + Each excited configuration has a sign, :math:`+1` or :math:`-1`, that is obtained by reordering + the creation operators. Args: electrons (int): number of electrons @@ -427,7 +430,8 @@ def _excited_configurations(electrons, orbitals, excitation): excitation (int): number of excited electrons Returns: - tuple(array, array): arrays of excited states and signs obtained by reordering of the creation operators + tuple(array, array): arrays of excited configurations and signs obtained by reordering the + creation operators **Example** @@ -464,19 +468,18 @@ def _excited_configurations(electrons, orbitals, excitation): def _ucisd_state(cisd_solver, tol=1e-15): - r"""Construct a wavefunction from PySCF's `UCISD` Solver object. - - The wavefunction is represented as a dictionary where the keys are tuples representing a - configuration, which corresponds to a Slater determinant, and the values are the CI coefficient - corresponding to that Slater determinant. Each dictionary key is a tuple of two integers. The - binary representation of these integers correspond to a specific configuration, or Slater - determinant. The first number represents the configuration of alpha electrons and the second - number represents the beta electrons. For instance, the Hartree-Fock state - :math:`|1 1 0 0 \rangle` will be represented by the binary string `0011`. This string can be - splited to `01` and `01` for the alpha and beta electrons. The integer corresponding to `01` is - `1`. Then the dictionary representation of a state with `0.99` contribution of the Hartree-Fock - state and `0.01` contribution from the doubly-excited state will be - `{(1, 1): 0.99, (2, 2): 0.01}`. + r"""Construct a wavefunction from PySCF's `UCISD` solver object. + + The generated wavefunction is a dictionary where the keys represent a configuration, which + corresponds to a Slater determinant, and the values are the CI coefficients of the that Slater + determinant. Each dictionary key is a tuple of two integers. The binary representation of these + integers correspond to a specific configuration such that the first number represents the + configuration of the alpha electrons and the second number represents the configuration of the + beta electrons. For instance, the Hartree-Fock state :math:`|1 1 0 0 \rangle` will be + represented by the flipped binary string `0011`. This string can be splitted to `01` and `01` for + the alpha and beta electrons. The integer corresponding to `01` is `1`. Then the dictionary + representation of a state with `0.99` contribution of the Hartree-Fock state and `0.01` + contribution from the doubly-excited state will be `{(1, 1): 0.99, (2, 2): 0.01}`. Args: cisd_solver (PySCF UCISD Class instance): the class object representing the UCISD calculation in PySCF @@ -485,7 +488,7 @@ def _ucisd_state(cisd_solver, tol=1e-15): Returns: dict: Dictionary of the form `{(int_a, int_b) :coeff}`, with integers `int_a, int_b` - having binary represention corresponding to the Fock occupation vector in alpha and beta + having binary representation corresponding to the Fock occupation vector in alpha and beta spin sectors, respectively, and coeff being the CI coefficients of those configurations. **Example** @@ -553,20 +556,17 @@ def _ucisd_state(cisd_solver, tol=1e-15): def import_state(solver, method, tol=1e-15): - r"""Construct a wavefunction from PySCF output. + r"""Convert an external wavefunction to a PennyLane state vector. Args: solver: external wavefunction object that will be converted - method (str): keyword specifying the method of calculation + method (str): keyword specifying the calculation method tol (float): the tolerance to which the wavefunction is being built -- Slater - determinants with coefficients smaller than this are discarded. Default - is 1e-15 (all coefficients are included). + determinants with coefficients smaller than this tolerance are discarded. Returns: - dict: Dictionary of the form `{(int_a, int_b) :coeff}`, with integers `int_a, int_b` - having binary represention corresponding to the Fock occupation vector in alpha and beta - spin sectors, respectively, and coeff being the CI coefficients of those configurations. + statevector (array): normalized state vector of length 2**len(number_of_spinorbitals) **Example** @@ -602,10 +602,9 @@ def _wfdict_to_statevector(wf_dict, norbs): norbs (int): number of molecular orbitals in the system Returns: - statevector (list): normalized statevector of length 2^(number_of_spinorbitals), standard in Pennylane + statevector (array): normalized state vector of length 2^(number_of_spinorbitals) """ - statevector = np.zeros(2 ** (2 * norbs), dtype=complex) for (int_a, int_b), coeff in wf_dict.items(): diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 8b678c5f7f6..c332d3ede45 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -975,7 +975,7 @@ def test_import_state(molecule, basis, symm, method, wf_ref): wf_comp = qchem.convert.import_state(myci, method) - # overall sign could be +/-1 different + # overall sign could be different in each PySCF run assert np.allclose(wf_comp, wf_ref) or np.allclose(wf_comp, -wf_ref) From b1867362221d9d196b5537aef36b07b501fd6fd3 Mon Sep 17 00:00:00 2001 From: soranjh Date: Thu, 3 Aug 2023 17:30:39 -0400 Subject: [PATCH 33/41] add review comments --- pennylane/qchem/convert.py | 24 +++++++++++++++--------- tests/qchem/of_tests/test_convert.py | 25 +++++++++++++++++-------- 2 files changed, 32 insertions(+), 17 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 9fa3c779086..6986c5f4faa 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -391,7 +391,7 @@ def _excitations(electrons, orbitals): >>> electrons = 2 >>> orbitals = 4 >>> _excitations(electrons, orbitals) - ([[[0, 2]], [[0, 3]], [[1, 2]], [[1, 3]]], [[0, 1, 2, 3]]) + ([[0, 2], [0, 3], [1, 2], [1, 3]], [[0, 1, 2, 3]]) """ singles_p, singles_q = [], [] doubles_pq, doubles_rs = [], [] @@ -441,6 +441,11 @@ def _excited_configurations(electrons, orbitals, excitation): >>> _excitated_states(electrons, orbitals, excitation) (array([28, 26, 25]), array([ 1, -1, 1])) """ + if excitation not in [1, 2]: + raise ValueError( + "Only single (excitation = 1) and double (excitation = 2) excitations are supported." + ) + hf_state = qml.qchem.hf_state(electrons, orbitals) singles, doubles = _excitations(electrons, orbitals) states, signs = [], [] @@ -448,14 +453,14 @@ def _excited_configurations(electrons, orbitals, excitation): if excitation == 1: for s in singles: state = hf_state.copy() - state[s] = state[s[::-1]] + state[s] = state[[s[1], s[0]]] # apply single excitation states += [state] signs.append((-1) ** len(range(s[0], electrons - 1))) if excitation == 2: for d in doubles: state = hf_state.copy() - state[d] = state[[d[2], d[3], d[0], d[1]]] + state[d] = state[[d[2], d[3], d[0], d[1]]] # apply double excitation states += [state] order_pq = len(range(d[0], electrons - 1)) order_rs = len(range(d[1], electrons - 1)) @@ -560,7 +565,7 @@ def import_state(solver, method, tol=1e-15): Args: solver: external wavefunction object that will be converted - method (str): keyword specifying the calculation method + method (str): keyword specifying the calculation method, possible option is 'ucisd' tol (float): the tolerance to which the wavefunction is being built -- Slater determinants with coefficients smaller than this tolerance are discarded. @@ -583,8 +588,7 @@ def import_state(solver, method, tol=1e-15): wf_dict = _ucisd_state(solver, tol=tol) else: raise ValueError( - "The supported method options are 'rcisd' for restricted and 'ucisd' for unrestricted" - " CISD calculations." + "The supported method option is 'ucisd' for unrestricted CISD calculations." ) wf = _wfdict_to_statevector(wf_dict, solver.mol.nao) @@ -608,10 +612,12 @@ def _wfdict_to_statevector(wf_dict, norbs): statevector = np.zeros(2 ** (2 * norbs), dtype=complex) for (int_a, int_b), coeff in wf_dict.items(): - bin_a, bin_b = bin(int_a)[2:][::-1], bin(int_b)[2:][::-1] + bin_a = bin(int_a)[2:][::-1] + bin_b = bin(int_b)[2:][::-1] + bin_a += "0" * (norbs - len(bin_a)) + bin_b += "0" * (norbs - len(bin_b)) bin_ab = "".join(i + j for i, j in zip(bin_a, bin_b)) - bin_ab_full = bin_ab + "0" * (2 * norbs - len(bin_ab)) - statevector[int(bin_ab_full, 2)] += coeff + statevector[int(bin_ab, 2)] += coeff statevector = statevector / np.sqrt(np.sum(statevector**2)) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index c332d3ede45..0c12d0b1829 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -903,24 +903,24 @@ def test_ucisd_state(molecule, basis, symm, tol, wf_ref): @pytest.mark.parametrize( ("wf_dict", "n_orbitals", "wf_ref"), [ - ( # -0.99 |1100> + 0.11 |0011> - {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602179}, + ( + {(1, 1): 0.87006284, (1, 2): 0.3866946, (2, 1): 0.29002095, (2, 2): 0.09667365}, 2, np.array( [ 0.0, 0.0, 0.0, - 0.1066467, - 0.0, - 0.0, + 0.09667365, 0.0, 0.0, + 0.29002095, 0.0, 0.0, + 0.3866946, 0.0, 0.0, - -0.99429698, + 0.87006284, 0.0, 0.0, 0.0, @@ -932,6 +932,8 @@ def test_ucisd_state(molecule, basis, symm, tol, wf_ref): def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): r"""Test that _wfdict_to_statevector returns the correct statevector.""" wf_comp = qchem.convert._wfdict_to_statevector(wf_dict, n_orbitals) + print(wf_comp) + print(wf_ref) assert np.allclose(wf_comp, wf_ref) @@ -980,10 +982,17 @@ def test_import_state(molecule, basis, symm, method, wf_ref): def test_import_state_error(): - r"""Test that an error is raised if a wrong/not-supported method symbol is entered.""" + r"""Test that an error is raised by import_state if a wrong method symbol is entered.""" myci = pyscf.ci.UCISD method = "wrongmethod" - with pytest.raises(ValueError, match="The supported method options are"): + with pytest.raises(ValueError, match="The supported method"): _ = qchem.convert.import_state(myci, method) + + +@pytest.mark.parametrize(("excitation"), [-1, 0, 3]) +def test_excited_configurations_error(excitation): + r"""Test that an error is raised by _excited_configurations if a wrong excitation is entered.""" + with pytest.raises(ValueError, match="excitations are supported"): + _ = qchem.convert._excited_configurations(2, 4, excitation) From d0f3707f2b84143300f78f3858b3ed4804d489a5 Mon Sep 17 00:00:00 2001 From: soranjh Date: Thu, 3 Aug 2023 17:43:14 -0400 Subject: [PATCH 34/41] add error to doctsring --- pennylane/qchem/convert.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 6986c5f4faa..d04cf78d293 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -570,6 +570,9 @@ def import_state(solver, method, tol=1e-15): tol (float): the tolerance to which the wavefunction is being built -- Slater determinants with coefficients smaller than this tolerance are discarded. + Raises: + ValueError: if ``method`` is not supported + Returns: statevector (array): normalized state vector of length 2**len(number_of_spinorbitals) From 6f47d29ab9b6a256bc48351c83d4fa173289da7f Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 09:46:24 -0400 Subject: [PATCH 35/41] modify test --- tests/qchem/of_tests/test_convert.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 0c12d0b1829..1d2f3def4df 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -903,7 +903,7 @@ def test_ucisd_state(molecule, basis, symm, tol, wf_ref): @pytest.mark.parametrize( ("wf_dict", "n_orbitals", "wf_ref"), [ - ( + ( # 0.87006284 |1100> + 0.3866946 |1001> + 0.29002095 |0110> + 0.09667365 |0011> {(1, 1): 0.87006284, (1, 2): 0.3866946, (2, 1): 0.29002095, (2, 2): 0.09667365}, 2, np.array( From 188104ab539bf64819c76368895ef16b0d4cbbc3 Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 10:10:03 -0400 Subject: [PATCH 36/41] add more tests --- tests/qchem/of_tests/test_convert.py | 41 +++++++++++----------------- 1 file changed, 16 insertions(+), 25 deletions(-) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index 1d2f3def4df..dcd7a83098d 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -901,39 +901,30 @@ def test_ucisd_state(molecule, basis, symm, tol, wf_ref): @pytest.mark.parametrize( - ("wf_dict", "n_orbitals", "wf_ref"), - [ + ("wf_dict", "n_orbitals", "string_ref", "coeff_ref"), + [ # reference data were obtained manually ( # 0.87006284 |1100> + 0.3866946 |1001> + 0.29002095 |0110> + 0.09667365 |0011> {(1, 1): 0.87006284, (1, 2): 0.3866946, (2, 1): 0.29002095, (2, 2): 0.09667365}, 2, - np.array( - [ - 0.0, - 0.0, - 0.0, - 0.09667365, - 0.0, - 0.0, - 0.29002095, - 0.0, - 0.0, - 0.3866946, - 0.0, - 0.0, - 0.87006284, - 0.0, - 0.0, - 0.0, - ] - ), + ["1100", "1001", "0110", "0011"], + [0.87006284, 0.3866946, 0.29002095, 0.09667365], + ), + ( # 0.80448616 |110000> + 0.53976564 |001100> + 0.22350293 |000011> + 0.10724511 |100100> + {(1, 1): 0.80448616, (2, 2): 0.53976564, (4, 4): 0.22350293, (1, 2): 0.10724511}, + 3, + ["110000", "001100", "000011", "100100"], + [0.80448616, 0.53976564, 0.22350293, 0.10724511], ), ], ) -def test_wfdict_to_statevector(wf_dict, n_orbitals, wf_ref): +def test_wfdict_to_statevector(wf_dict, n_orbitals, string_ref, coeff_ref): r"""Test that _wfdict_to_statevector returns the correct statevector.""" + wf_ref = np.zeros(2 ** (n_orbitals * 2)) + idx_nonzero = [int(s, 2) for s in string_ref] + wf_ref[idx_nonzero] = coeff_ref + wf_comp = qchem.convert._wfdict_to_statevector(wf_dict, n_orbitals) - print(wf_comp) - print(wf_ref) + assert np.allclose(wf_comp, wf_ref) From c4020cc1a3ff2c458e1766e1b8f5d974a9604b13 Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 10:16:27 -0400 Subject: [PATCH 37/41] update changelog --- doc/releases/changelog-dev.md | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/doc/releases/changelog-dev.md b/doc/releases/changelog-dev.md index 8c40590fefc..ce09034bb36 100644 --- a/doc/releases/changelog-dev.md +++ b/doc/releases/changelog-dev.md @@ -37,6 +37,9 @@ def circuit(): >>> circuit(shots=1) array([False, False]) +* Functions added to convert wavefunctions obtained from `PySCF` to a state vector. + [(#4427)](https://github.com/PennyLaneAI/pennylane/pull/4427) +

Improvements 🛠

* Transform Programs, `qml.transforms.core.TransformProgram`, can now be called on a batch of circuits @@ -278,6 +281,7 @@ array([False, False]) This release contains contributions from (in alphabetical order): Isaac De Vlugt, +Stepan Fomichev, Lillian M. A. Frederiksen, Soran Jahangiri, Edward Jiang, From 6c5219af3a540ce273d7580aa43bd39852545adb Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 13:46:57 -0400 Subject: [PATCH 38/41] modify docstring and add import --- pennylane/qchem/__init__.py | 2 +- pennylane/qchem/convert.py | 17 ++++++++--------- 2 files changed, 9 insertions(+), 10 deletions(-) diff --git a/pennylane/qchem/__init__.py b/pennylane/qchem/__init__.py index fe8ce45927e..9af19dedd25 100644 --- a/pennylane/qchem/__init__.py +++ b/pennylane/qchem/__init__.py @@ -25,7 +25,7 @@ ) from .basis_data import load_basisset from .basis_set import BasisFunction, atom_basis_data, mol_basis_data -from .convert import import_operator +from .convert import import_operator, import_state from .dipole import dipole_integrals, fermionic_dipole, dipole_moment from .factorization import basis_rotation, factorize from .hamiltonian import electron_integrals, fermionic_hamiltonian, diff_hamiltonian diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index d04cf78d293..c244b83520a 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -566,9 +566,7 @@ def import_state(solver, method, tol=1e-15): Args: solver: external wavefunction object that will be converted method (str): keyword specifying the calculation method, possible option is 'ucisd' - - tol (float): the tolerance to which the wavefunction is being built -- Slater - determinants with coefficients smaller than this tolerance are discarded. + tol (float): the tolerance for discarding Slater determinants with small coefficients Raises: ValueError: if ``method`` is not supported @@ -599,14 +597,15 @@ def import_state(solver, method, tol=1e-15): def _wfdict_to_statevector(wf_dict, norbs): - r"""Convert a wavefunction in sparce dictionary format to a PennyLane statevector. + r"""Convert a wavefunction in sparse dictionary format to a PennyLane statevector. - Args: - wf_dict (wf_dict format): dictionary with keys (int_a,int_b) being integers whose binary representation - shows the Fock occupation vector for alpha and beta electrons; and values being the CI - coefficients. + In the sparse dictionary format, the keys (int_a, int_b) are integers whose binary + representation shows the Fock occupation vector for alpha and beta electrons and values are the + CI coefficients. - norbs (int): number of molecular orbitals in the system + Args: + wf_dict (dict): the sparse dictionary format of a wavefunction + norbs (int): number of molecular orbitals Returns: statevector (array): normalized state vector of length 2^(number_of_spinorbitals) From 9465660792ff2634276c9a0a215e75753aed681c Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 14:17:26 -0400 Subject: [PATCH 39/41] add pyscf import --- pennylane/qchem/convert.py | 24 ++++++++++++---------- tests/qchem/of_tests/test_convert.py | 30 +++++++++++++++++++--------- 2 files changed, 34 insertions(+), 20 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index c244b83520a..e5da9e389f0 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -435,11 +435,11 @@ def _excited_configurations(electrons, orbitals, excitation): **Example** - >>> electrons = 2 + >>> electrons = 3 >>> orbitals = 5 >>> excitation = 2 - >>> _excitated_states(electrons, orbitals, excitation) - (array([28, 26, 25]), array([ 1, -1, 1])) + >>> _excited_configurations(electrons, orbitals, excitation) + ([28, 26, 25], [1, -1, 1]) """ if excitation not in [1, 2]: raise ValueError( @@ -502,7 +502,7 @@ def _ucisd_state(cisd_solver, tol=1e-15): >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g') >>> myhf = scf.UHF(mol).run() >>> myci = ci.UCISD(myhf).run() - >>> wf_cisd = cisd_state(myci, tol=1e-1) + >>> wf_cisd = _ucisd_state(myci, tol=1e-1) >>> print(wf_cisd) {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} """ @@ -560,12 +560,11 @@ def _ucisd_state(cisd_solver, tol=1e-15): return dict_fcimatr -def import_state(solver, method, tol=1e-15): +def import_state(solver, tol=1e-15): r"""Convert an external wavefunction to a PennyLane state vector. Args: solver: external wavefunction object that will be converted - method (str): keyword specifying the calculation method, possible option is 'ucisd' tol (float): the tolerance for discarding Slater determinants with small coefficients Raises: @@ -584,13 +583,17 @@ def import_state(solver, method, tol=1e-15): >>> print(wf_cisd) {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} """ + try: + import pyscf + except ImportError as Error: + raise ImportError( + "This feature requires pyscf. " "It can be installed with: pip install pyscf" + ) from Error - if method == "ucisd": + if isinstance(solver, pyscf.ci.ucisd.UCISD): wf_dict = _ucisd_state(solver, tol=tol) else: - raise ValueError( - "The supported method option is 'ucisd' for unrestricted CISD calculations." - ) + raise ValueError("The supported option is 'ucisd' for unrestricted CISD calculations.") wf = _wfdict_to_statevector(wf_dict, solver.mol.nao) return wf @@ -609,7 +612,6 @@ def _wfdict_to_statevector(wf_dict, norbs): Returns: statevector (array): normalized state vector of length 2^(number_of_spinorbitals) - """ statevector = np.zeros(2 ** (2 * norbs), dtype=complex) diff --git a/tests/qchem/of_tests/test_convert.py b/tests/qchem/of_tests/test_convert.py index dcd7a83098d..ccc6c7e05ab 100644 --- a/tests/qchem/of_tests/test_convert.py +++ b/tests/qchem/of_tests/test_convert.py @@ -929,13 +929,12 @@ def test_wfdict_to_statevector(wf_dict, n_orbitals, string_ref, coeff_ref): @pytest.mark.parametrize( - ("molecule", "basis", "symm", "method", "wf_ref"), + ("molecule", "basis", "symm", "wf_ref"), [ ( [["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], "sto6g", "d2h", - "ucisd", np.array( [ 0.0, @@ -959,27 +958,26 @@ def test_wfdict_to_statevector(wf_dict, n_orbitals, string_ref, coeff_ref): ), ], ) -def test_import_state(molecule, basis, symm, method, wf_ref): +def test_import_state(molecule, basis, symm, wf_ref): r"""Test that cisd_state returns the correct wavefunction.""" mol = pyscf.gto.M(atom=molecule, basis=basis, symmetry=symm) myhf = pyscf.scf.UHF(mol).run() myci = pyscf.ci.UCISD(myhf).run() - wf_comp = qchem.convert.import_state(myci, method) + wf_comp = qchem.convert.import_state(myci) # overall sign could be different in each PySCF run assert np.allclose(wf_comp, wf_ref) or np.allclose(wf_comp, -wf_ref) def test_import_state_error(): - r"""Test that an error is raised by import_state if a wrong method symbol is entered.""" + r"""Test that an error is raised by import_state if a wrong object is entered.""" - myci = pyscf.ci.UCISD - method = "wrongmethod" + myci = "wrongobject" - with pytest.raises(ValueError, match="The supported method"): - _ = qchem.convert.import_state(myci, method) + with pytest.raises(ValueError, match="The supported option"): + _ = qchem.convert.import_state(myci) @pytest.mark.parametrize(("excitation"), [-1, 0, 3]) @@ -987,3 +985,17 @@ def test_excited_configurations_error(excitation): r"""Test that an error is raised by _excited_configurations if a wrong excitation is entered.""" with pytest.raises(ValueError, match="excitations are supported"): _ = qchem.convert._excited_configurations(2, 4, excitation) + + +def test_fail_import_pyscf(monkeypatch): + """Test if an ImportError is raised when pyscf is requested but not installed.""" + + mol = pyscf.gto.M(atom=[["H", (0, 0, 0)], ["H", (0, 0, 0.71)]], basis="sto6g") + myhf = pyscf.scf.UHF(mol).run() + myci = pyscf.ci.UCISD(myhf).run() + + with monkeypatch.context() as m: + m.setitem(sys.modules, "pyscf", None) + + with pytest.raises(ImportError, match="This feature requires pyscf"): + qml.qchem.convert.import_state(myci) From 92064007bdc53133ebd838fde3438bf024c55c3f Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 14:20:11 -0400 Subject: [PATCH 40/41] correct docstring example --- pennylane/qchem/convert.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index e5da9e389f0..73d09d22ae1 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -579,9 +579,12 @@ def import_state(solver, tol=1e-15): >>> mol = gto.M(atom=[['H', (0, 0, 0)], ['H', (0,0,0.71)]], basis='sto6g') >>> myhf = scf.UHF(mol).run() >>> myci = ci.UCISD(myhf).run() - >>> wf_cisd = cisd_state(myci, tol=1e-1) + >>> wf_cisd = qml.qchem.convert.import_state(myci, tol=1e-1) >>> print(wf_cisd) - {(1, 1): -0.9942969785398778, (2, 2): 0.10664669927602159} + [ 0. +0.j 0. +0.j 0. +0.j 0.1066467 +0.j + 0. +0.j 0. +0.j 0. +0.j 0. +0.j + 0. +0.j 0. +0.j 0. +0.j 0. +0.j + -0.99429698+0.j 0. +0.j 0. +0.j 0. +0.j] """ try: import pyscf From c2edd1bab1a79fc9d45945a4068f82b715f8d54c Mon Sep 17 00:00:00 2001 From: soranjh Date: Fri, 4 Aug 2023 14:25:15 -0400 Subject: [PATCH 41/41] correct string concatenation error --- pennylane/qchem/convert.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pennylane/qchem/convert.py b/pennylane/qchem/convert.py index 73d09d22ae1..37fde41a611 100644 --- a/pennylane/qchem/convert.py +++ b/pennylane/qchem/convert.py @@ -590,7 +590,7 @@ def import_state(solver, tol=1e-15): import pyscf except ImportError as Error: raise ImportError( - "This feature requires pyscf. " "It can be installed with: pip install pyscf" + "This feature requires pyscf. It can be installed with: pip install pyscf" ) from Error if isinstance(solver, pyscf.ci.ucisd.UCISD):