AME 527

Department of Aerospace and Mechanical Engineering

University of Southern California

Homework #4

Due: 6:30pm on October 23, 2019

Fall 2019

USC AME 527 Homework #4 Fall 2019

Problem 1: Sequential Quadratic Programming [40 points]

NOTE: Computations and plots shall be generated using Matlab or Python.

A well studied multidisciplinary design optimization (MDO) problem is Golinski’s Speed

Reducer Problem, which is a nonlinear constrained minimization problem. The problem aims

to minimize the weight of a gear-box while satisfying many constraints imposed by gear and

shaft design approaches. The definition of the MDO problem is shown in this reference:

Ray, Tapabrata, Golinski’s Speed Reducer Problem Revisited, AIAA Journal, Vol. 41, No. 3,

2002.

(a) [20 points] Solve Golinski’s Speed Reducer Problem using a sequential quadratic

programming tool (e.g., fmincon in Matlab) and continuous design variables. Report

your starting design variable values and the optimal design values found, including

the optimal final weight. Report the value of the constraint equations at the optimal

design point. Is your final design feasible? Which constraints are active? How does

your optimal design compare to that of the reference above?

(b) [20 points] Using sensitivity analysis, plot the relative sensitivity of the objective

function f(x∗) to each optimal design variable as a horizontal bar plot. Identify which

design variable the objective function is most sensitive to. With a first-order estimate,

how much would the objective function change if the most sensitive design variable

were increased by 1%, 5% and 15%? How much would it change if that variable were

decreased by those same percentages? Compare these estimates with the true objective

function value that is computed when this sensitive design variable is perturbed by the

percentages above. Which first-order estimate is a poor approximation?

1

USC AME 527 Homework #4 Fall 2019

Problem 2: 1000 mph Car Design [100 points]

NOTE: Computations and plots shall be generated using Matlab or Python.

As the chief designer of a new concept car expected to reach speeds of 1000 mph, you are

responsible for the mission modeling and safety of the driver while sizing the car. During

the conceptual phase of the design program you choose to rely on empirical and low-fidelity

models for aerodynamics, propulsion and mission profile. You are responsible for prepar-

ing a simulation that incorporates these models and determine if the vehicle can attain

1000 mph performance at the Black-Rock Desert in Nevada, the site of previous land-speed

records. The problem setup and models are described below, followed by the specific tasks

and deliverables to complete this assignment.

Problem Setup

The design variables for sizing the vehicle are summarized in Table 2.1.

Variable Units Initial Value Domain Description

x1 = h meters 3.0 [2, 3] Maximum vehicle height

x2 = w meters 2.5 [2, 3] Maximum vehicle width

x3 = L meters 12.9 [8, 12.9] Maximum vehicle length

x4 = Ln meters 5.0 [0.5, 5.0 ] Nose length

x5 = θ deg 10.0 [1, 90] Nose wedge angle

x6 = trocket1,on seconds 5.0 [1, 120] Rocket #1 ignition time

Table 2.1: Description of baseline design variables.

The parameters for this design problem are summarized in Table 2.2 and reflect condi-

tions at the Black-Rock Desert in Nevada (standard atmosphere, elevation 3907 ft). The

optimization objective function is J = −Vmax.

Mission Model

A non-linear ODE is used to model the horizontal motion of the vehicle (assuming a planar

desert surface) as derived in Lecture #5:

f [V (t)] = m

∂V (t)

∂t

− T [V (t), t] +mgµfriction + 1

2

ρV (t)2SrefCD[V (t)] = 0, (2.1)

where V (t) is the vehicle velocity, t is time, T [V (t)] is thrust, Sref is a reference area (here

we will use the vehicle maximum cross-sectional area projected forward), and CD[V (t)] is

the drag coefficient. After initializing with V (0) = 0.0, the ODE can be solved numerically

in time-steps of ∆t = 0.1 seconds. An ODE solver can be used (as found in Matlab), or

individual time-steps can be converged using a Newton Method while advancing the solution

from the initial condition.

2

USC AME 527 Homework #4 Fall 2019

Parameter Units Value Description

ρ kg/m3 1.0910 Air density

T Kelvin 280.4109 Air temperature

µair N-s/m

2 1.789E-5 Air viscosity

g m/s2 9.81 Gravitational acceleration

µfriction 0.01 Rolling friction coefficient

m kg 6422.0 Vehicle total mass

R J/kg-K 287.1 Gas constant

Tjet N 90000 Jet engine max thrust

Trocket# N 40000 Rocket motor max thrust for motor #

γ 1.4 Specific-heat ratio

∆tfull-power sec 60.0 Full-power duration

∆tdelay sec 10.0 Ignition/shutdown time interval

Table 2.2: Description of parameters.

For the Newton Method approach at time-step i, the velocity at the initial Newton

iteration k = 0 is Vk = 0 and the equation (2.1) is used in the residual equation F [Vk(ti)]k =

V (ti−1) + Vk(ti) + ∆tf [Vk(ti)]k = 0. A velocity perturbation ∆Vk(ti) is made using the

Jacobian ∂F [Vk(ti)]k/∂Vk(ti) with the equation

∆Vk = −F [Vk(ti)]k/(∂F [Vk(ti)]k/∂Vk(ti)), (2.2)

after which the velocity is updated as Vk+1(ti) = Vk(ti) + ∆Vk. This is repeated until a given

tolerance is reached (e.g. 1.0E-12) for ∆Vk or F [Vk(ti)]k, yielding the value of V (ti) for the

ith time-step. This process is repeated for additional time-steps ti+1 = ti + ∆t unil ti = 300

seconds. A time history of velocity V (t) can be stored and analyzed, as well as other vehicle

data that depends on it, such as drag coefficient CD[V (t)], thrust T [V (t), t], Mach M [V (t)],

distance x[V (t)] and acceleration a[V (t)], after exercising the Newton Method from t = 0 to

t = 300 seconds.

Propulsion Model

A simple propulsion model is used to schedule thrust as a function of velocity and time. For

the acceleration phase the ignition intervals are

trocket2,on = trocket1,on + ∆tdelay (2.3)

trocket3,on = trocket2,on + ∆tdelay (2.4)

and the thrust schedule is

T [V (t), t] =

0, V (t) = 0

Tjet, V (t) > 0

Tjet + Trocket, t ≥ trocket1,on

Tjet + 2Trocket, t ≥ trocket2,on

Tjet + 3Trocket, t ≥ trocket3,on

. (2.5)

3

USC AME 527 Homework #4 Fall 2019

For the deceleration phase, full-power is reduced via the shutdown interval

trocket1,off = trocket3,on + ∆tfull-power (2.6)

trocket2,off = trocket1,off + ∆tdelay (2.7)

trocket3,off = trocket2,off + ∆tdelay (2.8)

tjet,off = trocket3,off + ∆tdelay (2.9)

and thrust is scheduled to reduce from Tmax = Tjet + 3Trocket as

T [V (t), t] =

Tmax − Trocket, t ≥ trocket1,off

Tmax − 2Trocket, t ≥ trocket2,off

Tmax − 3Trocket, t ≥ trocket3,off

Tmax − 3Trocket − Tjet, t ≥ tjet,off

. (2.10)

Aerodynamics Model

Empirical models typically used for axisymmetric bodies were selected and modified slightly

to estimate aerodynamic drag coefficient as a function of Mach number M = V (t)/

√

γRT .

These are expressed in Python syntax in the accompanying file named aero model.py. The

syntax can be modified to work in a Matlab environment if desired.

Geometry Model

Geometric relationships required for the aerodynamics model are:

Sref = hw (2.11)

La = L− Ln (2.12)

d = w (2.13)

FR = L/d (2.14)

Ss = wL+ 2(0.5Ln0.4h+ Lah) (2.15)

dequiv = d (2.16)

db = 0.2w (2.17)

fn = Ln/dequiv (2.18)

fa = La/dequiv (2.19)

Deliverables

For each task below, optimize the vehicle design variables for maximum speed using a SQP

method (e.g. fmincon in Matlab). Report the information in Table 2.5 and plot the perfor-

mance time-history as in Figure 2.1 for each task. The time-history plots should include a

dashed-line for the initial design result overlayed with a solid-line for the final design result.

4

USC AME 527 Homework #4 Fall 2019

(a) [25 Points] Assuming the vehicle mass is unchanged (i.e. the consumed fuel is ig-

nored), program the given models above in order to simulate a car mission started at

V (0) = 0 and run the simulation from t = 0 to t = 300 seconds. The mission time-

history should be similar to what is shown in Figure 2.1 for the initial design point.

Optimize the design vector for maximum speed.

(b) [25 Points] To avoid unsafe control conditions for the driver beyond 1000 mph, aug-

ment and repeat the design optimization problem with a new constraint in units of

mph:

g1(x;p) = Vmax − 1050 ≤ 0. (2.20)

(c) [25 Points] The full desert expanse may become unsafe for driving due to poor wind

conditions, thus augment and repeat the design optimization problem again with a new

constraint and utilize a new design vector that includes aero-braking and wheel-braking

capabilities, described in Table 2.3 as new design variables.

Variable Units Initial Value Domain Description

x7 = Taerobrake Newtons -10000 [-50000, 0] Aero-braking force

x8 = Twheelbrake Newtons -10000 [-25000, 0] Wheel braking force

Table 2.3: Description of new design variable vector including aero-braking and wheel-

braking.

The new constraints in units of miles and mph, respectively, are:

g2(x;p) = xmax − 50 ≤ 0 (2.21)

g3(x;p) = V (t = 300) ≤ 0. (2.22)

The new braking capability will require modifing the deceleration schedule in the

Thrust model to the following:

T [V (t), t] =

Tmax − Trocket, t ≥ trocket1,off

Tmax − 2Trocket, t ≥ trocket2,off

Tmax − 3Trocket, t ≥ trocket3,off

Tmax − 3Trocket − Tjet, t ≥ tjet,off

Tmax − 3Trocket − Tjet + Taerobrake, t ≥ taerobrake,on

Tmax − 3Trocket − Tjet + Taerobrake + Twheelbrake, t ≥ twheelbrake,on

(2.23)

with additional time intervals

taerobrake,on = tjet,off + ∆tdelay (2.24)

twheelbrake,on = taerobrake,on + ∆tdelay. (2.25)

(d) [25 Points] In order to protect the driver from sudden accelerations that may cause

injury, augment and repeat the design optimization problem one more time with an

5

USC AME 527 Homework #4 Fall 2019

Variable Units Initial Value Domain Description

x9 = ∆tdelay seconds 10.0 [1, 60] Ignition/shutdown time interval

x10 = ∆tfull-power seconds 60.0 [60, 120] Full-power duration

Table 2.4: Description of new design variables.

additional constraint and changing two parameters into new design variables. The

design vector is augmented with the variables shown in Table 2.4. The new constraint

is non-dimensional in number of g-forces:

g4(x;p) =

[√

a2

]

max

/g − 2 ≤ 0. (2.26)

Value

x1

...

xn

Vmax

xmax

amax

Table 2.5: Reporting table summary of results.

6

USC AME 527 Homework #4 Fall 2019

0

50

10

0

15

0

20

0

25

0

30

0

0

20

0

40

0

60

0

80

0

10

00

Velocity [mph]

0

50

10

0

15

0

20

0

25

0

30

0

0.

0

0.

2

0.

4

0.

6

0.

8

Mach

0

50

10

0

15

0

20

0

25

0

30

0

0.

0

0.

1

0.

2

0.

3

0.

4

0.

5

D r a g C o e f f i c i e n t , C

D

0

50

10

0

15

0

20

0

25

0

30

0

05010

0

15

0

20

0

Thrust [kN]

0

50

10

0

15

0

20

0

25

0

30

0

Ti

m

e

[s

]

010203040

Distance [miles]

0

50

10

0

15

0

20

0

25

0

30

0

Ti

m

e

[s

]

1.

5

1.

0

0.

5

0.

0

0.

5

1.

0

1.

5

2.

0

Acceleration [#gs]

F

ig

u

re

2.

1:

T

im

e

h

is

to

ry

re

su

lt

s

fo

r

th

e

in

it

ia

l

d

es

ig

n

.

7