diff --git a/docs/make.jl b/docs/make.jl
index 47fc79f..bc0f9bb 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -48,7 +48,7 @@ repo_url = "github.com/control-toolbox/Tutorials.jl"
 Draft = false
 ```
 =#
-draft = false  # Draft mode: if true, @example blocks in markdown are not executed
+draft = true  # Draft mode: if true, @example blocks in markdown are not executed
 
 # ═══════════════════════════════════════════════════════════════════════════════
 # Build documentation
diff --git a/docs/src/assets/Manifest.toml b/docs/src/assets/Manifest.toml
index c65e39b..ced7cf3 100644
--- a/docs/src/assets/Manifest.toml
+++ b/docs/src/assets/Manifest.toml
@@ -2,7 +2,7 @@
 
 julia_version = "1.12.5"
 manifest_format = "2.0"
-project_hash = "a947cf68efaf3e1dbdfceb29223045209d5bb275"
+project_hash = "1d73fc8b7306dbb42eb2a1fdd2ae0c8be50a6d19"
 
 [[deps.ADNLPModels]]
 deps = ["ADTypes", "ForwardDiff", "LinearAlgebra", "NLPModels", "Requires", "ReverseDiff", "SparseArrays", "SparseConnectivityTracer", "SparseMatrixColorings"]
@@ -38,6 +38,12 @@ git-tree-sha1 = "6252039f98492252f9e47c312c8ffda0e3b9e78d"
 uuid = "ae81ac8f-d209-56e5-92de-9978fef736f9"
 version = "0.1.3+0"
 
+[[deps.AbstractPlutoDingetjes]]
+deps = ["Pkg"]
+git-tree-sha1 = "6e1d2a35f2f90a4bc7c2ed98079b2ba09c35b83a"
+uuid = "6e696c72-6542-2067-7265-42206c756150"
+version = "1.3.2"
+
 [[deps.AbstractTrees]]
 git-tree-sha1 = "2d9c9a55f9c93e8887ad391fbae72f8ef55e1177"
 uuid = "1520ce14-60c1-5f80-bbc7-55ef81b5835c"
@@ -88,6 +94,12 @@ version = "1.1.3"
 uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"
 version = "1.1.2"
 
+[[deps.ArnoldiMethod]]
+deps = ["LinearAlgebra", "Random", "StaticArrays"]
+git-tree-sha1 = "d57bd3762d308bded22c3b82d033bff85f6195c6"
+uuid = "ec485272-7323-5ecc-a04f-4719b315124d"
+version = "0.4.0"
+
 [[deps.ArrayInterface]]
 deps = ["Adapt", "LinearAlgebra"]
 git-tree-sha1 = "78b3a7a536b4b0a747a0f296ea77091ca0a9f9a3"
@@ -151,10 +163,15 @@ uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"
 version = "1.11.0"
 
 [[deps.BenchmarkTools]]
-deps = ["Compat", "JSON", "Logging", "Printf", "Profile", "Statistics", "UUIDs"]
-git-tree-sha1 = "6876e30dc02dc69f0613cb6ece242144f2ca9e56"
+deps = ["Compat", "JSON", "Logging", "PrecompileTools", "Printf", "Profile", "Statistics", "UUIDs"]
+git-tree-sha1 = "9670d3febc2b6da60a0ae57846ba74670290653f"
 uuid = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
-version = "1.7.0"
+version = "1.8.0"
+
+[[deps.Bijections]]
+git-tree-sha1 = "a2d308fcd4c2fb90e943cf9cd2fbfa9c32b69733"
+uuid = "e2ed5e7c-b2de-5872-ae92-c73ca462fb04"
+version = "0.2.2"
 
 [[deps.BitFlags]]
 git-tree-sha1 = "0691e34b3bb8be9307330f88d1a3c3f25466c24d"
@@ -248,9 +265,9 @@ version = "0.8.13"
 
 [[deps.CTSolvers]]
 deps = ["ADNLPModels", "CTBase", "CTModels", "CommonSolve", "DocStringExtensions", "ExaModels", "KernelAbstractions", "NLPModels", "SolverCore"]
-git-tree-sha1 = "37a9c44ca8fcc6948854fb796b75e545a4efbb9e"
+git-tree-sha1 = "22c283a24bd1b51cf2795074a6ec3fe0bff78adb"
 uuid = "d3e8d392-8e4b-4d9b-8e92-d7d4e3650ef6"
-version = "0.4.14"
+version = "0.4.15"
 
     [deps.CTSolvers.extensions]
     CTSolversCUDA = "CUDA"
@@ -334,6 +351,11 @@ git-tree-sha1 = "37ea44092930b1811e666c3bc38065d7d87fcc74"
 uuid = "5ae59095-9a9b-59fe-a467-6f913c188581"
 version = "0.13.1"
 
+[[deps.Combinatorics]]
+git-tree-sha1 = "08c8b6831dc00bfea825826be0bc8336fc369860"
+uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa"
+version = "1.0.2"
+
 [[deps.CommonSolve]]
 git-tree-sha1 = "78ea4ddbcf9c241827e7035c3a03e2e456711470"
 uuid = "38540f10-b2f7-11e9-35d8-d573e4eb0ff2"
@@ -365,6 +387,11 @@ deps = ["Artifacts", "Libdl"]
 uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae"
 version = "1.3.0+1"
 
+[[deps.CompositeTypes]]
+git-tree-sha1 = "bce26c3dab336582805503bed209faab1c279768"
+uuid = "b152e2b5-7a66-4b01-a709-34e65c35f657"
+version = "0.1.4"
+
 [[deps.CompositionsBase]]
 git-tree-sha1 = "802bb88cd69dfd1509f6670416bd4434015693ad"
 uuid = "a33af91c-f02d-484b-be07-31d278c5ca2b"
@@ -389,17 +416,13 @@ version = "2.5.1"
 git-tree-sha1 = "b4b092499347b18a015186eae3042f72267106cb"
 uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9"
 version = "1.6.0"
+weakdeps = ["IntervalSets", "LinearAlgebra", "StaticArrays"]
 
     [deps.ConstructionBase.extensions]
     ConstructionBaseIntervalSetsExt = "IntervalSets"
     ConstructionBaseLinearAlgebraExt = "LinearAlgebra"
     ConstructionBaseStaticArraysExt = "StaticArrays"
 
-    [deps.ConstructionBase.weakdeps]
-    IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953"
-    LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-    StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
-
 [[deps.Contour]]
 git-tree-sha1 = "439e35b0b36e2e5881738abc8857bd92ad6ff9a8"
 uuid = "d38c429a-6771-53c6-b99e-75d170b6e991"
@@ -455,6 +478,12 @@ git-tree-sha1 = "9e2f36d3c96a820c678f2f1f1782582fcf685bae"
 uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab"
 version = "1.9.1"
 
+[[deps.Dictionaries]]
+deps = ["Indexing", "Random", "Serialization"]
+git-tree-sha1 = "a55766a9c8f66cf19ffcdbdb1444e249bb4ace33"
+uuid = "85a47980-9c8c-11e8-2b9f-f7ca1fa99fb4"
+version = "0.4.6"
+
 [[deps.DiffEqBase]]
 deps = ["ArrayInterface", "BracketingNonlinearSolve", "ConcreteStructs", "DocStringExtensions", "FastBroadcast", "FastClosures", "FastPower", "FunctionWrappers", "FunctionWrappersWrappers", "LinearAlgebra", "Logging", "Markdown", "MuladdMacro", "PrecompileTools", "Printf", "RecursiveArrayTools", "Reexport", "SciMLBase", "SciMLOperators", "SciMLStructures", "Setfield", "Static", "StaticArraysCore", "SymbolicIndexingInterface", "TruncatedStacktraces"]
 git-tree-sha1 = "87e2ad6d4ae98505218e2f97cafcfa296dc97d37"
@@ -588,11 +617,31 @@ git-tree-sha1 = "d8a8cb2d5b0181fbbd41861016b221b0202c9bc5"
 uuid = "d12716ef-a0f6-4df4-a9f1-a5a34e75c656"
 version = "1.1.0"
 
+[[deps.DomainSets]]
+deps = ["CompositeTypes", "FunctionMaps", "IntervalSets", "LinearAlgebra", "StaticArrays"]
+git-tree-sha1 = "1af39efbaf76fb648432b5efaac0d73af6760407"
+uuid = "5b8099bc-c8ec-5219-889f-1d9e522a28bf"
+version = "0.8.0"
+
+    [deps.DomainSets.extensions]
+    DomainSetsMakieExt = "Makie"
+    DomainSetsRandomExt = "Random"
+
+    [deps.DomainSets.weakdeps]
+    Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
+    Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+
 [[deps.Downloads]]
 deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"]
 uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"
 version = "1.7.0"
 
+[[deps.DynamicPolynomials]]
+deps = ["LinearAlgebra", "MultivariatePolynomials", "MutableArithmetics", "Reexport", "StarAlgebras", "Test"]
+git-tree-sha1 = "5bfabc3827dfdd164359bad6800c115a81280c00"
+uuid = "7c1d4256-1411-5781-91ec-d7bc3513ac07"
+version = "0.6.6"
+
 [[deps.EnumX]]
 git-tree-sha1 = "c49898e8438c828577f04b92fc9368c388ac783c"
 uuid = "4e289a0a-7415-4d19-859d-a7e5c4648b56"
@@ -815,6 +864,12 @@ git-tree-sha1 = "7a214fdac5ed5f59a22c2d9a885a16da1c74bbc7"
 uuid = "559328eb-81f9-559d-9380-de523a88c83c"
 version = "1.0.17+0"
 
+[[deps.FunctionMaps]]
+deps = ["CompositeTypes", "LinearAlgebra", "StaticArrays"]
+git-tree-sha1 = "31bd99a57edf98990d1c21486032963955450e8d"
+uuid = "a85aefff-f8ca-4649-a888-c8e5398bc76c"
+version = "0.1.2"
+
 [[deps.FunctionWrappers]]
 git-tree-sha1 = "d62485945ce5ae9c0c48f124a84998d755bae00e"
 uuid = "069b7b12-0de2-55c6-9aab-29f3d0a68a2e"
@@ -918,6 +973,19 @@ git-tree-sha1 = "8a6dbda1fd736d60cc477d99f2e7a042acfa46e8"
 uuid = "3b182d85-2403-5c21-9c21-1e1f0cc25472"
 version = "1.3.15+0"
 
+[[deps.Graphs]]
+deps = ["ArnoldiMethod", "DataStructures", "Inflate", "LinearAlgebra", "Random", "SimpleTraits", "SparseArrays", "Statistics"]
+git-tree-sha1 = "7eb45fe833a5b7c51cf6d89c5a841d5967e44be3"
+uuid = "86223c79-3864-5bf0-83f7-82e725a168b6"
+version = "1.14.0"
+
+    [deps.Graphs.extensions]
+    GraphsSharedArraysExt = "SharedArrays"
+
+    [deps.Graphs.weakdeps]
+    Distributed = "8ba89e20-285c-5b6f-9357-94700520ee1b"
+    SharedArrays = "1a1011a3-84de-559e-8e89-a11a2f7dc383"
+
 [[deps.Grisu]]
 git-tree-sha1 = "53bb909d1151e57e2484c3d1b53e19552b887fb2"
 uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe"
@@ -952,6 +1020,16 @@ git-tree-sha1 = "debdd00ffef04665ccbb3e150747a77560e8fad1"
 uuid = "615f187c-cbe4-4ef1-ba3b-2fcf58d6d173"
 version = "0.1.1"
 
+[[deps.Indexing]]
+git-tree-sha1 = "ce1566720fd6b19ff3411404d4b977acd4814f9f"
+uuid = "313cdc1a-70c2-5d6a-ae34-0150d3930a38"
+version = "1.1.1"
+
+[[deps.Inflate]]
+git-tree-sha1 = "d1b1b796e47d94588b3757fe84fbf65a5ec4a80d"
+uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9"
+version = "0.1.5"
+
 [[deps.InlineStrings]]
 git-tree-sha1 = "8f3d257792a522b4601c24a577954b0a8cd7334d"
 uuid = "842dd82b-1e85-43dc-bf29-5d0ee9dffc48"
@@ -965,6 +1043,11 @@ version = "1.4.5"
     ArrowTypes = "31f734f8-188a-4ce0-8406-c8a06bd891cd"
     Parsers = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
 
+[[deps.IntegerMathUtils]]
+git-tree-sha1 = "4c1acff2dc6b6967e7e750633c50bc3b8d83e617"
+uuid = "18e54dd8-cb9d-406c-a71d-865a43cbb235"
+version = "0.1.3"
+
 [[deps.IntelOpenMP_jll]]
 deps = ["Artifacts", "JLLWrappers", "LazyArtifacts", "Libdl"]
 git-tree-sha1 = "ec1debd61c300961f98064cfb21287613ad7f303"
@@ -976,6 +1059,17 @@ deps = ["Markdown"]
 uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
 version = "1.11.0"
 
+[[deps.IntervalSets]]
+git-tree-sha1 = "79d6bd28c8d9bccc2229784f1bd637689b256377"
+uuid = "8197267c-284f-5f27-9208-e0e47529a953"
+version = "0.7.14"
+weakdeps = ["Random", "RecipesBase", "Statistics"]
+
+    [deps.IntervalSets.extensions]
+    IntervalSetsRandomExt = "Random"
+    IntervalSetsRecipesBaseExt = "RecipesBase"
+    IntervalSetsStatisticsExt = "Statistics"
+
 [[deps.InverseFunctions]]
 git-tree-sha1 = "a779299d77cd080bf77b97535acecd73e1c5e5cb"
 uuid = "3587e190-3f89-42d0-90ee-14403ec27112"
@@ -1218,9 +1312,9 @@ weakdeps = ["LineSearches"]
 
 [[deps.LineSearches]]
 deps = ["LinearAlgebra", "NLSolversBase", "NaNMath", "Printf"]
-git-tree-sha1 = "738bdcacfef25b3a9e4a39c28613717a6b23751e"
+git-tree-sha1 = "a666999510c794fe1d9dfa629ef33366f11103aa"
 uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255"
-version = "7.6.0"
+version = "7.6.1"
 
 [[deps.LinearAlgebra]]
 deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"]
@@ -1460,6 +1554,22 @@ git-tree-sha1 = "cac9cc5499c25554cba55cd3c30543cff5ca4fab"
 uuid = "46d2c3a1-f734-5fdb-9937-b9b9aeba4221"
 version = "0.2.4"
 
+[[deps.MultivariatePolynomials]]
+deps = ["DataStructures", "LinearAlgebra", "MutableArithmetics", "StarAlgebras"]
+git-tree-sha1 = "b45f1ed8448ea20885cb4c5029c2a462fe2682bf"
+uuid = "102ac46a-7ee4-5c85-9060-abc95bfdeaa3"
+version = "0.5.17"
+weakdeps = ["ChainRulesCore"]
+
+    [deps.MultivariatePolynomials.extensions]
+    MultivariatePolynomialsChainRulesCoreExt = "ChainRulesCore"
+
+[[deps.MutableArithmetics]]
+deps = ["LinearAlgebra", "SparseArrays", "Test"]
+git-tree-sha1 = "7c25249fc13a070f5ba433c50e21e22bb33c6fb0"
+uuid = "d8a4904e-b15c-11e9-3269-09a3773c0cb0"
+version = "1.7.1"
+
 [[deps.NLPModels]]
 deps = ["FastClosures", "LinearAlgebra", "LinearOperators", "Printf", "SparseArrays"]
 git-tree-sha1 = "5919b1be82ca6ebfc9754cb22ddfb54f56f1d9ac"
@@ -1969,6 +2079,12 @@ version = "3.3.2"
     [deps.PrettyTables.weakdeps]
     Typstry = "f0ed7684-a786-439e-b1e3-3b82803b501e"
 
+[[deps.Primes]]
+deps = ["IntegerMathUtils"]
+git-tree-sha1 = "25cdd1d20cd005b52fc12cb6be3f75faaf59bb9b"
+uuid = "27ebfcd6-29c5-5fa9-bf4b-fb8fc14df3ae"
+version = "0.5.7"
+
 [[deps.Printf]]
 deps = ["Unicode"]
 uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7"
@@ -2024,6 +2140,11 @@ deps = ["SHA"]
 uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 version = "1.11.0"
 
+[[deps.ReadOnlyArrays]]
+git-tree-sha1 = "e6f7ddf48cf141cb312b078ca21cb2d29d0dc11d"
+uuid = "988b38a3-91fc-5605-94a2-ee2116b3bd83"
+version = "0.2.0"
+
 [[deps.RecipesBase]]
 deps = ["PrecompileTools"]
 git-tree-sha1 = "5c3d09cc4f31f5fc6af001c250bf1278733100ff"
@@ -2256,6 +2377,12 @@ version = "2.11.1"
     ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267"
     Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c"
 
+[[deps.SimpleTraits]]
+deps = ["InteractiveUtils", "MacroTools"]
+git-tree-sha1 = "be8eeac05ec97d379347584fa9fe2f5f76795bcb"
+uuid = "699a6c99-e7fa-54fc-8d76-47d257e15c1d"
+version = "0.9.5"
+
 [[deps.Sockets]]
 uuid = "6462fe0b-24de-5631-8697-dd941f90decc"
 version = "1.11.0"
@@ -2336,6 +2463,12 @@ git-tree-sha1 = "4f96c596b8c8258cc7d3b19797854d368f243ddc"
 uuid = "860ef19b-820b-49d6-a774-d7a799459cd3"
 version = "1.0.4"
 
+[[deps.StarAlgebras]]
+deps = ["LinearAlgebra", "MutableArithmetics", "SparseArrays"]
+git-tree-sha1 = "235b1f9d287bbf34083b3d0829343a7942c0ad1c"
+uuid = "0c0c59c1-dc5f-42e9-9a8b-b5dc384a6cd1"
+version = "0.3.0"
+
 [[deps.Static]]
 deps = ["CommonWorldInvalidations", "IfElse", "PrecompileTools", "SciMLPublic"]
 git-tree-sha1 = "49440414711eddc7227724ae6e570c7d5559a086"
@@ -2408,9 +2541,9 @@ version = "0.4.4"
 
 [[deps.StructUtils]]
 deps = ["Dates", "UUIDs"]
-git-tree-sha1 = "fa95b3b097bcef5845c142ea2e085f1b2591e92c"
+git-tree-sha1 = "aab80fbf866600f3299dd7f6656d80e7be177cfe"
 uuid = "ec057cc2-7a8d-4b58-b3b3-92acb9f63b42"
-version = "2.7.1"
+version = "2.7.2"
 
     [deps.StructUtils.extensions]
     StructUtilsMeasurementsExt = ["Measurements"]
@@ -2451,6 +2584,63 @@ weakdeps = ["PrettyTables"]
     [deps.SymbolicIndexingInterface.extensions]
     SymbolicIndexingInterfacePrettyTablesExt = "PrettyTables"
 
+[[deps.SymbolicLimits]]
+deps = ["SymbolicUtils", "TermInterface"]
+git-tree-sha1 = "5085671d2cba1eb02136a3d6661c583e801984c1"
+uuid = "19f23fe9-fdab-4a78-91af-e7b7767979c3"
+version = "1.1.0"
+
+[[deps.SymbolicUtils]]
+deps = ["AbstractTrees", "ArrayInterface", "Combinatorics", "ConstructionBase", "DataStructures", "Dictionaries", "DocStringExtensions", "DynamicPolynomials", "EnumX", "ExproniconLite", "Graphs", "LinearAlgebra", "MacroTools", "Moshi", "MultivariatePolynomials", "MutableArithmetics", "NaNMath", "PrecompileTools", "ReadOnlyArrays", "Setfield", "SparseArrays", "SpecialFunctions", "StaticArraysCore", "SymbolicIndexingInterface", "TaskLocalValues", "TermInterface", "WeakCacheSets"]
+git-tree-sha1 = "552eaeac659f5802d45c4d936272a712ca13a969"
+uuid = "d1185830-fcd6-423d-90d6-eec64667417b"
+version = "4.24.1"
+
+    [deps.SymbolicUtils.extensions]
+    SymbolicUtilsChainRulesCoreExt = "ChainRulesCore"
+    SymbolicUtilsDistributionsExt = "Distributions"
+    SymbolicUtilsLabelledArraysExt = "LabelledArrays"
+    SymbolicUtilsReverseDiffExt = "ReverseDiff"
+
+    [deps.SymbolicUtils.weakdeps]
+    ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
+    Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
+    LabelledArrays = "2ee39098-c373-598a-b85f-a56591580800"
+    ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267"
+
+[[deps.Symbolics]]
+deps = ["ADTypes", "AbstractPlutoDingetjes", "ArrayInterface", "Bijections", "CommonWorldInvalidations", "ConstructionBase", "DataStructures", "DiffRules", "DocStringExtensions", "DomainSets", "DynamicPolynomials", "Libdl", "LinearAlgebra", "LogExpFunctions", "MacroTools", "Markdown", "Moshi", "MultivariatePolynomials", "MutableArithmetics", "NaNMath", "PrecompileTools", "Preferences", "Primes", "RecipesBase", "Reexport", "RuntimeGeneratedFunctions", "SciMLPublic", "Setfield", "SparseArrays", "SpecialFunctions", "StaticArraysCore", "SymbolicIndexingInterface", "SymbolicLimits", "SymbolicUtils", "TermInterface"]
+git-tree-sha1 = "5c5db34512d5349b5e68d75cd707a3190fb26c81"
+uuid = "0c5d862f-8b57-4792-8d23-62f2024744c7"
+version = "7.18.1"
+
+    [deps.Symbolics.extensions]
+    SymbolicsD3TreesExt = "D3Trees"
+    SymbolicsDistributionsExt = "Distributions"
+    SymbolicsForwardDiffExt = "ForwardDiff"
+    SymbolicsGroebnerExt = "Groebner"
+    SymbolicsHypergeometricFunctionsExt = "HypergeometricFunctions"
+    SymbolicsLatexifyExt = ["Latexify", "LaTeXStrings"]
+    SymbolicsLuxExt = "Lux"
+    SymbolicsNemoExt = "Nemo"
+    SymbolicsPreallocationToolsExt = ["PreallocationTools", "ForwardDiff"]
+    SymbolicsSymPyExt = "SymPy"
+    SymbolicsSymPyPythonCallExt = "SymPyPythonCall"
+
+    [deps.Symbolics.weakdeps]
+    D3Trees = "e3df1716-f71e-5df9-9e2d-98e193103c45"
+    Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
+    ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
+    Groebner = "0b43b601-686d-58a3-8a1c-6623616c7cd4"
+    HypergeometricFunctions = "34004b35-14d8-5ef3-9330-4cdb6864b03a"
+    LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
+    Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
+    Lux = "b2108857-7c20-44ae-9111-449ecde12c47"
+    Nemo = "2edaba10-b0f1-5616-af89-8c11ac63239a"
+    PreallocationTools = "d236fae5-4411-538c-8e31-a6e3d9e00b46"
+    SymPy = "24249f21-da20-56a4-8eb1-6a02cf4ae2e6"
+    SymPyPythonCall = "bc8888f7-b21e-4b7c-a06a-5d9c9496438c"
+
 [[deps.TOML]]
 deps = ["Dates"]
 uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
@@ -2473,12 +2663,22 @@ deps = ["ArgTools", "SHA"]
 uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"
 version = "1.10.0"
 
+[[deps.TaskLocalValues]]
+git-tree-sha1 = "67e469338d9ce74fc578f7db1736a74d93a49eb8"
+uuid = "ed4db957-447d-4319-bfb6-7fa9ae7ecf34"
+version = "0.1.3"
+
 [[deps.TensorCore]]
 deps = ["LinearAlgebra"]
 git-tree-sha1 = "1feb45f88d133a655e001435632f019a9a1bcdb6"
 uuid = "62fd8b95-f654-4bbd-a8a5-9c27f68ccd50"
 version = "0.1.1"
 
+[[deps.TermInterface]]
+git-tree-sha1 = "d673e0aca9e46a2f63720201f55cc7b3e7169b16"
+uuid = "8ea1fca8-c5ef-4a55-8b96-4e9afe9c9a3c"
+version = "2.0.0"
+
 [[deps.Test]]
 deps = ["InteractiveUtils", "Logging", "Random", "Serialization"]
 uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
@@ -2566,6 +2766,11 @@ git-tree-sha1 = "96478df35bbc2f3e1e791bc7a3d0eeee559e60e9"
 uuid = "a2964d1f-97da-50d4-b82a-358c7fce9d89"
 version = "1.24.0+0"
 
+[[deps.WeakCacheSets]]
+git-tree-sha1 = "386050ae4353310d8ff9c228f83b1affca2f7f38"
+uuid = "d30d5f5c-d141-4870-aa07-aabb0f5fe7d5"
+version = "0.1.0"
+
 [[deps.XML2_jll]]
 deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Zlib_jll"]
 git-tree-sha1 = "80d3930c6347cfce7ccf96bd3bafdf079d9c0390"
diff --git a/docs/src/tutorial-continuation.md b/docs/src/tutorial-continuation.md
index b9b085a..acabcb7 100644
--- a/docs/src/tutorial-continuation.md
+++ b/docs/src/tutorial-continuation.md
@@ -1,5 +1,9 @@
 # [Discrete continuation](@id tutorial-continuation)
 
+```@meta
+Draft = false
+```
+
 By using the warm start option, it is easy to implement a basic discrete continuation method, in which a sequence of problems is solved by using each solution as the initial guess for the next problem. This approach typically leads to faster and more reliable convergence than solving each problem with the same initial guess and is particularly useful for problems that require a good initial guess to converge.
 
 ## Continuation on parametric OCP
@@ -128,20 +132,18 @@ function goddard_problem(Tmax)
 
         tf ∈ R, variable
         t ∈ [0, tf], time
-        x ∈ R^3, state
+        x = (r, v, m) ∈ R^3, state
         u ∈ R, control
 
         0.01 ≤ tf ≤ Inf
 
-        r = x[1]
-        v = x[2]
-        m = x[3]
         x(0) == x0
         m(tf) == mf
         r0 ≤ r(t) ≤ r0 + 0.1
         v0 ≤ v(t) ≤ vmax
         mf ≤ m(t) ≤ m0
         0 ≤ u(t) ≤ 1
+
         ẋ(t) == F0(x(t)) + u(t) * F1(x(t), Tmax)
 
         r(tf) → max
@@ -212,6 +214,348 @@ end
 plot(plt_sol; size=(800, 800), leftmargin=5mm)
 ```
 
