diff --git a/rrompy/reduction_methods/trained_model/trained_model_data.py b/rrompy/reduction_methods/trained_model/trained_model_data.py
index 766f236..00dd780 100644
--- a/rrompy/reduction_methods/trained_model/trained_model_data.py
+++ b/rrompy/reduction_methods/trained_model/trained_model_data.py
@@ -1,37 +1,36 @@
# Copyright (C) 2018 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 .
#
from copy import deepcopy as copy
from rrompy.utilities.base.types import Np2D, List, paramVal
from rrompy.utilities.exception_manager import RROMPyAssert
__all__ = ['TrainedModelData']
class TrainedModelData:
"""ROM approximant evaluation data (must be pickle-able)."""
- def __init__(self, name:str, mu0:paramVal, projMat:Np2D,
+ def __init__(self, mu0:paramVal, projMat:Np2D,
scaleFactor : List[float] = [1.],
rescalingExp : List[float] = [1.]):
self.npar = len(rescalingExp)
RROMPyAssert(mu0.shape[1], self.npar, "Number of parameters")
- self.name = name
self.mu0 = mu0
self.projMat = copy(projMat)
self.scaleFactor = scaleFactor
self.rescalingExp = rescalingExp
diff --git a/rrompy/reduction_methods/trained_model/trained_model_pivoted_data.py b/rrompy/reduction_methods/trained_model/trained_model_pivoted_data.py
index 79b8803..27b7a58 100644
--- a/rrompy/reduction_methods/trained_model/trained_model_pivoted_data.py
+++ b/rrompy/reduction_methods/trained_model/trained_model_pivoted_data.py
@@ -1,72 +1,72 @@
# Copyright (C) 2018 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 .
#
from .trained_model_data import TrainedModelData
from rrompy.utilities.base.types import Np2D, List, ListAny, paramVal
from rrompy.parameter import checkParameterList
__all__ = ['TrainedModelPivotedData']
class TrainedModelPivotedData(TrainedModelData):
"""ROM approximant evaluation data (must be pickle-able)."""
- def __init__(self, name:str, mu0:paramVal, projMat:Np2D,
+ def __init__(self, mu0:paramVal, projMat:Np2D,
scaleFactor : ListAny = [1.],
rescalingExp : List[float] = [1.],
directionPivot : ListAny = [0]):
- super().__init__(name, mu0, projMat, scaleFactor, rescalingExp)
+ super().__init__(mu0, projMat, scaleFactor, rescalingExp)
self.directionPivot = directionPivot
@property
def directionMarginal(self):
return tuple([x for x in range(self.npar) \
if x not in self.directionPivot])
@property
def mu0Pivot(self):
return checkParameterList(self.mu0(0, self.directionPivot),
self.nparPivot)[0]
@property
def mu0Marginal(self):
return checkParameterList(self.mu0(0, self.directionMarginal),
self.nparMarginal)[0]
@property
def nparPivot(self):
return len(self.directionPivot)
@property
def nparMarginal(self):
return self.npar - self.nparPivot
@property
def rescalingExpPivot(self):
return [self.rescalingExp[x] for x in self.directionPivot]
@property
def rescalingExpMarginal(self):
return [self.rescalingExp[x] for x in self.directionMarginal]
@property
def scaleFactorPivot(self):
return [self.scaleFactor[x] for x in self.directionPivot]
@property
def scaleFactorMarginal(self):
return [self.scaleFactor[x] for x in self.directionMarginal]
diff --git a/rrompy/reduction_methods/trained_model/trained_model_pivoted_rational.py b/rrompy/reduction_methods/trained_model/trained_model_pivoted_rational.py
index 98f8117..0e22cbe 100644
--- a/rrompy/reduction_methods/trained_model/trained_model_pivoted_rational.py
+++ b/rrompy/reduction_methods/trained_model/trained_model_pivoted_rational.py
@@ -1,107 +1,107 @@
# Copyright (C) 2018 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 .
#
from .trained_model_pivoted_general import TrainedModelPivotedGeneral
from .trained_model_rational import TrainedModelRational
from rrompy.utilities.base.types import Np1D, List, paramList, sampList, HFEng
from rrompy.utilities.base import verbosityManager as vbMng
from rrompy.utilities.poly_fitting.heaviside import rational2heaviside
from rrompy.utilities.exception_manager import RROMPyAssert
from rrompy.parameter import checkParameterList
from rrompy.sampling import sampleList
__all__ = ['TrainedModelPivotedRational']
class TrainedModelPivotedRational(TrainedModelPivotedGeneral,
TrainedModelRational):
"""
ROM approximant evaluation for pivoted Rational approximant (with pole
matching).
Attributes:
Data: dictionary with all that can be pickled.
"""
def initializeFromRational(self, HFEngine:HFEng,
matchingWeight : float = 1., POD : bool = True,
is_state : bool = True):
"""Initialize Heaviside representation."""
RROMPyAssert(self.data.nparPivot, 1, "Number of pivot parameters")
poles, coeffs = [], []
for Q, P in zip(self.data.Qs, self.data.Ps):
cfs, pls, basis = rational2heaviside(P, Q)
poles += [pls]
coeffs += [cfs]
self.initializeFromLists(poles, coeffs, basis, HFEngine,
matchingWeight, POD, is_state)
- self.data._temporary = False
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")
- if not self.data._temporary: return super().getPVal(mu)
+ if not hasattr(self.data, "_temporary"): return super().getPVal(mu)
mu = checkParameterList(mu, self.data.npar)[0]
muP = checkParameterList(mu.data[:, self.data.directionPivot],
self.data.nparPivot)[0]
muM = checkParameterList(mu.data[:, self.data.directionMarginal],
self.data.nparMarginal)[0]
muP = checkParameterList(muP, self.data.nparPivot)[0]
vbMng(self, "INIT", "Evaluating numerator at mu = {}.".format(muP), 17)
muPC = self.centerNormalizePivot(muP)
pP = [sampleList(P(muPC)) for P in self.data.Ps]
vbMng(self, "DEL", "Done evaluating numerator.", 17)
return self.interpolateMarginal(muM, pP)
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")
- if not self.data._temporary: return super().getQVal(mu, der, scl)
+ if not hasattr(self.data, "_temporary"):
+ return super().getQVal(mu, der, scl)
mu = checkParameterList(mu, self.data.npar)[0]
muP = checkParameterList(mu.data[:, self.data.directionPivot],
self.data.nparPivot)[0]
muM = checkParameterList(mu.data[:, self.data.directionMarginal],
self.data.nparMarginal)[0]
if der is None:
derP, derM = None, None
else:
derP = [der[x] for x in self.data.directionPivot]
derM = [der[x] for x in self.data.directionMarginal]
if scl is None:
sclP, sclM = None, None
else:
sclP = [scl[x] for x in self.data.directionPivot]
sclM = [scl[x] for x in self.data.directionMarginal]
muP = checkParameterList(muP, self.data.nparPivot)[0]
vbMng(self, "INIT", "Evaluating denominator at mu = {}.".format(muP),
17)
muPC = self.centerNormalizePivot(muP)
qP = [Q(muPC, derP, sclP).reshape(1, -1) for Q in self.data.Qs]
vbMng(self, "DEL", "Done evaluating denominator.", 17)
return self.interpolateMarginal(muM, qP, derM, sclM)
diff --git a/rrompy/reduction_methods/trained_model/trained_model_rational_mls.py b/rrompy/reduction_methods/trained_model/trained_model_rational_mls.py
index f012a55..180b232 100644
--- a/rrompy/reduction_methods/trained_model/trained_model_rational_mls.py
+++ b/rrompy/reduction_methods/trained_model/trained_model_rational_mls.py
@@ -1,174 +1,175 @@
# Copyright (C) 2018 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 .
#
import numpy as np
from .trained_model_rational import TrainedModelRational
from rrompy.utilities.base.types import Np1D, paramVal, paramList, sampList
from rrompy.utilities.base import verbosityManager as vbMng
from rrompy.utilities.poly_fitting.moving_least_squares import mlsweights
from rrompy.utilities.poly_fitting.polynomial import (
PolynomialInterpolator as PI)
from rrompy.utilities.numerical import customPInv, degreeTotalToFull
from rrompy.parameter import checkParameterList
from rrompy.sampling import emptySampleList
__all__ = ['TrainedModelRationalMLS']
class TrainedModelRationalMLS(TrainedModelRational):
"""
ROM approximant evaluation for rational moving least squares approximant.
Attributes:
Data: dictionary with all that can be pickled.
"""
def assembleReducedModel(self, mu:paramVal):
- if not hasattr(self, "lastSetupMu") or self.lastSetupMu != mu:
+ if not (hasattr(self.data, "lastSetupMu")
+ and self.data.lastSetupMu == mu):
vbMng(self, "INIT", "Assembling reduced model for mu = {}."\
.format(mu), 17)
vbMng(self, "INIT", "Starting computation of denominator.", 35)
muC = self.centerNormalize(mu)
muSC = self.centerNormalize(self.data.mus)
wQ = mlsweights(muC, muSC, self.data.radialBasisDen,
directionalWeights = self.data.radialWeightsDen,
nNearestNeighbor = self.data.nNearestNeighborDen)
if self.data.N > self.data.M:
PQVan = self.data.QVan
else:
PQVan = self.data.PVan
VQAdjW = PQVan.conj().T * wQ
VQAdjWVQ = VQAdjW.dot(PQVan)
interpPseudoInverse, info = customPInv(VQAdjWVQ, full = True,
rcond = self.data.interpRcond)
interpPseudoInverse = interpPseudoInverse.dot(VQAdjW).dot(
self.data.QBlocks)
if info[0] < interpPseudoInverse.shape[-1]:
q = np.zeros(interpPseudoInverse.shape[-1], dtype = np.complex)
q[0] = 1.
else:
halfGram = interpPseudoInverse[self.data.domQIdxs]
if self.data.POD:
Rstack = halfGram.reshape(-1, halfGram.shape[-1])
vbMng(self, "INIT",
"Solving svd for square root of gramian matrix.", 67)
_, s, eV = np.linalg.svd(Rstack, full_matrices = False)
condN = s[0] / s[-1]
q = eV[-1, :].T.conj()
vbMng(self, "MAIN",
("Solved svd problem of size {} x {} with condition "
"number {:.4e}.").format(*Rstack.shape, condN), 55)
vbMng(self, "DEL", "Done solving svd.", 67)
else:
RRstack = np.tensordot(self.trainedModel.gramian, halfGram,
1).reshape(-1, halfGram.shape[-1])
RLstack = halfGram.reshape(-1, halfGram.shape[-1])
gram = RLstack.T.conj().dot(RRstack)
vbMng(self, "INIT",
"Solving eigenvalue problem for gramian matrix.", 67)
ev, eV = np.linalg.eigh(gram)
condN = ev[-1] / ev[0]
q = eV[:, 0]
vbMng(self, "MAIN",
("Solved eigenvalue problem of size {} with "
"condition number {:.4e}.").format(gram.shape[0],
condN), 55)
vbMng(self, "DEL", "Done solving eigenvalue problem.", 67)
self.data.Q = PI()
self.data.Q.npar = self.npar
self.data.Q.polybasis = self.data.polybasis
if self.data.polydegreetype == "TOTAL":
self.data.Q.coeffs = degreeTotalToFull(
(self.data.N + 1,) * self.npar,
self.npar, q)
else:
self.data.Q.coeffs = q.reshape((self.data.N + 1,) * self.npar)
vbMng(self, "DEL", "Done computing denominator.", 35)
vbMng(self, "INIT", "Starting computation of numerator.", 35)
self.data.P = PI()
self.data.P.npar = self.npar
self.data.P.polybasis = self.data.polybasis
wP = mlsweights(muC, muSC, self.data.radialBasis,
directionalWeights = self.data.radialWeights,
nNearestNeighbor = self.data.nNearestNeighbor)
VAdjW = self.data.PVan.conj().T * wP
VAdjWV = VAdjW.dot(self.data.PVan)
interpPPseudoInverse = customPInv(VAdjWV, self.data.interpRcond)
Pcoeffs = np.tensordot(interpPPseudoInverse.dot(VAdjW),
self.data.QBlocks.dot(q), ([1], [1]))
if self.data.polydegreetype == "TOTAL":
self.data.P.coeffs = degreeTotalToFull(
(self.data.M + 1,) * self.npar
+ (self.data.QBlocks.shape[0],),
self.npar, Pcoeffs)
else:
self.data.P.coeffs = Pcoeffs.reshape(
(self.data.M + 1,) * self.npar
+ (self.data.QBlocks.shape[0],))
vbMng(self, "DEL", "Done computing numerator.", 35)
vbMng(self, "DEL", "Done assembling reduced model.", 17)
- self.lastSetupMu = mu
+ self.data.lastSetupMu = mu
def getApproxReduced(self, mu : paramList = []) -> sampList:
"""
Evaluate reduced representation of approximant at arbitrary parameter.
Args:
mu: Target parameter.
"""
mu = checkParameterList(mu, self.data.npar)[0]
if (not hasattr(self, "lastSolvedApproxReduced")
or self.lastSolvedApproxReduced != mu):
vbMng(self, "INIT",
"Evaluating approximant at mu = {}.".format(mu), 12)
self.uApproxReduced = emptySampleList()
for i in range(len(mu)):
self.assembleReducedModel(mu[i])
vbMng(self, "INIT",
"Solving reduced model for mu = {}.".format(mu[i]), 15)
uAppR = self.getPVal(mu[i]) / self.getQVal(mu[i])
if i == 0:
#self.data.P.shape[-1], len(mu)
self.uApproxReduced.reset((len(uAppR), len(mu)),
dtype = uAppR.dtype)
self.uApproxReduced[i] = uAppR
vbMng(self, "DEL", "Done solving reduced model.", 15)
vbMng(self, "DEL", "Done evaluating approximant.", 12)
self.lastSolvedApproxReduced = mu
return self.uApproxReduced
def getPoles(self, *args, mu : paramVal = None, **kwargs) -> Np1D:
"""
Obtain approximant poles.
Returns:
Numpy complex vector of poles.
"""
if mu is None: mu = self.data.mu0
self.assembleReducedModel(mu)
return super().getPoles(*args, **kwargs)
def getResidues(self, *args, mu : paramVal = None, **kwargs) -> Np1D:
"""
Obtain approximant residues.
Returns:
Numpy matrix with residues as columns.
"""
if mu is None: mu = self.data.mu0
self.assembleReducedModel(mu)
return super().getResidues(*args, **kwargs)
diff --git a/rrompy/reduction_methods/trained_model/trained_model_reduced_basis.py b/rrompy/reduction_methods/trained_model/trained_model_reduced_basis.py
index 33e8364..975987e 100644
--- a/rrompy/reduction_methods/trained_model/trained_model_reduced_basis.py
+++ b/rrompy/reduction_methods/trained_model/trained_model_reduced_basis.py
@@ -1,115 +1,116 @@
# Copyright (C) 2018 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 .
#
import numpy as np
from .trained_model import TrainedModel
from rrompy.utilities.base.types import (Np1D, ListAny, paramVal, paramList,
sampList)
from rrompy.utilities.base import verbosityManager as vbMng, freepar as fp
from rrompy.utilities.numerical import (eigvalsNonlinearDense,
marginalizePolyList)
from rrompy.utilities.expression import expressionEvaluator
from rrompy.utilities.exception_manager import RROMPyException, RROMPyWarning
from rrompy.parameter import checkParameter, checkParameterList
from rrompy.sampling import emptySampleList
__all__ = ['TrainedModelReducedBasis']
class TrainedModelReducedBasis(TrainedModel):
"""
ROM approximant evaluation for RB approximant.
Attributes:
Data: dictionary with all that can be pickled.
"""
def assembleReducedModel(self, mu:paramVal):
mu = checkParameter(mu, self.data.npar)
- if not hasattr(self, "lastSetupMu") or self.lastSetupMu != mu:
+ if not (hasattr(self.data, "lastSetupMu")
+ and self.data.lastSetupMu == mu):
vbMng(self, "INIT", "Assembling reduced model for mu = {}."\
.format(mu), 17)
muEff = mu ** self.data.rescalingExp
self.data.ARBmu, self.data.bRBmu = 0., 0.
for thA, ARB in zip(self.data.thAs, self.data.ARBs):
self.data.ARBmu = (expressionEvaluator(thA[0], muEff) * ARB
+ self.data.ARBmu)
for thb, bRB in zip(self.data.thbs, self.data.bRBs):
self.data.bRBmu = (expressionEvaluator(thb[0], muEff) * bRB
+ self.data.bRBmu)
vbMng(self, "DEL", "Done assembling reduced model.", 17)
- self.lastSetupMu = mu
+ self.data.lastSetupMu = mu
def getApproxReduced(self, mu : paramList = []) -> sampList:
"""
Evaluate reduced representation of approximant at arbitrary parameter.
Args:
mu: Target parameter.
"""
mu = checkParameterList(mu, self.data.npar)[0]
if (not hasattr(self, "lastSolvedApproxReduced")
or self.lastSolvedApproxReduced != mu):
vbMng(self, "INIT",
"Computing RB solution at mu = {}.".format(mu), 12)
self.uApproxReduced = emptySampleList()
for i in range(len(mu)):
self.assembleReducedModel(mu[i])
vbMng(self, "INIT",
"Solving reduced model for mu = {}.".format(mu[i]), 15)
uAppR = np.linalg.solve(self.data.ARBmu, self.data.bRBmu)
if i == 0:
#self.data.ARBs[0].shape[-1], len(mu)
self.uApproxReduced.reset((len(uAppR), len(mu)),
dtype = uAppR.dtype)
self.uApproxReduced[i] = uAppR
vbMng(self, "DEL", "Done solving reduced model.", 15)
vbMng(self, "DEL", "Done computing RB solution.", 12)
self.lastSolvedApproxReduced = mu
return self.uApproxReduced
def getPoles(self, marginalVals : ListAny = [fp], jSupp : int = 1,
**kwargs) -> Np1D:
"""
Obtain approximant poles.
Returns:
Numpy complex vector of poles.
"""
if not self.data.affinePoly:
RROMPyWarning(("Unable to compute approximate poles due "
"to parametric dependence (detected non-"
"polynomial). Change HFEngine.affinePoly to True "
"if necessary."))
return
if not hasattr(marginalVals, "__len__"): marginalVals = [marginalVals]
mVals = list(marginalVals)
try:
rDim = mVals.index(fp)
if rDim < len(mVals) - 1 and fp in mVals[rDim + 1 :]:
raise
except:
raise RROMPyException(("Exactly 1 'freepar' entry in "
"marginalVals must be provided."))
ARBs = self.data.ARBs
if self.data.npar > 1:
mVals[rDim] = self.data.mu0(rDim)
mVals = checkParameter(mVals).data.flatten()
mVals[rDim] = fp
ARBs = marginalizePolyList(ARBs, mVals, "auto")
ev = eigvalsNonlinearDense(ARBs, jSupp = jSupp, **kwargs)
return np.sort(np.power(ev, 1. / self.data.rescalingExp[rDim]))