Note

Go to the end to download the full example code.

Make a monodisciplinary optimization problem multidisciplinary#

Introduction#

The OptAsMDOScenario is a monodisciplinary optimization scenario made multidisciplinary. The only requirement is that the discipline has at least three scalar inputs and at least one output. This scenario can be used to enrich a catalog of benchmark MDO problems, based on the observation that MDO benchmark problems are far less numerous than optimization problems.

This example illustrates it in the case of the minimization of the 3-dimensional Rosenbrock function

\[f(z) = 100(z_2-z_1^2)^2 + (1-z_1)^2 + 100(z_1-z_0^2)^2 + (1-z_0)^2\]

over the hypercube \([-1,1]^3\).

The unique solution of this minimization problem is the design point \(z^*=(1,1,1)\) at which \(f\) is zero.

from __future__ import annotations

from numpy import array

from gemseo import configure_logger
from gemseo import create_design_space
from gemseo import create_discipline
from gemseo import create_scenario
from gemseo import generate_coupling_graph
from gemseo.problems.mdo.opt_as_mdo_scenario import OptAsMDOScenario

configure_logger()

<RootLogger root (INFO)>

Material#

First, we create the discipline implementing the Rosenbrock function:

discipline = create_discipline(
    "AnalyticDiscipline",
    expressions={"f": "100*(z_2-z_1**2)**2+(1-z_1)**2+100*(z_1-z_0**2)**2+(1-z_0)**2"},
    name="Rosenbrock",
)

as well as the design space:

design_space = create_design_space()
design_space.add_variable("z_0", lower_bound=-1, upper_bound=1)
design_space.add_variable("z_1", lower_bound=-1, upper_bound=1)
design_space.add_variable("z_2", lower_bound=-1, upper_bound=1)

and choose \(x^{(0)}=(-0.25, 0.75, -0.9)\) as the starting point of the optimization:

initial_point = array([-0.25, 0.75, -0.9])
design_space.set_current_value(initial_point)

Optimization problem#

Then, we define the optimization problem:

opt_scenario = create_scenario(
    [discipline], "f", design_space, formulation_name="DisciplinaryOpt"
)

and solve it using the SLSQP algorithm:

opt_scenario.execute(algo_name="NLOPT_SLSQP", max_iter=100)