+### Animation
+
+The animation below shows three rockets launching simultaneously with different maximum thrust values (Tmax = 3, 2, 1). The rockets start at the same time (t=0) but reach their final altitudes at different times according to their optimal trajectories. Each rocket displays:
+
+- **Altitude** (above the rocket)
+- **Velocity gauge** (horizontal bar, green→red)
+- **Fuel gauge** (vertical bar, showing remaining fuel m(t) - mf)
+- **Flame** (visible when thrust u(t) > 0)
+
+```@setup main-cont
+# Extract data for the 3 selected solutions (Tmax ≈ 3, 2, 1)
+selected_sols = [sols[idx] for idx in indices]
+selected_Tmax = [round(data.Tmax[idx], digits=2) for idx in indices]
+
+# Extract time grids and state/control data for each solution
+solutions_data = []
+for (sol, Tmax_val) in zip(selected_sols, selected_Tmax)
+    t_grid = time_grid(sol)
+    X = state(sol).(t_grid)
+    u = control(sol).(t_grid)
+    
+    # Extract r, v, m from state
+    r = [x[1] for x in X]
+    v = [x[2] for x in X]
+    m = [x[3] for x in X]
+    
+    push!(solutions_data, (t=t_grid, r=r, v=v, m=m, u=u))
+end
+
+# Serialize to JSON for JavaScript injection
+json_t1 = "[" * join(string.(solutions_data[1].t), ",") * "]"
+json_r1 = "[" * join(string.(solutions_data[1].r), ",") * "]"
+json_v1 = "[" * join(string.(solutions_data[1].v), ",") * "]"
+json_m1 = "[" * join(string.(solutions_data[1].m), ",") * "]"
+json_u1 = "[" * join(string.(solutions_data[1].u), ",") * "]"
+
+json_t2 = "[" * join(string.(solutions_data[2].t), ",") * "]"
+json_r2 = "[" * join(string.(solutions_data[2].r), ",") * "]"
+json_v2 = "[" * join(string.(solutions_data[2].v), ",") * "]"
+json_m2 = "[" * join(string.(solutions_data[2].m), ",") * "]"
+json_u2 = "[" * join(string.(solutions_data[2].u), ",") * "]"
+
+json_t3 = "[" * join(string.(solutions_data[3].t), ",") * "]"
+json_r3 = "[" * join(string.(solutions_data[3].r), ",") * "]"
+json_v3 = "[" * join(string.(solutions_data[3].v), ",") * "]"
+json_m3 = "[" * join(string.(solutions_data[3].m), ",") * "]"
+json_u3 = "[" * join(string.(solutions_data[3].u), ",") * "]"
+
+# Inject real Tmax values for labels
+json_Tmax1 = string(selected_Tmax[1])
+json_Tmax2 = string(selected_Tmax[2])
+json_Tmax3 = string(selected_Tmax[3])
+
+# Inject m0 and mf values
+json_m0 = string(m0)
+json_mf = string(mf)
+
+# Inject r0 value
+json_r0 = string(r0)
+
+# Calculate global max time for animation loop
+t_max_global = max(solutions_data[1].t[end], solutions_data[2].t[end], solutions_data[3].t[end])
+
+# Define RawHTML wrapper for Documenter
+struct RawHTML
+    raw::String
+end
+Base.show(io::IO, ::MIME"text/html", h::RawHTML) = print(io, h.raw)
+
+html_anim = """
+<div style="display: flex; justify-content: center; margin: 20px 0;">
+    <canvas id="goddardCanvas" width="900" height="650"
+            style="border:1px solid #ddd; border-radius: 8px; max-width: 100%;">
+    </canvas>
+</div>
+
+<script>
+(function() {
+    // Data injected from Julia
+    const t1 = $json_t1, r1 = $json_r1, v1 = $json_v1, m1 = $json_m1, u1 = $json_u1;
+    const t2 = $json_t2, r2 = $json_r2, v2 = $json_v2, m2 = $json_m2, u2 = $json_u2;
+    const t3 = $json_t3, r3 = $json_r3, v3 = $json_v3, m3 = $json_m3, u3 = $json_u3;
+    
+    const Tmax1 = $json_Tmax1, Tmax2 = $json_Tmax2, Tmax3 = $json_Tmax3;
+    
+    const t_max_global = $t_max_global;
+    const m0 = $json_m0, mf = $json_mf;
+    const r0 = $json_r0;
+    
+    const canvas = document.getElementById('goddardCanvas');
+    const ctx = canvas.getContext('2d');
+    
+    // Detect dark/light mode
+    const isDark = window.matchMedia && window.matchMedia('(prefers-color-scheme: dark)').matches;
+    
+    // Color palettes
+    const colors = isDark ? {
+        canvas_bg: '#1e1e2e',
+        ground: '#4a4a4a',
+        gauge_bg: '#3a3a3a',
+        gauge_high: '#2ecc71',
+        gauge_mid: '#f39c12',
+        gauge_low: '#e74c3c',
+        text: '#ecf0f1',
+        text_secondary: '#bdc3c7',
+        flame_outer: '#e67e22',
+        flame_inner: '#f1c40f',
+        window_fill: 'rgba(255,255,255,0.3)',
+        window_stroke: 'rgba(255,255,255,0.6)',
+        progress: '#4063D8'
+    } : {
+        canvas_bg: '#fafafa',
+        ground: '#ccc',
+        gauge_bg: '#ecf0f1',
+        gauge_high: '#2ecc71',
+        gauge_mid: '#f39c12',
+        gauge_low: '#e74c3c',
+        text: '#2c3e50',
+        text_secondary: '#7f8c8d',
+        flame_outer: '#e67e22',
+        flame_inner: '#f1c40f',
+        window_fill: 'rgba(255,255,255,0.55)',
+        window_stroke: 'rgba(255,255,255,0.85)',
+        progress: '#4063D8'
+    };
+    
+    // Animation duration in seconds (independent of real problem time)
+    const animation_duration = 10.0;
+    
+    const zone_width = canvas.width / 3;
+    const top_margin = 50;   // Space for time label at top
+    const bottom_margin = 120; // Space for ground line + v gauge + Tmax label
+    const baseline = canvas.height - bottom_margin;
+    
+    // Scale for altitude (r goes from r0 to ~r0+0.1)
+    const r_min = r0;
+    const r_max_all = Math.max(...r1, ...r2, ...r3);
+    const r_range = r_max_all - r_min;
+    const altitude_scale = (baseline - top_margin) / r_range;
+    
+    // Scale for velocity (v goes from 0 to ~0.1)
+    const v_max_all = Math.max(...v1, ...v2, ...v3);
+    
+    let start_time = null;
+
+    function interpolate(t, data_t, data_y) {
+        // Find interval
+        let i = 0;
+        while(i < data_t.length - 1 && data_t[i + 1] < t) { i++; }
+        
+        if (i >= data_t.length - 1) {
+            return data_y[data_y.length - 1];  // Return final value if past end
+        }
+        
+        const dt = data_t[i+1] - data_t[i];
+        const alpha = (dt > 0) ? (t - data_t[i]) / dt : 0;
+        return data_y[i] + alpha * (data_y[i+1] - data_y[i]);
+    }
+
+    function drawRocket(zone_idx, t, r, v, m, u, color, Tmax) {
+        const zone_x = zone_idx * zone_width;
+        const center_x = zone_x + zone_width / 2;
+        
+        // Rocket position
+        const rocket_y = baseline - altitude_scale * (r - r_min);
+        
+        // Draw ground line
+        ctx.beginPath();
+        ctx.moveTo(zone_x + 20, baseline + 20);
+        ctx.lineTo(zone_x + zone_width - 20, baseline + 20);
+        ctx.lineWidth = 2;
+        ctx.strokeStyle = colors.ground;
+        ctx.stroke();
+        
+        // Draw fuel gauge (vertical bar to the right of rocket)
+        const fuel_x = center_x + 80;
+        const fuel_height = 200;
+        const fuel_y = baseline - fuel_height;
+        const fuel_level = (m - mf) / (m0 - mf);  // 0 to 1
+        
+        // Fuel gauge background
+        ctx.fillStyle = colors.gauge_bg;
+        ctx.fillRect(fuel_x - 10, fuel_y, 20, fuel_height);
+        
+        // Fuel gauge fill
+        const fuel_fill_height = fuel_level * fuel_height;
+        ctx.fillStyle = fuel_level > 0.5 ? colors.gauge_high : (fuel_level > 0.2 ? colors.gauge_mid : colors.gauge_low);
+        ctx.fillRect(fuel_x - 10, fuel_y + (fuel_height - fuel_fill_height), 20, fuel_fill_height);
+        
+        // Fuel gauge label
+        ctx.fillStyle = colors.text_secondary;
+        ctx.font = '14px Arial';
+        ctx.textAlign = 'center';
+        ctx.fillText('Fuel', fuel_x, fuel_y - 8);
+        
+        // Draw rocket: nose cone + body + fins + window
+        const rocket_bottom = rocket_y + 18;
+        ctx.fillStyle = color;
+
+        // Nose cone
+        ctx.beginPath();
+        ctx.moveTo(center_x, rocket_y - 38);
+        ctx.lineTo(center_x - 12, rocket_y - 10);
+        ctx.lineTo(center_x + 12, rocket_y - 10);
+        ctx.closePath();
+        ctx.fill();
+
+        // Body
+        ctx.fillRect(center_x - 12, rocket_y - 10, 24, 28);
+
+        // Left fin
+        ctx.beginPath();
+        ctx.moveTo(center_x - 12, rocket_y + 6);
+        ctx.lineTo(center_x - 24, rocket_y + 22);
+        ctx.lineTo(center_x - 12, rocket_y + 18);
+        ctx.closePath();
+        ctx.fill();
+
+        // Right fin
+        ctx.beginPath();
+        ctx.moveTo(center_x + 12, rocket_y + 6);
+        ctx.lineTo(center_x + 24, rocket_y + 22);
+        ctx.lineTo(center_x + 12, rocket_y + 18);
+        ctx.closePath();
+        ctx.fill();
+
+        // Window
+        ctx.beginPath();
+        ctx.arc(center_x, rocket_y + 2, 5, 0, 2 * Math.PI);
+        ctx.fillStyle = colors.window_fill;
+        ctx.fill();
+        ctx.strokeStyle = colors.window_stroke;
+        ctx.lineWidth = 1.5;
+        ctx.stroke();
+
+        // Draw flame if thrust > 0.1
+        if (u > 0.1) {
+            const flame_size = 10 + 20 * u;
+            ctx.fillStyle = colors.flame_outer;
+            ctx.beginPath();
+            ctx.moveTo(center_x - 10, rocket_bottom);
+            ctx.lineTo(center_x, rocket_bottom + flame_size);
+            ctx.lineTo(center_x + 10, rocket_bottom);
+            ctx.closePath();
+            ctx.fill();
+
+            // Inner flame
+            ctx.fillStyle = colors.flame_inner;
+            ctx.beginPath();
+            ctx.moveTo(center_x - 5, rocket_bottom);
+            ctx.lineTo(center_x, rocket_bottom + flame_size * 0.6);
+            ctx.lineTo(center_x + 5, rocket_bottom);
+            ctx.closePath();
+            ctx.fill();
+        }
+        
+        // Draw altitude text
+        ctx.fillStyle = colors.text;
+        ctx.font = 'bold 14px Arial';
+        ctx.textAlign = 'center';
+        ctx.fillText('r = ' + r.toFixed(4), center_x, rocket_y - 40);
+        
+        // Draw velocity gauge (horizontal bar)
+        const v_gauge_width = 60;
+        const v_gauge_height = 8;
+        const v_gauge_x = center_x - v_gauge_width / 2;
+        const v_gauge_y = rocket_y + 50;
+        const v_level = Math.abs(v) / v_max_all;
+        
+        // Velocity gauge background
+        ctx.fillStyle = colors.gauge_bg;
+        ctx.fillRect(v_gauge_x, v_gauge_y, v_gauge_width, v_gauge_height);
+        
+        // Velocity gauge fill
+        ctx.fillStyle = v_level > 0.7 ? colors.gauge_low : (v_level > 0.4 ? colors.gauge_mid : colors.gauge_high);
+        ctx.fillRect(v_gauge_x, v_gauge_y, v_level * v_gauge_width, v_gauge_height);
+        
+        // Velocity label
+        ctx.fillStyle = colors.text_secondary;
+        ctx.font = '12px Arial';
+        ctx.fillText('v', center_x, v_gauge_y + 20);
+        
+        // Draw Tmax label at bottom
+        ctx.fillStyle = colors.text;
+        ctx.font = 'bold 16px Arial';
+        ctx.fillText('Tmax = ' + Tmax, center_x, baseline + 90);
+    }
+
+    function draw(time) {
+        if (!start_time) start_time = time;
+        
+        const elapsed = (time - start_time) / 1000.0;
+        const progress = (elapsed % animation_duration) / animation_duration;
+        const t_anim = progress * t_max_global;
+        
+        ctx.fillStyle = colors.canvas_bg;
+        ctx.fillRect(0, 0, canvas.width, canvas.height);
+        
+        // Draw time indicator at top
+        ctx.fillStyle = colors.text;
+        ctx.font = 'bold 18px Arial';
+        ctx.textAlign = 'center';
+        ctx.fillText('t = ' + t_anim.toFixed(3) + ' s', canvas.width / 2, 30);
+        
+        // Interpolate and draw each rocket
+        const r1_cur = interpolate(t_anim, t1, r1);
+        const v1_cur = interpolate(t_anim, t1, v1);
+        const m1_cur = interpolate(t_anim, t1, m1);
+        const u1_cur = interpolate(t_anim, t1, u1);
+        drawRocket(0, t_anim, r1_cur, v1_cur, m1_cur, u1_cur, '#CB3C33', Tmax1);
+        
+        const r2_cur = interpolate(t_anim, t2, r2);
+        const v2_cur = interpolate(t_anim, t2, v2);
+        const m2_cur = interpolate(t_anim, t2, m2);
+        const u2_cur = interpolate(t_anim, t2, u2);
+        drawRocket(1, t_anim, r2_cur, v2_cur, m2_cur, u2_cur, '#389826', Tmax2);
+        
+        const r3_cur = interpolate(t_anim, t3, r3);
+        const v3_cur = interpolate(t_anim, t3, v3);
+        const m3_cur = interpolate(t_anim, t3, m3);
+        const u3_cur = interpolate(t_anim, t3, u3);
+        drawRocket(2, t_anim, r3_cur, v3_cur, m3_cur, u3_cur, '#9558B2', Tmax3);
+        
+        // Draw progress bar at bottom
+        const bar_height = 4;
+        ctx.fillStyle = colors.gauge_bg;
+        ctx.fillRect(0, canvas.height - bar_height, canvas.width, bar_height);
+        ctx.fillStyle = colors.progress;
+        ctx.fillRect(0, canvas.height - bar_height, canvas.width * progress, bar_height);
+        
+        requestAnimationFrame(draw);
+    }
+    requestAnimationFrame(draw);
+})();
+</script>
+"""
+```
+
+```@example main-cont
+RawHTML(html_anim)  # hide
+```
+
 ## Best practices and limitations
 
 The examples above illustrate a basic discrete continuation method. Here are some important considerations: