From 285f0f3d7da8494ad9955d8f39d9b7f0bc294d4f Mon Sep 17 00:00:00 2001 From: Jean-Philippe ARGAUD Date: Thu, 27 Jun 2024 11:50:33 +0200 Subject: [PATCH] Minor documentation and source update for compatibilities --- doc/en/snippets/ModuleCompatibility.rst | 2 +- doc/fr/snippets/ModuleCompatibility.rst | 2 +- src/daComposant/daAlgorithms/Atoms/ecwnlls.py | 2 + .../daAlgorithms/Atoms/lbfgsb114hlt.py | 553 ++++++++++++++++++ .../daAlgorithms/Atoms/psas3dvar.py | 2 + .../daAlgorithms/Atoms/std3dvar.py | 2 + .../daAlgorithms/Atoms/std4dvar.py | 2 + .../daAlgorithms/Atoms/van3dvar.py | 2 + .../pst4mod/modelica_calibration/case.py | 90 ++- .../modelica_libraries/around_simulation.py | 10 +- 10 files changed, 613 insertions(+), 54 deletions(-) create mode 100644 src/daComposant/daAlgorithms/Atoms/lbfgsb114hlt.py diff --git a/doc/en/snippets/ModuleCompatibility.rst b/doc/en/snippets/ModuleCompatibility.rst index 60cf2a5..dc99912 100644 --- a/doc/en/snippets/ModuleCompatibility.rst +++ b/doc/en/snippets/ModuleCompatibility.rst @@ -16,7 +16,7 @@ versions within the range described below. Python, 3.6.5, 3.12.3 Numpy, 1.14.3, 1.26.4 - Scipy, 0.19.1, 1.13.1 + Scipy, 0.19.1, 1.14.0 MatplotLib, 2.2.2, 3.8.4 GnuplotPy, 1.8, 1.8 NLopt, 2.4.2, 2.7.1 diff --git a/doc/fr/snippets/ModuleCompatibility.rst b/doc/fr/snippets/ModuleCompatibility.rst index 0eecbfd..04b5748 100644 --- a/doc/fr/snippets/ModuleCompatibility.rst +++ b/doc/fr/snippets/ModuleCompatibility.rst @@ -17,7 +17,7 @@ l'étendue décrite ci-dessous. Python, 3.6.5, 3.12.3 Numpy, 1.14.3, 1.26.4 - Scipy, 0.19.1, 1.13.1 + Scipy, 0.19.1, 1.14.0 MatplotLib, 2.2.2, 3.8.4 GnuplotPy, 1.8, 1.8 NLopt, 2.4.2, 2.7.1 diff --git a/src/daComposant/daAlgorithms/Atoms/ecwnlls.py b/src/daComposant/daAlgorithms/Atoms/ecwnlls.py index b6ed31a..2762c02 100644 --- a/src/daComposant/daAlgorithms/Atoms/ecwnlls.py +++ b/src/daComposant/daAlgorithms/Atoms/ecwnlls.py @@ -148,6 +148,8 @@ def ecwnlls(selfA, Xb, Y, U, HO, CM, R, B, __storeState = False): import daAlgorithms.Atoms.lbfgsb112hlt as optimiseur elif vt("1.13.0") <= vt(scipy.version.version) <= vt("1.13.99"): import daAlgorithms.Atoms.lbfgsb113hlt as optimiseur + elif vt("1.14.0") <= vt(scipy.version.version) <= vt("1.14.99"): + import daAlgorithms.Atoms.lbfgsb114hlt as optimiseur else: import scipy.optimize as optimiseur Minimum, J_optimal, Informations = optimiseur.fmin_l_bfgs_b( diff --git a/src/daComposant/daAlgorithms/Atoms/lbfgsb114hlt.py b/src/daComposant/daAlgorithms/Atoms/lbfgsb114hlt.py new file mode 100644 index 0000000..9e4f5be --- /dev/null +++ b/src/daComposant/daAlgorithms/Atoms/lbfgsb114hlt.py @@ -0,0 +1,553 @@ +# Modification de la version 1.14.0 +# flake8: noqa +""" +Functions +--------- +.. autosummary:: + :toctree: generated/ + + fmin_l_bfgs_b + +""" + +## License for the Python wrapper +## ============================== + +## Copyright (c) 2004 David M. Cooke + +## Permission is hereby granted, free of charge, to any person obtaining a +## copy of this software and associated documentation files (the "Software"), +## to deal in the Software without restriction, including without limitation +## the rights to use, copy, modify, merge, publish, distribute, sublicense, +## and/or sell copies of the Software, and to permit persons to whom the +## Software is furnished to do so, subject to the following conditions: + +## The above copyright notice and this permission notice shall be included in +## all copies or substantial portions of the Software. + +## THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +## IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +## FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +## AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +## LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +## FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER +## DEALINGS IN THE SOFTWARE. + +## Modifications by Travis Oliphant and Enthought, Inc. for inclusion in SciPy + +import numpy as np +from numpy import array, asarray, float64, zeros +from scipy.optimize import _lbfgsb +from scipy.optimize._optimize import (MemoizeJac, OptimizeResult, _call_callback_maybe_halt, + _wrap_callback, _check_unknown_options, + _prepare_scalar_function) +from scipy.optimize._constraints import old_bound_to_new + +from scipy.sparse.linalg import LinearOperator + +__all__ = ['fmin_l_bfgs_b', 'LbfgsInvHessProduct'] + + +def fmin_l_bfgs_b(func, x0, fprime=None, args=(), + approx_grad=0, + bounds=None, m=10, factr=1e7, pgtol=1e-5, + epsilon=1e-8, + iprint=-1, maxfun=15000, maxiter=15000, disp=None, + callback=None, maxls=20): + """ + Minimize a function func using the L-BFGS-B algorithm. + + Parameters + ---------- + func : callable f(x,*args) + Function to minimize. + x0 : ndarray + Initial guess. + fprime : callable fprime(x,*args), optional + The gradient of `func`. If None, then `func` returns the function + value and the gradient (``f, g = func(x, *args)``), unless + `approx_grad` is True in which case `func` returns only ``f``. + args : sequence, optional + Arguments to pass to `func` and `fprime`. + approx_grad : bool, optional + Whether to approximate the gradient numerically (in which case + `func` returns only the function value). + bounds : list, optional + ``(min, max)`` pairs for each element in ``x``, defining + the bounds on that parameter. Use None or +-inf for one of ``min`` or + ``max`` when there is no bound in that direction. + m : int, optional + The maximum number of variable metric corrections + used to define the limited memory matrix. (The limited memory BFGS + method does not store the full hessian but uses this many terms in an + approximation to it.) + factr : float, optional + The iteration stops when + ``(f^k - f^{k+1})/max{|f^k|,|f^{k+1}|,1} <= factr * eps``, + where ``eps`` is the machine precision, which is automatically + generated by the code. Typical values for `factr` are: 1e12 for + low accuracy; 1e7 for moderate accuracy; 10.0 for extremely + high accuracy. See Notes for relationship to `ftol`, which is exposed + (instead of `factr`) by the `scipy.optimize.minimize` interface to + L-BFGS-B. + pgtol : float, optional + The iteration will stop when + ``max{|proj g_i | i = 1, ..., n} <= pgtol`` + where ``proj g_i`` is the i-th component of the projected gradient. + epsilon : float, optional + Step size used when `approx_grad` is True, for numerically + calculating the gradient + iprint : int, optional + Controls the frequency of output. ``iprint < 0`` means no output; + ``iprint = 0`` print only one line at the last iteration; + ``0 < iprint < 99`` print also f and ``|proj g|`` every iprint iterations; + ``iprint = 99`` print details of every iteration except n-vectors; + ``iprint = 100`` print also the changes of active set and final x; + ``iprint > 100`` print details of every iteration including x and g. + disp : int, optional + If zero, then no output. If a positive number, then this over-rides + `iprint` (i.e., `iprint` gets the value of `disp`). + maxfun : int, optional + Maximum number of function evaluations. Note that this function + may violate the limit because of evaluating gradients by numerical + differentiation. + maxiter : int, optional + Maximum number of iterations. + callback : callable, optional + Called after each iteration, as ``callback(xk)``, where ``xk`` is the + current parameter vector. + maxls : int, optional + Maximum number of line search steps (per iteration). Default is 20. + + Returns + ------- + x : array_like + Estimated position of the minimum. + f : float + Value of `func` at the minimum. + d : dict + Information dictionary. + + * d['warnflag'] is + + - 0 if converged, + - 1 if too many function evaluations or too many iterations, + - 2 if stopped for another reason, given in d['task'] + + * d['grad'] is the gradient at the minimum (should be 0 ish) + * d['funcalls'] is the number of function calls made. + * d['nit'] is the number of iterations. + + See also + -------- + minimize: Interface to minimization algorithms for multivariate + functions. See the 'L-BFGS-B' `method` in particular. Note that the + `ftol` option is made available via that interface, while `factr` is + provided via this interface, where `factr` is the factor multiplying + the default machine floating-point precision to arrive at `ftol`: + ``ftol = factr * numpy.finfo(float).eps``. + + Notes + ----- + License of L-BFGS-B (FORTRAN code): + + The version included here (in fortran code) is 3.0 + (released April 25, 2011). It was written by Ciyou Zhu, Richard Byrd, + and Jorge Nocedal . It carries the following + condition for use: + + This software is freely available, but we expect that all publications + describing work using this software, or all commercial products using it, + quote at least one of the references given below. This software is released + under the BSD License. + + References + ---------- + * R. H. Byrd, P. Lu and J. Nocedal. A Limited Memory Algorithm for Bound + Constrained Optimization, (1995), SIAM Journal on Scientific and + Statistical Computing, 16, 5, pp. 1190-1208. + * C. Zhu, R. H. Byrd and J. Nocedal. L-BFGS-B: Algorithm 778: L-BFGS-B, + FORTRAN routines for large scale bound constrained optimization (1997), + ACM Transactions on Mathematical Software, 23, 4, pp. 550 - 560. + * J.L. Morales and J. Nocedal. L-BFGS-B: Remark on Algorithm 778: L-BFGS-B, + FORTRAN routines for large scale bound constrained optimization (2011), + ACM Transactions on Mathematical Software, 38, 1. + + Examples + -------- + Solve a linear regression problem via `fmin_l_bfgs_b`. To do this, first we define + an objective function ``f(m, b) = (y - y_model)**2``, where `y` describes the + observations and `y_model` the prediction of the linear model as + ``y_model = m*x + b``. The bounds for the parameters, ``m`` and ``b``, are arbitrarily + chosen as ``(0,5)`` and ``(5,10)`` for this example. + + >>> import numpy as np + >>> from scipy.optimize import fmin_l_bfgs_b + >>> X = np.arange(0, 10, 1) + >>> M = 2 + >>> B = 3 + >>> Y = M * X + B + >>> def func(parameters, *args): + ... x = args[0] + ... y = args[1] + ... m, b = parameters + ... y_model = m*x + b + ... error = sum(np.power((y - y_model), 2)) + ... return error + + >>> initial_values = np.array([0.0, 1.0]) + + >>> x_opt, f_opt, info = fmin_l_bfgs_b(func, x0=initial_values, args=(X, Y), + ... approx_grad=True) + >>> x_opt, f_opt + array([1.99999999, 3.00000006]), 1.7746231151323805e-14 # may vary + + The optimized parameters in ``x_opt`` agree with the ground truth parameters + ``m`` and ``b``. Next, let us perform a bound contrained optimization using the `bounds` + parameter. + + >>> bounds = [(0, 5), (5, 10)] + >>> x_opt, f_op, info = fmin_l_bfgs_b(func, x0=initial_values, args=(X, Y), + ... approx_grad=True, bounds=bounds) + >>> x_opt, f_opt + array([1.65990508, 5.31649385]), 15.721334516453945 # may vary + """ + # handle fprime/approx_grad + if approx_grad: + fun = func + jac = None + elif fprime is None: + fun = MemoizeJac(func) + jac = fun.derivative + else: + fun = func + jac = fprime + + # build options + if disp is None: + disp = iprint + callback = _wrap_callback(callback) + opts = {'disp': disp, + 'iprint': iprint, + 'maxcor': m, + 'ftol': factr * np.finfo(float).eps, + 'gtol': pgtol, + 'eps': epsilon, + 'maxfun': maxfun, + 'maxiter': maxiter, + 'callback': callback, + 'maxls': maxls} + + res = _minimize_lbfgsb(fun, x0, args=args, jac=jac, bounds=bounds, + **opts) + d = {'grad': res['jac'], + 'task': res['message'], + 'funcalls': res['nfev'], + 'nit': res['nit'], + 'warnflag': res['status']} + f = res['fun'] + x = res['x'] + + return x, f, d + + +def _minimize_lbfgsb(fun, x0, args=(), jac=None, bounds=None, + disp=None, maxcor=10, ftol=2.2204460492503131e-09, + gtol=1e-5, eps=1e-8, maxfun=15000, maxiter=15000, + iprint=-1, callback=None, maxls=20, + finite_diff_rel_step=None, **unknown_options): + """ + Minimize a scalar function of one or more variables using the L-BFGS-B + algorithm. + + Options + ------- + disp : None or int + If `disp is None` (the default), then the supplied version of `iprint` + is used. If `disp is not None`, then it overrides the supplied version + of `iprint` with the behaviour you outlined. + maxcor : int + The maximum number of variable metric corrections used to + define the limited memory matrix. (The limited memory BFGS + method does not store the full hessian but uses this many terms + in an approximation to it.) + ftol : float + The iteration stops when ``(f^k - + f^{k+1})/max{|f^k|,|f^{k+1}|,1} <= ftol``. + gtol : float + The iteration will stop when ``max{|proj g_i | i = 1, ..., n} + <= gtol`` where ``proj g_i`` is the i-th component of the + projected gradient. + eps : float or ndarray + If `jac is None` the absolute step size used for numerical + approximation of the jacobian via forward differences. + maxfun : int + Maximum number of function evaluations. Note that this function + may violate the limit because of evaluating gradients by numerical + differentiation. + maxiter : int + Maximum number of iterations. + iprint : int, optional + Controls the frequency of output. ``iprint < 0`` means no output; + ``iprint = 0`` print only one line at the last iteration; + ``0 < iprint < 99`` print also f and ``|proj g|`` every iprint iterations; + ``iprint = 99`` print details of every iteration except n-vectors; + ``iprint = 100`` print also the changes of active set and final x; + ``iprint > 100`` print details of every iteration including x and g. + callback : callable, optional + Called after each iteration, as ``callback(xk)``, where ``xk`` is the + current parameter vector. + maxls : int, optional + Maximum number of line search steps (per iteration). Default is 20. + finite_diff_rel_step : None or array_like, optional + If `jac in ['2-point', '3-point', 'cs']` the relative step size to + use for numerical approximation of the jacobian. The absolute step + size is computed as ``h = rel_step * sign(x) * max(1, abs(x))``, + possibly adjusted to fit into the bounds. For ``method='3-point'`` + the sign of `h` is ignored. If None (default) then step is selected + automatically. + + Notes + ----- + The option `ftol` is exposed via the `scipy.optimize.minimize` interface, + but calling `scipy.optimize.fmin_l_bfgs_b` directly exposes `factr`. The + relationship between the two is ``ftol = factr * numpy.finfo(float).eps``. + I.e., `factr` multiplies the default machine floating-point precision to + arrive at `ftol`. + + """ + _check_unknown_options(unknown_options) + m = maxcor + pgtol = gtol + factr = ftol / np.finfo(float).eps + + x0 = asarray(x0).ravel() + n, = x0.shape + + # historically old-style bounds were/are expected by lbfgsb. + # That's still the case but we'll deal with new-style from here on, + # it's easier + if bounds is None: + pass + elif len(bounds) != n: + raise ValueError('length of x0 != length of bounds') + else: + bounds = np.array(old_bound_to_new(bounds)) + + # check bounds + if (bounds[0] > bounds[1]).any(): + raise ValueError( + "LBFGSB - one of the lower bounds is greater than an upper bound." + ) + + # initial vector must lie within the bounds. Otherwise ScalarFunction and + # approx_derivative will cause problems + x0 = np.clip(x0, bounds[0], bounds[1]) + + if disp is not None: + if disp == 0: + iprint = -1 + else: + iprint = disp + + # _prepare_scalar_function can use bounds=None to represent no bounds + sf = _prepare_scalar_function(fun, x0, jac=jac, args=args, epsilon=eps, + bounds=bounds, + finite_diff_rel_step=finite_diff_rel_step) + + func_and_grad = sf.fun_and_grad + + fortran_int = _lbfgsb.types.intvar.dtype + + nbd = zeros(n, fortran_int) + low_bnd = zeros(n, float64) + upper_bnd = zeros(n, float64) + bounds_map = {(-np.inf, np.inf): 0, + (1, np.inf): 1, + (1, 1): 2, + (-np.inf, 1): 3} + + if bounds is not None: + for i in range(0, n): + l, u = bounds[0, i], bounds[1, i] + if not np.isinf(l): + low_bnd[i] = l + l = 1 + if not np.isinf(u): + upper_bnd[i] = u + u = 1 + nbd[i] = bounds_map[l, u] + + if not maxls > 0: + raise ValueError('maxls must be positive.') + + x = array(x0, float64) + f = array(0.0, float64) + g = zeros((n,), float64) + wa = zeros(2*m*n + 5*n + 11*m*m + 8*m, float64) + iwa = zeros(3*n, fortran_int) + task = zeros(1, 'S60') + csave = zeros(1, 'S60') + lsave = zeros(4, fortran_int) + isave = zeros(44, fortran_int) + dsave = zeros(29, float64) + + task[:] = 'START' + + n_iterations = 0 + + while 1: + # g may become float32 if a user provides a function that calculates + # the Jacobian in float32 (see gh-18730). The underlying Fortran code + # expects float64, so upcast it + g = g.astype(np.float64) + # x, f, g, wa, iwa, task, csave, lsave, isave, dsave = \ + _lbfgsb.setulb(m, x, low_bnd, upper_bnd, nbd, f, g, factr, + pgtol, wa, iwa, task, iprint, csave, lsave, + isave, dsave, maxls) + task_str = task.tobytes() + if task_str.startswith(b'FG'): + # The minimization routine wants f and g at the current x. + # Note that interruptions due to maxfun are postponed + # until the completion of the current minimization iteration. + # Overwrite f and g: + f, g = func_and_grad(x) + if sf.nfev > maxfun: + task[:] = ('STOP: TOTAL NO. of f AND g EVALUATIONS ' + 'EXCEEDS LIMIT') + elif task_str.startswith(b'NEW_X'): + # new iteration + n_iterations += 1 + + intermediate_result = OptimizeResult(x=x, fun=f) + if _call_callback_maybe_halt(callback, intermediate_result): + task[:] = 'STOP: CALLBACK REQUESTED HALT' + if n_iterations >= maxiter: + task[:] = 'STOP: TOTAL NO. of ITERATIONS REACHED LIMIT' + elif sf.nfev > maxfun: + task[:] = ('STOP: TOTAL NO. of f AND g EVALUATIONS ' + 'EXCEEDS LIMIT') + else: + break + + task_str = task.tobytes().strip(b'\x00').strip() + if task_str.startswith(b'CONV'): + warnflag = 0 + elif sf.nfev > maxfun or n_iterations >= maxiter: + warnflag = 1 + else: + warnflag = 2 + + # These two portions of the workspace are described in the mainlb + # subroutine in lbfgsb.f. See line 363. + s = wa[0: m*n].reshape(m, n) + y = wa[m*n: 2*m*n].reshape(m, n) + + # See lbfgsb.f line 160 for this portion of the workspace. + # isave(31) = the total number of BFGS updates prior the current iteration; + n_bfgs_updates = isave[30] + + n_corrs = min(n_bfgs_updates, maxcor) + hess_inv = LbfgsInvHessProduct(s[:n_corrs], y[:n_corrs]) + + task_str = task_str.decode() + return OptimizeResult(fun=f, jac=g, nfev=sf.nfev, + njev=sf.ngev, + nit=n_iterations, status=warnflag, message=task_str, + x=x, success=(warnflag == 0), hess_inv=hess_inv) + + +class LbfgsInvHessProduct(LinearOperator): + """Linear operator for the L-BFGS approximate inverse Hessian. + + This operator computes the product of a vector with the approximate inverse + of the Hessian of the objective function, using the L-BFGS limited + memory approximation to the inverse Hessian, accumulated during the + optimization. + + Objects of this class implement the ``scipy.sparse.linalg.LinearOperator`` + interface. + + Parameters + ---------- + sk : array_like, shape=(n_corr, n) + Array of `n_corr` most recent updates to the solution vector. + (See [1]). + yk : array_like, shape=(n_corr, n) + Array of `n_corr` most recent updates to the gradient. (See [1]). + + References + ---------- + .. [1] Nocedal, Jorge. "Updating quasi-Newton matrices with limited + storage." Mathematics of computation 35.151 (1980): 773-782. + + """ + + def __init__(self, sk, yk): + """Construct the operator.""" + if sk.shape != yk.shape or sk.ndim != 2: + raise ValueError('sk and yk must have matching shape, (n_corrs, n)') + n_corrs, n = sk.shape + + super().__init__(dtype=np.float64, shape=(n, n)) + + self.sk = sk + self.yk = yk + self.n_corrs = n_corrs + self.rho = 1 / np.einsum('ij,ij->i', sk, yk) + + def _matvec(self, x): + """Efficient matrix-vector multiply with the BFGS matrices. + + This calculation is described in Section (4) of [1]. + + Parameters + ---------- + x : ndarray + An array with shape (n,) or (n,1). + + Returns + ------- + y : ndarray + The matrix-vector product + + """ + s, y, n_corrs, rho = self.sk, self.yk, self.n_corrs, self.rho + q = np.array(x, dtype=self.dtype, copy=True) + if q.ndim == 2 and q.shape[1] == 1: + q = q.reshape(-1) + + alpha = np.empty(n_corrs) + + for i in range(n_corrs-1, -1, -1): + alpha[i] = rho[i] * np.dot(s[i], q) + q = q - alpha[i]*y[i] + + r = q + for i in range(n_corrs): + beta = rho[i] * np.dot(y[i], r) + r = r + s[i] * (alpha[i] - beta) + + return r + + def todense(self): + """Return a dense array representation of this operator. + + Returns + ------- + arr : ndarray, shape=(n, n) + An array with the same shape and containing + the same data represented by this `LinearOperator`. + + """ + s, y, n_corrs, rho = self.sk, self.yk, self.n_corrs, self.rho + I = np.eye(*self.shape, dtype=self.dtype) + Hk = I + + for i in range(n_corrs): + A1 = I - s[i][:, np.newaxis] * y[i][np.newaxis, :] * rho[i] + A2 = I - y[i][:, np.newaxis] * s[i][np.newaxis, :] * rho[i] + + Hk = np.dot(A1, np.dot(Hk, A2)) + (rho[i] * s[i][:, np.newaxis] * + s[i][np.newaxis, :]) + return Hk diff --git a/src/daComposant/daAlgorithms/Atoms/psas3dvar.py b/src/daComposant/daAlgorithms/Atoms/psas3dvar.py index 2e31e58..9fc496e 100644 --- a/src/daComposant/daAlgorithms/Atoms/psas3dvar.py +++ b/src/daComposant/daAlgorithms/Atoms/psas3dvar.py @@ -130,6 +130,8 @@ def psas3dvar(selfA, Xb, Y, U, HO, CM, R, B, __storeState = False): import daAlgorithms.Atoms.lbfgsb112hlt as optimiseur elif vt("1.13.0") <= vt(scipy.version.version) <= vt("1.13.99"): import daAlgorithms.Atoms.lbfgsb113hlt as optimiseur + elif vt("1.14.0") <= vt(scipy.version.version) <= vt("1.14.99"): + import daAlgorithms.Atoms.lbfgsb114hlt as optimiseur else: import scipy.optimize as optimiseur Minimum, J_optimal, Informations = optimiseur.fmin_l_bfgs_b( diff --git a/src/daComposant/daAlgorithms/Atoms/std3dvar.py b/src/daComposant/daAlgorithms/Atoms/std3dvar.py index 303845d..db77218 100644 --- a/src/daComposant/daAlgorithms/Atoms/std3dvar.py +++ b/src/daComposant/daAlgorithms/Atoms/std3dvar.py @@ -133,6 +133,8 @@ def std3dvar(selfA, Xb, Y, U, HO, CM, R, B, __storeState = False): import daAlgorithms.Atoms.lbfgsb112hlt as optimiseur elif vt("1.13.0") <= vt(scipy.version.version) <= vt("1.13.99"): import daAlgorithms.Atoms.lbfgsb113hlt as optimiseur + elif vt("1.14.0") <= vt(scipy.version.version) <= vt("1.14.99"): + import daAlgorithms.Atoms.lbfgsb114hlt as optimiseur else: import scipy.optimize as optimiseur Minimum, J_optimal, Informations = optimiseur.fmin_l_bfgs_b( diff --git a/src/daComposant/daAlgorithms/Atoms/std4dvar.py b/src/daComposant/daAlgorithms/Atoms/std4dvar.py index 6af5464..0c4cbdd 100644 --- a/src/daComposant/daAlgorithms/Atoms/std4dvar.py +++ b/src/daComposant/daAlgorithms/Atoms/std4dvar.py @@ -194,6 +194,8 @@ def std4dvar(selfA, Xb, Y, U, HO, EM, CM, R, B, Q): import daAlgorithms.Atoms.lbfgsb112hlt as optimiseur elif vt("1.13.0") <= vt(scipy.version.version) <= vt("1.13.99"): import daAlgorithms.Atoms.lbfgsb113hlt as optimiseur + elif vt("1.14.0") <= vt(scipy.version.version) <= vt("1.14.99"): + import daAlgorithms.Atoms.lbfgsb114hlt as optimiseur else: import scipy.optimize as optimiseur Minimum, J_optimal, Informations = optimiseur.fmin_l_bfgs_b( diff --git a/src/daComposant/daAlgorithms/Atoms/van3dvar.py b/src/daComposant/daAlgorithms/Atoms/van3dvar.py index 34eb24d..b743fff 100644 --- a/src/daComposant/daAlgorithms/Atoms/van3dvar.py +++ b/src/daComposant/daAlgorithms/Atoms/van3dvar.py @@ -141,6 +141,8 @@ def van3dvar(selfA, Xb, Y, U, HO, CM, R, B, __storeState = False): import daAlgorithms.Atoms.lbfgsb112hlt as optimiseur elif vt("1.13.0") <= vt(scipy.version.version) <= vt("1.13.99"): import daAlgorithms.Atoms.lbfgsb113hlt as optimiseur + elif vt("1.14.0") <= vt(scipy.version.version) <= vt("1.14.99"): + import daAlgorithms.Atoms.lbfgsb114hlt as optimiseur else: import scipy.optimize as optimiseur Minimum, J_optimal, Informations = optimiseur.fmin_l_bfgs_b( diff --git a/src/daComposant/daNumerics/pst4mod/modelica_calibration/case.py b/src/daComposant/daNumerics/pst4mod/modelica_calibration/case.py index f605b9e..e2bd861 100644 --- a/src/daComposant/daNumerics/pst4mod/modelica_calibration/case.py +++ b/src/daComposant/daNumerics/pst4mod/modelica_calibration/case.py @@ -24,15 +24,15 @@ __all__ = ["Calibration", "name", "version", "year"] # ============================================================================== name = "Modelica and Dymola Calibration Tools" -version = "1.0.9.7.0" # "x.x"+"adao version" +version = "1.0.9.13.0" # "x.x"+"adao version" year = "2021" # ============================================================================== # Default configuration # --------------------- -import os, sys, io, shutil, numpy, scipy, tempfile, time, logging, pandas, subprocess, copy, csv -from datetime import timedelta -from datetime import datetime +import os, sys, io, shutil, time, logging, subprocess, copy, csv # noqa: E402 +import tempfile, warnings, numpy, scipy, pandas # noqa: E402 +from datetime import datetime # noqa: E402 try: import adao @@ -58,9 +58,10 @@ except ImportError: try: from fmpy import simulate_fmu from fmpy.fmi1 import FMICallException -# except ImportError: -# raise ImportError("fmpy library not found, please install it first") -except: +except ImportError: + __msg = "fmpy library not found, it will not be possible for you to simulate Functional Mock-up Units (FMUs) as models. Add it to Python installation it if necessary." + warnings.warn(__msg, ImportWarning, stacklevel=50) +except Exception: pass _configuration = { @@ -69,12 +70,7 @@ _configuration = { "DirectoryMethods" : "Methods", "DirectoryResults" : "Results", "Launcher" : "configuration.py", - } -__paths = [] -__paths.append( os.path.abspath(os.path.join(os.path.dirname(__file__), '..','..','buildingspy-1.5.0')) ) -__paths.append( os.path.abspath(os.path.join(os.path.dirname(__file__), '..','..','Packages')) ) -for __path in __paths: - if os.path.exists(__path): sys.path.append(__path) +} # ============================================================================== class Calibration(object): @@ -95,12 +91,12 @@ class Calibration(object): self.__verbose = bool(Verbose) self.__stdoutid = sys.stdout if SaveStdoutIn is not None: - sys.stdout = open(SaveStdoutIn,"w") + sys.stdout = open(SaveStdoutIn, "w") if self.__verbose: __msg = "[VERBOSE] %s"%self.__name print("") print(" %s"%__msg) - print(" %s"%("-"*len(__msg),)) + print(" %s"%("-" * len(__msg),)) print("") VersionsLogicielles() @@ -111,10 +107,10 @@ class Calibration(object): """ if Configuration is None or type(Configuration) is not dict: Configuration = _configuration - for __k,__v in Configuration.items(): + for __k, __v in Configuration.items(): setattr(self, __k, __v) - def _setModelTmpDir(self, __model_name="dsin.txt", __model_format="DYMOSIM", __model_dest = "dsin.txt", __before_tmp = None, __suffix=None): + def _setModelTmpDir(self, __model_name = "dsin.txt", __model_format = "DYMOSIM", __model_dest = "dsin.txt", __before_tmp = None, __suffix=None): if __model_name is None: raise ValueError('Model file or directory has to be set and can not be None') elif os.path.isfile(__model_name): @@ -128,21 +124,21 @@ class Calibration(object): else: raise ValueError('Model file or directory not found using %s'%str(__model_name)) # - if __before_tmp is not None: # Mettre "../.." si nécessaire + if __before_tmp is not None: # Mettre "../.." si nécessaire __mtmp = os.path.join(__model_dir, __before_tmp, "tmp") else: __mtmp = os.path.join(__model_dir, "tmp") if not os.path.exists(__mtmp): os.mkdir(__mtmp) - __prefix = time.strftime('%Y%m%d_%Hh%Mm%Ss_tmp_',time.localtime()) + __prefix = time.strftime('%Y%m%d_%Hh%Mm%Ss_tmp_', time.localtime()) if __suffix is not None: __prefix = __prefix + __suffix + "_" __ltmp = tempfile.mkdtemp( prefix=__prefix, dir=__mtmp ) # for sim, pmf, dst in ( - (__model_nam,__model_format,__model_dst), - ("dymosim.exe","DYMOSIM","dymosim.exe"), - ("pythonsim.exe","PYSIM","pythonsim.exe"), + (__model_nam, __model_format, __model_dst), + ("dymosim.exe", "DYMOSIM", "dymosim.exe"), + ("pythonsim.exe", "PYSIM", "pythonsim.exe"), ): # Recherche un fichier de données ou un simulateur dans cet ordre @@ -150,7 +146,7 @@ class Calibration(object): shutil.copy( os.path.join(__model_dir, sim), os.path.join(__ltmp, dst) - ) + ) # _secure_copy_textfile( # os.path.join(__model_dir, sim), @@ -710,7 +706,7 @@ class Calibration(object): try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs = dict(zip(LVariablesToChange, var_int)) @@ -734,7 +730,7 @@ class Calibration(object): debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt') try: shutil.copy(dslog_file,debug_model_file) - except: + except Exception: pass if auto_simul.success_code == 2 : @@ -818,7 +814,7 @@ class Calibration(object): try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs = dict(zip(LVariablesToChange, var_int)) @@ -859,7 +855,7 @@ class Calibration(object): debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt') try: shutil.copy(dslog_file,debug_model_file) - except: + except Exception: pass if not(os.path.exists(os.path.join(simudir, 'RES','1.0.mat'))): @@ -908,7 +904,7 @@ class Calibration(object): # try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs = dict(zip(LVariablesToChange, var_int)) writer_no_CWP = write_in_dsin.Write_in_dsin(dict_inputs = dict_inputs, filedir=simudir, dsin_name=str("dsin_" +LNames[etat]+ ".txt"), old_file_name='old_dsin.txt',new_file_name = str("dsin_" +LNames[etat]+ ".txt")) @@ -1282,7 +1278,7 @@ class Calibration(object): simulation_results = TOP_LEVEL_exeOpenModelicaMultiobs(x_whole_values_with_bg, KeepCalculationFolders = KeepCalculationFolders, VariablesToCalibrate=VariablesToCalibrate, OutputVariables=OutputVariables, LNames = LNames, LColumns = LColumns, LVariablesToChange=LVariablesToChange, ref_simudir = self._get_Name_ModelTmpDir_simple(ModelName), ModelName = ModelName, List_Multideltatime = List_Multideltatime, Linux = Linux, AdvancedDebugModel = AdvancedDebugModel, TimeoutModelExecution = TimeoutModelExecution) else: raise NotImplementedError("Not yet implemented for current model format: ", model_format) - except: + except Exception: list_failed_model_evaluation.append(x_value_bg_tested) dict_check_model_stability[x_bg_tested] = list_failed_model_evaluation @@ -1542,7 +1538,7 @@ class Calibration(object): try: os.remove(os.path.join(dir_for_dsin_opti,str("dsin_" +LNames[etat]+ "_optimal.txt"))) - except: + except Exception: pass shutil.copyfile( @@ -1568,14 +1564,14 @@ class Calibration(object): else: try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs_boundary_conditions_optimal_params = dict(zip(LVariablesToChange, var_int)) dict_inputs_boundary_conditions_optimal_params.update(dict_optimal_params) try: os.remove(os.path.join(dir_for_dsin_opti,str("dsin_" +LNames[etat]+ "_optimal.txt"))) - except: + except Exception: pass dir_for_dsin_opti = os.path.abspath(os.path.join(self._get_Name_ModelTmpDir_simple(ModelName),os.pardir)) shutil.copyfile( @@ -1628,7 +1624,7 @@ class Calibration(object): outwriter.writerow(row) try: os.remove(optimal_file_results_name_cl_tmp) - except: + except Exception: pass if InitialSimulation: #If initial simulation is required as well, the Initial results files are given as well @@ -1659,7 +1655,7 @@ class Calibration(object): outwriter.writerow(row) try: os.remove(initial_file_results_name_cl_tmp) - except: + except Exception: pass return __resultats @@ -1687,7 +1683,7 @@ def _readLink(__filename=None, __colnames=None, __indexname="Variable", __format try: #Solution provisoire pou gérer les cas sans CL précisées par l'utilisateur with ImportFromFile(__filename, __colnames, __indexname, __format, False) as reading: colnames, columns, indexname, variablestochange = reading.getvalue() - except: + except Exception: colnames =__colnames columns =[] variablestochange =[] @@ -2101,7 +2097,7 @@ def TOP_LEVEL_exefmuMultiobs( x_values_matrix , VariablesToCalibrate=None, Outpu try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs = dict(zip(LVariablesToChange, var_int)) @@ -2150,7 +2146,7 @@ def TOP_LEVEL_exefmuMultiobs( x_values_matrix , VariablesToCalibrate=None, Outpu 'log_debug_model.txt') try: #try toremove the previous file os.remove(log_file) - except: + except Exception: pass f=open(log_file,'a') @@ -2242,7 +2238,7 @@ def TOP_LEVEL_exeOpenModelicaMultiobs( x_values_matrix , KeepCalculationFolders try: #to handle situations in which there is only one CL (boundary conditions) file var_int = numpy.ravel(LColumns[:,etat]) - except : + except Exception: var_int = numpy.ravel(LColumns) dict_inputs = dict(zip(LVariablesToChange, var_int)) @@ -2259,7 +2255,7 @@ def TOP_LEVEL_exeOpenModelicaMultiobs( x_values_matrix , KeepCalculationFolders try: reader = Reader(results_file_name,'dymola') #dymola even if it is OpenModelica - except: + except Exception: raise ValueError("Simulation cannot be performed: reduce the number of parameters to calibrate and/or the range in which their optimal value should be found (or modify and simplify the model to make it easier to simulate)" ) y_whole = [reader.values(y_name) for y_name in OutputVariables] @@ -2360,7 +2356,7 @@ def TOPLEVEL_exedymosimMultiobs_simple( x_values_matrix , VariablesToCalibrate=N debug_model_file = os.path.join(ref_simudir,'log_debug_model.txt') try: shutil.copy(dslog_file,debug_model_file) - except: + except Exception: pass if auto_simul.success_code == 2 : @@ -2471,7 +2467,7 @@ def TOP_LEVEL_exedymosimMultiobs( x_values_matrix , VariablesToCalibrate=None, O debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt') try: shutil.copy(dslog_file,debug_model_file) - except: + except Exception: pass if auto_simul.success_code == 2 : @@ -2539,11 +2535,11 @@ def TOP_LEVEL_exedymosimMultiobs( x_values_matrix , VariablesToCalibrate=None, O debug_model_file = os.path.join(simudir, os.pardir, 'log_debug_model.txt') try: shutil.copy(dslog_file,debug_model_file) - except: + except Exception: pass try: reader = Reader(os.path.join(simudir, 'RES','1.0.mat'),'dymola') - except: + except Exception: raise ValueError("Simulation cannot be performed: reduce the number of parameters to calibrate and/or the range in which their optimal value should be found (or modify and simplify the model to make it easier to simulate)" ) @@ -2632,7 +2628,7 @@ def run_OM_model(xml_file, #req arg, file path to _init.xml file try: os.remove(log_file) - except: + except Exception: pass try: f=open(log_file,'a') @@ -2646,7 +2642,7 @@ def run_OM_model(xml_file, #req arg, file path to _init.xml file else: try: proc.communicate(timeout = TimeoutModelExecution) - except: + except Exception: raise ValueError("Timeout for simulation reached, please increase it in order to be able to simulate your model and/or check if your model is correct (use TimeoutModelExecution option in configuration.py file)") # proc.communicate() @@ -2694,7 +2690,7 @@ def readObsnamesfile(__filenames=None): try: df = pandas.read_csv(__filename, sep = ";", header =0) #so that repeated names are not modified obsnames_infile = df.loc[0, :].values.tolist()[2:] - except: + except Exception: df_tmp = pandas.read_csv(__filename, sep = ";", header =1) with open(__filename, 'r') as infile: readie=csv.reader(infile, delimiter=';') @@ -2711,7 +2707,7 @@ def readObsnamesfile(__filenames=None): try: os.remove('tmp_obs_file.csv') - except: + except Exception: pass return obsnames_infile diff --git a/src/daComposant/daNumerics/pst4mod/modelica_libraries/around_simulation.py b/src/daComposant/daNumerics/pst4mod/modelica_libraries/around_simulation.py index 8b605e4..6f33504 100644 --- a/src/daComposant/daNumerics/pst4mod/modelica_libraries/around_simulation.py +++ b/src/daComposant/daNumerics/pst4mod/modelica_libraries/around_simulation.py @@ -94,7 +94,7 @@ class Around_Simulation(object): self.__dict_iter_var = {} self.__list_complex_iter_var = [] self.__dict_complex_iter_var = {} - self.__dict_iter_var_from_ini = {} # pour récupérer la valeur du fichier ini.txt + self.__dict_iter_var_from_ini = {} # pour récupérer la valeur du fichier ini.txt self.computed_variables = [] #List to store nale of variable computed by user defined function (if any) @@ -203,7 +203,7 @@ class Around_Simulation(object): return self.__dict_iter_var @property - def dict_iter_var_from_ini(self): # pour récupérer la valeur du fichier ini.txt + def dict_iter_var_from_ini(self): # pour récupérer la valeur du fichier ini.txt """ Dictionary which associates the name of an iteration variable to its value in ini.txt file """ @@ -532,7 +532,7 @@ class Around_Simulation(object): out.append(words[0].split(' ',1)[1]) #OpenModelica : Errors lines copied from the compilation log contains the 'unit' as first word (Integer, Real...) else: #Dymola, default case out.append(words[0]) - self.__dict_iter_var_from_ini[words[0]] =[word.strip() for word in line.split('=')][1][:-1] # valeur du ini.txt brute_init + self.__dict_iter_var_from_ini[words[0]] = [word.strip() for word in line.split('=')][1][:-1] # valeur du ini.txt brute_init #Detection of the iteration variables that are given the values of other variables @@ -628,11 +628,11 @@ class Around_Simulation(object): except KeyError: #The variable is not found try: - node_pos = re.search(r"\[\d+\]",res_name_var).span() #Position of the node number if any (in case the error is due to a remesh of the model) + node_pos = re.search(r"\[\d+\]",res_name_var).span() # Position of the node number if any (in case the error is due to a remesh of the model) except AttributeError: #No node number found try: - node_pos = re.search(r"\[\]",res_name_var).span() #Position of empty brakets (all array values are due) + node_pos = re.search(r"\[\]",res_name_var).span() # Position of empty brakets (all array values are due) except AttributeError : pass else: -- 2.30.2