Eigenvalue problems

26. Eigenvalue problems#

The Helmholtz equation (with hard boundary conditions) is singular when \(\omega^2\) is an eigenvalue of the negative Laplace operator, i.e.

\[ -\Delta u = \lambda u \]

The weak form is to find eigen pairs \((\lambda, u) \in {\mathbf R} \times H^1\) such that

\[ \int \nabla u \nabla v \, dx = \lambda \int u v \, dx \]

Finite element discretization leads to the generalized matrix eigenvalue problem

\[ A u = \lambda M u. \]
from ngsolve import *
from ngsolve.webgui import Draw

mesh = Mesh(unit_square.GenerateMesh(maxh=0.1))

fes = H1(mesh, order=3, dirichlet=".*")
u,v = fes.TnT()
a = BilinearForm(grad(u)*grad(v)*dx).Assemble()
m = BilinearForm(u*v*dx).Assemble()

Eigenvectorsare stationary points of the Rayleigh quotient

\[ \frac{\left< A u, u \right> }{ \left< M u, u \right> } \]

This characterization is used in the LOBPCG (Locally Optimal Block Preconditioned Conjugate Gradient) method.

import ngsolve.solvers

pre = a.mat.Inverse(freedofs=fes.FreeDofs())
evals, evecs = solvers.LOBPCG(a.mat, m.mat, pre=pre, num=12, maxit=20, printrates=True)

for lam in evals:
    print (lam/pi**2)
0 : [19.811240988942, 54.338371257028705, 59.2512871605776, 128.1268173298293, 140.08299253955752, 170.8241639335056, 199.30981821142643, 322.82217832369565, 356.67867106662027, 461.590313503806, 517.387088963402, 625.7205006210877]
1 : [19.739301603250873, 49.38806350192722, 49.39964851873114, 80.26947296223456, 99.91354180553509, 101.16717245415818, 130.94799280504282, 150.73757278119194, 180.39373250563222, 199.2970669790016, 243.37890518670426, 275.00690221583096]
2 : [19.7392097361078, 49.34817916499002, 49.348239248633334, 78.97556274229088, 98.71753152603867, 98.77541632451623, 128.51404525146344, 129.49669541839535, 168.96924037508037, 173.53824643625023, 184.9648307575386, 204.37620613281666]
3 : [19.739209703758846, 49.348050576168006, 49.34805663276662, 78.95714216946482, 98.69661976740971, 98.69738490103187, 128.31140649384957, 128.32786041915872, 167.8358599192328, 168.23479011439866, 177.96647880522892, 197.8093379364253]
4 : [19.739209703752923, 49.34805036204214, 49.348056346966274, 78.95704924286466, 98.69648901545288, 98.69653276278288, 128.3063553990993, 128.30659054351054, 167.78848149696745, 167.8214223840422, 177.68669578321638, 197.43802071565338]
5 : [19.739209703752955, 49.348050361663745, 49.348056346575156, 78.95704877856993, 98.6964879939008, 98.69652648656134, 128.30613642309675, 128.30637471269, 167.7871194967775, 167.79291728480277, 177.66540362579084, 197.40227381534777]
6 : [19.739209703752934, 49.34805036166209, 49.34805634657469, 78.9570487766552, 98.69648798064453, 98.69652639235422, 128.30613002138543, 128.30636997448593, 167.78703931304648, 167.78762793676438, 177.65937008686828, 197.39974063572484]
7 : [19.739209703752945, 49.34805036166322, 49.348056346574566, 78.95704877664484, 98.6964879804717, 98.69652639140124, 128.30612985593964, 128.30636985368807, 167.78693710315147, 167.78714727828526, 177.65812002873284, 197.39944830695546]
8 : [19.739209703752945, 49.34805036166305, 49.3480563465747, 78.95704877664548, 98.69648798046988, 98.69652639139198, 128.30612985146632, 128.30636985120063, 167.78690759595779, 167.7871302491841, 177.65792228063992, 197.39942301576815]
9 : [19.739209703752934, 49.348050361662885, 49.34805634657521, 78.95704877664546, 98.6964879804699, 98.69652639139233, 128.3061298513787, 128.30636985114035, 167.78690570033996, 167.787129405907, 177.65790627722052, 197.39942161697064]
10 : [19.73920970375293, 49.348050361663425, 49.34805634657493, 78.95704877664546, 98.69648798047024, 98.69652639139207, 128.30612985137583, 128.3063698511375, 167.78690552992728, 167.78712932978502, 177.65790470624842, 197.39942139987554]
11 : [19.739209703752923, 49.34805036166256, 49.348056346574886, 78.95704877664511, 98.69648798047004, 98.69652639139213, 128.306129851376, 128.306369851137, 167.78690551927116, 167.78712932508057, 177.65790458107364, 197.39942136463657]
12 : [19.739209703752937, 49.348050361662935, 49.34805634657501, 78.9570487766455, 98.69648798046983, 98.6965263913923, 128.3061298513753, 128.30636985113705, 167.7869055183964, 167.78712932468494, 177.65790456981242, 197.39942133489967]
13 : [19.73920970375293, 49.34805036166266, 49.34805634657486, 78.95704877664535, 98.69648798047004, 98.69652639139244, 128.3061298513754, 128.30636985113728, 167.78690551829754, 167.78712932464003, 177.65790456854336, 197.39942121051288]
14 : [19.73920970375294, 49.34805036166319, 49.348056346575035, 78.9570487766453, 98.69648798046983, 98.69652639139196, 128.30612985137603, 128.30636985113728, 167.78690551829226, 167.78712932463745, 177.65790456847049, 197.39942031804236]
15 : [19.739209703752934, 49.34805036166255, 49.348056346574204, 78.95704877664524, 98.69648798047, 98.69652639139181, 128.30612985137566, 128.30636985113708, 167.7869055182919, 167.78712932463756, 177.65790456846702, 197.3994137246542]
16 : [19.739209703752934, 49.348050361664114, 49.34805634657491, 78.95704877664554, 98.6964879804701, 98.69652639139242, 128.30612985137563, 128.30636985113756, 167.78690551829263, 167.78712932463742, 177.65790456846696, 197.3993670689742]
17
 : [19.739209703752934, 49.34805036166221, 49.34805634657463, 78.95704877664478, 98.69648798047034, 98.6965263913921, 128.3061298513762, 128.30636985113733, 167.78690551829217, 167.78712932463753, 177.65790456846668, 197.399157647413]
18
 : [19.739209703752934, 49.3480503616639, 49.34805634657417, 78.95704877664504, 98.69648798046975, 98.69652639139161, 128.306129851376, 128.3063698511373, 167.786905518292, 167.78712932463765, 177.6579045684665, 197.39905618813304]
19
 : [19.739209703752937, 49.34805036166317, 49.34805634657462, 78.95704877664522, 98.69648798047, 98.69652639139208, 128.30612985137557, 128.30636985113694, 167.78690551829186, 167.78712932463762, 177.65790456846617, 197.39904610944936]
2.0000000913485674
5.000002873085407
5.00000347948372
8.000021638955491
10.000044983522983
10.000048875363074
13.000128945108862
13.000153262168757
17.000367866798427
17.000390543122286
18.00050917429448
20.000705001677822
gfu = GridFunction(fes, multidim=0)
for evec in evecs:
    gfu.AddMultiDimComponent(evec)

Draw (gfu, min=-2, max=2);