From d146e80b6a16bc75c82ee90b349d39d188905826 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 17 May 2024 21:06:42 -0400 Subject: [PATCH 1/2] Fix doc builds --- docs/src/showcase/symbolic_analysis.md | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/docs/src/showcase/symbolic_analysis.md b/docs/src/showcase/symbolic_analysis.md index 5c9d2b59d7f..e5463ad97df 100644 --- a/docs/src/showcase/symbolic_analysis.md +++ b/docs/src/showcase/symbolic_analysis.md @@ -57,8 +57,8 @@ tspan = (0, 10.0) pendulum_prob = ODEProblem(pendulum_fun!, u0, tspan, p) traced_sys = modelingtoolkitize(pendulum_prob) pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) -prob = ODAEProblem(pendulum_sys, [], tspan) -sol = solve(prob, Tsit5(), abstol = 1e-8, reltol = 1e-8) +prob = ODEProblem(pendulum_sys, [], tspan) +sol = solve(prob, Rodas5P(), abstol = 1e-8, reltol = 1e-8) plot(sol, vars = states(traced_sys)) ``` @@ -98,7 +98,7 @@ u0 = [1.0, 0, 0, 0, 0]; p = [9.8, 1]; tspan = (0, 10.0); pendulum_prob = ODEProblem(pendulum_fun!, u0, tspan, p) -solve(pendulum_prob, Rodas4()) +solve(pendulum_prob, Rodas5P()) ``` However, one will quickly be greeted with the unfortunate message: @@ -178,7 +178,7 @@ numerical solver. Let's try that out: traced_sys = modelingtoolkitize(pendulum_prob) pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) prob = ODEProblem(pendulum_sys, Pair[], tspan) -sol = solve(prob, Rodas4()) +sol = solve(prob, Rodas5P()) using Plots plot(sol, vars = states(traced_sys)) @@ -198,8 +198,8 @@ be solved via an explicit Runge-Kutta method: ```@example indexred traced_sys = modelingtoolkitize(pendulum_prob) pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) -prob = ODAEProblem(pendulum_sys, Pair[], tspan) -sol = solve(prob, Tsit5(), abstol = 1e-8, reltol = 1e-8) +prob = ODEProblem(pendulum_sys, Pair[], tspan) +sol = solve(prob, Rodas5P(), abstol = 1e-8, reltol = 1e-8) plot(sol, vars = states(traced_sys)) ``` From 62295f7aaee81fca48e5046b97a4cde361af53da Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 17 May 2024 21:09:55 -0400 Subject: [PATCH 2/2] Update symbolic_analysis.md --- docs/src/showcase/symbolic_analysis.md | 25 ++++++++++++------------- 1 file changed, 12 insertions(+), 13 deletions(-) diff --git a/docs/src/showcase/symbolic_analysis.md b/docs/src/showcase/symbolic_analysis.md index e5463ad97df..f21e2cdcbd4 100644 --- a/docs/src/showcase/symbolic_analysis.md +++ b/docs/src/showcase/symbolic_analysis.md @@ -59,7 +59,7 @@ traced_sys = modelingtoolkitize(pendulum_prob) pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) prob = ODEProblem(pendulum_sys, [], tspan) sol = solve(prob, Rodas5P(), abstol = 1e-8, reltol = 1e-8) -plot(sol, vars = states(traced_sys)) +plot(sol, vars = unknowns(traced_sys)) ``` ## Explanation @@ -170,42 +170,41 @@ we can avoid this by using methods like [the Pantelides algorithm](https://ptolemy.berkeley.edu/projects/embedded/eecsx44/lectures/Spring2013/modelica-dae-part-2.pdf) for automatically performing this reduction to index 1. While this requires the ModelingToolkit symbolic form, we use `modelingtoolkitize` to transform -the numerical code into symbolic code, run `dae_index_lowering` lowering, +the numerical code into symbolic code, run `structural_simplify` to +simplify the system and lower the index, then transform back to numerical code with `ODEProblem`, and solve with a numerical solver. Let's try that out: ```@example indexred traced_sys = modelingtoolkitize(pendulum_prob) -pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) +pendulum_sys = structural_simplify(traced_sys) prob = ODEProblem(pendulum_sys, Pair[], tspan) sol = solve(prob, Rodas5P()) using Plots -plot(sol, vars = states(traced_sys)) +plot(sol, vars = unknowns(traced_sys)) ``` -Note that plotting using `states(traced_sys)` is done so that any +Note that plotting using `unknowns(traced_sys)` is done so that any variables which are symbolically eliminated, or any variable reordering done for enhanced parallelism/performance, still show up in the resulting plot and the plot is shown in the same order as the original numerical code. -Note that we can even go a bit further. If we use the `ODAEProblem` -constructor, we can remove the algebraic equations from the states of the -system and fully transform the index-3 DAE into an index-0 ODE, which can -be solved via an explicit Runge-Kutta method: +Note that we can even go a bit further. If we use the `ODEProblem` +constructor, we represent the mass matrix DAE of the index-reduced system, +which can be solved via: ```@example indexred traced_sys = modelingtoolkitize(pendulum_prob) pendulum_sys = structural_simplify(dae_index_lowering(traced_sys)) prob = ODEProblem(pendulum_sys, Pair[], tspan) sol = solve(prob, Rodas5P(), abstol = 1e-8, reltol = 1e-8) -plot(sol, vars = states(traced_sys)) +plot(sol, vars = unknowns(traced_sys)) ``` And there you go: this has transformed the model from being too hard to -solve with implicit DAE solvers, to something that is easily solved with -explicit Runge-Kutta methods for non-stiff equations. +solve with implicit DAE solvers, to something that is easily solved. # Parameter Identifiability in ODE Models @@ -285,7 +284,7 @@ local_id_all = assess_local_identifiability(de, measured_quantities = measured_q # 1 ``` -We can see that all states (except $x_7$) and all parameters are locally identifiable with probability 0.99. +We can see that all unknowns (except $x_7$) and all parameters are locally identifiable with probability 0.99. Let's try to check specific parameters and their combinations