User-Provided Mesh Example

Here is an example script with a custom mesh provided as an array of coordinates. This should help you understand how meshes are defined in OpenAeroStruct and how to create them for your own custom planform shapes. This is an alternative to the helper-function approach described in Aerodynamic Optimization.

The following shows the portion of the example script in which the user provides the coordinates for the mesh. This example is for a wing with a kink and two distinct trapezoidal segments.

# -----------------------------------------------------------------------------
# CUSTOM MESH: Example mesh for a 2-segment wing with sweep
# -----------------------------------------------------------------------------

# Planform specifications
half_span = 12.0  # wing half-span in m
kink_location = 4.0  # spanwise location of the kink in m

root_chord = 6.0  # root chord in m
kink_chord = 3.0  # kink chord in m
tip_chord = 2.0  # tip chord in m

inboard_LE_sweep = 10.0  # inboard leading-edge sweep angle in deg
outboard_LE_sweep = -10.0  # outboard leading-edge sweep angle in deg

# Mesh specifications
nx = 5  # number of chordwise nodal points (should be odd)
ny_outboard = 9  # number of spanwise nodal points for the outboard segment
ny_inboard = 7  # number of spanwise nodal points for the inboard segment

# Initialize the 3-D mesh object. Indexing: Chordwise, spanwise, then the 3-D coordinates.
# We use ny_inboard+ny_outboard-1 because the 2 segments share the nodes where they connect.
mesh = np.zeros((nx, ny_inboard + ny_outboard - 1, 3))

# The form of this 3-D array can be confusing initially.
# For each node, we are providing the x, y, and z coordinates.
# x is streamwise, y is spanwise, and z is up.
# For example, the node for the leading edge at the tip would be specified as mesh[0, 0, :] = np.array([x, y, z]).
# And the node at the trailing edge at the root would be mesh[nx-1, ny-1, :] = np.array([x, y, z]).
# We only provide the right half of the wing here because we use symmetry.
# Print elements of the mesh to better understand the form.

####### THE Z-COORDINATES ######
# Assume no dihedral, so set the z-coordinate for all the points to 0.
mesh[:, :, 2] = 0.0

####### THE Y-COORDINATES ######
# Using uniform spacing for the spanwise locations of all the nodes within each of the two trapezoidal segments:
# Outboard
mesh[:, :ny_outboard, 1] = np.linspace(half_span, kink_location, ny_outboard)
# Inboard
mesh[:, ny_outboard : ny_outboard + ny_inboard, 1] = np.linspace(kink_location, 0, ny_inboard)[1:]

###### THE X-COORDINATES ######
# Start with the leading edge and create some intermediate arrays that we will use
x_LE = np.zeros(ny_inboard + ny_outboard - 1)

array_for_inboard_leading_edge_x_coord = np.linspace(0, kink_location, ny_inboard) * np.tan(
    inboard_LE_sweep / 180.0 * np.pi
)

array_for_outboard_leading_edge_x_coord = (
    np.linspace(0, half_span - kink_location, ny_outboard) * np.tan(outboard_LE_sweep / 180.0 * np.pi)
    + np.ones(ny_outboard) * array_for_inboard_leading_edge_x_coord[-1]
)

x_LE[:ny_inboard] = array_for_inboard_leading_edge_x_coord
x_LE[ny_inboard : ny_inboard + ny_outboard] = array_for_outboard_leading_edge_x_coord[1:]

# Then the trailing edge
x_TE = np.zeros(ny_inboard + ny_outboard - 1)

array_for_inboard_trailing_edge_x_coord = np.linspace(
    array_for_inboard_leading_edge_x_coord[0] + root_chord,
    array_for_inboard_leading_edge_x_coord[-1] + kink_chord,
    ny_inboard,
)

array_for_outboard_trailing_edge_x_coord = np.linspace(
    array_for_outboard_leading_edge_x_coord[0] + kink_chord,
    array_for_outboard_leading_edge_x_coord[-1] + tip_chord,
    ny_outboard,
)

x_TE[:ny_inboard] = array_for_inboard_trailing_edge_x_coord
x_TE[ny_inboard : ny_inboard + ny_outboard] = array_for_outboard_trailing_edge_x_coord[1:]

# # Quick plot to check leading and trailing edge x-coords
# plt.plot(x_LE, np.arange(0, ny_inboard+ny_outboard-1), marker='*')
# plt.plot(x_TE, np.arange(0, ny_inboard+ny_outboard-1), marker='*')
# plt.show()
# exit()

for i in range(0, ny_inboard + ny_outboard - 1):
    mesh[:, i, 0] = np.linspace(np.flip(x_LE)[i], np.flip(x_TE)[i], nx)

# -----------------------------------------------------------------------------
# END MESH
# -----------------------------------------------------------------------------

The following shows a visualization of the mesh.

../_images/two_part_mesh.png

The complete script for the optimization is as follows. Make sure you go through the Aerostructural Optimization with Wingbox before trying to understand this setup.

import warnings
import matplotlib
warnings.filterwarnings('ignore')
matplotlib.use('Agg')
"""
This example script can be used to run a multipoint aerostructural (w/ wingbox) optimization for a custom user-provided mesh.

The fuel burn from the cruise flight-point is the objective function and a 2.5g
maneuver flight-point is used for the structural sizing.
After running the optimization, use the 'plot_wingbox.py' script in the utils/
directory (e.g., as 'python ../utils/plot_wingbox.py aerostruct.db' if running
from this directory) to visualize the results.
'plot_wingbox.py' is still under development and will probably not work as it is for other types of cases for now.
"""

import numpy as np

from openaerostruct.integration.aerostruct_groups import AerostructGeometry, AerostructPoint
from openaerostruct.structures.wingbox_fuel_vol_delta import WingboxFuelVolDelta
import openmdao.api as om


# docs checkpoint 0
# -----------------------------------------------------------------------------
# CUSTOM MESH: Example mesh for a 2-segment wing with sweep
# -----------------------------------------------------------------------------

# Planform specifications
half_span = 12.0  # wing half-span in m
kink_location = 4.0  # spanwise location of the kink in m

root_chord = 6.0  # root chord in m
kink_chord = 3.0  # kink chord in m
tip_chord = 2.0  # tip chord in m

inboard_LE_sweep = 10.0  # inboard leading-edge sweep angle in deg
outboard_LE_sweep = -10.0  # outboard leading-edge sweep angle in deg

# Mesh specifications
nx = 5  # number of chordwise nodal points (should be odd)
ny_outboard = 9  # number of spanwise nodal points for the outboard segment
ny_inboard = 7  # number of spanwise nodal points for the inboard segment

# Initialize the 3-D mesh object. Indexing: Chordwise, spanwise, then the 3-D coordinates.
# We use ny_inboard+ny_outboard-1 because the 2 segments share the nodes where they connect.
mesh = np.zeros((nx, ny_inboard + ny_outboard - 1, 3))

# The form of this 3-D array can be confusing initially.
# For each node, we are providing the x, y, and z coordinates.
# x is streamwise, y is spanwise, and z is up.
# For example, the node for the leading edge at the tip would be specified as mesh[0, 0, :] = np.array([x, y, z]).
# And the node at the trailing edge at the root would be mesh[nx-1, ny-1, :] = np.array([x, y, z]).
# We only provide the right half of the wing here because we use symmetry.
# Print elements of the mesh to better understand the form.

####### THE Z-COORDINATES ######
# Assume no dihedral, so set the z-coordinate for all the points to 0.
mesh[:, :, 2] = 0.0

####### THE Y-COORDINATES ######
# Using uniform spacing for the spanwise locations of all the nodes within each of the two trapezoidal segments:
# Outboard
mesh[:, :ny_outboard, 1] = np.linspace(half_span, kink_location, ny_outboard)
# Inboard
mesh[:, ny_outboard : ny_outboard + ny_inboard, 1] = np.linspace(kink_location, 0, ny_inboard)[1:]

###### THE X-COORDINATES ######
# Start with the leading edge and create some intermediate arrays that we will use
x_LE = np.zeros(ny_inboard + ny_outboard - 1)

array_for_inboard_leading_edge_x_coord = np.linspace(0, kink_location, ny_inboard) * np.tan(
    inboard_LE_sweep / 180.0 * np.pi
)

array_for_outboard_leading_edge_x_coord = (
    np.linspace(0, half_span - kink_location, ny_outboard) * np.tan(outboard_LE_sweep / 180.0 * np.pi)
    + np.ones(ny_outboard) * array_for_inboard_leading_edge_x_coord[-1]
)

x_LE[:ny_inboard] = array_for_inboard_leading_edge_x_coord
x_LE[ny_inboard : ny_inboard + ny_outboard] = array_for_outboard_leading_edge_x_coord[1:]

# Then the trailing edge
x_TE = np.zeros(ny_inboard + ny_outboard - 1)

array_for_inboard_trailing_edge_x_coord = np.linspace(
    array_for_inboard_leading_edge_x_coord[0] + root_chord,
    array_for_inboard_leading_edge_x_coord[-1] + kink_chord,
    ny_inboard,
)

array_for_outboard_trailing_edge_x_coord = np.linspace(
    array_for_outboard_leading_edge_x_coord[0] + kink_chord,
    array_for_outboard_leading_edge_x_coord[-1] + tip_chord,
    ny_outboard,
)

x_TE[:ny_inboard] = array_for_inboard_trailing_edge_x_coord
x_TE[ny_inboard : ny_inboard + ny_outboard] = array_for_outboard_trailing_edge_x_coord[1:]

# # Quick plot to check leading and trailing edge x-coords
# plt.plot(x_LE, np.arange(0, ny_inboard+ny_outboard-1), marker='*')
# plt.plot(x_TE, np.arange(0, ny_inboard+ny_outboard-1), marker='*')
# plt.show()
# exit()

for i in range(0, ny_inboard + ny_outboard - 1):
    mesh[:, i, 0] = np.linspace(np.flip(x_LE)[i], np.flip(x_TE)[i], nx)

# -----------------------------------------------------------------------------
# END MESH
# -----------------------------------------------------------------------------
# docs checkpoint 1


# -----------------------------------------------------------------------------
# On to the problem setup (this is the same setup used for the Q400 example)
# -----------------------------------------------------------------------------

# Provide coordinates for a portion of an airfoil for the wingbox cross-section as an nparray with dtype=complex (to work with the complex-step approximation for derivatives).
# These should be for an airfoil with the chord scaled to 1.
# We use the 10% to 60% portion of the NASA SC2-0612 airfoil for this case
# We use the coordinates available from airfoiltools.com. Using such a large number of coordinates is not necessary.
# The first and last x-coordinates of the upper and lower surfaces must be the same

# fmt: off
upper_x = np.array([0.1, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16, 0.17, 0.18, 0.19, 0.2, 0.21, 0.22, 0.23, 0.24, 0.25, 0.26, 0.27, 0.28, 0.29, 0.3, 0.31, 0.32, 0.33, 0.34, 0.35, 0.36, 0.37, 0.38, 0.39, 0.4, 0.41, 0.42, 0.43, 0.44, 0.45, 0.46, 0.47, 0.48, 0.49, 0.5, 0.51, 0.52, 0.53, 0.54, 0.55, 0.56, 0.57, 0.58, 0.59, 0.6], dtype="complex128")
lower_x = np.array([0.1, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16, 0.17, 0.18, 0.19, 0.2, 0.21, 0.22, 0.23, 0.24, 0.25, 0.26, 0.27, 0.28, 0.29, 0.3, 0.31, 0.32, 0.33, 0.34, 0.35, 0.36, 0.37, 0.38, 0.39, 0.4, 0.41, 0.42, 0.43, 0.44, 0.45, 0.46, 0.47, 0.48, 0.49, 0.5, 0.51, 0.52, 0.53, 0.54, 0.55, 0.56, 0.57, 0.58, 0.59, 0.6], dtype="complex128")
upper_y = np.array([ 0.0447,  0.046,  0.0472,  0.0484,  0.0495,  0.0505,  0.0514,  0.0523,  0.0531,  0.0538, 0.0545,  0.0551,  0.0557, 0.0563,  0.0568, 0.0573,  0.0577,  0.0581,  0.0585,  0.0588,  0.0591,  0.0593,  0.0595,  0.0597,  0.0599,  0.06,    0.0601,  0.0602,  0.0602,  0.0602,  0.0602,  0.0602,  0.0601,  0.06,    0.0599,  0.0598,  0.0596,  0.0594,  0.0592,  0.0589,  0.0586,  0.0583,  0.058,   0.0576,  0.0572,  0.0568,  0.0563,  0.0558,  0.0553,  0.0547,  0.0541], dtype="complex128")  # noqa: E201, E241
lower_y = np.array([-0.0447, -0.046, -0.0473, -0.0485, -0.0496, -0.0506, -0.0515, -0.0524, -0.0532, -0.054, -0.0547, -0.0554, -0.056, -0.0565, -0.057, -0.0575, -0.0579, -0.0583, -0.0586, -0.0589, -0.0592, -0.0594, -0.0595, -0.0596, -0.0597, -0.0598, -0.0598, -0.0598, -0.0598, -0.0597, -0.0596, -0.0594, -0.0592, -0.0589, -0.0586, -0.0582, -0.0578, -0.0573, -0.0567, -0.0561, -0.0554, -0.0546, -0.0538, -0.0529, -0.0519, -0.0509, -0.0497, -0.0485, -0.0472, -0.0458, -0.0444], dtype="complex128")
# fmt: on

surf_dict = {
    # Wing definition
    "name": "wing",  # name of the surface
    "symmetry": True,  # if true, model one half of wing
    "S_ref_type": "wetted",  # how we compute the wing area,
    # can be 'wetted' or 'projected'
    "mesh": mesh,
    "twist_cp": np.array([6.0, 7.0, 7.0, 7.0]),
    "fem_model_type": "wingbox",
    "data_x_upper": upper_x,
    "data_x_lower": lower_x,
    "data_y_upper": upper_y,
    "data_y_lower": lower_y,
    "spar_thickness_cp": np.array([0.004, 0.004, 0.004, 0.004]),  # [m]
    "skin_thickness_cp": np.array([0.003, 0.006, 0.010, 0.012]),  # [m]
    "original_wingbox_airfoil_t_over_c": 0.12,
    # Aerodynamic deltas.
    # These CL0 and CD0 values are added to the CL and CD
    # obtained from aerodynamic analysis of the surface to get
    # the total CL and CD.
    # These CL0 and CD0 values do not vary wrt alpha.
    # They can be used to account for things that are not included, such as contributions from the fuselage, nacelles, tail surfaces, etc.
    "CL0": 0.0,
    "CD0": 0.0142,
    "with_viscous": True,  # if true, compute viscous drag
    "with_wave": True,  # if true, compute wave drag
    # Airfoil properties for viscous drag calculation
    "k_lam": 0.05,  # percentage of chord with laminar
    # flow, used for viscous drag
    "c_max_t": 0.38,  # chordwise location of maximum thickness
    "t_over_c_cp": np.array([0.1, 0.1, 0.15, 0.15]),
    # Structural values are based on aluminum 7075
    "E": 73.1e9,  # [Pa] Young's modulus
    "G": (73.1e9 / 2 / 1.33),  # [Pa] shear modulus (calculated using E and the Poisson's ratio here)
    "yield": (420.0e6 / 1.5),  # [Pa] allowable yield stress
    "mrho": 2.78e3,  # [kg/m^3] material density
    "strength_factor_for_upper_skin": 1.0,  # the yield stress is multiplied by this factor for the upper skin
    "wing_weight_ratio": 1.25,
    "exact_failure_constraint": False,  # if false, use KS function
    "struct_weight_relief": True,
    "distributed_fuel_weight": True,
    "fuel_density": 803.0,  # [kg/m^3] fuel density (only needed if the fuel-in-wing volume constraint is used)
    "Wf_reserve": 500.0,  # [kg] reserve fuel mass
}

surfaces = [surf_dict]

# Create the problem and assign the model group
prob = om.Problem()

# Add problem information as an independent variables component
indep_var_comp = om.IndepVarComp()
indep_var_comp.add_output("v", val=np.array([0.5 * 310.95, 0.3 * 340.294]), units="m/s")
indep_var_comp.add_output("alpha", val=0.0, units="deg")
indep_var_comp.add_output("alpha_maneuver", val=0.0, units="deg")
indep_var_comp.add_output("Mach_number", val=np.array([0.5, 0.3]))
indep_var_comp.add_output(
    "re",
    val=np.array([0.569 * 310.95 * 0.5 * 1.0 / (1.56 * 1e-5), 1.225 * 340.294 * 0.3 * 1.0 / (1.81206 * 1e-5)]),
    units="1/m",
)
indep_var_comp.add_output("rho", val=np.array([0.569, 1.225]), units="kg/m**3")
indep_var_comp.add_output("CT", val=0.43 / 3600, units="1/s")
indep_var_comp.add_output("R", val=2e6, units="m")
indep_var_comp.add_output("W0", val=25400 + surf_dict["Wf_reserve"], units="kg")
indep_var_comp.add_output("speed_of_sound", val=np.array([310.95, 340.294]), units="m/s")
indep_var_comp.add_output("load_factor", val=np.array([1.0, 2.5]))
indep_var_comp.add_output("empty_cg", val=np.zeros((3)), units="m")
indep_var_comp.add_output("fuel_mass", val=3000.0, units="kg")

prob.model.add_subsystem("prob_vars", indep_var_comp, promotes=["*"])

# Loop over each surface in the surfaces list
for surface in surfaces:

    # Get the surface name and create a group to contain components
    # only for this surface
    name = surface["name"]

    aerostruct_group = AerostructGeometry(surface=surface)

    # Add group to the problem with the name of the surface.
    prob.model.add_subsystem(name, aerostruct_group)

# Loop through and add a certain number of aerostruct points
for i in range(2):

    point_name = "AS_point_{}".format(i)
    # Connect the parameters within the model for each aerostruct point

    # Create the aero point group and add it to the model
    AS_point = AerostructPoint(surfaces=surfaces, internally_connect_fuelburn=False)

    prob.model.add_subsystem(point_name, AS_point)

    # Connect flow properties to the analysis point
    prob.model.connect("v", point_name + ".v", src_indices=[i])
    prob.model.connect("Mach_number", point_name + ".Mach_number", src_indices=[i])
    prob.model.connect("re", point_name + ".re", src_indices=[i])
    prob.model.connect("rho", point_name + ".rho", src_indices=[i])
    prob.model.connect("CT", point_name + ".CT")
    prob.model.connect("R", point_name + ".R")
    prob.model.connect("W0", point_name + ".W0")
    prob.model.connect("speed_of_sound", point_name + ".speed_of_sound", src_indices=[i])
    prob.model.connect("empty_cg", point_name + ".empty_cg")
    prob.model.connect("load_factor", point_name + ".load_factor", src_indices=[i])
    prob.model.connect("fuel_mass", point_name + ".total_perf.L_equals_W.fuelburn")
    prob.model.connect("fuel_mass", point_name + ".total_perf.CG.fuelburn")

    for surface in surfaces:

        name = surface["name"]

        if surf_dict["distributed_fuel_weight"]:
            prob.model.connect("load_factor", point_name + ".coupled.load_factor", src_indices=[i])

        com_name = point_name + "." + name + "_perf."
        prob.model.connect(
            name + ".local_stiff_transformed", point_name + ".coupled." + name + ".local_stiff_transformed"
        )
        prob.model.connect(name + ".nodes", point_name + ".coupled." + name + ".nodes")

        # Connect aerodynamic mesh to coupled group mesh
        prob.model.connect(name + ".mesh", point_name + ".coupled." + name + ".mesh")
        if surf_dict["struct_weight_relief"]:
            prob.model.connect(name + ".element_mass", point_name + ".coupled." + name + ".element_mass")

        # Connect performance calculation variables
        prob.model.connect(name + ".nodes", com_name + "nodes")
        prob.model.connect(name + ".cg_location", point_name + "." + "total_perf." + name + "_cg_location")
        prob.model.connect(name + ".structural_mass", point_name + "." + "total_perf." + name + "_structural_mass")

        # Connect wingbox properties to von Mises stress calcs
        prob.model.connect(name + ".Qz", com_name + "Qz")
        prob.model.connect(name + ".J", com_name + "J")
        prob.model.connect(name + ".A_enc", com_name + "A_enc")
        prob.model.connect(name + ".htop", com_name + "htop")
        prob.model.connect(name + ".hbottom", com_name + "hbottom")
        prob.model.connect(name + ".hfront", com_name + "hfront")
        prob.model.connect(name + ".hrear", com_name + "hrear")

        prob.model.connect(name + ".spar_thickness", com_name + "spar_thickness")
        prob.model.connect(name + ".t_over_c", com_name + "t_over_c")

prob.model.connect("alpha", "AS_point_0" + ".alpha")
prob.model.connect("alpha_maneuver", "AS_point_1" + ".alpha")

# Here we add the fuel volume constraint componenet to the model
prob.model.add_subsystem("fuel_vol_delta", WingboxFuelVolDelta(surface=surface))
prob.model.connect("wing.struct_setup.fuel_vols", "fuel_vol_delta.fuel_vols")
prob.model.connect("AS_point_0.fuelburn", "fuel_vol_delta.fuelburn")

if surf_dict["distributed_fuel_weight"]:
    prob.model.connect("wing.struct_setup.fuel_vols", "AS_point_0.coupled.wing.struct_states.fuel_vols")
    prob.model.connect("fuel_mass", "AS_point_0.coupled.wing.struct_states.fuel_mass")

    prob.model.connect("wing.struct_setup.fuel_vols", "AS_point_1.coupled.wing.struct_states.fuel_vols")
    prob.model.connect("fuel_mass", "AS_point_1.coupled.wing.struct_states.fuel_mass")

comp = om.ExecComp("fuel_diff = (fuel_mass - fuelburn) / fuelburn", units="kg")
prob.model.add_subsystem("fuel_diff", comp, promotes_inputs=["fuel_mass"], promotes_outputs=["fuel_diff"])
prob.model.connect("AS_point_0.fuelburn", "fuel_diff.fuelburn")


## Use these settings if you do not have pyOptSparse or SNOPT
prob.driver = om.ScipyOptimizeDriver()
prob.driver.options["optimizer"] = "SLSQP"
prob.driver.options["tol"] = 1e-4

recorder = om.SqliteRecorder("aerostruct.db")
prob.driver.add_recorder(recorder)

# We could also just use prob.driver.recording_options['includes']=['*'] here, but for large meshes the database file becomes extremely large. So we just select the variables we need.
prob.driver.recording_options["includes"] = [
    "alpha",
    "rho",
    "v",
    "cg",
    "AS_point_1.cg",
    "AS_point_0.cg",
    "AS_point_0.coupled.wing_loads.loads",
    "AS_point_1.coupled.wing_loads.loads",
    "AS_point_0.coupled.wing.normals",
    "AS_point_1.coupled.wing.normals",
    "AS_point_0.coupled.wing.widths",
    "AS_point_1.coupled.wing.widths",
    "AS_point_0.coupled.aero_states.wing_sec_forces",
    "AS_point_1.coupled.aero_states.wing_sec_forces",
    "AS_point_0.wing_perf.CL1",
    "AS_point_1.wing_perf.CL1",
    "AS_point_0.coupled.wing.S_ref",
    "AS_point_1.coupled.wing.S_ref",
    "wing.geometry.twist",
    "wing.mesh",
    "wing.skin_thickness",
    "wing.spar_thickness",
    "wing.t_over_c",
    "wing.structural_mass",
    "AS_point_0.wing_perf.vonmises",
    "AS_point_1.wing_perf.vonmises",
    "AS_point_0.coupled.wing.def_mesh",
    "AS_point_1.coupled.wing.def_mesh",
]

prob.driver.recording_options["record_objectives"] = True
prob.driver.recording_options["record_constraints"] = True
prob.driver.recording_options["record_desvars"] = True
prob.driver.recording_options["record_inputs"] = True

prob.model.add_objective("AS_point_0.fuelburn", scaler=1e-5)

