|
1 | | -@testset "Autodiff" verbose=true begin |
2 | | - @testset "sesolve" verbose=true begin |
3 | | - ψ0 = fock(2, 1) |
4 | | - t_max = 10 |
5 | | - tlist = range(0, t_max, 100) |
| 1 | +# ---- SESOLVE ---- |
| 2 | +const ψ0_sesolve = fock(2, 1) |
| 3 | +t_max = 10 |
| 4 | +const tlist_sesolve = range(0, t_max, 100) |
6 | 5 |
|
7 | | - # For direct Forward differentiation |
8 | | - function my_f_sesolve_direct(p) |
9 | | - H = p[1] * sigmax() |
10 | | - sol = sesolve(H, ψ0, tlist, progress_bar = Val(false)) |
| 6 | +# For direct Forward differentiation |
| 7 | +function my_f_sesolve_direct(p) |
| 8 | + H = p[1] * sigmax() |
| 9 | + sol = sesolve(H, ψ0_sesolve, tlist_sesolve, progress_bar = Val(false)) |
11 | 10 |
|
12 | | - return real(expect(projection(2, 0, 0), sol.states[end])) |
13 | | - end |
| 11 | + return real(expect(projection(2, 0, 0), sol.states[end])) |
| 12 | +end |
14 | 13 |
|
15 | | - # For SciMLSensitivity.jl |
16 | | - coef_Ω(p, t) = p[1] |
17 | | - H_evo = QobjEvo(sigmax(), coef_Ω) |
18 | | - |
19 | | - function my_f_sesolve(p) |
20 | | - sol = sesolve( |
21 | | - H_evo, |
22 | | - ψ0, |
23 | | - tlist, |
24 | | - progress_bar = Val(false), |
25 | | - params = p, |
26 | | - sensealg = BacksolveAdjoint(autojacvec = EnzymeVJP()), |
27 | | - ) |
28 | | - |
29 | | - return real(expect(projection(2, 0, 0), sol.states[end])) |
30 | | - end |
| 14 | +# For SciMLSensitivity.jl |
| 15 | +coef_Ω(p, t) = p[1] |
| 16 | +const H_evo = QobjEvo(sigmax(), coef_Ω) |
| 17 | + |
| 18 | +function my_f_sesolve(p) |
| 19 | + sol = sesolve( |
| 20 | + H_evo, |
| 21 | + ψ0_sesolve, |
| 22 | + tlist_sesolve, |
| 23 | + progress_bar = Val(false), |
| 24 | + params = p, |
| 25 | + sensealg = BacksolveAdjoint(autojacvec = EnzymeVJP()), |
| 26 | + ) |
| 27 | + |
| 28 | + return real(expect(projection(2, 0, 0), sol.states[end])) |
| 29 | +end |
31 | 30 |
|
32 | | - # Analytical solution |
33 | | - my_f_analytic(Ω) = abs2(sin(Ω * t_max)) |
34 | | - my_f_analytic_deriv(Ω) = 2 * t_max * sin(Ω * t_max) * cos(Ω * t_max) |
| 31 | +# Analytical solution |
| 32 | +my_f_analytic(Ω) = abs2(sin(Ω * t_max)) |
| 33 | +my_f_analytic_deriv(Ω) = 2 * t_max * sin(Ω * t_max) * cos(Ω * t_max) |
| 34 | + |
| 35 | +# ---- MESOLVE ---- |
| 36 | +const N = 20 |
| 37 | +const a = destroy(N) |
| 38 | +const ψ0_mesolve = fock(N, 0) |
| 39 | +const tlist_mesolve = range(0, 40, 100) |
| 40 | + |
| 41 | +# For direct Forward differentiation |
| 42 | +function my_f_mesolve_direct(p) |
| 43 | + H = p[1] * a' * a + p[2] * (a + a') |
| 44 | + c_ops = [sqrt(p[3]) * a] |
| 45 | + sol = mesolve(H, ψ0_mesolve, tlist_mesolve, c_ops, progress_bar = Val(false)) |
| 46 | + return real(expect(a' * a, sol.states[end])) |
| 47 | +end |
| 48 | + |
| 49 | +# For SciMLSensitivity.jl |
| 50 | +coef_Δ(p, t) = p[1] |
| 51 | +coef_F(p, t) = p[2] |
| 52 | +coef_γ(p, t) = sqrt(p[3]) |
| 53 | +H = QobjEvo(a' * a, coef_Δ) + QobjEvo(a + a', coef_F) |
| 54 | +c_ops = [QobjEvo(a, coef_γ)] |
| 55 | +const L = liouvillian(H, c_ops) |
| 56 | + |
| 57 | +function my_f_mesolve(p) |
| 58 | + sol = mesolve( |
| 59 | + L, |
| 60 | + ψ0_mesolve, |
| 61 | + tlist_mesolve, |
| 62 | + progress_bar = Val(false), |
| 63 | + params = p, |
| 64 | + sensealg = BacksolveAdjoint(autojacvec = EnzymeVJP()), |
| 65 | + ) |
| 66 | + |
| 67 | + return real(expect(a' * a, sol.states[end])) |
| 68 | +end |
35 | 69 |
|
| 70 | +# Analytical solution |
| 71 | +n_ss(Δ, F, γ) = abs2(F / (Δ + 1im * γ / 2)) |
| 72 | + |
| 73 | +@testset "Autodiff" verbose=true begin |
| 74 | + @testset "sesolve" verbose=true begin |
36 | 75 | Ω = 1.0 |
37 | 76 | params = [Ω] |
38 | 77 |
|
|
52 | 91 |
|
53 | 92 | @test grad_qt ≈ grad_exact atol=1e-6 |
54 | 93 | end |
55 | | - end |
56 | 94 |
|
57 | | - @testset "mesolve" verbose=true begin |
58 | | - N = 20 |
59 | | - a = destroy(N) |
60 | | - ψ0 = fock(N, 0) |
61 | | - tlist = range(0, 40, 100) |
62 | | - |
63 | | - # For direct Forward differentiation |
64 | | - function my_f_mesolve_direct(p) |
65 | | - H = p[1] * a' * a + p[2] * (a + a') |
66 | | - c_ops = [sqrt(p[3]) * a] |
67 | | - sol = mesolve(H, ψ0, tlist, c_ops, progress_bar = Val(false)) |
68 | | - return real(expect(a' * a, sol.states[end])) |
69 | | - end |
| 95 | + @testset "Enzyme.jl" begin |
| 96 | + dparams = Enzyme.make_zero(params) |
| 97 | + Enzyme.autodiff( |
| 98 | + Enzyme.set_runtime_activity(Enzyme.Reverse), |
| 99 | + my_f_sesolve, |
| 100 | + Active, |
| 101 | + Duplicated(params, dparams), |
| 102 | + )[1] |
70 | 103 |
|
71 | | - # For SciMLSensitivity.jl |
72 | | - coef_Δ(p, t) = p[1] |
73 | | - coef_F(p, t) = p[2] |
74 | | - coef_γ(p, t) = sqrt(p[3]) |
75 | | - H = QobjEvo(a' * a, coef_Δ) + QobjEvo(a + a', coef_F) |
76 | | - c_ops = [QobjEvo(a, coef_γ)] |
77 | | - L = liouvillian(H, c_ops) |
78 | | - |
79 | | - function my_f_mesolve(p) |
80 | | - sol = mesolve( |
81 | | - L, |
82 | | - ψ0, |
83 | | - tlist, |
84 | | - progress_bar = Val(false), |
85 | | - params = p, |
86 | | - sensealg = BacksolveAdjoint(autojacvec = EnzymeVJP()), |
87 | | - ) |
88 | | - |
89 | | - return real(expect(a' * a, sol.states[end])) |
| 104 | + @test dparams ≈ grad_exact atol=1e-6 |
90 | 105 | end |
| 106 | + end |
91 | 107 |
|
92 | | - # Analytical solution |
93 | | - n_ss(Δ, F, γ) = abs2(F / (Δ + 1im * γ / 2)) |
94 | | - |
| 108 | + @testset "mesolve" verbose=true begin |
95 | 109 | Δ = 1.0 |
96 | 110 | F = 1.0 |
97 | 111 | γ = 1.0 |
|
111 | 125 | grad_qt = Zygote.gradient(my_f_mesolve, params)[1] |
112 | 126 | @test grad_qt ≈ grad_exact atol=1e-6 |
113 | 127 | end |
| 128 | + |
| 129 | + @testset "Enzyme.jl" begin |
| 130 | + dparams = Enzyme.make_zero(params) |
| 131 | + Enzyme.autodiff( |
| 132 | + Enzyme.set_runtime_activity(Enzyme.Reverse), |
| 133 | + my_f_mesolve, |
| 134 | + Active, |
| 135 | + Duplicated(params, dparams), |
| 136 | + )[1] |
| 137 | + |
| 138 | + @test dparams ≈ grad_exact atol=1e-6 |
| 139 | + end |
114 | 140 | end |
115 | 141 | end |
0 commit comments