From 49e8a04405ec93f63528963b7aac5b58ea6aa89a Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Fri, 24 Jul 2026 07:01:09 +0000 Subject: [PATCH 1/2] docs: update MTK example for ModelingToolkit v11 The `@mtkmodel` DSL was removed in ModelingToolkit v11 in favor of the functional component form (`@component function ... System(eqs, t; name) end`). This updates the quad-tank parameter-fitting example accordingly: - `@mtkmodel QuadtankModel begin ... end` -> `@component function QuadtankModel(; name) ... end` - array variable defaults must match size: `h(t)[1:4] = 0` -> `= zeros(4)` - `ModelingToolkit.defaults(io_sys)[pi]` -> `ModelingToolkit.getdefault(pi)` Verified end-to-end on ModelingToolkit v11.36.0; the estimated parameters still match the true values. Output block, benchmark timing, and the tested-version note in the warning are updated to match. Co-Authored-By: Claude Opus 4.8 (1M context) --- docs/src/examples/modelingtoolkit.md | 29 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 15 deletions(-) diff --git a/docs/src/examples/modelingtoolkit.md b/docs/src/examples/modelingtoolkit.md index 4ecf725a..d08ce28d 100644 --- a/docs/src/examples/modelingtoolkit.md +++ b/docs/src/examples/modelingtoolkit.md @@ -17,7 +17,7 @@ t = ModelingToolkit.t_nounits ssqrt(x) = √(max(x, zero(x)) + 1e-3) # For numerical robustness at x = 0 @register_symbolic ssqrt(x) -@mtkmodel QuadtankModel begin +@component function QuadtankModel(; name) @parameters begin k1 = 1.4 k2 = 1.4 @@ -26,21 +26,20 @@ ssqrt(x) = √(max(x, zero(x)) + 1e-3) # For numerical robustness at x = 0 a = 0.03 γ = 0.25 end - begin - A1 = A2 = A3 = A4 = A - a1 = a3 = a2 = a4 = a - γ1 = γ2 = γ - end @variables begin - h(t)[1:4] = 0 - u(t)[1:2] = 0 + h(t)[1:4] = zeros(4) + u(t)[1:2] = zeros(2) end - @equations begin + A1 = A2 = A3 = A4 = A + a1 = a3 = a2 = a4 = a + γ1 = γ2 = γ + eqs = [ D(h[1]) ~ -a1/A1 * ssqrt(2g*h[1]) + a3/A1*ssqrt(2g*h[3]) + γ1*k1/A1 * u[1] D(h[2]) ~ -a2/A2 * ssqrt(2g*h[2]) + a4/A2*ssqrt(2g*h[4]) + γ2*k2/A2 * u[2] D(h[3]) ~ -a3/A3*ssqrt(2g*h[3]) + (1-γ2)*k2/A3 * u[2] D(h[4]) ~ -a4/A4*ssqrt(2g*h[4]) + (1-γ1)*k1/A4 * u[1] - end + ] + return System(eqs, t; name) end @named mtkmodel = QuadtankModel() @@ -77,7 +76,7 @@ function get_mtk_dynamics(mtkmodel, inputs, outputs, tunable_p) # A wrapper func continuous_dynamics = f_oop # This is ẋ = f(x, u, p, t) inner_discrete_dynamics = SeeToDee.Rk4(continuous_dynamics, Ts::Float64) # x⁺ = f(x, u, p, t) tunable_indices = [findfirst(isequal(pi), p) for pi in tunable_p] # Figure out what indices of the parameter array correspond to our tunable parameters - p0 = [ModelingToolkit.defaults(io_sys)[pi] for pi in p] + p0 = [ModelingToolkit.getdefault(pi) for pi in p] full_p = deepcopy(p0) output_indices = [findfirst(isequal(yi), statevars) for yi in outputs] # Figure out what indices of the state array correspond to our outputs @@ -141,8 +140,8 @@ model = ControlSystemIdentification.nonlinear_pem(d, discrete_dynamics, measurem ``` NonlinearPredictionErrorModel - p: [1.6130151977611773, 1.5995448472434575, 4.887899044534598, 0.20437506116084214] - x0: [2.5590156863624642, 1.674133802252665, 2.890730509103397, 2.114949939609547] + p: [1.611987373192685, 1.599185391084847, 4.877360360086102, 0.2044771473506988] + x0: [1.8477043095652341, 1.1766588459161196, 3.3770895080106467, 2.876924813784857] Ts: 1.0 ny = 2, nu = 2, nx = 4 ``` @@ -161,10 +160,10 @@ using BenchmarkTools @btime ControlSystemIdentification.nonlinear_pem(d, discrete_dynamics, measurement, p_guess, x0_guess, R1, R2, nu) ``` ``` -91.923 ms (876913 allocations: 97.23 MiB) +118.966 ms (876188 allocations: 106.92 MiB) ``` !!! warning - ModelingToolkit is a fast moving target that breaks frequently. The example below was tested with ModelingToolkit v10.0.1, but is not run as part of the build process for this documentation and is not to be considered a supported interface between ControlSystemIdentification and ModelingToolkit. \ No newline at end of file + ModelingToolkit is a fast moving target that breaks frequently. This example was tested with ModelingToolkit v11.36.0, but is not run as part of the build process for this documentation and is not to be considered a supported interface between ControlSystemIdentification and ModelingToolkit. \ No newline at end of file From 54e5a6f4f7006a2fd9d72267a44264bbb012e59f Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Fri, 24 Jul 2026 07:17:55 +0000 Subject: [PATCH 2/2] docs: fix wrapper name in MTK example prose The prose referred to `discrete_dynamics_mtk` but the code defines the wrapper as `discrete_dynamics_wrapper`. Co-Authored-By: Claude Opus 4.8 (1M context) --- docs/src/examples/modelingtoolkit.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/examples/modelingtoolkit.md b/docs/src/examples/modelingtoolkit.md index d08ce28d..d95ab1e3 100644 --- a/docs/src/examples/modelingtoolkit.md +++ b/docs/src/examples/modelingtoolkit.md @@ -58,7 +58,7 @@ We then specify the inputs and outputs of this model, since they are arrays in t ```math \dot x = f(x, u, p, t) ``` -from MTK, we call `ModelingToolkit.generate_control_function`. This example assumes that the system model has an external input, ``u``. If your example does not have this, you may leave this argument empty. This returns two versions of this function, where the second one operates in place (modifying its first argument). This function also returns the state variables chosen, the parameters of the model as well as a simplified system with inputs and outputs. The function returned from `ModelingToolkit.generate_control_function` expects all the parameters of the system to be provided, but we only want to optimize a few of them. We thus wrap this function in `discrete_dynamics_mtk` in order to insert the optimized parameters into a parameter array that contains also the non-optimized parameters. We make use of the function `similar` to ensure that the final parameter array has the correct type (Dual numbers for AD will be used). +from MTK, we call `ModelingToolkit.generate_control_function`. This example assumes that the system model has an external input, ``u``. If your example does not have this, you may leave this argument empty. This returns two versions of this function, where the second one operates in place (modifying its first argument). This function also returns the state variables chosen, the parameters of the model as well as a simplified system with inputs and outputs. The function returned from `ModelingToolkit.generate_control_function` expects all the parameters of the system to be provided, but we only want to optimize a few of them. We thus wrap this function in `discrete_dynamics_wrapper` in order to insert the optimized parameters into a parameter array that contains also the non-optimized parameters. We make use of the function `similar` to ensure that the final parameter array has the correct type (Dual numbers for AD will be used). For good performance, we wrap all the glue code in a function `get_mtk_dynamics` so that we avoid the use of too many global variables.