Skip to content
This repository was archived by the owner on Jul 19, 2023. It is now read-only.

Commit b9f5835

Browse files
add second order check
1 parent db4152a commit b9f5835

2 files changed

Lines changed: 52 additions & 0 deletions

File tree

test/2nd_order_check.jl

Lines changed: 51 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,51 @@
1+
using DiffEqOperators, Base.Test
2+
3+
ord = 2 # Order of approximation
4+
Δx = 0.025π
5+
N = Int(2*/Δx)) -2
6+
x = -π:Δx:π-Δx
7+
D1 = DerivativeOperator{Float64}(1,ord,Δx,N,:periodic,:periodic); # 1nd Derivative
8+
D2 = DerivativeOperator{Float64}(2,ord,Δx,N,:periodic,:periodic); # 2nd Derivative
9+
10+
u0 = cos(x)
11+
du_true = -cos(x)
12+
13+
M = full(Tridiagonal([1.0 for i in 1:N+1],[-2.0 for i in 1:N+2],[1.0 for i in 1:N+1]))
14+
# Do the reflections, different for x and y operators
15+
M[end,1] = 1.0
16+
M[1,end] = 1.0
17+
M = M/(Δx^2)
18+
@test M*u0 D2*u0
19+
20+
A = zeros(M)
21+
for i in 1:N+2, j in 1:N+2
22+
i == j && (A[i,j] = -0.5)
23+
abs(i - j) == 2 && (A[i,j] = 0.25)
24+
end
25+
A[end-1,1] = 0.25
26+
A[end,2] = 0.25
27+
A[1,end-1] = 0.25
28+
A[2,end] = 0.25
29+
A = A/(Δx^2)
30+
31+
@test A*u0 D1*(D1*u0)
32+
33+
κ(x) = 1.0
34+
const tmp1 = similar(x)
35+
const tmp2 = similar(x)
36+
function f(t,u,du)
37+
A_mul_B!(tmp1, D1, u)
38+
@. tmp2 = κ(x)*tmp1
39+
A_mul_B!(du,D1,tmp2)
40+
end
41+
42+
du = similar(u0)
43+
f(0,u0,du)
44+
du A*u0
45+
46+
error = du_true - du
47+
du2 = D2*u0
48+
error2 = du_true - du2
49+
50+
@test maximum(error) < 0.004
51+
@test maximum(error2) < 0.001

test/runtests.jl

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,6 +7,7 @@ tic()
77
@time @testset "Dirichlet BCs" begin include("dirichlet.jl") end
88
@time @testset "Periodic BCs" begin include("periodic.jl") end
99
@time @testset "Neumann BCs" begin include("neumann.jl") end
10+
@time @testset "2nd order check" begin include("2nd_order_check.jl") end
1011
@time @testset "None BCs" begin include("none.jl") end
1112
@time @testset "KdV" begin include("KdV.jl") end # KdV times out
1213
@time @testset "Heat Equation" begin include("heat_eqn.jl") end

0 commit comments

Comments
 (0)