sector_angle_engine.py
No OneTemporary
Actions

Subscribers

None

File Metadata

	import numpy as np
	import fenics as fen
	import mshr
	from rrompy.utilities.base.decorators import nonaffine_construct
	from rrompy.hfengines.fenics_engines import HelmholtzProblemEngine
	from rrompy.parameter import parameterMap as pMap

	class SectorAngleEngine(HelmholtzProblemEngine):
	def __init__(self, k0:float, t0:float, n:int):
	super().__init__(mu0 = [k0, t0])
	self._affinePoly = False
	self.npar = 2
	self.parameterMap = pMap([2., 1.])
	mesh = mshr.generate_mesh(
	mshr.Circle(fen.Point(0., 0.), 1.)
	- mshr.Rectangle(fen.Point(-1., -1.), fen.Point(0., 1.))
	- mshr.Rectangle(fen.Point(-1., -1.), fen.Point(1., 0.)), n)
	self.V = fen.FunctionSpace(mesh, "P", 1)
	x, y = fen.SpatialCoordinate(self.V.mesh())[:]
	self.forcingTerm = [fen.exp(x + y) * (1. - x 2. - y 2.),
	fen.exp(x - y) * (1. - x 2. - y 2.)]
	self._tBoundary = np.nan

	def setBoundary(self, t:float):
	while hasattr(t, "__len__"): t = t[0]
	if not np.isclose(t, self._tBoundary):
	eps = 1e-2
	self._tBoundary = t
	self.DirichletBoundary = lambda x, on_boundary: (
	on_boundary and x[0] >= eps
	and x[1] <= eps + np.sin(t * np.pi / 2.))
	self.NeumannBoundary = "REST"

	@nonaffine_construct
	def A(self, mu = [], der = 0):
	mu = self.checkParameter(mu)
	self.setBoundary(mu(1))
	return HelmholtzProblemEngine.A(self, mu, der)

	@nonaffine_construct
	def b(self, mu = [], der = 0):
	mu = self.checkParameter(mu)
	self.setBoundary(mu(1))
	return HelmholtzProblemEngine.b(self, mu, der)