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

Commit e102e06

Browse files
shivin9shivin9
committed
Old convolve functions from master
1 parent a77e6b7 commit e102e06

File tree

2 files changed

+81
-145
lines changed

2 files changed

+81
-145
lines changed

‎src/derivative_operators/convolutions.jl‎

Lines changed: 80 additions & 144 deletions
Original file line numberDiff line numberDiff line change
@@ -8,227 +8,163 @@ end
88
################################################
99

1010
# Against a standard vector, assume already padded and just apply the stencil
11-
function convolve_interior!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real}
11+
function convolve_interior!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real, N}
1212
@assert length(x_temp)+2 == length(x)
1313
stencil = A.stencil_coefs
1414
coeff = A.coefficients
15-
16-
# Upwind operators have a non-centred stencil
17-
if use_winding(A)
18-
mid = 1 + A.stencil_length%2
19-
else
20-
mid = div(A.stencil_length,2)
21-
end
22-
15+
mid = div(A.stencil_length,2)
2316
for i in (1+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count)
2417
xtempi = zero(T)
2518
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i-A.boundary_point_count] : stencil
2619
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
2720
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
2821
for idx in 1:A.stencil_length
29-
x_idx = use_winding(A) && cur_coeff < 0 ? x[i + mid - idx] : x[i - mid + idx]
30-
xtempi += cur_coeff * cur_stencil[idx] * x_idx
22+
xtempi += cur_coeff * cur_stencil[idx] * x[i - mid + idx]
3123
end
3224
x_temp[i] = xtempi + !overwrite*x_temp[i]
3325
end
3426
end
3527

36-
function convolve_BC_left!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real}
28+
function convolve_BC_left!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real, N}
3729
stencil = A.low_boundary_coefs
3830
coeff = A.coefficients
39-
40-
_bpc = A.boundary_point_count
41-
use_interior_stencil = false
42-
if isempty(stencil)
43-
_bpc = A.boundary_point_count + 1
44-
use_interior_stencil = true
45-
end
46-
47-
for i in 1 : _bpc
48-
if use_interior_stencil == true
49-
cur_stencil = A.stencil_coefs
50-
slen = length(A.stencil_coefs)
51-
else
52-
cur_stencil = stencil[i]
53-
slen = length(cur_stencil)
54-
end
55-
31+
for i in 1 : A.boundary_point_count
32+
cur_stencil = stencil[i]
5633
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
5734
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
58-
xtempi = cur_coeff*cur_stencil[1]*x[1]
59-
for idx in 2:slen
35+
xtempi = cur_coeff*stencil[i][1]*x[1]
36+
for idx in 2:A.boundary_stencil_length
6037
xtempi += cur_coeff * cur_stencil[idx] * x[idx]
6138
end
6239
x_temp[i] = xtempi + !overwrite*x_temp[i]
6340
end
6441
end
6542

66-
function convolve_BC_right!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real}
43+
function convolve_BC_right!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator; overwrite = true) where {T<:Real, N}
6744
stencil = A.high_boundary_coefs
6845
coeff = A.coefficients
69-
N = length(x)
70-
L = A.boundary_stencil_length
71-
72-
_bpc = A.boundary_point_count
73-
use_interior_stencil = false
74-
75-
if isempty(stencil)
76-
_bpc = 1
77-
use_interior_stencil = true
78-
end
79-
80-
for i in 1 : _bpc
81-
if use_interior_stencil == true
82-
cur_stencil = A.stencil_coefs
83-
slen = length(A.stencil_coefs)
84-
L = A.stencil_length
85-
else
86-
cur_stencil = stencil[i]
87-
slen = length(cur_stencil)
88-
end
89-
46+
for i in 1 : A.boundary_point_count
47+
cur_stencil = stencil[i]
9048
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
9149
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
92-
xtempi = zero(T)
93-
for idx in 1:slen
94-
xtempi += cur_coeff * cur_stencil[idx] * x[N-L+idx]
50+
xtempi = cur_coeff*stencil[i][end]*x[end]
51+
for idx in (A.boundary_stencil_length-1):-1:1
52+
xtempi += cur_coeff * cur_stencil[end-idx] * x[end-idx]
9553
end
96-
x_temp[end-_bpc+i] = xtempi + !overwrite*x_temp[end-_bpc+i]
54+
x_temp[end-A.boundary_point_count+i] = xtempi + !overwrite*x_temp[end-A.boundary_point_count+i]
9755
end
9856
end
9957

10058
###########################################
10159

10260
# Against A BC-padded vector, specialize the computation to explicitly use the left, right, and middle parts
103-
function convolve_interior!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,true}; overwrite = true) where {T<:Real,N}
104-
@assert length(x_temp) == length(_x.u)
61+
function convolve_interior!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real, N}
10562
stencil = A.stencil_coefs
10663
coeff = A.coefficients
107-
x_len = length(x_temp) + 2
108-
# Upwind operators have a non-centred stencil
109-
mid = 1 + A.stencil_length%2
11064
x = _x.u
65+
mid = div(A.stencil_length,2) + 1
11166
# Just do the middle parts
112-
for i in (1+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count)
67+
for i in (2+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count)-1
11368
xtempi = zero(T)
11469
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i-A.boundary_point_count] : stencil
11570
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_point_count] : coeff isa Number ? coeff : true
11671
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
11772
@inbounds for idx in 1:A.stencil_length
118-
if i-mid+idx == 1
119-
x_idx = _x.l
120-
elseif i-mid+idx == x_len
121-
x_idx = _x.r
122-
else
123-
x_idx = use_winding(A) && cur_coeff < 0 ? x[i + mid - idx - 1] : x[i - mid + idx - 1]
124-
end
125-
xtempi += cur_coeff * cur_stencil[idx] * x_idx
73+
xtempi += cur_coeff * cur_stencil[idx] * x[(i-1) - (mid-idx) + 1]
12674
end
12775
x_temp[i] = xtempi + !overwrite*x_temp[i]
12876
end
12977
end
13078

131-
function convolve_interior!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real,N}
132-
@assert length(x_temp) == length(_x.u)
133-
stencil = A.stencil_coefs
134-
coeff = A.coefficients
135-
x_len = length(x_temp) + 2
136-
137-
mid = div(A.stencil_length,2)
138-
x = _x.u
139-
140-
# Just do the middle parts
141-
for i in (2+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count-1)
142-
xtempi = zero(T)
143-
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i-A.boundary_point_count] : stencil
144-
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_point_count] : coeff isa Number ? coeff : true
145-
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
146-
@inbounds for idx in 1:A.stencil_length
147-
x_idx = use_winding(A) && cur_coeff < 0 ? x[i + mid - idx - 1] : x[i - mid + idx - 1]
148-
xtempi += cur_coeff * cur_stencil[idx] * x_idx
149-
end
150-
x_temp[i] = xtempi
151-
end
152-
end
153-
154-
155-
function convolve_BC_left!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real,N}
79+
function convolve_BC_left!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real, N}
15680
stencil = A.low_boundary_coefs
15781
coeff = A.coefficients
158-
159-
_bpc = A.boundary_point_count
160-
161-
for i in 1 : _bpc
82+
for i in 1 : A.boundary_point_count
16283
cur_stencil = stencil[i]
163-
slen = length(cur_stencil)
164-
16584
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
166-
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
167-
168-
# need to account for x.l in first interior
16985
xtempi = cur_coeff*cur_stencil[1]*_x.l
170-
@inbounds for idx in 2:slen
86+
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
87+
@inbounds for idx in 2:A.boundary_stencil_length
17188
xtempi += cur_coeff * cur_stencil[idx] * _x.u[idx-1]
17289
end
173-
x_temp[i] = xtempi
90+
x_temp[i] = xtempi + !overwrite*x_temp[i]
17491
end
175-
176-
# explicitely handling the last point which involves the ghost point in its calculations
177-
i = _bpc+1
178-
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
179-
cur_stencil = A.stencil_coefs
180-
slen = length(cur_stencil)
92+
# need to account for x.l in first interior
93+
mid = div(A.stencil_length,2) + 1
94+
x = _x.u
95+
i = 1 + A.boundary_point_count
96+
xtempi = zero(T)
97+
cur_stencil = eltype(A.stencil_coefs) <: AbstractVector ? A.stencil_coefs[i-A.boundary_point_count] : A.stencil_coefs
98+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_point_count] : coeff isa Number ? coeff : true
18199
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
182100
xtempi = cur_coeff*cur_stencil[1]*_x.l
183-
@inbounds for idx in 2:slen
184-
xtempi += cur_coeff * cur_stencil[idx] * _x.u[idx-1]
101+
@inbounds for idx in 2:A.stencil_length
102+
xtempi += cur_coeff * cur_stencil[idx] * x[(i-1) - (mid-idx) + 1]
185103
end
186-
x_temp[i] = xtempi
104+
x_temp[i] = xtempi + !overwrite*x_temp[i]
187105
end
188106

