1- using LinearAlgebra, SimpleDifferentialOperators, DiffEqOperators, Test, Parameters
1+ using LinearAlgebra, DiffEqOperators, Test, Parameters
2+
3+ # These tests are based on examples from SimpleDifferentialOperators.jl
24
35# parameters
46params = @with_kw (
@@ -7,7 +9,7 @@ params = @with_kw (
79 ρ = 0.05 ,
810 M = 3 , # size of grid (interior points)
911 x̄ = range (0.0 , 1.0 , length = (M+ 2 )),
10- x = interiornodes (x̄), # i.e., x̄ [2:end-1]
12+ x = x̄ [2 : end - 1 ],
1113 S = 3.0
1214)
1315p = params ()
@@ -17,14 +19,14 @@ p = params()
1719# ----------------------------------
1820# payoff function
1921pi_profit (x) = x^ 2
20- # SimpleDifferentialOperators setup
22+ #= # SimpleDifferentialOperators setup
2123function SDO_negative_drift(pi_profit, params)
2224 bc = (Reflecting(), Reflecting())
2325 Lₓ = params.μ*L₁₋bc(params.x̄, bc) + params.σ^2 / 2 * L₂bc(params.x̄, bc)
2426 L_bc = I * params.ρ - Lₓ
2527 v = L_bc \ pi_profit.(params.x);
2628 return (v = v, L_bc = L_bc, Lₓ = Lₓ, L₁₋bc = L₁₋bc(params.x̄, bc), L₂bc = L₂bc(params.x̄, bc))
27- end
29+ end=#
2830# DiffEqOperators Setup
2931function DEO_negative_drift (pi_profit, params)
3032 dx = params. x[2 ] - params. x[1 ]
@@ -45,25 +47,26 @@ function DEO_negative_drift(pi_profit, params)
4547 v = L_bc \ pi_profit .(params. x);
4648 return (v = v, L_bc = L_bc, Lₓ = Lₓ, L₁₋bc = L₁₋bc, L₂bc = Array (L2* Q)[1 ] )
4749end
50+
4851@testset " Constant Negative Drifts" begin
49- @test DEO_negative_drift (pi_profit, p). v ≈ SDO_negative_drift (pi_profit, p). v
50- @test DEO_negative_drift (pi_profit, p). L_bc ≈ SDO_negative_drift (pi_profit, p). L_bc
51- @test DEO_negative_drift (pi_profit, p). Lₓ ≈ SDO_negative_drift (pi_profit, p). Lₓ
52- @test DEO_negative_drift (pi_profit, p). L₁₋bc ≈ SDO_negative_drift (pi_profit, p). L₁₋bc
53- @test DEO_negative_drift (pi_profit, p). L₂bc ≈ SDO_negative_drift (pi_profit, p). L₂bc
52+ @test DEO_negative_drift (pi_profit, p). v ≈ [ 1.9182649086005406 , 2.3359304764758773 , 3.1768804315253277 ] # SDO_negative_drift(pi_profit, p).v
53+ @test DEO_negative_drift (pi_profit, p). L_bc ≈ [ 0.13 - 0.08 0. ; - 0.48 0.61 - 0.08 ; 0. - 0.48 0.53 ] # SDO_negative_drift(pi_profit, p).L_bc
54+ @test DEO_negative_drift (pi_profit, p). Lₓ ≈ [ - 0.08 0.08 0. ; . 48 - . 56 . 08 ; 0. . 48 - . 48 ] # SDO_negative_drift(pi_profit, p).Lₓ
55+ @test DEO_negative_drift (pi_profit, p). L₁₋bc ≈ [ 0. 0. 0. ; - 4. 4. 0. ; 0. - 4. 4. ] # SDO_negative_drift(pi_profit, p).L₁₋bc
56+ @test DEO_negative_drift (pi_profit, p). L₂bc ≈ [ - 16. 16. 0. ; 16. - 32. 16. ; 0. 16. - 16. ] # SDO_negative_drift(pi_profit, p).L₂bc
5457end
5558p = params (μ = 0.1 );
5659# ----------------------------------
5760# Testing For Positive Drifts
5861# ----------------------------------
5962# SimpleDifferentialOperators setup
60- function SDO_positive_drift (pi_profit, params)
63+ #= function SDO_positive_drift(pi_profit, params)
6164 bc = (Reflecting(), Reflecting())
6265 Lₓ = params.μ*L₁₊bc(params.x̄, bc) + params.σ^2 / 2 * L₂bc(params.x̄, bc)
6366 L_bc = I * params.ρ - Lₓ
6467 v = L_bc \ pi_profit.(params.x);
6568 return (v = v, L_bc = L_bc, Lₓ = Lₓ, L₁₊bc = L₁₊bc(params.x̄, bc), L₂bc = L₂bc(params.x̄, bc))
66- end
69+ end=#
6770# DiffEqOperators Setup
6871function DEO_positive_drift (pi_profit, params)
6972 dx = params. x[2 ] - params. x[1 ]
@@ -84,11 +87,11 @@ function DEO_positive_drift(pi_profit, params)
8487end
8588
8689@testset " Constant Positive Drifts" begin
87- @test DEO_positive_drift (pi_profit, p). v ≈ SDO_positive_drift (pi_profit, p). v
88- @test DEO_positive_drift (pi_profit, p). L_bc ≈ SDO_positive_drift (pi_profit, p). L_bc
89- @test DEO_positive_drift (pi_profit, p). Lₓ ≈ SDO_positive_drift (pi_profit, p). Lₓ
90- @test DEO_positive_drift (pi_profit, p). L₁₊bc ≈ SDO_positive_drift (pi_profit, p). L₁₊bc
91- @test DEO_positive_drift (pi_profit, p). L₂bc ≈ SDO_positive_drift (pi_profit, p). L₂bc
90+ @test DEO_positive_drift (pi_profit, p). v ≈ [ 8.855633802816902 , 9.647887323943664 , 10.264084507042257 ] # SDO_positive_drift(pi_profit, p).v
91+ @test DEO_positive_drift (pi_profit, p). L_bc ≈ [ 0.53 - 0.48 0. ; - 0.08 0.61 - 0.48 ; 0. - 0.08 0.13 ] # SDO_positive_drift(pi_profit, p).L_bc
92+ @test DEO_positive_drift (pi_profit, p). Lₓ ≈ [ - 0.48 0.48 0. ; 0.08 - 0.56 0.48 ; 0. 0.08 - 0.08 ] # SDO_positive_drift(pi_profit, p).Lₓ
93+ @test DEO_positive_drift (pi_profit, p). L₁₊bc ≈ [ - 4. 4. 0. ; 0. - 4. 4. ; 0. 0. 0. ] # SDO_positive_drift(pi_profit, p).L₁₊bc
94+ @test DEO_positive_drift (pi_profit, p). L₂bc ≈ [ - 16. 16. 0. ; 16. - 32. 16. ; 0. 16. - 16. ] # SDO_positive_drift(pi_profit, p).L₂bc
9295end
9396
9497
97100# ----------------------------------
98101μ (x) = - x
99102# SimpleDifferentialOperators setup
100- function SDO_state_dependent_drift (pi_profit, μ, params)
103+ #= function SDO_state_dependent_drift(pi_profit, μ, params)
101104 bc = (Reflecting(), Reflecting())
102105 L₁ = Diagonal(min.(μ.(params.x), 0.0)) * L₁₋bc(params.x̄, bc) + Diagonal(max.(μ.(params.x), 0.0)) * L₁₊bc(params.x̄, bc)
103106 Lₓ = L₁ - params.σ^2 / 2 * L₂bc(params.x̄, bc)
104107 L_bc = I * params.ρ - Lₓ
105108 v = L_bc \ pi_profit.(params.x)
106109 return (v = v, L_bc = L_bc, Lₓ = Lₓ, L₁ = L₁, L₂bc = L₂bc(params.x̄, bc))
107- end
110+ end=#
108111# DiffEqOperators Setup
109112function DEO_state_dependent_drift (pi_profit, μ, params)
110113 dx = params. x[2 ] - params. x[1 ]
@@ -127,11 +130,11 @@ function DEO_state_dependent_drift(pi_profit, μ, params)
127130end
128131
129132@testset " State Dependent Drifts" begin
130- @test DEO_state_dependent_drift (pi_profit, μ, p). v ≈ SDO_state_dependent_drift (pi_profit, μ, p). v
131- @test DEO_state_dependent_drift (pi_profit, μ, p). L_bc ≈ SDO_state_dependent_drift (pi_profit, μ, p). L_bc
132- @test DEO_state_dependent_drift (pi_profit, μ, p). Lₓ ≈ SDO_state_dependent_drift (pi_profit, μ, p). Lₓ
133- @test DEO_state_dependent_drift (pi_profit, μ, p). L₁ ≈ SDO_state_dependent_drift (pi_profit, μ, p). L₁
134- @test DEO_state_dependent_drift (pi_profit, μ, p). L₂bc ≈ SDO_state_dependent_drift (pi_profit, μ, p). L₂bc
133+ @test DEO_state_dependent_drift (pi_profit, μ, p). v ≈ [ 1.1027342984846056 , 1.194775361931727 , 1.3640552379934825 ] # SDO_state_dependent_drift(pi_profit, μ, p).v
134+ @test DEO_state_dependent_drift (pi_profit, μ, p). L_bc ≈ [ - 0.03 0.08 0. ; - 1.92 1.89 0.08 ; 0. - 2.92 2.97 ] # SDO_state_dependent_drift(pi_profit, μ, p).L_bc
135+ @test DEO_state_dependent_drift (pi_profit, μ, p). Lₓ ≈ [ 0.08 - 0.08 0. ; 1.92 - 1.84 - 0.08 ; 0. 2.92 - 2.92 ] # SDO_state_dependent_drift(pi_profit, μ, p).Lₓ
136+ @test DEO_state_dependent_drift (pi_profit, μ, p). L₁ ≈ [ 0.0 0.0 0. ; 2.0 - 2.0 0.0 ; 0. 3.0 - 3.0 ] # SDO_state_dependent_drift(pi_profit, μ, p).L₁
137+ @test DEO_state_dependent_drift (pi_profit, μ, p). L₂bc ≈ [ - 16. 16. 0. ; 16. - 32. 16. ; 0. 16. - 16. ] # SDO_state_dependent_drift(pi_profit, μ, p).L₂bc
135138end
136139
137140
144147# payoff function
145148pi_profit (x) = x^ 2
146149# SimpleDifferentialOperators setup
147- function SDO_absorbing_bc (pi_profit, params)
150+ #= function SDO_absorbing_bc(pi_profit, params)
148151 bc = (NonhomogeneousAbsorbing(params.S), Reflecting())
149152 Lₓbc = params.μ*L₁₋bc(params.x̄, bc) + params.σ^2 / 2 * L₂bc(params.x̄ , bc)
150153 L_bc = I * params.ρ - Lₓbc
@@ -157,8 +160,7 @@ function SDO_absorbing_bc(pi_profit, params)
157160
158161 return (v = v, L_bc = L_bc, Lₓbc = Lₓbc, L₁₋bc = L₁₋bc(params.x̄, bc), L₂bc = L₂bc(params.x̄, bc),
159162 pi_profit_star = pi_profit_star, pi_profit = pi_profit.(params.x))
160- end
161-
163+ end=#
162164
163165# DiffEqOperators Setup
164166function DEO_absorbing_bc (pi_profit, params)
@@ -186,13 +188,13 @@ end
186188p = params (x̄ = range (0.0 , 1.0 , length = (p. M+ 2 )))
187189
188190@testset " Absorbing BC" begin
189- @test DEO_absorbing_bc (pi_profit, p). v ≈ SDO_absorbing_bc (pi_profit, p). v
190- @test DEO_absorbing_bc (pi_profit, p). L_bc ≈ SDO_absorbing_bc (pi_profit, p). L_bc
191- @test DEO_absorbing_bc (pi_profit, p). Lₓbc ≈ SDO_absorbing_bc (pi_profit, p). Lₓbc
192- @test DEO_absorbing_bc (pi_profit, p). L₁₋bc ≈ SDO_absorbing_bc (pi_profit, p). L₁₋bc
193- @test DEO_absorbing_bc (pi_profit, p). L₂bc ≈ SDO_absorbing_bc (pi_profit, p). L₂bc
194- @test DEO_absorbing_bc (pi_profit, p). pi_profit_star ≈ SDO_absorbing_bc (pi_profit, p). pi_profit_star
195- @test DEO_absorbing_bc (pi_profit, p). pi_profit ≈ SDO_absorbing_bc (pi_profit, p). pi_profit
191+ @test DEO_absorbing_bc (pi_profit, p). v ≈ [ 0.2085953844248779 , 0.8092898062396943 , 1.7942624660284021 ] # SDO_absorbing_bc(pi_profit, p).v
192+ @test DEO_absorbing_bc (pi_profit, p). L_bc ≈ [ 0.61 - 0.08 0. ; - 0.48 0.61 - 0.08 ; 0. - 0.48 0.53 ] # SDO_absorbing_bc(pi_profit, p).L_bc
193+ @test DEO_absorbing_bc (pi_profit, p). Lₓbc ≈ [ - 0.56 0.08 0. ; 0.48 - 0.56 0.08 ; 0. 0.48 - 0.48 ] # SDO_absorbing_bc(pi_profit, p).Lₓbc
194+ @test DEO_absorbing_bc (pi_profit, p). L₁₋bc ≈ [ 4.0 0.0 0. ; - 4.0 4.0 0.0 ; 0. - 4.0 4.0 ] # SDO_absorbing_bc(pi_profit, p).L₁₋bc
195+ @test DEO_absorbing_bc (pi_profit, p). L₂bc ≈ [ - 32.0 16.0 0. ; 16.0 - 32.0 16.0 ; 0. 16.0 - 16.0 ] # SDO_absorbing_bc(pi_profit, p).L₂bc
196+ @test DEO_absorbing_bc (pi_profit, p). pi_profit_star ≈ [. 0625 , . 25 , . 5625 ] # SDO_absorbing_bc(pi_profit, p).pi_profit_star
197+ @test DEO_absorbing_bc (pi_profit, p). pi_profit ≈ [. 0625 , . 25 , . 5625 ] # SDO_absorbing_bc(pi_profit, p).pi_profit
196198end
197199
198200
202204# ----------------------------------
203205
204206
205- function SDO_Solve_KFE (params)
207+ #= function SDO_Solve_KFE(params)
206208 # ξ values for mixed boundary conditions
207209 ξ_lb = ξ_ub = -2 * params.μ / params.σ^2
208210 # define the corresponding mixed boundary conditions
@@ -214,9 +216,9 @@ function SDO_Solve_KFE(params)
214216 L_KFE_with_drift = Array(-params.μ*L₁₊bc(params.x̄, bc) + params.σ^2 / 2 * L₂bc(params.x̄, bc))
215217 L_KFE_without = params.σ^2 / 2 * L₂bc(params.x̄, bc)
216218 return (L_KFE_with_drift = L_KFE_with_drift, L_KFE_without = L_KFE_without)
217- end
219+ end=#
218220
219- function SDO_Solve_KFE_forward (params)
221+ #= function SDO_Solve_KFE_forward(params)
220222 # ξ values for mixed boundary conditions
221223 ξ_lb = ξ_ub = -2 * params.μ / params.σ^2
222224 # define the corresponding mixed boundary conditions
@@ -227,7 +229,7 @@ function SDO_Solve_KFE_forward(params)
227229 L_KFE_with_drift = Array(-params.μ*L₁₊bc(params.x̄, bc) + params.σ^2 / 2 * L₂bc(params.x̄, bc))
228230 L_KFE_without = params.σ^2 / 2 * L₂bc(params.x̄, bc)
229231 return (L_KFE_with_drift = L_KFE_with_drift, L_KFE_without = L_KFE_without)
230- end
232+ end=#
231233
232234
233235function DEO_Solve_KFE (params)
@@ -256,13 +258,12 @@ function DEO_Solve_KFE(params)
256258 return (L_KFE_with_drift = L_KFE_with_drift, L_KFE_without = L_KFE_without)
257259end
258260
259-
260261p = params (x̄ = range (0.0 , 1.0 , length = (p. M+ 2 )))
261262@testset " Solving KFE" begin
262- @test_broken DEO_Solve_KFE (p). L_KFE_with_drift ≈ SDO_Solve_KFE (p). L_KFE_with_drift # concretization of Robin conditions appears broken
263- @test_broken DEO_Solve_KFE (p). L_KFE_without ≈ SDO_Solve_KFE (p). L_KFE_without # Hard to check what's wrong since KFE uses forward/backward at boundaries rather than central differences when w/mixed bcs
264- @test DEO_Solve_KFE (p). L_KFE_with_drift ≈ SDO_Solve_KFE_forward (p). L_KFE_with_drift
265- @test DEO_Solve_KFE (p). L_KFE_without ≈ SDO_Solve_KFE_forward (p). L_KFE_without
263+ # @test_broken DEO_Solve_KFE(p).L_KFE_with_drift ≈ SDO_Solve_KFE(p).L_KFE_with_drift # concretization of Robin conditions appears broken
264+ # @test_broken DEO_Solve_KFE(p).L_KFE_without ≈ SDO_Solve_KFE(p).L_KFE_without # Hard to check what's wrong since KFE uses forward/backward at boundaries rather than central differences when w/mixed bcs
265+ @test DEO_Solve_KFE (p). L_KFE_with_drift ≈ [ - 0.58 0.48 0.0 ; 0.08 - 0.56 0.48 ; 0.0 0.08 - 0.48 ] # SDO_Solve_KFE_forward(p).L_KFE_with_drift
266+ @test DEO_Solve_KFE (p). L_KFE_without ≈ [ - 0.18 0.08 0. ; 0.08 - 0.16 0.08 ; 0. 0.08 - 0.1466666666666667 ] # SDO_Solve_KFE_forward(p).L_KFE_without
266267end
267268
268269nothing
0 commit comments