Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F87986390
rational_interpolant.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
Wed, Oct 16, 04:12
Size
23 KB
Mime Type
text/x-python
Expires
Fri, Oct 18, 04:12 (2 d)
Engine
blob
Format
Raw Data
Handle
21695452
Attached To
R6746 RationalROMPy
rational_interpolant.py
View Options
# 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 <http://www.gnu.org/licenses/>.
#
from
copy
import
copy
import
numpy
as
np
from
scipy.special
import
factorial
as
fact
from
rrompy.reduction_methods.base
import
checkRobustTolerance
from
.generic_distributed_approximant
import
GenericDistributedApproximant
from
rrompy.utilities.poly_fitting
import
(
polybases
,
polyvander
,
polyfitname
,
customFit
)
from
rrompy.reduction_methods.trained_model
import
TrainedModelPade
as
tModel
from
rrompy.reduction_methods.trained_model
import
TrainedModelData
from
rrompy.utilities.base.types
import
Np1D
,
Np2D
,
HFEng
,
DictAny
,
Tuple
from
rrompy.utilities.base
import
verbosityDepth
,
purgeDict
from
rrompy.utilities.exception_manager
import
(
RROMPyException
,
modeAssert
,
RROMPyWarning
)
__all__
=
[
'RationalInterpolant'
]
class
RationalInterpolant
(
GenericDistributedApproximant
):
"""
ROM rational interpolant computation for parametric problems.
Args:
HFEngine: HF problem solver.
mu0(optional): Default parameter. Defaults to 0.
approxParameters(optional): Dictionary containing values for main
parameters of approximant. Recognized keys are:
- 'POD': whether to compute POD of snapshots; defaults to True;
- 'muBounds': list of bounds for parameter values; defaults to
[0, 1];
- 'S': total number of samples current approximant relies upon;
defaults to 2;
- 'sampler': sample point generator; defaults to uniform sampler on
muBounds;
- 'polybasis': type of polynomial basis for interpolation; allowed
values include 'MONOMIAL', 'CHEBYSHEV' and 'LEGENDRE'; defaults
to 'MONOMIAL';
- 'E': coefficient of interpolant to be minimized; defaults to
min(S, M + 1);
- 'M': degree of Pade' interpolant numerator; defaults to 0;
- 'N': degree of Pade' interpolant denominator; defaults to 0;
- 'interpRcond': tolerance for interpolation; defaults to None;
- 'robustTol': tolerance for robust Pade' denominator management;
defaults to 0.
Defaults to empty dict.
homogeneized(optional): Whether to homogeneize Dirichlet BCs. Defaults
to False.
verbosity(optional): Verbosity level. Defaults to 10.
Attributes:
HFEngine: HF problem solver.
mu0: Default parameter.
mus: Array of snapshot parameters.
ws: Array of snapshot weigths.
homogeneized: Whether to homogeneize Dirichlet BCs.
approxParameters: Dictionary containing values for main parameters of
approximant. Recognized keys are in parameterList.
parameterList: Recognized keys of approximant parameters:
- 'POD': whether to compute POD of snapshots;
- 'muBounds': list of bounds for parameter values;
- 'S': total number of samples current approximant relies upon;
- 'sampler': sample point generator;
- 'polybasis': type of polynomial basis for interpolation;
- 'E': coefficient of interpolant to be minimized;
- 'M': degree of Pade' interpolant numerator;
- 'N': degree of Pade' interpolant denominator;
- 'interpRcond': tolerance for interpolation via numpy.polyfit;
- 'robustTol': tolerance for robust Pade' denominator management.
extraApproxParameters: List of approxParameters keys in addition to
mother class's.
POD: Whether to compute POD of snapshots.
muBounds: list of bounds for parameter values.
S: Number of solution snapshots over which current approximant is
based upon.
sampler: Sample point generator.
polybasis: type of polynomial basis for interpolation.
M: Numerator degree of approximant.
N: Denominator degree of approximant.
interpRcond: Tolerance for interpolation via numpy.polyfit.
robustTol: Tolerance for robust Pade' denominator management.
samplingEngine: Sampling engine.
uHF: High fidelity solution with wavenumber lastSolvedHF as numpy
complex vector.
lastSolvedHF: Wavenumber corresponding to last computed high fidelity
solution.
Q: Numpy 1D vector containing complex coefficients of approximant
denominator.
P: Numpy 2D vector whose columns are FE dofs of coefficients of
approximant numerator.
uApp: Last evaluated approximant as numpy complex vector.
"""
def
__init__
(
self
,
HFEngine
:
HFEng
,
mu0
:
complex
=
0.
,
approxParameters
:
DictAny
=
{},
homogeneized
:
bool
=
False
,
verbosity
:
int
=
10
,
timestamp
:
bool
=
True
):
self
.
_preInit
()
self
.
_addParametersToList
([
"polybasis"
,
"E"
,
"M"
,
"N"
,
"interpRcond"
,
"robustTol"
])
super
()
.
__init__
(
HFEngine
=
HFEngine
,
mu0
=
mu0
,
approxParameters
=
approxParameters
,
homogeneized
=
homogeneized
,
verbosity
=
verbosity
,
timestamp
=
timestamp
)
self
.
_postInit
()
@property
def
approxParameters
(
self
):
"""
Value of approximant parameters.
"""
return
self
.
_approxParameters
@approxParameters.setter
def
approxParameters
(
self
,
approxParams
):
approxParameters
=
purgeDict
(
approxParams
,
self
.
parameterList
,
dictname
=
self
.
name
()
+
".approxParameters"
,
baselevel
=
1
)
approxParametersCopy
=
purgeDict
(
approxParameters
,
[
"polybasis"
,
"E"
,
"M"
,
"N"
,
"interpRcond"
,
"robustTol"
],
True
,
True
,
baselevel
=
1
)
if
hasattr
(
self
,
"_M"
)
and
self
.
M
is
not
None
:
Mold
=
self
.
M
self
.
_M
=
0
if
hasattr
(
self
,
"_N"
)
and
self
.
N
is
not
None
:
Nold
=
self
.
N
self
.
_N
=
0
if
hasattr
(
self
,
"_E"
)
and
self
.
E
is
not
None
:
self
.
_E
=
0
GenericDistributedApproximant
.
approxParameters
.
fset
(
self
,
approxParametersCopy
)
keyList
=
list
(
approxParameters
.
keys
())
if
"polybasis"
in
keyList
:
self
.
polybasis
=
approxParameters
[
"polybasis"
]
elif
not
hasattr
(
self
,
"_polybasis"
)
or
self
.
_polybasis
is
None
:
self
.
polybasis
=
"MONOMIAL"
if
"interpRcond"
in
keyList
:
self
.
interpRcond
=
approxParameters
[
"interpRcond"
]
elif
not
hasattr
(
self
,
"interpRcond"
)
or
self
.
interpRcond
is
None
:
self
.
interpRcond
=
None
if
"robustTol"
in
keyList
:
self
.
robustTol
=
approxParameters
[
"robustTol"
]
elif
not
hasattr
(
self
,
"_robustTol"
)
or
self
.
_robustTol
is
None
:
self
.
robustTol
=
0
if
"M"
in
keyList
:
self
.
M
=
approxParameters
[
"M"
]
elif
hasattr
(
self
,
"_M"
)
and
self
.
M
is
not
None
:
self
.
M
=
Mold
else
:
self
.
M
=
0
if
"N"
in
keyList
:
self
.
N
=
approxParameters
[
"N"
]
elif
hasattr
(
self
,
"_N"
)
and
self
.
N
is
not
None
:
self
.
N
=
Nold
else
:
self
.
N
=
0
if
"E"
in
keyList
:
self
.
E
=
approxParameters
[
"E"
]
else
:
self
.
E
=
min
(
self
.
S
-
1
,
self
.
M
+
1
)
@property
def
polybasis
(
self
):
"""Value of polybasis."""
return
self
.
_polybasis
@polybasis.setter
def
polybasis
(
self
,
polybasis
):
try
:
polybasis
=
polybasis
.
upper
()
.
strip
()
.
replace
(
" "
,
""
)
if
polybasis
not
in
polybases
:
raise
RROMPyException
(
"Prescribed polybasis not recognized."
)
self
.
_polybasis
=
polybasis
except
:
RROMPyWarning
((
"Prescribed polybasis not recognized. Overriding "
"to 'MONOMIAL'."
))
self
.
_sampleType
=
"MONOMIAL"
self
.
_approxParameters
[
"polybasis"
]
=
self
.
polybasis
@property
def
interpRcond
(
self
):
"""Value of interpRcond."""
return
self
.
_interpRcond
@interpRcond.setter
def
interpRcond
(
self
,
interpRcond
):
self
.
_interpRcond
=
interpRcond
self
.
_approxParameters
[
"interpRcond"
]
=
self
.
interpRcond
@property
def
M
(
self
):
"""Value of M. Its assignment may change S."""
return
self
.
_M
@M.setter
def
M
(
self
,
M
):
if
M
<
0
:
raise
RROMPyException
(
"M must be non-negative."
)
self
.
_M
=
M
self
.
_approxParameters
[
"M"
]
=
self
.
M
if
hasattr
(
self
,
"_S"
)
and
self
.
S
<
self
.
M
+
1
:
RROMPyWarning
(
"Prescribed S is too small. Updating S to M + 1."
)
self
.
S
=
self
.
M
+
1
@property
def
N
(
self
):
"""Value of N. Its assignment may change S."""
return
self
.
_N
@N.setter
def
N
(
self
,
N
):
if
N
<
0
:
raise
RROMPyException
(
"N must be non-negative."
)
self
.
_N
=
N
self
.
_approxParameters
[
"N"
]
=
self
.
N
if
hasattr
(
self
,
"_S"
)
and
self
.
S
<
self
.
N
+
1
:
RROMPyWarning
(
"Prescribed S is too small. Updating S to N + 1."
)
self
.
S
=
self
.
N
+
1
@property
def
E
(
self
):
"""Value of E. Its assignment may change S."""
return
self
.
_E
@E.setter
def
E
(
self
,
E
):
if
E
<
0
:
raise
RROMPyException
(
"E must be non-negative."
)
self
.
_E
=
E
self
.
_approxParameters
[
"E"
]
=
self
.
E
if
hasattr
(
self
,
"_S"
)
and
self
.
S
<
self
.
E
+
1
:
RROMPyWarning
(
"Prescribed S is too small. Updating S to E + 1."
)
self
.
S
=
self
.
E
+
1
@property
def
robustTol
(
self
):
"""Value of tolerance for robust Pade' denominator management."""
return
self
.
_robustTol
@robustTol.setter
def
robustTol
(
self
,
robustTol
):
if
robustTol
<
0.
:
RROMPyWarning
((
"Overriding prescribed negative robustness "
"tolerance to 0."
))
robustTol
=
0.
self
.
_robustTol
=
robustTol
self
.
_approxParameters
[
"robustTol"
]
=
self
.
robustTol
@property
def
S
(
self
):
"""Value of S."""
return
self
.
_S
@S.setter
def
S
(
self
,
S
):
if
S
<=
0
:
raise
RROMPyException
(
"S must be positive."
)
if
hasattr
(
self
,
"_S"
):
Sold
=
self
.
S
else
:
Sold
=
-
1
vals
,
label
=
[
0
]
*
3
,
{
0
:
"M"
,
1
:
"N"
,
2
:
"E"
}
if
hasattr
(
self
,
"_M"
)
and
self
.
_M
is
not
None
:
vals
[
0
]
=
self
.
M
if
hasattr
(
self
,
"_N"
)
and
self
.
_N
is
not
None
:
vals
[
1
]
=
self
.
N
if
hasattr
(
self
,
"_E"
)
and
self
.
_E
is
not
None
:
vals
[
2
]
=
self
.
E
idxmax
=
np
.
argmax
(
vals
)
if
vals
[
idxmax
]
+
1
>
S
:
RROMPyWarning
((
"Prescribed S is too small. Updating S to {} + "
"1."
)
.
format
(
label
[
idxmax
]))
self
.
S
=
vals
[
idxmax
]
+
1
else
:
self
.
_S
=
S
self
.
_approxParameters
[
"S"
]
=
self
.
S
if
Sold
!=
self
.
S
:
self
.
resetSamples
()
def
_setupDenominator
(
self
):
"""Compute Pade' denominator."""
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"INIT"
,
"Starting computation of denominator."
,
timestamp
=
self
.
timestamp
)
while
self
.
N
>
0
:
TE
=
polyvander
[
self
.
polybasis
](
self
.
radiusPade
(
self
.
mus
),
self
.
E
,
scl
=
1.
/
self
.
scaleFactor
)
TE
=
(
TE
.
T
*
self
.
ws
)
.
T
RHS
=
np
.
zeros
(
self
.
E
+
1
)
RHS
[
-
1
]
=
1.
fitOut
=
customFit
(
TE
.
T
,
RHS
,
full
=
True
,
rcond
=
self
.
interpRcond
)
if
self
.
verbosity
>=
5
:
condfit
=
fitOut
[
1
][
2
][
0
]
/
fitOut
[
1
][
2
][
-
1
]
verbosityDepth
(
"MAIN"
,
(
"Fitting {} samples with degree {} "
"through {}... Conditioning of LS "
"system: {:.4e}."
)
.
format
(
self
.
S
,
self
.
E
,
polyfitname
[
self
.
polybasis
],
condfit
),
timestamp
=
self
.
timestamp
)
if
fitOut
[
1
][
1
]
<
self
.
E
+
1
:
Enew
=
fitOut
[
1
][
1
]
-
1
Nnew
=
min
(
self
.
N
,
Enew
)
Mnew
=
min
(
self
.
M
,
Enew
)
if
Nnew
==
self
.
N
:
strN
=
""
else
:
strN
=
"N from {} to {} and "
.
format
(
self
.
N
,
Nnew
)
if
Mnew
==
self
.
M
:
strM
=
""
else
:
strM
=
"M from {} to {} and "
.
format
(
self
.
M
,
Mnew
)
RROMPyWarning
((
"Polyfit is poorly conditioned.
\n
Reducing {}{}"
"E from {} to {}."
)
.
format
(
strN
,
strM
,
self
.
E
,
Enew
))
newParams
=
{
"N"
:
Nnew
,
"M"
:
Mnew
,
"E"
:
Enew
}
self
.
approxParameters
=
newParams
continue
mus_un
,
idx_un
,
cnt_un
=
np
.
unique
(
self
.
mus
,
return_inverse
=
True
,
return_counts
=
True
)
TE
=
polyvander
[
self
.
polybasis
](
self
.
radiusPade
(
self
.
mus
),
self
.
N
,
scl
=
1.
/
self
.
scaleFactor
)
TE
=
(
TE
.
T
*
self
.
ws
)
.
T
if
len
(
mus_un
)
==
len
(
self
.
mus
):
Ghalf
=
(
TE
.
T
*
fitOut
[
0
])
.
T
else
:
pseudoInv
=
np
.
zeros
((
len
(
self
.
mus
),
len
(
self
.
mus
)),
dtype
=
np
.
complex
)
for
j
in
range
(
len
(
mus_un
)):
pseudoInv_loc
=
np
.
zeros
((
cnt_un
[
j
],
cnt_un
[
j
]),
dtype
=
np
.
complex
)
mask
=
np
.
arange
(
len
(
self
.
mus
))[
idx_un
==
j
]
for
der
in
range
(
cnt_un
[
j
]):
fitderj
=
fitOut
[
0
][
mask
[
der
]]
pseudoInv_loc
=
(
pseudoInv_loc
+
fitderj
*
np
.
diag
(
np
.
ones
(
1
+
der
),
k
=
der
-
cnt_un
[
j
]
+
1
))
I
=
np
.
ix_
(
mask
,
mask
)
pseudoInv
[
I
]
=
np
.
flipud
(
pseudoInv_loc
)
Ghalf
=
pseudoInv
.
dot
(
TE
)
if
self
.
POD
:
self
.
Ghalf
=
self
.
samplingEngine
.
RPOD
.
dot
(
Ghalf
)
ev
,
eV
=
self
.
findeveVGQR
()
else
:
self
.
Ghalf
=
self
.
samplingEngine
.
samples
.
dot
(
Ghalf
)
ev
,
eV
=
self
.
findeveVGExplicit
()
newParams
=
checkRobustTolerance
(
ev
,
self
.
E
,
self
.
robustTol
)
if
not
newParams
:
break
self
.
approxParameters
=
newParams
if
self
.
N
<=
0
:
self
.
_N
=
0
eV
=
np
.
ones
((
1
,
1
))
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"DEL"
,
"Done computing denominator."
,
timestamp
=
self
.
timestamp
)
return
eV
[:,
0
]
def
_setupNumerator
(
self
):
"""Compute Pade' numerator."""
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"INIT"
,
"Starting computation of numerator."
,
timestamp
=
self
.
timestamp
)
Qevaldiag
=
np
.
diag
(
self
.
trainedModel
.
getQVal
(
self
.
mus
))
verb
=
self
.
trainedModel
.
verbosity
self
.
trainedModel
.
verbosity
=
0
mus_un
,
idx_un
,
cnt_un
=
np
.
unique
(
self
.
mus
,
return_inverse
=
True
,
return_counts
=
True
)
for
j
in
range
(
len
(
mus_un
)):
if
cnt_un
[
j
]
>
1
:
Q_loc
=
np
.
zeros
((
cnt_un
[
j
],
cnt_un
[
j
]),
dtype
=
np
.
complex
)
for
der
in
range
(
1
,
cnt_un
[
j
]):
Qderj
=
(
self
.
trainedModel
.
getQVal
(
mus_un
[
j
],
der
,
scl
=
1.
/
self
.
scaleFactor
)
/
fact
(
der
))
Q_loc
=
Q_loc
+
Qderj
*
np
.
diag
(
np
.
ones
(
cnt_un
[
j
]
-
der
),
k
=
-
der
)
I
=
np
.
ix_
(
idx_un
==
j
,
idx_un
==
j
)
Qevaldiag
[
I
]
=
Qevaldiag
[
I
]
+
Q_loc
self
.
trainedModel
.
verbosity
=
verb
while
self
.
M
>=
0
:
fitVander
=
polyvander
[
self
.
polybasis
](
self
.
radiusPade
(
self
.
mus
),
self
.
M
,
scl
=
1.
/
self
.
scaleFactor
)
fitOut
=
customFit
(
fitVander
,
Qevaldiag
,
w
=
self
.
ws
,
full
=
True
,
rcond
=
self
.
interpRcond
)
if
self
.
verbosity
>=
5
:
condfit
=
fitOut
[
1
][
2
][
0
]
/
fitOut
[
1
][
2
][
-
1
]
verbosityDepth
(
"MAIN"
,
(
"Fitting {} samples with degree {} "
"through {}... Conditioning of LS "
"system: {:.4e}."
)
.
format
(
self
.
S
,
self
.
M
,
polyfitname
[
self
.
polybasis
],
condfit
),
timestamp
=
self
.
timestamp
)
if
fitOut
[
1
][
1
]
==
self
.
M
+
1
:
P
=
fitOut
[
0
]
.
T
break
RROMPyWarning
((
"Polyfit is poorly conditioned. Reducing M from {} "
"to {}. Exact snapshot interpolation not "
"guaranteed."
)
.
format
(
self
.
M
,
fitOut
[
1
][
1
]
-
1
))
self
.
M
=
fitOut
[
1
][
1
]
-
1
if
self
.
M
<=
0
:
raise
RROMPyException
((
"Instability in computation of numerator. "
"Aborting."
))
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"DEL"
,
"Done computing numerator."
,
timestamp
=
self
.
timestamp
)
return
np
.
atleast_2d
(
P
)
def
setupApprox
(
self
):
"""
Compute Pade' interpolant.
SVD-based robust eigenvalue management.
"""
if
self
.
checkComputedApprox
():
return
modeAssert
(
self
.
_mode
,
message
=
"Cannot setup approximant."
)
if
self
.
verbosity
>=
5
:
verbosityDepth
(
"INIT"
,
"Setting up {}."
.
format
(
self
.
name
()),
timestamp
=
self
.
timestamp
)
self
.
computeScaleFactor
()
self
.
computeSnapshots
()
if
self
.
trainedModel
is
None
:
self
.
trainedModel
=
tModel
()
self
.
trainedModel
.
verbosity
=
self
.
verbosity
self
.
trainedModel
.
timestamp
=
self
.
timestamp
data
=
TrainedModelData
(
self
.
trainedModel
.
name
(),
self
.
mu0
,
np
.
copy
(
self
.
samplingEngine
.
samples
),
self
.
HFEngine
.
rescalingExp
)
data
.
polytype
=
self
.
polybasis
data
.
scaleFactor
=
self
.
scaleFactor
data
.
mus
=
np
.
copy
(
self
.
mus
)
self
.
trainedModel
.
data
=
data
else
:
self
.
trainedModel
.
data
.
projMat
=
np
.
copy
(
self
.
samplingEngine
.
samples
)
if
self
.
N
>
0
:
Q
=
self
.
_setupDenominator
()
else
:
Q
=
np
.
ones
(
1
,
dtype
=
np
.
complex
)
self
.
trainedModel
.
data
.
Q
=
np
.
copy
(
Q
)
P
=
self
.
_setupNumerator
()
if
self
.
POD
:
P
=
self
.
samplingEngine
.
RPOD
.
dot
(
P
)
self
.
trainedModel
.
data
.
P
=
np
.
copy
(
P
)
self
.
trainedModel
.
data
.
approxParameters
=
copy
(
self
.
approxParameters
)
if
self
.
verbosity
>=
5
:
verbosityDepth
(
"DEL"
,
"Done setting up approximant."
,
timestamp
=
self
.
timestamp
)
def
findeveVGExplicit
(
self
,
verbOutput
:
int
=
5
)
->
Tuple
[
Np1D
,
Np2D
]:
"""
Compute explicitly eigenvalues and eigenvectors of Pade' denominator
matrix.
"""
modeAssert
(
self
.
_mode
,
message
=
"Cannot solve eigenvalue problem."
)
if
self
.
verbosity
>=
10
:
verbosityDepth
(
"INIT"
,
"Building gramian matrix."
,
timestamp
=
self
.
timestamp
)
self
.
G
=
self
.
HFEngine
.
innerProduct
(
self
.
Ghalf
,
self
.
Ghalf
)
if
self
.
verbosity
>=
10
:
verbosityDepth
(
"DEL"
,
"Done building gramian."
,
timestamp
=
self
.
timestamp
)
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"INIT"
,
(
"Solving eigenvalue problem for gramian "
"matrix."
),
timestamp
=
self
.
timestamp
)
ev
,
eV
=
np
.
linalg
.
eigh
(
self
.
G
)
if
self
.
verbosity
>=
verbOutput
:
try
:
condev
=
ev
[
-
1
]
/
ev
[
0
]
except
:
condev
=
np
.
inf
verbosityDepth
(
"MAIN"
,
(
"Solved eigenvalue problem of size {} "
"with condition number {:.4e}."
)
.
format
(
self
.
N
+
1
,
condev
),
timestamp
=
self
.
timestamp
)
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"DEL"
,
"Done solving eigenvalue problem."
,
timestamp
=
self
.
timestamp
)
return
ev
,
eV
def
findeveVGQR
(
self
,
verbOutput
:
int
=
5
)
->
Tuple
[
Np1D
,
Np2D
]:
"""
Compute eigenvalues and eigenvectors of Pade' denominator matrix
through SVD of R factor.
"""
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"INIT"
,
(
"Solving svd for square root of gramian "
"matrix."
),
timestamp
=
self
.
timestamp
)
_
,
s
,
eV
=
np
.
linalg
.
svd
(
self
.
Ghalf
,
full_matrices
=
False
)
ev
=
s
[::
-
1
]
eV
=
eV
[::
-
1
,
:]
.
T
.
conj
()
if
self
.
verbosity
>=
verbOutput
:
try
:
condev
=
s
[
0
]
/
s
[
-
1
]
except
:
condev
=
np
.
inf
verbosityDepth
(
"MAIN"
,
(
"Solved svd problem of size {} x {} with "
"condition number {:.4e}."
)
.
format
(
self
.
S
,
self
.
N
+
1
,
condev
),
timestamp
=
self
.
timestamp
)
if
self
.
verbosity
>=
7
:
verbosityDepth
(
"DEL"
,
"Done solving eigenvalue problem."
,
timestamp
=
self
.
timestamp
)
return
ev
,
eV
def
radiusPade
(
self
,
mu
:
Np1D
,
mu0
:
float
=
None
)
->
float
:
"""
Compute translated radius to be plugged into Pade' approximant.
Args:
mu: Parameter(s) 1.
mu0: Parameter(s) 2. If None, set to self.mu0.
Returns:
Translated radius to be plugged into Pade' approximant.
"""
return
self
.
trainedModel
.
radiusPade
(
mu
,
mu0
)
def
getResidues
(
self
)
->
Np1D
:
"""
Obtain approximant residues.
Returns:
Matrix with residues as columns.
"""
return
self
.
trainedModel
.
getResidues
()
Event Timeline
Log In to Comment