189-
190-
function convolve_BC_right!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real,N}
107+
function convolve_BC_right!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,false}; overwrite = true) where {T<:Real, N}
191108
stencil = A.high_boundary_coefs
192109
coeff = A.coefficients
193-
u_len = length(_x.u)
194-
bpc = A.boundary_point_count
195-
110+
bc_start = length(_x.u) - A.boundary_point_count
196111
# need to account for _x.r in last interior convolution
112+
mid = div(A.stencil_length,2) + 1
197113
x = _x.u
198-
L = A.boundary_stencil_length
199-
_bpc = A.boundary_point_count
200-
201-
# explicitely handling the first point which involves the ghost point in its calculations
202-
i = u_len-_bpc
203-
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
204-
205-
cur_stencil = A.stencil_coefs
114+
i = length(x_temp)-A.boundary_point_count
115+
xtempi = zero(T)
116+
cur_stencil = eltype(A.stencil_coefs) <: AbstractVector ? A.stencil_coefs[i-A.boundary_point_count] : A.stencil_coefs
117+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_point_count] : coeff isa Number ? coeff : true
206118
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
207119
xtempi = cur_coeff*cur_stencil[end]*_x.r
208-
slen = length(cur_stencil)
209-
210-
@inbounds for idx in slen-1:-1:1
211-
xtempi += cur_coeff * cur_stencil[end-idx] * _x.u[end-idx+1]
120+
@inbounds for idx in 1:A.stencil_length-1
121+
xtempi += cur_coeff * cur_stencil[idx] * x[(i-1) - (mid-idx) + 1]
212122
end
213123
x_temp[i] = xtempi + !overwrite*x_temp[i]
214-
215-
for i in u_len-_bpc+1 : u_len
216-
cur_stencil = stencil[i+_bpc-u_len]
217-
slen = length(cur_stencil)
218-
219-
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i] : coeff isa Number ? coeff : true
124+
for i in 1 : A.boundary_point_count
125+
cur_stencil = stencil[i]
126+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[bc_start + i] : coeff isa Number ? coeff : true
220127
xtempi = cur_coeff*cur_stencil[end]*_x.r
221128
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
222-
223-
@inbounds for idx in slen-1:-1:1
129+
@inbounds for idx in A.stencil_length:-1:1
224130
xtempi += cur_coeff * cur_stencil[end-idx] * _x.u[end-idx+1]
225131
end
226-
x_temp[i] = xtempi + !overwrite*x_temp[i]
132+
x_temp[bc_start + i] = xtempi + !overwrite*x_temp[bc_start + i]
133+
end
134+
end
135+
136+
137+
# Against A BC-padded vector, specialize the computation to explicitly use the left, right, and middle parts
138+
function convolve_interior!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,true}; overwrite = true) where {T<:Real,N}
139+
@assert length(x_temp) == length(_x.u)
140+
stencil = A.stencil_coefs
141+
coeff = A.coefficients
142+
x_len = length(x_temp) + 2
143+
# Upwind operators have a non-centred stencil
144+
mid = 1 + A.stencil_length%2
145+
x = _x.u
146+
# Just do the middle parts
147+
for i in (1+A.boundary_point_count) : (length(x_temp)-A.boundary_point_count)
148+
xtempi = zero(T)
149+
cur_stencil = eltype(stencil) <: AbstractVector ? stencil[i-A.boundary_point_count] : stencil
150+
cur_coeff = typeof(coeff) <: AbstractVector ? coeff[i-A.boundary_point_count] : coeff isa Number ? coeff : true
151+
cur_stencil = use_winding(A) && cur_coeff < 0 ? reverse(cur_stencil) : cur_stencil
152+
@inbounds for idx in 1:A.stencil_length
153+
if i-mid+idx == 1
154+
x_idx = _x.l
155+
elseif i-mid+idx == x_len
156+
x_idx = _x.r
157+
else
158+
x_idx = use_winding(A) && cur_coeff < 0 ? x[i + mid - idx - 1] : x[i - mid + idx - 1]
159+
end
160+
xtempi += cur_coeff * cur_stencil[idx] * x_idx
161+
end
162+
x_temp[i] = xtempi
227163
end
228164
end
229165

230166

231-
function convolve_BC_left!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,true}; overwrite = true) where {T<:Real,N}
167+
function convolve_BC_left!(x_temp::AbstractVector{T}, _x::BoundaryPaddedVector, A::DerivativeOperator{T,N,true}) where {T<:Real,N}
232168
stencil = A.low_boundary_coefs
233169
coeff = A.coefficients
234170

@@ -304,7 +240,7 @@ end
304240

305241
#################################################################################
306242

307-
function convolve_interior_add_range!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator, offset::Int) where {T<:Real}
243+
function convolve_interior_add_range!(x_temp::AbstractVector{T}, x::AbstractVector{T}, A::DerivativeOperator, offset::Int) where {T<:Real, N}
308244
@assert length(x_temp)+2 == length(x)
309245
stencil = A.stencil_coefs
310246
coeff = A.coefficients

‎test/runtests.jl‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,6 @@ import Base: isapprox
1515
@time @safetestset "Matrix-Free Operators" begin include("matrixfree.jl") end
1616
@time @safetestset "Convolutions" begin include("convolutions.jl") end
1717
@time @safetestset "Differentiation Dimension" begin include("differentiation_dimension.jl") end
18-
# @time @safetestset "2D and 3D fast multiplication" begin include("2D_3D_fast_multiplication.jl") end
18+
@time @safetestset "2D and 3D fast multiplication" begin include("2D_3D_fast_multiplication.jl") end
1919
@time @safetestset "Higher Dimensional Concretization" begin include("concretization.jl") end
2020
@time @safetestset "Upwind Operator Interface" begin include("upwind_operators_interface.jl") end

0 commit comments

Comments
 (0)