Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F73600817
trained_model_pivoted_rational_polymatch.py
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Tue, Jul 23, 01:28
Size
15 KB
Mime Type
text/x-python
Expires
Thu, Jul 25, 01:28 (2 d)
Engine
blob
Format
Raw Data
Handle
19232123
Attached To
R6746 RationalROMPy
trained_model_pivoted_rational_polymatch.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
copy
import
deepcopy
as
copy
from
.trained_model_pivoted_rational_nomatch
import
(
TrainedModelPivotedRationalNoMatch
)
from
rrompy.utilities.base.types
import
(
Tuple
,
Np1D
,
Np2D
,
List
,
ListAny
,
paramVal
,
paramList
,
sampList
,
HFEng
)
from
rrompy.utilities.base
import
verbosityManager
as
vbMng
,
freepar
as
fp
from
rrompy.utilities.numerical
import
dot
from
rrompy.utilities.numerical.point_matching
import
polynomialMatching
from
rrompy.utilities.numerical.degree
import
reduceDegreeN
from
rrompy.utilities.poly_fitting.polynomial
import
(
polybases
as
ppb
,
PolynomialInterpolator
as
PI
,
PolynomialInterpolatorNodal
as
PIN
)
from
rrompy.utilities.poly_fitting.radial_basis
import
(
polybases
as
rbpb
,
RadialBasisInterpolator
as
RBI
)
from
rrompy.utilities.poly_fitting.heaviside
import
rational2heaviside
from
rrompy.utilities.poly_fitting.nearest_neighbor
import
(
NearestNeighborInterpolator
as
NNI
)
from
rrompy.utilities.poly_fitting.piecewise_linear
import
(
PiecewiseLinearInterpolator
as
PLI
)
from
rrompy.utilities.exception_manager
import
RROMPyException
,
RROMPyAssert
from
rrompy.sampling
import
sampleList
__all__
=
[
'TrainedModelPivotedRationalPolyMatch'
]
class
TrainedModelPivotedRationalPolyMatch
(
TrainedModelPivotedRationalNoMatch
):
"""
ROM approximant evaluation for pivoted approximants based on interpolation
of rational approximants (with polynomial matching).
Attributes:
Data: dictionary with all that can be pickled.
"""
def
centerNormalizeMarginal
(
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.mu0Marginal.
Returns:
Normalized parameter.
"""
mu
=
self
.
checkParameterListMarginal
(
mu
)
if
mu0
is
None
:
mu0
=
self
.
data
.
mu0
(
self
.
data
.
directionMarginal
)
return
(
self
.
mapParameterList
(
mu
,
idx
=
self
.
data
.
directionMarginal
)
-
self
.
mapParameterList
(
mu0
,
idx
=
self
.
data
.
directionMarginal
)
)
/
[
self
.
data
.
scaleFactor
[
x
]
for
x
in
self
.
data
.
directionMarginal
]
def
setupMarginalInterp
(
self
,
approx
,
interpPars
:
ListAny
,
extraPar
=
None
):
vbMng
(
self
,
"INIT"
,
"Starting computation of marginal interpolator."
,
12
)
musMCN
=
self
.
centerNormalizeMarginal
(
self
.
data
.
musMarginal
)
nM
,
pbM
=
len
(
musMCN
),
approx
.
polybasisMarginal
if
pbM
in
ppb
+
rbpb
:
if
extraPar
:
approx
.
_setMMarginalAuto
()
_MMarginalEff
=
approx
.
paramsMarginal
[
"MMarginal"
]
if
pbM
in
ppb
:
p
=
PI
()
elif
pbM
in
rbpb
:
p
=
RBI
()
else
:
# if pbM in sparsekinds + ["NEARESTNEIGHBOR"]:
if
pbM
==
"NEARESTNEIGHBOR"
:
p
=
NNI
()
else
:
# if pbM in sparsekinds:
pllims
=
[[
-
1.
]
*
self
.
data
.
nparMarginal
,
[
1.
]
*
self
.
data
.
nparMarginal
]
p
=
PLI
()
if
pbM
in
ppb
+
rbpb
:
if
not
extraPar
:
approx
.
paramsMarginal
[
"MMarginal"
]
=
reduceDegreeN
(
_MMarginalEff
,
len
(
musMCN
),
self
.
data
.
nparMarginal
,
approx
.
paramsMarginal
[
"polydegreetypeMarginal"
])
MMEff
=
approx
.
paramsMarginal
[
"MMarginal"
]
while
MMEff
>=
0
:
wellCond
,
msg
=
p
.
setupByInterpolation
(
musMCN
,
np
.
eye
(
nM
),
MMEff
,
*
interpPars
)
vbMng
(
self
,
"MAIN"
,
msg
,
30
)
if
wellCond
:
break
vbMng
(
self
,
"MAIN"
,
(
"Polyfit is poorly conditioned. Reducing "
"MMarginal by 1."
),
35
)
MMEff
-=
1
if
MMEff
<
0
:
raise
RROMPyException
((
"Instability in computation of "
"interpolant. Aborting."
))
if
(
pbM
in
rbpb
and
len
(
interpPars
)
>
4
and
"optimizeScalingBounds"
in
interpPars
[
4
]
.
keys
()):
interpPars
[
4
][
"optimizeScalingBounds"
]
=
[
-
1.
,
-
1.
]
elif
pbM
==
"NEARESTNEIGHBOR"
:
p
.
setupByInterpolation
(
musMCN
,
np
.
eye
(
nM
),
*
interpPars
)
else
:
# if pbM in sparsekinds:
p
.
setupByInterpolation
(
musMCN
,
np
.
eye
(
nM
),
pllims
,
extraPar
,
*
interpPars
)
self
.
data
.
marginalInterp
=
p
if
pbM
in
ppb
+
rbpb
:
approx
.
paramsMarginal
[
"MMarginal"
]
=
_MMarginalEff
vbMng
(
self
,
"DEL"
,
"Done computing marginal interpolator."
,
12
)
def
updateEffectiveSamples
(
self
,
exclude
:
List
[
int
],
*
args
,
**
kwargs
):
super
()
.
updateEffectiveSamples
(
exclude
)
self
.
initializeFromRational
(
*
args
,
**
kwargs
)
def
initializeFromRational
(
self
,
matchingWeight
:
float
,
matchingKind
:
str
,
HFEngine
:
HFEng
,
is_state
:
bool
):
"""Initialize rational representation."""
RROMPyAssert
(
self
.
data
.
nparPivot
,
1
,
"Number of pivot parameters"
)
Qs
=
[
Q
.
coeffs
for
Q
in
self
.
data
.
Qs
]
Ps
=
[
P
.
coeffs
for
P
in
self
.
data
.
Ps
]
N
=
len
(
Qs
)
degQ
=
np
.
max
([
Q
.
shape
[
0
]
for
Q
in
Qs
])
degP
=
np
.
max
([
P
.
shape
[
0
]
for
P
in
Ps
])
for
j
in
range
(
N
):
if
Qs
[
j
]
.
shape
[
0
]
<
degQ
:
Qs
[
j
]
=
np
.
pad
(
Qs
[
j
],
(
0
,
degQ
-
Qs
[
j
]
.
shape
[
0
]),
"constant"
)
if
Ps
[
j
]
.
shape
[
0
]
<
degP
:
Ps
[
j
]
=
np
.
pad
(
Ps
[
j
],
[(
0
,
degP
-
Ps
[
j
]
.
shape
[
0
]),
(
0
,
0
)],
"constant"
)
Qs
,
Ps
=
polynomialMatching
(
Qs
,
Ps
,
self
.
data
.
musMarginal
.
data
,
matchingWeight
,
self
.
data
.
Psupp
,
self
.
data
.
projMat
,
HFEngine
,
is_state
,
0
,
matchingKind
)
for
j
in
range
(
N
):
if
isinstance
(
self
.
data
.
Qs
[
j
],
PIN
):
q
=
PI
()
q
.
npar
=
self
.
data
.
Qs
[
j
]
.
npar
q
.
polybasis
=
self
.
data
.
Qs
[
j
]
.
polybasis
self
.
data
.
Qs
[
j
]
=
q
if
isinstance
(
self
.
data
.
Ps
[
j
],
PIN
):
p
=
PI
()
p
.
npar
=
self
.
data
.
Ps
[
j
]
.
npar
p
.
polybasis
=
self
.
data
.
Ps
[
j
]
.
polybasis
self
.
data
.
Ps
[
j
]
=
p
self
.
data
.
Qs
[
j
]
.
coeffs
,
self
.
data
.
Ps
[
j
]
.
coeffs
=
Qs
[
j
],
Ps
[
j
]
def
getApproxReduced
(
self
,
mu
:
paramList
=
[])
->
sampList
:
"""
Evaluate reduced representation of approximant at arbitrary parameter.
Args:
mu: Target parameter.
"""
RROMPyAssert
(
self
.
data
.
nparPivot
,
1
,
"Number of pivot parameters"
)
mu
=
self
.
checkParameterList
(
mu
)
if
(
not
hasattr
(
self
,
"lastSolvedApproxReduced"
)
or
self
.
lastSolvedApproxReduced
!=
mu
):
vbMng
(
self
,
"INIT"
,
"Evaluating approximant at mu = {}."
.
format
(
mu
),
12
)
muP
=
mu
(
self
.
data
.
directionPivot
)
muMC
=
self
.
centerNormalizeMarginal
(
mu
(
self
.
data
.
directionMarginal
))
mIvals
=
self
.
data
.
marginalInterp
(
muMC
)
QV
=
self
.
getQVal
(
muP
,
mIvals
=
mIvals
)
QVzero
=
np
.
where
(
QV
==
0.
)[
0
]
if
len
(
QVzero
)
>
0
:
QV
[
QVzero
]
=
np
.
finfo
(
np
.
complex
)
.
eps
/
(
1.
+
self
.
data
.
Qs
[
0
]
.
deg
[
0
])
self
.
uApproxReduced
=
self
.
getPVal
(
muP
,
mIvals
=
mIvals
)
/
QV
vbMng
(
self
,
"DEL"
,
"Done evaluating approximant."
,
12
)
self
.
lastSolvedApproxReduced
=
mu
return
self
.
uApproxReduced
def
interpolateMarginalP
(
self
,
mu
:
paramList
=
[],
mIvals
:
Np2D
=
None
)
->
ListAny
:
"""Obtain interpolated approximant numerator."""
mu
=
self
.
checkParameterListMarginal
(
mu
)
vbMng
(
self
,
"INIT"
,
"Interpolating marginal P at mu = {}."
.
format
(
mu
),
95
)
if
self
.
data
.
_collapsed
:
outShape
=
self
.
data
.
Ps
[
0
]
.
coeffs
.
shape
else
:
outShape
=
(
self
.
data
.
Ps
[
0
]
.
coeffs
.
shape
[
0
],
self
.
data
.
projMat
.
shape
[
1
])
intMP
=
np
.
zeros
((
len
(
mu
),)
+
outShape
,
dtype
=
np
.
complex
)
if
mIvals
is
None
:
muC
=
self
.
centerNormalizeMarginal
(
mu
)
mIvals
=
self
.
data
.
marginalInterp
(
muC
)
for
P
,
Psupp
,
mI
in
zip
(
self
.
data
.
Ps
,
self
.
data
.
Psupp
,
mIvals
):
iL
,
iR
=
Psupp
,
Psupp
+
P
.
shape
[
0
]
for
j
,
m
in
enumerate
(
mI
):
intMP
[
j
,
:,
iL
:
iR
]
+=
m
*
P
.
coeffs
vbMng
(
self
,
"DEL"
,
"Done interpolating marginal P."
,
95
)
return
intMP
def
interpolateMarginalQ
(
self
,
mu
:
paramList
=
[],
mIvals
:
Np2D
=
None
)
->
ListAny
:
"""Obtain interpolated approximant denominator."""
mu
=
self
.
checkParameterListMarginal
(
mu
)
vbMng
(
self
,
"INIT"
,
"Interpolating marginal Q at mu = {}."
.
format
(
mu
),
95
)
intMQ
=
np
.
zeros
((
len
(
mu
),)
+
self
.
data
.
Qs
[
0
]
.
coeffs
.
shape
,
dtype
=
np
.
complex
)
if
mIvals
is
None
:
muC
=
self
.
centerNormalizeMarginal
(
mu
)
mIvals
=
self
.
data
.
marginalInterp
(
muC
)
for
Q
,
mI
in
zip
(
self
.
data
.
Qs
,
mIvals
):
for
j
,
m
in
enumerate
(
mI
):
intMQ
[
j
]
+=
m
*
Q
.
coeffs
vbMng
(
self
,
"DEL"
,
"Done interpolating marginal Q."
,
95
)
return
intMQ
def
getPVal
(
self
,
mu
:
paramList
=
[],
mIvals
:
List
[
int
]
=
None
)
->
sampList
:
"""
Evaluate rational numerator at arbitrary parameter.
Args:
mu: Target parameter.
"""
RROMPyAssert
(
self
.
data
.
nparPivot
,
1
,
"Number of pivot parameters"
)
if
mIvals
is
None
:
mu
=
self
.
checkParameterList
(
mu
)
muP
=
self
.
centerNormalizePivot
(
mu
(
self
.
data
.
directionPivot
))
muM
=
self
.
centerNormalizeMarginal
(
mu
(
self
.
data
.
directionMarginal
))
mIvals
=
self
.
data
.
marginalInterp
(
muM
)
else
:
mu
=
self
.
checkParameterListPivot
(
mu
)
muP
=
self
.
centerNormalizePivot
(
mu
)
if
self
.
data
.
_collapsed
:
outShape
=
self
.
data
.
Ps
[
0
]
.
shape
else
:
outShape
=
(
self
.
data
.
projMat
.
shape
[
1
],)
p
=
np
.
zeros
(
outShape
+
(
len
(
mu
),),
dtype
=
np
.
complex
)
for
P
,
Psupp
,
mI
in
zip
(
self
.
data
.
Ps
,
self
.
data
.
Psupp
,
mIvals
):
iL
,
iR
=
Psupp
,
Psupp
+
P
.
shape
[
0
]
for
j
,
m
in
enumerate
(
mI
):
p
[
iL
:
iR
,
j
]
+=
m
*
P
(
muP
[
j
])[:,
0
]
return
sampleList
(
p
)
def
getQVal
(
self
,
mu
:
paramList
=
[],
der
:
List
[
int
]
=
None
,
scl
:
Np1D
=
None
,
mIvals
:
List
[
int
]
=
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
mIvals
is
None
:
mu
=
self
.
checkParameterList
(
mu
)
muP
=
self
.
centerNormalizePivot
(
mu
(
self
.
data
.
directionPivot
))
muM
=
self
.
centerNormalizeMarginal
(
mu
(
self
.
data
.
directionMarginal
))
mIvals
=
self
.
data
.
marginalInterp
(
muM
)
else
:
mu
=
self
.
checkParameterListPivot
(
mu
)
muP
=
self
.
centerNormalizePivot
(
mu
)
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
]]
q
=
np
.
zeros
(
len
(
mu
),
dtype
=
np
.
complex
)
for
Q
,
mI
in
zip
(
self
.
data
.
Qs
,
mIvals
):
for
j
,
m
in
enumerate
(
mI
):
q
[
j
]
=
q
[
j
]
+
m
*
Q
(
muP
[
j
],
derP
,
sclP
)
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"
)
RROMPyAssert
(
self
.
data
.
npar
,
len
(
marginalVals
),
"Number of 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
(
0
,
rDim
)
mMarg
=
[
mVals
[
j
]
for
j
in
range
(
len
(
mVals
))
if
j
!=
rDim
]
Q
=
copy
(
self
.
data
.
Qs
[
0
])
Q
.
coeffs
=
self
.
interpolateMarginalQ
(
mMarg
)[
0
]
roots
=
self
.
data
.
scaleFactor
[
rDim
]
*
Q
.
roots
()
return
self
.
mapParameterList
(
self
.
mapParameterList
(
self
.
data
.
mu0
(
rDim
),
idx
=
[
rDim
])(
0
,
0
)
+
roots
,
"B"
,
[
rDim
])(
0
)
def
getResidues
(
self
,
marginalVals
:
ListAny
=
[
fp
])
->
Tuple
[
paramList
,
Np2D
]:
"""
Obtain approximant residues.
Returns:
Numpy matrix with residues as columns.
"""
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
(
0
,
rDim
)
mMarg
=
[
mVals
[
j
]
for
j
in
range
(
len
(
mVals
))
if
j
!=
rDim
]
P
,
Q
=
copy
(
self
.
data
.
Ps
[
0
]),
copy
(
self
.
data
.
Qs
[
0
])
P
.
coeffs
=
self
.
interpolateMarginalP
(
mMarg
)[
0
]
Q
.
coeffs
=
self
.
interpolateMarginalQ
(
mMarg
)[
0
]
res
,
pls
,
_
=
rational2heaviside
(
P
,
Q
)
if
len
(
pls
)
==
0
:
return
pls
,
np
.
empty
((
0
,
0
),
dtype
=
np
.
complex
)
res
=
(
res
[:
len
(
pls
),
:]
.
T
*
self
.
data
.
scaleFactor
[
self
.
data
.
directionPivot
[
0
]])
pls
=
self
.
mapParameterList
(
self
.
mapParameterList
(
self
.
data
.
mu0
(
rDim
),
idx
=
[
rDim
])(
0
,
0
)
+
self
.
data
.
scaleFactor
[
rDim
]
*
pls
,
"B"
,
[
rDim
])(
0
)
if
not
self
.
data
.
_collapsed
:
res
=
dot
(
self
.
data
.
projMat
,
res
)
return
pls
,
res
.
T
Event Timeline
Log In to Comment