MDF-based DOE on the Sobieski SSBJ test case

from os import name as os_name

from gemseo.api import configure_logger
from gemseo.api import create_discipline
from gemseo.api import create_scenario
from gemseo.problems.sobieski.core.problem import SobieskiProblem
from matplotlib import pyplot as plt

configure_logger()

Out:

<RootLogger root (INFO)>

Instantiate the disciplines

First, we instantiate the four disciplines of the use case: SobieskiPropulsion, SobieskiAerodynamics, SobieskiMission and SobieskiStructure.

disciplines = create_discipline(
    [
        "SobieskiPropulsion",
        "SobieskiAerodynamics",
        "SobieskiMission",
        "SobieskiStructure",
    ]
)

Build, execute and post-process the scenario

Then, we build the scenario which links the disciplines with the formulation and the optimization algorithm. Here, we use the BiLevel formulation. We tell the scenario to minimize -y_4 instead of minimizing y_4 (range), which is the default option.

We need to define the design space.

design_space = SobieskiProblem().design_space

Instantiate the scenario

scenario = create_scenario(
    disciplines,
    formulation="MDF",
    objective_name="y_4",
    design_space=design_space,
    maximize_objective=True,
    scenario_type="DOE",
)

Set the design constraints

for constraint in ["g_1", "g_2", "g_3"]:
    scenario.add_constraint(constraint, "ineq")

Execute the scenario

Use provided analytic derivatives

scenario.set_differentiation_method("user")

Multiprocessing

It is possible to run a DOE in parallel using multiprocessing, in order to do this, we specify the number of processes to be used for the computation of the samples.

#    The multiprocessing option has some limitations on Windows.
#    Due to problems with sphinx, we disable it in this example.
#    For Python versions < 3.7 and Numpy < 1.20.0, subprocesses may get hung
#    randomly during execution. It is strongly recommended to update your
#    environment to avoid this problem.
#    The features :class:`.MemoryFullCache` and :class:`.HDF5Cache` are not
#    available for multiprocessing on Windows.
#    As an alternative, we recommend the method
#    :meth:`.DOEScenario.set_optimization_history_backup`.

We define the algorithm options. Here the criterion = center option of pyDOE centers the points within the sampling intervals.

algo_options = {
    "criterion": "center",
    # Evaluate gradient of the MDA
    # with coupled adjoint
    "eval_jac": True,
    # Run in parallel on 1 or 4 processors
    "n_processes": 1 if os_name == "nt" else 4,
}

scenario.execute({"n_samples": 30, "algo": "lhs", "algo_options": algo_options})

Out:

    INFO - 07:17:46:
    INFO - 07:17:46: *** Start DOEScenario execution ***
    INFO - 07:17:46: DOEScenario
    INFO - 07:17:46:    Disciplines: SobieskiPropulsion SobieskiAerodynamics SobieskiMission SobieskiStructure
    INFO - 07:17:46:    MDO formulation: MDF
    INFO - 07:17:46: Optimization problem:
    INFO - 07:17:46:    minimize -y_4(x_shared, x_1, x_2, x_3)
    INFO - 07:17:46:    with respect to x_1, x_2, x_3, x_shared
    INFO - 07:17:46:    subject to constraints:
    INFO - 07:17:46:       g_1(x_shared, x_1, x_2, x_3) <= 0.0
    INFO - 07:17:46:       g_2(x_shared, x_1, x_2, x_3) <= 0.0
    INFO - 07:17:46:       g_3(x_shared, x_1, x_2, x_3) <= 0.0
    INFO - 07:17:46:    over the design space:
    INFO - 07:17:46:    +----------+-------------+-------+-------------+-------+
    INFO - 07:17:46:    | name     | lower_bound | value | upper_bound | type  |
    INFO - 07:17:46:    +----------+-------------+-------+-------------+-------+
    INFO - 07:17:46:    | x_shared |     0.01    |  0.05 |     0.09    | float |
    INFO - 07:17:46:    | x_shared |    30000    | 45000 |    60000    | float |
    INFO - 07:17:46:    | x_shared |     1.4     |  1.6  |     1.8     | float |
    INFO - 07:17:46:    | x_shared |     2.5     |  5.5  |     8.5     | float |
    INFO - 07:17:46:    | x_shared |      40     |   55  |      70     | float |
    INFO - 07:17:46:    | x_shared |     500     |  1000 |     1500    | float |
    INFO - 07:17:46:    | x_1      |     0.1     |  0.25 |     0.4     | float |
    INFO - 07:17:46:    | x_1      |     0.75    |   1   |     1.25    | float |
    INFO - 07:17:46:    | x_2      |     0.75    |   1   |     1.25    | float |
    INFO - 07:17:46:    | x_3      |     0.1     |  0.5  |      1      | float |
    INFO - 07:17:46:    +----------+-------------+-------+-------------+-------+
    INFO - 07:17:46: Solving optimization problem with algorithm lhs:
    INFO - 07:17:46: ...   0%|          | 0/30 [00:00<?, ?it]
    INFO - 07:17:46: Running DOE in parallel on n_processes = 4
    INFO - 07:17:46: ...   3%|▎         | 1/30 [00:00<00:00, 47.40 it/sec, obj=-567]
    INFO - 07:17:46: ...   7%|▋         | 2/30 [00:00<00:00, 39.51 it/sec, obj=-222]
    INFO - 07:17:47: ...  17%|█▋        | 5/30 [00:01<00:00, 25.76 it/sec, obj=-519]
    INFO - 07:17:47: ...  30%|███       | 9/30 [00:01<00:01, 18.96 it/sec, obj=-556]
    INFO - 07:17:47: ...  40%|████      | 12/30 [00:01<00:01, 16.67 it/sec, obj=-538]
    INFO - 07:17:47: ...  43%|████▎     | 13/30 [00:01<00:01, 15.35 it/sec, obj=-433]
    INFO - 07:17:48: ...  50%|█████     | 15/30 [00:02<00:01, 14.11 it/sec, obj=-574]
    INFO - 07:17:48: ...  57%|█████▋    | 17/30 [00:02<00:01, 12.45 it/sec, obj=-1.07e+3]
    INFO - 07:17:48: ...  63%|██████▎   | 19/30 [00:02<00:00, 11.52 it/sec, obj=-287]
    INFO - 07:17:48: ...  70%|███████   | 21/30 [00:02<00:00, 10.67 it/sec, obj=-247]
    INFO - 07:17:49: ...  73%|███████▎  | 22/30 [00:02<00:00, 10.07 it/sec, obj=-414]
    INFO - 07:17:49: ...  80%|████████  | 24/30 [00:03<00:00,  9.57 it/sec, obj=-273]
    INFO - 07:17:49: ...  87%|████████▋ | 26/30 [00:03<00:00,  8.73 it/sec, obj=-383]
    INFO - 07:17:49: ...  90%|█████████ | 27/30 [00:03<00:00,  8.46 it/sec, obj=-606]
    INFO - 07:17:49: ...  97%|█████████▋| 29/30 [00:03<00:00,  8.18 it/sec, obj=-485]
    INFO - 07:17:49: ... 100%|██████████| 30/30 [00:03<00:00,  7.99 it/sec, obj=-405]
    INFO - 07:17:49: Optimization result:
    INFO - 07:17:49:    Optimizer info:
    INFO - 07:17:49:       Status: None
    INFO - 07:17:49:       Message: None
    INFO - 07:17:49:       Number of calls to the objective function by the optimizer: 30
    INFO - 07:17:49:    Solution:
    INFO - 07:17:49:       The solution is feasible.
    INFO - 07:17:49:       Objective: -485.49229968733584
    INFO - 07:17:49:       Standardized constraints:
    INFO - 07:17:49:          g_1 = [-0.11350951 -0.10812292 -0.1045109  -0.10204971 -0.10028641 -0.01838903
    INFO - 07:17:49:  -0.22161097]
    INFO - 07:17:49:          g_2 = -0.02400000000000002
    INFO - 07:17:49:          g_3 = [-0.33063134 -0.66936866 -0.73821755 -0.07789536]
    INFO - 07:17:49:       Design space:
    INFO - 07:17:49:       +----------+-------------+---------------------+-------------+-------+
    INFO - 07:17:49:       | name     | lower_bound |        value        | upper_bound | type  |
    INFO - 07:17:49:       +----------+-------------+---------------------+-------------+-------+
    INFO - 07:17:49:       | x_shared |     0.01    | 0.05400000000000001 |     0.09    | float |
    INFO - 07:17:49:       | x_shared |    30000    |        46500        |    60000    | float |
    INFO - 07:17:49:       | x_shared |     1.4     |  1.686666666666667  |     1.8     | float |
    INFO - 07:17:49:       | x_shared |     2.5     |         5.2         |     8.5     | float |
    INFO - 07:17:49:       | x_shared |      40     |         66.5        |      70     | float |
    INFO - 07:17:49:       | x_shared |     500     |  583.3333333333334  |     1500    | float |
    INFO - 07:17:49:       | x_1      |     0.1     |        0.185        |     0.4     | float |
    INFO - 07:17:49:       | x_1      |     0.75    |  0.9416666666666667 |     1.25    | float |
    INFO - 07:17:49:       | x_2      |     0.75    |        0.775        |     1.25    | float |
    INFO - 07:17:49:       | x_3      |     0.1     |        0.115        |      1      | float |
    INFO - 07:17:49:       +----------+-------------+---------------------+-------------+-------+
    INFO - 07:17:49: *** End DOEScenario execution (time: 0:00:03.768810) ***

{'eval_jac': False, 'algo': 'lhs', 'algo_options': {'criterion': 'center', 'eval_jac': True, 'n_processes': 4, 'seed': 1, 'n_samples': 30}, 'n_samples': 30}

Warning

On Windows, the progress bar may show duplicated instances during the initialization of each subprocess. In some cases it may also print the conclusion of an iteration ahead of another one that was concluded first. This is a consequence of the pickling process and does not affect the computations of the scenario.

Plot the optimization history view

scenario.post_process("OptHistoryView", show=False, save=False)

• 
• 
• 
• 
• 

Out:

<gemseo.post.opt_history_view.OptHistoryView object at 0x7f28d84cc160>

Tip

Each post-processing method requires different inputs and offers a variety of customization options. Use the API function get_post_processing_options_schema() to print a table with the attributes for any post-processing algo. Or refer to our dedicated page: Post-processing algorithms.

Plot the scatter matrix

scenario.post_process(
    "ScatterPlotMatrix", show=False, save=False, variable_names=["y_4", "x_shared"]
)

Out:

<gemseo.post.scatter_mat.ScatterPlotMatrix object at 0x7f28d724c370>

Plot correlations

scenario.post_process("Correlations", show=False, save=False)
# Workaround for HTML rendering, instead of ``show=True``
plt.show()

Out:

INFO - 07:17:52: Detected 10 correlations > 0.95