import numpy as np
from dataclasses import dataclass
import warnings
import time
import os
from typing import Callable
from smt_optim.core import Problem, State
from smt_optim.surrogate_models import Surrogate
from smt_optim.acquisition_strategies import AcquisitionStrategy
from smt_optim.core import Evaluator
from smt_optim.utils.initial_design import generate_initial_design
from smt_optim.utils.stop_criteria import check_stop_criteria
from smt_optim.utils.logger import ConsoleLogger, JsonLogger
[docs]
def wrap_func(func: Callable, factor: float = 1, step: float = 0) -> Callable:
"""
Wrap function to return factor * (func - step).
:param func: Function to wrap.
:type func: Callable
:param factor: Multiplicative factor.
:type factor: float
:param step: Additive factor.
:type step: float
:return: Wrapped function.
:rtype: Callable
"""
def wrapped(x, f=func):
return factor * (f(x) - step)
return wrapped
[docs]
def wrap_array(
array: np.ndarray, factor: float | np.ndarray = 1.0, step: float | np.ndarray = 0.0
) -> np.ndarray:
return factor * (array - step)
[docs]
def check_bounds(x: np.ndarray, bounds: np.ndarray) -> np.ndarray:
"""
Apply L1 correction to x point to make sure it's within the problem's bounds.
:param x: Infill point.
:type x: np.ndarray
:param bounds: Problem boundaries.
:type bounds: np.ndarray
:return: The bounds corrected infill point.
:rtype: np.ndarray
"""
x_corrected = np.where(x < bounds[:, 0], bounds[:, 0], x)
x_corrected = np.where(x_corrected > bounds[:, 1], bounds[:, 1], x_corrected)
if np.any(x != x_corrected):
warnings.warn(
f"Infill point was outside of the bounds. L1 correction was applied: (initial = {x}; corrected = {x_corrected})."
)
return x_corrected
# def compute_rscv(self, cstr_array: np.ndarray, cstr_config: list[ConstraintConfig], g_tol: float = 0., h_tol: float = 0.) -> np.ndarray:
[docs]
def compute_rscv(
cstr_array: np.ndarray, cstr_config: list, g_tol: float = 0.0, h_tol: float = 0.0
) -> np.ndarray:
scv = np.full_like(cstr_array, 0.0) # Square Constraint Violation
for c_id, c_config in enumerate(cstr_config):
if c_config.type in ["less", "greater"]:
valid_mask = cstr_array[:, c_id] <= g_tol
scv[~valid_mask, c_id] = cstr_array[~valid_mask, c_id] ** 2
elif c_config.type == "equal":
valid_mask = np.abs(cstr_array[:, c_id]) <= h_tol
scv[~valid_mask, c_id] = cstr_array[~valid_mask, c_id] ** 2
else:
raise Exception(
f"{c_config.type} is not a valid constraint type. It must be 'less', 'greater' or 'equal'."
)
rscv = np.sqrt(scv.sum(axis=1))
return rscv
[docs]
@dataclass
class ObjectiveConfig:
"""
Configuration of the objective function used in the optimization problem.
This class stores an objective callable(s) together with surrogate
modeling information used to approximate the objective during optimization.
Attributes
----------
objective : list[Callable]
List of an objective functions. Each callable must accept a decision
variable vector ``x`` and return a scalar objective value. The functions
must be ordered in increasing level of fidelity.
type : {"minimize", "maximize"}, default="minimize"
Specifies whether the objective should be minimized or maximized.
surrogate : Surrogate or None, default=None
Surrogate model used to approximate the objective function.
surrogate_kwargs : dict or None, default=None
Optional keyword arguments passed to the surrogate model.
"""
objective: list[Callable]
surrogate: type[Surrogate]
type: str = "minimize" # problem's type -> "minimize" or "maximize"
surrogate_kwargs: dict | None = None
[docs]
@dataclass
class ConstraintConfig:
"""
Configuration of a constraint function used in the optimization problem.
This class stores a constraint callable(s) together with surrogate
modeling information used to approximate the constraint during optimization.
Attributes
----------
constraint : list[Callable]
List of constraint functions. Each callable must accept a decision
variable vector ``x`` and return a scalar constraint value. The functions
must be ordered in increasing level of fidelity.
lower: float | None
Lower bound of the constraint. If not specified, the constraint is considered
unconstrained in this direction.
upper: float | None
Upper bound of the constraint. If not specified, the constraint is considered
unconstrained in this direction.
equal: float | None
Equality constraint.
surrogate : Surrogate or None, default=None
Surrogate model used to approximate the constraint function.
surrogate_kwargs : dict or None, default=None
Optional keyword arguments passed to the surrogate model.
Notes
-----
A constraint must either be an inequality constraint or an equality constraint. For inequality
constraints, it is possible to define a lower and upper bound.
"""
constraint: list[Callable]
lower: float = None
upper: float = None
equal: float = None
surrogate: type[Surrogate] = None
surrogate_kwargs: dict | None = None
[docs]
@dataclass
class DriverConfig:
"""
Optimization driver configuration
Attributes
----------
max_iter: int
Maximum number of iteration
max_budget: float, default=inf
Maximum budget before termination of the optimization process
max_time: float, default=inf
Maximum time before termination of the optimization process
nt_init: int
Number of samples in the initial DoE
xt_init: list[np.ndarray]
Initial DoE to use. The Numpy array must be of shape (num_sample, num_dimension).
By providing an initial DoE, the driver will not generate an initial DoE. Cannot be
used with `nt_init`
results_dir: str or None, default=None
Name of the logging directory
verbose: bool, default=False
Print optimization information.
log_doe: bool, default=False
Log the value of the function values as soon as they are sampled. The values are
stored in a .csv file.
log_stats: bool, default=False
Log optimization statistics at the end of each iteration. The stats are store in
a .jsonl file.
scaling: bool, default=True
Scale the data. The objective is standardized. The constraints are divided by
their standard deviation. The design variables are normalized between 0 and 1.
seed: default=None
Seed for experiment reproducibility
"""
ctol: float = 1e-4 # tolerance for all constraints
max_iter: int | None = None # max number of BO iterations
max_budget: float = float("inf") # max BO budget
max_time: float = float("inf") # max BO elapsed time
nt_init: int | None = None # number of samples in initial DOE (with LHS)
xt_init: np.ndarray | None = (
None # initial training data [np.ndarray(nt, dim), np.ndarray(nt, dim)]
)
results_dir: str | None = "bo_result" # name for the results directory
verbose: bool = False # True/False print each iteration informations
log_doe: bool = False
log_stats: bool = False
callback_func: list[Callable] | Callable | None = (
None # additional method to call at the end of each iteration
)
scaling: bool = True # standardize the training data
seed: None = None
[docs]
class Driver:
def __init__(
self,
problem: Problem,
config: DriverConfig,
strategy: AcquisitionStrategy,
strategy_kwargs: dict = {},
):
"""
Initializes the object with the given parameters.
Parameters
----------
problem : Problem
The problem instance to be optimized.
config : DriverConfig
The configuration settings for the optimization process.
strategy : AcquisitionStrategy
The infill sampling strategy to use.
strategy_kwargs : dict, optional
Additional keyword arguments passed to the strategy constructor (default is an empty dictionary).
Notes
-----
This method sets up the object's internal state and initializes its components.
"""
self.problem = problem
self.config = config
self.state = State(problem)
self.state.dataset.log_data = True
self.strategy_kwargs = strategy_kwargs
self.strategy_kwargs["seed"] = config.seed
# = {"seed": config.seed if config.seed is not None else None}
self.strategy = strategy(self.state, **self.strategy_kwargs)
if self.config.log_doe or self.config.log_stats:
self.config.results_dir = self.make_res_dir(self.config.results_dir)
else:
self.config.results_dir = None
self.evaluator = Evaluator(problem, config.results_dir)
# setup loggers
self.loggers = []
if self.config.verbose:
self.loggers.append(ConsoleLogger(self.config))
if self.config.log_stats:
self.loggers.append(JsonLogger(self.config))
[docs]
def iteration(self, state):
"""
Performs an optimization iteration on the given state.
The iteration process involves the following steps:
1. Scaling the training data according to the specified scaling configuration.
2. Building surrogate models for approximating the expensive-to-evaluate functions.
3. Acquiring points to sample and their associated fidelity levels using the infill strategy.
4. Sampling the original, unmodified functions with the acquired points.
Parameters
----------
state : State
The optimization state on which to perform an iteration.
Returns
-------
State:
The updated optimization state after completing the iteration.
"""
state.iter += 1
# scale data
state.scale_dataset(self.config.scaling)
# build models
state.build_models()
# get infill
t0 = time.perf_counter()
infill = self.strategy.get_infill(state)
t1 = time.perf_counter()
state.iter_log["acq_opt_time"] = t1 - t0
# safe descale, clip infill to the bounds, verify design space
for i in range(len(infill)):
if infill[i] is not None:
infill[i] = safe_descale(infill[i], state)
state.iter_log["fidelity"] = i + 1
infill_not_in_xt(infill, state)
# evaluate infill points
self.evaluator.sample_func(infill, state)
# log iteration data
self.call_loggers(state)
return state
[docs]
def optimize(self):
"""
Performs an optimization process on the current state.
The process consists of two stages:
1. If the initial dataset is empty, it generates a Design of Experiment (DoE).
2. Iteratively performs optimization iterations until termination criteria are met.
Returns
-------
State:
The updated optimization state after completing the optimization process.
"""
self.start_optim()
# loop - check stop criteria
while check_stop_criteria(self.state, self.config):
# iteration
self.iteration(self.state)
return self.state
[docs]
def start_optim(self):
"""
Initializes the optimization process by creating an initial Design of Experiment (DoE) if necessary.
If the State dataset is empty (i.e., contains no samples), a new DoE will be generated.
Otherwise, no action is taken to avoid modifying existing sampling points.
Returns
-------
None
"""
# generate initial design
if len(self.state.dataset.samples) == 0:
generate_initial_design(self.state, self.evaluator, self.config)
self.call_loggers(self.state)
[docs]
def call_loggers(self, state):
# if self.loggers is not None:
for logger in self.loggers:
try:
logger.on_iter_end(state)
except Exception as e:
print(f"Error while logging: {e}")
[docs]
def make_res_dir(self, res_dir: str | None) -> str | None:
"""
Creates a unique results directory path based on the provided input.
If a directory with the same name already exists, the method will append an incrementing index to it
(e.g., "results" -> "results_1", "results_2", etc.) until a unique name is found. No directory will be
created if `res_dir` is set to `None`.
Parameters
----------
res_dir : Optional[str]
The desired results directory path. If `None`, no directory will be created.
Returns
-------
str or None:
The unique results directory path, or `None` if the input was `None`.
"""
og_res_dir = res_dir
if res_dir is not None:
# if results_dir already exists, append '_idx' to it to avoid overwriting existing data
idx = 1
while os.path.exists(res_dir):
res_dir = og_res_dir + f"_{idx}"
idx += 1
# create results_dir directory
os.makedirs(res_dir, exist_ok=False)
return res_dir
else:
return None
[docs]
def safe_descale(x_scaled: np.ndarray, state: State) -> np.ndarray:
"""
Descale point X and verify that it satisfies the design space. If necessary, correct it.
Parameters
----------
x_scaled: np.ndarray
state: State
Returns
-------
x_clip
"""
ds = state.problem.design_space
ds_bounds = ds.get_num_bounds()
# descale infill point
# descaling only applies to cont. variables -> int. and cat. variables are not scaled
x_raw = x_scaled * state.x_factor + state.x_step
# rounds cat. and int. variables
# only clips bounds for cat. variables
x_corr, _ = ds.correct_get_acting(x_raw)
# clips bounds for all variables
x_clip = np.clip(x_corr, ds_bounds[:, 0], ds_bounds[:, 1])
if np.linalg.norm(x_clip - x_raw) > 1e-10: # TODO: add customizable tolerance
warnings.warn("Infill point did not respect the design space.")
return x_clip
[docs]
def infill_not_in_xt(infills: list[np.ndarray], state: State) -> None:
"""
Raise exception if an infill point is already in the training data.
Parameters
----------
infills: list[np.ndarray]
state: State
Returns
-------
"""
dataset = state.dataset.export_as_dict()
fidelity = dataset["fidelity"].ravel()
xt = dataset["x"]
for lvl in range(len(infills)):
if infills[lvl] is not None:
fid_mask = fidelity == lvl
xt_lvl = xt[fid_mask, :]
for idx in range(infills[lvl].shape[0]):
diff = xt_lvl - infills[lvl][idx, :]
l2_norms = np.linalg.norm(diff, axis=1)
if np.min(l2_norms) < 1e-8: # TODO: add customizable tolerance
# raise Exception("Infill point already in training data.")
warnings.warn(
f"Infill point {infills[lvl][idx, :]} already in training data. L2 = {np.min(l2_norms)}"
)
return None
# # ======= old =======
# class OptimizerOld():
#
# def __init__(self, obj_config: ObjectiveConfig, config: OptimizerConfig, strategy: AcquisitionStrategy, strategy_params: dict | None = None):
#
# # initialize print setting
# self.verbose = config.verbose
#
# self.opt_data = {}
#
# self.results_dir = config.results_dir
# self.log_types = ["doe", "json"]
#
# # initialize objective function
# self.obj_config = obj_config
# self.obj_func = obj_config.objective
# self.domain = obj_config.design_space
# self.obj_type = obj_config.type
# self.obj_models = obj_config.surrogate
# self.costs = obj_config.costs
#
# # get constraint configurations
# self.cstr_config = config.constraints
# self.cstr_funcs = []
# self.ctol = config.ctol
#
# self.iter = 0
#
# # get stopping criteria
# self.max_iter = config.max_iter
# self.max_budget = config.max_budget
# self.max_time = config.max_time
#
# # get initial training data / setup
# self.nt_init = config.nt_init
# self.xt_init = config.xt_init
#
# # get misc configuration
# self.callback_func = config.callback_func
# self.scaling = config.scaling
# self.dynamic_costs = config.dynamic_costs
#
# self._check_optimizer_config()
# self._setup_logging()
#
# #
# self.strategy = strategy
# self.strategy_params = strategy_params
#
# self.num_dim = 0
# self.num_fidelity = 0
# self.num_obj = 0
# self.num_cstr = 0
#
# self.yt_factor = None
# self.ct_factor = None
# self.ct_step = None
#
# self._check_objective()
# self._check_constraints()
# self._setup_stopping_criteria()
# self._check_init_points()
#
# self.obj_models = []
# self.cstr_models = []
# self.g_models = []
# self.h_models = []
#
# self._initialize_surrogates()
#
# self.dataset = OptimizationDataset()
# self.scaled_dataset = None
# self.data = []
# self.xt = []
# self.yt = []
# self.ct = []
# self.f_min = np.inf
# self.rscv_min = np.inf
# self.x_min = None
# self.c_min = None
# self.samples_time = []
#
# self.xt_scaled = []
# self.yt_scaled = []
# self.ct_scaled = []
# self.f_min_scaled = np.inf
#
# # self._check_init_points()
# # self._gen_init_train_data()
# # self.update_f_min()
#
# self.acq_strategy = None
# self._initialize_acq_strategy()
#
# self.opt_data = {}
# self.iter_data = {}
#
# def _check_optimizer_config(self) -> None:
#
# if self.dynamic_costs is not None and self.dynamic_costs not in ["samples"]:
# raise Exception(f"Dynamic costs '{self.dynamic_costs}' is not supported.")
#
# if callable(self.callback_func):
# self.callback_func = [self.callback_func]
#
#
# def _setup_logging(self):
#
# results_dir = self.results_dir
#
# if results_dir is not None:
# # if results_dir already exists, append '_idx' to it to avoid overwriting existing data
# idx = 1
# while os.path.exists(results_dir):
# results_dir = self.results_dir + f"_{idx}"
# idx += 1
# self.results_dir = results_dir
#
# # create results_dir directory
# os.makedirs(self.results_dir)
#
# # create the DOE subdirectory (used to save the DOEs from each level)
# doe_path = os.path.join(self.results_dir, "DOE/")
# os.makedirs(doe_path, exist_ok=True)
#
#
# def _check_objective(self) -> None:
#
# if callable(self.obj_func):
# self.obj_func = [self.obj_func]
# elif type(self.obj_func) is list:
# pass
# else:
# raise Exception("ObjectiveConfig.objective must be of type list.")
#
# if self.obj_type == "minimize":
# maximize = False
# self.yt_factor = 1.
# elif self.obj_type == "maximize":
# maximize = True
# self.yt_factor = -1.
# else:
# raise Exception("ObjectiveConfig.type must be 'minimize' or 'maximize'.")
#
# # self.obj_func = self._wrap_objectives(self.obj_func, maximize=maximize)
#
# self.num_dim = self.domain.shape[0]
# self.num_obj += 1
#
# self.num_fidelity = len(self.obj_func)
#
# self.obj_func = [self.obj_func]
#
# if len(self.costs) != self.num_fidelity:
# raise Exception("ObjectiveConfig.costs must have the same number of levels as the objective.")
#
# # TODO: check costs are in ascending order
#
# def _check_constraints(self) -> None:
#
# if self.cstr_config is None:
# self.cstr_config = []
#
# self.num_cstr = len(self.cstr_config)
#
# if self.num_cstr > 0:
#
# self.ct_factor = np.empty(self.num_cstr)
# self.ct_step = np.empty(self.num_cstr)
#
# for c_id, c_config in enumerate(self.cstr_config):
#
# # self.cstr_funcs.append([])
#
# if c_config.type == "greater":
# self.ct_factor[c_id] = -1.
# else:
# self.ct_factor[c_id] = 1.
#
# self.ct_step[c_id] = c_config.value
#
# if callable(c_config.constraint):
# c_config.constraint = [c_config.constraint]
# elif type(c_config.constraint) is list:
# pass
# else:
# raise Exception("ConstraintConfig.constraint must be of type list.")
#
# if len(c_config.constraint) != self.num_fidelity:
# raise Exception("ConstraintConfig.constraint must have the same number of levels as the objective.")
#
# self.cstr_funcs.append( c_config.constraint )
# # self._wrap_constraints(c_config.constraint, c_config.type, c_config.value)
# # )
#
# def _setup_stopping_criteria(self):
#
# if self.num_dim == 0:
# raise Exception("Problem must have at least one dimension.")
#
# if self.max_iter is None:
# warnings.warn("Max number of iterations not specified. Set to 100.")
# self.max_iter = 10*self.num_dim
#
# if self.max_budget is None:
# self.max_budget = np.inf
#
# if self.max_time is None:
# self.max_time = np.inf
#
# def _check_init_points(self):
#
# if self.nt_init is not None and self.xt_init is not None:
# raise Exception("Define nt_init or xt_init, but not both.")
#
# elif self.nt_init is None:
# self.nt_init = 3*self.num_dim
#
# if self.xt_init is None:
# sampler = stats.qmc.LatinHypercube(self.num_dim)
# xt = sampler.random(self.nt_init)
# xt = stats.qmc.scale(xt, self.domain[:, 0], self.domain[:, 1])
# self.xt_init = [xt for _ in range(self.num_fidelity)]
#
# def _initialize_surrogates(self):
#
# self.obj_models = [self.obj_config.surrogate(optimizer=self)]
#
# for c_id, c_config in enumerate(self.cstr_config):
# self.cstr_models.append(
# c_config.surrogate(optimizer=self)
# )
#
# if c_config.type in ["less", "greater"]:
# self.g_models.append(self.cstr_models[c_id])
# elif c_config.type == "equal":
# self.h_models.append(self.cstr_models[c_id])
#
# def _gen_init_train_data(self):
#
# for lvl in range(self.num_fidelity):
# for idx in range(self.xt_init[lvl].shape[0]):
# self.sample_point(self.xt_init[lvl][idx, :], lvl)
#
# def _initialize_acq_strategy(self):
#
# acq_context = self.generate_acq_context()
#
# if type(self.strategy_params) is dict:
# self.acq_strategy = self.strategy(acq_context, **self.strategy_params)
# else:
# self.acq_strategy = self.strategy(acq_context)
#
# self.acq_strategy.validate_config(acq_context)
#
# # def update_f_min(self):
# # # feasible_mask = np.any(self.ct[-1] <= 1e-4, axis=1) # use cstr_tol in ConstraintConfig
# # # self.f_min = np.min(np.where(feasible_mask == True, self.yt[-1], np.inf))
# # # print(f"f_min = {self.f_min}")
# #
# # previous_f_min = self.f_min
# #
# # #feas_mask = np.all(self.ct[-1] <= self.ctol, axis=1)
# # valid_cstr = np.full((self.ct[-1].shape[0], self.num_cstr), False)
# #
# # for c_id, c_config in enumerate(self.cstr_config):
# # if c_config.type in ["less", "greater"]:
# # valid_cstr[:, c_id] = np.where(self.ct[-1][:, c_id] <= self.ctol, True, False)
# # elif c_config.type == "equal":
# # valid_cstr[:, c_id] = np.where(np.abs(self.ct[-1][:, c_id]) <= self.ctol, True, False)
# #
# # feas_mask = np.all(valid_cstr, axis=1)
# #
# # if np.any(feas_mask):
# # local_min_id = self.yt[-1][feas_mask].argmin()
# # # global_min_id = feas_mask[local_min_id]
# #
# # self.f_min = self.yt[-1][feas_mask][local_min_id].item()
# # self.x_min = self.xt[-1][feas_mask][local_min_id]
# # self.c_min = self.ct[-1][feas_mask][local_min_id]
# #
# # else:
# # self.f_min = np.inf
# # self.x_min = None
# # self.c_min = None
# #
# # self._check_f_min_decreasing(self.f_min, previous_f_min)
# #
# # def update_rscv_min(self):
# # rscv = compute_rscv(self.ct[-1], self.cstr_config, g_tol=0.0, h_tol=0.0)
# # self.rscv_min = rscv.min()
# #
# #
# # def _check_f_min_decreasing(self, next_f_min, previous_f_min):
# # if previous_f_min < next_f_min:
# # warnings.warn("f_min is increasing.")
#
# def sample_point(self, x_new: np.ndarray, level: int) -> None:
#
# x_new = check_bounds(x_new, self.domain)
#
# obj_values = np.empty(self.num_obj)
# cstr_values = np.empty(self.num_cstr)
# times = np.empty(self.num_obj+self.num_cstr)
#
# def sample_func(x_new: np.ndarray, func: Callable) -> tuple[float, float]:
# t0 = time.perf_counter()
# value = func(x_new)
# t1 = time.perf_counter()
# return value, t1-t0
#
# for obj_idx in range(self.num_obj):
# obj_values[obj_idx], times[obj_idx] = sample_func(x_new, self.obj_func[obj_idx][level])
#
# for cstr_idx in range(self.num_cstr):
# cstr_values[cstr_idx], times[self.num_obj+cstr_idx] = sample_func(x_new, self.cstr_funcs[cstr_idx][level])
#
# sample = Sample(
# x=x_new,
# fidelity=level,
# obj=obj_values,
# cstr=cstr_values,
# eval_time=times,
# metadata={"iter": self.iter}
# )
#
# self.dataset.add(sample)
#
# if "doe" in self.log_types:
# self.save_sample(sample)
#
# def _standardize_data(self, data: np.ndarray) -> tuple[np.ndarray | float]:
#
# mean = data.mean()
# std = data.std()
# std_data = (data - mean)/std
#
# return std_data, mean, std
#
# def scale_data(self):
#
# num_qoi = self.num_obj + self.num_cstr
# qoi_factor = np.empty(num_qoi)
# qoi_step = np.empty(num_qoi)
#
# for obj_idx in range(self.num_obj):
# if self.obj_config.type == "minimize":
# qoi_factor[obj_idx] = 1
# elif self.obj_config.type == "maximize":
# qoi_factor[obj_idx] = -1
#
# qoi_step[obj_idx] = 0
#
# for cstr_idx in range(self.num_cstr):
# c_config = self.cstr_config[cstr_idx]
#
# if c_config.type in ["less", "equal"]:
# qoi_factor[self.num_obj+cstr_idx] = 1
# elif c_config.type in ["greater"]:
# qoi_factor[self.num_obj+cstr_idx] = -1
#
# qoi_step[self.num_obj+cstr_idx] = c_config.value
#
# self.qoi_factor = qoi_factor
# self.qoi_step = qoi_step
#
# self.scaled_dataset = OptimizationDataset()
#
# for sample in self.dataset.samples:
#
# scaled_sample = copy.deepcopy(sample)
#
# # should only normalize real variables
# scaled_sample.x -= self.domain[:, 0]
# scaled_sample.x /= (self.domain[:, 1] - self.domain[:, 0])
# scaled_sample.obj[:] *= self.qoi_factor[:self.num_obj]
# scaled_sample.cstr[:] *= self.qoi_factor[self.num_obj:self.num_obj+self.num_cstr]
#
# self.scaled_dataset.add(scaled_sample)
#
#
# def build_models(self):
#
# def group_by_fidelity() -> tuple[list[np.ndarray], list[np.ndarray]]:
#
# x = []
# qoi = []
#
# for lvl in range(self.num_fidelity):
#
# samples = self.scaled_dataset.get_by_fidelity(lvl)
#
# x_lvl = np.empty((len(samples), self.num_dim))
# qoi_lvl = np.empty((len(samples), self.num_obj + self.num_cstr))
#
# for idx, sample in enumerate(samples):
# x_lvl[idx, :] = sample.x
# qoi_lvl[idx, :self.num_obj] = sample.obj
# qoi_lvl[idx, self.num_obj:] = sample.cstr
#
# x.append(x_lvl)
# qoi.append(qoi_lvl)
#
# return x, qoi
#
# x_train, qoi_train = group_by_fidelity()
#
# qoi_models = self.obj_models + self.cstr_models
#
# for qoi_idx in range(self.num_obj+self.num_cstr):
#
# idx_train = []
#
# for lvl in range(self.num_fidelity):
# idx_train.append(qoi_train[lvl][:, qoi_idx].reshape(-1, 1))
#
# qoi_models[qoi_idx].train(x_train, idx_train)
#
#
# def generate_acq_context(self):
#
# cstr_types = []
# for c_id, c_config in enumerate(self.cstr_config):
# if c_config.type in ["less", "greater"]:
# cstr_types.append("less")
# elif c_config.type == "equal":
# cstr_types.append("equal")
#
# acq_context = State(
# num_dim=self.num_dim,
# num_obj=1,
# num_cstr=self.num_cstr,
# num_fidelity=self.num_fidelity,
# design_space=np.array([[0, 1]] * self.num_dim),
# obj_models=self.obj_models,
# cstr_models=self.cstr_models,
# cstr_types=cstr_types,
# dataset=self.scaled_dataset,
# )
#
# return acq_context
#
# def sample_infills(self, next_x: list[np.ndarray]) -> None:
#
# # Convert the single fidelity acquisition function output to a list as
# # to make it compatible with the multi-fidelity approach
# if type(next_x) is not list:
# next_x = [next_x]
#
# if len(next_x) != self.num_fidelity:
# warnings.warn("")
#
# # Sample each fidelity level sequentially
# for k in range(self.num_fidelity):
#
# # if None -> the fidelity k is not to be sampled
# if self.next_x[k] is None:
# continue
#
# # unscale infill point
# if self.scaling:
# self.next_x[k] *= (self.domain[:, 1] - self.domain[:, 0])
# self.next_x[k] += self.domain[:, 0]
#
# # Check if the infill point is already in the training data
# # TODO: what to do if the next infill location is already in the training data?
# # if np.any(np.all(self.xt[k] == self.next_x[k], axis=1)):
# # warnings.warn("Infill point is already in the training data.")
# # continue
#
# # sample objective function and constraints
# self.sample_point(self.next_x[k], k)
#
#
# def perform_iteration(self):
#
# self.iter_data = dict() # reset iteration data dictionary
#
# # ------- Wrap training data -------
# # self._wrap_training_data()
#
# # ------- Scale training data -------
# # self.scale_training_data()
# self.scale_data()
#
# # ------- Update cost ratio -------
# self.update_costs()
#
# # ------- Surrogate models training -------
# gp_t0 = time.perf_counter()
#
# self.build_models()
#
# gp_t1 = time.perf_counter()
# gp_time = gp_t1 - gp_t0 # elapsed time for training the models
#
# # log gp training elapsed time
# self.iter_data["gp_training_time"] = gp_time
#
# # # ------- Acquisition function optimization -------
# acq_context = self.generate_acq_context()
#
# acq_t0 = time.perf_counter()
#
# # Find enrichment location
# self.next_x, acq_data = self.acq_strategy.get_infill(acq_context)
# self.iter_data["acq_data"] = acq_data
#
# acq_t1 = time.perf_counter()
# acq_time = acq_t1 - acq_t0 # elapsed time for finding the next acquisition point
#
# # log acquisition function maximization time
# self.iter_data["acq_opt_time"] = acq_time
#
# # ------- Sample infill location -------
# self.sample_infills(self.next_x)
#
# if self.callback_func is not None:
# for func in self.callback_func:
# func(self)
#
#
# def optimize(self):
#
# bo_start = time.perf_counter()
#
# # generate initial doe
# self._gen_init_train_data()
#
# # self.scale_training_data()
# self.scale_data()
#
# # update f_min
# # self.update_f_min()
# # self.update_rscv_min()
#
# iter_id = 0
# self.iter = iter_id
#
# # for l in range(self.num_fidelity):
# # self.iter_data[f"n{l + 1}"] = len(self.yt[l])
# #
# # self.update_costs()
# # self.budget = self.compute_used_budget()
# # self.iter_data["budget"] = self.compute_used_budget()
# # self.iter_data["f_min"] = self.f_min
# # self.iter_data["x_min"] = self.x_min
# # self.iter_data["c_min"] = self.c_min
#
# self.opt_data[0] = self.iter_data
# self.save_data()
#
# self.continue_bo = True
#
# # if self.verbose: print(
# # f"| iter= {iter_id}/{self.max_iter} | budget={self.budget:.3f}/{self.max_budget:.3f} | f_min={self.f_min:.3e} | rscv_min={self.rscv_min:.3e} |"
# # )
#
# while self.continue_bo:
#
# iter_id += 1
# self.iter = iter_id
#
# self.perform_iteration()
#
#
# # update f_min
# # self.update_f_min()
# # self.update_rscv_min()
#
# # for l in range(self.num_fidelity):
# # self.iter_data[f"n{l + 1}"] = len(self.yt[l])
# #
# # self.iter_data["f_min"] = self.f_min
# # self.iter_data["x_min"] = self.x_min
# # self.iter_data["c_min"] = self.c_min
#
# self.budget = self.compute_used_budget()
# self.iter_data["budget"] = self.budget
#
# self.iter_data["costs"] = self.costs
#
# # elapsed time since optimization start
# self.bo_time = time.perf_counter() - bo_start
# self.iter_data["bo_time"] = self.bo_time
#
# self.continue_bo = self.check_stop_criteria()
# self.iter_data["continue"] = self.continue_bo
#
# # add iteration data to the optimization data dictionary
# self.opt_data[iter_id] = self.iter_data
#
# # self.dump_pikle_log()
# self.save_data()
#
# # Display the iteration number, best feasible objective and fidelity level sampled
# # if self.verbose : print(
# # f"| iter= {iter_id}/{self.max_iter} | budget={self.budget:.3f}/{self.max_budget:.3f} | f_min={self.f_min:.3e} | rscv_min={self.rscv_min:.3e} | lvl={max_level}/{self.num_fidelity - 1} | gp_time={gp_time:.3f} | acq_time={acq_time:.3f}")
#
# if self.verbose : print(
# f"| iter= {iter_id}/{self.max_iter} | budget={self.budget:.3f}/{self.max_budget:.3f}")
#
# # ------- End of optimization loop -------
#
# return self.opt_data
#
# # def dump_pikle_log(self):
# # try:
# #
# # path = os.path.join(self.results_dir, "opt_data.pkl")
# #
# # with open(path, 'wb') as file:
# # pickle.dump(self.opt_data, file)
# #
# # except Exception as e:
# # # TODO: use warnings.warn
# # warnings.warn(f"Error while saving optimization data: {e}")
#
#
# def dump_json_log(self):
# try:
#
# path = os.path.join(self.results_dir, "opt_data.json")
#
# with open(path, 'w') as file:
# safe_data = json_safe(self.opt_data)
# json.dump(safe_data, file, indent=2)
#
# except Exception as e:
# # TODO: use warnings.warn
# warnings.warn(f"Error while saving optimization data: {e}")
#
#
# def save_data(self):
#
# if "json" in self.log_types:
# self.dump_json_log()
#
#
#
#
# def save_sample(self, sample: Sample) -> None:
#
# try:
# row = dict()
#
# row["iter"] = self.iter
# row["budget"] = self.compute_used_budget() # self.budget
#
# # save variable
# for i in range(self.num_dim):
# row[f"x{i}"] = sample.x[i]
#
# # save objectives
# for i in range(self.num_obj):
# row[f"{i}"] = sample.obj[i]
#
# # save constraints
# for i in range(self.num_cstr):
# row[f"c{i}"] = sample.cstr[i]
#
# row["time"] = np.sum(sample.eval_time)
#
# path = os.path.join(self.results_dir, "DOE", f"doe_fidelity_{sample.fidelity}.csv")
# file_exists = os.path.isfile(path)
#
# # possibly does not work on Windows -> to be tested
# with open(path, 'a') as file:
# writer = csv.DictWriter(file, fieldnames=row.keys())
#
# if not file_exists:
# writer.writeheader()
#
# writer.writerow(row)
#
# except Exception as e:
# print(f"Error while saving the DoE: {e}")
#
#
# def dump_csv_doe(self, level):
# pass
#
# try:
# row = dict()
#
# row["iter"] = self.iter
# row["budget"] = self.budget
#
# x = self.xt[level][-1, :]
# print(f"x = {x}")
# for i in range(len(x)):
# row[f"x{i}"] = x[i]
#
# row["f"] = self.data[level][-1, 0]
#
# for i in range(self.num_cstr):
# row[f"c{i}"] = self.data[level][-1, i+1]
#
# path = f"DoE/{self.log_filename}_{level}.csv"
# file_exists = os.path.isfile(path)
# print(f"file exists = {file_exists}")
#
# with open(path, 'w') as file:
# writer = csv.DictWriter(file, fieldnames=row.keys())
#
# if not file_exists:
# writer.writeheader()
#
# writer.writerow(row)
#
# except Exception as e:
# print(f"Error while logging the DoE: {e}")
#
# def compute_used_budget(self):
#
# budget = 0
#
# for k in range(self.num_fidelity):
# samples = self.dataset.get_by_fidelity(k)
# budget += self.costs[k] * len(samples)
#
# return budget
#
# def check_stop_criteria(self):
#
# if self.iter >= self.max_iter:
# return False
# elif self.budget >= self.max_budget:
# return False
# elif self.bo_time >= self.max_time:
# return False
# else:
# return True
#
# def update_costs(self) -> None:
# """
# Update the costs of each level.
# - If set to None, the costs are never updated.
# - If set to 'samples', the cost of each level corresponds to it's average time to be sampled.
#
# :return: None
# """
#
# if self.dynamic_costs == "samples":
# for lvl in range(self.num_fidelity):
# # average sampling time per level
# self.costs[lvl] = self.samples_time[lvl].sum(axis=1).mean().item()
#
# # def scale_training_data(self):
# #
# # self.xt_scaled = []
# # self.yt_scaled = []
# # self.ct_scaled = []
# #
# # for lvl in range(self.num_fidelity):
# #
# # self.xt_scaled.append(self.xt[lvl].copy())
# #
# # # transform the objective into a minimization problem
# # self.yt_scaled.append(wrap_array(self.yt[lvl], factor=self.yt_factor))
# #
# # if self.num_cstr >= 1:
# # # transform the constraints to define the feasible domain as: g <= 0 and h == 0
# # self.ct_scaled.append(wrap_array(self.ct[lvl], factor=self.ct_factor, step=self.ct_step))
# #
# # if self.scaling:
# # # scale xt between 0 and 1
# # self.xt_scaled[lvl][:] -= self.domain[:, 0]
# # self.xt_scaled[lvl][:] /= (self.domain[:, 1] - self.domain[:, 0])
# #
# # # update scaled domain boundaries
# # self.domain_scaled = np.empty((self.num_dim, 2))
# # self.domain_scaled[:, 0] = 0.
# # self.domain_scaled[:, 1] = 1.
# #
# #
# # # scaled objective to unit std
# # yt_scaled, yt_mean, yt_std = self._standardize_data(self.yt[lvl])
# # self.yt_scaled[lvl] = yt_scaled
# #
# # # update minimum objective
# # self.f_min_scaled = (self.f_min - yt_mean) / yt_std
# #
# # # scaled constraints to unit std
# # if self.num_cstr >= 1:
# # for c_id in range(self.ct[lvl].shape[1]):
# # self.ct_scaled[lvl][:, c_id] /= self.ct[lvl][:, c_id].std()
# # else:
# # self.domain_scaled = self.domain.copy()
if __name__ == "__main__":
pass
