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

Commit b169937

Browse files
ChrisRackauckasChrisRackauckas
committed
fix concretization for coefficients
1 parent 6f1bc5b commit b169937

File tree

2 files changed

+39
-15
lines changed

2 files changed

+39
-15
lines changed

‎src/derivative_operators/concretization.jl‎

Lines changed: 33 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -3,15 +3,23 @@ function LinearAlgebra.Array(A::DerivativeOperator{T}, N::Int=A.len) where T
33
bl = A.boundary_point_count
44
stl = A.stencil_length
55
bstl = A.boundary_stencil_length
6+
coeff = A.coefficients
67
stl_2 = div(stl,2)
78
for i in 1:A.boundary_point_count
8-
L[i,1:bstl] = A.low_boundary_coefs[i]
9+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
10+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.low_boundary_coefs[i]) : A.low_boundary_coefs[i]
11+
L[i,1:bstl] = cur_coeff * cur_stencil
912
end
1013
for i in bl+1:N-bl
11-
L[i,i+1-stl_2:i+1+stl_2] = A.stencil_coefs
14+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
15+
stencil = eltype(A.stencil_coefs) <: AbstractVector ? A.stencil_coefs[i] : A.stencil_coefs
16+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(stencil) : stencil
17+
L[i,i+1-stl_2:i+1+stl_2] = cur_coeff * cur_stencil
1218
end
1319
for i in N-bl+1:N
14-
L[i,N-bstl+3:N+2] = A.high_boundary_coefs[i-N+bl]
20+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
21+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.high_boundary_coefs[i-N+bl]) : A.high_boundary_coefs[i-N+bl]
22+
L[i,N-bstl+3:N+2] = cur_coeff * cur_stencil
1523
end
1624
return L / A.dx^A.derivative_order
1725
end
@@ -22,14 +30,22 @@ function SparseArrays.SparseMatrixCSC(A::DerivativeOperator{T}, N::Int=A.len) wh
2230
stl = A.stencil_length
2331
stl_2 = div(stl,2)
2432
bstl = A.boundary_stencil_length
33+
coeff = A.coefficients
2534
for i in 1:A.boundary_point_count
26-
L[i,1:bstl] = A.low_boundary_coefs[i]
35+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
36+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.low_boundary_coefs[i]) : A.low_boundary_coefs[i]
37+
L[i,1:bstl] = cur_coeff * cur_stencil
2738
end
2839
for i in bl+1:N-bl
29-
L[i,i+1-stl_2:i+1+stl_2] = A.stencil_coefs
40+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
41+
stencil = eltype(A.stencil_coefs) <: AbstractVector ? A.stencil_coefs[i] : A.stencil_coefs
42+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(stencil) : stencil
43+
L[i,i+1-stl_2:i+1+stl_2] = cur_coeff * cur_stencil
3044
end
3145
for i in N-bl+1:N
32-
L[i,N-bstl+3:N+2] = A.high_boundary_coefs[i-N+bl]
46+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
47+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.high_boundary_coefs[i-N+bl]) : A.high_boundary_coefs[i-N+bl]
48+
L[i,N-bstl+3:N+2] = cur_coeff * cur_stencil
3349
end
3450
return L / A.dx^A.derivative_order
3551
end
@@ -42,16 +58,24 @@ function BandedMatrices.BandedMatrix(A::DerivativeOperator{T}, N::Int=A.len) whe
4258
bl = A.boundary_point_count
4359
stl = A.stencil_length
4460
bstl = A.boundary_stencil_length
61+
coeff = A.coefficients
4562
stl_2 = div(stl,2)
4663
L = BandedMatrix{T}(Zeros(N, N+2), (max(stl-3,0,bstl),max(stl-1,0,bstl)))
4764
for i in 1:A.boundary_point_count
48-
L[i,1:bstl] = A.low_boundary_coefs[i]
65+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
66+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.low_boundary_coefs[i]) : A.low_boundary_coefs[i]
67+
L[i,1:bstl] = cur_coeff * cur_stencil
4968
end
5069
for i in bl+1:N-bl
51-
L[i,i+1-stl_2:i+1+stl_2] = A.stencil_coefs
70+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
71+
stencil = eltype(A.stencil_coefs) <: AbstractVector ? A.stencil_coefs[i] : A.stencil_coefs
72+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(stencil) : stencil
73+
L[i,i+1-stl_2:i+1+stl_2] = cur_coeff * cur_stencil
5274
end
5375
for i in N-bl+1:N
54-
L[i,N-bstl+3:N+2] = A.high_boundary_coefs[i-N+bl]
76+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
77+
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(A.high_boundary_coefs[i-N+bl]) : A.high_boundary_coefs[i-N+bl]
78+
L[i,N-bstl+3:N+2] = cur_coeff * cur_stencil
5579
end
5680
return L / A.dx^A.derivative_order
5781
end

‎src/derivative_operators/convolutions.jl‎

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -17,7 +17,7 @@ function convolve_interior!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::
1717
for i in (1+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count)
1818
xtempi = zero(T)
1919
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i] : stencil
20-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i] : true
20+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
2121
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
2222
for idx in 1:A.stencil_length
2323
xtempi += cur_coeff * cur_stencil[idx] * x[i - mid + idx]
@@ -32,7 +32,7 @@ function convolve_BC_left!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::D
3232
for i in 1 : A.boundary_point_count
3333
xtempi = stencil[i][1]*x[1]
3434
cur_stencil = stencil[i]
35-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i] : true
35+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
3636
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
3737
for idx in 2:A.boundary_stencil_length
3838
xtempi += cur_coeff * cur_stencil[idx] * x[idx]
@@ -47,7 +47,7 @@ function convolve_BC_right!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::
4747
for i in 1 : A.boundary_point_count
4848
xtempi = stencil[i][end]*x[end]
4949
cur_stencil = stencil[i]
50-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i] : true
50+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : true
5151
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
5252
for idx in (A.boundary_stencil_length-1):-1:1
5353
xtempi += cur_coeff * cur_stencil[end-idx] * x[end-idx]
@@ -68,7 +68,7 @@ function convolve_interior!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector,
6868
for i in (1+A.boundary_length) : (length(x_temp)-A.boundary_length)
6969
xtempi = zero(T)
7070
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i-A.boundary_length] : stencil
71-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
71+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
7272
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
7373
@inbounds for idx in 1:A.stencil_length
7474
xtempi += cur_coeff * cur_stencil[idx] * x[i - (mid-idx) + 1]
@@ -83,7 +83,7 @@ function convolve_BC_left!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector,
8383
for i in 1 : A.boundary_length
8484
xtempi = stencil[i][1]*x.l
8585
cur_stencil = stencil[i]
86-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
86+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
8787
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
8888
@inbounds for idx in 2:A.stencil_length
8989
xtempi += cur_coeff * cur_stencil[idx] * x.u[idx-1]
@@ -99,7 +99,7 @@ function convolve_BC_right!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector,
9999
for i in 1 : A.boundary_length
100100
xtempi = stencil[i][end]*x.r
101101
cur_stencil = stencil[i]
102-
cur_coeff = eltype(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
102+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_length] : true
103103
cur_stencil = cur_coeff < 0 && A.use_winding ? reverse(cur_stencil) : cur_stencil
104104
@inbounds for idx in A.stencil_length:-1:2
105105
xtempi += cur_coeff * cur_stencil[end-idx] * x.u[end-idx+1]

0 commit comments

Comments
 (0)