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.22909326272806, 22.636447500631803, 38.75912680375153, 56.956883730514164, 64.68904674603476, 71.61510643862434, 76.893112927557, 86.36950976780795, 90.92730267472626, 99.90890224834206, 108.07709801704047, 123.53078334767719]
1 : [3.231012675422761, 6.086234263790819, 6.406899224055859, 12.054533473094823, 12.67242102912832, 13.138599690189622, 14.490032246485464, 17.723681705664507, 18.662584395827814, 19.346120603700296, 23.999160961477475, 25.45185223799197]
2 : [3.220379348186411, 5.8866238226086445, 5.905548429297472, 10.85374739069025, 10.901287970136169, 11.033338300674648, 12.490608609918196, 13.713193966612167, 14.184813087593689, 15.348674546792484, 16.529828074312725, 17.96454034937476]
3 : [3.220253450726736, 5.880768708677629, 5.881908488028274, 10.709133420843397, 10.726052065962977, 10.778656219130982, 12.34596039697391, 12.76067604279276, 13.760281355775717, 14.018326043114628, 14.607718515366853, 16.59221414991472]
4 : [3.2202514669005953, 5.88049427998199, 5.880556917690783, 10.689894953486473, 10.698983267278935, 10.718281105556004, 12.323621310870177, 12.515134574948759, 13.530915636210617, 13.646352490574202, 14.013638403963862, 15.844042101998179]
5 : [3.2202514327806315, 5.880477114537648, 5.880480023263929, 10.686851036556432, 10.694700891024912, 10.700890407059273, 12.31910979366708, 12.41624050121676, 13.409477889922728, 13.552567887176908, 13.749216311270198, 15.246452467294823]
6 : [3.2202514321849813, 5.880473890511268, 5.880477843318853, 10.686161008730938, 10.694032876757856, 10.695820538591446, 12.318002400653882, 12.365108987042445, 13.375722155069228, 13.488497017386807, 13.613981824847672, 14.813746364763482]
7 : [3.2202514321746225, 5.880473676484473, 5.880477751399417, 10.68596928719647, 10.693915072917063, 10.694415397844784, 12.317696219798007, 12.338610207909346, 13.365112608796792, 13.452949540218249, 13.531662728919565, 14.54811834572785]
8 : [3.220251432174445, 5.880473664772754, 5.88047774601528, 10.685917731763416, 10.693892108137469, 10.694043457650784, 12.317604665085359, 12.326098068779144, 13.360878384970665, 13.435960447014432, 13.480566199315502, 14.39803652519751]
9 : [3.2202514321744404, 5.880473664135108, 5.88047774570595, 10.685904872914357, 10.693887275312889, 10.693949433865198, 12.317575537840627, 12.320797131169963, 13.358877986650509, 13.428258104719617, 13.45152542094424, 14.315735104930638]
10 : [3.2202514321744378, 5.880473664100327, 5.8804777456885, 10.685901792193624, 10.693886175009977, 10.693926453072802, 12.317565160998518, 12.318715753653372, 13.357892966632805, 13.424793975975147, 13.43629869333067, 14.2705435457338]
11 : [3.2202514321744404, 5.880473664098432, 5.8804777456875295, 10.685901066925112, 10.69388591658357, 10.693920952760008, 12.317558831679758, 12.317939833433579, 13.357414569606497, 13.423228034983753, 13.428747499055484, 14.245626480388964]
12 : [3.2202514321744413, 5.880473664098308, 5.8804777456874895, 10.68590089734632, 10.693885856844705, 10.693919648676124, 12.317545562864954, 12.317668127012109, 13.357188433639278, 13.42251587644623, 13.425120084202169, 14.231830179212022]
13 : [3.2202514321744427, 5.88047366409828, 5.880477745687438, 10.685900857785045, 10.693885843151763, 10.693919341026662, 12.317512856021644, 12.317598521356036, 13.357083761564395, 13.422189471755654, 13.423410552878961, 14.224177003267274]
14 : [3.2202514321744373, 5.8804736640983135, 5.880477745687476, 10.685900848563069, 10.693885840014035, 10.693919268627818, 12.317488499090619, 12.317586423603554, 13.357035961684623, 13.422036746198499, 13.422616269662484, 14.21992967521665]
15 : [3.2202514321744413, 5.8804736640984006, 5.8804777456875525, 10.685900846411023, 10.693885839291946, 10.693919251608746, 12.31747850407984, 12.31758347613697, 13.357014272869206, 13.42195981973702, 13.42225428299073, 14.21757048542279]
16 : [3.2202514321744427, 5.880473664098261, 5.880477745687563, 10.68590084590883, 10.693885839125066, 10.693919247607093, 12.317474828481311, 12.3175825593431, 13.357004467050531, 13.42191394677502, 13.422097042634293, 14.216258757362445]
17 : [3.2202514321744435, 5.880473664098342, 5.880477745687482, 10.685900845791176, 10.693885839086537, 10.693919246667692, 12.317473512275216, 12.317582249416791, 13.357000037428264, 13.4218837439372, 13.422032654246767, 14.215526772936776]
18 : [3.220251432174441, 5.880473664098205, 5.880477745687407, 10.685900845763685, 10.693885839077517, 10.693919246445967, 12.317473040163728, 12.317582140216805, 13.356998040345367, 13.421865728823438, 13.422006458781166, 14.215115240932263]
19 : [3.2202514321744418, 5.880473664098221, 5.88047774568738, 10.685900845757372, 10.693885839075815, 10.693919246393081, 12.31747287318861, 12.317582102235683, 13.356997120679289, 13.421856296565002, 13.421996295943519, 14.214883947332961]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.2202514321744418
5.880473664098221
5.88047774568738
10.685900845757372
10.693885839075815
10.693919246393081
12.31747287318861
12.317582102235683
13.356997120679289
13.421856296565002
13.421996295943519
14.214883947332961
[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);
[ ]:
[ ]: