Source code for gemseo.post.robustness

# Copyright 2021 IRT Saint Exupéry, https://www.irt-saintexupery.com
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU Lesser General Public
# License version 3 as published by the Free Software Foundation.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
# Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with this program; if not, write to the Free Software Foundation,
# Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
# Contributors:
#    INITIAL AUTHORS - initial API and implementation and/or initial
#                           documentation
#        :author: Damien Guenot
#        :author: Francois Gallard
#    OTHER AUTHORS   - MACROSCOPIC CHANGES
"""Boxplots to quantify the robustness of the optimum."""

from __future__ import annotations

from math import sqrt
from typing import TYPE_CHECKING

import matplotlib.pyplot as plt
from numpy import zeros
from numpy.random import default_rng

from gemseo.post.core.robustness_quantifier import RobustnessQuantifier
from gemseo.post.opt_post_processor import OptPostProcessor
from gemseo.utils.seeder import SEED
from gemseo.utils.string_tools import repr_variable

if TYPE_CHECKING:
    from matplotlib.figure import Figure


[docs] class Robustness(OptPostProcessor): """Uncertainty quantification at the optimum. Compute the quadratic approximations of all the output functions, propagate analytically a normal distribution centered on the optimal design variables with a standard deviation which is a percentage of the mean passed in option (default: 1%) and plot the corresponding output boxplot. """ DEFAULT_FIG_SIZE = (8.0, 5.0) SR1_APPROX = "SR1" def _plot( self, stddev: float = 0.01, ) -> None: """ Args: stddev: The standard deviation of the normal uncertain variable to be added to the optimal design value; expressed as a fraction of the bounds of the design variables. """ # noqa: D205, D212, D415 self._add_figure(self.__boxplot(stddev)) def __boxplot( self, standard_deviation: float = 0.01, ) -> Figure: """Plot the Hessian of the function. Args: standard_deviation: The standard deviation of the normal uncertain variable to be added to the optimal design value; expressed as a fraction of the bounds of the design variables. Returns: A plot of the Hessian of the function. """ problem = self.opt_problem design_space = problem.design_space bounds_range = design_space.get_upper_bounds() - design_space.get_lower_bounds() n_x = problem.get_dimension() cov = zeros((n_x, n_x)) cov[range(n_x), range(n_x)] = (standard_deviation * bounds_range) ** 2 robustness = RobustnessQuantifier(self.database) function_samples = [] function_names = [] for func in self.opt_problem.get_all_functions(): func_name = database_func_name = func.name if self._change_obj and func_name == self._neg_obj_name: func_name = self._obj_name dim = func.dim at_most_niter = int(1.5 * n_x) for func_index in range(dim): robustness.compute_approximation( funcname=database_func_name, at_most_niter=at_most_niter, func_index=func_index, b0_mat=zeros((n_x, n_x)), ) x_ref = robustness.x_ref mean = robustness.compute_expected_value(x_ref, cov) if self._change_obj: mean = -mean variance = robustness.compute_variance(x_ref, cov) if variance > 0: # Otherwise normal doesn't work function_samples.append( default_rng(SEED).normal(mean, sqrt(variance), 500) ) function_names.append(repr_variable(func_name, func_index, dim)) fig = plt.figure(figsize=self.DEFAULT_FIG_SIZE) fig.suptitle( "Boxplot of the optimization functions " f"with normalized stddev {standard_deviation}" ) plt.boxplot(function_samples, showfliers=False, labels=function_names) return fig