Problem updates

available_problems_cache.txt

@@ -0,0 +1,20 @@
+beam
+bioreactor
+cart_pendulum
+chain
+dielectrophoretic_particle
+double_oscillator
+ducted_fan
+electric_vehicle
+glider
+insurance
+jackson
+moonlander
+quadrotor
+robbins
+robot
+rocket
+space_shuttle
+steering
+truck_trailer
+vanderpol

docs/src/list_of_problems.md

@@ -33,7 +33,7 @@ The table below summarizes the names and status of the each problem:
 | dielectrophoretic_particle    | ✅        | ✅    | 
 | double_oscillator             | ✅        | ✅    |
 | ducted_fan                    | ✅        | ✅    |
-| electrical_vehicle            | ✅        | ✅    |
+| electric_vehicle              | ✅        | ✅    |
 | glider                        | ✅        | ❌    |
 | insurance                     | ✅        | ✅    |
 | jackson                       | ✅        | ✅    |

ext/JuMPModels/beam.jl

@@ -39,7 +39,7 @@ function OptimalControlProblems.beam(::JuMPBackend; nh::Int=100)
         end
     )
 
-    @objective(model, Min, sum(u[t]^2 for t in 0:nh))
+    @objective(model, Min, (sum(u[t]^2 for t in 0:nh-1)) / nh)
 
     return model
 end

ext/JuMPModels/bioreactor.jl

@@ -2,7 +2,7 @@
 The Bioreactor Problem:
     The problem is formulated as a JuMP model and can be found [here](https://github.com/control-toolbox/bocop/tree/main/bocop)
 """
-function OptimalControlProblems.bioreactor(::JuMPBackend; nh::Int=100, N::Int=30)
+function OptimalControlProblems.bioreactor(::JuMPBackend; nh::Int=250, N::Int=30)
     # Parameters
     beta = 1
     c = 2
@@ -63,7 +63,7 @@ function OptimalControlProblems.bioreactor(::JuMPBackend; nh::Int=100, N::Int=30
         end
     )
 
-    @objective(model, Max, sum(b[t] / (beta + c) for t in 0:nh))
+    @objective(model, Max, step * sum(mu2[t] * b[t] / (beta + c) for t in 0:nh-1))
 
     return model
 end

ext/JuMPModels/cart_pendulum.jl

@@ -4,7 +4,7 @@ The Cart-Pendulum Problem:
     The objective is to swing the pendulum from the downward position to the upright position in the shortest time possible.      
     The problem is formulated as a JuMP model, and can be found [here](https://arxiv.org/pdf/2303.16746).
 """
-function OptimalControlProblems.cart_pendulum(::JuMPBackend; nh::Int64=100)
+function OptimalControlProblems.cart_pendulum(::JuMPBackend; nh::Int64=250)
     ## parameters
     g = 9.81      # gravitation [m/s^2]
     L = 1.0       # pendulum length [m]
@@ -21,13 +21,13 @@ function OptimalControlProblems.cart_pendulum(::JuMPBackend; nh::Int64=100)
     @variables(
         model,
         begin
-            0.0 <= tf
-            ddx
-            -max_x <= x[0:nh] <= max_x
-            -max_v <= dx[0:nh] <= max_v
-            theta[0:nh]
-            omega[0:nh]
-            -max_f <= Fex[0:nh] <= max_f
+            0.01 <= tf, (start = 0.1)
+            ddx, (start = 0.1)
+            -max_x <= x[0:nh] <= max_x, (start = 0.1)
+            -max_v <= dx[0:nh] <= max_v, (start = 0.1)
+            theta[0:nh], (start = 0.1)
+            omega[0:nh], (start = 0.1)
+            -max_f <= Fex[0:nh] <= max_f, (start = 0.1)
         end
     )
 
@@ -46,14 +46,16 @@ function OptimalControlProblems.cart_pendulum(::JuMPBackend; nh::Int64=100)
         model,
         begin
             step, tf / nh
-            COG[i=0:nh], L / 2 * [sin(theta[i]) + x[i], -cos(theta[i])] + [x[i], 0]
+            COG_1[i=0:nh], L / 2 * sin(theta[i]) + x[i] + x[i]
+            COG_2[i=0:nh], - L / 2 * -cos(theta[i])
             alpha[i=0:nh], 1.0 / (I + 0.25 * m * L^2) * 0.5 * L * m * (-ddx * cos(theta[i]) - g * sin(theta[i]))
-            ddCOG[i=0:nh], L * [-sin(theta[i]), cos(theta[i])] * omega[i] + L / 2 * [cos(theta[i]), sin(theta[i])] * alpha[i] + [ddx, 0]
-            FXFY[i=0:nh], m * ddCOG[i] + [0, m * g]
-            eq[i=0:nh], -FXFY[i][1] + Fex[i] - m_cart * ddx
-            J[i=0:nh], m_cart
-            c[i=0:nh], eq[i] - J[i] * ddx
-            ddx_[i=0:nh], -1.0 / J[i] * c[i]
+            ddCOG_1[i=0:nh], - L * sin(theta[i]) * omega[i] + L / 2 * cos(theta[i]) * alpha[i] + ddx
+            ddCOG_2[i=0:nh], L * cos(theta[i]) * omega[i] + L / 2 * sin(theta[i]) * alpha[i]
+            FXFY_1[i=0:nh], m * ddCOG_1[i]
+            FXFY_2[i=0:nh], m * ddCOG_2[i] + m * g
+            eq[i=0:nh], -FXFY_1[i] + Fex[i] - m_cart * ddx
+            c[i=0:nh], eq[i] - m_cart * ddx
+            ddx_[i=0:nh], -1.0 / m_cart * c[i]
             alpha_[i=0:nh], 1.0 / (I + 0.25 * m * L^2) * 0.5 * L * m * (-ddx_[i] * cos(theta[i]) - g * sin(theta[i]))
         end
     )

ext/JuMPModels/double_oscillator.jl

@@ -19,16 +19,16 @@ function OptimalControlProblems.double_oscillator(::JuMPBackend; nh::Int=100)
     @variables(
         model,
         begin
-            x1[0:nh]
-            x2[0:nh]
-            x3[0:nh]
-            x4[0:nh]
-            -1.0 <= u[0:nh] <= 1.0
+            x1[0:nh], (start = 0.1)
+            x2[0:nh], (start = 0.1)
+            x3[0:nh], (start = 0.1)
+            x4[0:nh], (start = 0.1)
+            -1.0 <= u[0:nh] <= 1.0, (start = 0.1)
         end
     )
 
     # Objective
-    @objective(model, Min, 0.5 * sum(x1[t]^2 + x2[t]^2 + u[t]^2 for t in 0:nh))
+    @objective(model, Min, 0.5 * step * sum(x1[t]^2 + x2[t]^2 + u[t]^2 for t in 0:nh-1))
 
     # Dynamics
     @expressions(

ext/JuMPModels/ducted_fan.jl

@@ -5,7 +5,7 @@ The Ducted Fan Problem:
     The problem is formulated as a JuMP model.
 Ref: Graichen, K., & Petit, N. (2009). Incorporating a class of constraints into the dynamics of optimal control problems. Optimal Control Applications and Methods, 30(6), 537-561.
 """
-function OptimalControlProblems.ducted_fan(::JuMPBackend; nh::Int=100)
+function OptimalControlProblems.ducted_fan(::JuMPBackend; nh::Int=250)
     r = 0.2         # [m]
     J = 0.05        # [kg.m2]
     m = 2.2         # [kg]
@@ -14,18 +14,16 @@ function OptimalControlProblems.ducted_fan(::JuMPBackend; nh::Int=100)
 
     model = Model()
 
-    @variable(model, x1[0:nh])
-    @variable(model, v1[0:nh])
-    @variable(model, x2[0:nh])
-    @variable(model, v2[0:nh])
-    @variable(model, -deg2rad(30.0) <= α[0:nh] <= deg2rad(30.0)) # radian
-    @variable(model, vα[0:nh])
-    @variable(model, -5.0 <= u1[0:nh] <= 5.0) # [nh]
-    @variable(model, 0.0 <= u2[0:nh] <= 17.0) # [nh]
+    @variable(model, x1[0:nh], start = 0.1)
+    @variable(model, v1[0:nh], start = 0.1)
+    @variable(model, x2[0:nh], start = 0.1)
+    @variable(model, v2[0:nh], start = 0.1)
+    @variable(model, -deg2rad(30.0) <= α[0:nh] <= deg2rad(30.0), start = 0.1) # radian
+    @variable(model, vα[0:nh], start = 0.1)
+    @variable(model, -5.0 <= u1[0:nh] <= 5.0, start = 0.1) # [nh]
+    @variable(model, 0.0 <= u2[0:nh] <= 17.0, start = 0.1) # [nh]
     @variable(model, 0 <= tf, start = 1.0)
 
-    @objective(model, Min, tf / nh * sum(2 * u1[t]^2 + u2[t]^2 for t in 0:nh) + μ * tf)
-
     # Dynamics
     @expressions(
         model,
@@ -70,5 +68,7 @@ function OptimalControlProblems.ducted_fan(::JuMPBackend; nh::Int=100)
         end
     )
 
+    @objective(model, Min, 1. / tf * step * sum(2. * u1[t]^2 + u2[t]^2 for t in 0:nh-1) + μ * tf)
+
     return model
 end

ext/JuMPModels/electrical_vehicle.jl → ext/JuMPModels/electric_vehicle.jl

@@ -1,10 +1,10 @@
 """
-The Electrical Vehicle Problem
-    Implement optimal control of an electrical vehicle.
+The electric Vehicle Problem
+    Implement optimal control of an electric vehicle.
     The problem is formulated as a JuMP model.
 Ref: [PS2011] Nicolas Petit and Antonio Sciarretta. "Optimal drive of electric vehicles using an inversion-based trajectory generation approach." IFAC Proceedings Volumes 44, no. 1 (2011): 14519-14526.
 """
-function OptimalControlProblems.electrical_vehicle(::JuMPBackend; nh::Int=100)
+function OptimalControlProblems.electric_vehicle(::JuMPBackend; nh::Int=250)
     D = 10.0
     T = 1.0
     b1 = 1e3
@@ -18,11 +18,11 @@ function OptimalControlProblems.electrical_vehicle(::JuMPBackend; nh::Int=100)
 
     model = Model()
 
-    @variable(model, x[0:nh])
-    @variable(model, v[0:nh])
-    @variable(model, u[0:nh])
+    @variable(model, x[0:nh], start = 0.1)
+    @variable(model, v[0:nh], start = 0.1)
+    @variable(model, u[0:nh], start = 0.1)
 
-    @objective(model, Min, sum(b1 * u[t] * v[t] + b2 * u[t]^2 for t in 0:nh))
+    @objective(model, Min, step * sum(b1 * u[t] * v[t] + b2 * u[t]^2 for t in 0:nh-1))
 
     # Dynamics
     @expressions(

ext/JuMPModels/glider.jl

@@ -30,7 +30,7 @@ function OptimalControlProblems.glider(::JuMPBackend; nh::Int64=100)
     @variables(
         model,
         begin
-            0 <= t_f, (start = 1.0)
+            0 <= tf, (start = 1.0)
             0.0 <= x[k=0:nh], (start = x_0 + vx_0 * (k / nh))
             y[k=0:nh], (start = y_0 + (k / nh) * (y_f - y_0))
             0.0 <= vx[k=0:nh], (start = vx_0)
@@ -44,7 +44,7 @@ function OptimalControlProblems.glider(::JuMPBackend; nh::Int64=100)
     @expressions(
         model,
         begin
-            step, t_f / nh
+            step, tf / nh
             r[i=0:nh], (x[i] / r_0 - 2.5)^2
             u[i=0:nh], u_c * (1 - r[i]) * exp(-r[i])
             w[i=0:nh], vy[i] - u[i]

ext/JuMPModels/insurance.jl

@@ -2,7 +2,7 @@
 The Insurance Problem:
     The problem is formulated as a JuMP model and can be found [here](https://github.com/control-toolbox/bocop/tree/main/bocop)
 """
-function OptimalControlProblems.insurance(::JuMPBackend; nh::Int=100, N::Int=30)
+function OptimalControlProblems.insurance(::JuMPBackend; nh::Int=100)
     # constants
     gamma = 0.2
     lambda = 0.25
@@ -82,7 +82,7 @@ function OptimalControlProblems.insurance(::JuMPBackend; nh::Int=100, N::Int=30)
     )
 
     # Objective
-    @objective(model, Max, sum(U[t] * fx[t] for t in 0:nh))
+    @objective(model, Max, step * sum(U[t] * fx[t] for t in 0:nh-1))
 
     return model
 end

ext/JuMPModels/moonlander.jl

@@ -5,7 +5,7 @@ The Moonlander Problem:
     The problem is formulated as a JuMP model, and can be found [here](https://arxiv.org/pdf/2303.16746)
 """
 function OptimalControlProblems.moonlander(
-    ::JuMPBackend; target::Array{Float64}=[5.0, 5.0], nh::Int64=100
+    ::JuMPBackend; target::Array{Float64}=[5.0, 5.0], nh::Int64=1000
 )
     ## parameters
     if size(target) != (2,)
@@ -23,14 +23,14 @@ function OptimalControlProblems.moonlander(
     @variables(
         model,
         begin
-            0.0 <= step
+            0.01 <= tf, (start = 0.1)
             # state variables
-            p1[k=0:nh]
-            p2[k=0:nh]
-            dp1[k=0:nh]
-            dp2[k=0:nh]
-            theta[k=0:nh]
-            dtheta[k=0:nh]
+            p1[k=0:nh], (start = 0.1)
+            p2[k=0:nh], (start = 0.1)
+            dp1[k=0:nh], (start = 0.1)
+            dp2[k=0:nh], (start = 0.1)
+            theta[k=0:nh], (start = 0.1)
+            dtheta[k=0:nh], (start = 0.1)
             # control variables
             0 <= F1[k=0:nh] <= max_thrust, (start = 5.0)
             0 <= F2[k=0:nh] <= max_thrust, (start = 5.0)
@@ -55,6 +55,13 @@ function OptimalControlProblems.moonlander(
     )
 
     #dynamics
+    @expressions(
+        model,
+        begin
+            step, tf / nh
+        end
+    )
+
     @expressions(
         model,
         begin
@@ -95,7 +102,7 @@ function OptimalControlProblems.moonlander(
         end
     )
 
-    @objective(model, Min, step * nh)
+    @objective(model, Min, tf)
 
     return model
 end

ext/JuMPModels/quadrotor.jl

@@ -22,17 +22,17 @@ function OptimalControlProblems.quadrotor(::JuMPBackend; nh::Int64=60)
     @variables(
         model,
         begin
-            0.0 <= tf
-            p1[0:nh]
-            p2[0:nh]
-            p3[0:nh]
-            v1[0:nh]
-            v2[0:nh]
-            v3[0:nh]
+            0.0 <= tf, (start = 0.1)
+            p1[0:nh], (start = 0.1)
+            p2[0:nh], (start = 0.1)
+            p3[0:nh], (start = 0.1)
+            v1[0:nh], (start = 0.1)
+            v2[0:nh], (start = 0.1)
+            v3[0:nh], (start = 0.1)
             atmin <= at[0:nh] <= atmax, (start = 10.0)
-            -pi / 2 <= ϕ[0:nh] <= pi / 2
-            -pi / 2 <= θ[0:nh] <= pi / 2
-            ψ[0:nh]
+            -pi / 2 <= ϕ[0:nh] <= pi / 2, (start = 0.1)
+            -pi / 2 <= θ[0:nh] <= pi / 2, (start = 0.1)
+            ψ[0:nh], (start = 0.1)
         end
     )
 

ext/JuMPModels/robbins.jl

@@ -2,7 +2,7 @@
 The Robbins Problem:
     The problem is formulated as a JuMP model and can be found [here](https://github.com/control-toolbox/bocop/tree/main/bocop)
 """
-function OptimalControlProblems.robbins(::JuMPBackend; nh::Int=100, N::Int=30)
+function OptimalControlProblems.robbins(::JuMPBackend; nh::Int=250)
     # constants
     alpha = 3
     beta = 0
@@ -48,7 +48,7 @@ function OptimalControlProblems.robbins(::JuMPBackend; nh::Int=100, N::Int=30)
 
     # Objective
     @objective(
-        model, Min, sum(alpha * x1[t] + beta * x1[t]^2 + gamma * u[t]^2 for t in 0:nh)
+        model, Min, step * sum(alpha * x1[t] + beta * x1[t]^2 + gamma * u[t]^2 for t in 0:nh-1)
     )
 
     return model

ext/JuMPModels/space_Shuttle.jl → ext/JuMPModels/space_shuttle.jl

@@ -1,11 +1,12 @@
 """
 Space Shuttle Reentry Trajectory Problem:
     We want to find the optimal trajectory of a space shuttle reentry.
-    The objective is to minimize the angle of attack at the terminal point.
+    The objective is to maximize the latitude (cross range) at the terminal point.
     The problem is formulated as a JuMP model, and can be found [here](https://jump.dev/JuMP.jl/stable/tutorials/nonlinear/space_shuttle_reentry_trajectory/)
+    Note: no heating limit path constraint
 """
 function OptimalControlProblems.space_shuttle(
-    ::JuMPBackend; integration_rule::String="rectangular", nh::Int64=503
+    ::JuMPBackend; integration_rule::String="trapezoidal", nh::Int64=503
 )
     ## Global variables
     w = 203000.0  # weight (lb)