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.2290932597990825, 22.636447495057702, 38.759126804149616, 56.95688371479104, 64.68904673318777, 71.61510645856109, 76.89311293359025, 86.36950974553933, 90.92730267650916, 99.90890224160074, 108.07709799915821, 123.53078333601272]
1 : [3.2310126754138575, 6.08623426370834, 6.406899224083785, 12.054533473884993, 12.672421028372735, 13.138599686951272, 14.49003224620985, 17.72368170638937, 18.662584395930168, 19.346120606901422, 23.999160958716793, 25.45185223766388]
2 : [3.220379348186336, 5.886623822605638, 5.905548429308459, 10.853747390971815, 10.90128796982568, 11.03333830065748, 12.490608609859404, 13.713193966580635, 14.184813087318906, 15.348674547815579, 16.52982807769545, 17.964540347986482]
3 : [3.2202534507267275, 5.880768708677494, 5.881908488028898, 10.709133420829584, 10.726052065988616, 10.778656219195971, 12.345960396947204, 12.76067604262261, 13.760281355653893, 14.018326044004773, 14.607718515489994, 16.59221414958171]
4 : [3.2202514669005926, 5.880494279982013, 5.880556917690805, 10.68989495348247, 10.69898326728357, 10.71828110558335, 12.323621310857753, 12.515134574836173, 13.530915636346348, 13.64635249063558, 14.013638403787326, 15.844042101752738]
5 : [3.220251432780625, 5.880477114537647, 5.880480023263914, 10.686851036556146, 10.694700891025288, 10.700890407067513, 12.319109793661797, 12.41624050114092, 13.409477889948343, 13.552567887241976, 13.749216310994187, 15.246452467122802]
6 : [3.2202514321849827, 5.880473890511294, 5.880477843318854, 10.686161008731094, 10.694032876757841, 10.695820538593546, 12.318002400651816, 12.365108986996558, 13.37572215506594, 13.488497017445612, 13.613981824599188, 14.813746364705365]
7 : [3.220251432174622, 5.880473676484467, 5.880477751399424, 10.685969287196569, 10.693915072917028, 10.694415397845303, 12.31769621979728, 12.338610207885845, 13.365112608790442, 13.452949540255998, 13.53166272876122, 14.548118345728923]
8 : [3.220251432174444, 5.88047366477277, 5.880477746015268, 10.68591773176345, 10.693892108137469, 10.694043457650842, 12.31760466508513, 12.326098068775995, 13.360878384967355, 13.435960447028384, 13.480566199320787, 14.398036525382338]
9 : [3.220251432174438, 5.880473664135078, 5.88047774570596, 10.685904872915662, 10.693887275314308, 10.693949433866518, 12.317575537838206, 12.320797131368685, 13.358877986635894, 13.428258103358283, 13.45152542398028, 14.315735120248586]
10 : [3.220251432174439, 5.880473664100285, 5.880477745688488, 10.685901792194251, 10.693886175010649, 10.693926453073553, 12.317565160998592, 12.318715753997143, 13.357892966463336, 13.424793976264317, 13.436298704567632, 14.270543567924637]
11 : [3.2202514321744404, 5.880473664098421, 5.8804777456875215, 10.685901066925029, 10.693885916583877, 10.693920952758948, 12.317558831676285, 12.317939833399313, 13.357414569351556, 13.423228035619706, 13.428747506526332, 14.245626482489008]
12 : [3.220251432174442, 5.880473664098299, 5.880477745687489, 10.685900897346341, 10.69388585684475, 10.69391964867532, 12.317545562865641, 12.31766812707094, 13.357188433361182, 13.42251587653664, 13.425120099377907, 14.23183019805073]
13 : [3.2202514321744387, 5.880473664098321, 5.880477745687445, 10.685900857789344, 10.693885843151833, 10.693919341041056, 12.31751286367068, 12.317598526649194, 13.357083769917889, 13.42218946839027, 13.423410983683144, 14.224178955127895]
14 : [3.2202514321744378, 5.880473664098316, 5.8804777456874575, 10.685900848563273, 10.693885840012488, 10.693919268633545, 12.317488501642039, 12.317586424464327, 13.357035961200681, 13.422036762081625, 13.422616130634463, 14.219929727511293]
15 : [3.220251432174439, 5.880473664098237, 5.880477745687068, 10.685900846410881, 10.693885839291257, 10.693919251610193, 12.317478504660908, 12.317583476448924, 13.35701427352718, 13.421959836607755, 13.4222540958462, 14.217569832757551]
16 : [3.2202514321744395, 5.880473664098373, 5.880477745687224, 10.685900845908503, 10.693885839124952, 10.693919247607951, 12.317474829013609, 12.317582559465293, 13.357004464638342, 13.421913948718169, 13.422096879886185, 14.216258236112452]
17 : [3.2202514321744378, 5.880473664098298, 5.8804777456872435, 10.685900845790915, 10.693885839086708, 10.693919246666882, 12.317473510157791, 12.317582249069735, 13.357000022353269, 13.421883645275964, 13.422032627691188, 14.215523974588015]
18 : [3.2202514321744395, 5.880473664098347, 5.880477745687406, 10.685900845763642, 10.693885839077526, 10.69391924644599, 12.317473044991312, 12.317582141728746, 13.356998027156385, 13.421865904557777, 13.422007215030018, 14.215119952198785]
19 : [3.2202514321744355, 5.880473664098193, 5.88047774568781, 10.685900845757363, 10.693885839075698, 10.69391924639401, 12.317472879727099, 12.317582103851871, 13.356997124878026, 13.421856650257704, 13.421996494091124, 14.214893790820598]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.2202514321744355
5.880473664098193
5.88047774568781
10.685900845757363
10.693885839075698
10.69391924639401
12.317472879727099
12.317582103851871
13.356997124878026
13.421856650257704
13.421996494091124
14.214893790820598
[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);
[ ]:
[ ]: