Skip to content

SSP IMEX and Linearized SSP IMEX #34

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 21 commits into
base: ma/dirk
Choose a base branch
from
Open

SSP IMEX and Linearized SSP IMEX #34

Changes from 1 commit
Commits
Show all changes
21 commits
Select commit Hold shift + click to select a range
cd7974a
Add Trixi example
vchuravy Mar 3, 2025
e5016aa
implement ODE interface
vchuravy May 19, 2025
37e9f91
fix TR-BDF2 implementation
vchuravy May 27, 2025
800d79d
remove hallucinated coefficients
vchuravy Jun 6, 2025
772c66f
add SSP Rosenbrock
MarcoArtiano Jun 26, 2025
d0805ca
some cleanup
MarcoArtiano Jun 26, 2025
f8ea589
some more cleanup
MarcoArtiano Jun 26, 2025
b0aa1b5
add struct for tableau
MarcoArtiano Jun 26, 2025
bea7d6a
delete some useless files
MarcoArtiano Jun 26, 2025
6e81d6f
add DIRK
MarcoArtiano Jun 26, 2025
2ee5590
Add three-stage DIRK
MarcoArtiano Jun 26, 2025
8735a38
add TRBDF2 as DIRK
MarcoArtiano Jun 27, 2025
05eb94f
add quick convergence analysis for TRBDF2 as DIRK
MarcoArtiano Jun 27, 2025
db867b6
add linear imex and imex
MarcoArtiano Jul 17, 2025
0a04a09
remove out files
MarcoArtiano Jul 17, 2025
dd53d2f
fix linear implicit
MarcoArtiano Jul 17, 2025
d9e9d75
minor changes
MarcoArtiano Jul 17, 2025
17228e0
fix lhs for linear implicit
MarcoArtiano Jul 18, 2025
530df4a
fix time integrator for linear implicit
MarcoArtiano Jul 18, 2025
cad182b
implement linear ssp imex
MarcoArtiano Jul 18, 2025
edcafb5
add rtol and atol to krylov_solve!
MarcoArtiano Jul 23, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
38 changes: 17 additions & 21 deletions libs/Implicit/src/linear_imex.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,17 +4,23 @@ abstract type RKLIMEX{N} <: SimpleLinearImplicitExplicitAlgorithm{N} end

struct RKLinearImplicitExplicitEuler <: RKLIMEX{1} end

function (::RKLinearImplicitExplicitEuler)(res, uₙ, Δt, f1!, f2!, du, du_tmp, u, p, t, stages, stage, RK, J, lin_du_tmp, lin_du_tmp1)
function (::RKLinearImplicitExplicitEuler)(res, uₙ, Δt, f1!, f2!, du, du_tmp, u, p, t, stages, stage, RK, M, lin_du_tmp, lin_du_tmp1, workspace)
if stage == 1
# Stage 1:
## f2 is the conservative part
## f1 is the parabolic part
mul!(lin_du_tmp, J, uₙ)
mul!(lin_du_tmp1, J, u)
f2!(du, u, p, t + RK.c[stage] * Δt)
f1!(du_tmp, u, p, t + RK.c[stage] * Δt)
return res .= u .- uₙ .- RK.a[stage, stage] * Δt .* (du .+ du_tmp .- lin_du_tmp1 .+ lin_du_tmp)
else
f2!(du, uₙ, p, t + RK.c[stage] * Δt)
f1!(du_tmp, uₙ, p, t + RK.c[stage] * Δt)

res .= uₙ .+ RK.a[stage, stage] * Δt .* (du .+ du_tmp .- lin_du_tmp)
krylov_solve!(workspace, M, copy(res))
Copy link
Member

@vchuravy vchuravy Jul 17, 2025

Choose a reason for hiding this comment

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

Is this meant to be M\res?

Copy link
Collaborator Author

@MarcoArtiano MarcoArtiano Jul 17, 2025

Choose a reason for hiding this comment

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

Yes


f2!(du, workspace.x, p, t + RK.c[stage] * Δt)
f1!(du_tmp, workspace.x, p, t + RK.c[stage] * Δt)

stages[stage] .= du .+ du_tmp
@. u = uₙ + RK.b[1] * Δt * stages[1]
end

Expand Down Expand Up @@ -211,27 +217,17 @@ function step!(integrator::SimpleLinearImplicitExplicit)
end
end


function stage!(integrator, alg::RKLIMEX)
F!(du, u, p) = integrator.f1(du, u, p, integrator.t)
F!(du, u, p) = integrator.f1(du, u, p, integrator.t) ## parabolic
J = JacobianOperator(F!, integrator.du, integrator.u, integrator.p)
for stage in 1:stages(alg)
F! = nonlinear_problem(alg, integrator.f2)
# TODO: Pass in `stages[1:(stage-1)]` or full tuple?
_, stats = Ariadne.newton_krylov!(
F!, integrator.u_tmp, (integrator.u, integrator.dt, integrator.f1, integrator.du, integrator.du_tmp, integrator.p, integrator.t, integrator.stages, stage, integrator.RK, J, integrator.lin_du_tmp, integrator.lin_du_tmp1), integrator.res;
verbose=integrator.opts.verbose, krylov_kwargs=integrator.opts.krylov_kwargs,
algo=integrator.opts.algo, tol_abs=6.0e-6,
)
@assert stats.solved
M = MOperator(J, integrator.dt)
kc = KrylovConstructor(integrator.res)
workspace = krylov_workspace(:gmres, kc)
for stage in 1:stages(alg)
# Store the solution for each stage in stages
## For a split Problem we need to compute rhs_conservative and rhs_parabolic
integrator.f2(integrator.du, integrator.u_tmp, integrator.p, integrator.t + integrator.RK.c[stage] * integrator.dt)
integrator.stages[stage] .= integrator.du
integrator.f1(integrator.du, integrator.u_tmp, integrator.p, integrator.t + integrator.RK.c[stage] * integrator.dt)
integrator.stages[stage] .+= integrator.du
if stage == stages(alg)
alg(integrator.res, integrator.u, integrator.dt, integrator.f1, integrator.f2, integrator.du, integrator.du_tmp, integrator.u_tmp, integrator.p, integrator.t, integrator.stages, stage + 1, integrator.RK, J, integrator.lin_du_tmp, integrator.lin_du_tmp1)
alg(integrator.res, integrator.u, integrator.dt, integrator.f1, integrator.f2, integrator.du, integrator.du_tmp, integrator.u_tmp, integrator.p, integrator.t, integrator.stages, stage + 1, integrator.RK, M, integrator.lin_du_tmp, integrator.lin_du_tmp1, workspace)
end

end
Expand Down
Loading