diff --git a/docs/src/examples/Optimal Control/Rocket.jl b/docs/src/examples/Optimal Control/Rocket.jl new file mode 100644 index 000000000..8ca5906f1 --- /dev/null +++ b/docs/src/examples/Optimal Control/Rocket.jl @@ -0,0 +1,100 @@ +# # Rocket Launch +# We will solve an optimal control problem that aims to minimize final time. + +# # Problem Statement and Model +# In this problem, we seek to minimize the final time ``t_f`` of a rocket launch between 0 and 10 metres. +# Minimizing fuel use is also important for this problem. The force- or thrust- of the rocket can be between +# -1.1 and 1.1, while the velocity of the rocket cannot exceed 1.7. The drag resistance of the rocket is +# proportional to the square of the velocity. The mass of the rocket decreases as fuel is burned over time. +# Our optimization guidelines are subject to these constraints: +# ```math +# \begin{aligned} +# &&\min t_f \\ +# &&&\frac{ds}{dt}= v \\ +# &&&\frac{dv}{dt}= (u-0.2*(v^2))/2\\ +# &&&\frac{dm}{dt}= -0.01(u^2)\\ +# &&&0 \leq v(t) \leq 1.7 \\ +# &&&-1.1 \leq u(t) \leq 1.1 \\ +# &&&s(0) = 0 \\ +# &&&v(0) = 0 \\ +# &&&m(0) = 1 \\ +# &&&s(t_f) = 10 \\ +# &&&v(t_f) = 0 \\ +# \end{aligned} +# ``` +# where ``v(t)`` represents the rocket's velocity, ``s(t)`` represents the rocket's position and +# ``m(t)`` represents the rocket's mass. Additionally, ``u(t)`` represents the rocket's force (or thrust). + +# # Modeling in InfiniteOpt +# First we must import ``InfiniteOpt`` and other packages. +using JuMP, InfiniteOpt, Ipopt, Plots; + +# # Model Initialization +# We then initialize the infinite model with [`InfiniteModel`](@ref), selecting Ipopt as the desired optimizer to solve this model. +model = InfiniteModel(Ipopt.Optimizer); + +# # Infinite Parameter Definition +# We are defining the infinite parameter ``t \in [0, 1]``. The bounds give the optimizer parameters to iterate over. +@infinite_parameter(model, t in [0,1], num_supports = 50, derivative_method = FiniteDifference(Backward())); + +# # Infinite Variable Definition +# Moving on to specifying variables, we have v as velocity of the rocket, u as the force of the rocket. +# The position, p, and mass, m, vary over time. The final time, t_f, must be between 0.1 and 100 seconds. +@variable(model, 0 <= v <= 1.7, Infinite(t)) +@variable(model, -1.1 <= u <= 1.1, Infinite(t)) +@variable(model, 0.1 <= t_f <= 100) +@variable(model, p, Infinite(t)) +@variable(model, m, Infinite(t)); + +# # Initial and Final Condition Definition +# Now that the model is defined, we then want to specify both our initial and final conditions. +@constraint(model, v(0) == 0) +@constraint(model, v(1) == 0) +@constraint(model, p(0) == 0) +@constraint(model, p(1) == 10) +@constraint(model, m(0) == 1) +@constraint(model, m(1) >= 0.2); + +# # Objective Definition +# We now add the objective, using `@objective`, which is to minimize final time t_f. +@objective(model, Min, t_f); + +# # Constraint Definition +# We finally are adding our constraints, which are defining the ODEs for our system model. +@constraint(model, v * t_f == deriv(p, t)) +@constraint(model, ((u - (0.2 * (v ^ 2))) * t_f) == (deriv(v, t)) * m) +@constraint(model, (-0.01 * (u ^ 2)) * t_f == deriv(m, t)); + +# # Problem Solution +# Now that everything is set up, we can solve the model by using `@optimize!` +optimize!(model) + +# # Extract and Plot the Results +# We can now extract the optimized results and plot them to visualize four results graphically: +# Position over time, velocity over time, mass over time and force over time. In addition, we can print the final time value of t_f. +tf_data = value(t_f) +t_data = supports(t) .* tf_data +v_data = value(v) +p_data = value(p) +m_data = value(m) +u_data = value(u); + +# Printing the final time value. +println("Minimized Final Time: ", tf_data," seconds"); + +# Creating the four plots as described above. +plot1 = plot(t_data, v_data, title="Velocity Over Time", linecolor = :red, xlabel="Time", label = "Velocity", ylabel="Velocity", lw=2); +plot2 = plot(t_data, p_data, title="Position Over Time", linecolor = :orange, xlabel="Time", label = "Position", ylabel="Position", lw=2); +plot3 = plot(t_data, m_data, title="Mass Over Time", linecolor = :yellow, xlabel="Time", label = "Position", ylabel="Position", lw=2); +plot4 = plot(t_data, u_data, title="Force Over Time", linecolor = :limegreen, xlabel="Time", label = "Position", ylabel="Position", lw=2); + +# Visualizing all four plots on one display. +display(plot(plot1, plot2, plot3, plot4, layout=(4,1), size = (1000, 1000))) + +# ### Maintenance Tests +# These are here to ensure this example stays up to date. +using Test +tol = 1E-4 +@test termination_status(model) == MOI.LOCALLY_SOLVED +@test has_values(model) +@test isapprox(objective_value(model), 7.4482495165303710, atol=tol) \ No newline at end of file