Skip to content

Commit

Permalink
Update derivatives.jl
Browse files Browse the repository at this point in the history
  • Loading branch information
matthieugomez committed Jun 19, 2024
1 parent 364faab commit a1d9662
Showing 1 changed file with 10 additions and 10 deletions.
20 changes: 10 additions & 10 deletions src/derivatives.jl
Original file line number Diff line number Diff line change
Expand Up @@ -57,16 +57,16 @@ Base.IndexStyle(d::SecondDerivative) = IndexLinear()

function Base.getindex(d::SecondDerivative{T}, i::Int) where {T}
x, y, bc = d.x, d.y, d.bc
Δxp = x[min(i, length(x)-1)+1] - x[min(i, length(x)-1)]
Δxm = x[max(i-1, 1) + 1] - x[max(i-1, 1)]
Δx = (Δxm + Δxp) / 2
if i == 1
return convert(T, y[2] / (Δxp * Δx) + (y[1] - bc[1] * Δxm) / (Δxm * Δx) - 2 * y[1] / (Δxp * Δxm))
elseif i == length(x)
return convert(T, (y[end] + bc[end] * Δxp) / (Δxp * Δx) + y[end - 1] / (Δxm * Δx) - 2 * y[end] / (Δxp * Δxm))
else
return convert(T, y[i + 1] / (Δxp * Δx) + y[i - 1] / (Δxm * Δx) - 2 * y[i] / (Δxp * Δxm))
end
Δxp = x[min(i, length(x)-1)+1] - x[min(i, length(x)-1)]
Δxm = x[max(i-1, 1) + 1] - x[max(i-1, 1)]
Δx = (Δxm + Δxp) / 2
if i == 1
return convert(T, y[2] / (Δxp * Δx) + (y[1] - bc[1] * Δxm) / (Δxm * Δx) - 2 * y[1] / (Δxp * Δxm))
elseif i == length(x)
return convert(T, (y[end] + bc[end] * Δxp) / (Δxp * Δx) + y[end - 1] / (Δxm * Δx) - 2 * y[end] / (Δxp * Δxm))
else
return convert(T, y[i + 1] / (Δxp * Δx) + y[i - 1] / (Δxm * Δx) - 2 * y[i] / (Δxp * Δxm))
end
end


Expand Down

0 comments on commit a1d9662

Please sign in to comment.