prob.model.add_design_var("wing.twist_cp", lower=-15.0, upper=15.0, scaler=0.1)
prob.model.add_design_var("wing.spar_thickness_cp", lower=0.003, upper=0.1, scaler=1e2)
prob.model.add_design_var("wing.skin_thickness_cp", lower=0.003, upper=0.1, scaler=1e2)
prob.model.add_design_var("wing.geometry.t_over_c_cp", lower=0.07, upper=0.2, scaler=10.0)
prob.model.add_design_var("fuel_mass", lower=0.0, upper=2e5, scaler=1e-5)
prob.model.add_design_var("alpha_maneuver", lower=-15.0, upper=15)

prob.model.add_constraint("AS_point_0.CL", equals=0.6)
prob.model.add_constraint("AS_point_1.L_equals_W", equals=0.0)
prob.model.add_constraint("AS_point_1.wing_perf.failure", upper=0.0)

prob.model.add_constraint("fuel_vol_delta.fuel_vol_delta", lower=0.0)
prob.model.add_constraint("fuel_diff", equals=0.0)

# Set up the problem
prob.setup()

# om.view_model(prob)

# prob.check_partials(form='central', compact_print=True)

prob.run_driver()
==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 59474.4786 1
NL: NLBGS 2 ; 58887.5619 0.990131622
NL: NLBGS 3 ; 2197.7072 0.0369521054
NL: NLBGS 4 ; 67.9021348 0.00114170206
NL: NLBGS 5 ; 0.185879282 3.12536212e-06
NL: NLBGS 6 ; 0.00905010794 1.52167924e-07
NL: NLBGS 7 ; 7.21557777e-05 1.21322254e-09
NL: NLBGS 8 ; 1.54661908e-07 2.60047523e-12
NL: NLBGS 9 ; 2.80110505e-09 4.70975975e-14
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 63459.5261 1
NL: NLBGS 2 ; 54460.64 0.858194874
NL: NLBGS 3 ; 1868.49116 0.0294438246
NL: NLBGS 4 ; 53.4050567 0.000841560913
NL: NLBGS 5 ; 0.136184357 2.14600337e-06
NL: NLBGS 6 ; 0.00600916275 9.469284e-08
NL: NLBGS 7 ; 5.31043859e-05 8.3682292e-10
NL: NLBGS 8 ; 8.32399117e-08 1.31170081e-12
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 1.76433335e-10 1
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 1.8599746e-09 1
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 12750.4121 1
NL: NLBGS 2 ; 11862.2719 0.930344201
NL: NLBGS 3 ; 1908.36166 0.149670587
NL: NLBGS 4 ; 88.3873596 0.00693211788
NL: NLBGS 5 ; 0.96287727 7.55173452e-05
NL: NLBGS 6 ; 0.0841092157 6.59658803e-06
NL: NLBGS 7 ; 0.00233803599 1.83369446e-07
NL: NLBGS 8 ; 2.70095336e-05 2.11832633e-09
NL: NLBGS 9 ; 1.64718835e-06 1.29187068e-10
NL: NLBGS 10 ; 1.53681498e-07 1.20530612e-11
NL: NLBGS 11 ; 5.87751086e-09 4.60966345e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 99441.5133 1
NL: NLBGS 2 ; 97674.107 0.982226675
NL: NLBGS 3 ; 16423.2767 0.165155137
NL: NLBGS 4 ; 936.87811 0.00942139836
NL: NLBGS 5 ; 35.1249448 0.000353222146
NL: NLBGS 6 ; 1.740572 1.75034746e-05
NL: NLBGS 7 ; 0.0617139255 6.20605252e-07
NL: NLBGS 8 ; 0.00344667647 3.46603381e-08
NL: NLBGS 9 ; 6.4077225e-05 6.44370977e-10
NL: NLBGS 10 ; 1.18287202e-05 1.1895153e-10
NL: NLBGS 11 ; 2.12569352e-06 2.13763191e-11
NL: NLBGS 12 ; 5.23896522e-08 5.26838847e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 5539.38445 1
NL: NLBGS 2 ; 5477.27013 0.988786784
NL: NLBGS 3 ; 590.741732 0.106643931
NL: NLBGS 4 ; 34.2111635 0.00617598649
NL: NLBGS 5 ; 0.168616946 3.04396541e-05
NL: NLBGS 6 ; 0.0113631236 2.05133328e-06
NL: NLBGS 7 ; 0.000326404857 5.89243913e-08
NL: NLBGS 8 ; 1.57198278e-06 2.83782936e-10
NL: NLBGS 9 ; 7.4455816e-08 1.34411714e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 11777.1844 1
NL: NLBGS 2 ; 10952.5376 0.929979285
NL: NLBGS 3 ; 1788.44443 0.151856706
NL: NLBGS 4 ; 60.9773026 0.00517757897
NL: NLBGS 5 ; 0.437788331 3.717258e-05
NL: NLBGS 6 ; 0.0508793958 4.32016634e-06
NL: NLBGS 7 ; 0.000849896845 7.21646883e-08
NL: NLBGS 8 ; 1.59615921e-05 1.3552978e-09
NL: NLBGS 9 ; 7.85299277e-07 6.66797128e-11
NL: NLBGS 10 ; 3.41265488e-08 2.89768314e-12
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 8614.67151 1
NL: NLBGS 2 ; 8654.39792 1.00461148
NL: NLBGS 3 ; 663.934336 0.0770701861
NL: NLBGS 4 ; 53.7962891 0.00624472902
NL: NLBGS 5 ; 0.203597221 2.36337765e-05
NL: NLBGS 6 ; 0.00470142446 5.45746225e-07
NL: NLBGS 7 ; 0.000316492563 3.67387848e-08
NL: NLBGS 8 ; 9.10961111e-07 1.0574531e-10
NL: NLBGS 9 ; 1.62422942e-08 1.88542236e-12
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 8614.53461 1
NL: NLBGS 2 ; 8050.13429 0.934482785
NL: NLBGS 3 ; 913.897946 0.106087907
NL: NLBGS 4 ; 63.1593862 0.00733172355
NL: NLBGS 5 ; 0.444320326 5.15779838e-05
NL: NLBGS 6 ; 0.0142120316 1.64977359e-06
NL: NLBGS 7 ; 0.000765748307 8.889027e-08
NL: NLBGS 8 ; 4.09141381e-06 4.74943105e-10
NL: NLBGS 9 ; 1.12308452e-07 1.30370887e-11
NL: NLBGS 10 ; 8.73013894e-09 1.01341968e-12
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 7540.79319 1
NL: NLBGS 2 ; 7560.41376 1.00260192
NL: NLBGS 3 ; 531.841463 0.0705285836
NL: NLBGS 4 ; 44.7595975 0.00593566172
NL: NLBGS 5 ; 0.177697437 2.35648204e-05
NL: NLBGS 6 ; 0.00266036312 3.52796192e-07
NL: NLBGS 7 ; 0.000195485334 2.59237097e-08
NL: NLBGS 8 ; 4.46564206e-07 5.92197922e-11
NL: NLBGS 9 ; 6.78440399e-09 8.99693682e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 7725.02638 1
NL: NLBGS 2 ; 7176.13027 0.92894573
NL: NLBGS 3 ; 551.565217 0.0713997842
NL: NLBGS 4 ; 53.136538 0.00687849275
NL: NLBGS 5 ; 0.323731268 4.19068171e-05
NL: NLBGS 6 ; 0.00299241284 3.87366035e-07
NL: NLBGS 7 ; 0.00023467972 3.0379148e-08
NL: NLBGS 8 ; 3.75173794e-07 4.85660211e-11
NL: NLBGS 9 ; 2.41306498e-09 3.12369805e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 2017.48797 1
NL: NLBGS 2 ; 1981.26916 0.98204757
NL: NLBGS 3 ; 139.420836 0.0691061544
NL: NLBGS 4 ; 11.5023529 0.00570132416
NL: NLBGS 5 ; 0.0461678481 2.28838281e-05
NL: NLBGS 6 ; 0.000751236633 3.72362385e-07
NL: NLBGS 7 ; 5.39772018e-05 2.67546585e-08
NL: NLBGS 8 ; 1.14137391e-07 5.65740132e-11
NL: NLBGS 9 ; 1.79325304e-09 8.88854387e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 2315.02802 1
NL: NLBGS 2 ; 1947.17547 0.841102334
NL: NLBGS 3 ; 129.690426 0.0560211041
NL: NLBGS 4 ; 12.9194059 0.00558066936
NL: NLBGS 5 ; 0.084312961 3.64198448e-05
NL: NLBGS 6 ; 0.00106888305 4.61714952e-07
NL: NLBGS 7 ; 8.63545087e-05 3.73017122e-08
NL: NLBGS 8 ; 8.47044777e-08 3.65889644e-11
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 182.356897 1
NL: NLBGS 2 ; 180.701446 0.990921916
NL: NLBGS 3 ; 11.6449018 0.0638577534
NL: NLBGS 4 ; 1.02824367 0.00563863327
NL: NLBGS 5 ; 0.00408511219 2.24017421e-05
NL: NLBGS 6 ; 2.74852894e-05 1.50722511e-07
NL: NLBGS 7 ; 2.21944067e-06 1.21708622e-08
NL: NLBGS 8 ; 4.34092368e-09 2.38045489e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 194.682609 1
NL: NLBGS 2 ; 167.911289 0.862487358
NL: NLBGS 3 ; 12.7173453 0.065323479
NL: NLBGS 4 ; 1.35751281 0.00697295363
NL: NLBGS 5 ; 0.00813143779 4.1767664e-05
NL: NLBGS 6 ; 0.000186070551 9.55763601e-07
NL: NLBGS 7 ; 1.61412033e-05 8.29103504e-08
NL: NLBGS 8 ; 4.3112499e-08 2.21450181e-10
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 1343.78717 1
NL: NLBGS 2 ; 1347.29176 1.002608
NL: NLBGS 3 ; 91.1709384 0.0678462637
NL: NLBGS 4 ; 7.68169212 0.00571644998
NL: NLBGS 5 ; 0.0304621912 2.26689107e-05
NL: NLBGS 6 ; 0.000360842003 2.68526157e-07
NL: NLBGS 7 ; 2.70284734e-05 2.01136565e-08
NL: NLBGS 8 ; 5.76507609e-08 4.29017052e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 1370.66385 1
NL: NLBGS 2 ; 1272.96048 0.928718213
NL: NLBGS 3 ; 92.4206633 0.0674276654
NL: NLBGS 4 ; 9.40435711 0.00686116959
NL: NLBGS 5 ; 0.0598106248 4.36362459e-05
NL: NLBGS 6 ; 0.000975568836 7.11749155e-07
NL: NLBGS 7 ; 8.11386081e-05 5.91965771e-08
NL: NLBGS 8 ; 1.65356052e-07 1.20639391e-10
NL: NLBGS 9 ; 8.54907231e-10 6.23717648e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 4796.16569 1
NL: NLBGS 2 ; 4815.89788 1.00411416
NL: NLBGS 3 ; 322.713386 0.0672857042
NL: NLBGS 4 ; 26.2718184 0.00547767116
NL: NLBGS 5 ; 0.0997039686 2.07882661e-05
NL: NLBGS 6 ; 0.00118106394 2.46251698e-07
NL: NLBGS 7 ; 8.5493814e-05 1.78254505e-08
NL: NLBGS 8 ; 1.86108046e-07 3.88035064e-11
NL: NLBGS 9 ; 2.16961447e-09 4.52364371e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 4833.53316 1
NL: NLBGS 2 ; 4509.32209 0.93292462
NL: NLBGS 3 ; 333.284798 0.0689526248
NL: NLBGS 4 ; 32.828267 0.00679177443
NL: NLBGS 5 ; 0.21309467 4.40867298e-05
NL: NLBGS 6 ; 0.00304251288 6.294594e-07
NL: NLBGS 7 ; 0.000246641359 5.10271371e-08
NL: NLBGS 8 ; 6.07710133e-07 1.25727933e-10
NL: NLBGS 9 ; 2.66538342e-09 5.5143584e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 3465.06458 1
NL: NLBGS 2 ; 3480.15596 1.0043553
NL: NLBGS 3 ; 234.283573 0.0676130464
NL: NLBGS 4 ; 19.5143219 0.00563173395
NL: NLBGS 5 ; 0.0720890735 2.08045397e-05
NL: NLBGS 6 ; 0.000902271391 2.60390931e-07
NL: NLBGS 7 ; 6.67075723e-05 1.92514658e-08
NL: NLBGS 8 ; 1.50432111e-07 4.34139415e-11
NL: NLBGS 9 ; 1.86640411e-09 5.38634725e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 3491.04651 1
NL: NLBGS 2 ; 3258.37708 0.933352525
NL: NLBGS 3 ; 241.553227 0.069192211
NL: NLBGS 4 ; 24.2372912 0.00694270074
NL: NLBGS 5 ; 0.146059838 4.18384106e-05
NL: NLBGS 6 ; 0.00231904887 6.64284723e-07
NL: NLBGS 7 ; 0.000189302484 5.42251393e-08
NL: NLBGS 8 ; 3.97976935e-07 1.13999322e-10
NL: NLBGS 9 ; 2.12332889e-09 6.08221312e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 662.11056 1
NL: NLBGS 2 ; 663.253152 1.00172568
NL: NLBGS 3 ; 44.7672993 0.0676130271
NL: NLBGS 4 ; 3.70985422 0.00560307363
NL: NLBGS 5 ; 0.0144471893 2.18199047e-05
NL: NLBGS 6 ; 0.000182792571 2.76075601e-07
NL: NLBGS 7 ; 1.34092446e-05 2.02522742e-08
NL: NLBGS 8 ; 2.87120403e-08 4.33644199e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 692.351233 1
NL: NLBGS 2 ; 639.808153 0.924109213
NL: NLBGS 3 ; 45.5138587 0.0657381059
NL: NLBGS 4 ; 4.40617384 0.00636407308
NL: NLBGS 5 ; 0.0279964497 4.04367731e-05
NL: NLBGS 6 ; 0.000310347281 4.48251215e-07
NL: NLBGS 7 ; 2.44746078e-05 3.5349988e-08
NL: NLBGS 8 ; 2.33769846e-08 3.37646321e-11
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 651.433414 1
NL: NLBGS 2 ; 655.374592 1.00605001
NL: NLBGS 3 ; 43.0208036 0.0660402164
NL: NLBGS 4 ; 3.54540835 0.00544247235
NL: NLBGS 5 ; 0.0128300647 1.96951284e-05
NL: NLBGS 6 ; 0.000147865682 2.26985105e-07
NL: NLBGS 7 ; 1.09060663e-05 1.6741644e-08
NL: NLBGS 8 ; 2.31624648e-08 3.55561509e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 631.057486 1
NL: NLBGS 2 ; 594.936795 0.942761648
NL: NLBGS 3 ; 44.5705552 0.0706283598
NL: NLBGS 4 ; 4.60981844 0.00730491048
NL: NLBGS 5 ; 0.0273817785 4.33903077e-05
NL: NLBGS 6 ; 0.000567852174 8.99842227e-07
NL: NLBGS 7 ; 4.80210952e-05 7.60962294e-08
NL: NLBGS 8 ; 1.42701966e-07 2.26131485e-10
NL: NLBGS 9 ; 7.47900754e-10 1.18515471e-12
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 2947.8912 1
NL: NLBGS 2 ; 2967.78558 1.00674868
NL: NLBGS 3 ; 186.960873 0.0634219042
NL: NLBGS 4 ; 14.5094538 0.00492197738
NL: NLBGS 5 ; 0.0498259355 1.69022302e-05
NL: NLBGS 6 ; 0.000564755775 1.91579586e-07
NL: NLBGS 7 ; 3.93348026e-05 1.33433699e-08
NL: NLBGS 8 ; 8.36832878e-08 2.8387509e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 2837.07981 1
NL: NLBGS 2 ; 2680.51337 0.944814228
NL: NLBGS 3 ; 195.05891 0.0687534094
NL: NLBGS 4 ; 19.0016168 0.00669759686
NL: NLBGS 5 ; 0.112361827 3.96047468e-05
NL: NLBGS 6 ; 0.00202122504 7.12431506e-07
NL: NLBGS 7 ; 0.000162759992 5.73688451e-08
NL: NLBGS 8 ; 5.01178953e-07 1.76653103e-10
NL: NLBGS 9 ; 2.56061251e-09 9.02552158e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 6498.02001 1
NL: NLBGS 2 ; 6555.96558 1.00891742
NL: NLBGS 3 ; 376.600411 0.0579561789
NL: NLBGS 4 ; 25.529461 0.00392880614
NL: NLBGS 5 ; 0.0772799339 1.18928433e-05
NL: NLBGS 6 ; 0.000856878032 1.31867558e-07
NL: NLBGS 7 ; 5.24697736e-05 8.07473253e-09
NL: NLBGS 8 ; 1.21891205e-07 1.8758207e-11
NL: NLBGS 9 ; 1.15609594e-09 1.77915109e-13
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 6112.88667 1
NL: NLBGS 2 ; 5822.07202 0.952425971
NL: NLBGS 3 ; 408.110629 0.0667623417
NL: NLBGS 4 ; 34.4532279 0.00563616337
NL: NLBGS 5 ; 0.197189588 3.22580146e-05
NL: NLBGS 6 ; 0.00244910071 4.00645528e-07
NL: NLBGS 7 ; 0.000174603643 2.85632064e-08
NL: NLBGS 8 ; 5.75131718e-07 9.40851268e-11
NL: NLBGS 9 ; 2.6876007e-09 4.39661463e-13
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 2049.18654 1
NL: NLBGS 2 ; 2069.74055 1.01003033
NL: NLBGS 3 ; 113.319699 0.0552998457
NL: NLBGS 4 ; 7.25235685 0.00353913942
NL: NLBGS 5 ; 0.0187400119 9.14509809e-06
NL: NLBGS 6 ; 0.000221917413 1.0829537e-07
NL: NLBGS 7 ; 1.28778168e-05 6.28435556e-09
NL: NLBGS 8 ; 3.44444643e-08 1.68088477e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 1869.07026 1
NL: NLBGS 2 ; 1770.62527 0.94732943
NL: NLBGS 3 ; 138.137521 0.0739070775
NL: NLBGS 4 ; 10.6039355 0.00567337445
NL: NLBGS 5 ; 0.0540334817 2.89092833e-05
NL: NLBGS 6 ; 0.000416869401 2.23035704e-07
NL: NLBGS 7 ; 2.70839805e-05 1.44906166e-08
NL: NLBGS 8 ; 8.98826069e-08 4.80894746e-11
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 1504.59529 1
NL: NLBGS 2 ; 1522.02926 1.01158715
NL: NLBGS 3 ; 81.299914 0.0540344066
NL: NLBGS 4 ; 5.04625033 0.00335389214
NL: NLBGS 5 ; 0.0128464693 8.53815597e-06
NL: NLBGS 6 ; 0.000147445086 9.79965088e-08
NL: NLBGS 7 ; 8.31852713e-06 5.52874728e-09
NL: NLBGS 8 ; 2.3431675e-08 1.55734071e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 1372.22594 1
NL: NLBGS 2 ; 1315.73474 0.958832439
NL: NLBGS 3 ; 94.7558345 0.069052648
NL: NLBGS 4 ; 7.10814432 0.00518001018
NL: NLBGS 5 ; 0.03706181 2.70085334e-05
NL: NLBGS 6 ; 0.000297080581 2.16495384e-07
NL: NLBGS 7 ; 1.88762542e-05 1.37559375e-08
NL: NLBGS 8 ; 6.43310138e-08 4.68807739e-11
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 744.092338 1
NL: NLBGS 2 ; 752.886104 1.01181811
NL: NLBGS 3 ; 39.6927074 0.0533437926
NL: NLBGS 4 ; 2.42026763 0.0032526442
NL: NLBGS 5 ; 0.00597073854 8.02419033e-06
NL: NLBGS 6 ; 6.95658397e-05 9.34908696e-08
NL: NLBGS 7 ; 3.85729484e-06 5.18389271e-09
NL: NLBGS 8 ; 1.11818993e-08 1.50275695e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 676.253092 1
NL: NLBGS 2 ; 647.979011 0.958190091
NL: NLBGS 3 ; 47.0007497 0.0695017151
NL: NLBGS 4 ; 3.42083727 0.00505851628
NL: NLBGS 5 ; 0.0176157993 2.60491219e-05
NL: NLBGS 6 ; 0.000119428947 1.7660392e-07
NL: NLBGS 7 ; 7.32626391e-06 1.08336124e-08
NL: NLBGS 8 ; 2.47355757e-08 3.65773938e-11
NL: NLBGS Converged

==================
AS_point_0.coupled
==================
NL: NLBGS 1 ; 25.0985467 1
NL: NLBGS 2 ; 25.4254066 1.01302306
NL: NLBGS 3 ; 1.28290816 0.0511148384
NL: NLBGS 4 ; 0.0780752779 0.00311074895
NL: NLBGS 5 ; 0.000168739364 6.72307308e-06
NL: NLBGS 6 ; 2.04724949e-06 8.15684476e-08
NL: NLBGS 7 ; 1.13754415e-07 4.53231081e-09
NL: NLBGS 8 ; 3.17113056e-10 1.26347178e-11
NL: NLBGS Converged

==================
AS_point_1.coupled
==================
NL: NLBGS 1 ; 22.9694892 1
NL: NLBGS 2 ; 21.832609 0.950504769
NL: NLBGS 3 ; 1.63563168 0.0712088834
NL: NLBGS 4 ; 0.105752815 0.00460405601
NL: NLBGS 5 ; 0.000529192185 2.30389183e-05
NL: NLBGS 6 ; 2.86582787e-06 1.24766722e-07
NL: NLBGS 7 ; 1.64082588e-07 7.14350182e-09
NL: NLBGS 8 ; 5.76755837e-10 2.51096501e-11
NL: NLBGS Converged
/home/docs/checkouts/readthedocs.org/user_builds/mdolab-openaerostruct/envs/latest/lib/python3.7/site-packages/openmdao/core/system.py:1961: PromotionWarning:'AS_point_0.coupled.wing.struct_states' <class SpatialBeamStates>: input variable 'nodes', promoted using 'nodes', was already promoted using 'nodes'.
/home/docs/checkouts/readthedocs.org/user_builds/mdolab-openaerostruct/envs/latest/lib/python3.7/site-packages/openmdao/core/system.py:1961: PromotionWarning:'AS_point_0.coupled.wing.struct_states' <class SpatialBeamStates>: input variable 'load_factor', promoted using 'load_factor', was already promoted using 'load_factor'.
/home/docs/checkouts/readthedocs.org/user_builds/mdolab-openaerostruct/envs/latest/lib/python3.7/site-packages/openmdao/core/system.py:1961: PromotionWarning:'AS_point_1.coupled.wing.struct_states' <class SpatialBeamStates>: input variable 'nodes', promoted using 'nodes', was already promoted using 'nodes'.
/home/docs/checkouts/readthedocs.org/user_builds/mdolab-openaerostruct/envs/latest/lib/python3.7/site-packages/openmdao/core/system.py:1961: PromotionWarning:'AS_point_1.coupled.wing.struct_states' <class SpatialBeamStates>: input variable 'load_factor', promoted using 'load_factor', was already promoted using 'load_factor'.
Optimization terminated successfully    (Exit mode 0)
            Current function value: [0.02670353]
            Iterations: 16
            Function evaluations: 18
            Gradient evaluations: 16
Optimization Complete
-----------------------------------
# prob.run_model()

print("The fuel burn value is", prob["AS_point_0.fuelburn"][0], "[kg]")
The fuel burn value is 2670.3528978501763 [kg]
print(
    "The wingbox mass (excluding the wing_weight_ratio) is",
    prob["wing.structural_mass"][0] / surf_dict["wing_weight_ratio"],
    "[kg]",
)
The wingbox mass (excluding the wing_weight_ratio) is 1123.6042237756928 [kg]

There is plenty of room for improvement. A finer mesh and a tighter optimization tolerance should be used.