# pml on waveguides

[attachment=undefined]schema.PNG[/attachment] [attachment=undefined]schema.PNG[/attachment]
https://ngsolve.org/media/kunena/attachments/1052/pulse2d_waveguide_pml.pyDear NGSolve developers,
I have trouble with pml when used for a waveguide problem. Consider a 2D geometry: \Omega = [0, Range] X [0, Depth] and the Helmholtz equation:

\Delta u + k^2 u = -f in \Omega

with

\frac{\partial u}{\partial n} = 0 on \Gamma_N (top and bottom boundaries).
To the left and right, we have \Gamma_R where we use the PML, so that \partial \Omega = \Gamma_N \cup \Gamma_R

I wish to prescribe f as a point source at (x0,y0) located in \Omega so I followed the pulse implementation given in the NGSolve’s the i-tutorials page (1.7 Complex-valued waves — NGS-Py 6.2.2302-87-ga5a5eff3b documentation)
For the geometry described above, I created the mesh externally in GMsh and specified the regions for PML (to left and right of the source) by using the domain keywords i used in GMsh. While I can see that the solution decays within the PML regions, I would have expected it to undergo an exponential decay due to geometrical spreading loss as it travels to left and right of the point source. Does the NGSolve pml implementation take into account the modal behaviour of waves in a confined cavity unlike fully exterior problem?
I have attached the script that I am using and a schematic of the problem to get an idea.

Not sure why the attachment is not visible in my earlier message. Here i am pasting the code below:

[code]from ngsolve import *
from netgen.geom2d import SplineGeometry
import numpy as np
import subprocess

Rmax = 3e3
Depth = 1e2
Lpml = 1e2

# wrap it into an NGSolve mesh

from ngsolve import *
mesh = Mesh(mesh)
Draw(mesh)
pcR = pml.HalfSpace((Rmax,0.),(1,0),2j)
pcL = pml.HalfSpace((0,0.),(-1,0),2j)

# pc2 = pml.Cartesian((Rmax,0.),(Rmax+Lpml,Depth),2j) # doesn’t absorb!

mesh.SetPML(pcR,“pmlR”)
mesh.SetPML(pcL,“pmlL”)

fes = H1(mesh, order=2, complex=True)
u, v = fes.TnT()

# Wavenumber & source

cspeed = 1500.
freq = 20.
x0=0.
y0=36.
omega = 2.np.pifreq / cspeed
ScaleFactor = 1e1

a = BilinearForm(fes, symmetric=True)
a.Assemble()

pulse = ScaleFactor exp(-(omega**2)((x-x0)(x-x0) + (y-y0)(y-y0)))
f = LinearForm(fes)
f += -pulse * v * dx
f.Assemble()
gfu = GridFunction(fes, name=“u”)
gfu.vec.data = a.mat.Inverse() * f.vec

vtk = VTKOutput(ma=mesh, coefs=[gfu.real], names = [“ureal”],filename=“wgfem_real”,subdivision=3)
vtk.Do()
vtk = VTKOutput(ma=mesh, coefs=[gfu.imag], names = [“uimag”],filename=“wgfem_imag”,subdivision=3)
vtk.Do()

[/code]

can you attach the mesh as well? In your example the left side pml starts at 0 and the right at Rmax, the cartesian at Rmax and Rmax+Lpml
Cartesian specifies a box, and outside of the box there are pml layers
Best
Christopher

Hi Christopher,