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.229093262810293, 22.6364475051225, 38.75912681055939, 56.95688372413629, 64.68904674421539, 71.61510643577651, 76.89311296541675, 86.36950975758471, 90.92730267324319, 99.90890224739357, 108.07709802016265, 123.53078335008216]
1 : [3.231012675426317, 6.0862342637755535, 6.406899224356592, 12.054533473709698, 12.672421029114952, 13.138599688915065, 14.49003224838112, 17.723681704687735, 18.66258439707, 19.3461206054422, 23.999160961162637, 25.45185223928206]
2 : [3.2203793481864573, 5.886623822608889, 5.905548429312819, 10.85374739081217, 10.901287970181603, 11.033338300478853, 12.49060860995466, 13.71319396709547, 14.184813087307623, 15.34867454866488, 16.529828077616358, 17.964540350031246]
3 : [3.2202534507267333, 5.880768708677691, 5.881908488029245, 10.709133420848218, 10.726052066013555, 10.778656219073557, 12.345960396981575, 12.760676043044072, 13.760281355560204, 14.018326044509678, 14.607718515898389, 16.592214150833055]
4 : [3.2202514669005975, 5.880494279982003, 5.880556917690833, 10.689894953487268, 10.698983267290268, 10.718281105538928, 12.323621310872804, 12.51513457509298, 13.530915636506077, 13.64635249050193, 14.013638404198998, 15.84404210325453]
5 : [3.220251432780632, 5.880477114537634, 5.880480023263915, 10.686851036556556, 10.694700891026889, 10.700890407054954, 12.319109793668042, 12.416240501305136, 13.409477890006137, 13.5525678871285, 13.74921631135577, 15.246452468516857]
6 : [3.2202514321849853, 5.880473890511285, 5.880477843318867, 10.686161008730885, 10.694032876758255, 10.695820538590603, 12.31800240065423, 12.365108987095137, 13.375722155094387, 13.48849701737709, 13.613981824881884, 14.813746365678831]
7 : [3.2202514321746216, 5.880473676484483, 5.880477751399422, 10.685969287196457, 10.693915072917129, 10.694415397844635, 12.317696219798139, 12.338610207937302, 13.365112608807758, 13.452949540217883, 13.531662728957844, 14.548118346316624]
8 : [3.220251432174444, 5.8804736647727704, 5.880477746015271, 10.685917731763398, 10.693892108137492, 10.694043457650656, 12.317604665085396, 12.326098068805, 13.3608783849769, 13.435960447006165, 13.480566199504791, 14.398036525811872]
9 : [3.2202514321744395, 5.880473664135098, 5.880477745705945, 10.685904872914715, 10.693887275313298, 10.69394943386472, 12.317575537839142, 12.320797131208757, 13.358877986662451, 13.42825810504428, 13.451525423084068, 14.315735109132918]
10 : [3.220251432174439, 5.880473664100304, 5.8804777456885, 10.685901792193867, 10.693886175010105, 10.693926453072581, 12.317565160998626, 12.318715754147926, 13.357892966568976, 13.4247939758549, 13.436298700684368, 14.270543561681373]
11 : [3.220251432174439, 5.8804736640984006, 5.880477745687526, 10.685901066925178, 10.69388591658353, 10.693920952759887, 12.317558831682375, 12.317939833810456, 13.357414569546277, 13.423228034855255, 13.428747506413684, 14.245626495227171]
12 : [3.2202514321744458, 5.880473664098324, 5.880477745687471, 10.685900897346514, 10.693885856844494, 10.693919648675587, 12.317545562913883, 12.31766812735211, 13.357188433338655, 13.422515874577, 13.425120101201918, 14.231830227382503]
13 : [3.2202514321744413, 5.880473664098325, 5.880477745687483, 10.685900857787113, 10.693885843152671, 10.693919341029755, 12.317512858984296, 12.31759852349252, 13.3570837648275, 13.422189469070833, 13.423410819523834, 14.22417801130031]
14 : [3.2202514321744413, 5.880473664098325, 5.880477745687477, 10.685900848563524, 10.693885840013818, 10.693919268629985, 12.317488500273726, 12.31758642389661, 13.357035962915184, 13.422036742638483, 13.422616246960317, 14.219930051110927]
15 : [3.220251432174443, 5.880473664098006, 5.8804777456872195, 10.685900846411087, 10.693885839291879, 10.693919251608369, 12.317478502915344, 12.317583475807117, 13.35701427462487, 13.421959821315893, 13.422254103177599, 14.217570023366102]
16 : [3.220251432174443, 5.880473664098331, 5.880477745687461, 10.685900845908144, 10.693885839123954, 10.693919247589909, 12.317474787811106, 12.317582547658269, 13.357004470637055, 13.42191372424843, 13.42209565673255, 14.216247270296387]
17 : [3.2202514321744395, 5.880473664098172, 5.880477745687558, 10.68590084579099, 10.693885839086494, 10.693919246660318, 12.317473488109693, 12.317582242666246, 13.357000040944076, 13.421883404756123, 13.422032261635998, 14.215517067734218]
18 : [3.2202514321744387, 5.880473664098511, 5.880477745687355, 10.685900845761152, 10.693885839077742, 10.693919246423855, 12.317472970427636, 12.317582125581, 13.356997740090485, 13.42186171287991, 13.422005483272814, 14.215012822715138]
19 : [3.220251432174441, 5.880473664098466, 5.880477745687607, 10.68590084575605, 10.693885839075625, 10.693919246387065, 12.317472843465227, 12.317582095997194, 13.356996926043278, 13.421853753854704, 13.421995470758452, 14.214807415116624]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.220251432174441
5.880473664098466
5.880477745687607
10.68590084575605
10.693885839075625
10.693919246387065
12.317472843465227
12.317582095997194
13.356996926043278
13.421853753854704
13.421995470758452
14.214807415116624
[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);
[ ]:
[ ]: