Dynamic Optimization

The focus of this course is on modeling, simulation, estimation, and optimization of dynamic systems. See the Engineering Design Optimization course for information on solving static design problems. This section of the course starts with dynamic modeling or methods to mathematically describe time-evolving systems, particularly for the purpose of dynamic optimization in engineering disciplines. The examples are particularly focused on engineering applications although much of the theory and applications can also be applied to other fields as well.

• Introduction to Dynamic Modeling (pdf)

The discussion starts with modeling because a reasonably accurate model of the system must first be created to approximate input to output relationships between values that can be adjusted (inputs) and those values that are used to judge the desirability of solution (outputs).

Dynamic Systems

Dynamic systems can often be described with differential equations. The equations can either be derived empirically from data or from fundamental relationships. One example of fundamental relationships are equations of motion for an object in a friction-less environment in one dimension. In this case velocity (v) is the equal to the time-derivative of position (y) and acceleration (a) is the time-derivative of velocity.

 dy/dt = v
 dv/dt = a

The above differential equations are expressed in semi-explicit form where the derivative terms are isolated on the left side of the equation and all other variables are on the right side of the expression. A more general notation for the semi-explicit form is dx/dt = f(x,p) where f(x,p) is any combination of variables (x) or parameters (p). Variables (x) are those values that are determined by the solution of the equations while parameters (p) are those quantities that are either specified by a user or determined by an optimizer. A more general form for differential equations is the open-equation format as f(dx/dt,x,p)=0 where all terms are brought to one side of the equation.

 dx/dt = f(x,p), Semi-Explicit Form

 f(dx/dt,x,p)=0, Open Equation Form

In the case of the object motion, the acceleration is not typically manipulated directly but may be adjusted by a force that is acting on that object. The relationship between acceleration (a) and force (F) includes another parameter, the mass of the object (m). There are now three equations that describe the motion of the system and relate force (F) to acceleration (a), velocity (v), and position (y).

 F = m a
 dy/dt = v
 dv/dt = a

The equation F = m a is an algebraic equation with no differential terms. While it would be simple to eliminate a from the equation by substituting for F/m, suppose that it is not possible or convenient to rearrange the equations to eliminate the algebraic expressions. This collection of equations is now referred to as Differential Algebraic Equations (DAEs). When there are no algebraic expressions, the system of equations are simplified to Ordinary Differential Equations (ODEs). If the differential equations have spatial and temporal derivatives, they become Partial Differential Equations (PDEs). For the purpose of this modeling discussion, any PDEs are discretized in space to return to ODE or DAE form where only time derivatives (not spatial derivatives) are present.

Dynamic versus Steady State

Dynamic systems are those where values that describe the system are expected to evolve over time and not necessarily remain at a steady state. When a dynamic system is at steady state, all time derivative values are either set to zero (dx/dt=0) or do not otherwise appear in the equations.

 Steady state: Values (x) do not change with time
 Dynamic: At least one value (x) changes with time

For physical systems the values of x may be measured at specific times. For virtual or simulated dynamic systems there can either be analytic solutions or a numerical solution.

Analytic Solution

An exact solution x(t) may exist from an analytic approach to solving the differential equations. In the in example case, the solution is found by integrating both sides of the differential equations. With the object initially at rest, this translates to zero initial conditions. An initial condition must be specified for each differential variable (those variables that appear in differential terms). In this case, the acceleration a is determined by the first equation but does not appear in differential form such as da/dt. Values of the parameters must be specified over the selected integration time horizon. The time horizon spans from the initial time where the initial conditions are specified to a future time t_f. Suppose there is a constant force (F=5) applied to the body starting at t₀=0. The analytic solution is then solved as a function of time.

 a = F/m
 v = a t + v0
 y = a t2/2 + v0 t + y0

Notice that the initial condition for acceleration (a₀) does not appear in the solution but initial conditions for velocity (v₀) and position (y₀) do influence the solution. If the force (F) changes at any time in the future, the time horizon is divided into separate segments where the force is constant. The first segment is computed and the final values become the initial values for the next segment.

Numerical Solution

Another method to solve dynamic system is with a numerical approach that approximates the exact solution x(t) at discrete time intervals, much like measuring a physical system at particular time points. The numerical approach is able to solve much larger and more complex systems of equations, especially when the exact analytic solution does not exist. One drawback to numerical solutions is that the accuracy of the solution is not exact and depends on factors such as the discretization methods employed and the error tolerance of the solver before convergence is confirmed.

Excel, MATLAB, Python, and Simulink are used in the following example to solve the differential equation that describes the velocity of a vehicle.

Simulate with Excel, MATLAB, Python, and Simulink

Mathematical Description of Dynamic Systems

Models are not just equations but are collections of assumptions, mathematical expressions, boundary conditions, initial conditions, and constraints. These mathematical expressions can range from a simple equation of motion for an object in a frictionless environment to complex systems that describe multi-body physics and interactions. A standard model form is shown below:

 0 = f(dx/dt,x,p)
 0 < g(dx/dt,x,p)

In this case, the differential and algebraic equations are condensed into f(dx/dt,x) and any algebraic expressions simply omit the derivative terms. Expressions may either be equality (f) or inequality (g) constraints.

Exercise

Objective: Provide a basic introduction to dynamic system modeling and simulation with the equations of motion. Create a MATLAB (ode23, ode15s, etc) or Python (ODEINT) script to simulate and display the results. Simulate with APM MATLAB or Python GEKKO. Compare the sequential (ODE integrators in MATLAB / Python) versus the simultaneous method (APMonitor/Gekko). Observe how the number of time points in APMonitor affects the solution accuracy. Estimated Time: 1 hour

Predict the position and velocity of a skydiver in two dimensions (horizontal and vertical) from the time of the initial jump through the first 90 seconds. At 60 seconds after the jump, the skydiver pulls the chute and the drag coefficient increases to slow the decent. The airplane is flying at a constant altitude of 5000 meters and 50 m/s when the skydiver jumps. The drag coefficient is 0.2 N-s²/m² while free-falling and 10 N-s²/m² with the parachute open. The gravitational constant is 9.8 m/s² and the mass is 80 kg for the skydiver and chute.

Solution

Skydiver Simulation in MATLAB and Python

This problem can also be solved with an ODE integrator such as Python's ODEINT function in SciPy.Integrate or MATLAB's ode23 function. The focus of this class is on solving dynamic optimization problems where ODE integrators are not suitable because of the inefficiencies of shooting methods. Regardless, it is valuable to know how to use ODE integrators for simulation. Below is the source code in MATLAB and Python.