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 = 43192
[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 : [8.213129810314612, 17.799813247101675, 37.21982566100833, 44.019300333753456, 62.18233337178291, 65.53333261929534, 73.98716308337256, 81.08231277335221, 88.06001748322421, 91.25055694259324, 113.61980594341678, 122.21673397598462]
1 : [3.238097539179341, 5.998650630588756, 6.161390072348573, 11.485860551390664, 12.037309679960355, 14.413520504343095, 15.234814499827536, 16.21826937528625, 17.346726952882108, 19.892283926110174, 22.63881322440006, 24.807560961479695]
2 : [3.2205426885757555, 5.88487670889949, 5.888226095534204, 10.785471632144075, 10.847156783156317, 11.381986378644582, 13.343238320210492, 13.372706506482231, 14.112393285842757, 14.98863956329978, 15.344080016156944, 17.68832957714287]
3 : [3.2202954145601486, 5.88071344655438, 5.88083952659158, 10.707466526321864, 10.720069178385552, 10.94477331109386, 12.527973954871857, 13.075036550562103, 13.54069791583285, 14.030695748836376, 14.367494218070775, 15.906741583465786]
4 : [3.2202898941021965, 5.880502253291038, 5.880518097321824, 10.695114363079757, 10.699129958351238, 10.78775181170959, 12.354098370721175, 12.871053038558092, 13.445683487781059, 13.649403154578406, 14.21385996032969, 14.60766544701946]
5 : [3.2202897664212062, 5.88049337525466, 5.880495963116373, 10.69277636760732, 10.69499950969464, 10.721906317359633, 12.324202578884352, 12.640658945302796, 13.417157954198926, 13.509762215045422, 13.816759276283394, 14.260863114611452]
6 : [3.2202897638070227, 5.880492978117826, 5.8804945104449295, 10.69109718181177, 10.694161153653544, 10.698643594372449, 12.318883281766514, 12.460746604943322, 13.402394959830636, 13.454795875825525, 13.563293263168166, 14.232709515973163]
7 : [3.220289763759514, 5.880492958816416, 5.880494425037908, 10.687861627679673, 10.693985793779893, 10.694552792295557, 12.3178995711059, 12.37236442651547, 13.38745084267182, 13.43393941283865, 13.473271158046472, 14.223431358499429]
8 : [3.2202897637586956, 5.88049295783442, 5.880494420312156, 10.686464368101749, 10.693942498545207, 10.694085266048548, 12.317711055329683, 12.337392167598846, 13.373735427610427, 13.426640065251366, 13.441852510571445, 14.219267826275654]
9 : [3.220289763758685, 5.8804929577828275, 5.880494420057018, 10.686108698963897, 10.693925371312716, 10.694004893999336, 12.317673740214728, 12.32468016497018, 13.36519418847121, 13.423921247791194, 13.430320729828258, 14.217191023776492]
10 : [3.220289763758678, 5.88049295778006, 5.880494420043315, 10.686025216524408, 10.69391907743595, 10.693989757203505, 12.317666113412718, 12.320210206752254, 13.360871917039423, 13.422816823154236, 13.425650479964027, 14.216103594764382]
11 : [3.2202897637586827, 5.880492957779939, 5.880494420042612, 10.686006011050935, 10.693917414577, 10.69398657096079, 12.317664492531192, 12.318659203411382, 13.358849295201624, 13.422344037656002, 13.423640055250557, 14.21552028525456]
12 : [3.220289763758682, 5.880492957779889, 5.880494420042576, 10.686001612157382, 10.693917021615375, 10.693985859212644, 12.317664119924542, 12.318123653043754, 13.357930692968063, 13.422136540517268, 13.422748130146685, 14.215203525656346]
13 : [3.220289763758685, 5.8804929577798815, 5.880494420042572, 10.686000605556455, 10.693916931076009, 10.693985697532778, 12.317664014899652, 12.317939023421149, 13.357518488365468, 13.422044242302649, 13.422347230645611, 14.215030370910645]
14 : [3.220289763758686, 5.8804929577799125, 5.880494420042557, 10.686000375002681, 10.693916910322711, 10.693985660616315, 12.317663977204242, 12.317875405900041, 13.357334534424274, 13.422002971318507, 13.422165839902437, 14.214935412065216]
15 : [3.2202897637586885, 5.880492957779939, 5.880494420042595, 10.686000322154042, 10.693916905559107, 10.69398565216939, 12.31766396300961, 12.317853449260866, 13.357252487631786, 13.421984348148433, 13.422083626493652, 14.214883204928634]
16 : [3.2202897637586854, 5.880492957779734, 5.880494420042597, 10.686000310003367, 10.693916904464501, 10.693985650229452, 12.31766395791089, 12.317845854377264, 13.357215848348714, 13.421975901398268, 13.422046255573406, 14.214854438756443]
17 : [3.2202897637586854, 5.880492957780129, 5.880494420042567, 10.686000307197357, 10.693916904212244, 10.693985649781714, 12.31766395612913, 12.317843214006668, 13.357199383179092, 13.421972043557554, 13.422029145822696, 14.214838415519848]
18 : [3.2202897637586845, 5.880492957779859, 5.880494420042721, 10.686000306534602, 10.69391690415065, 10.693985649677249, 12.31766395551425, 12.317842251539174, 13.357191712457794, 13.42197016232848, 13.422020927264638, 14.214829625609891]
19 : [3.220289763758685, 5.880492957780034, 5.880494420042466, 10.686000306404903, 10.693916904139869, 10.693985649656259, 12.317663955310232, 12.317841965205567, 13.357188590238186, 13.421969407533032, 13.422017545638544, 14.214824795594225]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.220289763758685
5.880492957780034
5.880494420042466
10.686000306404903
10.693916904139869
10.693985649656259
12.317663955310232
12.317841965205567
13.357188590238186
13.421969407533032
13.422017545638544
14.214824795594225
[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);
[ ]:
[ ]: