Page MenuHomec4science

LagrangeSweep.py
No OneTemporary

File Metadata

Created
Mon, Nov 4, 11:37

LagrangeSweep.py

from copy import copy
import numpy as np
#from rrompy.hfengines.scipy import HelmholtzBoxScatteringProblemEngine as \
# HBSPE
from rrompy.hfengines.scipy import HelmholtzCavityScatteringProblemEngine as \
HBSPE
from rrompy.reduction_methods.lagrange import ApproximantLagrangePade as Pade
from rrompy.reduction_methods.lagrange import ApproximantLagrangeRB as RB
from rrompy.utilities.parameter_sweeper import ParameterSweeper as Sweeper
from rrompy.utilities.parameter_sampling import QuadratureSampler as QS
from rrompy.utilities.parameter_sampling import WarpingFunction as WF
from rrompy.utilities.parameter_sampling import WarpingSampler as WS
verb = 0
npoints = 100
homog = True
homog = False
LSratio = 2. / 3.
sampling = "Uniform"
#sampling = "Cheby"
#sampling = "SinC"
assert LSratio <= 1. + np.finfo(float).eps
ks = [10 + 0.j, 14 + 0.j]
solver = HBSPE(kappa = 3, n = 15)
solver.omega = np.real(np.mean(ks))
mutars = np.linspace(9, 15, npoints)
homogMSG = "Homog"
if not homog: homogMSG = "Non" + homogMSG
filenamebase = '../data/output/ScatteringSquareLSweep'
#filenamebase = '../data/plots/LagrangeScatteringSquare1p5'
#filenamebase = filenamebase + sampling + "/HelmholtzBoxLSweep" + homogMSG
k0 = np.mean(ks)
shift = 7
shift = np.int(8 / LSratio - 1)
nsets = 3
stride = np.int(8 / LSratio)
Smax = stride * (nsets - 1) + shift + 2
rescaling = lambda x: np.power(x, 2.)
rescalingInv = lambda x: np.power(x, .5)
if sampling == "Uniform":
polyBasis = "MONOMIAL"
sampler = QS(ks, "UNIFORM", rescaling, rescalingInv)
if sampling == "Cheby":
polyBasis = "CHEBYSHEV"
sampler = QS(ks, "CHEBYSHEV", rescaling, rescalingInv)
if sampling == "SinC":
polyBasis = "MONOMIAL"
warping = WF(call = lambda x: (x - 2. * (1. - LSratio) / np.pi
* np.sin(np.pi * x)),
repr = "x-{}*sin(pi*x)".format(2. * (1. - LSratio) / np.pi))
sampler = WS(ks, warping, rescaling, rescalingInv)
paramsPade = {'S':Smax, 'POD':True, 'basis':polyBasis, 'sampler':sampler}
paramsRB = copy(paramsPade)
paramsPoly = copy(paramsPade)
paramsSetsPade = [None] * nsets
paramsSetsRB = [None] * nsets
paramsSetsPoly = [None] * nsets
for i in range(nsets):
paramsSetsPade[i] = {'N': np.int(LSratio * (stride * i + shift + 1)),
'M': np.int(LSratio * (stride * i + shift + 1)),
'S': stride * i + shift + 2}
paramsSetsRB[i] = {'R': np.int(LSratio * (stride * i + shift + 1)),
'S': stride * i + shift + 2}
paramsSetsPoly[i] = {'N': 0,
'M': np.int(LSratio * (stride * i + shift + 1)),
'S': stride * i + shift + 2}
appPade = Pade(solver, mu0 = k0, approxParameters = paramsPade,
verbosity = verb, homogeneized = homog)
appRB = RB(solver, mu0 = k0, approxParameters = paramsRB,
verbosity = verb, homogeneized = homog)
appPoly = Pade(solver, mu0 = k0, approxParameters = paramsPoly,
verbosity = verb, homogeneized = homog)
sweeper = Sweeper(mutars = mutars, mostExpensive = 'Approx')
sweeper.ROMEngine = appPade
sweeper.params = paramsSetsPade
filenamePade = sweeper.sweep(filenamebase + 'Pade.dat')
sweeper.ROMEngine = appRB
sweeper.params = paramsSetsRB
filenameRB = sweeper.sweep(filenamebase + 'RB.dat')
sweeper.ROMEngine = appPoly
sweeper.params = paramsSetsPoly
filenamePoly = sweeper.sweep(filenamebase + 'Poly.dat')
sweeper.plotCompare([filenamePade, filenameRB, filenamePoly], ['muRe'],
['normHF', 'normApp'], ['S'], onePlot = True,
save = filenamebase + 'Norm', ylims = {'top' : 1e1},
saveFormat = "png", labels = ["Pade'", "RB", "Poly"],
# figsize = (5, 3.75))
figsize = (10, 7.5))
sweeper.plotCompare([filenamePade, filenameRB, filenamePoly], ['muRe'],
['normResRel'], ['S'], save = filenamebase + 'Res',
ylims = {'top' : 1e1}, saveFormat = "png",
labels = ["Pade'", "RB", "Poly"],
# figsize = (5, 3.75))
figsize = (10, 7.5))
sweeper.plotCompare([filenamePade, filenameRB, filenamePoly], ['muRe'],
['normErrRel'], ['S'], save = filenamebase + 'Err',
ylims = {'top' : 1e1}, saveFormat = "png",
labels = ["Pade'", "RB", "Poly"],
# figsize = (5, 3.75))
figsize = (10, 7.5))

Event Timeline