From 68d40e1c087dde43dd55c013406242c0dd2ac918 Mon Sep 17 00:00:00 2001 From: Nicola Demo Date: Fri, 9 Jan 2026 11:31:25 +0100 Subject: [PATCH 1/2] Black formatter --- docs/source/installation.rst | 1 - ezyrb/__init__.py | 23 +- ezyrb/approximation/__init__.py | 11 +- ezyrb/approximation/ann.py | 9 + ezyrb/approximation/approximation.py | 1 + ezyrb/approximation/gpr.py | 35 +- ezyrb/approximation/kneighbors_regressor.py | 10 +- ezyrb/approximation/linear.py | 37 +- ezyrb/approximation/neighbors_regressor.py | 7 +- .../radius_neighbors_regressor.py | 7 +- ezyrb/approximation/rbf.py | 41 +- ezyrb/approximation/sklearn_approximation.py | 31 +- ezyrb/database.py | 89 +++-- ezyrb/parallel/__init__.py | 15 +- ezyrb/parallel/ae.py | 38 +- ezyrb/parallel/ae_eddl.py | 172 +++++---- ezyrb/parallel/ann.py | 7 +- ezyrb/parallel/approximation.py | 1 + ezyrb/parallel/gpr.py | 39 +- ezyrb/parallel/kneighbors_regressor.py | 1 + ezyrb/parallel/linear.py | 21 +- ezyrb/parallel/neighbors_regressor.py | 5 +- ezyrb/parallel/pod.py | 41 +- ezyrb/parallel/radius_neighbors_regressor.py | 1 + ezyrb/parallel/rbf.py | 31 +- ezyrb/parallel/reducedordermodel.py | 76 ++-- ezyrb/parallel/reduction.py | 1 + ezyrb/parameter.py | 38 +- ezyrb/plugin/__init__.py | 16 +- ezyrb/plugin/aggregation.py | 147 ++++--- ezyrb/plugin/automatic_shift.py | 125 ++++-- ezyrb/plugin/database_splitter.py | 64 +-- ezyrb/plugin/plugin.py | 34 +- ezyrb/plugin/scaler.py | 5 +- ezyrb/plugin/shift.py | 46 ++- ezyrb/reducedordermodel.py | 364 ++++++++++-------- ezyrb/reduction/__init__.py | 10 +- ezyrb/reduction/ae.py | 58 +-- ezyrb/reduction/pod.py | 21 +- ezyrb/reduction/pod_ae.py | 11 +- ezyrb/reduction/reduction.py | 3 +- ezyrb/reduction/sklearn_reduction.py | 40 +- ezyrb/regular_grid.py | 32 +- ezyrb/snapshot.py | 62 +-- pyproject.toml | 4 +- 45 files changed, 1065 insertions(+), 766 deletions(-) diff --git a/docs/source/installation.rst b/docs/source/installation.rst index 146f1057..bd89ab62 100644 --- a/docs/source/installation.rst +++ b/docs/source/installation.rst @@ -11,7 +11,6 @@ EZyRB requires: - scipy - matplotlib - scikit-learn -- torch (for neural network-based methods) Install via pip --------------- diff --git a/ezyrb/__init__.py b/ezyrb/__init__.py index e8ac9502..f57bb34c 100644 --- a/ezyrb/__init__.py +++ b/ezyrb/__init__.py @@ -1,11 +1,24 @@ """EZyRB package""" __all__ = [ - 'Database', 'Snapshot', 'Reduction', 'POD', 'Approximation', 'RBF', - 'Linear', 'GPR', 'ANN', 'KNeighborsRegressor', - 'RadiusNeighborsRegressor', 'AE', 'ReducedOrderModel', 'PODAE', - 'RegularGrid', 'MultiReducedOrderModel', 'SklearnApproximation', - 'SklearnReduction' + "Database", + "Snapshot", + "Reduction", + "POD", + "Approximation", + "RBF", + "Linear", + "GPR", + "ANN", + "KNeighborsRegressor", + "RadiusNeighborsRegressor", + "AE", + "ReducedOrderModel", + "PODAE", + "RegularGrid", + "MultiReducedOrderModel", + "SklearnApproximation", + "SklearnReduction", ] from .database import Database diff --git a/ezyrb/approximation/__init__.py b/ezyrb/approximation/__init__.py index 15aa07ce..69f06a84 100644 --- a/ezyrb/approximation/__init__.py +++ b/ezyrb/approximation/__init__.py @@ -1,9 +1,14 @@ """EZyRB package""" __all__ = [ - 'Approximation', 'RBF', 'Linear', 'GPR', - 'ANN', 'KNeighborsRegressor', 'RadiusNeighborsRegressor', - 'SklearnApproximation' + "Approximation", + "RBF", + "Linear", + "GPR", + "ANN", + "KNeighborsRegressor", + "RadiusNeighborsRegressor", + "SklearnApproximation", ] from .approximation import Approximation diff --git a/ezyrb/approximation/ann.py b/ezyrb/approximation/ann.py index 269ddef1..1deaa632 100755 --- a/ezyrb/approximation/ann.py +++ b/ezyrb/approximation/ann.py @@ -46,7 +46,16 @@ class ANN(Approximation): >>> print(y_pred) >>> print(len(ann.loss_trend)) >>> print(ann.loss_trend[-1]) + + .. note:: + This module provides a wrapper around sklearn's MLPRegressor for + multidimensional function approximation using feed-forward neural + networks. It is not intended for deep learning tasks, but rather for + approximating functions based on given data points. For more advanced + deep learning applications, consider using dedicated libraries such as + :ref:`PINA`. """ + def __init__( self, layers, diff --git a/ezyrb/approximation/approximation.py b/ezyrb/approximation/approximation.py index 862ad788..914df668 100644 --- a/ezyrb/approximation/approximation.py +++ b/ezyrb/approximation/approximation.py @@ -10,6 +10,7 @@ class Approximation(ABC): All the classes that implement the input-output mapping should be inherited from this class. """ + @abstractmethod def fit(self, points, values): """Abstract `fit`""" diff --git a/ezyrb/approximation/gpr.py b/ezyrb/approximation/gpr.py index 17c40153..d9fff86b 100644 --- a/ezyrb/approximation/gpr.py +++ b/ezyrb/approximation/gpr.py @@ -1,6 +1,7 @@ """ Module wrapper exploiting `GPy` for Gaussian Process Regression """ + import logging import numpy as np from scipy.optimize import minimize @@ -40,19 +41,24 @@ class GPR(Approximation): >>> print(np.allclose(y, y_pred)) """ + def __init__(self, kern=None, normalizer=True, optimization_restart=20): """ Initialize a Gaussian Process Regressor. - + :param kern: Kernel object from sklearn. Default is None. :param bool normalizer: Whether to normalize values. Default is True. :param int optimization_restart: Number of restarts for optimization. Default is 20. """ - logger.debug("Initializing GPR with kernel=%s, normalizer=%s, " - "optimization_restart=%d", - kern, normalizer, optimization_restart) + logger.debug( + "Initializing GPR with kernel=%s, normalizer=%s, " + "optimization_restart=%d", + kern, + normalizer, + optimization_restart, + ) self.X_sample = None self.Y_sample = None self.kern = kern @@ -67,8 +73,11 @@ def fit(self, points, values): :param array_like points: the coordinates of the points. :param array_like values: the values in the points. """ - logger.debug("Fitting GPR with points shape: %s, values shape: %s", - np.array(points).shape, np.array(values).shape) + logger.debug( + "Fitting GPR with points shape: %s, values shape: %s", + np.array(points).shape, + np.array(values).shape, + ) self.X_sample = np.array(points) self.Y_sample = np.array(values) if self.X_sample.ndim == 1: @@ -80,8 +89,10 @@ def fit(self, points, values): logger.debug("Creating GaussianProcessRegressor") self.model = GaussianProcessRegressor( - kernel=self.kern, n_restarts_optimizer=self.optimization_restart, - normalize_y=self.normalizer) + kernel=self.kern, + n_restarts_optimizer=self.optimization_restart, + normalize_y=self.normalizer, + ) self.model.fit(self.X_sample, self.Y_sample) logger.info("GPR fitted successfully") @@ -118,14 +129,14 @@ def optimal_mu(self, bounds, optimization_restart=10): def min_obj(X): return -1 * np.linalg.norm(self.predict(X.reshape(1, -1), True)[1]) - initial_starts = np.random.uniform(bounds[:, 0], - bounds[:, 1], - size=(optimization_restart, dim)) + initial_starts = np.random.uniform( + bounds[:, 0], bounds[:, 1], size=(optimization_restart, dim) + ) # Find the best optimum by starting from n_restart different random # points. for x0 in initial_starts: - res = minimize(min_obj, x0, bounds=bounds, method='L-BFGS-B') + res = minimize(min_obj, x0, bounds=bounds, method="L-BFGS-B") if res.fun < min_val: min_val = res.fun diff --git a/ezyrb/approximation/kneighbors_regressor.py b/ezyrb/approximation/kneighbors_regressor.py index 44c2e945..0776289e 100644 --- a/ezyrb/approximation/kneighbors_regressor.py +++ b/ezyrb/approximation/kneighbors_regressor.py @@ -14,9 +14,9 @@ class KNeighborsRegressor(NeighborsRegressor): :param kwargs: arguments passed to the internal instance of KNeighborsRegressor. - + :Example: - + >>> import numpy as np >>> from ezyrb import KNeighborsRegressor >>> x = np.random.uniform(-1, 1, size=(20, 2)) @@ -26,12 +26,12 @@ class KNeighborsRegressor(NeighborsRegressor): >>> new_x = np.array([[0.5, 0.5]]) >>> y_pred = knn.predict(new_x) """ + def __init__(self, **kwargs): """ Initialize a K-Neighbors Regressor. - + :param kwargs: Arguments passed to sklearn's KNeighborsRegressor. """ - logger.debug("Initializing KNeighborsRegressor with kwargs: %s", - kwargs) + logger.debug("Initializing KNeighborsRegressor with kwargs: %s", kwargs) self.regressor = Regressor(**kwargs) diff --git a/ezyrb/approximation/linear.py b/ezyrb/approximation/linear.py index ca27b02e..31b1f028 100644 --- a/ezyrb/approximation/linear.py +++ b/ezyrb/approximation/linear.py @@ -20,10 +20,11 @@ class Linear(Approximation): of the convex hull of the input points. If not provided, then the default is numpy.nan. """ + def __init__(self, fill_value=np.nan): """ Initialize a Linear interpolator. - + :param float fill_value: Value for points outside the convex hull. Default is numpy.nan. """ @@ -38,33 +39,41 @@ def fit(self, points, values): :param array_like points: the coordinates of the points. :param array_like values: the values in the points. """ - logger.debug("Fitting Linear with points shape: %s, " - "values shape: %s", - np.array(points).shape, np.array(values).shape) + logger.debug( + "Fitting Linear with points shape: %s, " "values shape: %s", + np.array(points).shape, + np.array(values).shape, + ) # the first dimension is the list of parameters, the second one is # the dimensionality of each tuple of parameters (we look for # parameters of dimensionality one) as_np_array = np.array(points) if not np.issubdtype(as_np_array.dtype, np.number): logger.error("Invalid points format/dimension") - raise ValueError('Invalid format or dimension for the argument' - '`points`.') + raise ValueError( + "Invalid format or dimension for the argument" "`points`." + ) if as_np_array.shape[-1] == 1: as_np_array = np.squeeze(as_np_array, axis=-1) logger.debug("Squeezed points array") - if as_np_array.ndim == 1 or (as_np_array.ndim == 2 - and as_np_array.shape[1] == 1): + if as_np_array.ndim == 1 or ( + as_np_array.ndim == 2 and as_np_array.shape[1] == 1 + ): logger.debug("Using 1D interpolation") - self.interpolator = interp1d(as_np_array, values, axis=0, - bounds_error=False, - fill_value=self.fill_value) + self.interpolator = interp1d( + as_np_array, + values, + axis=0, + bounds_error=False, + fill_value=self.fill_value, + ) else: logger.debug("Using ND interpolation") - self.interpolator = LinearNDInterp(points, - values, - fill_value=self.fill_value) + self.interpolator = LinearNDInterp( + points, values, fill_value=self.fill_value + ) logger.info("Linear fitted successfully") def predict(self, new_point): diff --git a/ezyrb/approximation/neighbors_regressor.py b/ezyrb/approximation/neighbors_regressor.py index 5fab2195..79494963 100644 --- a/ezyrb/approximation/neighbors_regressor.py +++ b/ezyrb/approximation/neighbors_regressor.py @@ -7,14 +7,14 @@ class NeighborsRegressor(Approximation): """ A generic superclass for wrappers of *NeighborsRegressor from sklearn. - + This class provides a common interface for neighbor-based regression methods. :param kwargs: Arguments passed to the internal instance of *NeighborsRegressor. - + :Example: - + >>> import numpy as np >>> from ezyrb import KNeighborsRegressor >>> x = np.random.uniform(-1, 1, size=(20, 2)) @@ -23,6 +23,7 @@ class NeighborsRegressor(Approximation): >>> knn.fit(x, y) >>> y_pred = knn.predict(x[:5]) """ + def fit(self, points, values): """ Construct the interpolator given `points` and `values`. diff --git a/ezyrb/approximation/radius_neighbors_regressor.py b/ezyrb/approximation/radius_neighbors_regressor.py index 26cfb100..edf172d4 100644 --- a/ezyrb/approximation/radius_neighbors_regressor.py +++ b/ezyrb/approximation/radius_neighbors_regressor.py @@ -11,9 +11,9 @@ class RadiusNeighborsRegressor(NeighborsRegressor): :param kwargs: arguments passed to the internal instance of RadiusNeighborsRegressor. - + :Example: - + >>> import numpy as np >>> from ezyrb import RadiusNeighborsRegressor >>> x = np.random.uniform(-1, 1, size=(20, 2)) @@ -23,10 +23,11 @@ class RadiusNeighborsRegressor(NeighborsRegressor): >>> new_x = np.array([[0.0, 0.0]]) >>> y_pred = rnn.predict(new_x) """ + def __init__(self, **kwargs): """ Initialize a Radius Neighbors Regressor. - + :param kwargs: Arguments passed to sklearn's RadiusNeighborsRegressor. """ self.regressor = Regressor(**kwargs) diff --git a/ezyrb/approximation/rbf.py b/ezyrb/approximation/rbf.py index a1d6f303..3c0b9986 100644 --- a/ezyrb/approximation/rbf.py +++ b/ezyrb/approximation/rbf.py @@ -42,15 +42,24 @@ class RBF(Approximation): >>> print(np.allclose(y, y_pred)) """ - def __init__(self, - kernel='thin_plate_spline', - smooth=0, - neighbors=None, - epsilon=None, - degree=None): - logger.debug("Initializing RBF with kernel=%s, smooth=%s, " - "neighbors=%s, epsilon=%s, degree=%s", - kernel, smooth, neighbors, epsilon, degree) + + def __init__( + self, + kernel="thin_plate_spline", + smooth=0, + neighbors=None, + epsilon=None, + degree=None, + ): + logger.debug( + "Initializing RBF with kernel=%s, smooth=%s, " + "neighbors=%s, epsilon=%s, degree=%s", + kernel, + smooth, + neighbors, + epsilon, + degree, + ) self.kernel = kernel self.smooth = smooth self.neighbors = neighbors @@ -66,8 +75,11 @@ def fit(self, points, values): :param array_like points: the coordinates of the points. :param array_like values: the values in the points. """ - logger.debug("Fitting RBF with points shape: %s, values shape: %s", - np.asarray(points).shape, np.asarray(values).shape) + logger.debug( + "Fitting RBF with points shape: %s, values shape: %s", + np.asarray(points).shape, + np.asarray(values).shape, + ) self.xi = np.asarray(points) if self.epsilon is None: @@ -78,8 +90,8 @@ def fit(self, points, values): ximin = np.amin(self.xi, axis=0) edges = ximax - ximin edges = edges[np.nonzero(edges)] - self.epsilon = np.power(np.prod(edges)/N, 1.0/edges.size) - if self.kernel in ['thin_plate_spline', 'cubic', 'quintic']: + self.epsilon = np.power(np.prod(edges) / N, 1.0 / edges.size) + if self.kernel in ["thin_plate_spline", "cubic", "quintic"]: self.epsilon = 1 logger.debug("Auto-computed epsilon: %f", self.epsilon) @@ -91,7 +103,8 @@ def fit(self, points, values): smoothing=self.smooth, kernel=self.kernel, epsilon=self.epsilon, - degree=self.degree) + degree=self.degree, + ) logger.info("RBF fitted successfully") def predict(self, new_point): diff --git a/ezyrb/approximation/sklearn_approximation.py b/ezyrb/approximation/sklearn_approximation.py index 0522a5af..0dcc715a 100644 --- a/ezyrb/approximation/sklearn_approximation.py +++ b/ezyrb/approximation/sklearn_approximation.py @@ -49,17 +49,13 @@ def __init__(self, sklearn_model, fit_params=None): """ logger.debug( "Initializing SklearnApproximation with model: %s", - type(sklearn_model).__name__ + type(sklearn_model).__name__, ) - if not hasattr(sklearn_model, 'fit'): - raise ValueError( - "sklearn_model must have a 'fit' method" - ) - if not hasattr(sklearn_model, 'predict'): - raise ValueError( - "sklearn_model must have a 'predict' method" - ) + if not hasattr(sklearn_model, "fit"): + raise ValueError("sklearn_model must have a 'fit' method") + if not hasattr(sklearn_model, "predict"): + raise ValueError("sklearn_model must have a 'predict' method") self.model = sklearn_model self.fit_params = fit_params if fit_params is not None else {} @@ -75,12 +71,10 @@ def fit(self, points, values): logger.info( "Fitting %s with %d samples", type(self.model).__name__, - points.shape[0] + points.shape[0], ) logger.debug( - "Input shape: %s, Output shape: %s", - points.shape, - values.shape + "Input shape: %s, Output shape: %s", points.shape, values.shape ) # Ensure 2D arrays @@ -103,13 +97,11 @@ def predict(self, new_points): :rtype: numpy.ndarray """ if not self._fitted: - raise RuntimeError( - "Model must be fitted before calling predict()" - ) + raise RuntimeError("Model must be fitted before calling predict()") logger.debug( "Predicting for %d new points", - new_points.shape[0] if new_points.ndim > 1 else 1 + new_points.shape[0] if new_points.ndim > 1 else 1, ) # Ensure 2D array @@ -118,7 +110,8 @@ def predict(self, new_points): predictions = self.model.predict(new_points) - logger.debug("Prediction completed, output shape: %s", - predictions.shape) + logger.debug( + "Prediction completed, output shape: %s", predictions.shape + ) return predictions diff --git a/ezyrb/database.py b/ezyrb/database.py index 85397670..f9810073 100644 --- a/ezyrb/database.py +++ b/ezyrb/database.py @@ -9,7 +9,7 @@ logger = logging.getLogger(__name__) -class Database(): +class Database: """ Database class for storing parameter-snapshot pairs. @@ -20,9 +20,9 @@ class Database(): :param Scale scaler_snapshots: the scaler for the snapshots. Default is None meaning no scaling. :param array_like space: the input spatial data - + :Example: - + >>> import numpy as np >>> from ezyrb import Database, Parameter, Snapshot >>> params = np.array([[1.0, 2.0], [3.0, 4.0], [5.0, 6.0]]) @@ -35,10 +35,15 @@ class Database(): >>> print(db.snapshots_matrix.shape) (3, 100) """ + def __init__(self, parameters=None, snapshots=None, space=None): - logger.debug("Initializing Database with parameters=%s, " - "snapshots=%s, space=%s", - type(parameters), type(snapshots), type(space)) + logger.debug( + "Initializing Database with parameters=%s, " + "snapshots=%s, space=%s", + type(parameters), + type(snapshots), + type(space), + ) self._pairs = [] if parameters is None and snapshots is None: @@ -48,19 +53,25 @@ def __init__(self, parameters=None, snapshots=None, space=None): if parameters is None: n_snaps = len(snapshots) if snapshots is not None else 0 parameters = [None] * n_snaps - logger.debug("Parameters were None, created %d None parameters", - n_snaps) + logger.debug( + "Parameters were None, created %d None parameters", n_snaps + ) elif snapshots is None: n_params = len(parameters) if parameters is not None else 0 snapshots = [None] * n_params - logger.debug("Snapshots were None, created %d None snapshots", - n_params) + logger.debug( + "Snapshots were None, created %d None snapshots", n_params + ) if len(parameters) != len(snapshots): - logger.error("Mismatch: %d parameters vs %d snapshots", - len(parameters), len(snapshots)) - raise ValueError('parameters and snapshots must have the same ' - 'length') + logger.error( + "Mismatch: %d parameters vs %d snapshots", + len(parameters), + len(snapshots), + ) + raise ValueError( + "parameters and snapshots must have the same " "length" + ) logger.debug("Adding %d parameter-snapshot pairs", len(parameters)) for param, snap in zip(parameters, snapshots): @@ -122,9 +133,10 @@ def __len__(self): return len(self._pairs) def __str__(self): - """ Print minimal info about the Database """ - s = 'Database with {} snapshots and {} parameters'.format( - self.snapshots_matrix.shape[1], self.parameters_matrix.shape[1]) + """Print minimal info about the Database""" + s = "Database with {} snapshots and {} parameters".format( + self.snapshots_matrix.shape[1], self.parameters_matrix.shape[1] + ) return s def add(self, parameter, snapshot): @@ -144,8 +156,9 @@ def add(self, parameter, snapshot): raise ValueError self._pairs.append((parameter, snapshot)) - logger.debug("Added parameter-snapshot pair. Total pairs: %d", - len(self._pairs)) + logger.debug( + "Added parameter-snapshot pair. Total pairs: %d", len(self._pairs) + ) return self @@ -157,8 +170,7 @@ def split(self, chunks, seed=None): >>> train, test = db.split([80, 20]) # n snapshots """ - logger.debug("Splitting database with chunks=%s, seed=%s", - chunks, seed) + logger.debug("Splitting database with chunks=%s, seed=%s", chunks, seed) if seed is not None: np.random.seed(seed) @@ -166,30 +178,31 @@ def split(self, chunks, seed=None): if all(isinstance(n, int) for n in chunks): if sum(chunks) != len(self): - logger.error("Sum of chunks %d != database size %d", - sum(chunks), len(self)) - raise ValueError('chunk elements are inconsistent') + logger.error( + "Sum of chunks %d != database size %d", + sum(chunks), + len(self), + ) + raise ValueError("chunk elements are inconsistent") logger.debug("Splitting by absolute numbers: %s", chunks) - ids = [ - j for j, chunk in enumerate(chunks) - for i in range(chunk) - ] + ids = [j for j, chunk in enumerate(chunks) for i in range(chunk)] np.random.shuffle(ids) elif all(isinstance(n, float) for n in chunks): - if not np.isclose(sum(chunks), 1.): + if not np.isclose(sum(chunks), 1.0): logger.error("Sum of chunk ratios %f != 1.0", sum(chunks)) - raise ValueError('chunk elements are inconsistent') + raise ValueError("chunk elements are inconsistent") logger.debug("Splitting by ratios: %s", chunks) cum_chunks = np.cumsum(chunks) cum_chunks = np.insert(cum_chunks, 0, 0.0) - ids = np.ones(len(self)) * -1. + ids = np.ones(len(self)) * -1.0 tmp = np.random.uniform(0, 1, size=len(self)) - for i in range(len(cum_chunks)-1): + for i in range(len(cum_chunks) - 1): is_between = np.logical_and( - tmp >= cum_chunks[i], tmp < cum_chunks[i+1]) + tmp >= cum_chunks[i], tmp < cum_chunks[i + 1] + ) ids[is_between] = i else: @@ -202,12 +215,14 @@ def split(self, chunks, seed=None): for p, s in np.asarray(self._pairs)[chunk_ids]: new_database[i].add(p, s) - logger.info("Database split into %d parts with sizes: %s", - len(new_database), - [len(db) for db in new_database]) + logger.info( + "Database split into %d parts with sizes: %s", + len(new_database), + [len(db) for db in new_database], + ) return new_database - + def get_snapshot_space(self, index): """ Get the space coordinates of a snapshot by its index. diff --git a/ezyrb/parallel/__init__.py b/ezyrb/parallel/__init__.py index 4b84a5b5..0d60b177 100644 --- a/ezyrb/parallel/__init__.py +++ b/ezyrb/parallel/__init__.py @@ -1,9 +1,18 @@ """EZyRB package""" __all__ = [ - 'Reduction', 'POD', 'Approximation', 'RBF', 'Linear', 'GPR', - 'ANN', 'KNeighborsRegressor', 'RadiusNeighborsRegressor', 'AE', 'AE_EDDL', - 'ReducedOrderModel' + "Reduction", + "POD", + "Approximation", + "RBF", + "Linear", + "GPR", + "ANN", + "KNeighborsRegressor", + "RadiusNeighborsRegressor", + "AE", + "AE_EDDL", + "ReducedOrderModel", ] from .reduction import Reduction diff --git a/ezyrb/parallel/ae.py b/ezyrb/parallel/ae.py index 89592d5d..2353a320 100644 --- a/ezyrb/parallel/ae.py +++ b/ezyrb/parallel/ae.py @@ -10,6 +10,7 @@ from pycompss.api.parameter import INOUT, IN from .reduction import Reduction + class AE(Reduction, ANN): """ Feed-Forward AutoEncoder class (AE) @@ -51,15 +52,18 @@ class AE(Reduction, ANN): >>> reduced_snapshots = ae.reduce(snapshots) >>> expanded_snapshots = ae.expand(reduced_snapshots) """ - def __init__(self, - layers_encoder, - layers_decoder, - function_encoder, - function_decoder, - stop_training, - loss=None, - optimizer=torch.optim.Adam, - lr=0.001): + + def __init__( + self, + layers_encoder, + layers_decoder, + function_encoder, + function_decoder, + stop_training, + loss=None, + optimizer=torch.optim.Adam, + lr=0.001, + ): if loss is None: loss = torch.nn.MSELoss() @@ -92,6 +96,7 @@ class InnerAE(torch.nn.Module): """ Autoencoder as a pytorch module """ + def __init__(self, outer_instance): super().__init__() self.encoder = outer_instance.encoder @@ -126,18 +131,22 @@ def _build_model(self, values): layers_encoder_torch = [] for i in range(len(layers_encoder) - 2): layers_encoder_torch.append( - nn.Linear(layers_encoder[i], layers_encoder[i + 1])) + nn.Linear(layers_encoder[i], layers_encoder[i + 1]) + ) layers_encoder_torch.append(self.function_encoder[i]) layers_encoder_torch.append( - nn.Linear(layers_encoder[-2], layers_encoder[-1])) + nn.Linear(layers_encoder[-2], layers_encoder[-1]) + ) layers_decoder_torch = [] for i in range(len(layers_decoder) - 2): layers_decoder_torch.append( - nn.Linear(layers_decoder[i], layers_decoder[i + 1])) + nn.Linear(layers_decoder[i], layers_decoder[i + 1]) + ) layers_decoder_torch.append(self.function_decoder[i]) layers_decoder_torch.append( - nn.Linear(layers_decoder[-2], layers_decoder[-1])) + nn.Linear(layers_decoder[-2], layers_decoder[-1]) + ) self.encoder = nn.Sequential(*layers_encoder_torch) self.decoder = nn.Sequential(*layers_decoder_torch) # Creating the model adding the encoder and the decoder @@ -215,7 +224,8 @@ def inverse_transform(self, g, database): if database and database.scaler_snapshots: predicted_sol = database.scaler_snapshots.inverse_transform( - predicted_sol.T).T + predicted_sol.T + ).T if 1 in predicted_sol.shape: predicted_sol = predicted_sol.ravel() diff --git a/ezyrb/parallel/ae_eddl.py b/ezyrb/parallel/ae_eddl.py index 38353715..9c47f2bd 100644 --- a/ezyrb/parallel/ae_eddl.py +++ b/ezyrb/parallel/ae_eddl.py @@ -2,8 +2,9 @@ Module for FNN-Autoencoders. """ -import sys # Path to find PyEDDL library in PyCOMPSs container -sys.path.append('/usr/local/miniconda3/lib/python3.8/site-packages') +import sys # Path to find PyEDDL library in PyCOMPSs container + +sys.path.append("/usr/local/miniconda3/lib/python3.8/site-packages") from pyeddl import eddl from pyeddl.tensor import Tensor import numpy as np @@ -11,6 +12,7 @@ from pycompss.api.parameter import INOUT, IN from .reduction import Reduction + class AE_EDDL(Reduction): """ Feed-Forward AutoEncoder class (AE) @@ -61,19 +63,22 @@ class AE_EDDL(Reduction): >>> reduced_snapshots = ae.reduce(snapshots) >>> expanded_snapshots = ae.expand(reduced_snapshots) """ - def __init__(self, - layers_encoder, - layers_decoder, - function_encoder, - function_decoder, - stop_training, - batch_size, - loss=None, - metric=None, - optimizer=eddl.adam, - lr=0.001, - cs = eddl.CS_CPU, - training_type = 2): + + def __init__( + self, + layers_encoder, + layers_decoder, + function_encoder, + function_decoder, + stop_training, + batch_size, + loss=None, + metric=None, + optimizer=eddl.adam, + lr=0.001, + cs=eddl.CS_CPU, + training_type=2, + ): if loss is None: loss = "mean_squared_error" @@ -131,11 +136,11 @@ def __getstate__(self): state = self.__dict__.copy() # remove unpicklable entries - del state['encoder'] - del state['decoder'] - del state['decoder2'] - del state['model_Autoencoder'] - del state['model_Decoder'] + del state["encoder"] + del state["decoder"] + del state["decoder2"] + del state["model_Autoencoder"] + del state["model_Decoder"] return state def __setstate__(self, state): @@ -154,12 +159,12 @@ def __setstate__(self, state): self.__dict__.update(state) # restore unpicklable entries - if self.fitted: # inputs for all tasks except the fit() task + if self.fitted: # inputs for all tasks except the fit() task self.imported = True # (n) hidden layers + (n-1) activation layers + (1) input layer n = len(self.layers_encoder) - encoder_layer_index = 2*n - 1 + encoder_layer_index = 2 * n - 1 self.model_Autoencoder = eddl.import_net_from_onnx_file(self.file_1) self.encoder = self.model_Autoencoder.layers[encoder_layer_index] self.decoder = self.model_Autoencoder.layers[-1] @@ -167,31 +172,31 @@ def __setstate__(self, state): eddl.build( self.model_Autoencoder, self.optimizer(self.lr), - [self.loss ], + [self.loss], [self.metric], self.cs(mem="low_mem"), - False # don't initialize random weights (using trained model) + False, # don't initialize random weights (using trained model) ) # resize manually since we don't use "fit" --> # size of each layer = (batch_size, layer_dof) self.model_Autoencoder.resize(self.batch_size) - #------------------------------------------------------------------- + # ------------------------------------------------------------------- self.model_Decoder = eddl.import_net_from_onnx_file(self.file_2) self.decoder2 = self.model_Decoder.layers[-1] eddl.build( self.model_Decoder, self.optimizer(self.lr), - [self.loss ], + [self.loss], [self.metric], self.cs(mem="low_mem"), - False # don't initialize random weights (using trained model) + False, # don't initialize random weights (using trained model) ) # resize manually since we don't use "fit" --> # size of each layer = (batch_size, layer_dof) self.model_Decoder.resize(self.batch_size) - else: # inputs for the fit() method + else: # inputs for the fit() method self.encoder = None self.decoder = None self.decoder2 = None @@ -217,19 +222,22 @@ def _build_model(self, values): layers_decoder.append(values.shape[1]) def EncoderBlock(layer): - for i in range(1, len(layers_encoder)-1): + for i in range(1, len(layers_encoder) - 1): layer = self.function_encoder[i]( - eddl.Dense(layer, layers_encoder[i])) + eddl.Dense(layer, layers_encoder[i]) + ) layer = eddl.Dense(layer, layers_encoder[-1]) return layer def DecoderBlock(layer): - for i in range(1, len(layers_decoder)-1): + for i in range(1, len(layers_decoder) - 1): layer = self.function_encoder[i]( - eddl.Dense(layer, layers_decoder[i])) + eddl.Dense(layer, layers_decoder[i]) + ) layer = eddl.Dense(layer, layers_decoder[-1]) - return layer - #----------------------------------------------------------------------- + return layer + + # ----------------------------------------------------------------------- # Define network1 in1 = eddl.Input([layers_encoder[0]]) self.encoder = EncoderBlock(in1) @@ -239,31 +247,31 @@ def DecoderBlock(layer): eddl.build( self.model_Autoencoder, self.optimizer(self.lr), - [self.loss ], + [self.loss], [self.metric], self.cs(mem="low_mem"), - True # initialize weights to random values + True, # initialize weights to random values ) - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- # Define network2 in2 = eddl.Input([layers_decoder[0]]) self.decoder2 = DecoderBlock(in2) self.model_Decoder = eddl.Model([in2], [self.decoder2]) - + eddl.build( self.model_Decoder, self.optimizer(self.lr), - [self.loss ], + [self.loss], [self.metric], self.cs(mem="low_mem"), - True # initialize weights to random values + True, # initialize weights to random values ) - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- eddl.summary(self.model_Autoencoder) eddl.plot(self.model_Autoencoder, "Autoencoder.pdf") eddl.summary(self.model_Decoder) eddl.plot(self.model_Decoder, "Decoder.pdf") - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- @task(target_direction=INOUT) def fit(self, values): @@ -285,17 +293,22 @@ def fit(self, values): :param numpy.ndarray values: the (training) values in the points. """ values = values.T - values = Tensor.fromarray(values) # Numpy array to EDDL.Tensor + values = Tensor.fromarray(values) # Numpy array to EDDL.Tensor self._build_model(values) - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- if self.training_type == 1: - print('Coarse training:') - n_epoch =1 + print("Coarse training:") + n_epoch = 1 flag = True while flag: - print('Iteration number: ',n_epoch) - eddl.fit(self.model_Autoencoder, [values], [values], - self.batch_size, 1) + print("Iteration number: ", n_epoch) + eddl.fit( + self.model_Autoencoder, + [values], + [values], + self.batch_size, + 1, + ) losses = eddl.get_losses(self.model_Autoencoder) metrics = eddl.get_metrics(self.model_Autoencoder) self.loss_trend.append(losses) @@ -311,13 +324,13 @@ def fit(self, values): if losses[0] < criteria: flag = False n_epoch += 1 - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- elif self.training_type == 2: - print('Fine training_1:') + print("Fine training_1:") s = values.shape num_batches = s[0] // self.batch_size xbatch = Tensor([self.batch_size, 3067]) - n_epoch =1 + n_epoch = 1 flag = True while flag: eddl.reset_loss(self.model_Autoencoder) @@ -335,7 +348,7 @@ def fit(self, values): # print("Epoch %d (%d batches)" % (n_epoch, num_batches)) # for l, m in zip(losses, metrics): # print("Loss: %.6f\tMetric: %.6f" % (l, m)) - + for criteria in self.stop_training: if isinstance(criteria, int): # stop criteria is an integer @@ -346,20 +359,21 @@ def fit(self, values): if losses[0] < criteria: flag = False n_epoch += 1 - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- elif self.training_type == 3: - print('Fine training_2:') + print("Fine training_2:") s = values.shape num_batches = s[0] // self.batch_size - n_epoch =1 + n_epoch = 1 flag = True while flag: eddl.reset_loss(self.model_Autoencoder) for j in range(num_batches): # 2) using samples indices indices = np.random.randint(0, s[0], self.batch_size) - eddl.train_batch(self.model_Autoencoder, [values], [values], - indices) + eddl.train_batch( + self.model_Autoencoder, [values], [values], indices + ) losses = eddl.get_losses(self.model_Autoencoder) metrics = eddl.get_metrics(self.model_Autoencoder) @@ -370,7 +384,7 @@ def fit(self, values): # print("Epoch {} ({} batches)".format(n_epoch, num_batches)) # for l, m in zip(losses, metrics): # print("Loss: %.6f\tMetric: %.6f" % (l, m)) - + for criteria in self.stop_training: if isinstance(criteria, int): # stop criteria is an integer @@ -381,20 +395,21 @@ def fit(self, values): if losses[0] < criteria: flag = False n_epoch += 1 - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- # Copy parameters from model_Autoencoder to model_Decoder # munber of Decoder layers = (n) hidden layers + (n-1) activation layers n = len(self.layers_encoder) - num_lay_decoder = 2*n - 1 - decoder_parameters = eddl.get_parameters(self.model_Autoencoder, - True)[-num_lay_decoder:] + num_lay_decoder = 2 * n - 1 + decoder_parameters = eddl.get_parameters(self.model_Autoencoder, True)[ + -num_lay_decoder: + ] # Insert empty parameter for the new input layer of decoder - decoder_parameters.insert(0,[]) + decoder_parameters.insert(0, []) eddl.set_parameters(self.model_Decoder, decoder_parameters) ## For debugging # for i in decoder_parameters: # print(len(i)) - #----------------------------------------------------------------------- + # ----------------------------------------------------------------------- # For PyCOMPSs: objects used as task parameters must be automatically # serializable/picklable so we save the trained model and delete all # unpicklables at the time of serialization @@ -409,19 +424,19 @@ def transform(self, X, scaler_red): :param numpy.ndarray X: the input snapshots matrix (stored by column). """ - #----------------------------------------------------------------------- - if self.imported: # Means PyCOMPSs is used. + # ----------------------------------------------------------------------- + if self.imported: # Means PyCOMPSs is used. # # For debugging # print("Encoder layer {} --> {}".format( # self.encoder.input.shape, self.encoder.output.shape)) print(f"Trained model imported from ({self.file_1})") eddl.summary(self.model_Autoencoder) - #----------------------------------------------------------------------- - X = Tensor.fromarray(X.T) # Numpy array to EDDL.Tensor + # ----------------------------------------------------------------------- + X = Tensor.fromarray(X.T) # Numpy array to EDDL.Tensor # One prediction for the fitted model(1 forward pass after training) eddl.predict(self.model_Autoencoder, [X]) g = eddl.getOutput(self.encoder) - reduced_output = (Tensor.getdata(g).T).T # EDDL.Tensor to Numpy array + reduced_output = (Tensor.getdata(g).T).T # EDDL.Tensor to Numpy array if scaler_red: reduced_output = scaler_red.fit_transform(reduced_output) @@ -441,12 +456,12 @@ def inverse_transform(self, g, database): :param: numpy.ndarray g the latent variables. """ - #----------------------------------------------------------------------- - if self.imported: # Means PyCOMPSs is used. + # ----------------------------------------------------------------------- + if self.imported: # Means PyCOMPSs is used. print(f"Trained model imported from ({self.file_2})") eddl.summary(self.model_Decoder) - #----------------------------------------------------------------------- - g = Tensor.fromarray(g.T) # Numpy array to EDDL.Tensor + # ----------------------------------------------------------------------- + g = Tensor.fromarray(g.T) # Numpy array to EDDL.Tensor # One forward pass for the new decoder (without training, parameters # copied from the trained autoencoder) eddl.forward(self.model_Decoder, [g]) @@ -456,17 +471,18 @@ def inverse_transform(self, g, database): # print('Before expansion', g.shape) # print('Before expansion', u.shape) - predicted_sol = Tensor.getdata(u).T # EDDL.Tensor to Numpy array - + predicted_sol = Tensor.getdata(u).T # EDDL.Tensor to Numpy array + if database and database.scaler_snapshots: predicted_sol = database.scaler_snapshots.inverse_transform( - predicted_sol.T).T + predicted_sol.T + ).T if 1 in predicted_sol.shape: predicted_sol = predicted_sol.ravel() else: predicted_sol = predicted_sol.T - + return predicted_sol def reduce(self, X, scaler_red): diff --git a/ezyrb/parallel/ann.py b/ezyrb/parallel/ann.py index fe2cfbb8..396d36dc 100644 --- a/ezyrb/parallel/ann.py +++ b/ezyrb/parallel/ann.py @@ -9,6 +9,7 @@ from pycompss.api.parameter import INOUT, IN from .approximation import Approximation + class ANN(Approximation): """ Feed-Forward Artifical Neural Network class (ANN). @@ -39,6 +40,7 @@ class ANN(Approximation): >>> print(len(ann.loss_trend)) >>> print(ann.loss_trend[-1]) """ + def __init__(self, layers, function, stop_training, loss=None): if loss is None: @@ -162,7 +164,6 @@ def predict(self, new_point, scaler_red): y_new = self.model(new_point) predicted_red_sol = np.atleast_2d(self._convert_torch_to_numpy(y_new)) if scaler_red: # rescale modal coefficients - predicted_red_sol = scaler_red.inverse_transform( - predicted_red_sol) + predicted_red_sol = scaler_red.inverse_transform(predicted_red_sol) predicted_red_sol = predicted_red_sol.T - return predicted_red_sol \ No newline at end of file + return predicted_red_sol diff --git a/ezyrb/parallel/approximation.py b/ezyrb/parallel/approximation.py index fd4c297e..ae87881e 100644 --- a/ezyrb/parallel/approximation.py +++ b/ezyrb/parallel/approximation.py @@ -10,6 +10,7 @@ class Approximation(ABC): All the classes that implement the input-output mapping should be inherited from this class. """ + @abstractmethod def fit(self, points, values): """Abstract `fit`""" diff --git a/ezyrb/parallel/gpr.py b/ezyrb/parallel/gpr.py index 61e9927f..460c59b2 100644 --- a/ezyrb/parallel/gpr.py +++ b/ezyrb/parallel/gpr.py @@ -1,6 +1,7 @@ """ Module wrapper exploiting `GPy` for Gaussian Process Regression """ + import numpy as np from scipy.optimize import minimize from sklearn.gaussian_process import GaussianProcessRegressor @@ -9,6 +10,7 @@ from pycompss.api.parameter import INOUT, IN from .approximation import Approximation + class GPR(Approximation): """ Multidimensional regression using Gaussian process. @@ -31,18 +33,21 @@ class GPR(Approximation): >>> print(np.allclose(y, y_pred)) """ + def __init__(self): self.X_sample = None self.Y_sample = None self.model = None @task(target_direction=INOUT) - def fit(self, - points, - values, - kern=None, - normalizer=True, - optimization_restart=20): + def fit( + self, + points, + values, + kern=None, + normalizer=True, + optimization_restart=20, + ): """ Construct the regression given `points` and `values`. @@ -63,8 +68,10 @@ def fit(self, self.Y_sample = self.Y_sample.reshape(-1, 1) self.model = GaussianProcessRegressor( - kernel=kern, n_restarts_optimizer=optimization_restart, - normalize_y=normalizer) + kernel=kern, + n_restarts_optimizer=optimization_restart, + normalize_y=normalizer, + ) self.model.fit(self.X_sample, self.Y_sample) @task(returns=np.ndarray, target_direction=IN) @@ -80,11 +87,11 @@ def predict(self, new_points, scaler_red, return_variance=False): :rtype: (numpy.ndarray, numpy.ndarray) """ new_points = np.atleast_2d(new_points) - predicted_red_sol = np.atleast_2d(self.model.predict(new_points, - return_std=return_variance)) + predicted_red_sol = np.atleast_2d( + self.model.predict(new_points, return_std=return_variance) + ) if scaler_red: # rescale modal coefficients - predicted_red_sol = scaler_red.inverse_transform( - predicted_red_sol) + predicted_red_sol = scaler_red.inverse_transform(predicted_red_sol) predicted_red_sol = predicted_red_sol.T return predicted_red_sol @@ -107,14 +114,14 @@ def optimal_mu(self, bounds, optimization_restart=10): def min_obj(X): return -1 * np.linalg.norm(self.predict(X.reshape(1, -1), True)[1]) - initial_starts = np.random.uniform(bounds[:, 0], - bounds[:, 1], - size=(optimization_restart, dim)) + initial_starts = np.random.uniform( + bounds[:, 0], bounds[:, 1], size=(optimization_restart, dim) + ) # Find the best optimum by starting from n_restart different random # points. for x0 in initial_starts: - res = minimize(min_obj, x0, bounds=bounds, method='L-BFGS-B') + res = minimize(min_obj, x0, bounds=bounds, method="L-BFGS-B") if res.fun < min_val: min_val = res.fun diff --git a/ezyrb/parallel/kneighbors_regressor.py b/ezyrb/parallel/kneighbors_regressor.py index 60afdcee..0b86668b 100644 --- a/ezyrb/parallel/kneighbors_regressor.py +++ b/ezyrb/parallel/kneighbors_regressor.py @@ -12,5 +12,6 @@ class KNeighborsRegressor(NeighborsRegressor): :param kwargs: arguments passed to the internal instance of KNeighborsRegressor. """ + def __init__(self, **kwargs): self.regressor = Regressor(**kwargs) diff --git a/ezyrb/parallel/linear.py b/ezyrb/parallel/linear.py index 487624a8..2400d229 100644 --- a/ezyrb/parallel/linear.py +++ b/ezyrb/parallel/linear.py @@ -10,6 +10,7 @@ from pycompss.api.parameter import INOUT, IN from .approximation import Approximation + class Linear(Approximation): """ Multidimensional linear interpolator. @@ -18,6 +19,7 @@ class Linear(Approximation): of the convex hull of the input points. If not provided, then the default is numpy.nan. """ + def __init__(self, fill_value=np.nan): self.fill_value = fill_value self.interpolator = None @@ -35,19 +37,21 @@ def fit(self, points, values): # parameters of dimensionality one) as_np_array = np.array(points) if not np.issubdtype(as_np_array.dtype, np.number): - raise ValueError('Invalid format or dimension for the argument' - '`points`.') + raise ValueError( + "Invalid format or dimension for the argument" "`points`." + ) if as_np_array.shape[-1] == 1: as_np_array = np.squeeze(as_np_array, axis=-1) - if as_np_array.ndim == 1 or (as_np_array.ndim == 2 - and as_np_array.shape[1] == 1): + if as_np_array.ndim == 1 or ( + as_np_array.ndim == 2 and as_np_array.shape[1] == 1 + ): self.interpolator = interp1d(as_np_array, values, axis=0) else: - self.interpolator = LinearNDInterp(points, - values, - fill_value=self.fill_value) + self.interpolator = LinearNDInterp( + points, values, fill_value=self.fill_value + ) @task(returns=np.ndarray, target_direction=IN) def predict(self, new_point, scaler_red): @@ -68,8 +72,7 @@ def predict(self, new_point, scaler_red): predicted_red_sol = np.atleast_2d(new_red_snap) if scaler_red: # rescale modal coefficients - predicted_red_sol = scaler_red.inverse_transform( - predicted_red_sol) + predicted_red_sol = scaler_red.inverse_transform(predicted_red_sol) predicted_red_sol = predicted_red_sol.T return predicted_red_sol diff --git a/ezyrb/parallel/neighbors_regressor.py b/ezyrb/parallel/neighbors_regressor.py index 0bf872ee..a91ca553 100644 --- a/ezyrb/parallel/neighbors_regressor.py +++ b/ezyrb/parallel/neighbors_regressor.py @@ -5,6 +5,7 @@ from pycompss.api.parameter import INOUT, IN from .approximation import Approximation + class NeighborsRegressor(Approximation): """ A generic superclass for wrappers of *NeighborsRegressor from sklearn. @@ -12,6 +13,7 @@ class NeighborsRegressor(Approximation): :param kwargs: arguments passed to the internal instance of *NeighborsRegressor. """ + @task(target_direction=INOUT) def fit(self, points, values): """ @@ -41,7 +43,6 @@ def predict(self, new_point, scaler_red): predicted_red_sol = np.atleast_2d(self.regressor.predict(new_point)) if scaler_red: # rescale modal coefficients - predicted_red_sol = scaler_red.inverse_transform( - predicted_red_sol) + predicted_red_sol = scaler_red.inverse_transform(predicted_red_sol) predicted_red_sol = predicted_red_sol.T return predicted_red_sol diff --git a/ezyrb/parallel/pod.py b/ezyrb/parallel/pod.py index 26901885..40ad8c62 100644 --- a/ezyrb/parallel/pod.py +++ b/ezyrb/parallel/pod.py @@ -4,6 +4,7 @@ Decomposition, Truncated Randomized Singular Value Decomposition, Truncated Singular Value Decomposition via correlation matrix. """ + try: from scipy.linalg import eigh except ImportError: @@ -14,6 +15,7 @@ from pycompss.api.parameter import INOUT, IN from .reduction import Reduction + class POD(Reduction): """ Perform the Proper Orthogonal Decomposition. @@ -53,21 +55,19 @@ class POD(Reduction): omega_rank=10) >>> pod = POD('correlation_matrix', rank=10, save_memory=False) """ - def __init__(self, method='svd', **kwargs): + + def __init__(self, method="svd", **kwargs): available_methods = { - 'svd': (self._svd, { - 'rank': -1 - }), - 'randomized_svd': (self._rsvd, { - 'rank': -1, - 'subspace_iteration': 1, - 'omega_rank': 0 - }), - 'correlation_matrix': (self._corrm, { - 'rank': -1, - 'save_memory': False - }), + "svd": (self._svd, {"rank": -1}), + "randomized_svd": ( + self._rsvd, + {"rank": -1, "subspace_iteration": 1, "omega_rank": 0}, + ), + "correlation_matrix": ( + self._corrm, + {"rank": -1, "save_memory": False}, + ), } self._modes = None @@ -77,7 +77,9 @@ def __init__(self, method='svd', **kwargs): if method is None: raise RuntimeError( "Invalid method for POD. Please chose one among {}".format( - ', '.join(available_methods))) + ", ".join(available_methods) + ) + ) self.__method, args = method args.update(kwargs) @@ -137,7 +139,8 @@ def inverse_transform(self, X, database): if database and database.scaler_snapshots: predicted_sol = database.scaler_snapshots.inverse_transform( - predicted_sol.T).T + predicted_sol.T + ).T if 1 in predicted_sol.shape: predicted_sol = predicted_sol.ravel() @@ -181,6 +184,7 @@ def _truncation(self, X, s): :return: the number of modes :rtype: int """ + def omega(x): return 0.56 * x**3 - 0.95 * x**2 + 1.82 * x + 1.43 @@ -231,8 +235,11 @@ def _rsvd(self, X): constructing approximate matrix decompositions. N. Halko, P. G. Martinsson, J. A. Tropp. """ - if (self.omega_rank == 0 and isinstance(self.rank, int) - and self.rank not in [0, -1]): + if ( + self.omega_rank == 0 + and isinstance(self.rank, int) + and self.rank not in [0, -1] + ): omega_rank = self.rank * 2 elif self.omega_rank == 0: omega_rank = X.shape[1] * 2 diff --git a/ezyrb/parallel/radius_neighbors_regressor.py b/ezyrb/parallel/radius_neighbors_regressor.py index a649b4db..c2379070 100644 --- a/ezyrb/parallel/radius_neighbors_regressor.py +++ b/ezyrb/parallel/radius_neighbors_regressor.py @@ -12,5 +12,6 @@ class RadiusNeighborsRegressor(NeighborsRegressor): :param kwargs: arguments passed to the internal instance of RadiusNeighborsRegressor. """ + def __init__(self, **kwargs): self.regressor = Regressor(**kwargs) diff --git a/ezyrb/parallel/rbf.py b/ezyrb/parallel/rbf.py index 168f112b..836de014 100644 --- a/ezyrb/parallel/rbf.py +++ b/ezyrb/parallel/rbf.py @@ -6,6 +6,7 @@ from pycompss.api.parameter import INOUT, IN from .approximation import Approximation + class RBF(Approximation): """ Multidimensional interpolator using Radial Basis Function. @@ -40,12 +41,15 @@ class RBF(Approximation): >>> print(np.allclose(y, y_pred)) """ - def __init__(self, - kernel='thin_plate_spline', - smooth=0, - neighbors=None, - epsilon=None, - degree=None): + + def __init__( + self, + kernel="thin_plate_spline", + smooth=0, + neighbors=None, + epsilon=None, + degree=None, + ): self.kernel = kernel self.smooth = smooth self.neighbors = neighbors @@ -72,8 +76,8 @@ def fit(self, points, values): ximin = np.amin(self.xi, axis=0) edges = ximax - ximin edges = edges[np.nonzero(edges)] - self.epsilon = np.power(np.prod(edges)/N, 1.0/edges.size) - if self.kernel in ['thin_plate_spline', 'cubic', 'quintic']: + self.epsilon = np.power(np.prod(edges) / N, 1.0 / edges.size) + if self.kernel in ["thin_plate_spline", "cubic", "quintic"]: self.epsilon = 1 self.interpolator = RBFInterpolator( @@ -83,7 +87,8 @@ def fit(self, points, values): smoothing=self.smooth, kernel=self.kernel, epsilon=self.epsilon, - degree=self.degree) + degree=self.degree, + ) @task(returns=np.ndarray, target_direction=IN) def predict(self, new_point, scaler_red): @@ -94,10 +99,10 @@ def predict(self, new_point, scaler_red): :return: the interpolated values. :rtype: numpy.ndarray """ - predicted_red_sol = np.atleast_2d(self.interpolator( - np.asarray(new_point))) + predicted_red_sol = np.atleast_2d( + self.interpolator(np.asarray(new_point)) + ) if scaler_red: # rescale modal coefficients - predicted_red_sol = scaler_red.inverse_transform( - predicted_red_sol) + predicted_red_sol = scaler_red.inverse_transform(predicted_red_sol) predicted_red_sol = predicted_red_sol.T return predicted_red_sol diff --git a/ezyrb/parallel/reducedordermodel.py b/ezyrb/parallel/reducedordermodel.py index c65b58a9..fe3bb4b6 100644 --- a/ezyrb/parallel/reducedordermodel.py +++ b/ezyrb/parallel/reducedordermodel.py @@ -8,7 +8,8 @@ from sklearn.model_selection import KFold from pycompss.api.api import compss_wait_on -class ReducedOrderModel(): + +class ReducedOrderModel: """ Reduced Order Model class. @@ -45,6 +46,7 @@ class ReducedOrderModel(): >>> rom.predict(new_param) """ + def __init__(self, database, reduction, approximation, scaler_red=None): self.database = database self.reduction = reduction @@ -60,13 +62,12 @@ def fit(self, *args, **kwargs): """ self.reduction.fit(self.database.snapshots.T) reduced_output = self.reduction.transform( - self.database.snapshots.T, self.scaler_red) + self.database.snapshots.T, self.scaler_red + ) self.approximation.fit( - self.database.parameters, - reduced_output, - *args, - **kwargs) + self.database.parameters, reduced_output, *args, **kwargs + ) return self @@ -75,13 +76,14 @@ def predict(self, mu): Calculate predicted solution for given mu """ mu = np.atleast_2d(mu) - if hasattr(self, 'database') and self.database.scaler_parameters: + if hasattr(self, "database") and self.database.scaler_parameters: mu = self.database.scaler_parameters.transform(mu) predicted_red_sol = self.approximation.predict(mu, self.scaler_red) predicted_sol = self.reduction.inverse_transform( - predicted_red_sol, self.database) + predicted_red_sol, self.database + ) return predicted_sol @@ -113,7 +115,7 @@ def save(self, fname, save_db=True, save_reduction=True, save_approx=True): if not save_approx: del rom_to_store.approximation - with open(fname, 'wb') as output: + with open(fname, "wb") as output: pickle.dump(rom_to_store, output, pickle.HIGHEST_PROTOCOL) @staticmethod @@ -129,7 +131,7 @@ def load(fname): >>> rom = ROM.load('ezyrb.rom') >>> rom.predict(new_param) """ - with open(fname, 'rb') as output: + with open(fname, "rb") as output: rom = pickle.load(output) return rom @@ -152,14 +154,16 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, **kwargs): :rtype: numpy.ndarray """ error = [] - predicted_test = [] # to save my future objects - original_test= [] + predicted_test = [] # to save my future objects + original_test = [] kf = KFold(n_splits=n_splits) for train_index, test_index in kf.split(self.database): new_db = self.database[train_index] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation)).fit( - *args, **kwargs) + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + ).fit(*args, **kwargs) test = self.database[test_index] predicted_test.append(rom.predict(test.parameters)) @@ -167,9 +171,12 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, **kwargs): predicted_test = compss_wait_on(predicted_test) for j in range(len(predicted_test)): - error.append(np.mean( - norm(predicted_test[j] - original_test[j], axis=1) / - norm(original_test[j], axis=1))) + error.append( + np.mean( + norm(predicted_test[j] - original_test[j], axis=1) + / norm(original_test[j], axis=1) + ) + ) return np.array(error) @@ -194,8 +201,8 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): """ error = np.zeros(len(self.database)) db_range = list(range(len(self.database))) - predicted_test = [] # to save my future objects - original_test= [] + predicted_test = [] # to save my future objects + original_test = [] for j in db_range: indeces = np.array([True] * len(self.database)) @@ -203,9 +210,11 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): new_db = self.database[indeces] test_db = self.database[~indeces] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation)).fit( - *args, **kwargs) + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + ).fit(*args, **kwargs) predicted_test.append(rom.predict(test_db.parameters)) original_test.append(test_db.snapshots) @@ -213,8 +222,9 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): predicted_test = compss_wait_on(predicted_test) for j in range(len(predicted_test)): error[j] = np.mean( - norm(predicted_test[j] - original_test[j],axis=1) / - norm(original_test[j], axis=1)) + norm(predicted_test[j] - original_test[j], axis=1) + / norm(original_test[j], axis=1) + ) return error @@ -239,10 +249,12 @@ def optimal_mu(self, error=None, k=1): mu = self.database.parameters tria = Delaunay(mu) - error_on_simplex = np.array([ - np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) - for smpx in tria.simplices - ]) + error_on_simplex = np.array( + [ + np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) + for smpx in tria.simplices + ] + ) barycentric_point = [] for index in np.argpartition(error_on_simplex, -k)[-k:]: @@ -250,7 +262,8 @@ def optimal_mu(self, error=None, k=1): worst_tria_err = error[tria.simplices[index]] barycentric_point.append( - np.average(worst_tria_pts, axis=0, weights=worst_tria_err)) + np.average(worst_tria_pts, axis=0, weights=worst_tria_err) + ) return np.asarray(barycentric_point) @@ -269,4 +282,5 @@ def _simplex_volume(self, vertices): """ distance = np.transpose([vertices[0] - vi for vi in vertices[1:]]) return np.abs( - np.linalg.det(distance) / math.factorial(vertices.shape[1])) + np.linalg.det(distance) / math.factorial(vertices.shape[1]) + ) diff --git a/ezyrb/parallel/reduction.py b/ezyrb/parallel/reduction.py index 200c88d2..2c117f26 100644 --- a/ezyrb/parallel/reduction.py +++ b/ezyrb/parallel/reduction.py @@ -10,6 +10,7 @@ class Reduction(ABC): All the classes that implement the input-output mapping should be inherited from this class. """ + @abstractmethod def fit(self): """Abstract `fit`""" diff --git a/ezyrb/parameter.py b/ezyrb/parameter.py index 6f2ba619..3205a202 100644 --- a/ezyrb/parameter.py +++ b/ezyrb/parameter.py @@ -1,4 +1,5 @@ -""" Module for parameter object """ +"""Module for parameter object""" + import logging import numpy as np @@ -8,14 +9,14 @@ class Parameter: """ Class for representing a parameter in the reduced order model. - + This class encapsulates parameter values and provides validation to ensure parameters are 1-dimensional arrays. - + :param array_like values: The parameter values as a 1D array. - + :Example: - + >>> import numpy as np >>> from ezyrb import Parameter >>> param = Parameter([1.0, 2.5, 3.0]) @@ -29,38 +30,41 @@ class Parameter: def __init__(self, values): """ Initialize a Parameter object. - + :param array_like values: The parameter values. Can be a Parameter instance or an array-like object that can be converted to a 1D numpy array. """ - logger.debug("Initializing Parameter with values type: %s", - type(values)) + logger.debug( + "Initializing Parameter with values type: %s", type(values) + ) if isinstance(values, Parameter): self.values = values.values logger.debug("Copied values from existing Parameter") else: self.values = values - logger.debug("Created Parameter with shape: %s", - np.asarray(values).shape) + logger.debug( + "Created Parameter with shape: %s", np.asarray(values).shape + ) @property def values(self): - """ Get the snapshot values. """ + """Get the snapshot values.""" return self._values @values.setter def values(self, new_values): """ Set the parameter values with validation. - + :param array_like new_values: The new parameter values. :raises ValueError: If the new values are not a 1D array. """ if np.asarray(new_values).ndim != 1: - logger.error("Invalid parameter dimension: %d (expected 1D)", - np.asarray(new_values).ndim) - raise ValueError('only 1D array are usable as parameter.') + logger.error( + "Invalid parameter dimension: %d (expected 1D)", + np.asarray(new_values).ndim, + ) + raise ValueError("only 1D array are usable as parameter.") self._values = np.asarray(new_values) - logger.debug("Parameter values set with shape: %s", - self._values.shape) + logger.debug("Parameter values set with shape: %s", self._values.shape) diff --git a/ezyrb/plugin/__init__.py b/ezyrb/plugin/__init__.py index acfc7068..83279ee3 100644 --- a/ezyrb/plugin/__init__.py +++ b/ezyrb/plugin/__init__.py @@ -1,13 +1,13 @@ -""" Plugins submodule """ +"""Plugins submodule""" __all__ = [ - 'Plugin', - 'DatabaseScaler', - 'ShiftSnapshots', - 'AutomaticShiftSnapshots', - 'Aggregation', - 'DatabaseSplitter', - 'DatabaseDictionarySplitter' + "Plugin", + "DatabaseScaler", + "ShiftSnapshots", + "AutomaticShiftSnapshots", + "Aggregation", + "DatabaseSplitter", + "DatabaseDictionarySplitter", ] from .scaler import DatabaseScaler diff --git a/ezyrb/plugin/aggregation.py b/ezyrb/plugin/aggregation.py index ef9507cf..1d5f0c40 100644 --- a/ezyrb/plugin/aggregation.py +++ b/ezyrb/plugin/aggregation.py @@ -47,17 +47,19 @@ class Aggregation(Plugin): def __init__(self, fit_function=None, predict_function=Linear()): """ Initialize the Aggregation plugin. - + :param fit_function: Regression model to fit weights in validation set. If None, uses standard space-dependent methods. Default is None. :param predict_function: Regression model to predict weights in test set. Default is Linear(). """ super().__init__() - logger.debug("Initializing Aggregation with fit_function=%s, " - "predict_function=%s", - type(fit_function).__name__ if fit_function else None, - type(predict_function).__name__) + logger.debug( + "Initializing Aggregation with fit_function=%s, " + "predict_function=%s", + type(fit_function).__name__ if fit_function else None, + type(predict_function).__name__, + ) self.fit_function = fit_function self.predict_function = predict_function @@ -74,10 +76,12 @@ def _check_sum_gaussians(self, mrom, sum_gaussians, gaussians): :return: the corrected Gaussian functions. """ - zero_values = sum_gaussians < 1e-5 # tolerance to avoid numerical issues - + zero_values = ( + sum_gaussians < 1e-5 + ) # tolerance to avoid numerical issues + if zero_values.any(): - gaussians[:, zero_values] = 1/len(mrom.roms) # equal weights + gaussians[:, zero_values] = 1 / len(mrom.roms) # equal weights return gaussians def _optimize_sigma(self, mrom, sigma_range=[1e-5, 1e-2]): @@ -90,20 +94,39 @@ def _optimize_sigma(self, mrom, sigma_range=[1e-5, 1e-2]): :return: the optimal sigma value. """ + def obj_func(sigma): # compute test error of multiROM from current sigma - g_sigma = self._compute_validation_weights(mrom, sigma, normalized=True) - mrom_prediction = np.zeros_like(mrom.validation_full_database.snapshots_matrix) + g_sigma = self._compute_validation_weights( + mrom, sigma, normalized=True + ) + mrom_prediction = np.zeros_like( + mrom.validation_full_database.snapshots_matrix + ) for i, rom in enumerate(mrom.roms): - mrom_prediction += g_sigma[i, ...] * mrom.roms[rom].predict(mrom.validation_full_database).snapshots_matrix + mrom_prediction += ( + g_sigma[i, ...] + * mrom.roms[rom] + .predict(mrom.validation_full_database) + .snapshots_matrix + ) test_error = np.mean( - np.linalg.norm(mrom_prediction - mrom.validation_full_database.snapshots_matrix, - axis=1) / np.linalg.norm(mrom.validation_full_database.snapshots_matrix, axis=1)) + np.linalg.norm( + mrom_prediction + - mrom.validation_full_database.snapshots_matrix, + axis=1, + ) + / np.linalg.norm( + mrom.validation_full_database.snapshots_matrix, axis=1 + ) + ) return test_error + # minimization procedure - res = minimize(obj_func, x0=sigma_range[0], - method="L-BFGS-B", bounds=[sigma_range]) - print('Optimal sigma value in weights: ', res.x) + res = minimize( + obj_func, x0=sigma_range[0], method="L-BFGS-B", bounds=[sigma_range] + ) + print("Optimal sigma value in weights: ", res.x) return res.x def _compute_validation_weights(self, mrom, sigma, normalized=False): @@ -120,15 +143,19 @@ def _compute_validation_weights(self, mrom, sigma, normalized=False): """ validation_predicted = dict() for name, rom in mrom.roms.items(): - validation_predicted[name] = rom.predict(rom.validation_full_database).snapshots_matrix + validation_predicted[name] = rom.predict( + rom.validation_full_database + ).snapshots_matrix g = {} for k, v in validation_predicted.items(): snap = rom.validation_full_database.snapshots_matrix - g[k] = np.exp(- ((v - snap)**2)/(2 * (sigma**2))) + g[k] = np.exp(-((v - snap) ** 2) / (2 * (sigma**2))) g_tensor = np.array([g[k] for k in g.keys()]) - g_tensor = self._check_sum_gaussians(mrom, g_tensor.sum(axis=0), g_tensor) + g_tensor = self._check_sum_gaussians( + mrom, g_tensor.sum(axis=0), g_tensor + ) if normalized: g_tensor /= np.sum(g_tensor, axis=0) @@ -148,7 +175,7 @@ def fit_postprocessing(self, mrom): :return: None """ rom = list(mrom.roms.values())[0] - + # concatenate params and space params = rom.validation_full_database.parameters_matrix input_list = [] @@ -159,10 +186,9 @@ def fit_postprocessing(self, mrom): # Get the parameters for the i-th snapshot param = rom.validation_full_database.parameters_matrix[i, :] # Create the input array for the i-th snapshot - snapshot_input = np.hstack([ - space, - np.tile(param, (space.shape[0], 1)) - ]) + snapshot_input = np.hstack( + [space, np.tile(param, (space.shape[0], 1))] + ) # Append the input array to the list input_list.append(snapshot_input) # Concatenate the input arrays for all snapshots @@ -171,28 +197,32 @@ def fit_postprocessing(self, mrom): # Fit the regression/interpolation that will be used to predict the # weights in the test database if self.fit_function is None: - + optimal_sigma = self._optimize_sigma(mrom) - g_tensor = self._compute_validation_weights(mrom, optimal_sigma, normalized=False) - + g_tensor = self._compute_validation_weights( + mrom, optimal_sigma, normalized=False + ) + self.predict_functions = [] for i, rom in enumerate(mrom.roms.values()): - g_ = g_tensor[i, ...].reshape(space.shape[0] * params.shape[0], -1) + g_ = g_tensor[i, ...].reshape( + space.shape[0] * params.shape[0], -1 + ) # do a copy of the function used to predict the weights, # otherwise we use the same for all the ROMs rom_func = copy.deepcopy(self.predict_function) # replace NaN with 0 # TODO: this is a temporary fix, we should handle NaNs in a better way - g_[np.isnan(g_)] = 1/len(mrom.roms) + g_[np.isnan(g_)] = 1 / len(mrom.roms) rom_func.fit(input_, g_) self.predict_functions.append(rom_func) - - - # directly fit a regression to minimize the discrepancy between the aggregation and the snapshots + + # directly fit a regression to minimize the discrepancy between the aggregation and the snapshots elif self.fit_function is not None: snaps = rom.validation_full_database.snapshots_matrix - self.fit_function.fit(input_, snaps.reshape(space.shape[0] * params.shape[0], 1)) - + self.fit_function.fit( + input_, snaps.reshape(space.shape[0] * params.shape[0], 1) + ) def predict_postprocessing(self, mrom): """ @@ -207,61 +237,66 @@ def predict_postprocessing(self, mrom): # prepare the input for the prediction predict_weights = {} - + db = list(mrom.multi_predict_database.values())[0] - + input_list = [] - + # Loop over the number of snapshots in the prediction database - for i in range(mrom.predict_full_database.parameters_matrix.shape[0]): + for i in range(mrom.predict_full_database.parameters_matrix.shape[0]): # Get the space coordinates for the i-th snapshot space = mrom.train_full_database.get_snapshot_space(i) # Get the parameters for the i-th snapshot param = mrom.predict_full_database.parameters_matrix[i, :] # Create the input array for the i-th snapshot - snapshot_input = np.hstack([ - space, - np.tile(param, (space.shape[0], 1)) - ]) + snapshot_input = np.hstack( + [space, np.tile(param, (space.shape[0], 1))] + ) # Append the input array to the list input_list.append(snapshot_input) # Concatenate the input arrays for all snapshots input_ = np.vstack(input_list) - + # initialize weights mrom.weights_predict = {} # predict weights with regression and normalize them # case without fit_function (where we use the standard space-dependent aggregation methods) if self.fit_function is None and self.predict_function is not None: - + gaussians_test = [] for i, rom in enumerate(mrom.roms.values()): g_ = self.predict_functions[i].predict(input_) gaussians_test.append(g_) gaussians_test = np.array(gaussians_test) - gaussians_test = gaussians_test.reshape(len(mrom.roms), - db.parameters_matrix.shape[0], space.shape[0]) + gaussians_test = gaussians_test.reshape( + len(mrom.roms), db.parameters_matrix.shape[0], space.shape[0] + ) # normalize weights - predict_weights = gaussians_test/np.sum(gaussians_test, axis=0) + predict_weights = gaussians_test / np.sum(gaussians_test, axis=0) # compute predicted solution as convex combination of the reduced solutions for i, rom in enumerate(mrom.roms): mrom.weights_predict[rom] = predict_weights[i, ...] db = mrom.multi_predict_database[rom] - + # case with fit_function (ANN) elif self.fit_function is not None: weights = self.fit_function.predict(input_) # compute predicted solution as convex combination of the reduced solutions for i, rom in enumerate(mrom.roms): mrom.weights_predict[rom] = weights[..., i].reshape( - db.parameters_matrix.shape[0], -1) - + db.parameters_matrix.shape[0], -1 + ) + # compute the prediction - prediction = np.sum([mrom.weights_predict[k] * - mrom.multi_predict_database[k].snapshots_matrix for k in - list(mrom.roms.keys())], axis=0) + prediction = np.sum( + [ + mrom.weights_predict[k] + * mrom.multi_predict_database[k].snapshots_matrix + for k in list(mrom.roms.keys()) + ], + axis=0, + ) mrom.predict_full_database = Database( - mrom.predict_full_database.parameters_matrix, - prediction - ) \ No newline at end of file + mrom.predict_full_database.parameters_matrix, prediction + ) diff --git a/ezyrb/plugin/automatic_shift.py b/ezyrb/plugin/automatic_shift.py index 07cbfebb..9dede5b4 100644 --- a/ezyrb/plugin/automatic_shift.py +++ b/ezyrb/plugin/automatic_shift.py @@ -1,4 +1,4 @@ -""" Module for Scaler plugin """ +"""Module for Scaler plugin""" import numpy as np @@ -47,11 +47,19 @@ class AutomaticShiftSnapshots(Plugin): >>> rom = ROM(db, pod, rbf, plugins=[nnspod]) >>> rom.fit() """ - def __init__(self, shift_network, interp_network, interpolator, - parameter_index=0, reference_index=0, barycenter_loss=0): + + def __init__( + self, + shift_network, + interp_network, + interpolator, + parameter_index=0, + reference_index=0, + barycenter_loss=0, + ): """ Initialize the AutomaticShiftSnapshots plugin. - + :param shift_network: Neural network for learning the shift function. :param interp_network: Neural network for interpolation. :param Approximation interpolator: Interpolator for shifted snapshots evaluation. @@ -74,32 +82,46 @@ def _train_interp_network(self): """ self.interp_network.fit( self.reference_snapshot.space.reshape(-1, 1), - self.reference_snapshot.values.reshape(-1, 1) + self.reference_snapshot.values.reshape(-1, 1), ) def _train_shift_network(self, db): """ Train the shift network using the database snapshots. - + :param Database db: The database containing snapshots. """ - ref_center = torch.tensor(np.average( - self.reference_snapshot.space * self.reference_snapshot.values)) - - input_ = torch.from_numpy(np.vstack([ - np.vstack([s.space, np.ones(shape=(s.space.shape[0],))*p.values]).T - for p, s in db._pairs - ])).float() - output_ = torch.from_numpy(np.vstack([ - self.reference_snapshot.space.reshape(-1, 1) - for _ in db._pairs - ])) + ref_center = torch.tensor( + np.average( + self.reference_snapshot.space * self.reference_snapshot.values + ) + ) + + input_ = torch.from_numpy( + np.vstack( + [ + np.vstack( + [s.space, np.ones(shape=(s.space.shape[0],)) * p.values] + ).T + for p, s in db._pairs + ] + ) + ).float() + output_ = torch.from_numpy( + np.vstack( + [ + self.reference_snapshot.space.reshape(-1, 1) + for _ in db._pairs + ] + ) + ) self.shift_network._build_model(input_, output_) optimizer = self.shift_network.optimizer( self.shift_network.model.parameters(), lr=self.shift_network.lr, - weight_decay=self.shift_network.l2_regularization) + weight_decay=self.shift_network.l2_regularization, + ) n_epoch = 1 flag = True @@ -112,16 +134,25 @@ def _train_shift_network(self, db): tensor_space = torch.from_numpy(snap.space) tensor_values = torch.from_numpy(snap.values) - translated_space = tensor_space - shift.reshape(snap.space.shape) + translated_space = tensor_space - shift.reshape( + snap.space.shape + ) translated_space = translated_space.float() - interpolated_reference_values = self.interp_network.model(translated_space.reshape(-1, 1)).float().flatten() + interpolated_reference_values = ( + self.interp_network.model(translated_space.reshape(-1, 1)) + .float() + .flatten() + ) diff = torch.mean( - (tensor_values - interpolated_reference_values)**2) + (tensor_values - interpolated_reference_values) ** 2 + ) if self.barycenter_loss: snap_center = torch.mean(translated_space * tensor_values) - diff += self.barycenter_loss*(ref_center - snap_center)**2 + diff += ( + self.barycenter_loss * (ref_center - snap_center) ** 2 + ) loss += diff @@ -140,9 +171,12 @@ def _train_shift_network(self, db): if scalar_loss < criteria: flag = False - if (flag is False or - n_epoch == 1 or n_epoch % self.shift_network.frequency_print == 0): - print(f'[epoch {n_epoch:6d}]\t{scalar_loss:e}') + if ( + flag is False + or n_epoch == 1 + or n_epoch % self.shift_network.frequency_print == 0 + ): + print(f"[epoch {n_epoch:6d}]\t{scalar_loss:e}") n_epoch += 1 @@ -156,35 +190,46 @@ def fit_preprocessing(self, rom): self._train_shift_network(db) for param, snap in db._pairs: - input_shift = np.hstack([ - snap.space.reshape(-1, 1), - np.ones(shape=(snap.space.shape[0], 1))*param.values]) + input_shift = np.hstack( + [ + snap.space.reshape(-1, 1), + np.ones(shape=(snap.space.shape[0], 1)) * param.values, + ] + ) shift = self.shift_network.predict(input_shift) - + self.interpolator.fit( - snap.space.reshape(-1, 1) - shift, - snap.values.reshape(-1, 1)) + snap.space.reshape(-1, 1) - shift, snap.values.reshape(-1, 1) + ) snap.values = self.interpolator.predict( - reference_snapshot.space.reshape(-1, 1)).flatten() # reconstructing shifted snapshots in physical space + reference_snapshot.space.reshape(-1, 1) + ).flatten() # reconstructing shifted snapshots in physical space def predict_postprocessing(self, rom): - + ref_space = self.reference_snapshot.space db = Database() for param, snap in rom.predict_full_database._pairs: - input_shift = np.hstack([ - ref_space.reshape(-1, 1), - np.ones(shape=(ref_space.shape[0], 1))*param.values]) + input_shift = np.hstack( + [ + ref_space.reshape(-1, 1), + np.ones(shape=(ref_space.shape[0], 1)) * param.values, + ] + ) shift = self.shift_network.predict(input_shift) - snap.space = ref_space + shift.flatten() # shifted space transports to correct physical frame + snap.space = ( + ref_space + shift.flatten() + ) # shifted space transports to correct physical frame snap.space = snap.space.flatten() self.interpolator.fit(snap.space, snap.values.reshape(-1, 1)) - snap.values = self.interpolator.predict(ref_space) # reconstruct snapshot in physical space + snap.values = self.interpolator.predict( + ref_space + ) # reconstruct snapshot in physical space - snaps = Snapshot(values = snap.values, space = ref_space) + snaps = Snapshot(values=snap.values, space=ref_space) db.add(Parameter(param.values), snaps) - rom._full_database = db \ No newline at end of file + rom._full_database = db diff --git a/ezyrb/plugin/database_splitter.py b/ezyrb/plugin/database_splitter.py index 0096c5b8..4f013f98 100644 --- a/ezyrb/plugin/database_splitter.py +++ b/ezyrb/plugin/database_splitter.py @@ -1,4 +1,3 @@ - import logging from .plugin import Plugin from ..database import Database @@ -9,18 +8,18 @@ class DatabaseSplitter(Plugin): """ Plugin for splitting the database into training, test, validation, and prediction sets. - + This plugin automatically splits the database according to specified ratios before the fitting process begins. - + :param float train: Ratio or number of samples for training set. Default is 0.9. :param float test: Ratio or number of samples for test set. Default is 0.1. :param float validation: Ratio or number of samples for validation set. Default is 0.0. :param float predict: Ratio or number of samples for prediction set. Default is 0.0. :param int seed: Random seed for reproducibility. Default is None. - + :Example: - + >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database >>> from ezyrb.plugin import DatabaseSplitter @@ -35,11 +34,12 @@ class DatabaseSplitter(Plugin): >>> rom.fit() """ - def __init__(self, train=0.9, test=0.1, validation=0.0, predict=0.0, - seed=None): + def __init__( + self, train=0.9, test=0.1, validation=0.0, predict=0.0, seed=None + ): """ Initialize the DatabaseSplitter plugin. - + :param float train: Ratio for training set. Default is 0.9. :param float test: Ratio for test set. Default is 0.1. :param float validation: Ratio for validation set. Default is 0.0. @@ -47,9 +47,15 @@ def __init__(self, train=0.9, test=0.1, validation=0.0, predict=0.0, :param int seed: Random seed. Default is None. """ super().__init__() - logger.debug("Initializing DatabaseSplitter with train=%f, test=%f, " - "validation=%f, predict=%f, seed=%s", - train, test, validation, predict, seed) + logger.debug( + "Initializing DatabaseSplitter with train=%f, test=%f, " + "validation=%f, predict=%f, seed=%s", + train, + test, + validation, + predict, + seed, + ) self.train = train self.test = test @@ -60,7 +66,7 @@ def __init__(self, train=0.9, test=0.1, validation=0.0, predict=0.0, def fit_preprocessing(self, rom): """ Split the database before fitting begins. - + :param ReducedOrderModel rom: The ROM instance. """ logger.debug("Splitting database for ROM") @@ -69,14 +75,14 @@ def fit_preprocessing(self, rom): logger.debug("Splitting single Database") train, test, validation, predict = db.split( [self.train, self.test, self.validation, self.predict], - seed=self.seed + seed=self.seed, ) elif isinstance(db, dict): logger.debug("Splitting Database dictionary") train, test, validation, predict = list(db.values())[0].split( [self.train, self.test, self.validation, self.predict], - seed=self.seed + seed=self.seed, ) # TODO improve this splitting if needed (now only reading the # database of the first ROM) @@ -85,9 +91,13 @@ def fit_preprocessing(self, rom): rom.test_full_database = test rom.validation_full_database = validation rom.predict_full_database = predict - logger.info("Database split: train=%d, test=%d, validation=%d, " - "predict=%d", - len(train), len(test), len(validation), len(predict)) + logger.info( + "Database split: train=%d, test=%d, validation=%d, " "predict=%d", + len(train), + len(test), + len(validation), + len(predict), + ) # print('train', train.snapshots_matrix.shape) # print('test', test.snapshots_matrix.shape) # print('validation', validation.snapshots_matrix.shape) @@ -101,9 +111,9 @@ class DatabaseDictionarySplitter(Plugin): predict are already database objects stored in a dictionary. Given the desired keys of the dictionary as input, the plugin will assign the corresponding database objects to the train, test, validation and predict attributes of the ROM. - + :Example: - + >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database >>> from ezyrb.plugin import DatabaseDictionarySplitter @@ -117,12 +127,17 @@ class DatabaseDictionarySplitter(Plugin): >>> rom = ROM(db_dict['train'], pod, rbf, plugins=[splitter]) >>> rom.fit() """ - - def __init__(self, train_key=None, test_key=None, validation_key=None, - predict_key=None): + + def __init__( + self, + train_key=None, + test_key=None, + validation_key=None, + predict_key=None, + ): """ Initialize the DatabaseDictionarySplitter plugin. - + :param str train_key: Dictionary key for training database. Default is None. :param str test_key: Dictionary key for test database. Default is None. :param str validation_key: Dictionary key for validation database. Default is None. @@ -137,7 +152,7 @@ def __init__(self, train_key=None, test_key=None, validation_key=None, def fit_preprocessing(self, rom): """ Assign the database splits from the dictionary before fitting. - + :param ReducedOrderModel rom: The ROM instance. :raises ValueError: If the database is not a dictionary. """ @@ -153,4 +168,3 @@ def fit_preprocessing(self, rom): rom.predict_full_database = db[self.predict_key] else: raise ValueError("The database must be a dictionary of databases.") - \ No newline at end of file diff --git a/ezyrb/plugin/plugin.py b/ezyrb/plugin/plugin.py index 1eb739cb..34e6c9b5 100644 --- a/ezyrb/plugin/plugin.py +++ b/ezyrb/plugin/plugin.py @@ -6,14 +6,15 @@ class Plugin(ABC): """ The abstract Plugin class for ROM preprocessing and postprocessing. - + All plugin classes should inherit from this class and override the methods corresponding to the stages where they need to intervene. """ + def fit_preprocessing(self, rom): """ Execute before the fit process begins. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -21,7 +22,7 @@ def fit_preprocessing(self, rom): def fit_before_reduction(self, rom): """ Execute before the reduction step during fit. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -29,15 +30,15 @@ def fit_before_reduction(self, rom): def fit_after_reduction(self, rom): """ Execute after the reduction step during fit. - + :param ReducedOrderModel rom: The ROM instance. """ pass - + def fit_before_approximation(self, rom): """ Execute before the approximation step during fit. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -45,7 +46,7 @@ def fit_before_approximation(self, rom): def fit_after_approximation(self, rom): """ Execute after the approximation step during fit. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -53,7 +54,7 @@ def fit_after_approximation(self, rom): def fit_postprocessing(self, rom): """ Execute after the fit process completes. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -61,7 +62,7 @@ def fit_postprocessing(self, rom): def predict_preprocessing(self, rom): """ Execute before the prediction process begins. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -69,7 +70,7 @@ def predict_preprocessing(self, rom): def predict_before_approximation(self, rom): """ Execute before the approximation step during prediction. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -77,7 +78,7 @@ def predict_before_approximation(self, rom): def predict_after_approximation(self, rom): """ Execute after the approximation step during prediction. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -85,15 +86,15 @@ def predict_after_approximation(self, rom): def predict_before_expansion(self, rom): """ Execute before the expansion step during prediction. - + :param ReducedOrderModel rom: The ROM instance. """ pass - + def predict_after_expansion(self, rom): """ Execute after the expansion step during prediction. - + :param ReducedOrderModel rom: The ROM instance. """ pass @@ -101,10 +102,7 @@ def predict_after_expansion(self, rom): def predict_postprocessing(self, rom): """ Execute after the prediction process completes. - + :param ReducedOrderModel rom: The ROM instance. """ pass - - - diff --git a/ezyrb/plugin/scaler.py b/ezyrb/plugin/scaler.py index 5892fce0..588f54e7 100644 --- a/ezyrb/plugin/scaler.py +++ b/ezyrb/plugin/scaler.py @@ -46,8 +46,9 @@ def __init__(self, scaler, mode, target) -> None: :param str target: 'parameters' or 'snapshots' - what to scale. """ super().__init__() - logger.debug("Initializing DatabaseScaler with mode=%s, target=%s", - mode, target) + logger.debug( + "Initializing DatabaseScaler with mode=%s, target=%s", mode, target + ) self.scaler = scaler self.mode = mode diff --git a/ezyrb/plugin/shift.py b/ezyrb/plugin/shift.py index 836445f7..47877a52 100644 --- a/ezyrb/plugin/shift.py +++ b/ezyrb/plugin/shift.py @@ -1,4 +1,4 @@ -""" Module for Scaler plugin """ +"""Module for Scaler plugin""" import logging from .plugin import Plugin @@ -47,19 +47,25 @@ class ShiftSnapshots(Plugin): >>> rom.fit() """ - def __init__(self, shift_function, interpolator, parameter_index=0, - reference_index=0): + + def __init__( + self, shift_function, interpolator, parameter_index=0, reference_index=0 + ): """ Initialize the ShiftSnapshots plugin. - + :param callable shift_function: Function that returns the shift for a parameter. :param Approximation interpolator: Interpolator for evaluating shifted snapshots. :param int parameter_index: Index of parameter component for shift. Default is 0. :param int reference_index: Index of reference snapshot. Default is 0. """ super().__init__() - logger.debug("Initializing ShiftSnapshots with parameter_index=%d, " - "reference_index=%d", parameter_index, reference_index) + logger.debug( + "Initializing ShiftSnapshots with parameter_index=%d, " + "reference_index=%d", + parameter_index, + reference_index, + ) self.__shift_function = shift_function self.interpolator = interpolator @@ -69,27 +75,30 @@ def __init__(self, shift_function, interpolator, parameter_index=0, def fit_preprocessing(self, rom): """ Shift snapshots to a reference space during fit preprocessing. - + :param ReducedOrderModel rom: The ROM instance. """ logger.debug("Applying shift preprocessing") db = rom.database reference_snapshot = db._pairs[self.reference_index][1] - logger.debug("Using reference snapshot at index %d", - self.reference_index) + logger.debug( + "Using reference snapshot at index %d", self.reference_index + ) for param, snap in db._pairs: snap.space = reference_snapshot.space shift = self.__shift_function(param.values[self.parameter_index]) - logger.debug("Computed shift: %f for param: %s", - shift, param.values) + logger.debug( + "Computed shift: %f for param: %s", shift, param.values + ) self.interpolator.fit( - snap.space.reshape(-1, 1) - shift, - snap.values.reshape(-1, 1)) + snap.space.reshape(-1, 1) - shift, snap.values.reshape(-1, 1) + ) snap.values = self.interpolator.predict( - reference_snapshot.space.reshape(-1, 1)).flatten() + reference_snapshot.space.reshape(-1, 1) + ).flatten() rom.database = db logger.info("Shift preprocessing completed") @@ -97,11 +106,10 @@ def fit_preprocessing(self, rom): def predict_postprocessing(self, rom): """ Shift predicted snapshots back to their original space. - + :param ReducedOrderModel rom: The ROM instance. """ for param, snap in rom.predicted_full_database._pairs: - snap.space = ( - rom.database._pairs[self.reference_index][1].space + - self.__shift_function(param.values) - ) + snap.space = rom.database._pairs[self.reference_index][ + 1 + ].space + self.__shift_function(param.values) diff --git a/ezyrb/reducedordermodel.py b/ezyrb/reducedordermodel.py index af6efc19..03f31f04 100644 --- a/ezyrb/reducedordermodel.py +++ b/ezyrb/reducedordermodel.py @@ -20,12 +20,12 @@ class ReducedOrderModelInterface(ABC): """ Abstract interface for Reduced Order Model classes. - + This class defines the common interface and plugin execution mechanism for all ROM implementations. - + :Example: - + >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database >>> import numpy as np @@ -58,8 +58,9 @@ def _execute_plugins(self, when): logger.debug("Executing plugins at stage: %s", when) for plugin in self.plugins: if hasattr(plugin, when): - logger.debug("Running plugin %s.%s", - plugin.__class__.__name__, when) + logger.debug( + "Running plugin %s.%s", plugin.__class__.__name__, when + ) getattr(plugin, when)(self) @@ -98,11 +99,11 @@ class ReducedOrderModel(ReducedOrderModelInterface): >>> rom.predict(new_param) """ - def __init__(self, database, reduction, approximation, - plugins=None): + + def __init__(self, database, reduction, approximation, plugins=None): """ Initialize a Reduced Order Model. - + :param Database database: The database for training. :param Reduction reduction: The reduction method. :param Approximation approximation: The approximation method. @@ -131,7 +132,7 @@ def __init__(self, database, reduction, approximation, def clean(self): """ Clean all internal databases used during training and prediction. - + This method resets all training, prediction, test, and validation databases to None. """ @@ -153,7 +154,8 @@ def database(self, value): if not isinstance(value, (Database, dict)): raise TypeError( - "The database has to be an instance of the Database class, or a dictionary of Database.") + "The database has to be an instance of the Database class, or a dictionary of Database." + ) self._database = value @@ -169,7 +171,8 @@ def reduction(self): def reduction(self, value): if not isinstance(value, Reduction): raise TypeError( - "The reduction has to be an instance of the Reduction class, or a dict of Reduction.") + "The reduction has to be an instance of the Reduction class, or a dict of Reduction." + ) self._reduction = value @@ -185,7 +188,8 @@ def approximation(self): def approximation(self, value): if not isinstance(value, Approximation): raise TypeError( - "The approximation has to be an instance of the Approximation class, or a dict of Approximation.") + "The approximation has to be an instance of the Approximation class, or a dict of Approximation." + ) self._approximation = value @@ -211,16 +215,16 @@ def n_approximation(self): def fit_reduction(self): """ Fit the reduction method on the training database. - + This method applies the reduction technique to the snapshots matrix of the training database. - + :raises RuntimeError: If the training database has not been set. """ # for k, rom_ in self.roms.items(): # rom_['reduction'].fit(rom_['database'].snapshots_matrix.T) - if not hasattr(self, 'train_full_database'): + if not hasattr(self, "train_full_database"): raise RuntimeError self.reduction.fit(self.train_full_database.snapshots_matrix.T) @@ -228,31 +232,33 @@ def fit_reduction(self): def _reduce_database(self, db): """ Reduce a database using the fitted reduction method. - + :param Database db: The database to reduce. :return: A new database with reduced snapshots. :rtype: Database """ return Database( db.parameters_matrix, - self.reduction.transform(db.snapshots_matrix.T).T + self.reduction.transform(db.snapshots_matrix.T).T, ) def fit_approximation(self): """ Fit the approximation method on the reduced training database. - + This method trains the approximation technique on the reduced space representation of the snapshots. - + :raises RuntimeError: If the reduced training database has not been created. """ - if not hasattr(self, 'train_reduced_database'): + if not hasattr(self, "train_reduced_database"): raise RuntimeError - self.approximation.fit(self.train_reduced_database.parameters_matrix, - self.train_reduced_database.snapshots_matrix) + self.approximation.fit( + self.train_reduced_database.parameters_matrix, + self.train_reduced_database.snapshots_matrix, + ) # if self.n_database == 1 and self.n_reduction == 1: # self.train_full_database = self.database @@ -284,21 +290,24 @@ def fit(self): """ logger.info("Starting ROM fit process") - self._execute_plugins('fit_preprocessing') + self._execute_plugins("fit_preprocessing") if self.train_full_database is None: self.train_full_database = copy.deepcopy(self.database) logger.debug("Copied database to train_full_database") - self._execute_plugins('fit_before_reduction') + self._execute_plugins("fit_before_reduction") logger.debug("Fitting reduction method") self.fit_reduction() self.train_reduced_database = self._reduce_database( - self.train_full_database) - logger.info("Reduction completed. Reduced dimension: %s", - self.train_reduced_database.snapshots_matrix.shape) + self.train_full_database + ) + logger.info( + "Reduction completed. Reduced dimension: %s", + self.train_reduced_database.snapshots_matrix.shape, + ) - self._execute_plugins('fit_after_reduction') + self._execute_plugins("fit_after_reduction") # FULL-ORDER PREPROCESSING here # for plugin in self.plugins: # plugin.fom_preprocessing(self) @@ -309,14 +318,14 @@ def fit(self): # reduced_snapshots = self.reduction.transform( # self._full_database.snapshots_matrix.T).T - self._execute_plugins('fit_before_approximation') + self._execute_plugins("fit_before_approximation") logger.debug("Fitting approximation method") self.fit_approximation() - self._execute_plugins('fit_after_approximation') + self._execute_plugins("fit_after_approximation") - self._execute_plugins('fit_postprocessing') + self._execute_plugins("fit_postprocessing") logger.info("ROM fit process completed successfully") return self @@ -334,58 +343,56 @@ def predict(self, parameters): """ logger.info("Starting ROM predict process") logger.debug("Parameters type: %s", type(parameters)) - self._execute_plugins('predict_preprocessing') + self._execute_plugins("predict_preprocessing") if isinstance(parameters, Database): self.predict_reduced_database = parameters - logger.debug("Using Database with %d parameters", - len(parameters)) + logger.debug("Using Database with %d parameters", len(parameters)) elif isinstance(parameters, (list, np.ndarray, tuple)): parameters = np.atleast_2d(parameters) logger.debug("Parameters shape: %s", parameters.shape) self.predict_reduced_database = Database( - parameters=parameters, - snapshots=[None] * len(parameters) + parameters=parameters, snapshots=[None] * len(parameters) ) else: logger.error("Invalid parameters type: %s", type(parameters)) raise TypeError - self._execute_plugins('predict_before_approximation') + self._execute_plugins("predict_before_approximation") logger.debug("Predicting with approximation method") self.predict_reduced_database = Database( self.predict_reduced_database.parameters_matrix, self.approximation.predict( - self.predict_reduced_database.parameters_matrix).reshape( - self.predict_reduced_database.parameters_matrix.shape[0], - -1 - ) + self.predict_reduced_database.parameters_matrix + ).reshape( + self.predict_reduced_database.parameters_matrix.shape[0], -1 + ), ) - self._execute_plugins('predict_after_approximation') + self._execute_plugins("predict_after_approximation") # mu = np.atleast_2d(mu) # for plugin in self.plugins: # if plugin.target == 'parameters': # mu = plugin.scaler.transform(mu) - # self._reduced_database = Database( # mu, np.atleast_2d(self.approximation.predict(mu))) - self._execute_plugins('predict_before_expansion') + self._execute_plugins("predict_before_expansion") logger.debug("Expanding to full space") self.predicted_full_database = Database( self.predict_reduced_database.parameters_matrix, self.reduction.inverse_transform( - self.predict_reduced_database.snapshots_matrix.T).T + self.predict_reduced_database.snapshots_matrix.T + ).T, ) - self._execute_plugins('predict_after_expansion') + self._execute_plugins("predict_after_expansion") - self._execute_plugins('predict_postprocessing') + self._execute_plugins("predict_postprocessing") logger.info("ROM predict process completed") if isinstance(parameters, Database): @@ -393,7 +400,6 @@ def predict(self, parameters): else: return self.predicted_full_database.snapshots_matrix - def save(self, fname, save_db=True, save_reduction=True, save_approx=True): """ Save the object to `fname` using the pickle module. @@ -422,7 +428,7 @@ def save(self, fname, save_db=True, save_reduction=True, save_approx=True): if not save_approx: del rom_to_store.approximation - with open(fname, 'wb') as output: + with open(fname, "wb") as output: pickle.dump(rom_to_store, output, pickle.HIGHEST_PROTOCOL) @staticmethod @@ -438,7 +444,7 @@ def load(fname): >>> rom = ROM.load('ezyrb.rom') >>> rom.predict(new_param) """ - with open(fname, 'rb') as output: + with open(fname, "rb") as output: rom = pickle.load(output) return rom @@ -452,8 +458,8 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): :param function norm: the function used to assign at the vector of errors a float number. It has to take as input a 'numpy.ndarray' and returns a float. Default value is the L2 norm. - :param relative: True if the error computed is relative. Default is - True. + :param relative: True if the error computed is relative. Default is + True. :return: the mean L2 norm of the relative errors of the estimated test snapshots. :rtype: numpy.ndarray @@ -461,15 +467,15 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): predicted_test = self.predict(test.parameters_matrix) if relative: return np.mean( - norm(predicted_test - test.snapshots_matrix, - axis=1) / norm(test.snapshots_matrix, axis=1)) + norm(predicted_test - test.snapshots_matrix, axis=1) + / norm(test.snapshots_matrix, axis=1) + ) else: - return np.mean( - norm(predicted_test - test.snapshots_matrix, - axis=1)) + return np.mean(norm(predicted_test - test.snapshots_matrix, axis=1)) - def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, - **kwargs): + def kfold_cv_error( + self, n_splits, *args, norm=np.linalg.norm, relative=True, **kwargs + ): r""" Split the database into k consecutive folds (no shuffling by default). Each fold is used once as a validation while the k - 1 remaining folds @@ -481,8 +487,8 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, :param function norm: function to apply to compute the relative error between the true snapshot and the predicted one. Default value is the L2 norm. - :param relative: True if the error computed is relative. Default is - True. + :param relative: True if the error computed is relative. Default is + True. :param \*args: additional parameters to pass to the `fit` method. :param \**kwargs: additional parameters to pass to the `fit` method. :return: the vector containing the errors corresponding to each fold. @@ -492,12 +498,16 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, kf = KFold(n_splits=n_splits) for train_index, test_index in kf.split(self.database): new_db = self.database[train_index] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation), - plugins=[copy.deepcopy(p) for p in self.plugins]).fit( - *args, **kwargs) - - error.append(rom.test_error(self.database[test_index], norm, relative)) + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + plugins=[copy.deepcopy(p) for p in self.plugins], + ).fit(*args, **kwargs) + + error.append( + rom.test_error(self.database[test_index], norm, relative) + ) return np.array(error) @@ -529,9 +539,12 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): new_db = self.database[indeces] test_db = self.database[~indeces] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation), - plugins=[copy.deepcopy(p) for p in self.plugins]).fit() + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + plugins=[copy.deepcopy(p) for p in self.plugins], + ).fit() error[j] = rom.test_error(test_db) @@ -558,10 +571,12 @@ def optimal_mu(self, error=None, k=1): mu = self.database.parameters_matrix tria = Delaunay(mu) - error_on_simplex = np.array([ - np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) - for smpx in tria.simplices - ]) + error_on_simplex = np.array( + [ + np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) + for smpx in tria.simplices + ] + ) barycentric_point = [] for index in np.argpartition(error_on_simplex, -k)[-k:]: @@ -569,7 +584,8 @@ def optimal_mu(self, error=None, k=1): worst_tria_err = error[tria.simplices[index]] barycentric_point.append( - np.average(worst_tria_pts, axis=0, weights=worst_tria_err)) + np.average(worst_tria_pts, axis=0, weights=worst_tria_err) + ) return np.asarray(barycentric_point) @@ -588,26 +604,27 @@ def _simplex_volume(self, vertices): """ distance = np.transpose([vertices[0] - vi for vi in vertices[1:]]) return np.abs( - np.linalg.det(distance) / math.factorial(vertices.shape[1])) + np.linalg.det(distance) / math.factorial(vertices.shape[1]) + ) def reduction_error(self, db=None, relative=True, eps=1e-12): """ Calculate the reconstruction error between the original snapshots and the ones reconstructed by the ROM. - + :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is + :param bool relative: True if the error computed is relative. Default is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. :return: the vector containing the reconstruction errors. - - Esempio: + + Esempio: >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database - >>> db = Database(param, snapshots) # param and snapshots are assumed + >>> db = Database(param, snapshots) # param and snapshots are assumed to be declared >>> db_train = db[:10] # training database >>> db_test = db[10:] # test database @@ -618,10 +635,10 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): >>> err_train_reduct = rom.reconstruction_error(relative=True) >>> err_test_reduct = rom.reconstruction_error(db_test, relative=True) """ - + errs = [] if db is None: - db = self.database + db = self.database snap = db.snapshots_matrix snap_red = self.reduction.transform(snap.T) snap_full = self.reduction.inverse_transform(snap_red).T @@ -631,33 +648,33 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): if relative: num = np.linalg.norm(E, axis=1) den = np.linalg.norm(snap, axis=1) + eps - - err = float(np.mean(num/den)) + + err = float(np.mean(num / den)) else: err = float(np.mean(np.linalg.norm(E, axis=1))) errs.append(err) - + return np.array(errs) def approximation_error(self, db=None, relative=True, eps=1e-12): """ Calculate the approximation error between the true modal coefficients and the approximated ones. - + :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is + :param bool relative: True if the error computed is relative. Default is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. - + :return: the vector containing the approximation errors. - + Esempio: >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database - >>> db = Database(param, snapshots) # param and snapshots are assumed + >>> db = Database(param, snapshots) # param and snapshots are assumed to be declared >>> db_train = db[:10] # training database >>> db_test = db[10:] # test database @@ -686,7 +703,7 @@ def approximation_error(self, db=None, relative=True, eps=1e-12): num = np.linalg.norm(E, axis=1) den = np.linalg.norm(params_true, axis=1) + eps - err = float(np.mean(num/den)) + err = float(np.mean(num / den)) else: err = float(np.mean(np.linalg.norm(E, axis=1))) errs.append(err) @@ -694,7 +711,6 @@ def approximation_error(self, db=None, relative=True, eps=1e-12): return np.array(errs) - class MultiReducedOrderModel(ReducedOrderModelInterface): """ Multiple Reduced Order Model class. @@ -730,15 +746,16 @@ class MultiReducedOrderModel(ReducedOrderModelInterface): >>> rom.predict(new_param) """ + def __init__(self, *args, plugins=None, rom_plugin=None): """ Initialize a Multi-ROM with multiple databases and methods. - + Supports multiple initialization signatures: - (database_dict, reduction_dict, approximation_dict) - (database, roms_dict) - (roms_dict,) - + :param args: Variable arguments for different initialization modes. :param list plugins: Global plugins for the Multi-ROM. Default is None. :param rom_plugin: Plugin to add to each individual ROM. Default is None. @@ -754,27 +771,31 @@ def __init__(self, *args, plugins=None, rom_plugin=None): element_keys = product( self.database.keys(), self.reduction.keys(), - self.approximation.keys() + self.approximation.keys(), ) self.roms = { - tuple(key): ReducedOrderModel( copy.deepcopy(self.database[key[0]]), copy.deepcopy(self.reduction[key[1]]), - copy.deepcopy(self.approximation[key[2]]) + copy.deepcopy(self.approximation[key[2]]), ) for key in element_keys } - elif (len(args) == 2 and isinstance(args[0], Database) - and isinstance(args[1], dict)): + elif ( + len(args) == 2 + and isinstance(args[0], Database) + and isinstance(args[1], dict) + ): self.database = args[0] self.roms = args[1] elif len(args) == 1 and isinstance(args[0], dict): self.roms = args[0] roms_keys = list(self.roms.keys()) - self.database = {roms_keys[i]: self.roms[roms_keys[i]].database - for i in range(len(self.roms))} + self.database = { + roms_keys[i]: self.roms[roms_keys[i]].database + for i in range(len(self.roms)) + } if plugins is None: plugins = [] @@ -784,7 +805,6 @@ def __init__(self, *args, plugins=None, rom_plugin=None): for rom_ in self.roms.values(): rom_.plugins.append(rom_plugin) - @property def database(self): return self._database @@ -794,7 +814,8 @@ def database(self, value): if not isinstance(value, (Database, dict)): raise TypeError( - "The database has to be an instance of the Database class, or a dictionary of Database.") + "The database has to be an instance of the Database class, or a dictionary of Database." + ) if isinstance(value, Database): self._database = {0: value} @@ -813,7 +834,8 @@ def reduction(self): def reduction(self, value): if not isinstance(value, (Reduction, dict)): raise TypeError( - "The reduction has to be an instance of the Reduction class, or a dict of Reduction.") + "The reduction has to be an instance of the Reduction class, or a dict of Reduction." + ) if isinstance(value, Reduction): self._reduction = {0: value} @@ -832,7 +854,8 @@ def approximation(self): def approximation(self, value): if not isinstance(value, (Approximation, dict)): raise TypeError( - "The approximation has to be an instance of the Approximation class, or a dict of Approximation.") + "The approximation has to be an instance of the Approximation class, or a dict of Approximation." + ) if isinstance(value, Approximation): self._approximation = {0: value} @@ -853,21 +876,19 @@ def n_reduction(self): def n_approximation(self): value_, class_ = self.approximation, Approximation return len(value_) if isinstance(value_, class_) else 1 - + def fit(self): r""" Calculate reduced space """ - self._execute_plugins('fit_preprocessing') + self._execute_plugins("fit_preprocessing") for rom_ in self.roms.values(): - if rom_.train_reduced_database==None: + if rom_.train_reduced_database == None: rom_.fit() - - self._execute_plugins('fit_postprocessing') - return self - + self._execute_plugins("fit_postprocessing") + return self def predict(self, parameters=None): """ @@ -880,7 +901,7 @@ def predict(self, parameters=None): :rtype: Database """ - self._execute_plugins('predict_preprocessing') + self._execute_plugins("predict_preprocessing") # convert parameters from Database to numpy array (if database) if isinstance(parameters, Database): @@ -889,7 +910,9 @@ def predict(self, parameters=None): elif isinstance(parameters, (list, np.ndarray, tuple)): print(parameters) parameters = np.atleast_2d(parameters) - self.predict_reduced_database = Database(parameters, [None]*len(parameters)) + self.predict_reduced_database = Database( + parameters, [None] * len(parameters) + ) elif parameters is None: if self.predict_full_database is None: raise RuntimeError @@ -898,18 +921,19 @@ def predict(self, parameters=None): self.multi_predict_database = {} for k, rom_ in self.roms.items(): - self.multi_predict_database[k] = rom_.predict(self.predict_reduced_database) + self.multi_predict_database[k] = rom_.predict( + self.predict_reduced_database + ) print(self.multi_predict_database) - self._execute_plugins('predict_postprocessing') + self._execute_plugins("predict_postprocessing") if isinstance(parameters, Database): return self.multi_predict_database else: return { - k:db.snapshots_matrix + k: db.snapshots_matrix for k, db in self.multi_predict_database.items() } - def save(self, fname, save_db=True, save_reduction=True, save_approx=True): """ @@ -939,7 +963,7 @@ def save(self, fname, save_db=True, save_reduction=True, save_approx=True): if not save_approx: del rom_to_store.approximation - with open(fname, 'wb') as output: + with open(fname, "wb") as output: pickle.dump(rom_to_store, output, pickle.HIGHEST_PROTOCOL) @staticmethod @@ -955,7 +979,7 @@ def load(fname): >>> rom = ROM.load('ezyrb.rom') >>> rom.predict(new_param) """ - with open(fname, 'rb') as output: + with open(fname, "rb") as output: rom = pickle.load(output) return rom @@ -969,8 +993,8 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): :param function norm: the function used to assign at the vector of errors a float number. It has to take as input a 'numpy.ndarray' and returns a float. Default value is the L2 norm. - :param relative: True if the error computed is relative. Default is - True. + :param relative: True if the error computed is relative. Default is + True. :return: the mean L2 norm of the relative errors of the estimated test snapshots. :rtype: numpy.ndarray @@ -978,15 +1002,23 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): predicted_test = self.predict(test.parameters_matrix) if relative: return np.mean( - norm(predicted_test.snapshots_matrix - test.snapshots_matrix, - axis=1) / norm(test.snapshots_matrix, axis=1)) + norm( + predicted_test.snapshots_matrix - test.snapshots_matrix, + axis=1, + ) + / norm(test.snapshots_matrix, axis=1) + ) else: return np.mean( - norm(predicted_test.snapshots_matrix - test.snapshots_matrix, - axis=1)) + norm( + predicted_test.snapshots_matrix - test.snapshots_matrix, + axis=1, + ) + ) - def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, - **kwargs): + def kfold_cv_error( + self, n_splits, *args, norm=np.linalg.norm, relative=True, **kwargs + ): r""" Split the database into k consecutive folds (no shuffling by default). Each fold is used once as a validation while the k - 1 remaining folds @@ -998,8 +1030,8 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, :param function norm: function to apply to compute the relative error between the true snapshot and the predicted one. Default value is the L2 norm. - :param relative: True if the error computed is relative. Default is - True. + :param relative: True if the error computed is relative. Default is + True. :param \*args: additional parameters to pass to the `fit` method. :param \**kwargs: additional parameters to pass to the `fit` method. :return: the vector containing the errors corresponding to each fold. @@ -1010,11 +1042,15 @@ def kfold_cv_error(self, n_splits, *args, norm=np.linalg.norm, relative=True, for train_index, test_index in kf.split(self.database): new_db = self.database[train_index] # TODO: Fix plugins handling - should pass: plugins=[copy.deepcopy(p) for p in self.plugins] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation)).fit( - *args, **kwargs) - - error.append(rom.test_error(self.database[test_index], norm, relative)) + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + ).fit(*args, **kwargs) + + error.append( + rom.test_error(self.database[test_index], norm, relative) + ) return np.array(error) @@ -1047,8 +1083,11 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): new_db = self.database[indeces] test_db = self.database[~indeces] # TODO: Fix plugins handling - should pass: plugins=[copy.deepcopy(p) for p in self.plugins] - rom = type(self)(new_db, copy.deepcopy(self.reduction), - copy.deepcopy(self.approximation)).fit() + rom = type(self)( + new_db, + copy.deepcopy(self.reduction), + copy.deepcopy(self.approximation), + ).fit() error[j] = rom.test_error(test_db) @@ -1075,10 +1114,12 @@ def optimal_mu(self, error=None, k=1): mu = self.database.parameters_matrix tria = Delaunay(mu) - error_on_simplex = np.array([ - np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) - for smpx in tria.simplices - ]) + error_on_simplex = np.array( + [ + np.sum(error[smpx]) * self._simplex_volume(mu[smpx]) + for smpx in tria.simplices + ] + ) barycentric_point = [] for index in np.argpartition(error_on_simplex, -k)[-k:]: @@ -1086,7 +1127,8 @@ def optimal_mu(self, error=None, k=1): worst_tria_err = error[tria.simplices[index]] barycentric_point.append( - np.average(worst_tria_pts, axis=0, weights=worst_tria_err)) + np.average(worst_tria_pts, axis=0, weights=worst_tria_err) + ) return np.asarray(barycentric_point) @@ -1105,26 +1147,27 @@ def _simplex_volume(self, vertices): """ distance = np.transpose([vertices[0] - vi for vi in vertices[1:]]) return np.abs( - np.linalg.det(distance) / math.factorial(vertices.shape[1])) + np.linalg.det(distance) / math.factorial(vertices.shape[1]) + ) def reduction_error(self, db=None, relative=True, eps=1e-12): """ Calculate the reconstruction error between the original snapshots and the ones reconstructed by the ROM. - + :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is + :param bool relative: True if the error computed is relative. Default is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. :return: the vector containing the reconstruction errors. - - Esempio: + + Esempio: >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database - >>> db = Database(param, snapshots) # param and snapshots are assumed + >>> db = Database(param, snapshots) # param and snapshots are assumed to be declared >>> db_train = db[:10] # training database >>> db_test = db[10:] # test database @@ -1135,10 +1178,10 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): >>> err_train_reduct = rom.reconstruction_error(relative=True) >>> err_test_reduct = rom.reconstruction_error(db_test, relative=True) """ - + errs = [] if db is None: - db = self.database + db = self.database snap = db.snapshots_matrix snap_red = self.reduction.transform(snap.T) snap_full = self.reduction.inverse_transform(snap_red).T @@ -1148,33 +1191,33 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): if relative: num = np.linalg.norm(E, axis=1) den = np.linalg.norm(snap, axis=1) + eps - - err = float(np.mean(num/den)) + + err = float(np.mean(num / den)) else: err = float(np.mean(np.linalg.norm(E, axis=1))) errs.append(err) - + return np.array(errs) def approximation_error(self, db=None, relative=True, eps=1e-12): """ Calculate the approximation error between the true modal coefficients and the approximated ones. - + :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is + :param bool relative: True if the error computed is relative. Default is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. - + :return: the vector containing the approximation errors. - + Esempio: >>> from ezyrb import ReducedOrderModel as ROM >>> from ezyrb import POD, RBF, Database - >>> db = Database(param, snapshots) # param and snapshots are assumed + >>> db = Database(param, snapshots) # param and snapshots are assumed to be declared >>> db_train = db[:10] # training database >>> db_test = db[10:] # test database @@ -1203,10 +1246,9 @@ def approximation_error(self, db=None, relative=True, eps=1e-12): num = np.linalg.norm(E, axis=1) den = np.linalg.norm(params_true, axis=1) + eps - err = float(np.mean(num/den)) + err = float(np.mean(num / den)) else: err = float(np.mean(np.linalg.norm(E, axis=1))) errs.append(err) return np.array(errs) - diff --git a/ezyrb/reduction/__init__.py b/ezyrb/reduction/__init__.py index 6db8d842..83b247a0 100644 --- a/ezyrb/reduction/__init__.py +++ b/ezyrb/reduction/__init__.py @@ -1,12 +1,6 @@ -""" Reduction submodule """ +"""Reduction submodule""" -__all__ = [ - 'Reduction', - 'POD', - 'AE', - 'PODAE', - 'SklearnReduction' -] +__all__ = ["Reduction", "POD", "AE", "PODAE", "SklearnReduction"] from .reduction import Reduction from .pod import POD diff --git a/ezyrb/reduction/ae.py b/ezyrb/reduction/ae.py index c47402fe..b3972686 100644 --- a/ezyrb/reduction/ae.py +++ b/ezyrb/reduction/ae.py @@ -44,9 +44,9 @@ def __init__( latent_dim, layers_encoder, layers_decoder, - activation='tanh', + activation="tanh", max_iter=200, - solver='adam', + solver="adam", learning_rate_init=0.001, alpha=0, frequency_print=10, @@ -68,7 +68,9 @@ def __init__( raise ValueError("Wrong dimension in encoder and decoder layers") if not isinstance(latent_dim, int): - logger.error("latent_dim should be an integer, got %s", type(latent_dim)) + logger.error( + "latent_dim should be an integer, got %s", type(latent_dim) + ) raise ValueError("latent_dim should be an integer") self.latent_dim = latent_dim @@ -82,7 +84,7 @@ def __init__( self.frequency_print = frequency_print self.loss_trend = [] self.extra_kwargs = kwargs - + # Models self._autoencoder = None # Full trained model self.encoder = None # For encoding @@ -96,14 +98,16 @@ def fit(self, values): """ logger.info("Starting AE training") values = values.T - + # Combine encoder and decoder layers - combined_layers = self.layers_encoder + [self.latent_dim] + self.layers_decoder - + combined_layers = ( + self.layers_encoder + [self.latent_dim] + self.layers_decoder + ) + logger.debug( "Training full autoencoder with layers: %s", combined_layers ) - + # Train full autoencoder: input -> latent -> reconstruction self._autoencoder = ANN( combined_layers, @@ -114,19 +118,19 @@ def fit(self, values): alpha=self.alpha, **self.extra_kwargs, ) - + # Train to reconstruct input self._autoencoder.fit(values, values) self.loss_trend = self._autoencoder.loss_trend - + # Now create encoder and decoder and copy weights logger.debug("Creating encoder and decoder from trained autoencoder") - + # Create encoder self.encoder = ANN( self.layers_encoder, activation=self.activation, - solver='adam', + solver="adam", max_iter=1, learning_rate_init=self.learning_rate_init, alpha=self.alpha, @@ -135,7 +139,7 @@ def fit(self, values): self.decoder = ANN( self.layers_decoder, activation=self.activation, - solver='adam', + solver="adam", max_iter=1, learning_rate_init=self.learning_rate_init, alpha=self.alpha, @@ -145,29 +149,29 @@ def fit(self, values): # Dummy fit to initialize structure self.decoder.fit(dummy_latent, values) self.encoder.fit(values, dummy_latent) - + # Copy encoder weights from autoencoder n_encoder_layers = len(self.encoder.model.coefs_) for i in range(n_encoder_layers): - self.encoder.model.coefs_[i] = ( - self._autoencoder.model.coefs_[i].copy() - ) + self.encoder.model.coefs_[i] = self._autoencoder.model.coefs_[ + i + ].copy() self.encoder.model.intercepts_[i] = ( self._autoencoder.model.intercepts_[i].copy() ) - + # Copy decoder weights from autoencoder n_decoder_layers = len(self.decoder.model.coefs_) for i in range(n_decoder_layers): src_idx = len(self.encoder.model.coefs_) + i - self.decoder.model.coefs_[i] = ( - self._autoencoder.model.coefs_[src_idx].copy() - ) + self.decoder.model.coefs_[i] = self._autoencoder.model.coefs_[ + src_idx + ].copy() self.decoder.model.intercepts_[i] = ( self._autoencoder.model.intercepts_[src_idx].copy() ) - - print('dentro fit') + + print("dentro fit") print(values.shape) print(self._autoencoder.predict(values)) red = self.encoder.predict(values) @@ -198,13 +202,13 @@ def inverse_transform(self, g): if self.decoder is None: raise RuntimeError("Autoencoder not fitted yet") - if self.activation == 'tanh': + if self.activation == "tanh": g = np.tanh(g) - elif self.activation == 'logistic': + elif self.activation == "logistic": g = 1 / (1 + np.exp(-g)) - elif self.activation == 'relu': + elif self.activation == "relu": g = np.maximum(0, g) - elif self.activation == 'identity': + elif self.activation == "identity": pass return self.decoder.predict(g.T).T diff --git a/ezyrb/reduction/pod.py b/ezyrb/reduction/pod.py index 157cf7c1..ef77f429 100755 --- a/ezyrb/reduction/pod.py +++ b/ezyrb/reduction/pod.py @@ -61,8 +61,7 @@ def __init__( >>> pod = POD('correlation_matrix', rank=10, save_memory=False) """ logger.debug("Initializing POD with method=%s, rank=%s", method, rank) - self.available_methods = ["svd", "randomized_svd", - "correlation_matrix"] + self.available_methods = ["svd", "randomized_svd", "correlation_matrix"] self.rank = rank if method == "svd": self._method = self._svd @@ -104,9 +103,9 @@ def singular_values(self): def fit(self, X): """ Create the reduced space for the given snapshots using POD. - + Computes the POD modes and singular values using the specified method. - + :param numpy.ndarray X: The input snapshots matrix (stored by column). :return: self """ @@ -120,9 +119,11 @@ def fit(self, X): ) self._modes, self._singular_values = self._method(X) logger.info("POD fitted: %d modes extracted", self._modes.shape[1]) - logger.debug("Singular values range: [%f, %f]", - self._singular_values.min(), - self._singular_values.max()) + logger.debug( + "Singular values range: [%f, %f]", + self._singular_values.min(), + self._singular_values.max(), + ) return self def transform(self, X): @@ -216,13 +217,13 @@ def _svd(self, X): def _rsvd(self, X): """ Truncated randomized Singular Value Decomposition. - + Computes an approximate SVD using randomized algorithms for efficiency. - + :param numpy.ndarray X: The matrix to decompose. :return: Tuple of (truncated left-singular vectors, truncated singular values). :rtype: tuple(numpy.ndarray, numpy.ndarray) - + References: Finding structure with randomness: probabilistic algorithms for constructing approximate matrix decompositions. N. Halko, P. G. diff --git a/ezyrb/reduction/pod_ae.py b/ezyrb/reduction/pod_ae.py index 00e4cca8..a01036bb 100644 --- a/ezyrb/reduction/pod_ae.py +++ b/ezyrb/reduction/pod_ae.py @@ -9,17 +9,18 @@ class PODAE(POD, AE): """ Combined POD and AutoEncoder reduction class. - + This class first applies POD to reduce the dimensionality, then uses an autoencoder for further reduction in the latent space. - + :param POD pod: The POD instance for initial reduction. :param AE ae: The AutoEncoder instance for latent space reduction. """ + def __init__(self, pod, ae): """ Initialize the PODAE reducer. - + :param POD pod: The POD instance. :param AE ae: The AutoEncoder instance. """ @@ -29,9 +30,9 @@ def __init__(self, pod, ae): def fit(self, X): """ Fit the PODAE on the snapshots. - + First applies POD, then trains the autoencoder on POD coefficients. - + :param numpy.ndarray X: The input snapshots matrix (stored by column). """ self.pod.fit(X) diff --git a/ezyrb/reduction/reduction.py b/ezyrb/reduction/reduction.py index a5bc15c6..e95095e8 100644 --- a/ezyrb/reduction/reduction.py +++ b/ezyrb/reduction/reduction.py @@ -10,6 +10,7 @@ class Reduction(ABC): All the classes that implement dimensionality reduction should be inherited from this class. """ + @abstractmethod def fit(self): """Abstract `fit`""" @@ -28,4 +29,4 @@ def expand(self): @abstractmethod def reduce(self): - """Abstract `reduce`""" \ No newline at end of file + """Abstract `reduce`""" diff --git a/ezyrb/reduction/sklearn_reduction.py b/ezyrb/reduction/sklearn_reduction.py index 842747d9..7e0c519c 100644 --- a/ezyrb/reduction/sklearn_reduction.py +++ b/ezyrb/reduction/sklearn_reduction.py @@ -57,28 +57,24 @@ def __init__(self, sklearn_model, fit_params=None): """ logger.debug( "Initializing SklearnReduction with model: %s", - type(sklearn_model).__name__ + type(sklearn_model).__name__, ) - if not hasattr(sklearn_model, 'fit'): - raise ValueError( - "sklearn_model must have a 'fit' method" - ) - if not hasattr(sklearn_model, 'transform'): - raise ValueError( - "sklearn_model must have a 'transform' method" - ) + if not hasattr(sklearn_model, "fit"): + raise ValueError("sklearn_model must have a 'fit' method") + if not hasattr(sklearn_model, "transform"): + raise ValueError("sklearn_model must have a 'transform' method") self.model = sklearn_model self.fit_params = fit_params if fit_params is not None else {} self._fitted = False - self._has_inverse = hasattr(sklearn_model, 'inverse_transform') + self._has_inverse = hasattr(sklearn_model, "inverse_transform") if not self._has_inverse: logger.warning( "%s does not have inverse_transform method. " "inverse_transform() will raise an error.", - type(sklearn_model).__name__ + type(sklearn_model).__name__, ) def fit(self, values): @@ -90,7 +86,7 @@ def fit(self, values): logger.info( "Fitting %s with %d snapshots", type(self.model).__name__, - values.shape[0] + values.shape[0], ) logger.debug("Input shape: %s", values.shape) @@ -104,12 +100,9 @@ def fit(self, values): logger.debug("Model fitting completed") # Log explained variance if available (e.g., PCA) - if hasattr(self.model, 'explained_variance_ratio_'): + if hasattr(self.model, "explained_variance_ratio_"): total_var = self.model.explained_variance_ratio_.sum() - logger.info( - "Explained variance ratio: %.4f", - total_var - ) + logger.info("Explained variance ratio: %.4f", total_var) def transform(self, values): """ @@ -124,10 +117,7 @@ def transform(self, values): "Model must be fitted before calling transform()" ) - logger.debug( - "Transforming %d snapshots", - values.shape[0] - ) + logger.debug("Transforming %d snapshots", values.shape[0]) # Transpose for scikit-learn values_T = values.T @@ -137,8 +127,7 @@ def transform(self, values): reduced = reduced_T.T logger.debug( - "Transformation completed, output shape: %s", - reduced.shape + "Transformation completed, output shape: %s", reduced.shape ) return reduced @@ -163,8 +152,7 @@ def inverse_transform(self, reduced_values): ) logger.debug( - "Inverse transforming %d reduced vectors", - reduced_values.shape[0] + "Inverse transforming %d reduced vectors", reduced_values.shape[0] ) # Transpose for scikit-learn @@ -176,7 +164,7 @@ def inverse_transform(self, reduced_values): logger.debug( "Inverse transformation completed, output shape: %s", - reconstructed.shape + reconstructed.shape, ) return reconstructed diff --git a/ezyrb/regular_grid.py b/ezyrb/regular_grid.py index 1036dee8..391d3cdc 100644 --- a/ezyrb/regular_grid.py +++ b/ezyrb/regular_grid.py @@ -1,6 +1,7 @@ """ Module for higher order interpolation on regular grids """ + import logging import numpy as np from scipy.interpolate import RegularGridInterpolator as RGI @@ -52,7 +53,7 @@ class RegularGrid(Approximation): def __init__(self): """ Initialize a RegularGrid interpolator. - + The interpolator will be configured during the fit method. """ self.interpolator = None @@ -75,15 +76,18 @@ def get_grid_axes(self, pts_scrmbld, vals_scrmbld): nN = [] # size of dimension N dim = pts_scrmbld.shape[1] for i in range(dim): - xn, unique_inverse_n = np.unique(pts_scrmbld[:, i], - return_inverse=True) + xn, unique_inverse_n = np.unique( + pts_scrmbld[:, i], return_inverse=True + ) grid_axes.append(xn) nN.append(len(xn)) iN.append(unique_inverse_n) if np.prod(nN) != len(vals_scrmbld): - msg = "Values don't match grid. Make sure to pass a list of "\ - "points that are on a regular grid! Be aware of floating "\ - "point precision." + msg = ( + "Values don't match grid. Make sure to pass a list of " + "points that are on a regular grid! Be aware of floating " + "point precision." + ) raise ValueError(msg) new_row_index = np.ravel_multi_index(iN, nN) reverse_scrambling = np.argsort(new_row_index) @@ -99,14 +103,17 @@ def fit(self, points, values, **kwargs): :param array_like points: the coordinates of the points. :param array_like values: the values in the points. """ - logger.debug("Fitting RegularGrid with points shape: %s, " - "values shape: %s", - np.asarray(points).shape, np.asarray(values).shape) + logger.debug( + "Fitting RegularGrid with points shape: %s, " "values shape: %s", + np.asarray(points).shape, + np.asarray(values).shape, + ) points = np.asarray(points) if not np.issubdtype(points.dtype, np.number): logger.error("Invalid points format/dimension") - raise ValueError('Invalid format or dimension for the argument' - '`points`.') + raise ValueError( + "Invalid format or dimension for the argument" "`points`." + ) if points.ndim == 1: points = points[:, None] logger.debug("Expanded points to 2D") @@ -116,8 +123,7 @@ def fit(self, points, values, **kwargs): shape = [len(ax) for ax in grid_axes] shape.append(-1) logger.debug("Grid shape: %s", shape) - self.interpolator = RGI(grid_axes, values_grd.reshape(shape), - **kwargs) + self.interpolator = RGI(grid_axes, values_grd.reshape(shape), **kwargs) logger.info("RegularGrid fitted successfully") def predict(self, new_point): diff --git a/ezyrb/snapshot.py b/ezyrb/snapshot.py index bc5fd55a..ece11d35 100644 --- a/ezyrb/snapshot.py +++ b/ezyrb/snapshot.py @@ -1,4 +1,4 @@ -""" Module for discretized solution object """ +"""Module for discretized solution object""" import logging import numpy as np @@ -10,16 +10,16 @@ class Snapshot: """ Class for representing a discretized solution snapshot. - + This class encapsulates solution values and their spatial coordinates, providing methods for manipulation and visualization. - + :param array_like values: The solution values. :param array_like space: The spatial coordinates corresponding to the values. Default is None. - + :Example: - + >>> import numpy as np >>> from ezyrb import Snapshot >>> space = np.linspace(0, 1, 50) @@ -34,27 +34,33 @@ class Snapshot: def __init__(self, values, space=None): """ Initialize a Snapshot object. - + :param values: The solution values. Can be a Snapshot instance or an array-like object. :param space: The spatial coordinates. Default is None. """ - logger.debug("Initializing Snapshot with values type: %s, " - "space type: %s", type(values), type(space)) + logger.debug( + "Initializing Snapshot with values type: %s, " "space type: %s", + type(values), + type(space), + ) if isinstance(values, Snapshot): self.values = values.values self.space = values.space - logger.debug("Copied Snapshot with shape: %s", - np.asarray(values.values).shape) + logger.debug( + "Copied Snapshot with shape: %s", + np.asarray(values.values).shape, + ) else: self.values = values self.space = space - logger.debug("Created Snapshot with shape: %s", - np.asarray(values).shape) + logger.debug( + "Created Snapshot with shape: %s", np.asarray(values).shape + ) @property def values(self): - """ + """ Get the snapshot values. """ return self._values @@ -63,26 +69,28 @@ def values(self): def values(self, new_values): """ Set the snapshot values with validation. - + :param array_like new_values: The new snapshot values. :raises ValueError: If the length of new values doesn't match the space. """ - if hasattr(self, 'space') and self.space is not None: + if hasattr(self, "space") and self.space is not None: if new_values is not None and len(self.space) != len(new_values): - logger.error("Space length %d != values length %d", - len(self.space), len(new_values)) - raise ValueError('invalid ndof for the current space.') + logger.error( + "Space length %d != values length %d", + len(self.space), + len(new_values), + ) + raise ValueError("invalid ndof for the current space.") self._values = new_values if new_values is not None: - logger.debug("Snapshot values set with length: %d", - len(new_values)) + logger.debug("Snapshot values set with length: %d", len(new_values)) else: logger.debug("Snapshot values set to None") @property def space(self): - """ + """ Get the snapshot space. """ return self._space @@ -91,26 +99,26 @@ def space(self): def space(self, new_space): """ Set the snapshot space with validation. - + :param array_like new_space: The new spatial coordinates. :raises ValueError: If the length of new space doesn't match the values. """ - if hasattr(self, 'values') and self.values is not None: + if hasattr(self, "values") and self.values is not None: if new_space is not None and len(self.values) != len(new_space): - raise ValueError('invalid ndof for the current space.') + raise ValueError("invalid ndof for the current space.") self._space = new_space @property def flattened(self): - """ return the values in 1D array """ + """return the values in 1D array""" return self.values.flatten() def plot(self): - """ Plot the snapshot, if possible. """ + """Plot the snapshot, if possible.""" if self.space is None: - print('No space set, unable to plot.') + print("No space set, unable to plot.") return if np.asarray(self.space).ndim == 1: diff --git a/pyproject.toml b/pyproject.toml index 6331b040..863409b5 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -24,7 +24,9 @@ requires-python = ">=3.8" [project.optional-dependencies] docs = [ "sphinx", - "sphinx_rtd_theme" + "sphinx_rtd_theme", + "nbsphinx", + "pandoc" ] test = [ "pytest", From 0b4f3a2ababd8c850071178c47f5bd1d374c5ef8 Mon Sep 17 00:00:00 2001 From: Nicola Demo Date: Fri, 9 Jan 2026 11:54:22 +0100 Subject: [PATCH 2/2] black --- ezyrb/approximation/kneighbors_regressor.py | 4 +- ezyrb/database.py | 4 +- ezyrb/parallel/ae_eddl.py | 8 +- ezyrb/plugin/aggregation.py | 5 +- ezyrb/plugin/automatic_shift.py | 5 +- ezyrb/plugin/shift.py | 6 +- ezyrb/reducedordermodel.py | 98 +++++++++------------ ezyrb/reduction/pod.py | 6 +- ezyrb/snapshot.py | 4 +- 9 files changed, 76 insertions(+), 64 deletions(-) diff --git a/ezyrb/approximation/kneighbors_regressor.py b/ezyrb/approximation/kneighbors_regressor.py index 0776289e..e3c64c5f 100644 --- a/ezyrb/approximation/kneighbors_regressor.py +++ b/ezyrb/approximation/kneighbors_regressor.py @@ -33,5 +33,7 @@ def __init__(self, **kwargs): :param kwargs: Arguments passed to sklearn's KNeighborsRegressor. """ - logger.debug("Initializing KNeighborsRegressor with kwargs: %s", kwargs) + logger.debug( + "Initializing KNeighborsRegressor with kwargs: %s", kwargs + ) self.regressor = Regressor(**kwargs) diff --git a/ezyrb/database.py b/ezyrb/database.py index f9810073..53bc15b9 100644 --- a/ezyrb/database.py +++ b/ezyrb/database.py @@ -170,7 +170,9 @@ def split(self, chunks, seed=None): >>> train, test = db.split([80, 20]) # n snapshots """ - logger.debug("Splitting database with chunks=%s, seed=%s", chunks, seed) + logger.debug( + "Splitting database with chunks=%s, seed=%s", chunks, seed + ) if seed is not None: np.random.seed(seed) diff --git a/ezyrb/parallel/ae_eddl.py b/ezyrb/parallel/ae_eddl.py index 9c47f2bd..673c6df2 100644 --- a/ezyrb/parallel/ae_eddl.py +++ b/ezyrb/parallel/ae_eddl.py @@ -165,7 +165,9 @@ def __setstate__(self, state): # (n) hidden layers + (n-1) activation layers + (1) input layer n = len(self.layers_encoder) encoder_layer_index = 2 * n - 1 - self.model_Autoencoder = eddl.import_net_from_onnx_file(self.file_1) + self.model_Autoencoder = eddl.import_net_from_onnx_file( + self.file_1 + ) self.encoder = self.model_Autoencoder.layers[encoder_layer_index] self.decoder = self.model_Autoencoder.layers[-1] @@ -337,7 +339,9 @@ def fit(self, values): for j in range(num_batches): # 1) using next_batch eddl.next_batch([values], [xbatch]) - eddl.train_batch(self.model_Autoencoder, [xbatch], [xbatch]) + eddl.train_batch( + self.model_Autoencoder, [xbatch], [xbatch] + ) losses = eddl.get_losses(self.model_Autoencoder) metrics = eddl.get_metrics(self.model_Autoencoder) diff --git a/ezyrb/plugin/aggregation.py b/ezyrb/plugin/aggregation.py index 1d5f0c40..3e23d0ce 100644 --- a/ezyrb/plugin/aggregation.py +++ b/ezyrb/plugin/aggregation.py @@ -124,7 +124,10 @@ def obj_func(sigma): # minimization procedure res = minimize( - obj_func, x0=sigma_range[0], method="L-BFGS-B", bounds=[sigma_range] + obj_func, + x0=sigma_range[0], + method="L-BFGS-B", + bounds=[sigma_range], ) print("Optimal sigma value in weights: ", res.x) return res.x diff --git a/ezyrb/plugin/automatic_shift.py b/ezyrb/plugin/automatic_shift.py index 9dede5b4..3c74685a 100644 --- a/ezyrb/plugin/automatic_shift.py +++ b/ezyrb/plugin/automatic_shift.py @@ -101,7 +101,10 @@ def _train_shift_network(self, db): np.vstack( [ np.vstack( - [s.space, np.ones(shape=(s.space.shape[0],)) * p.values] + [ + s.space, + np.ones(shape=(s.space.shape[0],)) * p.values, + ] ).T for p, s in db._pairs ] diff --git a/ezyrb/plugin/shift.py b/ezyrb/plugin/shift.py index 47877a52..da8dab6c 100644 --- a/ezyrb/plugin/shift.py +++ b/ezyrb/plugin/shift.py @@ -49,7 +49,11 @@ class ShiftSnapshots(Plugin): """ def __init__( - self, shift_function, interpolator, parameter_index=0, reference_index=0 + self, + shift_function, + interpolator, + parameter_index=0, + reference_index=0, ): """ Initialize the ShiftSnapshots plugin. diff --git a/ezyrb/reducedordermodel.py b/ezyrb/reducedordermodel.py index 03f31f04..6f6ef9b2 100644 --- a/ezyrb/reducedordermodel.py +++ b/ezyrb/reducedordermodel.py @@ -12,7 +12,7 @@ from .approximation import Approximation -from abc import ABC, abstractmethod +from abc import ABC logger = logging.getLogger(__name__) @@ -154,7 +154,8 @@ def database(self, value): if not isinstance(value, (Database, dict)): raise TypeError( - "The database has to be an instance of the Database class, or a dictionary of Database." + "The database has to be an instance of the Database class, or " + "a dictionary of Database." ) self._database = value @@ -171,7 +172,8 @@ def reduction(self): def reduction(self, value): if not isinstance(value, Reduction): raise TypeError( - "The reduction has to be an instance of the Reduction class, or a dict of Reduction." + "The reduction has to be an instance of the Reduction class, " + "or a dict of Reduction." ) self._reduction = value @@ -188,7 +190,8 @@ def approximation(self): def approximation(self, value): if not isinstance(value, Approximation): raise TypeError( - "The approximation has to be an instance of the Approximation class, or a dict of Approximation." + "The approximation has to be an instance of the Approximation " + "class, or a dict of Approximation." ) self._approximation = value @@ -249,7 +252,8 @@ def fit_approximation(self): This method trains the approximation technique on the reduced space representation of the snapshots. - :raises RuntimeError: If the reduced training database has not been created. + :raises RuntimeError: If the reduced training database has not been + created. """ if not hasattr(self, "train_reduced_database"): @@ -260,30 +264,6 @@ def fit_approximation(self): self.train_reduced_database.snapshots_matrix, ) - # if self.n_database == 1 and self.n_reduction == 1: - # self.train_full_database = self.database - # self.reduction.fit(self.database.snapshots_matrix.T) - - # elif self.n_database == 1 and self.n_reduction > 1: - # self.train_full_database = self.database - # for reduction in self.reduction: - # reduction.fit(self.database.snapshots_matrix.T) - - # elif self.n_database > 1 and self.n_reduction == 1: - # self.train_full_database = self.database - # self.reduction = [ - # (k, copy.deepcopy(self.reduction)) - # for k in self.train_full_database - # ] - # print(self.reduction) - # for reduction, database in zip(self.reduction, self.train_full_database): - # self.reduction[reduction].fit(self.train_full_database[database].snapshots_matrix.T) - - # elif self.n_database > 1 and self.n_reduction > 1: - # raise NotImplementedError - # else: - # raise RuntimeError - def fit(self): r""" Calculate reduced space @@ -471,7 +451,9 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): / norm(test.snapshots_matrix, axis=1) ) else: - return np.mean(norm(predicted_test - test.snapshots_matrix, axis=1)) + return np.mean( + norm(predicted_test - test.snapshots_matrix, axis=1) + ) def kfold_cv_error( self, n_splits, *args, norm=np.linalg.norm, relative=True, **kwargs @@ -615,8 +597,8 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is - True. + :param bool relative: True if the error computed is relative. Default + is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. :return: the vector containing the reconstruction errors. @@ -664,8 +646,8 @@ def approximation_error(self, db=None, relative=True, eps=1e-12): :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is - True. + :param bool relative: True if the error computed is relative. Default + is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. @@ -758,7 +740,8 @@ def __init__(self, *args, plugins=None, rom_plugin=None): :param args: Variable arguments for different initialization modes. :param list plugins: Global plugins for the Multi-ROM. Default is None. - :param rom_plugin: Plugin to add to each individual ROM. Default is None. + :param rom_plugin: Plugin to add to each individual ROM. Default is + None. """ if len(args) == 3: @@ -814,7 +797,8 @@ def database(self, value): if not isinstance(value, (Database, dict)): raise TypeError( - "The database has to be an instance of the Database class, or a dictionary of Database." + "The database has to be an instance of the Database class, " + "or a dictionary of Database." ) if isinstance(value, Database): @@ -834,7 +818,8 @@ def reduction(self): def reduction(self, value): if not isinstance(value, (Reduction, dict)): raise TypeError( - "The reduction has to be an instance of the Reduction class, or a dict of Reduction." + "The reduction has to be an instance of the Reduction class" + " or a dict of Reduction." ) if isinstance(value, Reduction): @@ -854,7 +839,8 @@ def approximation(self): def approximation(self, value): if not isinstance(value, (Approximation, dict)): raise TypeError( - "The approximation has to be an instance of the Approximation class, or a dict of Approximation." + "The approximation has to be an instance of the Approximation" + " class, or a dict of Approximation." ) if isinstance(value, Approximation): @@ -884,7 +870,7 @@ def fit(self): self._execute_plugins("fit_preprocessing") for rom_ in self.roms.values(): - if rom_.train_reduced_database == None: + if rom_.train_reduced_database is None: rom_.fit() self._execute_plugins("fit_postprocessing") @@ -892,9 +878,9 @@ def fit(self): def predict(self, parameters=None): """ - Calculate predicted solution for given mu - If gaussians is True, the standard aggregation technique (gaussians) is used, - otherwise the ANNs are used. + Calculate predicted solution for given mu. + If gaussians is True, the standard aggregation technique (gaussians) is + used, otherwise the ANNs are used. :return: the database containing all the predicted solution (with corresponding parameters). @@ -1008,13 +994,13 @@ def test_error(self, test, norm=np.linalg.norm, relative=True): ) / norm(test.snapshots_matrix, axis=1) ) - else: - return np.mean( - norm( - predicted_test.snapshots_matrix - test.snapshots_matrix, - axis=1, - ) + + return np.mean( + norm( + predicted_test.snapshots_matrix - test.snapshots_matrix, + axis=1, ) + ) def kfold_cv_error( self, n_splits, *args, norm=np.linalg.norm, relative=True, **kwargs @@ -1041,7 +1027,8 @@ def kfold_cv_error( kf = KFold(n_splits=n_splits) for train_index, test_index in kf.split(self.database): new_db = self.database[train_index] - # TODO: Fix plugins handling - should pass: plugins=[copy.deepcopy(p) for p in self.plugins] + # TODO: Fix plugins handling - should pass: + # plugins=[copy.deepcopy(p) for p in self.plugins] rom = type(self)( new_db, copy.deepcopy(self.reduction), @@ -1082,14 +1069,15 @@ def loo_error(self, *args, norm=np.linalg.norm, **kwargs): new_db = self.database[indeces] test_db = self.database[~indeces] - # TODO: Fix plugins handling - should pass: plugins=[copy.deepcopy(p) for p in self.plugins] + # TODO: Fix plugins handling - should pass: + # plugins=[copy.deepcopy(p) for p in self.plugins] rom = type(self)( new_db, copy.deepcopy(self.reduction), copy.deepcopy(self.approximation), ).fit() - error[j] = rom.test_error(test_db) + error[j] = rom.test_error(test_db, norm=norm) return error @@ -1158,8 +1146,8 @@ def reduction_error(self, db=None, relative=True, eps=1e-12): :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is - True. + :param bool relative: True if the error computed is relative. Default + is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. :return: the vector containing the reconstruction errors. @@ -1207,8 +1195,8 @@ def approximation_error(self, db=None, relative=True, eps=1e-12): :param database.Database db: the database to use to compute the error. If None, the error is computed on the training database. Default is None. - :param bool relative: True if the error computed is relative. Default is - True. + :param bool relative: True if the error computed is relative. Default + is True. :param float eps: small number to avoid division by zero in relative error computation. Default is 1e-12. diff --git a/ezyrb/reduction/pod.py b/ezyrb/reduction/pod.py index ef77f429..efa995c0 100755 --- a/ezyrb/reduction/pod.py +++ b/ezyrb/reduction/pod.py @@ -61,7 +61,11 @@ def __init__( >>> pod = POD('correlation_matrix', rank=10, save_memory=False) """ logger.debug("Initializing POD with method=%s, rank=%s", method, rank) - self.available_methods = ["svd", "randomized_svd", "correlation_matrix"] + self.available_methods = [ + "svd", + "randomized_svd", + "correlation_matrix", + ] self.rank = rank if method == "svd": self._method = self._svd diff --git a/ezyrb/snapshot.py b/ezyrb/snapshot.py index ece11d35..d1712941 100644 --- a/ezyrb/snapshot.py +++ b/ezyrb/snapshot.py @@ -84,7 +84,9 @@ def values(self, new_values): self._values = new_values if new_values is not None: - logger.debug("Snapshot values set with length: %d", len(new_values)) + logger.debug( + "Snapshot values set with length: %d", len(new_values) + ) else: logger.debug("Snapshot values set to None")