Skip to content

Commit 8a37fcf

Browse files
baggepinnenbaggepinnen
committed
import fewer solvers
1 parent eddfe06 commit 8a37fcf

30 files changed

+262
-268
lines changed

‎Project.toml

Lines changed: 8 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -17,13 +17,12 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
1717
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
1818

1919
[weakdeps]
20-
Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
20+
Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
2121
JuliaFormatter = "98e50ef6-434e-11e9-1051-2b60c6c9e899"
2222
LightXML = "9c8b4983-aa76-5018-a973-4c85ecc9e179"
23-
Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
23+
Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
2424
MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377"
2525

26-
2726
[extensions]
2827
Render = ["Makie"]
2928
URDF = ["LightXML", "Graphs", "MetaGraphsNext", "JuliaFormatter"]
@@ -47,12 +46,14 @@ StaticArrays = "1"
4746
julia = "1"
4847

4948
[extras]
50-
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
51-
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
49+
Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
5250
JuliaFormatter = "98e50ef6-434e-11e9-1051-2b60c6c9e899"
5351
LightXML = "9c8b4983-aa76-5018-a973-4c85ecc9e179"
54-
Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
5552
MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377"
53+
OrdinaryDiffEqBDF = "6ad6398a-0878-4a85-9266-38940aa047c8"
54+
OrdinaryDiffEqRosenbrock = "43230ef6-c299-4910-a778-202eb28ce4ce"
55+
OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
56+
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
5657

5758
[targets]
58-
test = ["OrdinaryDiffEq", "Test", "JuliaFormatter", "LightXML", "Graphs", "MetaGraphsNext"]
59+
test = ["OrdinaryDiffEqBDF", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5", "Test", "JuliaFormatter", "LightXML", "Graphs", "MetaGraphsNext"]

‎docs/Manifest.toml

Lines changed: 128 additions & 220 deletions
Large diffs are not rendered by default.

‎docs/Project.toml

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,8 @@ Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
88
ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
99
ModelingToolkitStandardLibrary = "16a59e39-deab-5bd0-87e4-056b12336739"
1010
Multibody = "e1cad5d1-98ef-44f9-a79a-9ca4547f95b9"
11-
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
11+
OrdinaryDiffEqBDF = "6ad6398a-0878-4a85-9266-38940aa047c8"
12+
OrdinaryDiffEqRosenbrock = "43230ef6-c299-4910-a778-202eb28ce4ce"
13+
OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
1214
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
1315
RobustAndOptimalControl = "21fd56a4-db03-40ee-82ee-a87907bee541"

‎docs/src/examples/free_motion.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@ using Multibody
66
using ModelingToolkit
77
using Plots
88
using JuliaSimCompiler
9-
using OrdinaryDiffEq
9+
using OrdinaryDiffEqRosenbrock
1010
1111
t = Multibody.t
1212
D = Differential(t)

‎docs/src/examples/gearbox.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -11,7 +11,7 @@ using Multibody
1111
using ModelingToolkit
1212
using Plots
1313
using JuliaSimCompiler
14-
using OrdinaryDiffEq
14+
using OrdinaryDiffEqRosenbrock
1515
1616
t = Multibody.t
1717
D = Differential(t)

‎docs/src/examples/gyroscopic_effects.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -12,7 +12,7 @@ using Multibody
1212
using ModelingToolkit
1313
using Plots
1414
using JuliaSimCompiler
15-
using OrdinaryDiffEq
15+
using OrdinaryDiffEqBDF
1616
1717
t = Multibody.t
1818
D = Differential(t)

‎docs/src/examples/kinematic_loops.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -27,7 +27,7 @@ using Multibody
2727
using ModelingToolkit
2828
import ModelingToolkitStandardLibrary.Mechanical.Rotational
2929
using Plots
30-
using OrdinaryDiffEq
30+
using OrdinaryDiffEqBDF
3131
using LinearAlgebra
3232
using JuliaSimCompiler
3333

‎docs/src/examples/pendulum.md

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,7 @@ To start, we load the required packages
55
```@example pendulum
66
using ModelingToolkit
77
using Multibody, JuliaSimCompiler
8-
using OrdinaryDiffEq # Contains the ODE solver we will use
8+
using OrdinaryDiffEqRosenbrock # Contains the ODE solver we will use
99
using Plots
1010
```
1111
We then access the world frame and time variable from the Multibody module
@@ -57,7 +57,7 @@ This results in a simplified model with the minimum required variables and equat
5757
multibody
5858
```
5959

60-
We are now ready to create an `ODEProblem` and simulate it. We use the `Rodas4` solver from OrdinaryDiffEq.jl, and pass a dictionary for the initial conditions. We specify only initial condition for some variables, for those variables where no initial condition is specified, the default initial condition defined the model will be used.
60+
We are now ready to create an `ODEProblem` and simulate it. We use the `Rodas4` solver from OrdinaryDiffEqRosenbrock.jl, and pass a dictionary for the initial conditions. We specify only initial condition for some variables, for those variables where no initial condition is specified, the default initial condition defined the model will be used.
6161
```@example pendulum
6262
D = Differential(t)
6363
defs = Dict() # We may specify the initial condition here
@@ -187,7 +187,7 @@ The systems we have modeled so far have all been _planar_ mechanisms. We now ext
187187
This pendulum, sometimes referred to as a _rotary pendulum_, has two joints, one in the "shoulder", which is typically configured to rotate around the gravitational axis, and one in the "elbow", which is typically configured to rotate around the axis of the upper arm. The upper arm is attached to the shoulder joint, and the lower arm is attached to the elbow joint. The tip of the pendulum is attached to the lower arm.
188188

189189
```@example pendulum
190-
using ModelingToolkit, Multibody, JuliaSimCompiler, OrdinaryDiffEq, Plots
190+
using ModelingToolkit, Multibody, JuliaSimCompiler, OrdinaryDiffEqRosenbrock, Plots
191191
import ModelingToolkitStandardLibrary.Mechanical.Rotational.Damper as RDamper
192192
import Multibody.Rotations
193193
@@ -326,6 +326,7 @@ We will continue the pendulum theme and design an inverted pendulum on cart. The
326326
We start by putting the model together and control it in open loop using a simple periodic input signal:
327327

328328
```@example pendulum
329+
using OrdinaryDiffEqTsit5
329330
import ModelingToolkitStandardLibrary.Mechanical.TranslationalModelica
330331
import ModelingToolkitStandardLibrary.Blocks
331332
using Plots

‎docs/src/examples/prescribed_pose.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,7 @@ using Multibody
2323
using Multibody.Rotations # To specify orientations using Euler angles
2424
using ModelingToolkit
2525
using Plots
26-
using OrdinaryDiffEq
26+
using OrdinaryDiffEqRosenbrock
2727
using LinearAlgebra
2828
using JuliaSimCompiler
2929
using Test

‎docs/src/examples/quad.md

Lines changed: 77 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,7 +13,7 @@ using ModelingToolkitStandardLibrary.Blocks
1313
using LinearAlgebra
1414
using Plots
1515
using JuliaSimCompiler
16-
using OrdinaryDiffEq
16+
using OrdinaryDiffEqBDF
1717
using Test
1818
t = Multibody.t
1919
D = Differential(t)
@@ -336,3 +336,79 @@ nothing
336336
While gain and phase margins appear to be reasonable, we have a large high-frequency gain in the transfer functions from measurement noise to control signal, ``C(s)S(s)``. For a rotor craft where the control signal manipulates the current through motor windings, this may lead to excessive heat generation in the motors if the sensor measurements are noisy.
337337

338338
The diskmargin at the plant output is small, luckily, the gain variation appears at the plant input where the diskmargin is significantly larger. The diskmargins are visualized in the figure titled "Stable region for combined gain and phase variation". See [Diskmargin example](https://juliacontrol.github.io/RobustAndOptimalControl.jl/dev/#Diskmargin-example) to learn more about diskmargins.
339+
340+
```@example QUAD
341+
using LowLevelParticleFilters, SeeToDee, StaticArrays
342+
343+
Ts = 0.02
344+
tv = range(0, 7, step=Ts)
345+
Y = [SVector(i...) .+ 0.01 .* randn.() for i in sol(tv, idxs=outputs)]
346+
Y0 = [SVector(i...) for i in sol(tv, idxs=outputs)]
347+
U = [SVector(i...) .+ 0.2 .* randn.() for i in sol(tv, idxs=inputs)]
348+
X = Matrix(sol(tv, idxs=x_sym))
349+
##
350+
351+
state_parameters = [quad.body.m]
352+
353+
function dyn_with_state_parameters(mmodel, inputs, state_parameters)
354+
fdyn, x_sym, p_sym, io_sys = ModelingToolkit.generate_control_function(mmodel, inputs)
355+
356+
p_inds = [findfirst(x -> string(x) == string(p), p_sym) for p in state_parameters]
357+
# @info "found" p_sym[p_inds]
358+
359+
nx_orig = length(x_sym)
360+
nx_extended = nx_orig + length(state_parameters)
361+
362+
p1 = copy(p)
363+
dx = zeros(nx_orig)
364+
out = zeros(nx_extended)
365+
366+
f_ext = function (x, u, p, t)
367+
for (i, p_ind) in enumerate(p_inds)
368+
p1[p_ind] = max(x[nx_orig + i], 0.1) # NOTE: temp hack
369+
end
370+
dx .= 0
371+
fdyn(dx, x, u, p1, t)
372+
@views out[1:nx_orig] .= dx
373+
# [dx; x[nx_orig+1:end]]
374+
@views out[nx_orig+1:end] .= x[nx_orig+1:end]
375+
out
376+
end
377+
# f_ext, [x_sym; p_sym[p_inds]], p_sym[1:end-length(inputs)], io_sys
378+
f_ext, [x_sym; p_sym[p_inds]], p_sym, io_sys
379+
end
380+
381+
fdyn, x_sym, p_sym, io_sys = dyn_with_state_parameters(multibody(quad), inputs, state_parameters)
382+
383+
g = JuliaSimCompiler.build_explicit_observed_function(io_sys, outputs; inputs)
384+
nx = length(x_sym)
385+
nu = length(inputs)
386+
ny = length(outputs)
387+
x = randn(nx)
388+
u = randn(nu)
389+
390+
op = [
391+
quad.body.r_0[2] => 1e-32
392+
quad.v_alt => 1e-32 # To avoid singularity in linearization
393+
quad.world.g => 9.81
394+
inputs .=> 13;
395+
quad.body.m => 2 # We pretend that the body mass is close to cable+load+body mass for better altitude estimation
396+
] |> Dict
397+
##
398+
x0, p = JuliaSimCompiler.initial_conditions(io_sys, op, op)
399+
y0 = g(x0, u, p, 0)
400+
R1 = collect(0.001I(nx))
401+
R1[end] = 0.01
402+
R2 = 0.01^2 * I(ny)
403+
404+
d0 = LowLevelParticleFilters.MvNormal([x0; 2], R1)
405+
discrete_dynamics = SeeToDee.Rk4(fdyn, 0.01)
406+
407+
ukf = UnscentedKalmanFilter(discrete_dynamics, g, R1, R2, d0; p, nu, ny)
408+
409+
410+
@time filtersol = forward_trajectory(ukf, U, Y)
411+
plot(filtersol, size=(1200, 1000), ploty=false, plotu=false)
412+
X[end,:] .= 5.3
413+
plot!(X', sp=(1:nx)', label=false, l=(:black, :dash))
414+
```

‎docs/src/examples/robot.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,7 @@ using Multibody
77
using ModelingToolkit
88
using Plots
99
using JuliaSimCompiler
10-
using OrdinaryDiffEq
10+
using OrdinaryDiffEqRosenbrock
1111
using Test
1212
1313
t = Multibody.t

‎docs/src/examples/ropes_and_cables.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -16,7 +16,7 @@ using Multibody
1616
using ModelingToolkit
1717
using Plots
1818
using JuliaSimCompiler
19-
using OrdinaryDiffEq
19+
using OrdinaryDiffEqRosenbrock
2020
using Test
2121
t = Multibody.t
2222

‎docs/src/examples/sensors.md

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@ using Multibody
88
using ModelingToolkit
99
using JuliaSimCompiler
1010
using Plots
11-
using OrdinaryDiffEq
11+
using OrdinaryDiffEqRosenbrock
1212
using LinearAlgebra
1313
1414
t = Multibody.t
@@ -37,7 +37,6 @@ ssys = structural_simplify(multibody(model))
3737
D = Differential(t)
3838
prob = ODEProblem(ssys, [], (0, 3))
3939
40-
using OrdinaryDiffEq
4140
sol = solve(prob, Rodas4())
4241
@assert SciMLBase.successful_retcode(sol)
4342

‎docs/src/examples/space.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,7 @@ using Multibody
1919
using ModelingToolkit
2020
using Plots
2121
using JuliaSimCompiler
22-
using OrdinaryDiffEq
22+
using OrdinaryDiffEqRosenbrock
2323
2424
t = Multibody.t
2525
D = Differential(t)

‎docs/src/examples/spherical_pendulum.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -10,7 +10,7 @@ using Multibody
1010
using ModelingToolkit
1111
using Plots
1212
using JuliaSimCompiler
13-
using OrdinaryDiffEq
13+
using OrdinaryDiffEqRosenbrock
1414
1515
t = Multibody.t
1616
D = Differential(t)

‎docs/src/examples/spring_damper_system.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,7 +13,7 @@ using Multibody
1313
using ModelingToolkit
1414
using Plots
1515
using JuliaSimCompiler
16-
using OrdinaryDiffEq
16+
using OrdinaryDiffEqRosenbrock
1717
1818
t = Multibody.t
1919
D = Differential(t)

‎docs/src/examples/spring_mass_system.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -14,7 +14,7 @@ using ModelingToolkit
1414
using Plots
1515
using JuliaSimCompiler
1616
using ModelingToolkitStandardLibrary.Mechanical.TranslationalModelica
17-
using OrdinaryDiffEq
17+
using OrdinaryDiffEqRosenbrock
1818
1919
t = Multibody.t
2020
D = Differential(t)

‎docs/src/examples/suspension.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -11,7 +11,7 @@ using Multibody
1111
using ModelingToolkit
1212
import ModelingToolkitStandardLibrary.Mechanical.TranslationalModelica as Translational
1313
using Plots
14-
using OrdinaryDiffEq
14+
using OrdinaryDiffEqBDF, OrdinaryDiffEqRosenbrock
1515
using LinearAlgebra
1616
using JuliaSimCompiler
1717
using Test

‎docs/src/examples/swing.md

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -12,7 +12,7 @@ using Multibody
1212
using ModelingToolkit
1313
using Plots
1414
using JuliaSimCompiler
15-
using OrdinaryDiffEq
15+
using OrdinaryDiffEqBDF
1616
1717
t = Multibody.t
1818
D = Differential(t)
@@ -64,8 +64,9 @@ ssys = structural_simplify(multibody(model))
6464
prob = ODEProblem(ssys, [
6565
collect(model.body.v_0) .=> 0;
6666
collect(model.body.w_a) .=> 0;
67+
collect(model.body.r_0) .=> [0.5, 2, 0.1]; # To make the animation interesting
6768
], (0, 4))
68-
sol = solve(prob, ImplicitEuler(autodiff=false), reltol=5e-3)
69+
sol = solve(prob, QNDF1(autodiff=false), reltol=5e-3)
6970
@assert SciMLBase.successful_retcode(sol)
7071
```
7172

@@ -149,7 +150,7 @@ prob = ODEProblem(ssys, [
149150
], (0.0, 6))
150151
151152
152-
@time sol = solve(prob, ImplicitEuler(autodiff=false), reltol=1e-2)
153+
@time sol = solve(prob, QNDF1(autodiff=false), reltol=1e-2)
153154
@assert SciMLBase.successful_retcode(sol)
154155
155156
Plots.plot(sol, idxs = [collect(model.body.r_0);])
@@ -162,5 +163,3 @@ nothing # hide
162163
```
163164

164165
![animation](swing.gif)
165-
166-
There is an initial transient in the beginning where the springs in the ropes are vibrating substantially, but after about a second this transient is damped out (thanks in part to the, in this case desired, numerical damping contributed by the implicit Euler solver).

‎docs/src/examples/three_springs.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -12,7 +12,7 @@ using Multibody
1212
using ModelingToolkit
1313
using Plots
1414
using JuliaSimCompiler
15-
using OrdinaryDiffEq
15+
using OrdinaryDiffEqRosenbrock
1616
1717
t = Multibody.t
1818
D = Differential(t)

‎docs/src/examples/wheel.md

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -24,7 +24,7 @@ using ModelingToolkit
2424
import ModelingToolkitStandardLibrary.Mechanical.Rotational
2525
import ModelingToolkitStandardLibrary.Blocks
2626
using Plots
27-
using OrdinaryDiffEq
27+
using OrdinaryDiffEqRosenbrock, OrdinaryDiffEqTsit5
2828
using LinearAlgebra
2929
using JuliaSimCompiler
3030
using Test
@@ -307,7 +307,7 @@ import ModelingToolkitStandardLibrary.Mechanical.Rotational
307307
import ModelingToolkitStandardLibrary.Blocks
308308
import Multibody.PlanarMechanics as Pl
309309
using Plots
310-
using OrdinaryDiffEq
310+
using OrdinaryDiffEqRosenbrock
311311
using LinearAlgebra
312312
using JuliaSimCompiler
313313
using Test

‎docs/src/index.md

Lines changed: 13 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -9,7 +9,7 @@ Documentation for [Multibody](https://github.com/JuliaComputing/Multibody.jl).
99
```@setup logo
1010
using ModelingToolkit
1111
using Multibody, JuliaSimCompiler
12-
using OrdinaryDiffEq # Contains the ODE solver we will use
12+
using OrdinaryDiffEqRosenbrock # Contains the ODE solver we will use
1313
using Plots
1414
t = Multibody.t
1515
@@ -161,6 +161,18 @@ import Pkg
161161
Pkg.add("Multibody")
162162
```
163163

164+
In addition to Multibody.jl, we recommend installing the following **solver packages** to follow along with the examples in the documentation
165+
```julia
166+
Pkg.add([
167+
"OrdinaryDiffEqBDF", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5"
168+
])
169+
```
170+
Alternatively, install all available DAE solvers at the same time using (incurs lengthy pre compilation)
171+
```julia
172+
Pkg.add("OrdinaryDiffEq")
173+
```
174+
175+
164176

165177
## Notable differences from Modelica
166178

0 commit comments

Comments
 (0)