Note
Go to the end to download the full example code
MDF-based DOE on the Sobieski SSBJ test case¶
from __future__ import annotations
from os import name as os_name
from gemseo import configure_logger
from gemseo import create_discipline
from gemseo import create_scenario
from gemseo import generate_n2_plot
from gemseo.problems.sobieski.core.design_space import SobieskiDesignSpace
configure_logger()
<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",
])
We can quickly access the most relevant information of any discipline (name, inputs,
and outputs) with Python’s print() function. Moreover, we can get the default
input values of a discipline with the attribute MDODiscipline.default_inputs
for discipline in disciplines:
print(discipline)
print(f"Default inputs: {discipline.default_inputs}")
SobieskiPropulsion
Default inputs: {'x_3': array([0.5]), 'y_23': array([12562.01206488]), 'x_shared': array([5.0e-02, 4.5e+04, 1.6e+00, 5.5e+00, 5.5e+01, 1.0e+03]), 'c_3': array([4360.])}
SobieskiAerodynamics
Default inputs: {'c_4': array([0.01375]), 'y_12': array([5.06069742e+04, 9.50000000e-01]), 'y_32': array([0.50279625]), 'x_shared': array([5.0e-02, 4.5e+04, 1.6e+00, 5.5e+00, 5.5e+01, 1.0e+03]), 'x_2': array([1.])}
SobieskiMission
Default inputs: {'y_14': array([50606.9741711 , 7306.20262124]), 'y_24': array([4.15006276]), 'x_shared': array([5.0e-02, 4.5e+04, 1.6e+00, 5.5e+00, 5.5e+01, 1.0e+03]), 'y_34': array([1.10754577])}
SobieskiStructure
Default inputs: {'y_21': array([50606.9741711]), 'x_1': array([0.25, 1. ]), 'x_shared': array([5.0e-02, 4.5e+04, 1.6e+00, 5.5e+00, 5.5e+01, 1.0e+03]), 'c_0': array([2000.]), 'c_1': array([25000.]), 'c_2': array([6.]), 'y_31': array([6354.32430691])}
You may also be interested in plotting the couplings of your disciplines.
A quick way of getting this information is the high-level function
generate_n2_plot(). A much more detailed explanation of coupling
visualization is available here.
generate_n2_plot(disciplines, save=False, show=True)
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 = SobieskiDesignSpace()
design_space
Instantiate the scenario¶
scenario = create_scenario(
disciplines,
"MDF",
"y_4",
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, constraint_type="ineq")
Visualize the XDSM¶
Generate the XDSM on the fly:
log_workflow_status=Truewill log the status of the workflow in the console,save_html(defaultTrue) will generate a self-contained HTML file, that can be automatically opened usingshow_html=True.
scenario.xdsmize(save_html=False)
Execute the scenario¶
Use provided analytic derivatives
scenario.set_differentiation_method()
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.
Warning
The multiprocessing option has some limitations on Windows.
Due to problems with sphinx, we disable it in this example.
The features MemoryFullCache and HDF5Cache are not
available for multiprocessing on Windows.
As an alternative, we recommend the method
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})
INFO - 13:09:41:
INFO - 13:09:41: *** Start DOEScenario execution ***
INFO - 13:09:41: DOEScenario
INFO - 13:09:41: Disciplines: SobieskiAerodynamics SobieskiMission SobieskiPropulsion SobieskiStructure
INFO - 13:09:41: MDO formulation: MDF
INFO - 13:09:41: Optimization problem:
INFO - 13:09:41: minimize -y_4(x_shared, x_1, x_2, x_3)
INFO - 13:09:41: with respect to x_1, x_2, x_3, x_shared
INFO - 13:09:41: subject to constraints:
INFO - 13:09:41: g_1(x_shared, x_1, x_2, x_3) <= 0.0
INFO - 13:09:41: g_2(x_shared, x_1, x_2, x_3) <= 0.0
INFO - 13:09:41: g_3(x_shared, x_1, x_2, x_3) <= 0.0
INFO - 13:09:41: over the design space:
INFO - 13:09:41: +-------------+-------------+-------+-------------+-------+
INFO - 13:09:41: | Name | Lower bound | Value | Upper bound | Type |
INFO - 13:09:41: +-------------+-------------+-------+-------------+-------+
INFO - 13:09:41: | x_shared[0] | 0.01 | 0.05 | 0.09 | float |
INFO - 13:09:41: | x_shared[1] | 30000 | 45000 | 60000 | float |
INFO - 13:09:41: | x_shared[2] | 1.4 | 1.6 | 1.8 | float |
INFO - 13:09:41: | x_shared[3] | 2.5 | 5.5 | 8.5 | float |
INFO - 13:09:41: | x_shared[4] | 40 | 55 | 70 | float |
INFO - 13:09:41: | x_shared[5] | 500 | 1000 | 1500 | float |
INFO - 13:09:41: | x_1[0] | 0.1 | 0.25 | 0.4 | float |
INFO - 13:09:41: | x_1[1] | 0.75 | 1 | 1.25 | float |
INFO - 13:09:41: | x_2 | 0.75 | 1 | 1.25 | float |
INFO - 13:09:41: | x_3 | 0.1 | 0.5 | 1 | float |
INFO - 13:09:41: +-------------+-------------+-------+-------------+-------+
INFO - 13:09:41: Solving optimization problem with algorithm lhs:
INFO - 13:09:41: Running DOE in parallel on n_processes = 4
INFO - 13:09:42: 3%|▎ | 1/30 [00:00<00:28, 1.00 it/sec, obj=-285]
INFO - 13:09:42: 7%|▋ | 2/30 [00:01<00:14, 1.95 it/sec, obj=-222]
INFO - 13:09:42: 10%|█ | 3/30 [00:01<00:09, 2.88 it/sec, obj=-344]
INFO - 13:09:43: 13%|█▎ | 4/30 [00:01<00:07, 3.64 it/sec, obj=-567]
INFO - 13:09:43: 17%|█▋ | 5/30 [00:01<00:08, 2.98 it/sec, obj=-519]
INFO - 13:09:43: 20%|██ | 6/30 [00:01<00:06, 3.50 it/sec, obj=-481]
INFO - 13:09:43: 23%|██▎ | 7/30 [00:01<00:05, 3.93 it/sec, obj=-429]
INFO - 13:09:43: 27%|██▋ | 8/30 [00:01<00:04, 4.48 it/sec, obj=-235]
INFO - 13:09:44: 30%|███ | 9/30 [00:02<00:05, 3.75 it/sec, obj=-556]
INFO - 13:09:44: 33%|███▎ | 10/30 [00:02<00:04, 4.16 it/sec, obj=-306]
INFO - 13:09:44: 37%|███▋ | 11/30 [00:02<00:04, 4.53 it/sec, obj=-538]
INFO - 13:09:44: 40%|████ | 12/30 [00:02<00:03, 4.94 it/sec, obj=-419]
INFO - 13:09:44: 43%|████▎ | 13/30 [00:03<00:03, 4.30 it/sec, obj=-433]
INFO - 13:09:44: 47%|████▋ | 14/30 [00:03<00:03, 4.61 it/sec, obj=-495]
INFO - 13:09:44: 50%|█████ | 15/30 [00:03<00:03, 4.90 it/sec, obj=-574]
INFO - 13:09:45: 53%|█████▎ | 16/30 [00:03<00:02, 5.18 it/sec, obj=-621]
INFO - 13:09:45: 57%|█████▋ | 17/30 [00:03<00:02, 4.67 it/sec, obj=-481]
INFO - 13:09:45: 60%|██████ | 18/30 [00:03<00:02, 4.86 it/sec, obj=-1.07e+3]
INFO - 13:09:45: 63%|██████▎ | 19/30 [00:03<00:02, 5.13 it/sec, obj=-602]
INFO - 13:09:45: 67%|██████▋ | 20/30 [00:03<00:01, 5.34 it/sec, obj=-287]
INFO - 13:09:46: 70%|███████ | 21/30 [00:04<00:01, 4.94 it/sec, obj=-247]
INFO - 13:09:46: 73%|███████▎ | 22/30 [00:04<00:01, 5.17 it/sec, obj=-414]
INFO - 13:09:46: 77%|███████▋ | 23/30 [00:04<00:01, 5.29 it/sec, obj=-1.18e+3]
INFO - 13:09:46: 80%|████████ | 24/30 [00:04<00:01, 5.38 it/sec, obj=-273]
INFO - 13:09:46: 83%|████████▎ | 25/30 [00:04<00:00, 5.13 it/sec, obj=-383]
INFO - 13:09:46: 87%|████████▋ | 26/30 [00:04<00:00, 5.32 it/sec, obj=-624]
INFO - 13:09:46: 90%|█████████ | 27/30 [00:05<00:00, 5.35 it/sec, obj=-606]
INFO - 13:09:47: 93%|█████████▎| 28/30 [00:05<00:00, 5.46 it/sec, obj=-351]
INFO - 13:09:47: 97%|█████████▋| 29/30 [00:05<00:00, 5.43 it/sec, obj=-405]
INFO - 13:09:47: 100%|██████████| 30/30 [00:05<00:00, 5.58 it/sec, obj=-485]
INFO - 13:09:47: Optimization result:
INFO - 13:09:47: Optimizer info:
INFO - 13:09:47: Status: None
INFO - 13:09:47: Message: None
INFO - 13:09:47: Number of calls to the objective function by the optimizer: 30
INFO - 13:09:47: Solution:
INFO - 13:09:47: The solution is feasible.
INFO - 13:09:47: Objective: -485.49213288049987
INFO - 13:09:47: Standardized constraints:
INFO - 13:09:47: g_1 = [-0.11350951 -0.10812292 -0.1045109 -0.10204971 -0.10028641 -0.01838903
INFO - 13:09:47: -0.22161097]
INFO - 13:09:47: g_2 = -0.02400000000000002
INFO - 13:09:47: g_3 = [-0.33063169 -0.66936831 -0.73821755 -0.07789536]
INFO - 13:09:47: Design space:
INFO - 13:09:47: +-------------+-------------+---------------------+-------------+-------+
INFO - 13:09:47: | Name | Lower bound | Value | Upper bound | Type |
INFO - 13:09:47: +-------------+-------------+---------------------+-------------+-------+
INFO - 13:09:47: | x_shared[0] | 0.01 | 0.05400000000000001 | 0.09 | float |
INFO - 13:09:47: | x_shared[1] | 30000 | 46500 | 60000 | float |
INFO - 13:09:47: | x_shared[2] | 1.4 | 1.686666666666667 | 1.8 | float |
INFO - 13:09:47: | x_shared[3] | 2.5 | 5.2 | 8.5 | float |
INFO - 13:09:47: | x_shared[4] | 40 | 66.5 | 70 | float |
INFO - 13:09:47: | x_shared[5] | 500 | 583.3333333333334 | 1500 | float |
INFO - 13:09:47: | x_1[0] | 0.1 | 0.185 | 0.4 | float |
INFO - 13:09:47: | x_1[1] | 0.75 | 0.9416666666666667 | 1.25 | float |
INFO - 13:09:47: | x_2 | 0.75 | 0.775 | 1.25 | float |
INFO - 13:09:47: | x_3 | 0.1 | 0.115 | 1 | float |
INFO - 13:09:47: +-------------+-------------+---------------------+-------------+-------+
INFO - 13:09:47: *** End DOEScenario execution (time: 0:00:05.464534) ***
{'eval_jac': False, 'algo_options': {'criterion': 'center', 'eval_jac': True, 'n_processes': 4}, 'n_samples': 30, 'algo': 'lhs'}
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.
Exporting the problem data.¶
After the execution of the scenario, you may want to export your data to use it
elsewhere. The method Scenario.to_dataset() will allow you to export
your results to a Dataset, the basic GEMSEO class to store data.
dataset = scenario.to_dataset("a_name_for_my_dataset")
Plot the optimization history view¶
scenario.post_process("OptHistoryView", save=False, show=True)
<gemseo.post.opt_history_view.OptHistoryView object at 0x7f8b9f965ac0>
Tip
Each post-processing method requires different inputs and offers a variety
of customization options. Use the high-level 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",
save=False,
show=True,
variable_names=["y_4", "x_shared"],
)
<gemseo.post.scatter_mat.ScatterPlotMatrix object at 0x7f8b9c6c1490>
Plot correlations¶
scenario.post_process("Correlations", save=False, show=True)
INFO - 13:09:50: Detected 10 correlations > 0.95
<gemseo.post.correlations.Correlations object at 0x7f8b9d83e3a0>
Total running time of the script: (0 minutes 9.998 seconds)