Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Minor fixes and updates for chemistry workflows #9

Merged
merged 1 commit into from
Jul 29, 2024
Merged

Minor fixes and updates for chemistry workflows #9

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion cudaqlib/kernels/uccsd.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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() +
Expand Down
122 changes: 87 additions & 35 deletions python/bindings/operators/py_operators.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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<std::complex<double> *>(hpqInfo.ptr);
auto *hpqrsData =
reinterpret_cast<std::complex<double> *>(hpqrsInfo.ptr);
std::vector<std::complex<double>> hpqVec(
hpqData, hpqData + hpqInfo.shape[0] * hpqInfo.shape[1]);
std::vector<std::complex<double>> 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<std::complex<double> *>(info.ptr);
std::size_t size = 1;
for (auto &s : info.shape)
size *= s;
std::vector<std::complex<double>> 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_<fermion_op>(operators, "FermionOperator", "")
.def_readonly("hpq", &fermion_op::hpq, "")
.def_readonly("hpqrs", &fermion_op::hpqrs, "");

py::class_<fermion_compiler>(operators, "fermion_compiler")
.def_static(
"get",
Expand Down Expand Up @@ -70,15 +112,13 @@ void bindOperators(py::module &mod) {
calculateStrides(m.shape));
});

py::class_<fermion_op>(operators, "FermionOperator", "")
.def_readonly("hpq", &fermion_op::hpq, "")
.def_readonly("hpqrs", &fermion_op::hpqrs, "");

py::class_<molecular_hamiltonian>(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",
Expand All @@ -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<std::string>(options, "type", "gas_phase");
std::optional<std::size_t> nele_cas =
getValueOr<std::size_t>(options, "nele_cas", -1);
inOptions.nele_cas = nele_cas == -1 ? std::nullopt : nele_cas;
std::optional<std::size_t> norb_cas =
getValueOr<std::size_t>(options, "norb_cas", -1);
inOptions.norb_cas = norb_cas == -1 ? std::nullopt : norb_cas;
inOptions.symmetry = getValueOr<bool>(options, "symmetry", false);
inOptions.memory = getValueOr<double>(options, "memory", 4000.);
inOptions.cycles = getValueOr<std::size_t>(options, "cycles", 100);
inOptions.initguess =
getValueOr<std::string>(options, "initguess", "minao");
inOptions.UR = getValueOr<bool>(options, "UR", false);
inOptions.MP2 = getValueOr<bool>(options, "MP2", false);
inOptions.natorb = getValueOr<bool>(options, "natorb", false);
inOptions.casci = getValueOr<bool>(options, "casci", false);
inOptions.ccsd = getValueOr<bool>(options, "ccsd", false);
inOptions.casscf = getValueOr<bool>(options, "casscf", false);
inOptions.integrals_natorb =
getValueOr<bool>(options, "integrals_natorb", false);
inOptions.integrals_casscf =
getValueOr<bool>(options, "integrals_casscf", false);
inOptions.verbose = getValueOr<bool>(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<atom> atoms;
for (auto el : geometry) {
if (!py::isinstance<py::tuple>(el))
Expand All @@ -110,36 +182,16 @@ void bindOperators(py::module &mod) {
}
molecular_geometry molGeom(atoms);

molecule_options inOptions;
inOptions.type = getValueOr<std::string>(options, "type", "gas_phase");
std::optional<std::size_t> nele_cas =
getValueOr<std::size_t>(options, "nele_cas", -1);
inOptions.nele_cas = nele_cas == -1 ? std::nullopt : nele_cas;
std::optional<std::size_t> norb_cas =
getValueOr<std::size_t>(options, "norb_cas", -1);
inOptions.norb_cas = norb_cas == -1 ? std::nullopt : norb_cas;
inOptions.symmetry = getValueOr<bool>(options, "symmetry", false);
inOptions.memory = getValueOr<double>(options, "memory", 4000.);
inOptions.cycles = getValueOr<std::size_t>(options, "cycles", 100);
inOptions.initguess =
getValueOr<std::string>(options, "initguess", "minao");
inOptions.UR = getValueOr<bool>(options, "UR", false);
inOptions.MP2 = getValueOr<bool>(options, "MP2", false);
inOptions.natorb = getValueOr<bool>(options, "natorb", false);
inOptions.casci = getValueOr<bool>(options, "casci", false);
inOptions.ccsd = getValueOr<bool>(options, "ccsd", false);
inOptions.casscf = getValueOr<bool>(options, "casscf", false);
inOptions.integrals_natorb =
getValueOr<bool>(options, "integrals_natorb", false);
inOptions.integrals_casscf =
getValueOr<bool>(options, "integrals_casscf", false);
inOptions.verbose = getValueOr<bool>(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
60 changes: 51 additions & 9 deletions python/tests/test_operators.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand All @@ -26,23 +33,58 @@ 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])
exp_pauli(thetas[0], q, 'XXXY')

# 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)
2 changes: 1 addition & 1 deletion tests/kernels/UCCSDTester.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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<double> init{-2., -2., -2.};

cudaq::optim::cobyla optimizer;
Expand Down