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.229093259408164, 22.63644749871427, 38.759126811031365, 56.95688371834956, 64.68904673494364, 71.61510645722878, 76.89311293960826, 86.36950975763345, 90.927302674657, 99.90890224082283, 108.0770980188762, 123.53078334016317]
1 : [3.231012675410567, 6.086234263822988, 6.406899224299162, 12.054533472914915, 12.6724210290655, 13.138599687886162, 14.490032247894778, 17.7236817049896, 18.662584393779284, 19.346120605161712, 23.999160961181655, 25.451852237676196]
2 : [3.2203793481862757, 5.886623822609913, 5.905548429311276, 10.853747390755373, 10.901287969994906, 11.033338300567413, 12.490608609915988, 13.713193966056817, 14.184813087457245, 15.348674545613056, 16.52982807581595, 17.964540349069836]
3 : [3.2202534507267275, 5.880768708677623, 5.881908488029124, 10.709133420839656, 10.72605206596796, 10.77865621906485, 12.345960396958173, 12.76067604256242, 13.760281355732566, 14.018326042924281, 14.607718515213694, 16.59221415089495]
4 : [3.2202514669005953, 5.880494279982022, 5.880556917690844, 10.689894953485773, 10.698983267279337, 10.718281105531359, 12.323621310862592, 12.51513457485661, 13.530915636138472, 13.646352490651916, 14.013638403663302, 15.844042103498385]
5 : [3.2202514327806284, 5.8804771145376336, 5.880480023263916, 10.686851036556083, 10.694700891024613, 10.70089040705254, 12.319109793663676, 12.416240501188843, 13.409477889901574, 13.552567887226605, 13.749216311101861, 15.246452468777944]
6 : [3.220251432184983, 5.880473890511302, 5.880477843318856, 10.686161008730757, 10.69403287675776, 10.695820538590015, 12.318002400652476, 12.365108987042, 13.375722155059078, 13.488497017411383, 13.613981824811725, 14.813746365870585]
7 : [3.2202514321746234, 5.880473676484481, 5.880477751399412, 10.685969287196427, 10.693915072917033, 10.694415397844528, 12.317696219797469, 12.338610207915208, 13.365112608791692, 13.452949540230735, 13.531662728948632, 14.548118346419352]
8 : [3.2202514321744458, 5.880473664772765, 5.88047774601527, 10.685917731763432, 10.693892108137478, 10.694043457650642, 12.31760466508515, 12.32609806879308, 13.360878384968734, 13.435960447010164, 13.4805661994864, 14.398036525843478]
9 : [3.220251432174437, 5.880473664135102, 5.880477745705956, 10.685904872915286, 10.693887275314635, 10.693949433865447, 12.317575537838728, 12.320797131376127, 13.358877986624973, 13.428258104476354, 13.451525430098645, 14.31573512303175]
10 : [3.220251432174444, 5.880473664100338, 5.880477745688493, 10.685901792194166, 10.693886175010654, 10.693926453073034, 12.317565160998459, 12.318715755252375, 13.35789296633689, 13.424793977067829, 13.436298726064182, 14.270543600543458]
11 : [3.2202514321744413, 5.880473664098391, 5.880477745687529, 10.685901066925114, 10.693885916583827, 10.693920952759273, 12.317558831681099, 12.317939834410154, 13.357414569291278, 13.423228036049139, 13.428747528426282, 14.245626526451288]
12 : [3.2202514321744427, 5.880473664098333, 5.880477745687467, 10.685900897346425, 10.693885856844487, 10.693919648675356, 12.317545562935976, 12.317668127624847, 13.357188433281062, 13.422515875800016, 13.425120115595622, 14.23183024486441]
13 : [3.220251432174441, 5.880473664098306, 5.880477745687492, 10.685900857783858, 10.693885843152062, 10.69391934102008, 12.317512856495869, 12.317598521965907, 13.357083759102022, 13.4221894607061, 13.423410686002796, 14.224177161855884]
14 : [3.220251432174439, 5.880473664098329, 5.880477745687465, 10.685900848562277, 10.69388584001321, 10.693919268625418, 12.317488499058712, 12.31758642377891, 13.35703595917838, 13.422036742901984, 13.422616271732494, 14.219929393959701]
15 : [3.220251432174444, 5.880473664098114, 5.880477745687598, 10.685900846410757, 10.693885839291434, 10.693919251604605, 12.317478500777579, 12.317583475416692, 13.357014272936526, 13.421959799565162, 13.422254151415107, 14.217569203483857]
16 : [3.220251432174445, 5.880473664098207, 5.880477745687538, 10.6859008459084, 10.693885839121236, 10.693919247580212, 12.317474764946377, 12.317582539180991, 13.357004441751164, 13.421913701531789, 13.422095508425018, 14.216246027469076]
17 : [3.220251432174442, 5.880473664098266, 5.880477745687443, 10.685900845790577, 10.693885839085107, 10.693919246656865, 12.31747347014127, 12.317582237014216, 13.356999997137473, 13.421883056986692, 13.422030983627705, 14.215506312850298]
18 : [3.220251432174439, 5.8804736640982895, 5.880477745687561, 10.685900845763648, 10.693885839077474, 10.693919246442865, 12.317473023613074, 12.317582135502942, 13.35699802105203, 13.421865261476619, 13.422006114481439, 14.21510264851809]
19 : [3.220251432174442, 5.88047366409846, 5.88047774568751, 10.685900845757171, 10.69388583907567, 10.693919246391687, 12.317472863865504, 12.317582099745646, 13.356997066450417, 13.421855760785219, 13.421995923468444, 14.21486853350137]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.220251432174442
5.88047366409846
5.88047774568751
10.685900845757171
10.69388583907567
10.693919246391687
12.317472863865504
12.317582099745646
13.356997066450417
13.421855760785219
13.421995923468444
14.21486853350137
[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);
[ ]:
[ ]: