Page MenuHomec4science
Files F61465890

trained_model_pivoted_rational_nomatch.py
No OneTemporary
Actions

File Metadata

Created
Mon, May 6, 20:55

trained_model_pivoted_rational_nomatch.py
View Options

# Copyright (C) 2018-2020 by the RROMPy authors
#
# This file is part of RROMPy.
#
# RROMPy is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# RROMPy is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with RROMPy. If not, see <http://www.gnu.org/licenses/>.
#
import numpy as np
from rrompy.reduction_methods.standard.trained_model.trained_model_rational \
import TrainedModelRational
from rrompy.utilities.base.types import (Np1D, Np2D, Tuple, List, ListAny,
paramVal, paramList, sampList)
from rrompy.utilities.base import verbosityManager as vbMng, freepar as fp
from rrompy.utilities.numerical import dot
from rrompy.utilities.numerical.compress_matrix import compressMatrix
from rrompy.utilities.poly_fitting.heaviside import rational2heaviside
from rrompy.utilities.poly_fitting.nearest_neighbor import (
NearestNeighborInterpolator as NNI)
from rrompy.utilities.exception_manager import (RROMPyException, RROMPyAssert,
RROMPyWarning)
from rrompy.parameter import checkParameterList
from rrompy.sampling import emptySampleList
__all__ = ['TrainedModelPivotedRationalNoMatch']
class TrainedModelPivotedRationalNoMatch(TrainedModelRational):
"""
ROM approximant evaluation for pivoted approximants based on interpolation
of rational approximants (without pole matching).
Attributes:
Data: dictionary with all that can be pickled.
"""
def checkParameterListPivot(self, mu:paramList,
check_if_single : bool = False) -> paramList:
return checkParameterList(mu, self.data.nparPivot, check_if_single)
def checkParameterListMarginal(self, mu:paramList,
check_if_single : bool = False) -> paramList:
return checkParameterList(mu, self.data.nparMarginal, check_if_single)
def compress(self, collapse : bool = False, tol : float = 0.,
returnRMat : bool = False, **compressMatrixkwargs):
if not collapse and tol <= 0.: return
RMat = self.data.projMat
if not collapse:
if hasattr(self.data, "_compressTol"):
RROMPyWarning(("Recompressing already compressed model is "
"ineffective. Aborting."))
return
self.data.projMat, RMat, _ = compressMatrix(RMat, tol,
**compressMatrixkwargs)
if hasattr(self.data, "Ps"):
for obj, suppj in zip(self.data.Ps, self.data.Psupp):
obj.postmultiplyTensorize(RMat.T[suppj : suppj + obj.shape[0]])
if hasattr(self, "_PsExcl"):
for obj, suppj in zip(self._PsExcl, self._PsuppExcl):
obj.postmultiplyTensorize(RMat.T[suppj : suppj + obj.shape[0]])
self._PsuppExcl = [0] * len(self._PsuppExcl)
self.data.Psupp = [0] * len(self.data.Psupp)
super(TrainedModelRational, self).compress(collapse, tol)
if returnRMat: return RMat
def centerNormalizePivot(self, mu : paramList = [],
mu0 : paramVal = None) -> paramList:
"""
Compute normalized parameter to be plugged into approximant.
Args:
mu: Parameter(s) 1.
mu0: Parameter(s) 2. If None, set to self.data.mu0Pivot.
Returns:
Normalized parameter.
"""
mu = self.checkParameterListPivot(mu)
if mu0 is None:
mu0 = self.checkParameterListPivot(
self.data.mu0(0, self.data.directionPivot))
return (self.mapParameterList(mu, idx = self.data.directionPivot)
- self.mapParameterList(mu0, idx = self.data.directionPivot)
) / [self.data.scaleFactor[x] for x in self.data.directionPivot]
def setupMarginalInterp(self, interpPars:ListAny):
self.data.marginalInterp = NNI()
self.data.marginalInterp.setupByInterpolation(self.data.musMarginal,
np.arange(len(self.data.musMarginal)),
1, *interpPars)
def updateEffectiveSamples(self, exclude:List[int]):
if hasattr(self, "_idxExcl"):
for j, excl in enumerate(self._idxExcl):
self.data.musMarginal.insert(self._musMExcl[j], excl)
self.data.Ps.insert(excl, self._PsExcl[j])
self.data.Qs.insert(excl, self._QsExcl[j])
self.data.Psupp.insert(excl, self._PsuppExcl[j])
self._idxExcl, self._musMExcl = list(np.sort(exclude)), []
self._PsExcl, self._QsExcl, self._PsuppExcl = [], [], []
for excl in self._idxExcl[::-1]:
self._musMExcl = [self.data.musMarginal[excl]] + self._musMExcl
self.data.musMarginal.pop(excl)
self._PsExcl = [self.data.Ps.pop(excl)] + self._PsExcl
self._QsExcl = [self.data.Qs.pop(excl)] + self._QsExcl
self._PsuppExcl = [self.data.Psupp.pop(excl)] + self._PsuppExcl
def getPVal(self, mu : paramList = []) -> sampList:
"""
Evaluate rational numerator at arbitrary parameter.
Args:
mu: Target parameter.
"""
RROMPyAssert(self.data.nparPivot, 1, "Number of pivot parameters")
mu = self.checkParameterList(mu)
muP = self.centerNormalizePivot(mu(self.data.directionPivot))
muM = self.checkParameterListMarginal(mu(self.data.directionMarginal))
idxMUnique, idxMmap = np.unique(self.data.marginalInterp(muM),
return_inverse = True)
idxMUnique = np.array(idxMUnique, dtype = int)
p = emptySampleList()
vbMng(self, "INIT", "Evaluating numerator at mu = {}.".format(mu), 17)
for i, iM in enumerate(idxMUnique):
idx = np.where(idxMmap == i)[0]
Pval, supp = self.data.Ps[iM](muP[idx]), self.data.Psupp[iM]
if i == 0:
p.reset((self.data.projMat.shape[1], len(mu)),
dtype = Pval.dtype)
p.data[:] = 0.
p.data[supp : supp + len(Pval), idx] = Pval
vbMng(self, "DEL", "Done evaluating numerator.", 17)
return p
def getQVal(self, mu:Np1D, der : List[int] = None,
scl : Np1D = None) -> Np1D:
"""
Evaluate rational denominator at arbitrary parameter.
Args:
mu: Target parameter.
der(optional): Derivatives to take before evaluation.
"""
RROMPyAssert(self.data.nparPivot, 1, "Number of pivot parameters")
mu = self.checkParameterList(mu)
muP = self.centerNormalizePivot(mu(self.data.directionPivot))
muM = self.checkParameterListMarginal(mu(self.data.directionMarginal))
if der is None:
derP, derM = 0, [0]
else:
derP = der[self.data.directionPivot[0]]
derM = [der[x] for x in self.data.directionMarginal]
if np.any(np.array(derM) != 0):
raise RROMPyException(("Derivatives of Q with respect to marginal "
"parameters not allowed."))
sclP = 1 if scl is None else scl[self.data.directionPivot[0]]
idxMUnique, idxMmap = np.unique(self.data.marginalInterp(muM),
return_inverse = True)
idxMUnique = np.array(idxMUnique, dtype = int)
vbMng(self, "INIT", "Evaluating denominator at mu = {}.".format(mu),
17)
for i, iM in enumerate(idxMUnique):
idx = np.where(idxMmap == i)[0]
Qval = self.data.Qs[iM](muP[idx], derP, sclP)
if i == 0: q = np.empty(len(mu), dtype = Qval.dtype)
q[idx] = Qval
vbMng(self, "DEL", "Done evaluating denominator.", 17)
return q
def getPoles(self, marginalVals : ListAny = [fp]) -> paramList:
"""
Obtain approximant poles.
Returns:
Numpy complex vector of poles.
"""
RROMPyAssert(self.data.nparPivot, 1, "Number of pivot parameters")
mVals = list(marginalVals)
rDim = mVals.index(fp)
if rDim < len(mVals) - 1 and fp in mVals[rDim + 1 :]:
raise RROMPyException(("Exactly 1 'freepar' entry in "
"marginalVals must be provided."))
if rDim != self.data.directionPivot[0]:
raise RROMPyException(("'freepar' entry in marginalVals must "
"coincide with pivot direction."))
mVals[rDim] = self.data.mu0(rDim)[0]
muM = self.checkParameterListMarginal([mVals[j]
for j in range(len(mVals)) if j != rDim])
iM = int(self.data.marginalInterp(muM))
roots = self.data.scaleFactor[rDim] * self.data.Qs[iM].roots()
return self.mapParameterList(self.mapParameterList(self.data.mu0(rDim),
idx = [rDim])(0, 0)
+ roots, "B", [rDim])(0)
def getResidues(self, *args, **kwargs) -> Tuple[paramList, Np2D]:
"""
Obtain approximant residues.
Returns:
Numpy matrix with residues as columns.
"""
pls = self.getPoles(*args, **kwargs)
if len(pls) == 0:
return pls, np.empty((0, 0), dtype = self.data.Ps[0].coeffs.dtype)
if len(args) == 1:
mVals = args[0]
elif len(args) == 0:
mVals = [fp]
else:
mVals = kwargs["marginalVals"]
rDim = mVals.index(fp)
mVals[rDim] = self.data.mu0(rDim)[0]
muM = self.checkParameterListMarginal([mVals[j]
for j in range(len(mVals)) if j != rDim])
iM = int(self.data.marginalInterp(muM))
res = rational2heaviside(self.data.Ps[iM], self.data.Qs[iM])[0]
res = res[: len(pls), :].T
if not self.data._collapsed: res = dot(self.data.projMat, res).T
return pls, res

Event Timeline

c4science · Help