diff --git a/cudaqlib/kernels/uccsd.h b/cudaqlib/kernels/uccsd.h index d1030dd..1544c32 100644 --- a/cudaqlib/kernels/uccsd.h +++ b/cudaqlib/kernels/uccsd.h @@ -94,7 +94,7 @@ get_uccsd_excitations(std::size_t numElectrons, std::size_t numQubits) { doublesBeta); } -auto get_num_uccsd_parameters(std::size_t numElectrons, std::size_t numQubits) { +auto uccsd_num_parameters(std::size_t numElectrons, std::size_t numQubits) { auto [singlesAlpha, singlesBeta, doublesMixed, doublesAlpha, doublesBeta] = get_uccsd_excitations(numElectrons, numQubits); return singlesAlpha.size() + singlesBeta.size() + doublesMixed.size() + diff --git a/python/bindings/operators/py_operators.cpp b/python/bindings/operators/py_operators.cpp index 631c0c4..1b99fe6 100644 --- a/python/bindings/operators/py_operators.cpp +++ b/python/bindings/operators/py_operators.cpp @@ -26,6 +26,48 @@ void bindOperators(py::module &mod) { return fermion_compiler::get("jordan_wigner")->generate(op); }); + operators.def( + "jordan_wigner", + [](py::buffer hpq, py::buffer hpqrs, double core_energy = 0.0) { + auto hpqInfo = hpq.request(); + auto hpqrsInfo = hpqrs.request(); + fermion_op op(hpqInfo.shape[0], core_energy); + auto *hpqData = reinterpret_cast *>(hpqInfo.ptr); + auto *hpqrsData = + reinterpret_cast *>(hpqrsInfo.ptr); + std::vector> hpqVec( + hpqData, hpqData + hpqInfo.shape[0] * hpqInfo.shape[1]); + std::vector> hpqrsVec( + hpqrsData, hpqrsData + hpqrsInfo.shape[0] * hpqrsInfo.shape[1] * + hpqrsInfo.shape[2] * hpqrsInfo.shape[3]); + op.hpq.set_data(hpqVec); + op.hpqrs.set_data(hpqrsVec); + return fermion_compiler::get("jordan_wigner")->generate(op); + }, + py::arg("hpq"), py::arg("hpqrs"), py::arg("core_energy") = 0.0); + + operators.def( + "jordan_wigner", + [](py::buffer buffer, double core_energy = 0.0) { + auto info = buffer.request(); + fermion_op op(info.shape[0], core_energy); + auto *data = reinterpret_cast *>(info.ptr); + std::size_t size = 1; + for (auto &s : info.shape) + size *= s; + std::vector> vec(data, data + size); + if (info.shape.size() == 2) + op.hpq.set_data(vec); + else + op.hpqrs.set_data(vec); + return fermion_compiler::get("jordan_wigner")->generate(op); + }, + py::arg("hpq"), py::arg("core_energy") = 0.0); + + py::class_(operators, "FermionOperator", "") + .def_readonly("hpq", &fermion_op::hpq, "") + .def_readonly("hpqrs", &fermion_op::hpqrs, ""); + py::class_(operators, "fermion_compiler") .def_static( "get", @@ -70,15 +112,13 @@ void bindOperators(py::module &mod) { calculateStrides(m.shape)); }); - py::class_(operators, "FermionOperator", "") - .def_readonly("hpq", &fermion_op::hpq, "") - .def_readonly("hpqrs", &fermion_op::hpqrs, ""); - py::class_(operators, "MolecularHamiltonian") .def_readonly("energies", &molecular_hamiltonian::energies) .def_readonly("hamiltonian", &molecular_hamiltonian::hamiltonian) .def_readonly("n_electrons", &molecular_hamiltonian::n_electrons) - .def_readonly("n_orbitals", &molecular_hamiltonian::n_orbitals); + .def_readonly("n_orbitals", &molecular_hamiltonian::n_orbitals) + .def_readonly("fermion_operator", + &molecular_hamiltonian::fermionOperator); operators.def( "one_particle_op", @@ -87,10 +127,42 @@ void bindOperators(py::module &mod) { }, ""); + auto creator = [](molecular_geometry &molGeom, const std::string basis, + int spin, int charge, py::kwargs options) { + molecule_options inOptions; + inOptions.type = getValueOr(options, "type", "gas_phase"); + std::optional nele_cas = + getValueOr(options, "nele_cas", -1); + inOptions.nele_cas = nele_cas == -1 ? std::nullopt : nele_cas; + std::optional norb_cas = + getValueOr(options, "norb_cas", -1); + inOptions.norb_cas = norb_cas == -1 ? std::nullopt : norb_cas; + inOptions.symmetry = getValueOr(options, "symmetry", false); + inOptions.memory = getValueOr(options, "memory", 4000.); + inOptions.cycles = getValueOr(options, "cycles", 100); + inOptions.initguess = + getValueOr(options, "initguess", "minao"); + inOptions.UR = getValueOr(options, "UR", false); + inOptions.MP2 = getValueOr(options, "MP2", false); + inOptions.natorb = getValueOr(options, "natorb", false); + inOptions.casci = getValueOr(options, "casci", false); + inOptions.ccsd = getValueOr(options, "ccsd", false); + inOptions.casscf = getValueOr(options, "casscf", false); + inOptions.integrals_natorb = + getValueOr(options, "integrals_natorb", false); + inOptions.integrals_casscf = + getValueOr(options, "integrals_casscf", false); + inOptions.verbose = getValueOr(options, "verbose", false); + + if (inOptions.verbose) + inOptions.dump(); + return create_molecule(molGeom, basis, spin, charge, inOptions); + }; + operators.def( "create_molecule", - [](py::list geometry, const std::string basis, int spin, int charge, - py::kwargs options) { + [&](py::list geometry, const std::string basis, int spin, int charge, + py::kwargs options) { std::vector atoms; for (auto el : geometry) { if (!py::isinstance(el)) @@ -110,36 +182,16 @@ void bindOperators(py::module &mod) { } molecular_geometry molGeom(atoms); - molecule_options inOptions; - inOptions.type = getValueOr(options, "type", "gas_phase"); - std::optional nele_cas = - getValueOr(options, "nele_cas", -1); - inOptions.nele_cas = nele_cas == -1 ? std::nullopt : nele_cas; - std::optional norb_cas = - getValueOr(options, "norb_cas", -1); - inOptions.norb_cas = norb_cas == -1 ? std::nullopt : norb_cas; - inOptions.symmetry = getValueOr(options, "symmetry", false); - inOptions.memory = getValueOr(options, "memory", 4000.); - inOptions.cycles = getValueOr(options, "cycles", 100); - inOptions.initguess = - getValueOr(options, "initguess", "minao"); - inOptions.UR = getValueOr(options, "UR", false); - inOptions.MP2 = getValueOr(options, "MP2", false); - inOptions.natorb = getValueOr(options, "natorb", false); - inOptions.casci = getValueOr(options, "casci", false); - inOptions.ccsd = getValueOr(options, "ccsd", false); - inOptions.casscf = getValueOr(options, "casscf", false); - inOptions.integrals_natorb = - getValueOr(options, "integrals_natorb", false); - inOptions.integrals_casscf = - getValueOr(options, "integrals_casscf", false); - inOptions.verbose = getValueOr(options, "verbose", false); - - if (inOptions.verbose) - inOptions.dump(); - return create_molecule(molGeom, basis, spin, charge, inOptions); + return creator(molGeom, basis, spin, charge, options); }, ""); + + operators.def("create_molecule", + [&](const std::string &xyzFileName, const std::string basis, + int spin, int charge, py::kwargs options) { + auto geom = molecular_geometry::from_xyz(xyzFileName); + return creator(geom, basis, spin, charge, options); + }); } } // namespace cudaq::operators \ No newline at end of file diff --git a/python/tests/test_operators.py b/python/tests/test_operators.py index e600cd1..227178e 100644 --- a/python/tests/test_operators.py +++ b/python/tests/test_operators.py @@ -8,15 +8,22 @@ import os -import pytest +import pytest, pathlib import numpy as np import cudaq, cudaqlib +currentPath = pathlib.Path(__file__).parent.resolve() + + def test_operators(): geometry = [('H', (0., 0., 0.)), ('H', (0., 0., .7474))] - molecule = cudaqlib.operators.create_molecule( - geometry, 'sto-3g', 0, 0, verbose=True, casci=True) + molecule = cudaqlib.operators.create_molecule(geometry, + 'sto-3g', + 0, + 0, + verbose=True, + casci=True) print(molecule.hamiltonian.to_string()) print(molecule.energies) assert np.isclose(-1.11, molecule.energies['hf_energy'], atol=1e-2) @@ -26,11 +33,43 @@ def test_operators(): +def test_from_xyz_filename(): + molecule = cudaqlib.operators.create_molecule(str(currentPath) + + '/resources/LiH.xyz', + 'sto-3g', + 0, + 0, + verbose=True) + print(molecule.energies) + print(molecule.n_orbitals) + print(molecule.n_electrons) + assert molecule.n_orbitals == 6 + assert molecule.n_electrons == 4 + +def test_jordan_wigner(): + + geometry = [('H', (0., 0., 0.)), ('H', (0., 0., .7474))] + molecule = cudaqlib.operators.create_molecule(geometry, + 'sto-3g', + 0, + 0, + verbose=True, + casci=True) + op = cudaqlib.operators.jordan_wigner(molecule.fermion_operator) + print(op.to_string()) + assert molecule.hamiltonian == op + hpq = np.array(molecule.fermion_operator.hpq) + hpqrs = np.array(molecule.fermion_operator.hpqrs) + hpqJw = cudaqlib.operators.jordan_wigner(hpq, molecule.energies['nuclear_energy']) + hpqrsJw = cudaqlib.operators.jordan_wigner(hpqrs) + op2 = hpqJw + hpqrsJw + assert op2 == molecule.hamiltonian + def test_chemistry_operators(): - @cudaq.kernel - def ansatz(thetas:list[float]): + @cudaq.kernel + def ansatz(thetas: list[float]): q = cudaq.qvector(4) x(q[0]) x(q[1]) @@ -38,11 +77,14 @@ def ansatz(thetas:list[float]): # Define the molecule geometry = [('H', (0., 0., 0.)), ('H', (0., 0., .7474))] - molecule = cudaqlib.operators.create_molecule( - geometry, 'sto-3g', 0, 0, verbose=True) + molecule = cudaqlib.operators.create_molecule(geometry, + 'sto-3g', + 0, + 0, + verbose=True) # Run VQE energy, params, all_data = cudaqlib.gse.vqe(ansatz, - molecule.hamiltonian, [0.], - optimizer='cobyla') + molecule.hamiltonian, [0.], + optimizer='cobyla') assert np.isclose(-1.137, energy, atol=1e-2) diff --git a/tests/kernels/UCCSDTester.cpp b/tests/kernels/UCCSDTester.cpp index a4ff72d..c506780 100644 --- a/tests/kernels/UCCSDTester.cpp +++ b/tests/kernels/UCCSDTester.cpp @@ -314,7 +314,7 @@ TEST(UCCSDTester, checkUCCSDAnsatz) { EXPECT_TRUE(singlesBeta.size() == 1); EXPECT_TRUE(doublesMixed.size() == 1); - auto numParams = cudaq::get_num_uccsd_parameters(numElectrons, numQubits); + auto numParams = cudaq::uccsd_num_parameters(numElectrons, numQubits); std::vector init{-2., -2., -2.}; cudaq::optim::cobyla optimizer;