pml on waveguides

Main Forum Old-Forum - do not post here
1

[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.
Many thanks for your help.

2

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
from netgen.read_gmsh import ReadGmsh

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

import the Gmsh file to a Netgen mesh object

mesh = ReadGmsh(‘…/…/…/meshes/waveguide2d_sym.msh’)

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 += grad(u)grad(v)dx - omegaomegauvdx
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]

3

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

4

Hi Christopher,
Thanks for your reply.

[quote]Cartesian specifies a box, and outside of the box there are pml layers
[/quote]
Yes, but i have commented that specific line for cartesian pml as i did not understand how to use that or why it did not show any reduction in wave amplitude in the pml region. Most likely i have got the halfspace pml wrong as well as by definition it should be applied on one side alone, right?
I have attached my geo file below.
Thanks again for your help.
https://ngsolve.org/media/kunena/attachments/1052/waveguide2d_sym.geo
https://ngsolve.org/media/kunena/attachments/1052/waveguide2d_sym.geo

Attachment: waveguide2d_sym_2019-12-20-5.geo