Skip to content

Commit

Permalink
Merge pull request #2293 from oscardssmith/os/fix-oop-BDF
Browse files Browse the repository at this point in the history
don't broadcast the left side for out of place
  • Loading branch information
oscardssmith authored Jul 25, 2024
2 parents 6774e5d + 367c34b commit 17bd407
Showing 1 changed file with 11 additions and 15 deletions.
26 changes: 11 additions & 15 deletions lib/OrdinaryDiffEqBDF/src/bdf_perform_step.jl
Original file line number Diff line number Diff line change
Expand Up @@ -47,12 +47,11 @@ end
mass_matrix = f.mass_matrix

if mass_matrix === I
nlsolver.tmp = @.. broadcast=false ((dtₙ * β₁) * fₙ₋₁ +
(α₁ * uₙ₋₁ + α₂ * uₙ₋₂))/(dtₙ *
β₀)
nlsolver.tmp = @.. ((dtₙ * β₁) * fₙ₋₁ +
(α₁ * uₙ₋₁ + α₂ * uₙ₋₂))/(dtₙ * β₀)
else
_tmp = mass_matrix * @.. broadcast=false (α₁ * uₙ₋₁+α₂ * uₙ₋₂)
nlsolver.tmp = @.. broadcast=false ((dtₙ * β₁) * fₙ₋₁ + _tmp)/(dtₙ * β₀)
_tmp = mass_matrix * @.. (α₁ * uₙ₋₁+α₂ * uₙ₋₂)
nlsolver.tmp = @.. ((dtₙ * β₁) * fₙ₋₁ + _tmp)/(dtₙ * β₀)
end
nlsolver.γ = β₀
nlsolver.α = α₀
Expand Down Expand Up @@ -772,17 +771,17 @@ function perform_step!(integrator, cache::QNDFConstantCache{max_order},
u₀ = reshape(sum(D[:, 1:k], dims = 2) .+ uprev, size(u))
ϕ = zero(u)
for i in 1:k
@.. broadcast=false ϕ+=γₖ[i] * D[:, i]
ϕ = @.. ϕ + γₖ[i] * D[:, i]
end
end
markfirststage!(nlsolver)
nlsolver.z = u₀
mass_matrix = f.mass_matrix

if mass_matrix === I
nlsolver.tmp = @.. broadcast=false (u₀ / β₀ - ϕ)/dt
nlsolver.tmp = @.. (u₀ / β₀ - ϕ)/dt
else
nlsolver.tmp = mass_matrix * @.. broadcast=false (u₀ / β₀ - ϕ)/dt
nlsolver.tmp = mass_matrix * @.. (u₀ / β₀ - ϕ)/dt
end

nlsolver.γ = β₀
Expand Down Expand Up @@ -1110,18 +1109,16 @@ function perform_step!(integrator, cache::FBDFConstantCache{max_order},
end
else
for i in 1:(k - 1)
@.. broadcast=false @views u_corrector[:, i] = $calc_Lagrange_interp(
@.. @views u_corrector[:, i] = $calc_Lagrange_interp(
k, weights,
equi_ts[i],
ts,
u_history,
u_corrector[:,
i])
u_corrector[:, i])
end
tmp = -uprev * bdf_coeffs[k, 2]
vc = _vec(tmp)
for i in 1:(k - 1)
@.. broadcast=false @views vc -= u_corrector[:, i] * bdf_coeffs[k, i + 2]
@views tmp = @.. tmp - u_corrector[:, i] * bdf_coeffs[k, i + 2]
end
end

Expand Down Expand Up @@ -1177,9 +1174,8 @@ function perform_step!(integrator, cache::FBDFConstantCache{max_order},
end
terk *= abs(dt^(k))
else
vc = _vec(terk)
for i in 2:(k + 1)
@.. broadcast=false @views vc += fd_weights[i, k + 1] * u_history[:, i - 1]
@views terk = @.. terk + fd_weights[i, k + 1] * u_history[:, i - 1]
end
terk *= abs(dt^(k))
end
Expand Down

0 comments on commit 17bd407

Please sign in to comment.