Page MenuHomec4science
Files F60225450

synthetic_pod.py
No OneTemporary
Actions

File Metadata

Created
Sun, Apr 28, 11:10

synthetic_pod.py
View Options

import numpy as np
from rrompy.hfengines.linear_problem.bidimensional import \
SyntheticBivariateEngine as SBE
from rrompy.reduction_methods.centered import RationalPade as RP
from rrompy.reduction_methods.distributed import RationalInterpolant as RI
from rrompy.reduction_methods.centered import RBCentered as RBC
from rrompy.reduction_methods.distributed import RBDistributed as RBD
from rrompy.parameter.parameter_sampling import (QuadratureSampler as QS,
QuadratureSamplerTotal as QST,
ManualSampler as MS,
RandomSampler as RS)
verb = 5
size = 1
show_sample = True
show_norm = True
clip = -1
#clip = .4
#clip = .6
Delta = 0
MN = 10
R = (MN + 2) * (MN + 1) // 2
S = [int(np.ceil(R ** .5))] * 2
PODTol = 1e-10
samples = "centered"
samples = "centered_fake"
samples = "distributed"
algo = "rational"
#algo = "RB"
sampling = "quadrature"
#sampling = "quadrature_total"
#sampling = "random"
if samples == "distributed":
radial = 0
# radial = "gaussian"
# radial = "thinplate"
# radial = "multiquadric"
rW0 = 10.
radialWeight = [r * rW0 for r in [5., 5.]]
if size == 1: # small
mu0 = [10. ** .5, 15. ** .5]
mutar = [12. ** .5, 14. ** .5]
murange = [[5. ** .5, 10. ** .5], [15 ** .5, 20 ** .5]]
if size == 2: # large
mu0 = [15. ** .5, 17.5 ** .5]
mutar = [18. ** .5, 22. ** .5]
murange = [[5. ** .5, 10. ** .5], [25 ** .5, 25 ** .5]]
if size == 3: # medium
mu0 = [17.5 ** .5, 15 ** .5]
mutar = [20. ** .5, 18. ** .5]
murange = [[10. ** .5, 10. ** .5], [25 ** .5, 20 ** .5]]
assert Delta <= 0
aEff = 1.#25
bEff = 1. - aEff
murangeEff = [[(aEff*murange[0][0]**2. + bEff*murange[1][0]**2.) ** .5,
aEff*murange[0][1] + bEff*murange[1][1]],
[(aEff*murange[1][0]**2. + bEff*murange[0][0]**2.) ** .5,
aEff*murange[1][1] + bEff*murange[0][1]]]
kappa = 20. ** .5
theta = - np.pi / 6.
n = 20
L = np.pi
solver = SBE(kappa = kappa, theta = theta, n = n, L = L,
mu0 = mu0, verbosity = verb)
rescaling = [lambda x: np.power(x, 2.)] * 2
rescalingInv = [lambda x: np.power(x, .5)] * 2
if algo == "rational":
params = {'N':MN, 'M':MN + Delta, 'S':S, 'POD':True}
if samples == "distributed":
params['polybasis'] = "CHEBYSHEV"
# params['polybasis'] = "LEGENDRE"
# params['polybasis'] = "MONOMIAL"
params['E'] = MN
params['radialBasis'] = radial
params['radialBasisWeights'] = radialWeight
method = RI
elif samples == "centered_fake":
params['polybasis'] = "MONOMIAL"
params['S'] = R
method = RI
else:
params['S'] = R
method = RP
else: #if algo == "RB":
params = {'R':(MN + 2 + Delta) * (MN + 1 + Delta) // 2, 'S':S,
'POD':True, 'PODTolerance':PODTol}
if samples == "distributed":
method = RBD
elif samples == "centered_fake":
params['S'] = R
method = RBD
else:
params['S'] = R
method = RBC
if samples == "distributed":
if sampling == "quadrature":
params['sampler'] = QS(murange, "CHEBYSHEV", scaling = rescaling,
scalingInv = rescalingInv)
# params['sampler'] = QS(murange, "GAUSSLEGENDRE", scaling = rescaling,
# scalingInv = rescalingInv)
# params['sampler'] = QS(murange, "UNIFORM", scaling = rescaling,
# scalingInv = rescalingInv)
params['S'] = [max(j, MN + 1) for j in params['S']]
elif sampling == "quadrature_total":
params['sampler'] = QST(murange, "CHEBYSHEV", scaling = rescaling,
scalingInv = rescalingInv)
params['S'] = R
else: # if sampling == "random":
params['sampler'] = RS(murange, "HALTON", scaling = rescaling,
scalingInv = rescalingInv)
params['S'] = R
elif samples == "centered_fake":
params['sampler'] = MS(murange, points = [mu0], scaling = rescaling,
scalingInv = rescalingInv)
approx = method(solver, mu0 = mu0, approxParameters = params, verbosity = verb)
if samples == "distributed": approx.samplingEngine.allowRepeatedSamples = False
approx.setupApprox()
if show_sample:
L = mutar[1]
approx.plotApprox(mutar, name = 'u_app', homogeneized = False,
what = "REAL")
approx.plotHF(mutar, name = 'u_HF', homogeneized = False, what = "REAL")
approx.plotErr(mutar, name = 'err', homogeneized = False, what = "REAL")
# approx.plotRes(mutar, name = 'res', homogeneized = False, what = "REAL")
appErr = approx.normErr(mutar)
solNorm = approx.normHF(mutar)
resNorm = approx.normRes(mutar)
RHSNorm = approx.normRHS(mutar)
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)))
if algo == "rational" and approx.N > 0:
from plot_zero_set import plotZeroSet2
muZeroVals, Qvals = plotZeroSet2(murange, murangeEff, approx, mu0,
200, [2., 2.], clip = clip)
if show_norm:
solver._solveBatchSize = 100
from plot_inf_set import plotInfSet2
muInfVals, normEx, normApp, normRes, normErr, beta = plotInfSet2(
murange, murangeEff, approx, mu0, 50, [2., 2.], clip = clip)

Event Timeline

c4science · Help