Note
Go to the end to download the full example code.
Change the seed of a DOE#
Latin hypercube sampling is an example of stochastic DOE algorithm: given an input dimension and a number of samples, running the algorithm twice will give two different DOEs.
For the sake of reproducibility, GEMSEO uses a random seed: given an input dimension, a number of samples and a random seed, running the algorithm twice will give the same DOE.
In this example, we will see how GEMSEO uses the random seed and how the user can change its value.
from __future__ import annotations
from gemseo import create_design_space
from gemseo import create_discipline
from gemseo import create_scenario
from gemseo import execute_algo
from gemseo.algos.doe.openturns.openturns import OpenTURNS
from gemseo.algos.optimization_problem import OptimizationProblem
from gemseo.core.mdo_functions.mdo_function import MDOFunction
At the scenario level#
First,
we illustrate the use of the random seed at the DOEScenario level
which is the appropriate level for most users.
Then,
we will illustrate this use at the OptimizationProblem level
which can be useful for developers.
Let us consider an Discipline representing the function \(y=x^2\):
discipline = create_discipline("AnalyticDiscipline", expressions={"y": "x**2"})
This function is defined over the interval \([-1,1]\):
design_space = create_design_space()
design_space.add_variable("x", lower_bound=-1, upper_bound=1)
We want to sample this discipline over this design space.
For that,
we express the sampling problem as a DOEScenario:
scenario = create_scenario(
[discipline],
"y",
design_space,
scenario_type="DOE",
formulation_name="DisciplinaryOpt",
)
and solve it:
scenario.execute(algo_name="OT_OPT_LHS", n_samples=2)
scenario.formulation.optimization_problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: *** Start DOEScenario execution ***
INFO - 16:25:00: DOEScenario
INFO - 16:25:00: Disciplines: AnalyticDiscipline
INFO - 16:25:00: MDO formulation: DisciplinaryOpt
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize y(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: | x | -1 | None | 1 | float |
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: 50%|█████ | 1/2 [00:00<00:00, 605.50 it/sec, feas=True, obj=0.802]
INFO - 16:25:00: 100%|██████████| 2/2 [00:00<00:00, 978.83 it/sec, feas=True, obj=0.54]
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: *** End DOEScenario execution (time: 0:00:00.005349) ***
[array([-0.89534931]), array([0.73475608])]
When using the same DOE algorithm, solving again this problem with the same scenario leads to a new result:
scenario.execute(algo_name="OT_OPT_LHS", n_samples=2)
scenario.formulation.optimization_problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: *** Start DOEScenario execution ***
INFO - 16:25:00: DOEScenario
INFO - 16:25:00: Disciplines: AnalyticDiscipline
INFO - 16:25:00: MDO formulation: DisciplinaryOpt
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize y(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: 50%|█████ | 1/2 [00:00<00:00, 4686.37 it/sec, feas=True, obj=0.00248]
INFO - 16:25:00: 100%|██████████| 2/2 [00:00<00:00, 3623.59 it/sec, feas=True, obj=0.237]
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.0024824697156153775
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: *** End DOEScenario execution (time: 0:00:00.003039) ***
[array([0.04982439]), array([-0.48653892])]
The value of the default seed was incremented
from 0 (at first execution) to 1 (at last execution),
as we can check it by setting the custom "seed" to 1:
scenario.execute(algo_name="OT_OPT_LHS", n_samples=2, seed=1)
scenario.formulation.optimization_problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: *** Start DOEScenario execution ***
INFO - 16:25:00: DOEScenario
INFO - 16:25:00: Disciplines: AnalyticDiscipline
INFO - 16:25:00: MDO formulation: DisciplinaryOpt
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize y(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.0024824697156153775
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:00: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:00: *** End DOEScenario execution (time: 0:00:00.002405) ***
[array([0.04982439]), array([-0.48653892])]
The default seed is incremented at each execution, whatever the custom one.
At the problem level#
Basic#
Let us consider an MDOFunction representing the function \(y=x^2\):
function = MDOFunction(lambda x: x**2, "f", input_names=["x"], output_names=["y"])
and defined over the unit interval \(x\in[0,1]\):
design_space = create_design_space()
design_space.add_variable("x", lower_bound=-1, upper_bound=1)
We want to sample this function over this design space.
For that,
we express the sampling problem as an OptimizationProblem:
problem = OptimizationProblem(design_space)
problem.objective = function
and solve it:
execute_algo(problem, algo_name="OT_OPT_LHS", algo_type="doe", n_samples=2)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: | x | -1 | None | 1 | float |
INFO - 16:25:00: +------+-------------+-------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: 50%|█████ | 1/2 [00:00<00:00, 10672.53 it/sec, feas=True, obj=0.802]
INFO - 16:25:00: 100%|██████████| 2/2 [00:00<00:00, 5911.63 it/sec, feas=True, obj=0.54]
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
[array([-0.89534931]), array([0.73475608])]
- Note:
We use the function
execute_algo()as theOptimizationProblemdoes not have a methodexecute()unlike theScenario.
Solving again this problem with the same configuration leads to the same result:
execute_algo(problem, algo_name="OT_OPT_LHS", algo_type="doe", n_samples=2)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
[array([-0.89534931]), array([0.73475608])]
and the result is still the same if we take 1 as random seed, as 1 is the default value of this seed:
execute_algo(problem, algo_name="OT_OPT_LHS", algo_type="doe", n_samples=2, seed=1)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
[array([-0.89534931]), array([0.73475608])]
If you want to use a different random seed,
you only have to change the value of seed:
execute_algo(problem, algo_name="OT_OPT_LHS", algo_type="doe", n_samples=2, seed=3)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: 50%|█████ | 1/2 [00:00<00:00, 11184.81 it/sec, feas=True, obj=0.867]
INFO - 16:25:00: 100%|██████████| 2/2 [00:00<00:00, 6316.72 it/sec, feas=True, obj=0.663]
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
[array([0.93088508]), array([-0.81399054])]
Advanced#
You can also solve your problem with a lower level API
by directly instantiating the BaseDOELibrary of interest.
A BaseDOELibrary has a default seed generated by a Seeder
that is incremented at the beginning of each execution:
algo = OpenTURNS("OT_OPT_LHS")
algo.execute(problem, n_samples=2)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:00: over the design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:00: Optimization result:
INFO - 16:25:00: Optimizer info:
INFO - 16:25:00: Status: None
INFO - 16:25:00: Message: None
INFO - 16:25:00: Solution:
INFO - 16:25:00: Objective: 0.5398664901859535
INFO - 16:25:00: Design space:
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:00: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:00: +------+-------------+--------------------+-------------+-------+
[array([0.93088508]), array([-0.81399054])]
Solving again the problem will give different samples:
algo.execute(problem, n_samples=2)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:00: Optimization problem:
INFO - 16:25:00: minimize f(x)
INFO - 16:25:00: with respect to x
INFO - 16:25:01: over the design space:
INFO - 16:25:01: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:01: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:01: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:01: | x | -1 | 0.7347560752970699 | 1 | float |
INFO - 16:25:01: +------+-------------+--------------------+-------------+-------+
INFO - 16:25:01: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:01: 50%|█████ | 1/2 [00:00<00:00, 10538.45 it/sec, feas=True, obj=0.00248]
INFO - 16:25:01: 100%|██████████| 2/2 [00:00<00:00, 5940.94 it/sec, feas=True, obj=0.237]
INFO - 16:25:01: Optimization result:
INFO - 16:25:01: Optimizer info:
INFO - 16:25:01: Status: None
INFO - 16:25:01: Message: None
INFO - 16:25:01: Solution:
INFO - 16:25:01: Objective: 0.0024824697156153775
INFO - 16:25:01: Design space:
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
[array([0.04982439]), array([-0.48653892])]
You can also use a specific seed instead of the default one:
algo.execute(problem, n_samples=2, seed=123)
problem.database.get_last_n_x_vect(2)
INFO - 16:25:01: Optimization problem:
INFO - 16:25:01: minimize f(x)
INFO - 16:25:01: with respect to x
INFO - 16:25:01: over the design space:
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: Solving optimization problem with algorithm OT_OPT_LHS:
INFO - 16:25:01: 50%|█████ | 1/2 [00:00<00:00, 10255.02 it/sec, feas=True, obj=0.00904]
INFO - 16:25:01: 100%|██████████| 2/2 [00:00<00:00, 5979.05 it/sec, feas=True, obj=0.968]
INFO - 16:25:01: Optimization result:
INFO - 16:25:01: Optimizer info:
INFO - 16:25:01: Status: None
INFO - 16:25:01: Message: None
INFO - 16:25:01: Solution:
INFO - 16:25:01: Objective: 0.0024824697156153775
INFO - 16:25:01: Design space:
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | Name | Lower bound | Value | Upper bound | Type |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
INFO - 16:25:01: | x | -1 | 0.04982438876308848 | 1 | float |
INFO - 16:25:01: +------+-------------+---------------------+-------------+-------+
[array([-0.09505332]), array([0.98400738])]
Total running time of the script: (0 minutes 0.038 seconds)