It’s been a long time since I’ve looked at waveguide problems. But here’s an example with assumed perfect conductor BCs
You need vtkplotter
if you want the plotting, installed via pip3 install vtkplotter --user
. Otherwise just remove the plotting lines. Just plotting the first 4 propagating modes.
Also I’m not sure how correct this is. You should construct a test with a solution known a priori to make sure you’re solving what you think you’re solving and don’t use this example as a black box.
from dolfin import *
width = 1.0
height = 0.5
mesh = RectangleMesh(Point(0, 0), Point(width, height), 8, 4, "left/right")
V = FunctionSpace(mesh, "Nedelec 1st kind H(curl)", 2)
v = TestFunction(V)
u = TrialFunction(V)
s = inner(curl(v), curl(u))*dx
t = inner(v, u)*dx
L = dot(Constant((0.0, 0.0)), v)*dx
# Hermitian assembly
perf_conductor_bc = DirichletBC(V, Constant((0, 0)), "on_boundary")
S = PETScMatrix()
T = PETScMatrix()
dummy = PETScVector()
assemble_system(s, L, A_tensor=S, b_tensor=dummy, bcs=perf_conductor_bc)
assemble_system(t, L, A_tensor=T, b_tensor=dummy, bcs=perf_conductor_bc)
esolver = SLEPcEigenSolver(S, T)
esolver.parameters["spectrum"] = "smallest real"
esolver.parameters["solver"] = "lapack"
esolver.solve()
from vtkplotter.dolfin import plot
cutoff = None
count = 0
for i in range(S.size(1)):
(lr, lc, rr, rc) = esolver.get_eigenpair(i)
u = Function(V)
if lr > 1 and lc == 0:
cutoff = sqrt(lr)
u.vector()[:] = rr
plot(u, shape=((2,2)), at=count)
count += 1
if count > 4:
break
if cutoff is None:
print("Unable to find dominant mode")
else:
print("Cutoff frequency:", cutoff)