diff --git a/docs/src/tutorials/custom_component.md b/docs/src/tutorials/custom_component.md index c7d77880d..f259848e8 100644 --- a/docs/src/tutorials/custom_component.md +++ b/docs/src/tutorials/custom_component.md @@ -114,14 +114,14 @@ Since the initial voltage of the capacitors was already specified via `v` and th ```@example components @mtkbuild sys = ChaoticAttractor() -prob = ODEProblem(sys, Pair[], (0, 5e4), saveat = 0.01) -sol = solve(prob) +prob = ODEProblem(sys, Pair[], (0, 5e4)) +sol = solve(prob; saveat = 1.0) -plot(sol[capacitor1.v], sol[capacitor2.v], title = "Chaotic Attractor", label = "", +plot(sol[sys.capacitor1.v], sol[sys.capacitor2.v], title = "Chaotic Attractor", label = "", ylabel = "C1 Voltage in V", xlabel = "C2 Voltage in V") ``` ```@example components -plot(sol; idxs = [capacitor1.v, capacitor2.v, inductor.i], +plot(sol; idxs = [sys.capacitor1.v, sys.capacitor2.v, sys.inductor.i], labels = ["C1 Voltage in V" "C2 Voltage in V" "Inductor Current in A"]) ``` diff --git a/docs/src/tutorials/dc_motor_pi.md b/docs/src/tutorials/dc_motor_pi.md index dbcc22f0b..cb712653d 100644 --- a/docs/src/tutorials/dc_motor_pi.md +++ b/docs/src/tutorials/dc_motor_pi.md @@ -20,51 +20,54 @@ using ModelingToolkitStandardLibrary.Mechanical.Rotational using ModelingToolkitStandardLibrary.Blocks using OrdinaryDiffEq using Plots - -R = 0.5 # [Ohm] armature resistance -L = 4.5e-3 # [H] armature inductance -k = 0.5 # [N.m/A] motor constant -J = 0.02 # [kg.m²] inertia -f = 0.01 # [N.m.s/rad] friction factor -tau_L_step = -0.3 # [N.m] amplitude of the load torque step -nothing # hide ``` The actual model can now be composed. ```@example dc_motor_pi -systems = @named begin - ground = Ground() - source = Voltage() - ref = Blocks.Step(height = 1, start_time = 0) - pi_controller = Blocks.LimPI(k = 1.1, T = 0.035, u_max = 10, Ta = 0.035) - feedback = Blocks.Feedback() - R1 = Resistor(R = R) - L1 = Inductor(L = L) - emf = EMF(k = k) - fixed = Fixed() - load = Torque() - load_step = Blocks.Step(height = tau_L_step, start_time = 3) - inertia = Inertia(J = J) - friction = Damper(d = f) - speed_sensor = SpeedSensor() +@mtkmodel DCMotor begin + @parameters begin + R = 0.5, [description = "Armature resistance"] # Ohm + L = 4.5e-3, [description = "Armature inductance"] # H + k = 0.5, [description = "Motor constant"] # N.m/A + J = 0.02, [description = "Inertia"] # kg.m² + f = 0.01, [description = "Friction factor"] # N.m.s/rad + tau_L_step = -0.3, [description = "Amplitude of the load torque step"] # N.m + end + @components begin + ground = Ground() + source = Voltage() + ref = Blocks.Step(height = 1, start_time = 0) + pi_controller = Blocks.LimPI(k = 1.1, T = 0.035, u_max = 10, Ta = 0.035) + feedback = Blocks.Feedback() + R1 = Resistor(R = R) + L1 = Inductor(L = L) + emf = EMF(k = k) + fixed = Fixed() + load = Torque() + load_step = Blocks.Step(height = tau_L_step, start_time = 3) + inertia = Inertia(J = J) + friction = Damper(d = f) + speed_sensor = SpeedSensor() + end + @equations begin + connect(fixed.flange, emf.support, friction.flange_b) + connect(emf.flange, friction.flange_a, inertia.flange_a) + connect(inertia.flange_b, load.flange) + connect(inertia.flange_b, speed_sensor.flange) + connect(load_step.output, load.tau) + connect(ref.output, feedback.input1) + connect(speed_sensor.w, :y, feedback.input2) + connect(feedback.output, pi_controller.err_input) + connect(pi_controller.ctr_output, :u, source.V) + connect(source.p, R1.p) + connect(R1.n, L1.p) + connect(L1.n, emf.p) + connect(emf.n, source.n, ground.g) + end end -connections = [connect(fixed.flange, emf.support, friction.flange_b) - connect(emf.flange, friction.flange_a, inertia.flange_a) - connect(inertia.flange_b, load.flange) - connect(inertia.flange_b, speed_sensor.flange) - connect(load_step.output, load.tau) - connect(ref.output, feedback.input1) - connect(speed_sensor.w, :y, feedback.input2) - connect(feedback.output, pi_controller.err_input) - connect(pi_controller.ctr_output, :u, source.V) - connect(source.p, R1.p) - connect(R1.n, L1.p) - connect(L1.n, emf.p) - connect(emf.n, source.n, ground.g)] - -@named model = ODESystem(connections, t; systems) +@named model = DCMotor() nothing # hide ``` @@ -75,12 +78,12 @@ so that it can be represented as a system of `ODEs` (ordinary differential equat ```@example dc_motor_pi sys = structural_simplify(model) prob = ODEProblem(sys, unknowns(sys) .=> 0.0, (0, 6.0)) -sol = solve(prob, Rodas4()) +sol = solve(prob) -p1 = plot(sol.t, sol[inertia.w], ylabel = "Angular Vel. in rad/s", +p1 = plot(sol.t, sol[sys.inertia.w], ylabel = "Angular Vel. in rad/s", label = "Measurement", title = "DC Motor with Speed Controller") -plot!(sol.t, sol[ref.output.u], label = "Reference") -p2 = plot(sol.t, sol[load.tau.u], ylabel = "Disturbance in Nm", label = "") +plot!(sol.t, sol[sys.ref.output.u], label = "Reference") +p2 = plot(sol.t, sol[sys.load.tau.u], ylabel = "Disturbance in Nm", label = "") plot(p1, p2, layout = (2, 1)) ``` @@ -89,8 +92,8 @@ plot(p1, p2, layout = (2, 1)) When implementing and tuning a control system in simulation, it is a good practice to analyze the closed-loop properties and verify robustness of the closed-loop with respect to, e.g., modeling errors. To facilitate this, we added two analysis points to the set of connections above, more specifically, we added the analysis points named `:y` and `:u` to the connections (for more details on analysis points, see [Linear Analysis](@ref)) ```julia -connect(speed_sensor.w, :y, feedback.input2) -connect(pi_controller.ctr_output, :u, source.V) +connect(sys.speed_sensor.w, :y, feedback.input2) +connect(sys.pi_controller.ctr_output, :u, source.V) ``` one at the plant output (`:y`) and one at the plant input (`:u`). We may use these analysis points to calculate, e.g., sensitivity functions, illustrated below. Here, we calculate the sensitivity function $S(s)$ and the complimentary sensitivity function $T(s) = I - S(s)$, defined as @@ -108,7 +111,7 @@ matrices_S, simplified_sys = Blocks.get_sensitivity( model, :y, op = Dict(unknowns(sys) .=> 0.0)) So = ss(matrices_S...) |> minreal # The output-sensitivity function as a StateSpace system matrices_T, simplified_sys = Blocks.get_comp_sensitivity( - model, :y, op = Dict(inertia.phi => 0.0, inertia.w => 0.0)) + model, :y, op = Dict(model.inertia.phi => 0.0, model.inertia.w => 0.0)) To = ss(matrices_T...)# The output complementary sensitivity function as a StateSpace system bodeplot([So, To], label = ["S" "T"], plot_title = "Sensitivity functions", plotphase = false) diff --git a/docs/src/tutorials/rc_circuit.md b/docs/src/tutorials/rc_circuit.md index 3ee518d44..1c05cb1a8 100644 --- a/docs/src/tutorials/rc_circuit.md +++ b/docs/src/tutorials/rc_circuit.md @@ -12,27 +12,32 @@ using ModelingToolkitStandardLibrary.Electrical using ModelingToolkitStandardLibrary.Blocks: Constant using ModelingToolkit: t_nounits as t -R = 1.0 -C = 1.0 -V = 1.0 -systems = @named begin - resistor = Resistor(R = R) - capacitor = Capacitor(C = C, v = 0.0) - source = Voltage() - constant = Constant(k = V) - ground = Ground() +@mtkmodel RC begin + @parameters begin + R = 1.0 + C = 1.0 + V = 1.0 + end + @components begin + resistor = Resistor(R = R) + capacitor = Capacitor(C = C, v = 0.0) + source = Voltage() + constant = Constant(k = V) + ground = Ground() + end + @equations begin + connect(constant.output, source.V) + connect(source.p, resistor.p) + connect(resistor.n, capacitor.p) + connect(capacitor.n, source.n, ground.g) + end end -rc_eqs = [connect(constant.output, source.V) - connect(source.p, resistor.p) - connect(resistor.n, capacitor.p) - connect(capacitor.n, source.n, ground.g)] - -@named rc_model = ODESystem(rc_eqs, t; systems) -sys = structural_simplify(rc_model) +@mtkbuild sys = RC() prob = ODEProblem(sys, Pair[], (0, 10.0)) -sol = solve(prob, Tsit5()) -plot(sol, idxs = [capacitor.v, resistor.i], +sol = solve(prob) + +plot(sol, idxs = [sys.capacitor.v, sys.resistor.i], title = "RC Circuit Demonstration", labels = ["Capacitor Voltage" "Resistor Current"]) ``` diff --git a/docs/src/tutorials/thermal_model.md b/docs/src/tutorials/thermal_model.md index a3a11a8d3..b52c9741a 100644 --- a/docs/src/tutorials/thermal_model.md +++ b/docs/src/tutorials/thermal_model.md @@ -10,31 +10,34 @@ from dividing the total initial energy in the system by the sum of the heat capa using ModelingToolkitStandardLibrary.Thermal, ModelingToolkit, OrdinaryDiffEq, Plots using ModelingToolkit: t_nounits as t -C1 = 15 -C2 = 15 -systems = @named begin - mass1 = HeatCapacitor(C = C1, T = 373.15) - mass2 = HeatCapacitor(C = C2, T = 273.15) - conduction = ThermalConductor(G = 10) - Tsensor1 = TemperatureSensor() - Tsensor2 = TemperatureSensor() +@mtkmodel HeatConductionModel begin + @parameters begin + C1 = 15 + C2 = 15 + end + @components begin + mass1 = HeatCapacitor(C = C1, T = 373.15) + mass2 = HeatCapacitor(C = C2, T = 273.15) + conduction = ThermalConductor(G = 10) + Tsensor1 = TemperatureSensor() + Tsensor2 = TemperatureSensor() + end + @equations begin + connect(mass1.port, conduction.port_a) + connect(conduction.port_b, mass2.port) + connect(mass1.port, Tsensor1.port) + connect(mass2.port, Tsensor2.port) + end end -connections = [ - connect(mass1.port, conduction.port_a), - connect(conduction.port_b, mass2.port), - connect(mass1.port, Tsensor1.port), - connect(mass2.port, Tsensor2.port) -] - -@named model = ODESystem(connections, t; systems) -sys = structural_simplify(model) +@mtkbuild sys = HeatConductionModel() prob = ODEProblem(sys, Pair[], (0, 5.0)) -sol = solve(prob, Tsit5()) +sol = solve(prob) -T_final_K = sol[(mass1.T * C1 + mass2.T * C2) / (C1 + C2)] +T_final_K = sol[(sys.mass1.T * sys.C1 + sys.mass2.T * sys.C2) / (sys.C1 + sys.C2)] plot(title = "Thermal Conduction Demonstration") -plot!(sol, idxs = [mass1.T, mass2.T], labels = ["Mass 1 Temperature" "Mass 2 Temperature"]) +plot!(sol, idxs = [sys.mass1.T, sys.mass2.T], + labels = ["Mass 1 Temperature" "Mass 2 Temperature"]) plot!(sol.t, T_final_K, label = "Steady-State Temperature") ```