INFO - 10:10:34: *** Start MDOScenario execution ***
INFO - 10:10:34: MDOScenario
INFO - 10:10:34:    Disciplines: Rosenbrock
INFO - 10:10:34:    MDO formulation: DisciplinaryOpt
INFO - 10:10:34: Optimization problem:
INFO - 10:10:34:    minimize f(z_0, z_1, z_2)
INFO - 10:10:34:    with respect to z_0, z_1, z_2
INFO - 10:10:34:    over the design space:
INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
INFO - 10:10:34:       | Name | Lower bound | Value | Upper bound | Type  |
INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
INFO - 10:10:34:       | z_0  |      -1     | -0.25 |      1      | float |
INFO - 10:10:34:       | z_1  |      -1     |  0.75 |      1      | float |
INFO - 10:10:34:       | z_2  |      -1     |  -0.9 |      1      | float |
INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
INFO - 10:10:34: Solving optimization problem with algorithm NLOPT_SLSQP:
INFO - 10:10:34:      1%|          | 1/100 [00:00<00:00, 292.76 it/sec, obj=263]
INFO - 10:10:34:      2%|▏         | 2/100 [00:00<00:00, 390.39 it/sec, obj=408]
INFO - 10:10:34:      3%|▎         | 3/100 [00:00<00:00, 550.77 it/sec, obj=20.2]
INFO - 10:10:34:      4%|▍         | 4/100 [00:00<00:00, 542.20 it/sec, obj=168]
INFO - 10:10:34:      5%|▌         | 5/100 [00:00<00:00, 654.42 it/sec, obj=8.91]
INFO - 10:10:34:      6%|▌         | 6/100 [00:00<00:00, 647.04 it/sec, obj=395]
INFO - 10:10:34:      7%|▋         | 7/100 [00:00<00:00, 732.34 it/sec, obj=4.06]
INFO - 10:10:34:      8%|▊         | 8/100 [00:00<00:00, 711.86 it/sec, obj=151]
INFO - 10:10:34:      9%|▉         | 9/100 [00:00<00:00, 778.50 it/sec, obj=3.31]
INFO - 10:10:34:     10%|█         | 10/100 [00:00<00:00, 752.36 it/sec, obj=2.09]
INFO - 10:10:34:     11%|█         | 11/100 [00:00<00:00, 775.04 it/sec, obj=2.02]
INFO - 10:10:34:     12%|█▏        | 12/100 [00:00<00:00, 798.46 it/sec, obj=1.97]
INFO - 10:10:34:     13%|█▎        | 13/100 [00:00<00:00, 817.98 it/sec, obj=1.92]
INFO - 10:10:34:     14%|█▍        | 14/100 [00:00<00:00, 838.25 it/sec, obj=1.85]
INFO - 10:10:34:     15%|█▌        | 15/100 [00:00<00:00, 857.19 it/sec, obj=1.73]
INFO - 10:10:34:     16%|█▌        | 16/100 [00:00<00:00, 876.39 it/sec, obj=1.3]
INFO - 10:10:34:     17%|█▋        | 17/100 [00:00<00:00, 890.40 it/sec, obj=0.995]
INFO - 10:10:34:     18%|█▊        | 18/100 [00:00<00:00, 901.57 it/sec, obj=0.721]
INFO - 10:10:34:     19%|█▉        | 19/100 [00:00<00:00, 914.74 it/sec, obj=0.625]
INFO - 10:10:34:     20%|██        | 20/100 [00:00<00:00, 928.11 it/sec, obj=6.19]
INFO - 10:10:34:     21%|██        | 21/100 [00:00<00:00, 966.67 it/sec, obj=0.491]
INFO - 10:10:34:     22%|██▏       | 22/100 [00:00<00:00, 942.66 it/sec, obj=5.86]
INFO - 10:10:34:     23%|██▎       | 23/100 [00:00<00:00, 975.91 it/sec, obj=0.403]
INFO - 10:10:34:     24%|██▍       | 24/100 [00:00<00:00, 954.17 it/sec, obj=0.87]
INFO - 10:10:34:     25%|██▌       | 25/100 [00:00<00:00, 984.47 it/sec, obj=0.279]
INFO - 10:10:34:     26%|██▌       | 26/100 [00:00<00:00, 961.39 it/sec, obj=0.243]
INFO - 10:10:34:     27%|██▋       | 27/100 [00:00<00:00, 970.33 it/sec, obj=0.244]
INFO - 10:10:34:     28%|██▊       | 28/100 [00:00<00:00, 999.61 it/sec, obj=0.224]
INFO - 10:10:34:     29%|██▉       | 29/100 [00:00<00:00, 977.15 it/sec, obj=0.171]
INFO - 10:10:34:     30%|███       | 30/100 [00:00<00:00, 985.29 it/sec, obj=0.276]
INFO - 10:10:34:     31%|███       | 31/100 [00:00<00:00, 1012.12 it/sec, obj=0.117]
INFO - 10:10:34:     32%|███▏      | 32/100 [00:00<00:00, 987.73 it/sec, obj=0.0792]
INFO - 10:10:34:     33%|███▎      | 33/100 [00:00<00:00, 991.70 it/sec, obj=0.0624]
INFO - 10:10:34:     34%|███▍      | 34/100 [00:00<00:00, 996.50 it/sec, obj=0.0366]
INFO - 10:10:34:     35%|███▌      | 35/100 [00:00<00:00, 1001.49 it/sec, obj=0.0298]
INFO - 10:10:34:     36%|███▌      | 36/100 [00:00<00:00, 1004.64 it/sec, obj=0.0164]
INFO - 10:10:34:     37%|███▋      | 37/100 [00:00<00:00, 1009.11 it/sec, obj=0.0103]
INFO - 10:10:34:     38%|███▊      | 38/100 [00:00<00:00, 1013.42 it/sec, obj=0.00429]
INFO - 10:10:34:     39%|███▉      | 39/100 [00:00<00:00, 1018.81 it/sec, obj=0.00128]
INFO - 10:10:34:     40%|████      | 40/100 [00:00<00:00, 1023.66 it/sec, obj=0.000368]
INFO - 10:10:34:     41%|████      | 41/100 [00:00<00:00, 1026.94 it/sec, obj=2.77e-5]
INFO - 10:10:34:     42%|████▏     | 42/100 [00:00<00:00, 1031.25 it/sec, obj=6.38e-7]
INFO - 10:10:34:     43%|████▎     | 43/100 [00:00<00:00, 1034.00 it/sec, obj=1.32e-9]
INFO - 10:10:34:     44%|████▍     | 44/100 [00:00<00:00, 1037.72 it/sec, obj=9.1e-12]
INFO - 10:10:34:     45%|████▌     | 45/100 [00:00<00:00, 1042.56 it/sec, obj=6.58e-14]
INFO - 10:10:34:     46%|████▌     | 46/100 [00:00<00:00, 1045.56 it/sec, obj=1.58e-21]
INFO - 10:10:34:     47%|████▋     | 47/100 [00:00<00:00, 1050.72 it/sec, obj=Not evaluated]
INFO - 10:10:34: Optimization result:
INFO - 10:10:34:    Optimizer info:
INFO - 10:10:34:       Status: None
INFO - 10:10:34:       Message: Successive iterates of the design variables are closer than xtol_rel or xtol_abs. GEMSEO stopped the driver.
INFO - 10:10:34:       Number of calls to the objective function by the optimizer: 56
INFO - 10:10:34:    Solution:
INFO - 10:10:34:       Objective: 1.584054020051718e-21
INFO - 10:10:34:       Design space:
INFO - 10:10:34:          +------+-------------+--------------------+-------------+-------+
INFO - 10:10:34:          | Name | Lower bound |       Value        | Upper bound | Type  |
INFO - 10:10:34:          +------+-------------+--------------------+-------------+-------+
INFO - 10:10:34:          | z_0  |      -1     |         1          |      1      | float |
INFO - 10:10:34:          | z_1  |      -1     | 0.9999999999996698 |      1      | float |
INFO - 10:10:34:          | z_2  |      -1     | 0.9999999999953735 |      1      | float |
INFO - 10:10:34:          +------+-------------+--------------------+-------------+-------+
INFO - 10:10:34: *** End MDOScenario execution (time: 0:00:00.047760) ***

We can see that the numerical solution corresponds to the analytical one.

MDO problem#

Now, we use the OptAsMDOScenario to rewrite this optimization problem as an MDO problem with two strongly coupled disciplines.

First, we reset the design space to the initial point:

design_space.set_current_value(initial_point)

and create the OptAsMDOScenario, orchestrated by an MDF formulation:

mdo_scenario = OptAsMDOScenario(discipline, "f", design_space, formulation_name="MDF")

Then, we can see that the design variables have been renamed:

design_space

Design space:

Name	Lower bound	Value	Upper bound	Type
x_0	-1	-0.25	1	float
x_1	-1	0.75	1	float
x_2	-1	-0.9	1	float

This renaming is based on the convention:

the first design variable is the global design variable and is named \(x_0\),
the \((1+i)\)-th design variable is the local design variable specific to the \(i\)-th strongly coupled discipline and is named \(x_{1+i}\).

We can also have a look to the coupling graph:

generate_coupling_graph(mdo_scenario.disciplines, file_path="")

and see that there are two strongly coupled disciplines \(D_1\) and \(D_2\), connected by the coupling variables \(y_1\) and \(y_2\). These disciplines are weakly coupled to a downstream link discipline \(L\), which is weakly coupled to the downstream original discipline. Let us note that the link discipline computes the values of the design variables in the original optimization problem from the values of the design and coupling variables in the MDO problem.

Lastly, we solve this scenario using the SLSQP algorithm:

mdo_scenario.set_differentiation_method(method="finite_differences")
mdo_scenario.execute(algo_name="NLOPT_SLSQP", max_iter=100)

   INFO - 10:10:34: *** Start OptAsMDOScenario execution ***
   INFO - 10:10:34: OptAsMDOScenario
   INFO - 10:10:34:    Disciplines: D1 D2 L Rosenbrock
   INFO - 10:10:34:    MDO formulation: MDF
   INFO - 10:10:34: Optimization problem:
   INFO - 10:10:34:    minimize f(x_0, x_1, x_2)
   INFO - 10:10:34:    with respect to x_0, x_1, x_2
   INFO - 10:10:34:    over the design space:
   INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
   INFO - 10:10:34:       | Name | Lower bound | Value | Upper bound | Type  |
   INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
   INFO - 10:10:34:       | x_0  |      -1     | -0.25 |      1      | float |
   INFO - 10:10:34:       | x_1  |      -1     |  0.75 |      1      | float |
   INFO - 10:10:34:       | x_2  |      -1     |  -0.9 |      1      | float |
   INFO - 10:10:34:       +------+-------------+-------+-------------+-------+
   INFO - 10:10:34: Solving optimization problem with algorithm NLOPT_SLSQP:
   INFO - 10:10:34:      1%|          | 1/100 [00:00<00:01, 67.29 it/sec, obj=263]
   INFO - 10:10:34:      2%|▏         | 2/100 [00:00<00:02, 36.42 it/sec, obj=408]
   INFO - 10:10:34:      3%|▎         | 3/100 [00:00<00:01, 54.33 it/sec, obj=20.2]
WARNING - 10:10:34: All components of the normalized vector should be between 0 and 1; upper bounds violated: [1.000001].
WARNING - 10:10:34: All components of the normalized vector should be between 0 and 1; upper bounds violated: [1.000001].
   INFO - 10:10:34:      4%|▍         | 4/100 [00:00<00:02, 42.05 it/sec, obj=168]
   INFO - 10:10:34:      5%|▌         | 5/100 [00:00<00:01, 52.40 it/sec, obj=8.91]
   INFO - 10:10:34:      6%|▌         | 6/100 [00:00<00:02, 44.79 it/sec, obj=395]
   INFO - 10:10:34:      7%|▋         | 7/100 [00:00<00:01, 52.14 it/sec, obj=4.06]
WARNING - 10:10:34: All components of the normalized vector should be between 0 and 1; upper bounds violated: [1.000001].
   INFO - 10:10:34:      8%|▊         | 8/100 [00:00<00:01, 46.50 it/sec, obj=151]
   INFO - 10:10:34:      9%|▉         | 9/100 [00:00<00:01, 52.23 it/sec, obj=3.31]
   INFO - 10:10:34:     10%|█         | 10/100 [00:00<00:02, 44.80 it/sec, obj=2.09]
   INFO - 10:10:34:     11%|█         | 11/100 [00:00<00:01, 44.90 it/sec, obj=2.02]
   INFO - 10:10:34:     12%|█▏        | 12/100 [00:00<00:01, 45.06 it/sec, obj=1.97]
   INFO - 10:10:34:     13%|█▎        | 13/100 [00:00<00:01, 45.06 it/sec, obj=1.92]
   INFO - 10:10:35:     14%|█▍        | 14/100 [00:00<00:01, 45.12 it/sec, obj=1.85]
   INFO - 10:10:35:     15%|█▌        | 15/100 [00:00<00:01, 45.22 it/sec, obj=1.73]
   INFO - 10:10:35:     16%|█▌        | 16/100 [00:00<00:01, 45.29 it/sec, obj=1.3]
   INFO - 10:10:35:     17%|█▋        | 17/100 [00:00<00:01, 45.41 it/sec, obj=0.996]
   INFO - 10:10:35:     18%|█▊        | 18/100 [00:00<00:01, 45.47 it/sec, obj=0.721]
WARNING - 10:10:35: All components of the normalized vector should be between 0 and 1; upper bounds violated: [1.000001].
   INFO - 10:10:35:     19%|█▉        | 19/100 [00:00<00:01, 45.49 it/sec, obj=0.6]
   INFO - 10:10:35:     20%|██        | 20/100 [00:00<00:01, 46.73 it/sec, obj=6.26]
   INFO - 10:10:35:     21%|██        | 21/100 [00:00<00:01, 49.03 it/sec, obj=0.473]
   INFO - 10:10:35:     22%|██▏       | 22/100 [00:00<00:01, 47.08 it/sec, obj=6.1]
   INFO - 10:10:35:     23%|██▎       | 23/100 [00:00<00:01, 49.18 it/sec, obj=0.394]
   INFO - 10:10:35:     24%|██▍       | 24/100 [00:00<00:01, 47.24 it/sec, obj=0.665]
   INFO - 10:10:35:     25%|██▌       | 25/100 [00:00<00:01, 49.17 it/sec, obj=0.28]
   INFO - 10:10:35:     26%|██▌       | 26/100 [00:00<00:01, 46.34 it/sec, obj=0.242]
   INFO - 10:10:35:     27%|██▋       | 27/100 [00:00<00:01, 46.33 it/sec, obj=0.207]
   INFO - 10:10:35:     28%|██▊       | 28/100 [00:00<00:01, 46.33 it/sec, obj=0.142]
   INFO - 10:10:35:     29%|██▉       | 29/100 [00:00<00:01, 46.35 it/sec, obj=0.091]
   INFO - 10:10:35:     30%|███       | 30/100 [00:00<00:01, 46.37 it/sec, obj=0.0673]
   INFO - 10:10:35:     31%|███       | 31/100 [00:00<00:01, 46.37 it/sec, obj=0.0512]
   INFO - 10:10:35:     32%|███▏      | 32/100 [00:00<00:01, 46.42 it/sec, obj=0.029]
   INFO - 10:10:35:     33%|███▎      | 33/100 [00:00<00:01, 46.41 it/sec, obj=0.0164]
   INFO - 10:10:35:     34%|███▍      | 34/100 [00:00<00:01, 46.48 it/sec, obj=0.00786]
   INFO - 10:10:35:     35%|███▌      | 35/100 [00:00<00:01, 46.48 it/sec, obj=0.00314]
   INFO - 10:10:35:     36%|███▌      | 36/100 [00:00<00:01, 46.44 it/sec, obj=0.00103]
   INFO - 10:10:35:     37%|███▋      | 37/100 [00:00<00:01, 46.46 it/sec, obj=0.000147]
   INFO - 10:10:35:     38%|███▊      | 38/100 [00:00<00:01, 46.44 it/sec, obj=1.19e-5]
   INFO - 10:10:35:     39%|███▉      | 39/100 [00:00<00:01, 46.43 it/sec, obj=9.5e-7]
   INFO - 10:10:35:     40%|████      | 40/100 [00:00<00:01, 46.44 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     41%|████      | 41/100 [00:00<00:01, 47.10 it/sec, obj=6.45e-7]
   INFO - 10:10:35:     42%|████▏     | 42/100 [00:00<00:01, 48.23 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     43%|████▎     | 43/100 [00:00<00:01, 49.11 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     44%|████▍     | 44/100 [00:00<00:01, 49.95 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     45%|████▌     | 45/100 [00:00<00:01, 50.74 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     46%|████▌     | 46/100 [00:00<00:01, 51.57 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     47%|████▋     | 47/100 [00:00<00:01, 52.38 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     48%|████▊     | 48/100 [00:00<00:00, 53.17 it/sec, obj=6.41e-7]
   INFO - 10:10:35:     49%|████▉     | 49/100 [00:00<00:00, 53.96 it/sec, obj=6.41e-7]
   INFO - 10:10:35: Optimization result:
   INFO - 10:10:35:    Optimizer info:
   INFO - 10:10:35:       Status: None
   INFO - 10:10:35:       Message: Successive iterates of the design variables are closer than xtol_rel or xtol_abs. GEMSEO stopped the driver.
   INFO - 10:10:35:       Number of calls to the objective function by the optimizer: 57
   INFO - 10:10:35:    Solution:
   INFO - 10:10:35:       Objective: 6.40943911289712e-07
   INFO - 10:10:35:       Design space:
   INFO - 10:10:35:          +------+-------------+--------------------+-------------+-------+
   INFO - 10:10:35:          | Name | Lower bound |       Value        | Upper bound | Type  |
   INFO - 10:10:35:          +------+-------------+--------------------+-------------+-------+
   INFO - 10:10:35:          | x_0  |      -1     | 0.9996445136827841 |      1      | float |
   INFO - 10:10:35:          | x_1  |      -1     | 0.999284303839403  |      1      | float |
   INFO - 10:10:35:          | x_2  |      -1     | 0.9985690703536887 |      1      | float |
   INFO - 10:10:35:          +------+-------------+--------------------+-------------+-------+
   INFO - 10:10:35: *** End OptAsMDOScenario execution (time: 0:00:00.912933) ***

We can see that the numerical solution corresponds to the analytical one.

Total running time of the script: (0 minutes 1.412 seconds)

Gallery generated by Sphinx-Gallery