Page MenuHomec4science
Files F61702546

symmetric_disk.py
No OneTemporary
Actions

File Metadata

Created
Wed, May 8, 10:41

symmetric_disk.py
View Options

import numpy as np
from symmetric_disk_engine import SymmetricDiskEngine as engine
from rrompy.reduction_methods.standard import (RationalInterpolant as RI,
ReducedBasis as RB)
from rrompy.reduction_methods.greedy import (RationalInterpolantGreedy as RIG,
ReducedBasisGreedy as RBG)
from rrompy.parameter.parameter_sampling import QuadratureSampler as QS
ks = [10., 20.]
k0, n = np.mean(np.power(ks, 2.)) ** .5, 150
solver = engine(k0, n)
k = 12.
method = "RI"
method = "RB"
method = "RI_GREEDY"
method = "RB_GREEDY"
if method == "RI":
params = {'N':39, 'M':39, 'S':40, 'POD':True, 'polybasis':"CHEBYSHEV",
'sampler':QS(ks, "CHEBYSHEV", scalingExp = 2.)}
algo = RI
if method == "RB":
params = {'R':40, 'S':40, 'POD':True,
'sampler':QS(ks, "CHEBYSHEV", scalingExp = 2.)}
algo = RB
if method == "RI_GREEDY":
params = {'S':10, 'POD':True, 'polybasis':"LEGENDRE", 'greedyTol':1e-2,
'sampler':QS(ks, "UNIFORM", scalingExp = 2.),
'errorEstimatorKind':"DISCREPANCY"}
algo = RIG
if method == "RB_GREEDY":
params = {'S':10, 'POD':True, 'greedyTol':1e-2,
'sampler':QS(ks, "UNIFORM", scalingExp = 2.)}
algo = RBG
approx = algo(solver, mu0 = k0, approx_state = True,
approxParameters = params, verbosity = 20)
if len(method) == 2:
approx.setupApprox()
else:
approx.setupApprox(True)
approx.plotApprox(k, name = 'u_app')
approx.plotHF(k, name = 'u_HF')
approx.plotErr(k, name = 'err')
approx.plotRes(k, name = 'res')
normErr = approx.normErr(k)[0]
normSol = approx.normHF(k)[0]
normRes = approx.normRes(k)[0]
normRHS = approx.normRHS(k)[0]
print("SolNorm:\t{:.5e}\nErr:\t{:.5e}\nErrRel:\t{:.5e}".format(
normSol, normErr, normErr / normSol))
print("RHSNorm:\t{:.5e}\nRes:\t{:.5e}\nResRel:\t{:.5e}".format(
normRHS, normRes, normRes / normRHS))
if method[:2] == "RI":
poles, residues = approx.getResidues()
if method[:2] == "RB":
poles = approx.getPoles()
print("Poles:\n{}".format(poles))
if method[:2] == "RI":
for pol, res in zip(poles, residues.T):
solver.plot(res)
print("pole = {:.5e}".format(pol))

Event Timeline

c4science · Help