Switch solver from Rodas5P to FBDF and log resolved solver settings

Replace Rodas5P with FBDF throughout run_simulate. Before solving, call
init() to resolve all default kwargs (abstol, reltol, adaptive, dense,
saveat) and write them to the simulation log. Fix saveat splatting bug
where a NamedTuple was spread as positional args instead of kwargs,
causing a convert(Vector{Float64}, ::Float64) error when the ODEProblem
carries a scalar saveat in its own kwargs.

Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>

Author: AnHeuermann

Project.toml

@@ -10,6 +10,7 @@ DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
 Logging = "56ddb016-857b-54e1-b83d-db4d58db5568"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 OMJulia = "0f4fe800-344e-11e9-2949-fb537ad918e1"
+OrdinaryDiffEqBDF = "6ad6398a-0878-4a85-9266-38940aa047c8"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
 Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

src/BaseModelicaLibraryTesting.jl

@@ -4,7 +4,8 @@ import Pkg
 import OMJulia
 import OMJulia: sendExpression
 import BaseModelica
-import DifferentialEquations: solve, Rodas5P, ReturnCode
+import DifferentialEquations: init, solve, ReturnCode
+import OrdinaryDiffEqBDF
 import ModelingToolkit
 import Dates: now
 import Printf: @sprintf

src/simulate.jl

@@ -1,14 +1,15 @@
 # ── Phase 3: ODE simulation with DifferentialEquations / MTK ──────────────────
 
-import DifferentialEquations: solve, Rodas5P, ReturnCode
+import DifferentialEquations: init, solve, ReturnCode
+import OrdinaryDiffEqBDF
 import Logging
 import ModelingToolkit
 import Printf: @sprintf
 
 """
     run_simulate(ode_prob, model_dir, model; cmp_signals, csv_max_size_mb) → (success, time, error, sol)
 
-Solve `ode_prob` with Rodas5P (stiff solver).  On success, also writes the
+Solve `ode_prob` with FBDF (stiff solver).  On success, also writes the
 solution as a CSV file `<Short>_sim.csv` in `model_dir`.
 Writes a `<model>_sim.log` file in `model_dir`.
 Returns `nothing` as the fourth element on failure.
@@ -24,17 +25,18 @@ function run_simulate(ode_prob, model_dir::String,
                       model::String;
                       cmp_signals    ::Vector{String} = String[],
                       csv_max_size_mb::Int            = CSV_MAX_SIZE_MB)::Tuple{Bool,Float64,String,Any}
-    sim_success = false
-    sim_time    = 0.0
-    sim_error   = ""
-    sol         = nothing
+    sim_success           = false
+    sim_time              = 0.0
+    sim_error             = ""
+    sol                   = nothing
+    solver_settings_string = ""
 
     log_file = open(joinpath(model_dir, "$(model)_sim.log"), "w")
     println(log_file, "Model:   $model")
     logger = Logging.SimpleLogger(log_file, Logging.Debug)
     t0 = time()
     try
-        # Rodas5P handles stiff DAE-like systems well.
+        # FBDF handles stiff DAE-like systems and purely algebraic systems well.
         # Redirect all library log output (including Symbolics/MTK warnings)
         # to the log file so they don't clutter stdout.
         sol = Logging.with_logger(logger) do
@@ -46,7 +48,22 @@ function run_simulate(ode_prob, model_dir::String,
             saveat = isempty(ModelingToolkit.unknowns(sys)) ?
                          collect(range(ode_prob.tspan[1], ode_prob.tspan[end]; length = 500)) :
                          Float64[]
-            solve(ode_prob, Rodas5P(); saveat = saveat, dense = true)
+            kwargs = (saveat = saveat, dense = true)
+
+            # Log solver settings
+            initializedSolver = init(ode_prob, OrdinaryDiffEqBDF.FBDF(); kwargs...)
+            solver_settings_string =
+            """
+            OrdinaryDiffEqBDF.FBDF()
+                saveat:   $(let sv = initializedSolver.opts.saveat; isempty(sv) ? "[]" : "$(length(sv)) points in [$(first(sv)), $(last(sv))]" end)
+                abstol:   $(@sprintf("%.2e", initializedSolver.opts.abstol))
+                reltol:   $(@sprintf("%.2e", initializedSolver.opts.reltol))
+                adaptive: $(initializedSolver.opts.adaptive)
+                dense:    $(initializedSolver.opts.dense)
+            """
+
+            # Solve
+            solve(ode_prob, OrdinaryDiffEqBDF.FBDF(); kwargs...)
         end
         sim_time = time() - t0
         if sol.retcode == ReturnCode.Success
@@ -67,7 +84,8 @@ function run_simulate(ode_prob, model_dir::String,
         sim_time  = time() - t0
         sim_error = sprint(showerror, e, catch_backtrace())
     end
-    println(log_file, "Time:    $(round(sim_time; digits=3)) s")
+    println(log_file, "Solver settings: $solver_settings_string")
+    println(log_file, "Time: $(round(sim_time; digits=3)) s")
     println(log_file, "Success: $sim_success")
     isempty(sim_error) || println(log_file, "\n--- Error ---\n$sim_error")
     close(log_file)