diff --git a/Project.toml b/Project.toml index 91a9be87..9574031e 100644 --- a/Project.toml +++ b/Project.toml @@ -61,6 +61,7 @@ Mocking = "78c3b35d-d492-501b-9361-3d52fe80e533" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StructTypes = "856f2bd8-1eba-4b0a-8007-ebc267875bd4" Tar = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/PyBraket/CondaPkg.toml b/PyBraket/CondaPkg.toml index 617d3cdc..4477a755 100644 --- a/PyBraket/CondaPkg.toml +++ b/PyBraket/CondaPkg.toml @@ -5,4 +5,4 @@ scipy = "" numpy = "" [pip.deps] -amazon-braket-sdk = ">=1.33.0.dev0" +amazon-braket-sdk = ">=1.35.0" diff --git a/PyBraket/Project.toml b/PyBraket/Project.toml index aef98400..02f6f6f6 100644 --- a/PyBraket/Project.toml +++ b/PyBraket/Project.toml @@ -11,6 +11,7 @@ DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MicroMamba = "0b3b1443-0f03-428d-bdfb-f27f9c1191ea" PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StructTypes = "856f2bd8-1eba-4b0a-8007-ebc267875bd4" [compat] @@ -30,4 +31,5 @@ Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" CondaPkg = "992eb4ea-22a4-4c89-a5bb-47a3300528ab" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/PyBraket/src/pycircuit.jl b/PyBraket/src/pycircuit.jl index f0ecdb78..ed349cc3 100644 --- a/PyBraket/src/pycircuit.jl +++ b/PyBraket/src/pycircuit.jl @@ -29,18 +29,20 @@ end PyCircuit(c::Circuit) = convert(PyCircuit, c) Py(o::Braket.Observables.HermitianObservable) = braketobs.Hermitian(Py(o.matrix)) -Py(o::Braket.Observables.TensorProduct) = braketobs.TensorProduct(pylist(map(Py, o.factors))) -Py(o::Expectation) = circuit.result_types.Expectation(Py(o.observable), pylist(o.targets)) -Py(o::Variance) = circuit.result_types.Variance(Py(o.observable), pylist(o.targets)) -Py(o::Sample) = circuit.result_types.Sample(Py(o.observable), pylist(o.targets)) -Py(o::Probability) = circuit.result_types.Probability(pylist(o.targets)) -Py(o::DensityMatrix) = circuit.result_types.DensityMatrix(pylist(o.targets)) -Py(o::Amplitude) = circuit.result_types.Amplitude(pylist(map(s->Py(s), o.states))) -Py(o::StateVector) = circuit.result_types.StateVector() +Py(o::Braket.Observables.TensorProduct) = o.coefficient * braketobs.TensorProduct(pylist(map(Py, o.factors))) +Py(o::Braket.Observables.Sum) = braketobs.Sum(pylist(map(Py, o.summands))) +Py(o::AdjointGradient) = circuit.result_types.AdjointGradient(Py(o.observable), pylist(o.target), pylist(o.parameters)) +Py(o::Expectation) = circuit.result_types.Expectation(Py(o.observable), pylist(o.targets)) +Py(o::Variance) = circuit.result_types.Variance(Py(o.observable), pylist(o.targets)) +Py(o::Sample) = circuit.result_types.Sample(Py(o.observable), pylist(o.targets)) +Py(o::Probability) = circuit.result_types.Probability(pylist(o.targets)) +Py(o::DensityMatrix) = circuit.result_types.DensityMatrix(pylist(o.targets)) +Py(o::Amplitude) = circuit.result_types.Amplitude(pylist(map(s->Py(s), o.states))) +Py(o::StateVector) = circuit.result_types.StateVector() for (typ, py_typ, label) in ((:(Braket.Observables.H), :H, "h"), (:(Braket.Observables.X), :X, "x"), (:(Braket.Observables.Y), :Y, "y"), (:(Braket.Observables.Z), :Z, "z"), (:(Braket.Observables.I), :I, "i")) @eval begin - Py(o::$typ) = braketobs.$py_typ() + Py(o::$typ) = o.coefficient * braketobs.$py_typ() end end diff --git a/PyBraket/src/pyschema.jl b/PyBraket/src/pyschema.jl index 117aec97..c7918222 100644 --- a/PyBraket/src/pyschema.jl +++ b/PyBraket/src/pyschema.jl @@ -2,7 +2,7 @@ module PySchema using PythonCall, Braket, Braket.IR import Braket: DwaveTiming, GateModelQpuParadigmProperties, ResultTypeValue, OqcDeviceCapabilities, GateFidelity2Q, OneQubitProperties, AdditionalMetadata, DwaveAdvantageDeviceLevelParameters, IonqDeviceCapabilities, TwoQubitProperties, IonqDeviceParameters, NativeQuilMetadata, QueraDeviceCapabilities, ExecutionDay, StandardizedGateModelQpuDeviceProperties, XanaduDeviceCapabilities, DeviceExecutionWindow, DwaveDeviceParameters, CoherenceTime, DwaveMetadata, Frame, DeviceCost, Dwave2000QDeviceLevelParameters, FidelityType, QueraMetadata, OqcProviderProperties, DeviceServiceProperties, XanaduDeviceParameters, BlackbirdProgram, Geometry, Fidelity1Q, AnnealingTaskResult, AnalogHamiltonianSimulationShotResult, RigettiDeviceCapabilities, DwaveProviderProperties, DeviceActionProperties, DwaveProviderLevelParameters, JaqcdDeviceActionProperties, Direction, PulseDeviceActionProperties, BlackbirdDeviceActionProperties, PerformanceLattice, PerformanceRydberg, DeviceActionType, IonqProviderProperties, PersistedJobDataFormat, PhotonicModelTaskResult, DeviceDocumentation, XanaduProviderProperties, QubitDirection, ContinuousVariableQpuParadigmProperties, PulseFunctionArgument, DwaveAdvantageDeviceParameters, RigettiDeviceParameters, AnalogHamiltonianSimulationTaskResult, TaskMetadata, GateModelSimulatorParadigmProperties, GateModelTaskResult, RigettiProviderProperties, OqcDeviceParameters, DeviceConnectivity, PulseFunction, Rydberg, PerformanceRydbergGlobal, RydbergGlobal, GateModelSimulatorDeviceCapabilities, XanaduMetadata, GateModelParameters, Performance, GateModelSimulatorDeviceParameters, AnalogHamiltonianSimulationShotMeasurement, OqcMetadata, QueraAhsParadigmProperties, Lattice, PersistedJobData, DwaveDeviceCapabilities, ProblemType, RigettiMetadata, SimulatorMetadata, Area, Problem, ResultType, PostProcessingType, ResultFormat, Dwave2000QDeviceParameters, AnalogHamiltonianSimulationShotMetadata, braketSchemaHeader, OpenQasmProgram, Port, OpenQASMDeviceActionProperties, AbstractProgram -import Braket.IR: Z, Sample, CPhaseShift01, PhaseDamping, Rz, GeneralizedAmplitudeDamping, XX, ZZ, PhaseFlip, Vi, Depolarizing, Variance, TwoQubitDepolarizing, DensityMatrix, CPhaseShift00, ECR, CompilerDirective, CCNot, Unitary, BitFlip, Y, Swap, CZ, EndVerbatimBox, Program, CNot, CSwap, Ry, I, Si, AmplitudeDamping, StateVector, ISwap, H, XY, YY, T, TwoQubitDephasing, X, Ti, CV, StartVerbatimBox, PauliChannel, PSwap, Expectation, Probability, PhaseShift, V, CPhaseShift, S, Rx, Kraus, Amplitude, CPhaseShift10, MultiQubitPauliChannel, CY, Setup, Hamiltonian, ShiftingField, AtomArrangement, TimeSeries, PhysicalField, AHSProgram, DrivingField, AbstractProgramResult +import Braket.IR: Z, Sample, CPhaseShift01, PhaseDamping, Rz, GeneralizedAmplitudeDamping, XX, ZZ, PhaseFlip, Vi, Depolarizing, Variance, TwoQubitDepolarizing, DensityMatrix, CPhaseShift00, ECR, CompilerDirective, CCNot, Unitary, BitFlip, Y, Swap, CZ, EndVerbatimBox, Program, CNot, AdjointGradient, CSwap, Ry, I, Si, AmplitudeDamping, StateVector, ISwap, H, XY, YY, T, TwoQubitDephasing, X, Ti, CV, StartVerbatimBox, PauliChannel, PSwap, Expectation, Probability, PhaseShift, V, CPhaseShift, S, Rx, Kraus, Amplitude, CPhaseShift10, MultiQubitPauliChannel, CY, Setup, Hamiltonian, ShiftingField, AtomArrangement, TimeSeries, PhysicalField, AHSProgram, DrivingField, AbstractProgramResult, IRObservable const instructions = PythonCall.pynew() const dwave_metadata_v1 = PythonCall.pynew() @@ -71,6 +71,12 @@ const shared_models = PythonCall.pynew() const schema_header = PythonCall.pynew() const openqasm_device_action_properties = PythonCall.pynew() +function union_convert(::Type{IRObservable}, x) + PythonCall.pyisinstance(x, PythonCall.pybuiltins.str) && return pyconvert(String, x) + x_vec = Union{String, Vector{Vector{Vector{Float64}}}}[PythonCall.pyisinstance(x_, PythonCall.pybuiltins.str) ? pyconvert(String, x_) : pyconvert(Vector{Vector{Vector{Float64}}}, x_) for x_ in x] + return x_vec +end + function union_convert(union_type, x) union_ts = union_type isa Union ? Type[] : [union_type] union_t = union_type @@ -117,11 +123,8 @@ function jl_convert_attr(n, t, attr) return pyconvert(t, attr) end else - if !PythonCall.pyisnone(attr) - return union_convert(t, attr) - else - return nothing - end + PythonCall.pyisnone(attr) && return nothing + return union_convert(t, attr) end end @@ -148,16 +151,28 @@ function jl_convert(::Type{T}, x::Py) where {T<:AbstractIR} end PythonCall.pyconvert_return(T(args..., pyconvert(String, pygetattr(x, "type")))) end - function jl_convert(::Type{AbstractProgram}, x::Py) T = Braket.lookup_type(pyconvert(braketSchemaHeader, pygetattr(x, "braketSchemaHeader"))) return PythonCall.pyconvert_return(pyconvert(T, x)) end function jl_convert(::Type{AbstractProgramResult}, x::Py) T_ = pyconvert(String, x.type) - T = Dict("expectation"=>Expectation, "variance"=>Variance, "statevector"=>StateVector, "densitymatrix"=>DensityMatrix, "sample"=>Sample, "amplitude"=>Amplitude, "probability"=>Probability) + T = Dict("adjoint_gradient"=>AdjointGradient, "expectation"=>Expectation, "variance"=>Variance, "statevector"=>StateVector, "densitymatrix"=>DensityMatrix, "sample"=>Sample, "amplitude"=>Amplitude, "probability"=>Probability) return PythonCall.pyconvert_return(pyconvert(T[T_], x)) end +#=for (irT, pyT) in ((:(Braket.IR.Expectation), :(pyjaqcd.Expectation)), + (:(Braket.IR.Variance), :(pyjaqcd.Variance)), + (:(Braket.IR.Sample), :(pyjaqcd.Sample)), + (:(Braket.IR.Amplitude), :(pyjaqcd.Amplitude)), + (:(Braket.IR.StateVector), :(pyjaqcd.StateVector)), + (:(Braket.IR.Probability), :(pyjaqcd.Probability)), + (:(Braket.IR.DensityMatrix), :(pyjaqcd.DensityMatrix)), + (:(Braket.IR.AdjointGradient), :(pyjaqcd.AdjointGradient))) + @eval begin + Py(o::$irT) = $pyT(;arg_gen(o, fieldnames($irT))...) + end +end +=# function __init__() PythonCall.pycopy!(instructions, pyimport("braket.ir.jaqcd.instructions")) @@ -363,6 +378,7 @@ function __init__() PythonCall.pyconvert_add_rule("braket.ir.annealing.problem_v1:ProblemType", ProblemType, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.instructions:ECR", ECR, jl_convert) PythonCall.pyconvert_add_rule("braket.task_result.rigetti_metadata_v1:RigettiMetadata", RigettiMetadata, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:AdjointGradient", AdjointGradient, jl_convert) PythonCall.pyconvert_add_rule("braket.task_result.simulator_metadata_v1:SimulatorMetadata", SimulatorMetadata, jl_convert) PythonCall.pyconvert_add_rule("braket.device_schema.quera.quera_ahs_paradigm_properties_v1:Area", Area, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.annealing.problem_v1:Problem", Problem, jl_convert) @@ -386,11 +402,12 @@ function __init__() PythonCall.pyconvert_add_rule("braket.device_schema.openqasm_device_action_properties:OpenQASMDeviceActionProperties", OpenQASMDeviceActionProperties, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.program_v1:Program", Program, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.program_v1:Program", AbstractProgram, jl_convert) - PythonCall.pyconvert_add_rule("braket.ir.openqasm.program_v1:Program", OpenQasmProgram, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.jaqcd.program_v1:Program", AbstractProgram, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.openqasm.program_v1:OpenQasmProgram", AbstractProgram, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.openqasm.program_v1:Program", AbstractProgram, jl_convert) - PythonCall.pyconvert_add_rule("braket.ir.blackbird.program_v1:Program", BlackbirdProgram, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.blackbird.program_v1:BlackbirdProgram", AbstractProgram, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.blackbird.program_v1:Program", AbstractProgram, jl_convert) - PythonCall.pyconvert_add_rule("braket.ir.ahs.program_v1:Program", AHSProgram, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.ahs.program_v1:AHSProgram", AbstractProgram, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.ahs.program_v1:Program", AbstractProgram, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:Amplitude", AbstractProgramResult, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:Expectation", AbstractProgramResult, jl_convert) @@ -399,6 +416,7 @@ function __init__() PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:StateVector", AbstractProgramResult, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:DensityMatrix", AbstractProgramResult, jl_convert) PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:Variance", AbstractProgramResult, jl_convert) + PythonCall.pyconvert_add_rule("braket.ir.jaqcd.results:AdjointGradient", AbstractProgramResult, jl_convert) end diff --git a/PyBraket/test/CondaPkg.toml b/PyBraket/test/CondaPkg.toml index c7d46f42..f6f99c84 100644 --- a/PyBraket/test/CondaPkg.toml +++ b/PyBraket/test/CondaPkg.toml @@ -5,4 +5,4 @@ scipy = "" numpy = "" [pip.deps] -amazon-braket-sdk = ">=1.33.0.dev0" +amazon-braket-sdk = ">=1.35.0" diff --git a/PyBraket/test/Project.toml b/PyBraket/test/Project.toml index 9ec11c0f..18cb87d7 100644 --- a/PyBraket/test/Project.toml +++ b/PyBraket/test/Project.toml @@ -6,6 +6,7 @@ CondaPkg = "992eb4ea-22a4-4c89-a5bb-47a3300528ab" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" PyBraket = "e85266a6-1825-490b-a80e-9b9469c53660" PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [compat] diff --git a/PyBraket/test/circuits.jl b/PyBraket/test/circuits.jl index dd31184a..67102166 100644 --- a/PyBraket/test/circuits.jl +++ b/PyBraket/test/circuits.jl @@ -1,163 +1,185 @@ using Braket, PyBraket, Test, LinearAlgebra, PythonCall using PythonCall: pyconvert, Py, pyisTrue, pyisinstance +@testset "PyBraket circuits" begin + @testset for ir_type in (:JAQCD, :OpenQASM) + Braket.IRType[] = ir_type + @testset "Expectation" begin + c = Circuit([(H, 0), (CNot, 0, 1), (Expectation, Braket.Observables.Z(), 0)]) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + @test res.values[1] ≈ 0.0 atol=1e-12 + end -@testset "Expectation" begin - c = Circuit([(H, 0), (CNot, 0, 1), (Expectation, Braket.Observables.Z(), 0)]) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - @test res.values[1] ≈ 0.0 atol=1e-12 -end - -@testset "Variance" begin - c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Variance(ci, Braket.Observables.Z(), 0)) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - @test res.values[1] ≈ 1.0 atol=1e-12 -end - -@testset "Probability" begin - c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Probability(ci)) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - @test res.values[1] ≈ [0.5, 0.0, 0.0, 0.5] atol=1e-12 -end + @testset "Variance" begin + c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Variance(ci, Braket.Observables.Z(), 0)) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + @test res.values[1] ≈ 1.0 atol=1e-12 + end -@testset "DensityMatrix" begin - c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->DensityMatrix(ci)) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - ρ = zeros(ComplexF64, 4, 4) - ρ[1, 1] = 0.5 - ρ[1, 4] = 0.5 - ρ[4, 1] = 0.5 - ρ[4, 4] = 0.5 - for i in 1:4, j in 1:4 - @test res.values[1][i][j] ≈ ρ[i,j] atol=1e-12 - end -end - -@testset "Results type" begin - c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Expectation(ci, Braket.Observables.Z(), 0)) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - @test res.values[1] ≈ 0.0 atol=1e-12 - @test sprint(show, res) == "GateModelQuantumTaskResult\n" -end + @testset "Probability" begin + c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Probability(ci)) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + @test res.values[1] ≈ [0.5, 0.0, 0.0, 0.5] atol=1e-12 + end -@testset "Observable on all qubits" begin - c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Expectation(ci, Braket.Observables.Z())) - dev = LocalSimulator() - res = result(run(dev, c, shots=0)) - @test pyconvert(Vector{Float64}, res.values[1]) ≈ [0.0, 0.0] atol=1e-12 -end + @testset "DensityMatrix" begin + c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->DensityMatrix(ci)) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + ρ = zeros(ComplexF64, 4, 4) + ρ[1, 1] = 0.5 + ρ[1, 4] = 0.5 + ρ[4, 1] = 0.5 + ρ[4, 4] = 0.5 + for i in 1:4, j in 1:4 + @test res.values[1][i][j] ≈ ρ[i,j] atol=1e-12 + end + end -@testset "Unitary" begin - c = Circuit() - ccz_mat = Matrix(Diagonal(ones(ComplexF64, 2^3))) - ccz_mat[end,end] = -one(ComplexF64) - H(c, [0, 1, 2]) - Unitary(c, [0, 1, 2], ccz_mat) - H(c, [0, 1, 2]) - StateVector(c) - dev = LocalSimulator() - braket_task = run(dev, c, shots=0) - res = result(braket_task) - @test res.values[1] ≈ ComplexF64[0.75, 0.25, 0.25, -0.25, 0.25, -0.25, -0.25, 0.25] atol=1e-12 -end + @testset "Results type" begin + c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Expectation(ci, Braket.Observables.Z(), 0)) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + @test res.values[1] ≈ 0.0 atol=1e-12 + @test sprint(show, res) == "GateModelQuantumTaskResult\n" + end -@testset "PyCircuit" begin - c = Circuit() - ccz_mat = Matrix(Diagonal(ones(ComplexF64, 2^3))) - ccz_mat[end,end] = -one(ComplexF64) - H(c, [0, 1, 2]) - Unitary(c, [0, 1, 2], ccz_mat) - H(c, [0, 1, 2]) - StateVector(c) - pc = PyCircuit(c) - @test qubit_count(c) == PythonCall.pyconvert(Int, pc.qubit_count) - @test length(c.result_types) == PythonCall.pyconvert(Int, length(getproperty(pc, "_result_types"))) -end + @testset "Observable on all qubits" begin + c = Circuit() |> (ci->H(ci, 0)) |> (ci->CNot(ci, 0, 1)) |> (ci->Expectation(ci, Braket.Observables.Z())) + dev = LocalSimulator() + res = result(run(dev, c, shots=0)) + @test pyconvert(Vector{Float64}, res.values[1]) ≈ [0.0, 0.0] atol=1e-12 + end -@testset "Conversion of result types" begin - @testset for (rt, ir_rt) in ((Braket.Expectation, Braket.IR.Expectation), - (Braket.Variance, Braket.IR.Variance), - (Braket.Sample, Braket.IR.Sample)) - o = Braket.Observables.H() - obs = rt(o, [0]) - py_obs = Py(obs) - ir_obs = ir(obs) - @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) - @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs - end - @testset for (rt, ir_rt) in ((Braket.Probability, Braket.IR.Probability), - (Braket.DensityMatrix, Braket.IR.DensityMatrix)) - obs = rt([0]) - py_obs = Py(obs) - ir_obs = ir(obs) - @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) - @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs + @testset "Unitary" begin + c = Circuit() + ccz_mat = Matrix(Diagonal(ones(ComplexF64, 2^3))) + ccz_mat[end,end] = -one(ComplexF64) + H(c, [0, 1, 2]) + Unitary(c, [0, 1, 2], ccz_mat) + H(c, [0, 1, 2]) + StateVector(c) + dev = LocalSimulator() + braket_task = run(dev, c, shots=0) + res = result(braket_task) + @test res.values[1] ≈ ComplexF64[0.75, 0.25, 0.25, -0.25, 0.25, -0.25, -0.25, 0.25] atol=1e-12 + end - obs = rt() - py_obs = Py(obs) - ir_obs = ir(obs) - @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) - @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs - end - @testset "(rt, ir_rt) = (Amplitude, IR.Amplitude)" begin - obs = Braket.Amplitude(["0000"]) - py_obs = Py(obs) - ir_obs = ir(obs) - @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) - @test pyconvert(Braket.IR.Amplitude, Py(ir_obs)) == ir_obs - end - @testset "(rt, ir_rt) = (StateVector, IR.StateVector)" begin - obs = Braket.StateVector() - py_obs = Py(obs) - ir_obs = ir(obs) - @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) - @test pyconvert(Braket.IR.StateVector, Py(ir_obs)) == ir_obs - end - @testset "HermitianObservable and TensorProduct" begin - m = [1. -im; im -1.] - ho = Braket.Observables.HermitianObservable(kron(m, m)) - tp = Braket.Observables.TensorProduct([ho, ho]) - py_obs = Py(tp) - @test pyisinstance(py_obs, PyBraket.braketobs.TensorProduct) - @test pyisinstance(py_obs.factors[0], PyBraket.braketobs.Hermitian) - end -end + @testset "PyCircuit" begin + c = Circuit() + ccz_mat = Matrix(Diagonal(ones(ComplexF64, 2^3))) + ccz_mat[end,end] = -one(ComplexF64) + H(c, [0, 1, 2]) + Unitary(c, [0, 1, 2], ccz_mat) + H(c, [0, 1, 2]) + StateVector(c) + pc = PyCircuit(c) + @test qubit_count(c) == PythonCall.pyconvert(Int, pc.qubit_count) + @test length(c.result_types) == PythonCall.pyconvert(Int, length(getproperty(pc, "_result_types"))) + end + if ir_type == :JAQCD + @testset "Conversion of result types" begin + @testset for (rt, ir_rt) in ((Braket.Expectation, Braket.IR.Expectation), + (Braket.Variance, Braket.IR.Variance), + (Braket.Sample, Braket.IR.Sample)) + o = Braket.Observables.H() + obs = rt(o, [0]) + py_obs = Py(obs) + ir_obs = ir(obs) + @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) + @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs + end + @testset for (rt, ir_rt) in ((Braket.Probability, Braket.IR.Probability), + (Braket.DensityMatrix, Braket.IR.DensityMatrix)) + obs = rt([0]) + py_obs = Py(obs) + ir_obs = ir(obs) + @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) + @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs -@testset "FreeParameter" begin - α = FreeParameter(:alpha) - θ = FreeParameter(:theta) - circ = Circuit() - circ = H(circ, 0) - circ = Rx(circ, 1, α) - circ = Ry(circ, 0, θ) - circ = Probability(circ) - new_circ = circ(theta=2.0, alpha=1.0) - non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 2.0)) |> Probability - py_c1 = PyCircuit(circ) - py_c2 = PyCircuit(non_para_circ) - @test py_c2 == py_c1(theta = 2.0, alpha=1.0) - - non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 1.0)) |> Probability - py_c2 = PyCircuit(non_para_circ) - @test py_c2 == py_c1(1.0) - @testset "running with inputs" begin - dev = LocalSimulator() - oq3_circ = ir(circ, Val(:OpenQASM)) - braket_task = run(dev, oq3_circ, shots=0, inputs=Dict(string(α)=>1.0, string(θ)=>2.0)) - res = result(braket_task) - non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 2.0)) |> Probability - oq3_circ2 = ir(non_para_circ, Val(:OpenQASM)) - braket_task2 = run(dev, oq3_circ2, shots=0) - res2 = result(braket_task2) - @test res.values == res2.values + obs = rt() + py_obs = Py(obs) + ir_obs = ir(obs) + @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) + @test pyconvert(ir_rt, Py(ir_obs)) == ir_obs + end + @testset "(rt, ir_rt) = (Amplitude, IR.Amplitude)" begin + obs = Braket.Amplitude(["0000"]) + py_obs = Py(obs) + ir_obs = ir(obs) + @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) + @test pyconvert(Braket.IR.Amplitude, Py(ir_obs)) == ir_obs + end + @testset "(rt, ir_rt) = (StateVector, IR.StateVector)" begin + obs = Braket.StateVector() + py_obs = Py(obs) + ir_obs = ir(obs) + @test pyisTrue(py_obs.to_ir() == Py(ir_obs)) + @test pyconvert(Braket.IR.StateVector, Py(ir_obs)) == ir_obs + end + @testset "HermitianObservable and TensorProduct" begin + m = [1. -im; im -1.] + ho = Braket.Observables.HermitianObservable(kron(m, m)) + tp = Braket.Observables.TensorProduct([ho, ho]) + py_obs = Py(tp) + @test pyisinstance(py_obs, PyBraket.braketobs.TensorProduct) + @test pyisinstance(py_obs.factors[0], PyBraket.braketobs.Hermitian) + end + @testset "Coefficient handling" begin + o = 3.0 * Braket.Observables.Z() + py_obs = Py(o) + @test pyisinstance(py_obs, PyBraket.braketobs.Z) + @test pyisTrue(py_obs.coefficient == 3.0) + end + @testset "Sum" begin + m = [1. -im; im -1.] + ho = Braket.Observables.HermitianObservable(kron(m, m)) + tp = Braket.Observables.TensorProduct([Braket.Observables.HermitianObservable(m), Braket.Observables.HermitianObservable(m)]) + py_obs = Py(2.0*tp + ho) + @test pyisinstance(py_obs, PyBraket.braketobs.Sum) + @test pyisinstance(py_obs.summands[0], PyBraket.braketobs.TensorProduct) + @show py_obs.summands[0].coefficient + @test pyconvert(Float64, py_obs.summands[0].coefficient) == 2.0 + @test pyisinstance(py_obs.summands[1], PyBraket.braketobs.Hermitian) + end + end + end + @testset "FreeParameter" begin + α = FreeParameter(:alpha) + θ = FreeParameter(:theta) + circ = Circuit() + circ = H(circ, 0) + circ = Rx(circ, 1, α) + circ = Ry(circ, 0, θ) + circ = Probability(circ) + new_circ = circ(theta=2.0, alpha=1.0) + non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 2.0)) |> Probability + py_c1 = PyCircuit(circ) + py_c2 = PyCircuit(non_para_circ) + @test py_c2 == py_c1(theta = 2.0, alpha=1.0) + + non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 1.0)) |> Probability + py_c2 = PyCircuit(non_para_circ) + @test py_c2 == py_c1(1.0) + @testset "running with inputs" begin + dev = LocalSimulator() + oq3_circ = ir(circ, Val(:OpenQASM)) + braket_task = run(dev, oq3_circ, shots=0, inputs=Dict(string(α)=>1.0, string(θ)=>2.0)) + res = result(braket_task) + non_para_circ = Circuit() |> (ci->H(ci, 0)) |> (ci->Rx(ci, 1, 1.0)) |> (ci->Ry(ci, 0, 2.0)) |> Probability + oq3_circ2 = ir(non_para_circ, Val(:OpenQASM)) + braket_task2 = run(dev, oq3_circ2, shots=0) + res2 = result(braket_task2) + @test res.values == res2.values + end + end end + Braket.IRType[] = :OpenQASM end diff --git a/PyBraket/test/gates.jl b/PyBraket/test/gates.jl index ac201968..76bb6667 100644 --- a/PyBraket/test/gates.jl +++ b/PyBraket/test/gates.jl @@ -2,6 +2,7 @@ using Test, PyBraket, Braket, Braket.IR, PythonCall using PythonCall: Py, pyconvert, pyisTrue @testset "Gates" begin + Braket.IRType[] = :JAQCD @testset for (gate, ir_gate) in ((H, IR.H), (Braket.I, IR.I), (X, IR.X), @@ -77,4 +78,5 @@ using PythonCall: Py, pyconvert, pyisTrue @test pyisTrue(py_n.to_ir([0]) == Py(ir_n)) @test pyconvert(IR.Unitary, Py(ir_n)) == ir_n end -end \ No newline at end of file + Braket.IRType[] = :OpenQASM +end diff --git a/PyBraket/test/integ_tests/create_local_quantum_job.jl b/PyBraket/test/integ_tests/create_local_quantum_job.jl index f5f5c4d1..5eae8cc8 100644 --- a/PyBraket/test/integ_tests/create_local_quantum_job.jl +++ b/PyBraket/test/integ_tests/create_local_quantum_job.jl @@ -1,4 +1,4 @@ -using Braket, JSON3, PyBraket, PythonCall, Test +using Braket, Braket.JSON3, PyBraket, PythonCall, Test using PythonCall: pybuiltins @@ -127,4 +127,4 @@ using PythonCall: pybuiltins end end end -end \ No newline at end of file +end diff --git a/PyBraket/test/noise.jl b/PyBraket/test/noise.jl index dbf378ed..e4423483 100644 --- a/PyBraket/test/noise.jl +++ b/PyBraket/test/noise.jl @@ -28,6 +28,7 @@ using PythonCall: Py, pyconvert, pyisTrue @test length(circ.moments) == orig_len @test length(ts) == orig_len @testset "Conversions to/from Python" begin + Braket.IRType[] = :JAQCD @testset for (noise, ir_noise) in ( (BitFlip, IR.BitFlip), (PhaseFlip, IR.PhaseFlip), (AmplitudeDamping, IR.AmplitudeDamping), @@ -81,5 +82,6 @@ using PythonCall: Py, pyconvert, pyisTrue @test pyisTrue(py_n.to_ir([0]) == Py(ir_n)) @test pyconvert(IR.Kraus, Py(ir_n)) == ir_n end + Braket.IRType[] = :OpenQASM end end diff --git a/docs/src/index.md b/docs/src/index.md index a6ab5c25..e549e659 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -10,7 +10,7 @@ isn't answered here, it may be at one of those two locations. This is *experimental* software, and support may be discontinued in the future. For a fully supported SDK, please use the [Python SDK](https://github.com/aws/amazon-braket-sdk-python). We may change, remove, or deprecate parts of the API when making new releases. -Please review the [CHANGELOG](CHANGELOG.md) for information about changes in each release. +Please review the [CHANGELOG](https://github.com/awslabs/Braket.jl/blob/main/CHANGELOG.md) for information about changes in each release. ## What is Amazon Braket? diff --git a/docs/src/observables.md b/docs/src/observables.md index cde4f764..05128678 100644 --- a/docs/src/observables.md +++ b/docs/src/observables.md @@ -9,4 +9,5 @@ Braket.Observables.H Braket.Observables.I Braket.Observables.TensorProduct Braket.Observables.HermitianObservable +Braket.Observables.Sum ``` diff --git a/docs/src/results.md b/docs/src/results.md index b84d27ba..48bc2a77 100644 --- a/docs/src/results.md +++ b/docs/src/results.md @@ -1,6 +1,9 @@ # Results ```@docs Braket.Result +AdjointGradient +AdjointGradient(::Observable, ::Vector{QubitSet}, ::Vector) +AdjointGradient(::Circuit, ::Observable, ::Vector{QubitSet}, ::Vector{String}) Expectation Expectation(::Any, ::Any) Expectation(::Circuit, ::Any, ::Any) diff --git a/examples/adjoint_gradient.jl b/examples/adjoint_gradient.jl new file mode 100644 index 00000000..c1945f70 --- /dev/null +++ b/examples/adjoint_gradient.jl @@ -0,0 +1,574 @@ +### A Pluto.jl notebook ### +# v0.19.16 + +using Markdown +using InteractiveUtils + +# ╔═╡ 3d3c1591-9bbc-4a2d-b375-049f7d91fe21 +begin + using Pkg + Pkg.develop("Braket") +end + +# ╔═╡ 492616b4-0cea-41fe-8a8e-11c6d60200cb +using Braket, Plots, Optim, Graphs, SimpleWeightedGraphs, Distributions, GraphPlot, SparseArrays, LinearAlgebra + +# ╔═╡ 471998e8-763b-11ed-0be4-e5293a2b29d4 +md" +# Gradient computations using the adjoint differentiation method on Amazon Braket. + +In this notebook, we'll introduce the adjoint differentiation method and show how you can use it to accelerate variational workflows running on Braket. +" + +# ╔═╡ 5a45f736-f603-4f01-b3ca-1d415d7864c4 +md""" +## Introduction to gradient computation + +First we need to develop the concept of a \"gradient\". A [gradient](https://en.wikipedia.org/wiki/Gradient) refers to a vector derivative of a scalar-valued function of multiple variables. If we have a function $f(x)$, which depends only on $x$ and maps $x \to \mathbb{R}$ (maps $x$ to a single real number), then $f$'s gradient is just its derivative with respect to $x$: $\frac{df}{dx}$. The gradient is denoted $\nabla$, so that $\nabla f(x) = \frac{df}{dx}$. However, if $f$ is a function of multiple variables, mapping $\mathbb{R}^n \to \mathbb{R}$, we must take partial derivatives with respect to each variable. For example: + +```math +\nabla f(x, y, z) = \left[\frac{\partial f}{\partial x}, \frac{\partial f}{\partial y}, \frac{\partial f}{\partial z}\right] +``` + +The gradient of \$f\$ is itself a function and can be evaluated on specific values of $x$, $y$, and $z$. In general, for a function $f$ of $n$ independent real variables, $\nabla f$ is a length $n$ vector. + +Gradients are of interest to us because many quantum algorithms -- hybrid classical-quantum algorithms such as the quantum approximate optimization algorithm (QAOA) or the variational quantum eigensolver (VQE) especially -- can be formulated as a problem of optimizing parameters (i.e. variables) in a quantum circuit with respect to some cost function, for example an expectation value of a Hamiltonian. To efficiently perform this optimization it's common to use a gradient based optimization method, such as gradient descent (stochastic or not). An efficient means of computing gradients allows us to arrive at a good solution to the optimization problem in fewer circuit evaluations, and thus less cost. +""" + +# ╔═╡ dfba9752-1545-4d51-94b6-739954d63afc +md""" +## Computing gradients of parameters in a quantum circuit + +Let's make this a little more concrete. Suppose we have a quantum circuit which depends on a set of parameters $\vec{p}$. We can compactly represent this circuit as $U(\vec{p})$, where $U$ is the unitary that represents the action of all the gates in the circuit. Further suppose that after running this circuit, we will compute the expectation value of some operator $\hat{O}$ (for example, a Hamiltonian) and use the result to determine how good our choice of parameters $\vec{p}$ was. This situation arises often when running hybrid algorithms or quantum machine learning workflows. + + +**Note:** Although, for the sake of simplicity, we will only discuss measuring expectation values to generate the function to differentiate, one can equally well compute variances or any other scalar valued function. + +We can express this whole procedure as: +```math +f(\vec{p}) = \left\langle \psi \right| \hat{O} \left| \psi \right\rangle = \left\langle 0 U^\dagger(\vec{p}) \right| \hat{O} \left| U(\vec{p}) 0 \right\rangle +``` + + +$f(\vec{p})$ is a scalar valued function and we can compute its gradient -- all its partial derivatives with respect to the parameters $\vec{p}$ -- in order to optimize those parameters. There are a variety of methods available to compute these derivatives, three of which will be discussed below. Of the three, the adjoint differentiation method is the fastest and most frugal in circuit executions and should be preferred when available, that is, when running on a state vector simulator in exact (`shots=0`) mode. We'll introduce these other two common approaches to better understand the benefit of using the adjoint differentiation method. +""" + +# ╔═╡ a516549b-0418-4123-a654-fdd17359493e +md""" +### Finite differences + +The [finite difference method](https://en.wikipedia.org/wiki/Finite_difference) is a common technique used to approximate derivatives. Suppose we have a function $f(\vec{p})$ and we want to compute the $i$-th partial derivative of $f$, $\frac{\partial f}{\partial p_i}$. We can do so by approximating: + +```math +\frac{\partial f}{\partial p_i} \approx \frac{f(p_1, p_2, ..., p_i + h, ..., p_n) - f(p_1, p_2, ..., p_i, ..., p_n)}{h} +``` + +Where $h$ is some small real number. This formula might seem familiar from introductory calculus. The smaller $h$ is the better approximated the derivative is. By keeping all other parameters fixed, we can approximate the partial derivative with respect to $p_i$, but as we can see, computing *each* partial derivative would require *two* full circuit executions (one to compute each value of $f$). Thus, the total number of circuit executions needed to compute the gradient of $f$ for *one* set of values $\vec{p}$ would be $2n$, if the length of $p$ is $n$. + +For a quantum circuit there can be additional problems. On a real quantum device, we can't compute the exact expectation value (or variance) of a circuit. We can only run many shots, each of which is a full circuit execution, and approximate the expectation value from the measurement statistics that result. This means that for very small $h$, it may be very difficult to approximate the gradient accurately. +""" + +# ╔═╡ 202ca719-465e-4f4f-b5b2-40f953857288 +md""" +### Parameter shift rules + +Let's return to our original formula for the gradient of $f$: + +```math +\nabla f(\vec{p}) = \left(\frac{\partial f}{\partial p_1}, \ldots , \frac{\partial f}{\partial p_n}\right) +``` + +and examine one of the vector elements a little more closely: + +```math +\frac{\partial f}{\partial p_i} = \frac{\partial}{\partial p_i} \left\langle 0 U^\dagger(\vec{p}) \right| \hat{O} \left| U(\vec{p}) 0 \right\rangle = \frac{\partial}{\partial p_i} \left\langle 0 \right| U^\dagger(\vec{p}) \hat{O} U(\vec{p}) \left| 0 \right\rangle +``` + +We can pull the derivative operator inside the expectation value so that: + +```math +\frac{\partial f}{\partial p_i} = \left\langle 0 \left|\frac{\partial}{\partial p_i} \left( U^\dagger(\vec{p}) \hat{O} U(\vec{p})\right) \right| 0 \right\rangle +``` + +We'll assume that each gate depends on at most one parameter, and each parameter appears in only one gate. What if we have repeated parameters? We can write down a mapping of each repeated parameter to a unique copy and sum the derivatives of those copies at the end using the [product rule](https://en.wikipedia.org/wiki/Product_rule). But for now, for simplicity we will assume that each parameter appears only once and each gate has at most one parameter. Further we'll state that gate $i$ is associated with the $i$-th parameter (every gate has a parameter). If non-parametrized gates are present, we can contract them into parametrized gates to achieve this, or assign them constant parameters (remember, the derivative of a constant is always 0). + +We can write that the overall circuit unitary $U$ is a product of individual gates: + +```math +U(\vec{p}) = \otimes_{i=1}^N U_{i}(p_i) +``` + +if there are $N$ gates in the circuit. + +Then, using the product rule, we can write: + +```math +\frac{\partial f}{\partial p_i} = \left\langle 0 \left| \otimes_{j=1}^{i-1} U^\dagger_j \otimes \frac{\partial U^\dagger_i(p_i)}{\partial p_i} \otimes_{j=i+1}^{n} U^\dagger_j \hat{O} U(\vec{p}) + U^\dagger(\vec{p}) \hat{O}\otimes_{j=i+1}^{n} U_j \otimes \frac{\partial U_i(p_i)}{\partial p_i}\otimes_{j=1}^{i-1} U_j \right| 0 \right\rangle +``` + +We can absorb the non-differentiated products so that: + +```math +\frac{\partial f}{\partial p_i} = \left\langle \phi \left| \frac{\partial U^\dagger_i(p_i)}{\partial p_i} \hat{\mathcal{O}}U_i(p_i) + U^\dagger_i(p_i)\hat{\mathcal{O}}\frac{\partial U_i(p_i)}{\partial p_i}\right| \phi \right\rangle +``` + +where + +```math +\left|\phi\right\rangle = \otimes_{j=1}^{i-1} U_j \left| 0 \right\rangle +``` + +and + +```math +\hat{\mathcal{O}} = \otimes_{j=i+1}^{n} U^\dagger_j \hat{O} \otimes_{j=i+1}^{n} U_j +``` + +Now we can see that + +```math +\frac{\partial U^\dagger_i(p_i)}{\partial p_i} \hat{\mathcal{O}}U_i(p_i) + U^\dagger_i(p_i)\hat{\mathcal{O}}\frac{\partial U_i(p_i)}{\partial p_i} = \frac{\partial}{\partial p_i} \left( U_i^\dagger(p_i) \hat{\mathcal{O}} U_i(p_i)\right) +``` + +so, in sum: + +```math +\frac{\partial f}{\partial p_i} = \left\langle \phi \left|\frac{\partial}{\partial p_i} \left( U_i^\dagger(p_i) \hat{\mathcal{O}} U_i(p_i)\right)\right| \phi \right\rangle +``` + +and in many cases (but not all!) we can define a *shift* $s$ such that: + +```math +\frac{\partial}{\partial p_i} \left( U_i^\dagger(p_i) \hat{\mathcal{O}} U_i(p_i)\right) = U_i^\dagger(p_i + s) \hat{\mathcal{O}} U_i(p_i + s) - U_i^\dagger(p_i - s) \hat{\mathcal{O}} U_i(p_i - s) +``` + +Thus the name "parameter shift". What makes this different from the finite differences method is that $s$ is not necessarily small. Detailed guides to choosing shifts and identifying which gates support the method can be found in Refs. [1](https://arxiv.org/abs/1803.00745) and [2](https://arxiv.org/abs/1811.11184). If a gate does *not* support a parameter shift rule, we can always fall back to the finite difference method. + +We can see that the parameter shift method *also* requires two circuit executions to compute the partial derivative of each parametrized gate. The advantage over finite difference is in numerical accuracy. Parameter shift can be used both when `shots=0` or when `shots>0`. Since the introduction of the method, many extensions and generalizations have been published, including Refs. [3](https://arxiv.org/abs/2107.12390), [4](https://arxiv.org/abs/2005.10299), and many more. +""" + +# ╔═╡ 42822b26-269e-4d9e-92a2-cd34ef2ca4be +md""" +### Adjoint differentiation + +The two methods we've examined so far, finite differences and parameter shift, both require two full circuit executions per parameter to compute the gradient. This can become very expensive, in both time and charges, for deep circuits and/or circuits with many parameters. Is there a way to compute gradients in a more "execution-frugal" way? For `shots=0`, the answer is yes. First introduced in Ref. [5](https://arxiv.org/abs/2009.02823), the adjoint differentiation method allows us to compute all partial derivatives in "1+1" circuit executions. How does it work? Recall that: + +```math +\frac{\partial f}{\partial p_i} = \left\langle 0 \left| \otimes_{j=1}^{i-1} U^\dagger_j \otimes \frac{\partial U^\dagger_i(p_i)}{\partial p_i} \otimes_{j=i+1}^{n} U^\dagger_j \hat{O} U(\vec{p}) + U^\dagger(\vec{p}) \hat{O}\otimes_{j=i+1}^{n} U_j \otimes \frac{\partial U_i(p_i)}{\partial p_i}\otimes_{j=1}^{i-1} U_j \right| 0 \right\rangle +``` + +In the adjoint method, we take a different approach to computing this derivative. We realize that: + +```math +\left( \left \langle 0 \left| \otimes_{j=1}^{i-1} U^\dagger_j \otimes \frac{\partial U^\dagger_i(p_i)}{\partial p_i} \otimes_{j=i+1}^{n} U^\dagger_j \hat{O} U(\vec{p}) \right| 0 \right\rangle \right)^\dagger = \left \langle 0 \left| U^\dagger(\vec{p}) \hat{O}\otimes_{j=i+1}^{n} U_j \otimes \frac{\partial U_i(p_i)}{\partial p_i}\otimes_{j=1}^{i-1} U_j \right| 0 \right\rangle +``` + +and thus, because all the gates are unitaries and the operator $\hat{O}$ is Hermitian, + +```math +\frac{\partial f}{\partial p_i} = 2\Re \left\langle 0 \left| U^\dagger(\vec{p}) \hat{O}\otimes_{j=i+1}^{n} U_j \otimes \frac{\partial U_i(p_i)}{\partial p_i}\otimes_{j=1}^{i-1} U_j \right| 0 \right\rangle +``` + +where $\Re$ denotes the real part. Now we can absorb some factors so that: + +```math +\frac{\partial f}{\partial p_i} = 2\Re \left\langle b_i \left| \frac{\partial U_i(p_i)}{\partial p_i} \right| k_i \right\rangle +``` + +where + +```math +\left\langle b_i \right| = \left\langle 0 \right| U^\dagger(\vec{p}) \hat{O}\otimes_{j=i+1}^{n} U_j(p_j) +``` + +and + +```math +\left | k_i \right\rangle = \otimes_{j=1}^{i-1} U_j(p_j) \left| 0 \right\rangle +``` + +The basis of the adjoint method is realizing that we can iteratively compute each partial derivative by "back stepping" through the circuit after having applied all its gates once. This is very similar to classical back propagation, if you're familiar with that technique from classical machine learning. We first apply all gates to compute $\left| k_n \right\rangle = \otimes_{j=1}^{n} U_j \left| 0 \right\rangle$, copy the state and apply $\hat{O}$ to acquire: + +```math +\left\langle b_n \right| = \left\langle 0 \right| U^\dagger(\vec{p}) \hat{O} +``` + +then compute $\frac{\partial f}{\partial p_n}$: + +```math +\frac{\partial f}{\partial p_n} = \left\langle b_n \left|\frac{\partial U_n(p_n)}{\partial p_n}\right| k_n \right\rangle +``` + +In a moment we'll address how to find $\frac{\partial U_n(p_n)}{\partial p_n}$. Once we've computed the first partial derivative, we update $\left\langle b_n\right|$ and $\left| k_n \right\rangle$ to generate: + +```math +\left\langle b_{n-1} \right| = \left\langle b_n \right| U_n(p_n) \\ +\left | k_{n-1} \right\rangle = U^\dagger_{n-1} \left| k_n \right\rangle +``` + +By iteratively updating these two states, we can compute all partial derivatives with only one circuit execution plus one "back step" execution, significantly less than what is required by finite differences or parameter shift. The cost is that there is additional memory overhead, as we have to store an additional state vector and compute a third in the expectation value $\left\langle b_i \left|\frac{\partial U_i(p_i)}{\partial p_i}\right| k_i \right\rangle$. + +How do we compute the derivative $\frac{\partial U_i(p_i)}{\partial p_i}$? In many cases, if $U_i(p_i)$ is continually differentiable with respect to $p_i$, we can simply take a matrix derivative. In particular, many parametrizable gates can be written as exponentials of Paulis, so that: + +```math +\frac{\partial U_i(p_i)}{\partial p_i} = \frac{\partial}{\partial p_i}\exp\left\{ i c p_i \hat{P}\right\} = i c \hat{P} \exp\left\{ i c p_i \hat{P}\right\} +``` + +where $c$ is some constant, $\hat{P}$ is some Pauli gate, and $i$ is the imaginary number $\sqrt{-1}$. This is easily generalizable to exponents of sums of Paulis through the [chain rule](https://en.wikipedia.org/wiki/Chain_rule). In cases where $U(p_i)$ is **not** continuously differentiable, the derivative can be computed numerically, e.g. through finite differences as discussed above. + +**Note:** Because it is formulated **only** for exact computations, the adjoint method can only be used on simulators, such as SV1, when running with `shots=0`. + + +The adjoint differentiation method is available through the `AdjointGradient` result type on Amazon Braket, which we'll introduce in the next section. With this result type, all derivatives are computed using the adjoint differentiation method. +""" + +# ╔═╡ df640675-a439-4074-a6d9-df4f8289e446 +md""" +## The `AdjointGradient` result type + +Amazon Braket now supports a result type, `AdjointGradient`, which allows you to conveniently compute gradients of free parameters with respect to the expectation value of some observable on your circuits. + +**Note:** Currently, the `AdjointGradient` result type is **only** supported on SV1 when running in `shots=0` mode. All derivatives are computed using the adjoint differentiation method. + +Let's see an example of this result type in action: +""" + +# ╔═╡ a86cabc6-77ba-4991-a02f-613219e8984c +device_arn = "arn:aws:braket:::device/quantum-simulator/amazon/sv1" + +# ╔═╡ 9a9f4565-f094-4352-b54d-fafd391a8e88 +device = AwsDevice(device_arn) + +# ╔═╡ d57a24ed-1b8d-4195-bdd6-7ab8969627bd +md""" +We can prepare a simple parametrized circuit and compute its gradient with respect to some observable. Note that you supply the observable to the `AdjointGradient` result type. Supported observables are: + - Any of `Observables.Z()`, `Observables.X()`, `Observables.Y()`, `Observables.H()`, or `Observable.I()` + - A `TensorProduct` + - A `HermitianObservable` + - A `Sum` + +You can also supply the list of parameters to compute partial derivatives with respect to. If a parameter is present in the circuit, but not in the `parameters` argument to `adjoint_gradient`, its corresponding partial derivative will not be computed. If the list of `parameters` is empty, the gradient will be computed with respect to all free parameters present in the circuit. +""" + +# ╔═╡ 8eea449b-2935-4aa4-9eaf-70d0c0af5c3d +begin + θ = FreeParameter(:theta) + γ = FreeParameter(:gamma) + ag_circuit = Circuit([(H, 0), (CNot, 0, 1), (Rx, 0, θ), (Rx, 1, θ), (XX, 0, 1, γ)]) + # add the adjoint gradient result type + + # we can either explicitly provide a list of parameters + #ag_circuit(AdjointGradient, Braket.Observables.Z() * Braket.Observables.Z(), [0, 1], [θ, γ]) + + # or use the default of all + ag_circuit(AdjointGradient, Braket.Observables.Z() * Braket.Observables.Z(), [0, 1], []) +end + +# ╔═╡ 7621e93a-a14a-47e0-9d49-9ce24ce0b9e0 +md""" +Now we can compute the gradient of the circuit with respect to our two free parameters for a given set of parameter values, which we supply to the `device` functor with the `inputs` argument: +""" + +# ╔═╡ 59178e2e-707c-4e1f-b2ea-54c3aa718956 +grad_1 = result(device(ag_circuit, shots=0, inputs = Dict("theta"=>0.1, "gamma"=>0.05))).values[1] + +# ╔═╡ f1279ca2-26a3-4a3d-b50a-0ef6666eca5f +grad_2 = result(device(ag_circuit, shots=0, inputs = Dict("theta"=>0.2, "gamma"=>0.1))).values[1] + +# ╔═╡ 12dc5ee9-1636-405b-b095-2a965d0ccaea +md""" +We can immediately see that although `θ` appears twice in the circuit (in two `Rx` gates), it only appears once in the result. `AdjointGradient` computes gradients **per parameter**, and **not** per-gate. We can also see that if `parameters` is empty, derivatives with respect to all free parameters will be computed. This is useful in cases when your circuit has a large number of free parameters. +""" + +# ╔═╡ 4ab8a07d-da9b-4ccc-9bce-644f7da6de21 +md""" +## Accelerating QAOA with `AdjointGradient` + +Now we can see how using the `AdjointGradient` result type can improve performance for a hybrid algorithm such as QAOA. For an introduction to QAOA, see [its example notebook](https://github.com/aws/amazon-braket-examples/blob/main/examples/hybrid_quantum_algorithms/QAOA/QAOA_braket.ipynb). We'll create a `train` function to use `AdjointGradient` and determine a Jacobian, and compare this approach with the Jacobian-free method used in the QAOA notebook. Much of the code here is further explained in that notebook, so we strongly suggest you review it before proceeding. We'll run the entire QAOA workflow in `shots=0` mode so that we can compare with `AdjointGradient`, which means we can directly compute the cost (energy). First, we set up the problem and import the circuit generator and training functions: +""" + +# ╔═╡ 4a29e5f3-5b28-41cd-af15-1a3e05ea6e83 +begin + # setup Erdos Renyi graph + n = 10 # number of nodes/vertices + m = 20 # number of edges + seed = 2 + + # define graph object + G = SimpleWeightedGraph(erdos_renyi(n, m, seed=seed)) + + # choose random weights + weights_d = Uniform(0, 1) + for edge in edges(G) + G.weights[edge.src, edge.dst] = rand(weights_d) + G.weights[edge.dst, edge.src] = rand(weights_d) + end +end + +# ╔═╡ c3edb0ce-ecef-49a7-a799-b288a3c53d91 +gplot(G) + +# ╔═╡ ba3d585e-5118-480d-b09f-bbd2616022d3 +# build Ising matrix corresponding to G +Jfull = triu(adjacency_matrix(G), 1) + +# ╔═╡ f3afc1e5-b9ca-4938-9fd6-33c3a6a1fc87 +n_qubits = size(Jfull, 1) + +# ╔═╡ 8f0e6be9-440b-44af-8a4f-f8c4aefc2120 +md""" +We can now specify hyperparameters for the training. We can choose the number of QAOA layers to use, the maximum iteration number, and the classical optimization algorithm to use. +""" + +# ╔═╡ 29039bc9-f52f-4c3f-8457-8ad79a2ebc7a +DEPTH = 2 # circuit depth for QAOA + +# ╔═╡ c4a122e2-a94e-4eaa-8c9b-44719b3be8d5 +OPT_METHOD = BFGS() + +# ╔═╡ 75d93548-30be-4ad7-8e65-de3a8556c224 +optim_opts = Optim.Options(iterations = 30, show_trace=true, store_trace=true) + +# ╔═╡ 3fe9dbbe-3ef2-44e2-9de8-a460eb06cdb8 +md""" +We can also initialize the parameters and the initial guess for the energy (cost): +""" + +# ╔═╡ 18d6f499-4c43-4448-ba1f-ccc70761c4cf +begin + params_d = Uniform(0, 2π) + γ_initial = rand(params_d, DEPTH) + β_initial = rand(params_d, DEPTH) + params_initial = vcat(γ_initial, β_initial) + # initialize reference solution (simple guess) + energy_init = 0.0 +end + +# ╔═╡ 301001e8-f48b-4198-aa70-d75b475dda5b +md" +Now we build up the QAOA circuit itself. For an introduction to QAOA, check out the main [Amazon Braket example notebook](https://github.com/aws/amazon-braket-examples/blob/main/examples/hybrid_quantum_algorithms/QAOA/QAOA_braket.ipynb). +" + +# ╔═╡ a948c2fe-4c00-4b0d-a3fa-ae8ec0f02cff +""" +Returns circuit for driver Hamiltonian U(Hb, β) +""" +function driver(β, n_qubits) + # instantiate circuit object + circ = Circuit() + # apply parametrized rotation around x to every qubit + for qubit in 0:n_qubits-1 + circ(Rx, qubit, β) + end + return circ +end + +# ╔═╡ fba6c5c9-d6af-4a04-8bfa-dd83a1ca99dc +""" +Returns circuit for evolution with cost Hamiltonian. +""" +function cost_circuit(γ, n_qubits, ising) + # instantiate circuit object + circ = Circuit() + + # get all non-zero entries (edges) from Ising matrix + Is, Js, Vs = findnz(ising) + edges = zip(Is, Js) + + # apply ZZ gate for every edge (with corresponding interaction strength) + for (ii, qubit_pair) in enumerate(edges) + circ(ZZ, qubit_pair[1], qubit_pair[2], γ[ii]) + end + return circ +end + +# ╔═╡ 8d011c5e-7924-455b-aac0-63afa674a562 +""" +Returns full QAOA circuit of provided depth (inferred from length of `params` argument). +""" +function qaoa_builder(params, n_qubits, ising) + # initialize qaoa circuit with first Hadamard layer: + # for minimization start in |-> + circ = Circuit([(H, collect(0:n_qubits-1))]) + + # setup two parameter families + circuit_length = div(length(params), 2) + γs = params[1:circuit_length] + βs = params[circuit_length+1:end] + + # add QAOA circuit layer blocks + for mm in 1:circuit_length + append!(circ, cost_circuit(γs[mm], n_qubits, ising)) + append!(circ, driver(βs[mm], n_qubits)) + end + return circ +end + +# ╔═╡ f941586e-2b88-4f7e-8f5f-a0f93ebe8eee +""" +Builds a `Sum` observable representing the cost Hamiltonian. +""" +function cost_H(ising) + Is, Js, Vs = findnz(ising) + edges = collect(zip(Is, Js)) + + H_terms = [] + # apply ZZ gate for every edge (with corresponding interaction strength) + for qubit_pair in edges[2:end] + # get interaction strength from Ising matrix + int_strength = ising[qubit_pair[1], qubit_pair[2]] + push!(H_terms, 2*int_strength * Braket.Observables.Z() * Braket.Observables.Z()) + end + targets = [QubitSet([edge[1], edge[2]]) for edge in edges] + H_0 = 2*ising[edges[1][1], edges[1][2]] * Braket.Observables.Z() * Braket.Observables.Z() + sum_H = sum(H_terms, init=H_0) + return sum_H, targets +end + +# ╔═╡ 0656f232-f00e-4da8-a6a1-2c7156d9d840 +""" +Applies the scaling terms to the raw derivatives provided by `AdjointGradient` in order to generate the correct Jacobian. +""" +function form_jacobian(n_params, gradient, ising) + jac = zeros(Float64, n_params) + Is, Js, Vs = findnz(ising) + edges = collect(zip(Is, Js)) + split = div(n_params, 2) + for i in 1:split + # handle βs + jac[split + i] += 2 * gradient[Symbol("beta_$i")] + # handle γs + for j in 1:length(edges) + int_strength = ising[edges[j][1], edges[j][2]] + jac[i] += 2 * int_strength * gradient[Symbol("gamma_$(i)_$(j)")] + end + end + return jac +end + +# ╔═╡ d9505af6-18c8-4381-9574-3c03645a8a11 +""" +Creates the dictionary of parameter values to use when running the task. +""" +function form_inputs_dict(params, ising) + n_params = length(params) + params_dict = Dict{String, Float64}() + Is, Js, Vs = findnz(ising) + edges = collect(zip(Is, Js)) + split = div(n_params, 2) + for i in 1:split + params_dict["beta_$i"] = 2 * params[split + i] + for j in 1:length(edges) + int_strength = ising[edges[j][1], edges[j][2]] + params_dict["gamma_$(i)_$(j)"] = 2 * int_strength * params[i] + end + end + return params_dict +end + +# ╔═╡ 14fa1f6d-2262-49f0-bc74-05e3974a29ca +md" +Now we have everything we need to define the objective function to pass to `optimize`. +" + +# ╔═╡ a277af58-6d4b-4389-8821-5824cc6e71b0 +begin + mutable struct AdjointClosure + qaoa_circuit::Circuit + ising::SparseMatrixCSC + device::AwsDevice + end + """ + Computes the cost and gradient in one call to SV1. + """ + function (oc::AdjointClosure)(F, G, params) + # create parameter dict + params_dict = form_inputs_dict(params, oc.ising) + # classically simulate the circuit + # set the parameter values using the inputs argument + timeout = 3 * 24 * 60 * 60 + task = oc.device(oc.qaoa_circuit, shots=0, inputs=params_dict, poll_timeout_seconds=timeout) + + # get result for this task + res = result(task) + gradient = res.values[1] + energy = gradient[:expectation] + jac = form_jacobian(length(params), gradient[:gradient], oc.ising) + if G != nothing + copyto!(G, jac) + end + if F != nothing + return energy + end + end + Base.show(io::IO, ac::AdjointClosure) = print(io, "AdjointClosure") +end + +# ╔═╡ 2c99922f-86cb-462d-84c4-a4d7793f0716 +md" +We build a QAOA circuit which will use the `AdjointGradient` result type. +" + +# ╔═╡ 22ebd728-520f-4c85-9198-1c15e1366058 +begin + γ_params = [[FreeParameter(Symbol("gamma_$(i)_$(j)")) for j in 1:nnz(Jfull)] for i in 1:DEPTH] + β_params = [FreeParameter(Symbol("beta_$i")) for i in 1:DEPTH] + params = vcat(γ_params, β_params) + qaoa_circ = qaoa_builder(params, n_qubits, Jfull) + + Hamiltonian, targets = cost_H(Jfull) + qaoa_circ(AdjointGradient, Hamiltonian, targets, []) +end + +# ╔═╡ c9d8c32f-4980-4e8c-942a-dce5136b4ad2 +ac = AdjointClosure(qaoa_circ, Jfull, device) + +# ╔═╡ 591ae7b2-3fda-4301-b1d1-a97817246481 +md""" +We can let `Optim.jl` know that we are able to compute the cost and gradient in one fell swoop by using [`only_fg!`](https://julianlsolvers.github.io/Optim.jl/stable/#user/tipsandtricks/#avoid-repeating-computations). Note that running the code below will submit tasks to an Amazon Braket on-demand simulator, and you may incur charges to your account. For that reason, the cell is commented out by default. Un-comment it if you wish to run the cell. +""" + +# ╔═╡ Cell order: +# ╟─471998e8-763b-11ed-0be4-e5293a2b29d4 +# ╟─5a45f736-f603-4f01-b3ca-1d415d7864c4 +# ╟─dfba9752-1545-4d51-94b6-739954d63afc +# ╟─a516549b-0418-4123-a654-fdd17359493e +# ╟─202ca719-465e-4f4f-b5b2-40f953857288 +# ╟─42822b26-269e-4d9e-92a2-cd34ef2ca4be +# ╟─df640675-a439-4074-a6d9-df4f8289e446 +# ╠═3d3c1591-9bbc-4a2d-b375-049f7d91fe21 +# ╠═492616b4-0cea-41fe-8a8e-11c6d60200cb +# ╠═a86cabc6-77ba-4991-a02f-613219e8984c +# ╠═9a9f4565-f094-4352-b54d-fafd391a8e88 +# ╟─d57a24ed-1b8d-4195-bdd6-7ab8969627bd +# ╠═8eea449b-2935-4aa4-9eaf-70d0c0af5c3d +# ╟─7621e93a-a14a-47e0-9d49-9ce24ce0b9e0 +# ╠═59178e2e-707c-4e1f-b2ea-54c3aa718956 +# ╠═f1279ca2-26a3-4a3d-b50a-0ef6666eca5f +# ╟─12dc5ee9-1636-405b-b095-2a965d0ccaea +# ╟─4ab8a07d-da9b-4ccc-9bce-644f7da6de21 +# ╠═4a29e5f3-5b28-41cd-af15-1a3e05ea6e83 +# ╠═c3edb0ce-ecef-49a7-a799-b288a3c53d91 +# ╠═ba3d585e-5118-480d-b09f-bbd2616022d3 +# ╠═f3afc1e5-b9ca-4938-9fd6-33c3a6a1fc87 +# ╟─8f0e6be9-440b-44af-8a4f-f8c4aefc2120 +# ╠═29039bc9-f52f-4c3f-8457-8ad79a2ebc7a +# ╠═c4a122e2-a94e-4eaa-8c9b-44719b3be8d5 +# ╠═75d93548-30be-4ad7-8e65-de3a8556c224 +# ╟─3fe9dbbe-3ef2-44e2-9de8-a460eb06cdb8 +# ╠═18d6f499-4c43-4448-ba1f-ccc70761c4cf +# ╟─301001e8-f48b-4198-aa70-d75b475dda5b +# ╠═a948c2fe-4c00-4b0d-a3fa-ae8ec0f02cff +# ╠═fba6c5c9-d6af-4a04-8bfa-dd83a1ca99dc +# ╠═8d011c5e-7924-455b-aac0-63afa674a562 +# ╠═f941586e-2b88-4f7e-8f5f-a0f93ebe8eee +# ╠═0656f232-f00e-4da8-a6a1-2c7156d9d840 +# ╠═d9505af6-18c8-4381-9574-3c03645a8a11 +# ╟─14fa1f6d-2262-49f0-bc74-05e3974a29ca +# ╠═a277af58-6d4b-4389-8821-5824cc6e71b0 +# ╟─2c99922f-86cb-462d-84c4-a4d7793f0716 +# ╠═22ebd728-520f-4c85-9198-1c15e1366058 +# ╠═c9d8c32f-4980-4e8c-942a-dce5136b4ad2 +# ╟─591ae7b2-3fda-4301-b1d1-a97817246481 +# ╠═71da8a3a-9632-4bf9-8b45-09b6a99dec60 +# ╠═26408581-1132-4705-b401-c73343cd5c7b diff --git a/src/Braket.jl b/src/Braket.jl index 89afc303..b1ca4ee7 100644 --- a/src/Braket.jl +++ b/src/Braket.jl @@ -9,7 +9,7 @@ export logs, log_metric, metrics export depth, qubit_count, qubits, ir, IRType, OpenQASMSerializationProperties export OpenQasmProgram -export Expectation, Sample, Variance, Amplitude, Probability, StateVector, DensityMatrix, Result +export AdjointGradient, Expectation, Sample, Variance, Amplitude, Probability, StateVector, DensityMatrix, Result using AWSS3 using AWS @@ -41,21 +41,21 @@ Currently, two formats are supported: - `:JAQCD`, the Amazon Braket IR - `:OpenQASM`, the OpenQASM3 representation -By default, `IRType` is initialized to use `:JAQCD`, although this may change +By default, `IRType` is initialized to use `:OpenQASM`, although this may change in the future. The current default value can be checked by calling `IRType[]`. To change the default IR format, set `IRType[]`. # Examples ```jldoctest julia> IRType[] -:JAQCD +:OpenQASM -julia> IRType[] = :OpenQASM; +julia> IRType[] = :JAQCD; julia> IRType[] -:OpenQASM +:JAQCD -julia> IRType[] = :JAQCD; +julia> IRType[] = :OpenQASM; ``` """ const IRType = Ref{Symbol}() @@ -65,7 +65,7 @@ const GlobalTrackerContext = Ref{TrackerContext}() function __init__() AWS.DEFAULT_BACKEND[] = AWS.DownloadsBackend() - IRType[] = :JAQCD + IRType[] = :OpenQASM Prices[] = Pricing([]) GlobalTrackerContext[] = TrackerContext() end @@ -82,6 +82,7 @@ qubit_count(t) = throw(MethodError(qubit_count, t)) include("qubit_set.jl") include("raw_schema.jl") +ir(p::Union{AHSProgram, BlackbirdProgram, OpenQasmProgram}) = p using .IR include("raw_jobs_config.jl") include("raw_task_result_types.jl") diff --git a/src/circuit.jl b/src/circuit.jl index ea61d800..d4a434b5 100644 --- a/src/circuit.jl +++ b/src/circuit.jl @@ -61,6 +61,7 @@ julia> Circuit(v); ``` """ Circuit(v::AbstractVector) = (c = Circuit(); c(v); return c) +Base.show(io::IO, circ::Circuit) = print(io, "Circuit($(qubit_count(circ)) qubits)") function make_bound_circuit(c::Circuit, param_values::Dict{Symbol, Number}) clamped_circ = Circuit() @@ -179,11 +180,12 @@ julia> qubit_count(c) ``` """ qubit_count(c::Circuit) = length(qubits(c)) +qubit_count(p::Program) = length(mapreduce(ix->ix.target, union, p.instructions)) (rt::Result)(c::Circuit) = add_result_type!(c, rt) Base.convert(::Type{Circuit}, p::Program) = Circuit(Moments(p.instructions), p.instructions, Result[StructTypes.constructfrom(Result, r) for r in p.results], p.basis_rotation_instructions) -Base.convert(::Type{Program}, c::Circuit) = (basis_rotation_instructions!(c); return Program(braketSchemaHeader("braket.ir.jaqcd.program" ,"1"), c.instructions, ir.(c.result_types), c.basis_rotation_instructions)) +Base.convert(::Type{Program}, c::Circuit) = (basis_rotation_instructions!(c); return Program(braketSchemaHeader("braket.ir.jaqcd.program" ,"1"), c.instructions, ir.(c.result_types, Val(:JAQCD)), c.basis_rotation_instructions)) Circuit(p::Program) = convert(Circuit, p) Program(c::Circuit) = convert(Program, c) @@ -221,9 +223,9 @@ function ir(c::Circuit, ::Val{:OpenQASM}; serialization_properties::Serializatio end return OpenQasmProgram(header_dict[OpenQasmProgram], join(vcat(header, ixs, rts), "\n"), nothing) end -ir(c::Circuit, ::Val{:JAQCD}; kwargs...) = JSON3.write(convert(Program, c)) -ir(p::Program, ::Val{:JAQCD}; kwargs...) = JSON3.write(p) -ir(p::Program; kwargs...) = ir(p, Val(IRType[]); kwargs...) +ir(c::Circuit, ::Val{:JAQCD}; kwargs...) = convert(Program, c) +ir(p::Program, ::Val{:JAQCD}; kwargs...) = p +ir(p::Program; kwargs...) = ir(p, Val(:JAQCD); kwargs...) function Base.append!(c1::Circuit, c2::Circuit) foreach(ix->add_instruction!(c1, ix), c2.instructions) @@ -371,13 +373,17 @@ function add_to_qubit_observable_mapping!(c::Circuit, o::Observables.TensorProdu end add_to_qubit_observable_mapping!(c::Circuit, o::Observables.Observable, obs_targets::Vector) = add_to_qubit_observable_mapping!(c, o, QubitSet(obs_targets)) add_to_qubit_observable_set!(c::Circuit, rt::ObservableResult) = union!(c.qubit_observable_set, Set(rt.targets)) +add_to_qubit_observable_set!(c::Circuit, rt::AdjointGradient) = union!(c.qubit_observable_set, Set(reduce(union, rt.targets))) add_to_qubit_observable_set!(c::Circuit, rt::Result) = c.qubit_observable_set function add_result_type!(c::Circuit, rt::Result) rt_to_add = rt if rt_to_add ∉ c.result_types + if rt_to_add isa AdjointGradient && any(rt_ isa AdjointGradient for rt_ in c.result_types) + throw(ArgumentError("only one AdjointGradient result can be present on a circuit.")) + end obs = extract_observable(rt_to_add) - if !isnothing(obs) && c.observables_simultaneously_measureable + if !isnothing(obs) && c.observables_simultaneously_measureable && !(rt isa AdjointGradient) add_to_qubit_observable_mapping!(c, obs, rt_to_add.targets) end add_to_qubit_observable_set!(c, rt_to_add) @@ -388,7 +394,6 @@ end add_result_type!(c::Circuit, rt::Result, target) = add_result_type!(c, remap(rt, target)) add_result_type!(c::Circuit, rt::Result, target_mapping::Dict{<:Integer, <:Integer}) = add_result_type!(c, remap(rt, target_mapping)) - function add_instruction!(c::Circuit, ix::Instruction) to_add = [ix] if ix.operator isa QuantumOperator && Parametrizable(ix.operator) == Parametrized() @@ -709,6 +714,36 @@ julia> c = StateVector(c); ``` """ StateVector(c::Circuit) = (sv = StateVector(); sv(c)) + +""" + AdjointGradient(c::Circuit, o::Observable, targets, parameters) -> Circuit + +Constructs an `AdjointGradient` computation with respect to the expectation value of an observable `o` +on qubits `targets`, computing partial derivatives of `parameters`, and adds it as a result to [`Circuit`](@ref) `c`. + +`o` may be any `Observable`. `targets` must be a `Vector` of `QubitSet`s (or a single `QubitSet`, if `o` is not a `Sum`), +each of which is the same length as the qubit count of the corresponding term in `o`. +`parameters` can have elements which are [`FreeParameter`](@ref)s or `String`s, or `["all"]`, +in which case the gradient is computed with respect to all parameters in the circuit. + +# Examples +```jldoctest +julia> c = Circuit(); + +julia> α = FreeParameter("alpha"); + +julia> c = H(c, collect(0:10)); + +julia> c = Rx(c, collect(0:10), α); + +julia> c = AdjointGradient(c, Braket.Observables.Z(), 0, [α]); +``` +""" +AdjointGradient(c::Circuit, o::Observable, target::Vector{QubitSet}, parameters) = (ag = AdjointGradient(o, target, parameters); ag(c)) +AdjointGradient(c::Circuit, o::Observable, target::Vector{Vector{T}}, parameters) where {T} = (ag = AdjointGradient(o, target, parameters); ag(c)) +AdjointGradient(c::Circuit, o::Observable, target::Vector{<:IntOrQubit}, parameters) = (ag = AdjointGradient(o, target, parameters); ag(c)) +AdjointGradient(c::Circuit, o::Observable, target::QubitSet, parameters) = (ag = AdjointGradient(o, [target], parameters); ag(c)) +AdjointGradient(c::Circuit, o::Observable, target::IntOrQubit, parameters) = (ag = AdjointGradient(o, [target], parameters); ag(c)) function validate_circuit_and_shots(c::Circuit, shots::Int) isempty(c.instructions) && throw(ErrorException("Circuit must have instructions to run on a device.")) @@ -716,7 +751,7 @@ function validate_circuit_and_shots(c::Circuit, shots::Int) if shots > 0 && !isempty(c.result_types) !c.observables_simultaneously_measureable && throw(ErrorException("Observables cannot be measured simultaneously.")) for rt in c.result_types - (rt isa StateVector || rt isa Amplitude) && throw(ErrorException("StateVector or Amplitude cannot be specified when shots>0")) + (rt isa StateVector || rt isa Amplitude || rt isa AdjointGradient) && throw(ErrorException("StateVector or Amplitude cannot be specified when shots>0")) if rt isa Probability num_qubits = length(rt.targets) == 0 ? qubit_count(c) : length(rt.targets) num_qubits > 40 && throw(ErrorException("Probability target must be less than or equal to 40 qubits.")) @@ -728,5 +763,5 @@ end Base.:(==)(c1::Circuit, c2::Circuit) = (c1.instructions == c2.instructions && c1.result_types == c2.result_types) function Base.show(io::IO, program::OpenQasmProgram) - println(io, "OpenQASM program") + print(io, "OpenQASM program") end diff --git a/src/device.jl b/src/device.jl index e19fd896..3e025db2 100644 --- a/src/device.jl +++ b/src/device.jl @@ -33,15 +33,15 @@ properties(d::AwsDevice) = d._properties Base.convert(::Type{String}, d::AwsDevice) = d._arn Base.show(io::IO, d::AwsDevice) = print(io, "AwsDevice(arn="*d._arn*")") -function (d::AwsDevice)(task_spec::Union{Circuit, AnalogHamiltonianSimulation, AbstractProgram}; s3_destination_folder=default_task_bucket(), shots=nothing, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds::Int=DEFAULT_RESULTS_POLL_INTERVAL, kwargs...) +function (d::AwsDevice)(task_spec::Union{Circuit, AnalogHamiltonianSimulation, AbstractProgram}; s3_destination_folder=default_task_bucket(), shots=nothing, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds::Int=DEFAULT_RESULTS_POLL_INTERVAL, inputs=Dict{String, Float64}(), kwargs...) shots_ = isnothing(shots) ? d._default_shots : shots - return AwsQuantumTask(d._arn, task_spec, s3_destination_folder=s3_destination_folder, shots=shots_, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds, kwargs...) + return AwsQuantumTask(d._arn, task_spec, s3_destination_folder=s3_destination_folder, shots=shots_, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds, inputs=inputs, kwargs...) end # currently no batch support for AHS -function (d::AwsDevice)(task_specs::Vector{<:Union{Circuit, AbstractProgram}}; s3_destination_folder=default_task_bucket(), shots=nothing, max_parallel=nothing, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds::Int=DEFAULT_RESULTS_POLL_INTERVAL, kwargs...) +function (d::AwsDevice)(task_specs::Vector{<:Union{Circuit, AbstractProgram}}; s3_destination_folder=default_task_bucket(), shots=nothing, max_parallel=nothing, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds::Int=DEFAULT_RESULTS_POLL_INTERVAL, inputs=Dict{String, Float64}(), kwargs...) shots_ = isnothing(shots) ? d._default_shots : shots - return AwsQuantumTaskBatch(d._arn, task_specs; s3_destination_folder=s3_destination_folder, shots=shots_, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds, kwargs...) + return AwsQuantumTaskBatch(d._arn, task_specs; s3_destination_folder=s3_destination_folder, shots=shots_, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds, inputs=inputs, kwargs...) end function _construct_topology_graph(d::AwsDevice) diff --git a/src/observables.jl b/src/observables.jl index 8a9720f8..78695225 100644 --- a/src/observables.jl +++ b/src/observables.jl @@ -2,8 +2,8 @@ module Observables using JSON3, StructTypes, LinearAlgebra -import ..Braket: ir, qubit_count, pauli_eigenvalues, complex_matrix_from_ir, complex_matrix_to_ir, Operator, SerializationProperties, format_qubits, format, format_matrix, OpenQASMSerializationProperties, IRType, QubitSet, Qubit -export Observable, TensorProduct, HermitianObservable +import ..Braket: ir, qubit_count, pauli_eigenvalues, complex_matrix_from_ir, complex_matrix_to_ir, Operator, SerializationProperties, format_qubits, format, format_matrix, OpenQASMSerializationProperties, IRType, QubitSet, Qubit, IntOrQubit, IRObservable +export Observable, TensorProduct, HermitianObservable, Sum """ Observable <: Operator @@ -14,30 +14,43 @@ have `eigvals` defined. See also: [`H`](@ref), [`I`](@ref), [`X`](@ref), [`Y`](@ref), [`Z`](@ref), [`TensorProduct`](@ref), [`HermitianObservable`](@ref). """ abstract type Observable <: Operator end - LinearAlgebra.ishermitian(o::Observable) = true +coef(o::O) where {O<:Observable} = o.coefficient + for (typ, label) in ((:H, "h"), (:X, "x"), (:Y, "y"), (:Z, "z"), (:I, "i")) @eval begin @doc """ $($typ) <: Observable - $($typ)() -> $($typ) + $($typ)([coeff::Float64]) -> $($typ) - Struct representing a `$($typ)` observable in a measurement. + Struct representing a `$($typ)` observable in a measurement. The observable + may be scaled by `coeff`. """ - struct $typ <: Observable end + struct $typ <: Observable + coefficient::Float64 + $typ(coef::Float64=1.0) = new(coef) + end StructTypes.StructType(::Type{$typ}) = StructTypes.CustomStruct() StructTypes.lower(x::$typ) = [$label] qubit_count(o::$typ) = 1 qubit_count(::Type{$typ}) = 1 - Base.copy(o::$typ) = $typ() - ir(o::$typ, target::QubitSet, ::Val{:OpenQASM}; serialization_properties::SerializationProperties=OpenQASMSerializationProperties()) = isempty(target) ? $label * " all" : $label * "(" * format_qubits(target, serialization_properties) * ")" + unscaled(o::$typ) = $typ() + Base.copy(o::$typ) = $typ(o.coefficient) + function ir(o::$typ, target::QubitSet, ::Val{:OpenQASM}; serialization_properties::SerializationProperties=OpenQASMSerializationProperties()) + target_str = isempty(target) ? " all" : "(" * format_qubits(target, serialization_properties) * ")" + coef_str = isone(o.coefficient) ? "" : string(o.coefficient) * " * " + return coef_str * $label * target_str + end ir(o::$typ, target::Nothing, ::Val{:OpenQASM}; kwargs...) = ir(o, QubitSet(), Val(:OpenQASM); kwargs...) + Base.:(*)(o::$typ, n::Real) = $typ(Float64(n*o.coefficient)) + Base.:(==)(o1::$typ, o2::$typ) = (o1.coefficient ≈ o2.coefficient) + Base.show(io::IO, o::$typ) = print(io, (isone(o.coefficient) ? "" : string(o.coefficient) * " * ") * uppercase($label)) end end const StandardObservable = Union{H, X, Y, Z} -LinearAlgebra.eigvals(o::I) = [1.0, 1.0] -LinearAlgebra.eigvals(o::StandardObservable) = [1.0, -1.0] +LinearAlgebra.eigvals(o::I) = o.coefficient .* [1.0, 1.0] +LinearAlgebra.eigvals(o::StandardObservable) = o.coefficient .* [1.0, -1.0] """ TensorProduct <: Observable @@ -59,26 +72,33 @@ Braket.Observables.TensorProduct(Braket.Observables.Observable[Braket.Observable """ struct TensorProduct <: Observable factors::Vector{<:Observable} - function TensorProduct(v::Vector{<:Observable}) + coefficient::Float64 + function TensorProduct(v::Vector{<:Observable}, coefficient::Float64=1.0) flattened_v = Observable[] for o in v if o isa TensorProduct foreach(o_->push!(flattened_v, o_), o.factors) + elseif o isa Sum + throw(ArgumentError("Sum observable not allowed in TensorProduct.")) else push!(flattened_v, o) end end - return new(flattened_v) + coeff = mapreduce(coef, *, flattened_v, init=1.0) * coefficient + unscaled_factors = Observable[unscaled(f) for f in flattened_v] + return new(unscaled_factors, coeff) end end -function TensorProduct(o::Vector{String}) +Base.:(*)(o::TensorProduct, n::Real) = TensorProduct(deepcopy(o.factors), Float64(n*o.coefficient)) +function TensorProduct(o::Vector{String}, coefficient::Float64=1.0) sts = StructTypes.subtypes(Observable) ops = Observable[sts[Symbol(lowercase(oi))]() for oi in o] - return TensorProduct(ops) + return TensorProduct(ops, coefficient) end function ir(tp::TensorProduct, target::QubitSet, ::Val{:OpenQASM}; serialization_properties::SerializationProperties=OpenQASMSerializationProperties()) factors = [] use_qubits = collect(target) + coef_str = isone(tp.coefficient) ? "" : "$(tp.coefficient) * " for obs in tp.factors obs_target = QubitSet() num_qubits = qubit_count(obs) @@ -88,22 +108,29 @@ function ir(tp::TensorProduct, target::QubitSet, ::Val{:OpenQASM}; serialization end push!(factors, ir(obs, obs_target, Val(:OpenQASM), serialization_properties=serialization_properties)) end - return join(factors, " @ ") + return coef_str * join(factors, " @ ") end ir(tp::TensorProduct; kwargs...) = ir(tp, Val(IRType[]); kwargs...) - +unscaled(o::TensorProduct) = TensorProduct(o.factors, 1.0) qubit_count(o::TensorProduct) = sum(qubit_count.(o.factors)) StructTypes.StructType(::Type{TensorProduct}) = StructTypes.CustomStruct() -StructTypes.lower(x::TensorProduct) = vcat(StructTypes.lower.(x.factors)...) -Base.:(==)(t1::TensorProduct, t2::TensorProduct) = t1.factors == t2.factors -Base.copy(t::TensorProduct) = TensorProduct(deepcopy(t.factors)) +StructTypes.lower(x::TensorProduct) = Union{String, Vector{Vector{Vector{Float64}}}}[convert(Union{String, Vector{Vector{Vector{Float64}}}}, o) for o in mapreduce(StructTypes.lower, vcat, x.factors)] +Base.:(==)(t1::TensorProduct, t2::TensorProduct) = t1.factors == t2.factors && t1.coefficient ≈ t2.coefficient +Base.copy(t::TensorProduct) = TensorProduct(deepcopy(t.factors), t.coefficient) function LinearAlgebra.eigvals(o::TensorProduct) - all(o_ isa StandardObservable for o_ in o.factors) && return pauli_eigenvalues(length(o.factors)) + all(o_ isa StandardObservable for o_ in o.factors) && return o.coefficient .* pauli_eigenvalues(length(o.factors)) evs = mapfoldl(eigvals, kron, o.factors, init=[1.0]) - #=for f in o.factors - evs = kron(evs, eigvals(f)) - end=# - return evs + return o.coefficient .* evs +end +function Base.show(io::IO, o::TensorProduct) + coef_str = isone(o.coefficient) ? "" : string(o.coefficient) * " * " + print(io, coef_str) + for f in o.factors[1:end-1] + print(io, f) + print(io, " @ ") + end + print(io, o.factors[end]) + return end """ @@ -120,18 +147,20 @@ Braket.Observables.HermitianObservable(Complex{Int64}[0 + 0im 1 + 0im; 1 + 0im 0 """ struct HermitianObservable <: Observable matrix::Matrix{<:Complex} + coefficient::Float64 function HermitianObservable(mat::Matrix{<:Number}) ishermitian(mat) || throw(ArgumentError("input matrix to HermitianObservable must be Hermitian.")) - new(complex(mat)) + new(complex(mat), 1.0) end end HermitianObservable(v::Vector{Vector{Vector{T}}}) where {T<:Number} = HermitianObservable(complex_matrix_from_ir(v)) Base.copy(o::HermitianObservable) = HermitianObservable(copy(o.matrix)) StructTypes.StructType(::Type{HermitianObservable}) = StructTypes.CustomStruct() -StructTypes.lower(x::HermitianObservable) = [complex_matrix_to_ir(ComplexF64.(x.matrix))] +StructTypes.lower(x::HermitianObservable) = Union{String, Vector{Vector{Vector{Float64}}}}[complex_matrix_to_ir(ComplexF64.(x.matrix))] Base.:(==)(h1::HermitianObservable, h2::HermitianObservable) = (size(h1.matrix) == size(h2.matrix) && h1.matrix ≈ h2.matrix) qubit_count(o::HermitianObservable) = convert(Int, log2(size(o.matrix, 1))) LinearAlgebra.eigvals(o::HermitianObservable) = eigvals(Hermitian(o.matrix)) +unscaled(o::HermitianObservable) = o function ir(ho::HermitianObservable, target::QubitSet, ::Val{:OpenQASM}; serialization_properties=OpenQASMSerializationProperties()) function c_str(x::Number) iszero(real(x)) && iszero(imag(x)) && return "0im" @@ -144,37 +173,91 @@ function ir(ho::HermitianObservable, target::QubitSet, ::Val{:OpenQASM}; seriali t = isempty(target) ? "all" : format_qubits(target, serialization_properties) return "hermitian($m) $t" end +Base.:(*)(o::HermitianObservable, n::Real) = HermitianObservable(Float64(n) .* o.matrix, 1.0) ir(ho::HermitianObservable, target::Nothing, ::Val{:OpenQASM}; kwargs...) = ir(ho, QubitSet(), Val(:OpenQASM); kwargs...) -StructTypes.StructType(::Type{Observable}) = StructTypes.AbstractType() -StructTypes.subtypes(::Type{Observable}) = (i=I, h=H, x=X, y=Y, z=Z) -StructTypes.constructfrom(::Type{Observable}, obj::String) = StructTypes.subtypes(Observable)[Symbol(obj)]() -function StructTypes.constructfrom(::Type{Observable}, obj::Vector{String}) - length(obj) == 1 && return StructTypes.constructfrom(Observable, obj[1]) - return TensorProduct(obj) -end -StructTypes.constructfrom(::Type{Observable}, obj::Vector{Vector{Vector{T}}}) where {T<:Number} = HermitianObservable(obj) -function StructTypes.constructfrom(::Type{Observable}, obj::Vector{Vector{Vector{Vector{T}}}}) where {T<:Number} - length(obj) == 1 && return HermitianObservable(obj[1]) - return TensorProduct([HermitianObservable(o) for o in obj]) -end +""" + Sum <: Observable + Sum(summands::Vector{<:Observable}) -> Sum + +Struct representing the sum of observables. + +# Examples +```jldoctest +julia> o1 = 2.0 * Observables.I() @ Observables.Z(); -function StructTypes.constructfrom(::Type{Observable}, obj::Vector{Union{String, Vector{Vector{Vector{T}}}}}) where {T<:Number} - all(o isa String for o in obj) && return StructTypes.constructfrom(Observable, convert(Vector{String}, obj)) - all(o isa Vector{Vector{Vector{T}}} for o in obj) && return StructTypes.constructfrom(Observable, convert(Vector{Vector{Vector{Vector{T}}}}, obj)) - o_list = Observable[] - for o in obj - if o isa String - push!(o_list, StructTypes.constructfrom(Observable, o)) - elseif o isa Vector{Vector{Vector{T}}} - push!(o_list, HermitianObservable(o)) - else - throw(ErrorException("invalid observable $o")) +julia> o2 = 3.0 * Observables.X() @ Observables.X(); + +julia> o = o1 + o2 +Braket.Observables.Sum() +``` +""" +struct Sum <: Observable + summands::Vector{Observable} + coefficient::Float64 + function Sum(observables) + flattened_observables = Observable[] + for obs in observables + if obs isa Sum + append!(flattened_observables, obs.summands) + else + push!(flattened_observables, obs) + end end + new(flattened_observables, 1.0) + end +end +Base.length(s::Sum) = length(s.summands) +Base.:(*)(s::Sum, n::Real) = Sum(Float64(n) .* deepcopy(s.summands)) +function Base.:(==)(s1::Sum, s2::Sum) + length(s1) != length(s2) && return false + are_eq = true + for (summand1, summand2) in zip(s1.summands, s2.summands) + are_eq &= (summand1 == summand2) + are_eq || return false + end + return true +end +function Base.:(==)(s::Sum, o::Observable) + length(s) == 1 || return false + return first(s.summands) == o +end +Base.:(==)(o::Observable, s::Sum) = s == o +function Base.show(io::IO, s::Sum) + print(io, "Sum(") + for summand in s.summands[1:end-1] + print(io, summand) + print(io, ", ") + end + print(io, s.summands[end]) + print(io, ")") + return +end +ir(s::Sum, target::Vector{QubitSet}, ::Val{:JAQCD}; kwargs...) = throw(ErrorException("Sum observables are not supported in JAQCD.")) +ir(s::Sum, target::Vector{<:IntOrQubit}, ::Val{:JAQCD}; kwargs...) = throw(ErrorException("Sum observables are not supported in JAQCD.")) +function ir(s::Sum, target::Vector{QubitSet}, ::Val{:OpenQASM}; kwargs...) + length(s.summands) == length(target) || throw(DimensionMismatch("number of summands ($(length(s.summands))) must match length of targets vector ($(length(targets))).")) + for (ii, (term, term_target)) in enumerate(zip(s.summands, target)) + length(term_target) == qubit_count(term) || throw(DimensionMismatch("qubit count of term $ii ($(qubit_count(term))) must match number of targets ($(length(term_target))).")) end - return TensorProduct(o_list) + summands_irs = [ir(obs, targ, Val(:OpenQASM); kwargs...) for (obs, targ) in zip(s.summands, target)] + sum_str = replace(join(summands_irs, " + "), "+ -"=>"- ") + return sum_str +end +ir(s::Sum, target::Vector{Vector{T}}, ::Val{:OpenQASM}; kwargs...) where {T} = ir(s, [QubitSet(t) for t in target], Val(:OpenQASM); kwargs...) + +StructTypes.StructType(::Type{Observable}) = StructTypes.AbstractType() +StructTypes.subtypes(::Type{Observable}) = (i=I, h=H, x=X, y=Y, z=Z) +StructTypes.lowertype(::Type{O}) where {O<:Observable} = IRObservable +function StructTypes.constructfrom(::Type{Observable}, obj::IRObservable) + obj isa String && return StructTypes.subtypes(Observable)[Symbol(obj)]() + length(obj) == 1 && obj[1] isa Vector{Vector{Vector{Float64}}} && return HermitianObservable(obj[1]) + return TensorProduct([StructTypes.constructfrom(Observable, (o isa String ? o : convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, [o]))) for o in obj]) end Base.:(*)(o1::Observable, o2::Observable) = TensorProduct([o1, o2]) +Base.:(*)(n::Real, o::Observable) = o*n +Base.:(+)(o1::Observable, o2::Observable) = Sum([o1, o2]) +Base.:(-)(o1::Observable, o2::Observable) = Sum([o1, -1.0 * o2]) end diff --git a/src/operators.jl b/src/operators.jl index c944f705..4225927b 100644 --- a/src/operators.jl +++ b/src/operators.jl @@ -14,8 +14,8 @@ Abstract type representing *quantum* operations that can be applied to a [`Circu Subtypes include [`Gate`](@ref) and [`Noise`](@ref). """ abstract type QuantumOperator <: Operator end -ir(o::Operator, t::Vector{T}, ::Val{V}; kwargs...) where {T, V} = ir(o, QubitSet(t), Val(V); kwargs...) -ir(o::Operator, t::Vector{T}; kwargs...) where {T} = ir(o, QubitSet(t), Val(IRType[]); kwargs...) +ir(o::Operator, t::Vector{T}, ::Val{V}; kwargs...) where {T<:IntOrQubit, V} = ir(o, QubitSet(t), Val(V); kwargs...) +ir(o::Operator, t::Vector{T}; kwargs...) where {T<:IntOrQubit} = ir(o, QubitSet(t), Val(IRType[]); kwargs...) ir(o::Operator, t::QubitSet; kwargs...) = ir(o, t, Val(IRType[]); kwargs...) ir(o::Operator, t::IntOrQubit; kwargs...) = ir(o, QubitSet(t), Val(IRType[]); kwargs...) ir(o::Operator, t::IntOrQubit, args...; kwargs...) = ir(o, QubitSet(t), args...; kwargs...) @@ -23,4 +23,4 @@ ir(o::Operator, t::IntOrQubit, args...; kwargs...) = ir(o, QubitSet(t), args...; function pauli_eigenvalues(n::Int) n == 1 && return [1.0, -1.0] return reduce(vcat, [pauli_eigenvalues(n - 1), -1 .* pauli_eigenvalues(n - 1)]) -end \ No newline at end of file +end diff --git a/src/raw_schema.jl b/src/raw_schema.jl index f27e1bde..f2fe0ea8 100644 --- a/src/raw_schema.jl +++ b/src/raw_schema.jl @@ -114,12 +114,18 @@ module IR import ..Braket: AbstractProgram, braketSchemaHeader, BraketSchemaBase using DecFP, StructTypes -export Program, AHSProgram, AbstractIR, AbstractProgramResult, CompilerDirective +export Program, AHSProgram, AbstractIR, AbstractProgramResult, CompilerDirective, IRObservable abstract type AbstractIR end StructTypes.StructType(::Type{AbstractIR}) = StructTypes.AbstractType() StructTypes.subtypekey(::Type{AbstractIR}) = :type -StructTypes.subtypes(::Type{AbstractIR}) = (z=Z, sample=Sample, cphaseshift01=CPhaseShift01, phase_damping=PhaseDamping, rz=Rz, generalized_amplitude_damping=GeneralizedAmplitudeDamping, xx=XX, zz=ZZ, phase_flip=PhaseFlip, vi=Vi, depolarizing=Depolarizing, variance=Variance, two_qubit_depolarizing=TwoQubitDepolarizing, densitymatrix=DensityMatrix, cphaseshift00=CPhaseShift00, ecr=ECR, compilerdirective=CompilerDirective, ccnot=CCNot, unitary=Unitary, bit_flip=BitFlip, y=Y, swap=Swap, cz=CZ, cnot=CNot, cswap=CSwap, ry=Ry, i=I, si=Si, amplitude_damping=AmplitudeDamping, statevector=StateVector, iswap=ISwap, h=H, xy=XY, yy=YY, t=T, two_qubit_dephasing=TwoQubitDephasing, x=X, ti=Ti, cv=CV, pauli_channel=PauliChannel, pswap=PSwap, expectation=Expectation, probability=Probability, phaseshift=PhaseShift, v=V, cphaseshift=CPhaseShift, s=S, rx=Rx, kraus=Kraus, amplitude=Amplitude, cphaseshift10=CPhaseShift10, multi_qubit_pauli_channel=MultiQubitPauliChannel, cy=CY, setup=Setup, hamiltonian=Hamiltonian, shiftingfield=ShiftingField, atomarrangement=AtomArrangement, timeseries=TimeSeries, physicalfield=PhysicalField, ahsprogram=AHSProgram, drivingfield=DrivingField, abstractprogramresult=AbstractProgramResult) +StructTypes.subtypes(::Type{AbstractIR}) = (z=Z, sample=Sample, cphaseshift01=CPhaseShift01, phase_damping=PhaseDamping, rz=Rz, generalized_amplitude_damping=GeneralizedAmplitudeDamping, xx=XX, zz=ZZ, phase_flip=PhaseFlip, vi=Vi, depolarizing=Depolarizing, variance=Variance, two_qubit_depolarizing=TwoQubitDepolarizing, densitymatrix=DensityMatrix, cphaseshift00=CPhaseShift00, ecr=ECR, compilerdirective=CompilerDirective, ccnot=CCNot, unitary=Unitary, bit_flip=BitFlip, y=Y, swap=Swap, cz=CZ, cnot=CNot, adjoint_gradient=AdjointGradient, cswap=CSwap, ry=Ry, i=I, si=Si, amplitude_damping=AmplitudeDamping, statevector=StateVector, iswap=ISwap, h=H, xy=XY, yy=YY, t=T, two_qubit_dephasing=TwoQubitDephasing, x=X, ti=Ti, cv=CV, pauli_channel=PauliChannel, pswap=PSwap, expectation=Expectation, probability=Probability, phaseshift=PhaseShift, v=V, cphaseshift=CPhaseShift, s=S, rx=Rx, kraus=Kraus, amplitude=Amplitude, cphaseshift10=CPhaseShift10, multi_qubit_pauli_channel=MultiQubitPauliChannel, cy=CY, setup=Setup, hamiltonian=Hamiltonian, shiftingfield=ShiftingField, atomarrangement=AtomArrangement, timeseries=TimeSeries, physicalfield=PhysicalField, ahsprogram=AHSProgram, drivingfield=DrivingField, abstractprogramresult=AbstractProgramResult) + +const IRObservable = Union{Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, String} +Base.convert(::Type{IRObservable}, v::Vector{String}) = convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, v) +Base.convert(::Type{IRObservable}, s::String) = s +Base.convert(::Type{IRObservable}, v::Vector{Vector{Vector{Vector{Float64}}}}) = convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, v) + abstract type CompilerDirective <: AbstractIR end @@ -130,7 +136,7 @@ abstract type AbstractProgramResult <: AbstractIR end StructTypes.StructType(::Type{AbstractProgramResult}) = StructTypes.AbstractType() -StructTypes.subtypes(::Type{AbstractProgramResult}) = (amplitude=Amplitude, expectation=Expectation, probability=Probability, sample=Sample, statevector=StateVector, densitymatrix=DensityMatrix, variance=Variance) +StructTypes.subtypes(::Type{AbstractProgramResult}) = (amplitude=Amplitude, expectation=Expectation, probability=Probability, sample=Sample, statevector=StateVector, densitymatrix=DensityMatrix, variance=Variance, adjoint_gradient=AdjointGradient) struct AtomArrangement <: AbstractIR sites::Vector{Vector{Dec128}} filling::Vector{Int} @@ -182,7 +188,7 @@ StructTypes.StructType(::Type{Z}) = StructTypes.UnorderedStruct() StructTypes.defaults(::Type{Z}) = Dict{Symbol, Any}(:type => "z") struct Sample <: AbstractProgramResult - observable::Vector{Union{String, Vector{Vector{Vector{Float64}}}}} + observable::IRObservable targets::Union{Nothing, Vector{Int}} type::String end @@ -263,7 +269,7 @@ StructTypes.StructType(::Type{Depolarizing}) = StructTypes.UnorderedStruct() StructTypes.defaults(::Type{Depolarizing}) = Dict{Symbol, Any}(:type => "depolarizing") struct Variance <: AbstractProgramResult - observable::Vector{Union{String, Vector{Vector{Vector{Float64}}}}} + observable::IRObservable targets::Union{Nothing, Vector{Int}} type::String end @@ -371,6 +377,15 @@ end StructTypes.StructType(::Type{CNot}) = StructTypes.UnorderedStruct() StructTypes.defaults(::Type{CNot}) = Dict{Symbol, Any}(:type => "cnot") +struct AdjointGradient <: AbstractProgramResult + parameters::Union{Nothing, Vector{String}} + observable::IRObservable + targets::Union{Nothing, Vector{Vector{Int}}} + type::String +end +StructTypes.StructType(::Type{AdjointGradient}) = StructTypes.UnorderedStruct() + +StructTypes.defaults(::Type{AdjointGradient}) = Dict{Symbol, Any}(:type => "adjoint_gradient") struct CSwap <: AbstractIR control::Int targets::Vector{Int} @@ -508,7 +523,7 @@ StructTypes.StructType(::Type{PSwap}) = StructTypes.UnorderedStruct() StructTypes.defaults(::Type{PSwap}) = Dict{Symbol, Any}(:type => "pswap") struct Expectation <: AbstractProgramResult - observable::Vector{Union{String, Vector{Vector{Vector{Float64}}}}} + observable::IRObservable targets::Union{Nothing, Vector{Int}} type::String end @@ -637,6 +652,12 @@ for g in [:Expectation, :Sample, :Variance, :DensityMatrix, :Probability] Target(::Type{$g}) = OptionalMultiTarget() end end +struct OptionalNestedMultiTarget <: Target end +for g in [:AdjointGradient] + @eval begin + Target(::Type{$g}) = OptionalNestedMultiTarget() + end +end struct DoubleTarget <: Target end for g in [:Swap, :CSwap, :ISwap, :PSwap, :XY, :ECR, :XX, :YY, :ZZ, :TwoQubitDepolarizing, :TwoQubitDephasing] @eval begin @@ -691,10 +712,6 @@ end # module using .IR -@enum _CopyMode always if_needed never -_CopyModeDict = Dict(string(inst)=>inst for inst in instances(_CopyMode)) -Base.:(==)(x::_CopyMode, y::String) = string(x) == y -Base.:(==)(x::String, y::_CopyMode) = string(y) == x @enum SafeUUID safe unsafe unknown SafeUUIDDict = Dict(string(inst)=>inst for inst in instances(SafeUUID)) Base.:(==)(x::SafeUUID, y::String) = string(x) == y @@ -1027,6 +1044,7 @@ struct DeviceServiceProperties <: BraketSchemaBase deviceDocumentation::Union{Nothing, DeviceDocumentation} deviceLocation::Union{Nothing, String} updatedAt::Union{Nothing, String} + getTaskPollIntervalMillis::Union{Nothing, Int} end StructTypes.StructType(::Type{DeviceServiceProperties}) = StructTypes.UnorderedStruct() diff --git a/src/results.jl b/src/results.jl index e0eeab7d..622cdcd5 100644 --- a/src/results.jl +++ b/src/results.jl @@ -10,7 +10,7 @@ See also: [`Expectation`](@ref), [`Variance`](@ref), """ abstract type Result end -union_obs_typ = @NamedTuple{observable::Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, targets::Vector{Int}, type::String} +union_obs_typ = @NamedTuple{observable::IRObservable, targets::Vector{Int}, type::String} for (typ, ir_typ, label) in ((:Expectation, :(Braket.IR.Expectation), "expectation"), (:Variance, :(Braket.IR.Variance), "variance"), (:Sample, :(Braket.IR.Sample), "sample")) @eval begin @@ -46,7 +46,7 @@ for (typ, ir_typ, label) in ((:Expectation, :(Braket.IR.Expectation), "expectati $typ(o) = $typ(o, QubitSet()) Base.:(==)(e1::$typ, e2::$typ) = (e1.observable == e2.observable && e1.targets == e2.targets) StructTypes.StructType(::Type{$typ}) = StructTypes.CustomStruct() - StructTypes.lower(x::$typ) = $ir_typ(ir(x.observable), (isempty(x.targets) ? nothing : Int.(x.targets)), $label) + StructTypes.lower(x::$typ) = $ir_typ(StructTypes.lower(x.observable), (isempty(x.targets) ? nothing : Int.(x.targets)), $label) StructTypes.lowertype(::Type{$typ}) = union_obs_typ function ir(r::$typ, ::Val{:OpenQASM}; serialization_properties::SerializationProperties=OpenQASMSerializationProperties()) obs_ir = ir(r.observable, r.targets, Val(:OpenQASM); serialization_properties=serialization_properties) @@ -100,6 +100,104 @@ function ir(dm::DensityMatrix, ::Val{:OpenQASM}; serialization_properties::Seria t = format_qubits(dm.targets, serialization_properties) return "#pragma braket result density_matrix $t" end +""" + AdjointGradient <: Result + +Struct which represents a gradient computation using the adjoint differentiation method on a [`Circuit`](@ref). +""" +struct AdjointGradient <: Result + observable::Observable + targets::Vector{QubitSet} + parameters::Vector{String} + function AdjointGradient(observable::Observable, targets::Vector{QubitSet}, parameters::Vector{String}=["all"]) + if observable isa Sum + length(targets) == length(observable) || throw(DimensionMismatch("length of targets ($(length(targets))) must be the same as number of summands ($(length(observable))).")) + all(length(term_target) == qubit_count(summand) for (term_target, summand) in zip(targets, observable.summands)) || throw(DimensionMismatch("each target must be the same size as the qubit count of its corresponding term.")) + else + (length(targets) == 1 && length(targets[1]) == qubit_count(observable)) || throw(DimensionMismatch("targets $targets must have only one element if adjoint gradient observable is not a Sum.")) + end + new(observable, targets, parameters) + end +end + +""" + AdjointGradient(o::Observable, targets, parameters::Vector) -> AdjointGradient + AdjointGradient(o::Observable, targets) -> AdjointGradient + +Constructs an `AdjointGradient` with respect to the expectation value of +an observable `o` on qubits `targets`. The gradient will be calculated by +computing partial derivatives with respect to `parameters`. If `parameters` +is not present, is empty, or is `["all"]`, all parameters in the circuit +will be used. + +`targets` may be one of: + - A [`QubitSet`](@ref) + - A `Vector` of `Int`s and/or [`Qubit`](@ref)s + - An `Int` or `Qubit` + +`AdjointGradient` supports using [`Sum`](@ref) observables. If `o` is a `Sum`, +`targets` should be a nested vector of target qubits, such that the `n`-th term of +`targets` has the same length as the `n`-th term of `o`. + +# Examples +```jldoctest +julia> α = FreeParameter(:alpha); + +julia> c = Circuit([(H, 0), (H, 1), (Rx, 0, α), (Rx, 1, α)]); + +julia> op = 2.0 * Braket.Observables.X() * Braket.Observables.X(); + +julia> c = AdjointGradient(c, op, [QubitSet(0, 1)], [α]); +``` + +Using a `Sum`: + +```jldoctest +julia> α = FreeParameter(:alpha); + +julia> c = Circuit([(H, 0), (H, 1), (Rx, 0, α), (Rx, 1, α)]); + +julia> op1 = 2.0 * Braket.Observables.X() * Braket.Observables.X(); + +julia> op2 = -3.0 * Braket.Observables.Y() * Braket.Observables.Y(); + +julia> c = AdjointGradient(c, op1 + op2, [QubitSet(0, 1), QubitSet(0, 2)], [α]); +``` +""" +function AdjointGradient(observable::Observable, targets::Vector{QubitSet}, parameters::Vector) + isempty(parameters) && return AdjointGradient(observable, targets, ["all"]) + return AdjointGradient(observable, targets, string.(parameters)) +end +AdjointGradient(observable::Observable, targets::QubitSet, parameters) = AdjointGradient(observable, [targets], parameters) +AdjointGradient(observable::Observable, targets::Vector{Vector{T}}, args...) where {T} = AdjointGradient(observable, [QubitSet(t) for t in targets], args...) +AdjointGradient(observable::Observable, targets::Vector{<:IntOrQubit}, args...) = AdjointGradient(observable, [QubitSet(targets)], args...) +AdjointGradient(observable::Observable, targets::IntOrQubit, args...) = AdjointGradient(observable, [QubitSet(targets)], args...) + +function ir(ag::AdjointGradient, ::Val{:OpenQASM}; serialization_properties::SerializationProperties=OpenQASMSerializationProperties()) + isempty(ag.parameters) && (ag.parameters = ["all"]) + if ag.observable isa Sum + obs_str = ir(ag.observable, ag.targets, Val(:OpenQASM), serialization_properties=serialization_properties) + else + obs_str = ir(ag.observable, first(ag.targets), Val(:OpenQASM), serialization_properties=serialization_properties) + end + term_str = replace(obs_str, "+ -"=>"- ") + return "#pragma braket result adjoint_gradient expectation($term_str) $(join(ag.parameters, ", "))" +end +StructTypes.StructType(::Type{AdjointGradient}) = StructTypes.CustomStruct() +function StructTypes.lower(x::AdjointGradient) + lowered_obs = ir(x.observable, Val(:JAQCD)) + lowered_targets = (isempty(x.targets) ? nothing : convert(Vector{Vector{Int}}, x.targets)) + Braket.IR.AdjointGradient(x.parameters, lowered_obs, lowered_targets, "adjoint_gradient") +end +StructTypes.lowertype(::Type{AdjointGradient}) = @NamedTuple{parameters::Union{Nothing, Vector{String}}, observable::IRObservable, targets::Union{Nothing, Vector{Vector{Int}}}, type::String} +function AdjointGradient(x::@NamedTuple{parameters::Union{Nothing, Vector{String}}, observable::IRObservable, targets::Union{Nothing, Vector{Vector{Int}}}, type::String}) + params = isnothing(x.parameters) || isempty(x.parameters) ? ["all"] : x.parameters + obs = StructTypes.constructfrom(Observables.Observable, x.observable) + targets = isnothing(x.targets) || isempty(x.targets) ? QubitSet[] : [QubitSet(t) for t in x.targets] + return AdjointGradient(obs, targets, parameters) +end + + """ Amplitude <: Result @@ -135,6 +233,7 @@ ir(s::StateVector, ::Val{:OpenQASM}; kwargs...) = "#pragma braket result state_v Base.:(==)(sv1::StateVector, sv2::StateVector) = true const ObservableResult = Union{Expectation, Variance, Sample} +const ObservableParameterResult = Union{AdjointGradient,} remap(rt::R, mapping::Dict{<:Integer, <:Integer}) where {R<:ObservableResult} = R(copy(rt.observable), [mapping[q] for q in rt.targets]) remap(rt::R, mapping::Dict{<:Integer, <:Integer}) where {R<:Result} = R([mapping[q] for q in rt.targets]) @@ -171,8 +270,15 @@ function StructTypes.constructfrom(::Type{Result}, r::Braket.IR.AbstractProgramR StructTypes.constructfrom(typ, r) end +function Base.getindex(r::GateModelQuantumTaskResult, rt::AdjointGradient) + # only one AdjointGradient result per circuit + ix = findfirst(rt_->rt_.type isa IR.AdjointGradient, r.result_types) + isnothing(ix) && throw(KeyError(rt)) + return r.values[ix] +end + function Base.getindex(r::GateModelQuantumTaskResult, rt::Result) ix = findfirst(rt_->rt_.type==ir(rt, Val(:JAQCD)), r.result_types) isnothing(ix) && throw(KeyError(rt)) return r.values[ix] -end \ No newline at end of file +end diff --git a/src/schemas.jl b/src/schemas.jl index 2a1b90fd..fe60a463 100644 --- a/src/schemas.jl +++ b/src/schemas.jl @@ -72,3 +72,4 @@ function StructTypes.constructfrom(::Type{Program}, obj) new_obj[:instructions] = StructTypes.constructfrom(Vector{Instruction}, new_obj[:instructions]) return StructTypes.constructfrom(StructTypes.StructType(Program), Program, new_obj) end + diff --git a/src/task.jl b/src/task.jl index a21c63f4..e9e182bd 100644 --- a/src/task.jl +++ b/src/task.jl @@ -83,6 +83,7 @@ Valid `kwargs` are: - `poll_timeout_seconds::Int` - maximum number of seconds to wait while polling for results. Default: $DEFAULT_RESULTS_POLL_TIMEOUT - `poll_interval_seconds::Int` - default number of seconds to wait between attempts while polling for results. Default: $DEFAULT_RESULTS_POLL_INTERVAL - `tags::Dict{String, String}` - tags for the `AwsQuantumTask` + - `inputs::Dict{String, Float64}` - input values for any free parameters in the `task_spec` """ function AwsQuantumTask(device_arn::String, task_spec::Union{AbstractProgram, Circuit, AnalogHamiltonianSimulation}; @@ -92,8 +93,9 @@ function AwsQuantumTask(device_arn::String, disable_qubit_rewiring::Bool=false, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds::Int=DEFAULT_RESULTS_POLL_INTERVAL, - tags::Dict{String, String}=Dict{String, String}()) - args = prepare_task_input(task_spec, device_arn, s3_destination_folder, shots, device_params, disable_qubit_rewiring, tags=tags, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds) + tags::Dict{String, String}=Dict{String, String}(), + inputs::Dict{String,Float64}=Dict{String, Float64}()) + args = prepare_task_input(task_spec, device_arn, s3_destination_folder, shots, device_params, disable_qubit_rewiring, tags=tags, poll_timeout_seconds=poll_timeout_seconds, poll_interval_seconds=poll_interval_seconds, inputs=inputs) return AwsQuantumTask(args) end @@ -139,7 +141,26 @@ function prepare_task_input(ahs::AnalogHamiltonianSimulation, device_arn::String return prepare_task_input(ir(ahs), device_arn, s3_folder, shots, device_params, disable_qubit_rewiring; kwargs...) end -function prepare_task_input(program::Union{AHSProgram, BlackbirdProgram, OpenQasmProgram}, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) +function prepare_task_input(program::OpenQasmProgram, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) + common = _create_common_params(device_arn, s3_folder, shots; kwargs...) + client_token = string(uuid1()) + tags = get(kwargs, :tags, Dict{String,String}()) + device_parameters = Dict{String, Any}() # not currently used + dev_params = JSON3.write(device_parameters) + extra_opts = Dict("deviceParameters"=>dev_params, "tags"=>tags) + + inputs = get(kwargs, :inputs, Dict{String, Float64}()) + if !isempty(inputs) + inputs_merged = merge(program.inputs, inputs) + program_ = OpenQasmProgram(program.braketSchemaHeader, program.source, inputs_merged) + action = JSON3.write(program_) + else + action = JSON3.write(program) + end + return merge((action=action, client_token=client_token, extra_opts=extra_opts), common) +end + +function prepare_task_input(program::Union{AHSProgram, BlackbirdProgram}, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) device_parameters = Dict{String, Any}() # not currently used common = _create_common_params(device_arn, s3_folder, shots; kwargs...) client_token = string(uuid1()) @@ -150,7 +171,7 @@ function prepare_task_input(program::Union{AHSProgram, BlackbirdProgram, OpenQas return merge((action=action, client_token=client_token, extra_opts=extra_opts), common) end -function prepare_task_input(circuit::Union{Program, Circuit}, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) +function prepare_task_input(circuit::Circuit, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) validate_circuit_and_shots(circuit, shots) common = _create_common_params(device_arn, s3_folder, shots; kwargs...) paradigm_parameters = GateModelParameters(header_dict[GateModelParameters], qubit_count(circuit), disable_qubit_rewiring) @@ -162,9 +183,41 @@ function prepare_task_input(circuit::Union{Program, Circuit}, device_arn::String elseif occursin("oqc", device_arn) T = OqcDeviceParameters end + qubit_reference_type = VIRTUAL + if disable_qubit_rewiring || Instruction(StartVerbatimBox()) in circuit.instructions + #|| any(instruction.operator isa PulseGate for instruction in circuit.instructions) + qubit_reference_type = QubitReferenceType.PHYSICAL + end + serialization_properties = OpenQASMSerializationProperties(qubit_reference_type=qubit_reference_type) + oq3_program = ir(circuit, Val(:OpenQASM), serialization_properties=serialization_properties) + inputs = get(kwargs, :inputs, Dict{String, Float64}()) + program_inputs = isnothing(oq3_program.inputs) ? Dict{String, Float64}() : oq3_program.inputs + inputs_merged = !isempty(inputs) ? merge(program_inputs, inputs) : oq3_program.inputs + oq3_program = OpenQasmProgram(oq3_program.braketSchemaHeader, oq3_program.source, inputs_merged) + device_parameters = T(header_dict[T], paradigm_parameters) client_token = string(uuid1()) - action = ir(circuit) + action = JSON3.write(oq3_program) + dev_params = JSON3.write(device_parameters) + tags = get(kwargs, :tags, Dict{String,String}()) + extra_opts = Dict("deviceParameters"=>dev_params, "tags"=>tags) + return merge((action=action, client_token=client_token, extra_opts=extra_opts), common) +end + +function prepare_task_input(circuit::Program, device_arn::String, s3_folder::Tuple{String, String}, shots::Int, device_params::Dict{String, Any}, disable_qubit_rewiring::Bool=false; kwargs...) + common = _create_common_params(device_arn, s3_folder, shots; kwargs...) + paradigm_parameters = GateModelParameters(header_dict[GateModelParameters], qubit_count(circuit), disable_qubit_rewiring) + T = GateModelSimulatorDeviceParameters # default to use simulator + if occursin("ionq", device_arn) + T = IonqDeviceParameters + elseif occursin("rigetti", device_arn) + T = RigettiDeviceParameters + elseif occursin("oqc", device_arn) + T = OqcDeviceParameters + end + device_parameters = T(header_dict[T], paradigm_parameters) + client_token = string(uuid1()) + action = JSON3.write(circuit) dev_params = JSON3.write(device_parameters) tags = get(kwargs, :tags, Dict{String,String}()) extra_opts = Dict("deviceParameters"=>dev_params, "tags"=>tags) @@ -384,7 +437,7 @@ function calculate_result_types(ir_obj, measurements, measured_qubits::Vector{In !haskey(ir_obj, "results") && return result_types for result_type in ir_obj["results"] ir_obs = get(result_type, "observable", nothing) - rt = StructTypes.constructfrom(AbstractProgramResult, result_type) + rt = JSON3.read(JSON3.write(result_type), AbstractProgramResult) obs = isnothing(ir_obs) ? nothing : StructTypes.constructfrom(Observables.Observable, rt.observable) targs = rt.targets typ = result_type["type"] @@ -441,13 +494,12 @@ function from_dict(::Type{GateModelQuantumTaskResult}, r::GateModelTaskResult) task_mtd = r.taskMetadata addi_mtd = r.additionalMetadata rts = r.resultTypes - #=vals = map(rts) do rt + vals = map(rts) do rt !(rt.type isa IR.Amplitude) && return rt.value - # need to fix return dict - fixed_dict = Dict(string(k)=>complex(v[1], v[2]) for (k,v) in rt.value) - return fixed_dict - end=# - vals = [rt.value for rt in rts] + val_keys = string.(keys(rt.value)) + val_vals = [complex(v...) for v in values(rt.value)] + return Dict{String, ComplexF64}(zip(val_keys, val_vals)) + end return GateModelQuantumTaskResult(task_mtd, addi_mtd, rts, vals, nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing) end diff --git a/src/task_batch.jl b/src/task_batch.jl index 1544d966..8456666b 100644 --- a/src/task_batch.jl +++ b/src/task_batch.jl @@ -31,9 +31,9 @@ mutable struct AwsQuantumTaskBatch - `poll_timeout_seconds::Int` - maximum number of seconds to wait while polling for results. Default: $DEFAULT_RESULTS_POLL_TIMEOUT - `poll_interval_seconds::Int` - default number of seconds to wait between attempts while polling for results. Default: $DEFAULT_RESULTS_POLL_INTERVAL """ - function AwsQuantumTaskBatch(device_arn::String, task_specs::Vector{<:Union{AbstractProgram, Circuit}}; s3_destination_folder::Tuple{String, String}=default_task_bucket(), shots::Int=DEFAULT_SHOTS, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds=DEFAULT_RESULTS_POLL_INTERVAL) - tasks = launch_batch(device_arn, task_specs; s3_destination_folder=s3_destination_folder, shots=shots, poll_interval_seconds=poll_interval_seconds) - new(tasks, nothing, Set(), device_arn, task_specs, s3_destination_folder, shots, poll_timeout_seconds, poll_interval_seconds) + function AwsQuantumTaskBatch(device_arn::String, task_specs::Vector{<:Union{AbstractProgram, Circuit}}; s3_destination_folder::Tuple{String, String}=default_task_bucket(), shots::Int=DEFAULT_SHOTS, poll_timeout_seconds::Int=DEFAULT_RESULTS_POLL_TIMEOUT, poll_interval_seconds=DEFAULT_RESULTS_POLL_INTERVAL, inputs=Dict{String, Float64}()) + tasks = launch_batch(device_arn, task_specs; s3_destination_folder=s3_destination_folder, shots=shots, poll_interval_seconds=poll_interval_seconds, inputs=inputs) + new(tasks, nothing, Set(), device_arn, task_specs, s3_destination_folder, shots, poll_timeout_seconds, poll_interval_seconds) end end @@ -42,7 +42,20 @@ function launch_batch(device_arn::String, task_specs::Vector{<:Union{AbstractPro s3_folder = haskey(kwargs, :s3_destination_folder) ? kwargs[:s3_destination_folder] : default_task_bucket() shots = get(kwargs, :shots, DEFAULT_SHOTS) device_params = get(kwargs, :device_params, Dict{String, Any}()) - input = [prepare_task_input(ts, device_arn, s3_folder, shots, device_params, disable_qubit_rewiring; kwargs...) for ts in task_specs] + input_dicts = get(kwargs, :inputs, Dict{String, Float64}()) + if !isempty(input_dicts) + if input_dicts isa Dict + input = [prepare_task_input(ts, device_arn, s3_folder, shots, device_params, disable_qubit_rewiring; inputs=input_dicts, kwargs...) for ts in task_specs] + elseif input_dicts isa Vector + length(input_dicts) == length(task_specs) || throw(DimensionMismatch("number of input dictionaries $(length(input_dicts)) must equal number of task specifications $(length(task_specs)).")) + + input = [prepare_task_input(ts, device_arn, s3_folder, shots, device_params, disable_qubit_rewiring; inputs=ts_input, kwargs...) for (ts, ts_input) in zip(task_specs, input_dicts)] + else + throw(ArgumentError("invalid inputs $input_dicts.")) + end + else + input = [prepare_task_input(ts, device_arn, s3_folder, shots, device_params, disable_qubit_rewiring; kwargs...) for ts in task_specs] + end to_launch = Threads.Atomic{Int}(length(task_specs)) tasks = Vector{AwsQuantumTask}(undef, length(task_specs)) Threads.@threads for ii in 1:length(task_specs) diff --git a/test/Project.toml b/test/Project.toml index a7efe02d..a0316e82 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -11,6 +11,7 @@ Mocking = "78c3b35d-d492-501b-9361-3d52fe80e533" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StructTypes = "856f2bd8-1eba-4b0a-8007-ebc267875bd4" Tar = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/circuits.jl b/test/circuits.jl index 6752abcd..42d4e6a0 100644 --- a/test/circuits.jl +++ b/test/circuits.jl @@ -270,6 +270,31 @@ using Braket: Instruction, Result, VIRTUAL, PHYSICAL, OpenQASMSerializationPrope c = Probability(c) @test c.observables_simultaneously_measureable end + @testset "Adjoint gradient" begin + α = FreeParameter(:alpha) + c = Circuit([(H, 0), (H, 1), (Rx, 0, α), (Rx, 1, α)]) + op = 2.0 * Braket.Observables.X() * Braket.Observables.X() + op2 = 2.0 * Braket.Observables.Y() * Braket.Observables.Y() + c = AdjointGradient(c, op, [QubitSet(0, 1)], [α]) + @test length(c.result_types) == 1 + @test c.result_types[1] isa AdjointGradient + @test c.result_types[1].observable == op + @test c.result_types[1].targets == [QubitSet(0, 1)] + @test c.result_types[1].parameters == ["alpha"] + @test_throws ArgumentError AdjointGradient(c, op2, [QubitSet(0, 1)], [α]) + c = Circuit([(H, 0), (H, 1), (Rx, 0, α), (Rx, 1, α)]) + op = 2.0 * Braket.Observables.X() * Braket.Observables.X() + @test_throws DimensionMismatch AdjointGradient(c, op, [QubitSet(0)], [α]) + op3 = op + op2 + @test_throws DimensionMismatch AdjointGradient(c, op3, [QubitSet(0, 1)], [α]) + # make sure qubit count is correct + α = FreeParameter(:alpha) + c = Circuit([(H, 0), (H, 1), (Rx, 0, α), (Rx, 1, α)]) + op = 2.0 * Braket.Observables.X() * Braket.Observables.X() + @test qubit_count(c) == 2 + c = AdjointGradient(c, op, [QubitSet(1, 2)], [α]) + @test qubit_count(c) == 3 + end @testset "Basis rotation instructions" begin @testset "Basic" begin circ = Circuit([(H, 0), (CNot, 0, 1)]) @@ -664,6 +689,8 @@ using Braket: Instruction, Result, VIRTUAL, PHYSICAL, OpenQASMSerializationPrope end @testset "Result types & OpenQASM" begin + sum_obs = 2.0 * Observables.H() - 5.0 * Observables.Z() * Observables.X() + herm = Observables.HermitianObservable(diagm(ones(2))) @testset for ir_bolus in [ (Expectation(Braket.Observables.I(), 0), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result expectation i(q[0])",), (Expectation(Braket.Observables.I()), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result expectation i all",), @@ -672,6 +699,14 @@ using Braket: Instruction, Result, VIRTUAL, PHYSICAL, OpenQASMSerializationPrope (DensityMatrix([0, 2]), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result density_matrix q[0], q[2]",), (DensityMatrix(0), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), "#pragma braket result density_matrix \$0",), (Amplitude(["01", "10"]), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), "#pragma braket result amplitude \"01\", \"10\"",), + (AdjointGradient(Observables.H(), 0, ["alpha"]), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(h(q[0])) alpha",), + (AdjointGradient(Observables.H(), 0), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(h(q[0])) all",), + (AdjointGradient(Observables.X() * Observables.Y(), [0, 1], []), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(x(q[0]) @ y(q[1])) all",), + (AdjointGradient(Observables.H(), 0, []), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(h(q[0])) all",), + (AdjointGradient(sum_obs, [[0], [1, 2]], ["alpha"]), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(2.0 * h(q[0]) - 5.0 * z(q[1]) @ x(q[2])) alpha",), + (AdjointGradient(sum_obs, [[0], [1, 2]]), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(2.0 * h(q[0]) - 5.0 * z(q[1]) @ x(q[2])) all",), + (AdjointGradient(sum_obs, [[0], [1, 2]], []), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(2.0 * h(q[0]) - 5.0 * z(q[1]) @ x(q[2])) all",), + (AdjointGradient(herm, 0, []), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result adjoint_gradient expectation(hermitian([[1+0im, 0im], [0im, 1+0im]]) q[0]) all",), (Probability(), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result probability all",), (Probability([0, 2]), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), "#pragma braket result probability q[0], q[2]",), (Probability(0), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), "#pragma braket result probability \$0",), diff --git a/test/device.jl b/test/device.jl index 56392db0..7ff7b6cf 100644 --- a/test/device.jl +++ b/test/device.jl @@ -162,7 +162,7 @@ MOCK_DWAVE_QPU() = """{ @test sprint(show, dev) == "AwsDevice(arn=fake:arn)" end execution_window = Braket.DeviceExecutionWindow("everyday",Dates.Time("00:00:00"),Dates.Time("23:59:59")) - dsp = Braket.DeviceServiceProperties(Braket.header_dict[Braket.DeviceServiceProperties], [execution_window], (0, 1000), nothing, nothing, nothing, nothing) + dsp = Braket.DeviceServiceProperties(Braket.header_dict[Braket.DeviceServiceProperties], [execution_window], (0, 1000), nothing, nothing, nothing, nothing, nothing) paradigm = Braket.GateModelSimulatorParadigmProperties(Braket.header_dict[Braket.GateModelSimulatorParadigmProperties], 100) dev_capa = Braket.GateModelSimulatorDeviceCapabilities(dsp, Dict(), Dict(), Braket.header_dict[Braket.GateModelSimulatorDeviceCapabilities], paradigm) resp_dict = Dict{String, Any}("devices"=>[Dict{String, Any}( diff --git a/test/gate_model_task_result.jl b/test/gate_model_task_result.jl index a55ce337..f1ca6531 100644 --- a/test/gate_model_task_result.jl +++ b/test/gate_model_task_result.jl @@ -1,17 +1,19 @@ -using Braket, Braket.IR, Test, JSON3 +using Braket, Braket.Observables, Braket.IR, Test, JSON3 using Braket: ResultTypeValue + +exact_results = [(Probability(0), [0.5, 0.5]), + (StateVector(), [complex(0.70710678, 0), 0, 0, complex(0.70710678, 0)]), + (Expectation(Braket.Observables.Y(), 0), 0.0), + (Variance(Braket.Observables.Y(), 0), 0.1), + (Amplitude("00"), Dict("00"=>complex(0.70710678, 0))), + (AdjointGradient(2.0 * Observables.X() * Observables.X(), [0, 1], ["p_1", "p_2"]), Dict(:expectation=>0.1, :gradient=>Dict("p_1"=>0.2, "p_2"=>0.3)))] + zero_shots_result(task_mtd, add_mtd) = Braket.GateModelTaskResult( Braket.header_dict[Braket.GateModelTaskResult], nothing, nothing, - [ - ResultTypeValue(IR.Probability([0], "probability"), [0.5, 0.5]), - ResultTypeValue(IR.StateVector("statevector"), [complex(0.70710678, 0), 0, 0, complex(0.70710678, 0)]), - ResultTypeValue(IR.Expectation(["y"], [0], "expectation"), 0.0), - ResultTypeValue(IR.Variance(["y"], [0], "variance"), 0.1), - ResultTypeValue(IR.Amplitude(["00"], "amplitude"), Dict("00"=>complex(0.70710678, 0))) - ], + map(r->ResultTypeValue(ir(r[1], Val(:JAQCD)), r[2]), exact_results), [0,1], task_mtd, add_mtd, @@ -36,6 +38,11 @@ non_zero_shots_result(task_mtd, add_mtd) = Braket.GateModelTaskResult( g = Braket.format_result(r) @test g isa Braket.GateModelQuantumTaskResult @test sprint(show, g) == "GateModelQuantumTaskResult\n" + if shots == 0 + for er in exact_results + @test g[er[1]] == er[2] + end + end end task_metadata = Braket.TaskMetadata(Braket.header_dict[Braket.TaskMetadata], "task_arn", 0, "arn1", nothing, nothing, nothing, nothing, nothing) additional_metadata = Braket.AdditionalMetadata(action, nothing, nothing, nothing, nothing, nothing, nothing) @@ -58,8 +65,13 @@ end ] mat = [1 0; 0 -1] ho = Braket.Observables.HermitianObservable(mat) - samp = ir(Braket.Sample(ho, Int[2])) - action = Braket.Program(Braket.header_dict[Braket.Program], [Braket.IR.CNot(0, 1, "cnot"), Braket.IR.CNot(2, 3, "cnot")], [Braket.IR.Probability([1], "probability"), Braket.IR.Expectation(["z"], nothing, "expectation"), Braket.IR.Variance(["x", "x"], [0, 2], "variance"), samp], []) + samp = ir(Braket.Sample(ho, Int[2]), Val(:JAQCD)) + action = Braket.Program(Braket.header_dict[Braket.Program], + [Braket.IR.CNot(0, 1, "cnot"), Braket.IR.CNot(2, 3, "cnot")], + [Braket.IR.Probability([1], "probability"), + Braket.IR.Expectation(["z"], nothing, "expectation"), + Braket.IR.Variance(["x", "x"], [0, 2], "variance"), samp], + []) task_metadata_shots = Braket.TaskMetadata(Braket.header_dict[Braket.TaskMetadata], "task_arn", length(measurements), "arn1", nothing, nothing, nothing, nothing, nothing) additional_metadata = Braket.AdditionalMetadata(action, nothing, nothing, nothing, nothing, nothing, nothing) task_result = Braket.GateModelTaskResult(Braket.header_dict[Braket.GateModelTaskResult], diff --git a/test/integ_tests/adjoint_gradient.jl b/test/integ_tests/adjoint_gradient.jl new file mode 100644 index 00000000..993ef80c --- /dev/null +++ b/test/integ_tests/adjoint_gradient.jl @@ -0,0 +1,84 @@ +using Braket, Test + +sv1_arn() = "arn:aws:braket:::device/quantum-simulator/amazon/sv1" +@testset "Adjoint Gradient" begin + @testset "with nested targets" begin + θ = FreeParameter(:theta) + inputs = Dict("theta"=>0.2) + obs = (-2*Braket.Observables.Y()) * (3*Braket.Observables.I()) + 0.75 * Braket.Observables.Y() * Braket.Observables.Z() + circ = Circuit([(Rx, 0, θ), (AdjointGradient, obs, [[0, 1], [2, 3]], ["theta"])]) + + expected_openqasm = """ + OPENQASM 3.0; + input float theta; + qubit[4] q; + rx(theta) q[0]; + #pragma braket result adjoint_gradient expectation(-6.0 * y(q[0]) @ i(q[1]) + 0.75 * y(q[2]) @ z(q[3])) theta""" + gradient_task = Braket.AwsQuantumTask(sv1_arn(), circ, s3_destination_folder=s3_destination_folder, shots=0, inputs=inputs) + res = result(gradient_task) + @test res.additional_metadata.action.source == expected_openqasm + @test res.values == [Dict(:expectation=>1.1920159847703675, :gradient=>Dict(:theta=>5.880399467047451))] + @test first(res.result_types).type.observable == "-6.0 * y @ i + 0.75 * y @ z" + @test res.additional_metadata.action.inputs == inputs + end + @testset "with standard observable terms" begin + θ = FreeParameter(:theta) + inputs = Dict("theta"=>0.2) + obs = 2*Braket.Observables.X() + 3*Braket.Observables.Y() - Braket.Observables.Z() + circ = Circuit([(Rx, 0, θ), (AdjointGradient, obs, [[0], [1], [2]], ["theta"])]) + + expected_openqasm = """ + OPENQASM 3.0; + input float theta; + qubit[3] q; + rx(theta) q[0]; + #pragma braket result adjoint_gradient expectation(2.0 * x(q[0]) + 3.0 * y(q[1]) - 1.0 * z(q[2])) theta""" + + gradient_task = Braket.AwsQuantumTask(sv1_arn(), circ, s3_destination_folder=s3_destination_folder, shots=0, inputs=inputs) + res = result(gradient_task) + @test res.additional_metadata.action.source == expected_openqasm + @test res.values == [Dict(:expectation=>-1, :gradient=>Dict(:theta=>0))] + @test first(res.result_types).type.observable == "2.0 * x + 3.0 * y - 1.0 * z" + @test res.additional_metadata.action.inputs == inputs + end + @testset "with batch" begin + θ = FreeParameter(:theta) + inputs = Dict("theta"=>0.2) + obs1 = 2*Braket.Observables.Y() * 3 * Braket.Observables.I() + obs2 = -2*Braket.Observables.Y() * 3 * Braket.Observables.I() + 0.75 * Braket.Observables.Y() * Braket.Observables.Z() + circ1 = Circuit([(Rx, 0, θ), (AdjointGradient, obs1, [0, 1], ["theta"])]) + circ2 = Circuit([(Rx, 0, θ), (AdjointGradient, obs2, [[0, 1], [0, 1]], ["theta"])]) + expected_openqasm = [ + """OPENQASM 3.0; + input float theta; + qubit[2] q; + rx(theta) q[0]; + #pragma braket result adjoint_gradient expectation(6.0 * y(q[0]) @ i(q[1])) theta""", + """OPENQASM 3.0; + input float theta; + qubit[2] q; + rx(theta) q[0]; + #pragma braket result adjoint_gradient expectation(-6.0 * y(q[0]) @ i(q[1]) + 0.75 * y(q[0]) @ z(q[1])) theta""" + , + ] + expected_result_values = [ + [Dict(:expectation=>-1.1920159847703675, :gradient=>Dict(:theta=>-5.880399467047451))], + [Dict(:expectation=>1.0430139866740715, :gradient=>Dict(:theta=>5.145349533666519))], + ] + expected_observables = ["6.0 * y @ i", "-6.0 * y @ i + 0.75 * y @ z"] + + gradient_batch_tasks = Braket.AwsQuantumTaskBatch(sv1_arn(), + [circ1, circ2], + s3_destination_folder=s3_destination_folder, + shots=0, + inputs=inputs + ) + for i in 1:length(Braket.tasks(gradient_batch_tasks)) + res = result(Braket.tasks(gradient_batch_tasks)[i]) + @test res.additional_metadata.action.source == expected_openqasm[i] + @test res.values == expected_result_values[i] + @test first(res.result_types).type.observable == expected_observables[i] + @test res.additional_metadata.action.inputs == inputs + end + end +end diff --git a/test/integ_tests/runtests.jl b/test/integ_tests/runtests.jl index a23f0a75..9534dea5 100644 --- a/test/integ_tests/runtests.jl +++ b/test/integ_tests/runtests.jl @@ -3,9 +3,12 @@ using Braket, Test get_tol(shots::Int) = return (shots > 0 ? Dict("atol"=> 0.1, "rtol"=>0.15) : Dict("atol"=>0.01, "rtol"=>0)) s3_destination_folder = Braket.default_task_bucket() +include("adjoint_gradient.jl") +#= include("create_quantum_job.jl") include("cost_tracking.jl") include("device_creation.jl") include("density_matrix_simulator.jl") include("simulator_quantum_task.jl") -include("tensor_network_simulator.jl") \ No newline at end of file +include("tensor_network_simulator.jl") +=# diff --git a/test/integ_tests/simulator_quantum_task.jl b/test/integ_tests/simulator_quantum_task.jl index f2f25eff..ca074703 100644 --- a/test/integ_tests/simulator_quantum_task.jl +++ b/test/integ_tests/simulator_quantum_task.jl @@ -128,7 +128,7 @@ end device = AwsDevice(_arn=simulator_arn) res = result(device(task, shots=SHOTS, s3_destination_folder=s3_destination_folder)) @test length(res.result_types) == 1 - @test isapprox(res[Probability(qubits(circuit))], [0.5, 0.0, 0.0, 0.5], rtol=tol["rtol"], atol=tol["atol"]) + @test isapprox(res[Probability()], [0.5, 0.0, 0.0, 0.5], rtol=tol["rtol"], atol=tol["atol"]) end end @testset "Bell pair marginal probability" begin @@ -468,4 +468,4 @@ end @test length(res.measurements) == SHOTS end end -end \ No newline at end of file +end diff --git a/test/observables.jl b/test/observables.jl index f2457818..6404a680 100644 --- a/test/observables.jl +++ b/test/observables.jl @@ -1,8 +1,9 @@ using Braket, Braket.Observables, Test, JSON3, StructTypes, LinearAlgebra -using Braket: VIRTUAL, PHYSICAL, OpenQASMSerializationProperties, pauli_eigenvalues +using Braket: VIRTUAL, PHYSICAL, OpenQASMSerializationProperties, pauli_eigenvalues, IRObservable using LinearAlgebra: eigvals @testset "Observables" begin + Braket.IRType[] = :JAQCD @testset "pauli eigenvalues" begin z = [1.0 0.0; 0.0 -1.0] @test pauli_eigenvalues(1) == diag(z) @@ -15,11 +16,12 @@ using LinearAlgebra: eigvals @test qubit_count(o) == 1 rt = Expectation(o, [0]) @test JSON3.read(JSON3.write(rt), Braket.Result) == rt + @test ir(o) isa IRObservable end @testset "TensorProduct" begin tp = Observables.TensorProduct(["x", "y", "z"]) @test ir(tp) == ["x", "y", "z"] - @test JSON3.write(tp) == replace("""$(Braket.ir(tp))""", " "=>"") + @test JSON3.write(tp) == JSON3.write(ir(tp)) @test qubit_count(tp) == 3 rt = Expectation(tp, [0, 1, 2]) @test JSON3.read(JSON3.write(rt), Braket.Result) == rt @@ -33,6 +35,7 @@ using LinearAlgebra: eigvals @test qubit_count(tp) == 2 rt = Expectation(tp, [0, 2]) @test JSON3.read(JSON3.write(rt), Braket.Result) == rt + @test ir(tp) isa IRObservable end @testset "TensorProduct mixed types" begin m = [1. -im; im -1.] @@ -40,21 +43,36 @@ using LinearAlgebra: eigvals tp = Observables.TensorProduct([Observables.X(), o, Observables.Z()]) @test qubit_count(tp) == 3 rt = Expectation(tp, [0, 1, 2]) + @test ir(tp) isa IRObservable @test JSON3.read(JSON3.write(rt), Braket.Result) == rt end + @testset "TensorProduct doesn't accept Sum" begin + s = 2.0 * Observables.X() * Observables.X() + 3.0 * Observables.Z() * Observables.Z() + @test_throws ArgumentError Observables.TensorProduct([Observables.X(), s]) + end + @testset "Sum" begin + s = 2.0 * Observables.X() * Observables.X() + 3.0 * Observables.Z() * Observables.Z() + @test length(s) == 2 + @test s.summands[1].coefficient == 2.0 + @test s.summands[2].coefficient == 3.0 + s2 = -1.0 * s + @test length(s2) == 2 + @test s2.summands[1].coefficient == -2.0 + @test s2.summands[2].coefficient == -3.0 + end @test_throws ErrorException StructTypes.constructfrom(Observables.Observable, ["x", 1, "z"]) m = [1 -im; im -1] HO = Observables.HermitianObservable(m) - @test JSON3.write(HO) == replace("""$(Braket.ir(HO))""", " "=>"") - m_raw = first(ir(HO)) - @test StructTypes.constructfrom(Observables.Observable, m_raw) == HO + HO_ir = ir(HO) + @test JSON3.write(HO) == JSON3.write(ir(HO)) + @test StructTypes.constructfrom(Observables.Observable, convert(IRObservable, HO_ir)) == HO @test copy(HO) == HO for typ in (Observables.H, Observables.I, Observables.Z, Observables.X, Observables.Y) @test copy(typ()) == typ() - @test JSON3.write(typ()) == """$(Braket.ir(typ()))""" + @test JSON3.write(typ()) == JSON3.write(ir(typ())) @test ishermitian(typ()) end - + Braket.IRType[] = :OpenQASM @testset "eigenvalues" begin mat = [1.0 1.0 - im; 1.0 + 1im -1.0] h = Observables.HermitianObservable(mat) @@ -92,6 +110,8 @@ using LinearAlgebra: eigvals ( Observables.H(), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), [3], "h(q[3])"), ( Observables.H(), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), [3], "h(\$3)"), ( Observables.H(), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), nothing, "h all"), + ( 2.0 * Observables.H() - 5.0 * Observables.Z() * Observables.X(), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), [[2], [3, 4]], "2.0 * h(\$2) - 5.0 * z(\$3) @ x(\$4)"), + ( 2.0 * Observables.H() - 5.0 * Observables.Z() * Observables.X(), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), [[2], [3, 4]], "2.0 * h(q[2]) - 5.0 * z(q[3]) @ x(q[4])"), ( Observables.HermitianObservable(diagm(ones(Int64, 4))), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), [1, 2], "hermitian([[1+0im, 0im, 0im, 0im], [0im, 1+0im, 0im, 0im], [0im, 0im, 1+0im, 0im], [0im, 0im, 0im, 1+0im]]) q[1], q[2]"), ( Observables.HermitianObservable(diagm(ones(Int64, 4))), OpenQASMSerializationProperties(qubit_reference_type=PHYSICAL), [1, 2], "hermitian([[1+0im, 0im, 0im, 0im], [0im, 1+0im, 0im, 0im], [0im, 0im, 1+0im, 0im], [0im, 0im, 0im, 1+0im]]) \$1, \$2"), ( Observables.HermitianObservable(diagm(ones(Int64, 2))), OpenQASMSerializationProperties(qubit_reference_type=VIRTUAL), nothing, "hermitian([[1+0im, 0im], [0im, 1+0im]]) all"), diff --git a/test/schemas_misc.jl b/test/schemas_misc.jl index be1f9305..ac68db27 100644 --- a/test/schemas_misc.jl +++ b/test/schemas_misc.jl @@ -251,4 +251,4 @@ using Braket, Braket.IR, Test, JSON3, StructTypes read_in = JSON3.read(str, Braket.BraketSchemaBase) @test read_in isa Braket.GateModelParameters end -end \ No newline at end of file +end diff --git a/test/task.jl b/test/task.jl index d1269a78..f33fd888 100644 --- a/test/task.jl +++ b/test/task.jl @@ -119,7 +119,7 @@ zero_shots_result(task_mtd, add_mtd) = Braket.GateModelTaskResult( device_params = Dict("fake_param_1"=>2, "fake_param_2"=>"hello") s3_folder = ("fake_bucket", "fake_folder") task_args = Braket.prepare_task_input(program(), arn, s3_folder, shots, device_params) - @test task_args[:action] == ir(program()) + @test task_args[:action] == JSON3.write(ir(program())) @test task_args[:device_arn] == arn @test UUID(task_args[:client_token]) isa UUID @test task_args[:shots] == shots @@ -132,7 +132,7 @@ zero_shots_result(task_mtd, add_mtd) = Braket.GateModelTaskResult( device_params = Dict("fake_param_1"=>2, "fake_param_2"=>"hello") s3_folder = ("fake_bucket", "fake_folder") task_args = Braket.prepare_task_input(program(), arn, s3_folder, shots, device_params) - @test task_args[:action] == JSON3.write(program()) + @test task_args[:action] == JSON3.write(ir(program())) @test task_args[:device_arn] == arn @test UUID(task_args[:client_token]) isa UUID @test task_args[:shots] == shots diff --git a/test/translation.jl b/test/translation.jl index 62895cea..e8c9f1e5 100644 --- a/test/translation.jl +++ b/test/translation.jl @@ -55,7 +55,7 @@ using Braket: operator, target @testset "Results translation" begin json_str = """[{"observable":["i"],"targets":[0],"type":"expectation"},{"observable":["z"],"targets":[0],"type":"expectation"}]""" parsed = JSON3.read(json_str, Vector{Braket.AbstractProgramResult}) - raw = Braket.AbstractProgramResult[Braket.IR.Expectation(["i"], [0], "expectation"), Braket.IR.Expectation(["z"], [0], "expectation")] + raw = Braket.AbstractProgramResult[Braket.IR.Expectation(convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, ["i"]), [0], "expectation"), Braket.IR.Expectation(convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, ["z"]), [0], "expectation")] for (e_r, e_p) in zip(raw, parsed) @test e_r == e_p end @@ -70,14 +70,14 @@ using Braket: operator, target @testset "Program translation" begin json_str = """{"braketSchemaHeader": {"name": "braket.ir.jaqcd.program", "version": "1"}, "instructions": [{"target": 0, "type": "x"}, {"target": 0, "type": "z"}], "results": [{"observable":["i"],"targets":[0],"type":"expectation"},{"observable":["z"],"targets":[0],"type":"expectation"}], "basis_rotation_instructions": []}""" parsed = Braket.parse_raw_schema(json_str) - raw = Braket.Program(Braket.braketSchemaHeader("braket.ir.jaqcd.program", "1"), [Braket.Instruction(X(), 0), Braket.Instruction(Z(), 0)], Braket.AbstractProgramResult[Braket.IR.Expectation(["i"], [0], "expectation"), Braket.IR.Expectation(["z"], [0], "expectation")], Braket.Instruction[]) + raw = Braket.Program(Braket.braketSchemaHeader("braket.ir.jaqcd.program", "1"), [Braket.Instruction(X(), 0), Braket.Instruction(Z(), 0)], Braket.AbstractProgramResult[Braket.IR.Expectation(convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, ["i"]), [0], "expectation"), Braket.IR.Expectation(convert(Vector{Union{String, Vector{Vector{Vector{Float64}}}}}, ["z"]), [0], "expectation")], Braket.Instruction[]) @test parsed.braketSchemaHeader == raw.braketSchemaHeader @test parsed.instructions == raw.instructions @test parsed.results == raw.results @test parsed.basis_rotation_instructions == raw.basis_rotation_instructions @test Braket.parse_raw_schema(JSON3.write(raw)) == raw @test Braket.Program(Circuit(raw)) == raw - @test ir(raw) == ir(Circuit(raw)) + @test ir(raw) == ir(Circuit(raw), Val(:JAQCD)) noisy_json_str = """{"braketSchemaHeader": {"name": "braket.ir.jaqcd.program", "version": "1"}, "instructions": [{"target": 0, "type": "x"}, {"gamma": 0.1, "target": 0, "type": "amplitude_damping"}, {"target": 1, "type": "y"}, {"control": 0, "target": 2, "type": "cnot"}, {"gamma": 0.1, "target": 0, "type": "amplitude_damping"}, {"gamma": 0.1, "target": 2, "type": "amplitude_damping"}, {"target": 1, "type": "x"}, {"target": 2, "type": "z"}, {"gamma": 0.1, "target": 2, "type": "amplitude_damping"}], "results": [{"targets": null, "type": "probability"}, {"observable": ["z"], "targets": [0], "type": "expectation"}, {"targets": [0, 1], "type": "densitymatrix"}], "basis_rotation_instructions": []}""" @@ -95,12 +95,11 @@ using Braket: operator, target for (prt, rrt) in zip(p.results, raw.results) @test prt == rrt end - - @test JSON3.read(ir(raw), Dict) == JSON3.read(ir(Circuit(raw)), Dict) + @test JSON3.read(JSON3.write(raw), Dict) == JSON3.read(JSON3.write(ir(Circuit(raw), Val(:JAQCD))), Dict) end @testset "enum translation" begin - for e in [Braket._CopyMode, Braket.SafeUUID, Braket.PaymentCardBrand, + for e in [Braket.SafeUUID, Braket.PaymentCardBrand, Braket.Protocol, Braket.Extra, Braket.DeviceActionType, Braket.ExecutionDay, Braket.QubitDirection, Braket.PostProcessingType, Braket.ResultFormat,