Page MenuHomec4science
Files F85447699

HelmholtzPadeTaylorApproximant.py
No OneTemporary
Actions

File Metadata

Created
Sun, Sep 29, 07:24

HelmholtzPadeTaylorApproximant.py
View Options

#!/usr/bin/env python3
import numpy as np
from context import FreeFemHelmholtzEngine as HFEngine
from context import FreeFemHelmholtzScatteringEngine as HFSEngine
from context import FreeFemHelmholtzScatteringAugmentedEngine as HFSAEngine
from context import FreeFemHSEngine as HSEngine
from context import FreeFemHSAugmentedEngine as HSAEngine
from context import ROMApproximantTaylorPade as Pade
testNo = 1
if testNo == 1:
params = {'N':3, 'M':3, 'E':3, 'sampleType':'Arnoldi', 'POD':True}
z0 = 12+.5j
ztar = 11
from FreeFem_snippets import SquareHomogeneousBubble
bdr, V, f = SquareHomogeneousBubble(kappa = 12 ** .5, theta = np.pi / 3,
n = 40)
solver = HFEngine(V, z0**.5, forcingTerm = f, DirichletBoundary = bdr)
plotter = HSEngine(solver.V)
approx = Pade(solver, plotter, k0 = z0, w = np.real(z0**.5),
approxParameters = params)
approx.plotApp(ztar, name = 'u_Pade''')
approx.plotHF(ztar, name = 'u_HF')
approx.plotErr(ztar, name = 'err')
appErr, solNorm = approx.approxError(ztar), approx.HFNorm(ztar)
print(('SolNorm:\t{}\nErr:\t{}\nErrRel:\t{}').format(solNorm, appErr,
np.divide(appErr, solNorm)))
print('\nPoles Pade'':')
print(approx.getPoles(True))
############
elif testNo == 2:
params = {'N':7, 'M':7, 'E':7, 'sampleType':'Arnoldi', 'POD':True}
z0 = 16
ztar = 3.95 ** 2
from FreeFem_snippets import SquareTransmissionDirichlet
bdr, V, n, u0 = SquareTransmissionDirichlet(nT = 2, nB = 1,
theta = np.pi * 45 / 180,
kappa = 4., n = 50)
solver = HFEngine(V, z0**.5, refractionIndex = n,
DirichletBoundary = "1,2,3,4", DirichletDatum = u0)
plotter = HSEngine(solver.V)
approx = Pade(solver, plotter, k0 = z0, w = np.real(z0**.5),
approxParameters = params, plotDer = 'ALL')
approx.plotApp(ztar, name = 'u_RB')
approx.plotHF(ztar, name = 'u_HF')
approx.plotErr(ztar, name = 'err')
appErr, solNorm = approx.approxError(ztar), approx.HFNorm(ztar)
print(('SolNorm:\t{}\nErr:\t{}\nErrRel:\t{}').format(solNorm, appErr,
np.divide(appErr, solNorm)))
print('\nPoles Pade'':')
print(approx.getPoles(True))
############
elif testNo in [3, 4]:
k0 = 2.
ktar = 1.8 - .3j
from FreeFem_snippets import SquareScatteringTB
bdrD, bdrN, V, f = SquareScatteringTB(kappa = 2, theta = np.pi / 6, n = 50)
if testNo == 3:
params = {'N':8, 'M':7, 'E':8, 'sampleType':'Krylov', 'POD':False}
solver = HFSEngine(V, k0, forcingTerm = f, DirichletBoundary = bdrD,
RobinBoundary = bdrN)
plotter = HSEngine(solver.V)
else:
params = {'N':8, 'M':7, 'E':8, 'sampleType':'Arnoldi', 'POD':True}
solver = HFSAEngine(V, k0, forcingTerm = f, DirichletBoundary = bdrD,
RobinBoundary = bdrN)
plotter = HSAEngine(solver.V, 2)
approx = Pade(solver, plotter, k0 = k0, approxParameters = params,
plotDer = 'ALL')
approx.plotApp(ktar, name = 'u_Pade''')
approx.plotHF(ktar, name = 'u_HF')
approx.plotErr(ktar, name = 'err')
appErr, solNorm = approx.approxError(ktar), approx.HFNorm(ktar)
print(('SolNorm:\t{}\nErr:\t{}\nErrRel:\t{}').format(solNorm, appErr,
np.divide(appErr, solNorm)))
print('\nPoles Pade'':')
print(approx.getPoles(True))

Event Timeline

c4science · Help