I have the following code:
using ModelingToolkit, OrdinaryDiffEq, DataInterpolations, Plots Cp=[0.024384377390186646,0.03811829845979822,0.05731175742000688,0.0803120051369546,0.10526323380995135,0.1312349486430945, 0.1588590278697544,0.1882995422205645,0.21937359147848456,0.25269290703920616,0.2882162631147021,0.3240355349501659, 0.3607290350346774,0.39295261085207084,0.4141821251393675,0.4274331262027119,0.43636422894857274,0.4414767530303278, 0.443412806515125,0.44365385445693295,0.44317448744425714,0.44218641629727234,0.4407405317795181,0.43888039054963845, 0.4362875461540461,0.4325702155989231,0.4278606704093836,0.4224263790651815,0.4164272616252002,0.40987027258773945, 0.4027832291503407,0.39518785578723936,0.38719847687832065] TSR = 2.0:0.25:10.0 cp = CubicSpline(Cp, TSR) @variables t ω(t) Pr(t) λ(t) @parameters J R ρ A U Pg D = Differential(t) eqs = [J * D(ω) * ω ~ Pr - Pg, λ ~ ω * R/U, Pr ~ 0.5 * ρ * A * U^3 * cp(λ) ] @named sys = ODESystem(eqs, t) sys = structural_simplify(sys) u0 = [ω => 0.9, # this is a fixed initial value D(ω) => 0.0] # the value for the derivative is an inital estimate p = [ J => 4.0470e+07, R => 63, ρ => 1.225, A => 1.2469e+04, U => 8.1, Pg => 1771e3] tspan = (0.0, 1000.0) prob = ODEProblem(sys, u0, tspan, p) sol = solve(prob, Rosenbrock23()) #### plot the result #### X = sol.t Y1 = [vec for vec in sol.u] Y2 = [vec for vec in sol.u] println("ω at 1000s: $(Y1[end])") display(plot(X, Y1, label="ω")) ### How to plot λ ? I have already the equation: eqs
Because I simplify the system of equations the variable \lambda does not appear in the solution. But I want to plot \lambda(t).
How can I do this without writing the second equation again?