This page was generated from unit-1.7-helmholtz/helmholtz.ipynb.

1.7 Complex-valued waves

In NGSolve finite element spaces can be built over the complex field and the resulting complex matrix systems can be solved. This tutorial shows how to compute the solution of the Helmholtz equation with impedance boundary conditions in complex arithmetic. The boundary value problem is to find \(u\) satisfying

\[-\Delta u - \omega^2 u = f\qquad \text{ in } \Omega\]

together with the impedance (outgoing) boundary condition

\[\frac{\partial u }{ \partial n} - i \omega u = 0 \quad \text{ on } \partial \Omega\]

where \(i =\) 1j is the imaginary unit.

from ngsolve import *
from ngsolve.webgui import Draw
from netgen.occ import *
from netgen.geom2d import SplineGeometry
air = Circle((0.5, 0.5), 0.8).Face() = 'outer'
scatterer = MoveTo(0.7, 0.3).Rectangle(0.05, 0.4).Face() = 'scat'
geo = OCCGeometry(air - scatterer, dim=2)
mesh = Mesh(geo.GenerateMesh(maxh=0.05))

Declare a complex-valued finite element space

fes = H1(mesh, order=5, complex=True)
u, v = fes.TnT()
omega = 100
pulse = 5e4*exp(-(40**2)*((x-0.5)*(x-0.5) + (y-0.5)*(y-0.5)))
Draw(pulse, mesh, order=3)

Forming the system

The weak form for \(u \in H^1\):

\[\int_\Omega\big[ \nabla u \cdot \nabla \bar v - \omega^2 u \bar v \big] \, dx - i \,\omega\, \int_{\partial \Omega} u \bar v \, ds = \int_{\Omega} f \bar v\]

for all \(v\) in \(H^1\).

# Forms
a = BilinearForm(fes)
a += grad(u)*grad(v)*dx - omega**2*u*v*dx
a += -omega*1j*u*v * ds("outer")

f = LinearForm(pulse * v * dx).Assemble();


gfu = GridFunction(fes, name="u") = a.mat.Inverse() * f.vec
Draw(gfu, mesh, min=-1, max=1, order=3, animate_complex=True);

Open controls and explore webgui’s further viewing options such as viewing real and imaginary parts, viewing absolute value, and viewing with deformation turned on.