Skip to content

Commit

Permalink
fix issue 208 (#209)
Browse files Browse the repository at this point in the history 
* fix issue 208

* bump version

* adapt tests and fix typo in second-order fourth-derivative operator
  • Loading branch information
@ranocha
ranocha committed Aug 11, 2023
1 parent 02ac046 commit 108b1d3
Show file tree
Hide file tree
Showing 4 changed files with 42 additions and 35 deletions.
2 changes: 1 addition & 1 deletion Project.toml
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "SummationByPartsOperators"
uuid = "9f78cca6-572e-554e-b819-917d2f1cf240"
author = ["Hendrik Ranocha"]
version = "0.5.43"
version = "0.5.44"

[deps]
ArgCheck = "dce04be8-c92d-5529-be00-80e4d2c0e197"
Expand Down
2 changes: 1 addition & 1 deletion dev/Mattsson2014_dev.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@ R = vcat(r1, r2, r3, r4, r5, r6, r7)
n1 = [13//10 -12//5 9//10 1//5 0 0 0]
n2 = [n1[2] 26//5 -16//5 2//5 0 0 0]
n3 = [n1[3] n2[3] 47//10 -17//5 1 0 0]
n4 = [n1[4] n2[4] n3[4] -29//5 -4 1 0]
n4 = [n1[4] n2[4] n3[4] 29//5 -4 1 0]
n5 = [0 0 1//1 -4 6 -4 1]
n6 = [0 0 0 1//1 -4 6 -4]
n7 = [0 0 0 0 1//1 -4 6]
Expand Down
11 changes: 6 additions & 5 deletions src/SBP_coefficients/Mattsson2014.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
Coefficients of the SBP operators given in
- Mattsson (2014)
Diagonal-norm summation by parts operators for fiite difference approximations
Diagonal-norm summation by parts operators for finite difference approximations
of third and fourth derivatives.
Journal of Computational Physics 274, pp. 432-454.
"""
Expand All @@ -16,7 +16,7 @@ function Base.show(io::IO, source::Mattsson2014)
else
print(io,
"Mattsson (2014) \n",
" Diagonal-norm summation by parts operators for fiite difference approximations\n",
" Diagonal-norm summation by parts operators for finite difference approximations\n",
" of third and fourth derivatives. \n",
" Journal of Computational Physics 274, pp. 432-454.")
end
Expand Down Expand Up @@ -458,7 +458,8 @@ function third_derivative_coefficients(source::Mattsson2014, order::Int, T=Float
upper_coef = SVector(T(-1), T(1//2))
central_coef = T(0)
lower_coef = -upper_coef
left_weights = SVector(T(1//2))
# left_weights = SVector(T(1//2))
left_weights = SVector(T(1//2), T(1), T(1))
right_weights = left_weights
left_boundary_derivatives = (
# first derivative
Expand Down Expand Up @@ -712,15 +713,15 @@ function fourth_derivative_coefficients(source::Mattsson2014, order::Int, T=Floa
DerivativeCoefficientRow{T,1,6}(SVector(T(1//5),
T(2//5),
T(-17//5),
T(-29//5),
T(29//5),
T(-4),
T(1) )),
)
right_boundary = left_boundary
upper_coef = SVector(T(-4), T(1))
central_coef = T(6)
lower_coef = upper_coef
left_weights = SVector(T(1//2))
left_weights = SVector(T(1//2), T(1), T(1), T(1))
right_weights = left_weights
left_boundary_derivatives = (
# first derivative
Expand Down
62 changes: 34 additions & 28 deletions test/SBP_operators_test.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -467,7 +467,7 @@ for source in D_test_list, T in (Float32,Float64)
end

# Accuracy tests of third derivative operators.
for source in D_test_list, T in (Float32,Float64)
@testset "third-derivative operators" for source in D_test_list, T in (Float32,Float64)
xmin = -one(T)
xmax = 2*one(T)
N = 101
Expand Down Expand Up @@ -499,16 +499,12 @@ for source in D_test_list, T in (Float32,Float64)
mul!(res, D, x0)
@test all(i->abs(res[i]) < eps(T), eachindex(res))
mul!(res, D, x1)
@test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
mul!(res, D, x2)
@test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
@test all(i->abs(res[i]) < 50_000*eps(T), eachindex(res))
# only interior
mul!(res, D, x2)
@test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
@test all(i->abs(res[i]) < 160_000*eps(T), inner_indices)
mul!(res, D, x3)
@test all(i->res[i] ≈ 6*x0[i], inner_indices)
mul!(res, D, x4)
@test any(i->!(res[i] ≈ 24*x1[i]), inner_indices)
@test all(i->abs(res[i] - 6*x0[i]) < 340_000*eps(T), inner_indices)
# boundary: first derivative
@test abs(derivative_left( D, x0, Val{1}())) < eps(T)
@test abs(derivative_right(D, x0, Val{1}())) < eps(T)
Expand All @@ -519,8 +515,8 @@ for source in D_test_list, T in (Float32,Float64)
# boundary: second derivative
@test abs(derivative_left( D, x0, Val{2}())) < eps(T)
@test abs(derivative_right(D, x0, Val{2}())) < eps(T)
@test abs(derivative_left( D, x1, Val{2}())) < eps(T)
@test abs(derivative_right(D, x1, Val{2}())) < eps(T)
@test abs(derivative_left( D, x1, Val{2}())) < 600 * eps(T)
@test abs(derivative_right(D, x1, Val{2}())) < 600 * eps(T)
@test derivative_left( D, x2, Val{2}()) ≈ 2
@test derivative_right(D, x2, Val{2}()) ≈ 2
end
Expand Down Expand Up @@ -648,7 +644,7 @@ for source in D_test_list, T in (Float32,Float64)
end

# Accuracy tests of fourth derivative operators.
for source in D_test_list, T in (Float32,Float64)
@testset "fourth-derivative operators" for source in D_test_list, T in (Float32,Float64)
xmin = -one(T)
xmax = 2*one(T)
N = 101
Expand Down Expand Up @@ -678,18 +674,16 @@ for source in D_test_list, T in (Float32,Float64)
@test issymmetric(D) == false
# interior and boundary
mul!(res, D, x0)
@test all(i->abs(res[i]) < eps(T), eachindex(res))
@test all(i->abs(res[i]) < 10 * eps(T) / D.Δx^4, eachindex(res))
mul!(res, D, x1)
@test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
@test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, eachindex(res))
mul!(res, D, x2)
@test all(i->abs(res[i]) < 5000*eps(T), eachindex(res))
@test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, eachindex(res))
# only interior
mul!(res, D, x2)
@test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
@test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, inner_indices)
mul!(res, D, x3)
@test all(i->abs(res[i]) < 5000*eps(T), inner_indices)
mul!(res, D, x4)
@test any(i->!(res[i] ≈ 24*x0[i]), inner_indices)
@test all(i->abs(res[i]) < 100 * eps(T) / D.Δx^4, inner_indices)
# boundary: first derivative
@test abs(derivative_left( D, x0, Val{1}())) < eps(T)
@test abs(derivative_right(D, x0, Val{1}())) < eps(T)
Expand All @@ -700,17 +694,17 @@ for source in D_test_list, T in (Float32,Float64)
# boundary: second derivative
@test abs(derivative_left( D, x0, Val{2}())) < eps(T)
@test abs(derivative_right(D, x0, Val{2}())) < eps(T)
@test abs(derivative_left( D, x1, Val{2}())) < eps(T)
@test abs(derivative_right(D, x1, Val{2}())) < eps(T)
@test derivative_left( D, x2, Val{2}()) ≈ one(T)
@test derivative_right(D, x2, Val{2}()) ≈ one(T)
@test abs(derivative_left( D, x1, Val{2}())) < 600 * eps(T)
@test abs(derivative_right(D, x1, Val{2}())) < 600 * eps(T)
@test derivative_left( D, x2, Val{2}()) ≈ 2 * one(T)
@test derivative_right(D, x2, Val{2}()) ≈ 2 * one(T)
# boundary: third derivative
@test abs(derivative_left( D, x0, Val{3}())) < eps(T)
@test abs(derivative_right(D, x0, Val{3}())) < eps(T)
@test abs(derivative_left( D, x1, Val{3}())) < eps(T)
@test abs(derivative_right(D, x1, Val{3}())) < eps(T)
@test abs(derivative_left( D, x2, Val{3}())) < eps(T)
@test abs(derivative_right(D, x2, Val{3}())) < eps(T)
@test abs(derivative_left( D, x0, Val{3}())) < 10 * eps(T) / D.Δx^3
@test abs(derivative_right(D, x0, Val{3}())) < 10 * eps(T) / D.Δx^3
@test abs(derivative_left( D, x1, Val{3}())) < 10 * eps(T) / D.Δx^3
@test abs(derivative_right(D, x1, Val{3}())) < 10 * eps(T) / D.Δx^3
@test abs(derivative_left( D, x2, Val{3}())) < 10 * eps(T) / D.Δx^3
@test abs(derivative_right(D, x2, Val{3}())) < 10 * eps(T) / D.Δx^3
end

acc_order = 4
Expand Down Expand Up @@ -882,3 +876,15 @@ for T in (Float32, Float64), acc_order in (2,4,6,8)
dest3 = D_safe*u
@test all(i->dest1[i] ≈ dest3[i], eachindex(u))
end


# https://github.com/ranocha/SummationByPartsOperators.jl/issues/208
@testset "issue #208" begin
for acc_order in (2, 4, 6), der_order in 1:4
@test_nowarn derivative_operator(Mattsson2014(),
derivative_order = der_order,
accuracy_order = acc_order,
xmin = -1.0, xmax = 1.0,
N = 20)
end
end

2 comments on commit 108b1d3

@ranocha
Copy link
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JuliaRegistrator register

@JuliaRegistrator
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Registration pull request created: JuliaRegistries/General/89432

After the above pull request is merged, it is recommended that a tag is created on this repository for the registered package version.

This will be done automatically if the Julia TagBot GitHub Action is installed, or can be done manually through the github interface, or via:

git tag -a v0.5.44 -m "<description of version>" 108b1d3f171596d022df397b81b07c095196c050
git push origin v0.5.44

Please sign in to comment.