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)
/usr/lib/python3/dist-packages/scipy/__init__.py:146: UserWarning: A NumPy version >=1.17.3 and <1.25.0 is required for this version of SciPy (detected version 1.26.4
warnings.warn(f"A NumPy version >={np_minversion} and <{np_maxversion}"
0 : [8.213129806057713, 17.799813243338956, 37.219825641537795, 44.01930033044001, 62.18233336959265, 65.53333261977797, 73.98716308312365, 81.08231276750318, 88.06001747696237, 91.25055696454325, 113.61980593708506, 122.21673394805836]
1 : [3.238097539149225, 5.998650630544873, 6.161390071999765, 11.485860551447896, 12.037309680313156, 14.41352050412174, 15.234814500485424, 16.218269374433127, 17.346726951842793, 19.892283930144583, 22.63881321383744, 24.807560970367817]
2 : [3.2205426885753967, 5.884876708899002, 5.888226095526375, 10.785471632184338, 10.847156783241985, 11.381986378991414, 13.343238318750997, 13.372706506976911, 14.112393286178811, 14.988639564796163, 15.344080018300994, 17.688329579825503]
3 : [3.2202954145601432, 5.880713446554435, 5.880839526591644, 10.707466526331345, 10.720069178412029, 10.944773311402418, 12.527973954678302, 13.075036551362594, 13.540697916078342, 14.030695749435447, 14.367494218554844, 15.906741593107071]
4 : [3.220289894102196, 5.880502253291033, 5.880518097321875, 10.695114363082206, 10.699129958358926, 10.78775181198475, 12.354098370687973, 12.87105304011837, 13.44568348780674, 13.649403154985748, 14.213859960708305, 14.607665455646535]
5 : [3.220289766421212, 5.880493375254697, 5.880495963116384, 10.692776367607893, 10.694999509695974, 10.721906317531444, 12.324202578877406, 12.640658947229953, 13.417157954170529, 13.509762215271307, 13.816759279719056, 14.260863115474386]
6 : [3.220289763807021, 5.880492978117788, 5.880494510444957, 10.691097181821734, 10.69416115365383, 10.698643594428463, 12.318883281764892, 12.460746606131526, 13.402394959818235, 13.454795875961702, 13.56329326455004, 14.232709516228564]
7 : [3.220289763759517, 5.880492958816453, 5.880494425037921, 10.687861627691762, 10.69398579378006, 10.694552792302378, 12.31789957110548, 12.37236442703147, 13.387450842724002, 13.433939412910702, 13.473271158548004, 14.223431358615931]
8 : [3.220289763758697, 5.88049295783441, 5.880494420312183, 10.68646436810567, 10.693942498545345, 10.694085266049742, 12.31771105532959, 12.337392167807614, 13.37373542768452, 13.426640065283378, 13.441852510757519, 14.219267826336777]
9 : [3.220289763758683, 5.880492957782819, 5.880494420056998, 10.686108698958616, 10.69392537131257, 10.694004893998295, 12.317673740214369, 12.32468016472636, 13.365194188292238, 13.42392124783351, 13.430320729604144, 14.217191023820366]
10 : [3.220289763758678, 5.880492957780062, 5.880494420043336, 10.68602521651832, 10.693919077435579, 10.693989757202383, 12.317666113412589, 12.320210206243106, 13.36087191647814, 13.422816823165268, 13.425650479488423, 14.216103594724629]
11 : [3.220289763758686, 5.880492957779894, 5.880494420042619, 10.686006011056314, 10.693917414577276, 10.693986570961862, 12.317664492529202, 12.318659204070991, 13.358849296088803, 13.422344037493998, 13.423640056489242, 14.21552028537592]
12 : [3.220289763758682, 5.880492957779895, 5.880494420042561, 10.6860016121426, 10.693917021613988, 10.693985859209993, 12.317664119924187, 12.318123649810389, 13.357930686296474, 13.42213654047664, 13.42274812708441, 14.215203524568722]
13 : [3.2202897637586854, 5.880492957779926, 5.8804944200425675, 10.686000605564292, 10.693916931074813, 10.693985697533886, 12.317664014900116, 12.317939021098224, 13.357518483149597, 13.422044244556526, 13.422347245586174, 14.215030373383465]
14 : [3.2202897637586854, 5.880492957779915, 5.88049442004258, 10.68600037498212, 10.693916910321748, 10.69398566061076, 12.317663977193185, 12.317875385828804, 13.357334441951656, 13.422002974282956, 13.422165841794655, 14.21493536475173]
15 : [3.220289763758689, 5.880492957780002, 5.880494420042207, 10.686000322145082, 10.693916905558963, 10.6939856521666, 12.317663963001284, 12.317853436771408, 13.357252411301962, 13.421984350878834, 13.422083624470973, 14.214883162232478]
16 : [3.2202897637586867, 5.880492957779895, 5.880494420042461, 10.686000309999288, 10.693916904463967, 10.69398565022856, 12.317663957899498, 12.317845848468126, 13.357215804941431, 13.421975903543212, 13.422046246842399, 14.214854381145646]
17 : [3.2202897637586867, 5.880492957780343, 5.880494420042588, 10.68600030720389, 10.693916904212518, 10.693985649783382, 12.317663956133881, 12.31784321996404, 13.357199441615363, 13.421972059639721, 13.42202922762824, 14.214838505127304]
18 : [3.2202897637586845, 5.880492957779826, 5.880494420042571, 10.686000306561999, 10.693916904154477, 10.693985649680613, 12.317663955525218, 12.31784231274086, 13.35719215196121, 13.421970318945998, 13.422021485622963, 14.214829760210774]
19 : [3.220289763758688, 5.880492957780051, 5.880494420042437, 10.686000306414883, 10.693916904140877, 10.693985649657046, 12.31766395531741, 12.317841993632726, 13.357188856384836, 13.421969451974041, 13.42201798844545, 14.214824965302794]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.220289763758688
5.880492957780051
5.880494420042437
10.686000306414883
10.693916904140877
10.693985649657046
12.31766395531741
12.317841993632726
13.357188856384836
13.421969451974041
13.42201798844545
14.214824965302794
[7]:
gfu = GridFunction(fes, multidim=len(evecs))
for i in range(len(evecs)):
gfu.vecs[i].data = evecs[i]
Draw (Norm(gfu), mesh, order=2, min=0, max=2);
[ ]:
[ ]: