Page MenuHomec4science
Files F58996382

anisotropic_square.py
No OneTemporary
Actions

File Metadata

Created
Fri, Apr 19, 00:23

anisotropic_square.py
View Options

### example from Smetana, Zahm, Patera. Randomized residual-based error
### estimators for parametrized equations.
import numpy as np
import matplotlib.pyplot as plt
from anisotropic_square_engine import (AnisotropicSquareEngine as engine,
AnisotropicSquareEnginePoles as plsEx)
from rrompy.reduction_methods import (
RationalInterpolantGreedyPivotedGreedyPoleMatch as RIGPG)
from rrompy.parameter.parameter_sampling import (QuadratureSampler as QS,
SparseGridSampler as SGS)
zs, Ls = [10., 50.], [.2, 1.2]
z0, L0, n = np.mean(zs), np.mean(Ls), 50
murange = [[zs[0], Ls[0]], [zs[-1], Ls[-1]]]
np.random.seed(4020)
mu = [zs[0] + np.random.rand() * (zs[-1] - zs[0]),
Ls[0] + np.random.rand() * (Ls[-1] - Ls[0])]
solver = engine(z0, L0, n)
solver.cutOffPolesRMinRel, solver.cutOffPolesRMaxRel = -2, 2
fighandles = []
params = {"POD": True, "nTestPoints": 100, "greedyTol": 1e-4, "S": 3,
"polybasisMarginal": "PIECEWISE_LINEAR_UNIFORM",
"polybasis": "LEGENDRE", "samplerPivot":QS(zs, "UNIFORM"),
"samplerTrainSet":QS(zs, "UNIFORM"), "maxIterMarginal":50,
"errorEstimatorKind":"LOOK_AHEAD_RES",
"errorEstimatorKindMarginal":"LOOK_AHEAD", "SMarginal": 3,
"greedyTolMarginal": 1e-2, "samplerMarginal":SGS(Ls),
"radialDirectionalWeightsMarginal": [4.], "matchingWeight": 1.,
"badPoleCorrection": "POLYNOMIAL", "autoCollapse": 1}
for shared in [1., .5, 0.]:
print("Testing shared ratio {}.".format(shared))
params["matchingShared"] = shared
approx = RIGPG([0], solver, mu0 = [z0, L0], approxParameters = params,
verbosity = 5)
approx.setupApprox("ALL")
verb = approx.verbosity
approx.verbosity = 0
tspace = np.linspace(Ls[0], Ls[-1], 100)
for j, t in enumerate(tspace):
plsE = plsEx(t, 0., zs[-1])
pls = approx.getPoles([None, t])
pls[np.abs(np.imag(pls)) > 1e-5] = np.nan
if j == 0:
polesE = np.empty((len(tspace), len(plsE)))
poles = np.empty((len(tspace), len(pls)))
polesE[:] = np.nan
if len(plsE) > polesE.shape[1]:
nanR = np.empty((len(tspace), len(plsE) - polesE.shape[1]))
nanR[:] = np.nan
polesE = np.hstack((polesE, nanR))
polesE[j, : len(plsE)] = np.real(plsE)
poles[j] = np.real(pls)
approx.verbosity = verb
fighandles += [plt.figure(figsize = (17, 5))]
ax1 = fighandles[-1].add_subplot(1, 2, 1)
ax2 = fighandles[-1].add_subplot(1, 2, 2)
ax1.plot(poles, tspace)
ax1.set_ylim(Ls)
ax1.set_xlabel("mu_1")
ax1.set_ylabel("mu_2")
ax1.grid()
ax2.plot(polesE, tspace, "k-.", linewidth = 1)
ax2.plot(poles, tspace)
for mm in approx.musMarginal:
ax2.plot(zs, [mm[0, 0]] * 2, "k--", linewidth = 1)
ax2.set_xlim(zs)
ax2.set_ylim(Ls)
ax2.set_xlabel("mu_1")
ax2.set_ylabel("mu_2")
ax2.grid()
plt.show()
print("\n")

Event Timeline

c4science · Help