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.229093261714403, 22.63644749318042, 38.75912681180241, 56.95688372013053, 64.68904674753378, 71.61510644245911, 76.89311294002941, 86.36950975553901, 90.92730267786348, 99.90890224151633, 108.07709800778228, 123.53078334149133]
1 : [3.2310126754181994, 6.086234263687496, 6.406899224091792, 12.054533473686952, 12.672421029208772, 13.138599687852365, 14.490032247319215, 17.723681705952572, 18.66258439748134, 19.34612060564295, 23.9991609590725, 25.451852238015153]
2 : [3.2203793481863814, 5.886623822606728, 5.905548429299844, 10.853747390758606, 10.901287970125566, 11.033338300549701, 12.490608609824438, 13.713193967012526, 14.184813087392756, 15.348674547156051, 16.5298280759016, 17.964540349078803]
3 : [3.2202534507267324, 5.880768708677578, 5.881908488028324, 10.709133420837182, 10.726052065993828, 10.778656219086399, 12.34596039694463, 12.76067604275376, 13.760281355656769, 14.018326043437497, 14.607718515859561, 16.59221415029018]
4 : [3.2202514669005966, 5.880494279981975, 5.880556917690787, 10.689894953486107, 10.69898326728354, 10.718281105539235, 12.3236213108603, 12.51513457489917, 13.530915636186986, 13.64635249070268, 14.013638404084066, 15.844042102588517]
5 : [3.2202514327806284, 5.880477114537642, 5.880480023263921, 10.686851036555822, 10.694700891025521, 10.700890407053775, 12.319109793663554, 12.416240501186929, 13.409477889905824, 13.552567887266529, 13.749216311260872, 15.246452467873231]
6 : [3.220251432184985, 5.88047389051127, 5.880477843318863, 10.686161008730522, 10.694032876757944, 10.695820538589931, 12.318002400652595, 12.365108987029457, 13.375722155049521, 13.488497017451312, 13.613981824824304, 14.813746365209795]
7 : [3.2202514321746194, 5.8804736764844865, 5.880477751399419, 10.685969287196333, 10.693915072917054, 10.694415397844372, 12.317696219797526, 12.338610207905043, 13.365112608781715, 13.452949540256226, 13.531662728919692, 14.548118346023294]
8 : [3.220251432174442, 5.8804736647727625, 5.880477746015275, 10.685917731763388, 10.693892108137483, 10.694043457650613, 12.31760466508517, 12.326098068782427, 13.360878384961916, 13.43596044702398, 13.480566199356437, 14.398036525467688]
9 : [3.2202514321744387, 5.880473664135103, 5.880477745705924, 10.685904872915183, 10.6938872753128, 10.69394943386701, 12.317575537838808, 12.320797131569897, 13.358877986637934, 13.428258103102229, 13.451525421749281, 14.31573511847231]
10 : [3.2202514321744387, 5.8804736641003235, 5.880477745688508, 10.685901792194066, 10.6938861750099, 10.693926453073903, 12.31756516099861, 12.318715754440102, 13.357892966401675, 13.424793976227775, 13.436298705389925, 14.27054356984254]
11 : [3.2202514321744427, 5.8804736640984, 5.880477745687531, 10.685901066924878, 10.693885916583623, 10.693920952758935, 12.317558831674262, 12.317939833654432, 13.357414569261891, 13.423228035588174, 13.428747506893888, 14.245626481849511]
12 : [3.2202514321744413, 5.8804736640983135, 5.880477745687468, 10.685900897346238, 10.693885856844746, 10.693919648675427, 12.317545562866185, 12.317668127142397, 13.35718843333593, 13.422515877039396, 13.42512009705971, 14.23183018815208]
13 : [3.2202514321744418, 5.880473664098327, 5.880477745687495, 10.685900857789068, 10.693885843152731, 10.693919341041157, 12.317512863482362, 12.317598526549562, 13.35708376987382, 13.422189471968322, 13.423410997313438, 14.224178903257007]
14 : [3.2202514321744418, 5.8804736640983, 5.880477745687461, 10.685900848563929, 10.693885840013596, 10.693919268634385, 12.31748850193651, 12.317586424390697, 13.357035967273292, 13.422036762634443, 13.422616232099486, 14.219930217351619]
15 : [3.220251432174442, 5.880473664098286, 5.880477745687441, 10.685900846411293, 10.693885839291683, 10.693919251609591, 12.317478502156169, 12.317583475450512, 13.357014277643062, 13.421959838864616, 13.422253990831637, 14.217569680332081]
16 : [3.220251432174443, 5.880473664098444, 5.880477745687482, 10.685900845908437, 10.693885839124592, 10.693919247603763, 12.31747482013691, 12.31758255680931, 13.357004469668809, 13.421913911848836, 13.422096742194954, 14.216256617678267]
17 : [3.2202514321744404, 5.880473664098235, 5.880477745687429, 10.685900845790172, 10.693885839086388, 10.693919246653975, 12.317473485447065, 12.317582243315513, 13.356999986013696, 13.421883010315524, 13.422032629963255, 14.215509599166571]
18 : [3.2202514321744427, 5.880473664098224, 5.88047774568731, 10.68590084576344, 10.693885839077831, 10.693919246440913, 12.317473030746003, 12.317582138572753, 13.356997983343492, 13.42186531701566, 13.422007257029442, 14.21510519820315]
19 : [3.2202514321744435, 5.880473664098209, 5.88047774568747, 10.68590084575725, 10.693885839075804, 10.693919246392865, 12.317472874649583, 12.317582102739008, 13.35699710677541, 13.421856393967724, 13.421996692273082, 14.214886297040666]
[6]:
print ("Eigenvalues")
for lam in evals:
    print (lam)
Eigenvalues
3.2202514321744435
5.880473664098209
5.88047774568747
10.68590084575725
10.693885839075804
10.693919246392865
12.317472874649583
12.317582102739008
13.35699710677541
13.421856393967724
13.421996692273082
14.214886297040666
[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);
[ ]:

[ ]: