This page was generated from unit-2.4-Maxwell/Maxwellevp.ipynb.
2.4.1 Maxwell eigenvalue problem¶
We solve the Maxwell eigenvalue problem
\[\int \operatorname{curl} u \, \operatorname{curl} v
= \lambda \int u v\]
for \(u, v \; \bot \; \nabla H^1\) using a PINVIT solver from the ngsolve solvers module.
The orthogonality to gradient fields is important to eliminate the huge number of zero eigenvalues. The orthogonal sub-space is implemented using a Poisson projection:
\[P u = u - \nabla \Delta^{-1} \operatorname{div} u\]
The algorithm and example is take form the Phd thesis of Sabine Zaglmayr, p 145-150.
[1]:
from ngsolve import *
from ngsolve.webgui import Draw
from netgen.occ import *
[2]:
from netgen.occ import *
cube1 = Box( (-1,-1,-1), (1,1,1) )
cube2 = Box( (0,0,0), (2,2,2) )
cube2.edges.hpref=1 # mark edges for geometric refinement
fichera = cube1-cube2
Draw (fichera);
[3]:
mesh = Mesh(OCCGeometry(fichera).GenerateMesh(maxh=0.4))
mesh.RefineHP(levels=2, factor=0.2)
Draw (mesh);
[4]:
# SetHeapSize(100*1000*1000)
fes = HCurl(mesh, order=3)
print ("ndof =", fes.ndof)
u,v = fes.TnT()
a = BilinearForm(curl(u)*curl(v)*dx)
m = BilinearForm(u*v*dx)
apre = BilinearForm(curl(u)*curl(v)*dx + u*v*dx)
pre = Preconditioner(apre, "direct", inverse="sparsecholesky")
ndof = 42562
[5]:
with TaskManager():
a.Assemble()
m.Assemble()
apre.Assemble()
# build gradient matrix as sparse matrix (and corresponding scalar FESpace)
gradmat, fesh1 = fes.CreateGradient()
gradmattrans = gradmat.CreateTranspose() # transpose sparse matrix
math1 = gradmattrans @ m.mat @ gradmat # multiply matrices
math1[0,0] += 1 # fix the 1-dim kernel
invh1 = math1.Inverse(inverse="sparsecholesky")
# build the Poisson projector with operator Algebra:
proj = IdentityMatrix() - gradmat @ invh1 @ gradmattrans @ m.mat
projpre = proj @ pre.mat
evals, evecs = solvers.PINVIT(a.mat, m.mat, pre=projpre, num=12, maxit=20)
0 : [6.229093260802441, 22.636447504088302, 38.75912680627621, 56.956883714691465, 64.68904674861636, 71.61510644043861, 76.893112935674, 86.36950976505287, 90.92730267510284, 99.90890224270714, 108.07709801566953, 123.53078334715204]
1 : [3.23101267541736, 6.086234263830918, 6.40689922388243, 12.054533474392166, 12.672421028669005, 13.138599688337294, 14.49003224872494, 17.723681705711154, 18.66258439443327, 19.346120606254257, 23.99916096254276, 25.451852239184642]
2 : [3.2203793481863494, 5.886623822611208, 5.9055484292936296, 10.853747390857121, 10.901287970372485, 11.033338300559066, 12.490608610031462, 13.713193967098444, 14.18481308745683, 15.348674547325022, 16.52982807956572, 17.9645403493986]
3 : [3.220253450726731, 5.880768708677782, 5.881908488028408, 10.709133420856965, 10.726052066046002, 10.778656219154524, 12.345960396986666, 12.760676043221903, 13.760281355787631, 14.018326044196534, 14.607718515772051, 16.592214150570925]
4 : [3.2202514669005913, 5.880494279982009, 5.880556917690803, 10.689894953491189, 10.698983267292547, 10.718281105571734, 12.323621310871333, 12.515134575197996, 13.530915636489448, 13.64635249060352, 14.013638403865045, 15.844042102914896]
5 : [3.2202514327806293, 5.880477114537637, 5.880480023263932, 10.686851036557918, 10.69470089102668, 10.70089040706558, 12.319109793666804, 12.416240501357652, 13.409477889986652, 13.552567887220418, 13.74921631111496, 15.246452468181433]
6 : [3.220251432184981, 5.8804738905112846, 5.880477843318869, 10.68616100873137, 10.694032876758154, 10.695820538593537, 12.318002400653674, 12.365108987116837, 13.375722155086057, 13.488497017420366, 13.613981824738378, 14.813746365423762]
7 : [3.220251432174619, 5.880473676484498, 5.880477751399407, 10.68596928719661, 10.693915072917108, 10.694415397845358, 12.317696219797913, 12.338610207944908, 13.365112608805404, 13.452949540236407, 13.531662728872357, 14.548118346151323]
8 : [3.2202514321744444, 5.880473664772752, 5.880477746015261, 10.68591773176346, 10.693892108137472, 10.694043457650832, 12.317604665085298, 12.326098068801393, 13.360878384978147, 13.435960446995683, 13.480566199342151, 14.398036525644585]
9 : [3.2202514321744387, 5.880473664135108, 5.88047774570595, 10.685904872915353, 10.693887275313868, 10.69394943386634, 12.317575537838845, 12.320797131762573, 13.35887798664392, 13.428258104515187, 13.451525432590538, 14.315735127639284]
10 : [3.220251432174438, 5.880473664100337, 5.880477745688508, 10.685901792194137, 10.693886175010228, 10.69392645307352, 12.317565160998589, 12.318715754683822, 13.357892966462508, 13.424793977304303, 13.436298713900722, 14.270543579031608]
11 : [3.22025143217444, 5.880473664098397, 5.880477745687531, 10.685901066925132, 10.693885916583659, 10.693920952759672, 12.317558831681046, 12.317939833962908, 13.357414569383426, 13.423228036137857, 13.42874751450684, 14.245626502366115]
12 : [3.2202514321744435, 5.8804736640983055, 5.880477745687466, 10.685900897346416, 10.693885856844496, 10.693919648675605, 12.317545562910128, 12.317668127390549, 13.357188433399417, 13.422515876487354, 13.42512010482387, 14.231830220058612]
13 : [3.2202514321744395, 5.880473664098327, 5.88047774568747, 10.685900857785523, 10.693885843151753, 10.693919341024609, 12.31751285248995, 12.317598518189332, 13.357083765328747, 13.422189473527713, 13.423410327256816, 14.22417651114587]
14 : [3.220251432174437, 5.880473664098315, 5.8804777456874735, 10.685900848562289, 10.693885840013465, 10.693919268627546, 12.317488493795265, 12.317586421690708, 13.357035963363064, 13.422036772954359, 13.422615840135368, 14.21992807323386]
15 : [3.22025143217444, 5.880473664098218, 5.880477745687506, 10.685900846410991, 10.69388583929155, 10.69391925160843, 12.317478499925098, 12.317583475029808, 13.357014276548009, 13.421959836162374, 13.422253850312332, 14.217568533892093]
16 : [3.2202514321744444, 5.880473664098459, 5.880477745687512, 10.685900845908849, 10.693885839125029, 10.693919247606868, 12.317474824472727, 12.31758255808151, 13.357004463539068, 13.42191393052443, 13.42209669056108, 14.216257009550949]
17 : [3.220251432174442, 5.8804736640978685, 5.880477745687415, 10.685900845791144, 10.69388583908695, 10.693919246666914, 12.31747351176069, 12.31758224908956, 13.357000005827043, 13.42188374850517, 13.422033161599158, 14.215528130299358]
18 : [3.2202514321744418, 5.880473664098388, 5.880477745687522, 10.685900845762367, 10.693885839077437, 10.693919246433953, 12.317473028888404, 12.317582139365612, 13.35699775116248, 13.42186477833449, 13.42200751035509, 14.2150944164215]
19 : [3.2202514321744435, 5.880473664098384, 5.880477745687468, 10.685900845756862, 10.693885839075909, 10.693919246390264, 12.317472872381021, 12.317582102899559, 13.356996938831932, 13.421855965170014, 13.421996718321212, 14.214873770247992]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.2202514321744435
5.880473664098384
5.880477745687468
10.685900845756862
10.693885839075909
10.693919246390264
12.317472872381021
12.317582102899559
13.356996938831932
13.421855965170014
13.421996718321212
14.214873770247992
[7]:
gfu = GridFunction(fes, multidim=len(evecs))
for i in range(len(evecs)):
gfu.vecs[i].data = evecs[i]
Draw (Norm(gfu), mesh, order=4, min=0, max=2);
[ ]:
[ ]: