Chemistry library, SPSA sample, and notebook

Work in progress...

build.py

@@ -69,10 +69,10 @@
 )
 
 parser.add_argument(
-        "--ci-bench",
-        action=argparse.BooleanOptionalAction,
-        default=False,
-        help="Run the benchmarking script that is run in CI (default is --no-ci-bench)",
+    "--ci-bench",
+    action=argparse.BooleanOptionalAction,
+    default=False,
+    help="Run the benchmarking script that is run in CI (default is --no-ci-bench)",
 )
 
 args = parser.parse_args()
@@ -303,25 +303,46 @@ def run_python_integration_tests(cwd, interpreter):
 def run_ci_historic_benchmark():
     branch = "main"
     output = subprocess.check_output(
-        ["git", "rev-list", "--since=1 week ago", "--pretty=format:%ad__%h", "--date=short", branch]
+        [
+            "git",
+            "rev-list",
+            "--since=1 week ago",
+            "--pretty=format:%ad__%h",
+            "--date=short",
+            branch,
+        ]
     ).decode("utf-8")
-    print('\n'.join([line for i, line in enumerate(output.split('\n')) if i % 2 == 1]))
+    print("\n".join([line for i, line in enumerate(output.split("\n")) if i % 2 == 1]))
 
     output = subprocess.check_output(
-        ["git", "rev-list", "--since=1 week ago", "--pretty=format:%ad__%h", "--date=short", branch]
+        [
+            "git",
+            "rev-list",
+            "--since=1 week ago",
+            "--pretty=format:%ad__%h",
+            "--date=short",
+            branch,
+        ]
     ).decode("utf-8")
-    date_and_commits = [line for i, line in enumerate(output.split('\n')) if i % 2 == 1]
+    date_and_commits = [line for i, line in enumerate(output.split("\n")) if i % 2 == 1]
 
     for date_and_commit in date_and_commits:
         print("benching commit", date_and_commit)
         result = subprocess.run(
-            ["cargo", "criterion", "--message-format=json", "--history-id", date_and_commit],
+            [
+                "cargo",
+                "criterion",
+                "--message-format=json",
+                "--history-id",
+                date_and_commit,
+            ],
             capture_output=True,
-            text=True
+            text=True,
         )
         with open(f"{date_and_commit}.json", "w") as f:
             f.write(result.stdout)
 
+
 if build_pip:
     step_start("Building the pip package")
 
@@ -514,6 +535,7 @@ def run_ci_historic_benchmark():
         "ipykernel",
         "nbconvert",
         "pandas",
+        "qutip",
         "qiskit>=1.3.0,<2.0.0",
     ]
     subprocess.run(pip_install_args, check=True, text=True, cwd=root_dir, env=pip_env)

library/chemistry/qsharp.json

@@ -0,0 +1,42 @@
+{
+  "lints": [
+    {
+      "lint": "deprecatedNewtype",
+      "level": "error"
+    },
+    {
+      "lint": "deprecatedFunctionConstructor",
+      "level": "error"
+    },
+    {
+      "lint": "deprecatedWithOperator",
+      "level": "error"
+    },
+    {
+      "lint": "deprecatedDoubleColonOperator",
+      "level": "error"
+    },
+    {
+      "lint": "deprecatedSet",
+      "level": "error"
+    },
+    {
+      "lint": "needlessParens",
+      "level": "error"
+    }
+  ],
+  "files": [
+    "src/Tests.qs",
+    "src/Utils.qs",
+    "src/JordanWigner/Convenience.qs",
+    "src/JordanWigner/ConvenienceVQE.qs",
+    "src/JordanWigner/JordanWignerBlockEncoding.qs",
+    "src/JordanWigner/JordanWignerClusterOperatorEvolutionSet.qs",
+    "src/JordanWigner/JordanWignerEvolutionSet.qs",
+    "src/JordanWigner/JordanWignerOptimizedBlockEncoding.qs",
+    "src/JordanWigner/JordanWignerVQE.qs",
+    "src/JordanWigner/OptimizedBEOperator.qs",
+    "src/JordanWigner/StatePreparation.qs",
+    "src/JordanWigner/Utils.qs"
+  ]
+}

library/chemistry/src/JordanWigner/Convenience.qs

@@ -0,0 +1,112 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+export
+    TrotterStepOracle,
+    QubitizationOracle,
+    OptimizedQubitizationOracle;
+
+import JordanWigner.JordanWignerBlockEncoding.JordanWignerBlockEncodingGeneratorSystem;
+import JordanWigner.JordanWignerEvolutionSet.JordanWignerFermionEvolutionSet;
+import JordanWigner.JordanWignerEvolutionSet.JordanWignerGeneratorSystem;
+import JordanWigner.JordanWignerOptimizedBlockEncoding.JordanWignerOptimizedBlockEncoding;
+import JordanWigner.JordanWignerOptimizedBlockEncoding.PauliBlockEncoding;
+import JordanWigner.JordanWignerOptimizedBlockEncoding.QuantumWalkByQubitization;
+import JordanWigner.StatePreparation.PrepareArbitraryStateD;
+import JordanWigner.Utils.JordanWignerEncodingData;
+import JordanWigner.Utils.MultiplexOperationsFromGenerator;
+import JordanWigner.Utils.TrotterSimulationAlgorithm;
+import Std.Convert.IntAsDouble;
+import Std.Math.Ceiling;
+import Std.Math.Lg;
+import Utils.EvolutionGenerator;
+import Utils.GeneratorSystem;
+
+// Convenience functions for performing simulation.
+
+/// # Summary
+/// Returns Trotter step operation and the parameters necessary to run it.
+///
+/// # Input
+/// ## qSharpData
+/// Hamiltonian described by `JordanWignerEncodingData` format.
+/// ## trotterStepSize
+/// Step size of Trotter integrator.
+/// ## trotterOrder
+/// Order of Trotter integrator.
+///
+/// # Output
+/// A tuple where: `Int` is the number of qubits allocated,
+/// `Double` is `1.0/trotterStepSize`, and the operation
+/// is the Trotter step.
+function TrotterStepOracle(qSharpData : JordanWignerEncodingData, trotterStepSize : Double, trotterOrder : Int) : (Int, (Double, (Qubit[] => Unit is Adj + Ctl))) {
+    let (nSpinOrbitals, data, statePrepData, energyShift) = qSharpData!;
+    let generatorSystem = JordanWignerGeneratorSystem(data);
+    let evolutionGenerator = new EvolutionGenerator { EvolutionSet = JordanWignerFermionEvolutionSet(), System = generatorSystem };
+    let simulationAlgorithm = TrotterSimulationAlgorithm(trotterStepSize, trotterOrder);
+    let oracle = simulationAlgorithm(trotterStepSize, evolutionGenerator, _);
+    let nTargetRegisterQubits = nSpinOrbitals;
+    let rescaleFactor = 1.0 / trotterStepSize;
+    return (nTargetRegisterQubits, (rescaleFactor, oracle));
+}
+
+
+function QubitizationOracleSeperatedRegisters(qSharpData : JordanWignerEncodingData) : ((Int, Int), (Double, ((Qubit[], Qubit[]) => Unit is Adj + Ctl))) {
+    let (nSpinOrbitals, data, statePrepData, energyShift) = qSharpData!;
+    let generatorSystem = JordanWignerBlockEncodingGeneratorSystem(data);
+    let (nTerms, genIdxFunction) = generatorSystem!;
+    let (oneNorm, blockEncodingReflection) = PauliBlockEncoding(generatorSystem);
+    let nTargetRegisterQubits = nSpinOrbitals;
+    let nCtrlRegisterQubits = Ceiling(Lg(IntAsDouble(nTerms)));
+    return ((nCtrlRegisterQubits, nTargetRegisterQubits), (oneNorm, QuantumWalkByQubitization(blockEncodingReflection)));
+}
+
+/// # Summary
+/// Returns Qubitization operation and the parameters necessary to run it.
+///
+/// # Input
+/// ## qSharpData
+/// Hamiltonian described by `JordanWignerEncodingData` format.
+///
+/// # Output
+/// A tuple where: `Int` is the number of qubits allocated,
+/// `Double` is the one-norm of Hamiltonian coefficients, and the operation
+/// is the Quantum walk created by Qubitization.
+function QubitizationOracle(qSharpData : JordanWignerEncodingData) : (Int, (Double, (Qubit[] => Unit is Adj + Ctl))) {
+    let ((nCtrlRegisterQubits, nTargetRegisterQubits), (oneNorm, oracle)) = QubitizationOracleSeperatedRegisters(qSharpData);
+    let nQubits = nCtrlRegisterQubits + nTargetRegisterQubits;
+    return (nQubits, (oneNorm, MergeTwoRegisters(oracle, nTargetRegisterQubits, _)));
+}
+
+
+operation MergeTwoRegisters(oracle : ((Qubit[], Qubit[]) => Unit is Adj + Ctl), nSystemQubits : Int, allQubits : Qubit[]) : Unit is Adj + Ctl {
+    oracle(allQubits[nSystemQubits..Length(allQubits) - 1], allQubits[0..nSystemQubits - 1]);
+}
+
+
+function OptimizedQubitizationOracleSeperatedRegisters(qSharpData : JordanWignerEncodingData, targetError : Double) : ((Int, Int), (Double, ((Qubit[], Qubit[]) => Unit is Adj + Ctl))) {
+    let (nSpinOrbitals, data, statePrepData, energyShift) = qSharpData!;
+    let ((nCtrlRegisterQubits, nTargetRegisterQubits), (oneNorm, blockEncodingReflection)) = JordanWignerOptimizedBlockEncoding(targetError, data, nSpinOrbitals);
+    return ((nCtrlRegisterQubits, nTargetRegisterQubits), (oneNorm, QuantumWalkByQubitization(blockEncodingReflection)));
+}
+
+
+/// # Summary
+/// Returns T-count optimized Qubitization operation
+/// and the parameters necessary to run it.
+///
+/// # Input
+/// ## qSharpData
+/// Hamiltonian described by `JordanWignerEncodingData` format.
+/// ## targetError
+/// Error of auxillary state preparation step.
+///
+/// # Output
+/// A tuple where: `Int` is the number of qubits allocated,
+/// `Double` is the one-norm of Hamiltonian coefficients, and the operation
+/// is the Quantum walk created by Qubitization.
+function OptimizedQubitizationOracle(qSharpData : JordanWignerEncodingData, targetError : Double) : (Int, (Double, (Qubit[] => Unit is Adj + Ctl))) {
+    let ((nCtrlRegisterQubits, nTargetRegisterQubits), (oneNorm, oracle)) = OptimizedQubitizationOracleSeperatedRegisters(qSharpData, targetError);
+    let nQubits = nCtrlRegisterQubits + nTargetRegisterQubits;
+    return (nQubits, (oneNorm, MergeTwoRegisters(oracle, nTargetRegisterQubits, _)));
+}

library/chemistry/src/JordanWigner/ConvenienceVQE.qs

@@ -0,0 +1,93 @@
+// Copyright (c) Microsoft Corporation. All rights reserved.
+// Licensed under the MIT License.
+
+export
+    EstimateEnergyWrapper,
+    EstimateEnergy;
+
+import JordanWigner.JordanWignerEvolutionSet.JordanWignerGeneratorSystem;
+import JordanWigner.JordanWignerVQE.EstimateTermExpectation;
+import JordanWigner.JordanWignerVQE.ExpandedCoefficients;
+import JordanWigner.JordanWignerVQE.MeasurementOperators;
+import JordanWigner.StatePreparation.PrepareTrialState;
+import JordanWigner.Utils.JordanWignerEncodingData;
+import JordanWigner.Utils.JordanWignerInputState;
+import JordanWigner.Utils.JWOptimizedHTerms;
+import Std.Arrays.ForEach;
+
+/// # Summary
+/// Estimates the energy of the molecule by summing the energy contributed by the individual Jordan-Wigner terms.
+/// This convenience wrapper takes input in raw data types and converts them to the Jordan-Wigner encoding data type.
+///
+/// # Description
+/// This operation implicitly relies on the spin up-down indexing convention.
+///
+/// # Input
+/// ## jwHamiltonian
+/// The Jordan-Wigner Hamiltonian, represented in simple data types rather than a struct.
+/// ## nSamples
+/// The number of samples to use for the estimation of the term expectations.
+///
+/// # Output
+/// The estimated energy of the molecule
+operation EstimateEnergyWrapper(jwHamiltonian : (Int, ((Int[], Double[])[], (Int[], Double[])[], (Int[], Double[])[], (Int[], Double[])[]), (Int, ((Double, Double), Int[])[]), Double), nSamples : Int) : Double {
+    let (nQubits, jwTerms, inputState, energyOffset) = jwHamiltonian;
+    let (hterm0, hterm1, hterm2, hterm3) = jwTerms;
+    let jwTerms = new JWOptimizedHTerms { HTerm0 = hterm0, HTerm1 = hterm1, HTerm2 = hterm2, HTerm3 = hterm3 };
+    let (inputState1, inputState2) = inputState;
+    mutable jwInputState = [];
+    for entry in inputState2 {
+        let (amp, idicies) = entry;
+        jwInputState += [new JordanWignerInputState { Amplitude = amp, FermionIndices = idicies }];
+    }
+    let inputState = (inputState1, jwInputState);
+    let jwHamiltonian = new JordanWignerEncodingData {
+        NumQubits = nQubits,
+        Terms = jwTerms,
+        InputState = inputState,
+        EnergyOffset = energyOffset
+    };
+    return EstimateEnergy(jwHamiltonian, nSamples);
+}
+
+/// # Summary
+/// Estimates the energy of the molecule by summing the energy contributed by the individual Jordan-Wigner terms.
+///
+/// # Description
+/// This operation implicitly relies on the spin up-down indexing convention.
+///
+/// # Input
+/// ## jwHamiltonian
+/// The Jordan-Wigner Hamiltonian.
+/// ## nSamples
+/// The number of samples to use for the estimation of the term expectations.
+///
+/// # Output
+/// The estimated energy of the molecule
+operation EstimateEnergy(jwHamiltonian : JordanWignerEncodingData, nSamples : Int) : Double {
+    // Initialize return value
+    mutable energy = 0.;
+
+    // Unpack information and qubit Hamiltonian terms
+    let (nQubits, jwTerms, inputState, energyOffset) = jwHamiltonian!;
+
+    // Loop over all qubit Hamiltonian terms
+    let (nTerms, indexFunction) = (JordanWignerGeneratorSystem(jwTerms))!;
+
+    for idxTerm in 0..nTerms - 1 {
+        let term = indexFunction(idxTerm);
+        let ((idxTermType, coeff), idxFermions) = term!;
+        let termType = idxTermType[0];
+
+        let ops = MeasurementOperators(nQubits, idxFermions, termType);
+        let coeffs = ExpandedCoefficients(coeff, termType);
+
+        // The private wrapper enables fast emulation during expectation estimation
+        let inputStateUnitary = PrepareTrialState(inputState, _);
+
+        let jwTermEnergy = EstimateTermExpectation(inputStateUnitary, ops, coeffs, nQubits, nSamples);
+        energy += jwTermEnergy;
+    }
+
+    return energy + energyOffset;
+}