beuler_step Subroutine

private subroutine beuler_step(state)

Advance the solution by one time step using backward (implicit) Euler.

Solves Q^{n+1} = Q^n + dt · R(Q^{n+1}) by Newton-Raphson iteration:

Q^{0} = Q^n for k = 0 .. n_newton-1: compute R(Q^k) form FD Jacobian J = ∂R/∂Q |_{Q^k} via column-wise perturbation assemble A = I - dt·J (banded structure, kl=ku=neq·coupling_radius) solve A·ΔQ = dt·R(Q^k) - (Q^k - Q^n) Q^{k+1} = Q^k + ΔQ if max|ΔQ| < tol_newton: exit

WARNING: each Newton step performs (neq·n_pt + 1) residual evaluations. This is expensive; use modest grid sizes (n_cell <= 200) when testing. The main benefit is unconditional linear stability, which permits CFL > 1.

Also computes the global L2 residual norm for convergence monitoring.

References: LeVeque, "Finite Volume Methods for Hyperbolic Problems" (2002), Ch. 12.

Arguments

Type IntentOptional Attributes Name
type(solver_state_t), intent(inout) :: state