Page MenuHomec4science
Files F68548151

pod.py
No OneTemporary
Actions

File Metadata

Created
Thu, Jun 27, 22:07

pod.py
View Options

import numpy as np
from all_forcing_engine import AllForcingEngine
from rrompy.reduction_methods.standard import RationalInterpolant as RI
from rrompy.reduction_methods.standard import ReducedBasis as RB
from rrompy.parameter.parameter_sampling import QuadratureSampler as QS
verb = 100
sol = "single"
sol = "sweep"
algo = "RI"
#algo = "RB"
polyBasis = "LEGENDRE"
polyBasis = "CHEBYSHEV"
#polyBasis = "MONOMIAL"
ztar = 2.
z0s = [-3., 3.]
z0 = np.mean(z0s)
n = 30
solver = AllForcingEngine(mu0 = z0, n = n, degree_threshold = 8, verbosity = 0)
params = {'N':3, 'M':3, 'S':4, 'POD':True, 'polybasis':polyBasis,
'sampler':QS(z0s, "CHEBYSHEV")}
if algo == "RI":
approx = RI(solver, mu0 = z0, approxParameters = params, verbosity = verb)
else:
params.pop("N")
params.pop("M")
params.pop("polybasis")
approx = RB(solver, mu0 = z0, approxParameters = params, verbosity = verb)
approx.setupApprox()
if sol == "single":
approx.plotSamples(what = "REAL")
approx.plotApprox(ztar, what = "REAL", name = "uApp")
approx.plotHF(ztar, what = "REAL", name = "uHF")
approx.plotErr(ztar, what = "REAL", name = "err")
approx.plotRes(ztar, what = "REAL", name = "res")
appErr, solNorm = approx.normErr(ztar), approx.normHF(ztar)
resNorm, RHSNorm = approx.normRes(ztar), approx.normRHS(ztar)
print(('SolNorm:\t{}\nErr:\t{}\nErrRel:\t{}').format(solNorm, appErr,
np.divide(appErr, solNorm)))
print(('RHSNorm:\t{}\nRes:\t{}\nResRel:\t{}').format(RHSNorm, resNorm,
np.divide(resNorm, RHSNorm)))
poles = approx.getPoles()
print('Poles:', poles)
if sol == "sweep":
z0s = np.linspace(z0s[0], z0s[1], 100)
zl, zr = min(z0s), max(z0s)
approx.samplingEngine.verbosity = 0
approx.trainedModel.verbosity = 0
approx.verbosity = 0
from matplotlib import pyplot as plt
normRHS = approx.normRHS(z0s)
norm = approx.normHF(z0s)
normApp = approx.normApprox(z0s)
res = approx.normRes(z0s) / normRHS
err = approx.normErr(z0s) / norm
plt.figure()
plt.semilogy(z0s, norm)
plt.semilogy(z0s, normApp, '--')
plt.semilogy(np.real(approx.mus.data),
1.05*np.max(norm)*np.ones_like(approx.mus.data, dtype = float),
'rx')
plt.xlim([zl, zr])
plt.grid()
plt.show()
plt.close()
plt.figure()
plt.semilogy(z0s, res)
plt.xlim([zl, zr])
plt.grid()
plt.show()
plt.close()
plt.figure()
plt.semilogy(z0s, err)
# plt.semilogy(k0s, errApp)
plt.xlim([zl, zr])
plt.grid()
plt.show()
plt.close()

Event Timeline

c4science · Help