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.

from ngsolve import *
from ngsolve.webgui import Draw
from netgen.occ import *

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);

mesh = Mesh(OCCGeometry(fichera).GenerateMesh(maxh=0.4))
mesh.RefineHP(levels=2, factor=0.2)
Draw (mesh);

# 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

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.213129811011335, 17.799813246615486, 37.219825663351244, 44.01930032549319, 62.18233337285682, 65.53333261751895, 73.98716309631672, 81.08231275864784, 88.06001747674402, 91.25055695789207, 113.61980590015652, 122.21673395577942]
1 : [3.238097539161607, 5.998650630652925, 6.161390072384866, 11.485860551215653, 12.037309679416365, 14.413520503235686, 15.234814499988252, 16.2182693714381, 17.346726952352434, 19.892283926354168, 22.638813221713743, 24.807560947138814]
2 : [3.220542688575768, 5.884876708905918, 5.888226095534121, 10.78547163217801, 10.847156783091094, 11.381986378070312, 13.343238319924394, 13.372706506473122, 14.1123932850333, 14.98863956410815, 15.34408001640501, 17.688329581656248]
3 : [3.2202954145601637, 5.880713446554873, 5.880839526591903, 10.707466526334926, 10.72006917838689, 10.944773311149744, 12.527973954765736, 13.075036551741949, 13.540697915835546, 14.03069574860852, 14.367494218452597, 15.906741600122379]
4 : [3.2202898941021942, 5.880502253291083, 5.880518097321926, 10.69511436308063, 10.699129958355185, 10.787751812044556, 12.35409837070336, 12.871053041101332, 13.445683487740236, 13.649403154482927, 14.21385996188592, 14.6076654612774]
5 : [3.2202897664212125, 5.880493375254671, 5.880495963116384, 10.692776367606932, 10.694999509695593, 10.721906317620975, 12.324202578880916, 12.640658948554771, 13.417157954140654, 13.50976221505814, 13.816759282887416, 14.260863115546334]
6 : [3.220289763807022, 5.880492978117802, 5.88049451044493, 10.691097181828331, 10.694161153653688, 10.698643594462188, 12.318883281765832, 12.460746606980234, 13.40239495982329, 13.454795875902684, 13.56329326575988, 14.232709516208413]
7 : [3.220289763759513, 5.880492958816448, 5.880494425037923, 10.687861627699991, 10.69398579377996, 10.694552792306624, 12.317899571105713, 12.3723644274072, 13.387450842788077, 13.43393941290556, 13.473271158964511, 14.22343135860298]
8 : [3.220289763758694, 5.880492957834404, 5.880494420312163, 10.686464368108027, 10.693942498545283, 10.694085266050525, 12.317711055329621, 12.337392167943955, 13.373735427751464, 13.426640065283255, 13.441852510889378, 14.219267826329984]
9 : [3.220289763758684, 5.880492957782824, 5.880494420057021, 10.6861086989603, 10.693925371312693, 10.694004893998963, 12.317673740214804, 12.324680164956348, 13.365194188494293, 13.423921247803927, 13.430320729640801, 14.217191023838426]
10 : [3.2202897637586823, 5.880492957780036, 5.880494420043316, 10.686025216516871, 10.69391907743562, 10.693989757202164, 12.317666113412779, 12.320210206194142, 13.360871916368243, 13.422816823142306, 13.425650479320467, 14.216103594667338]
11 : [3.220289763758686, 5.880492957779917, 5.880494420042604, 10.68600601102599, 10.693917414576534, 10.693986570954896, 12.317664492533753, 12.318659200253846, 13.358849289023421, 13.422344037869118, 13.423640051856154, 14.215520282336078]
12 : [3.220289763758684, 5.880492957779903, 5.880494420042574, 10.686001612140867, 10.693917021614645, 10.693985859208603, 12.31766411992484, 12.318123648747553, 13.357930681858582, 13.42213653981108, 13.422748126480556, 14.21520352058704]
13 : [3.2202897637586885, 5.88049295777992, 5.880494420042557, 10.68600060556362, 10.69391693107647, 10.693985697530403, 12.317664014898892, 12.31793901909235, 13.357518463666642, 13.422044242716245, 13.4223472483057, 14.215030345945383]
14 : [3.220289763758688, 5.880492957779895, 5.88049442004254, 10.686000374956436, 10.693916910317306, 10.69398566060478, 12.317663977185587, 12.317875367628229, 13.35733437840781, 13.42200295986272, 13.42216582581799, 14.214935340602322]
15 : [3.220289763758688, 5.880492957779963, 5.880494420042638, 10.68600032212786, 10.693916905556677, 10.693985652163375, 12.317663962994972, 12.317853421396947, 13.357252345652647, 13.421984335828677, 13.42208359676174, 14.214883138190078]
16 : [3.2202897637586845, 5.880492957779913, 5.880494420042598, 10.686000309994208, 10.693916904463226, 10.693985650227889, 12.317663957895537, 12.317845840705457, 13.357215763011023, 13.421975890962896, 13.422046227831483, 14.21485436380243]
17 : [3.2202897637586867, 5.880492957780057, 5.880494420042594, 10.686000307200423, 10.693916904211937, 10.693985649782599, 12.317663956127937, 12.31784321447429, 13.357199402842385, 13.421972049953053, 13.422029196194313, 14.214838470664322]
18 : [3.2202897637586823, 5.880492957779927, 5.880494420042497, 10.686000306520048, 10.693916904152713, 10.693985649673136, 12.31766395545589, 12.317842225678033, 13.357191225142532, 13.421970265534913, 13.422020673176918, 14.214828805017413]
19 : [3.220289763758683, 5.880492957779885, 5.880494420042631, 10.686000306400198, 10.693916904140842, 10.693985649654861, 12.317663955277824, 12.317841950952001, 13.357188228407987, 13.421969474652757, 13.422017385323521, 14.214824068580635]

print ("Eigenvalues")
for lam in evals:
    print (lam)

Eigenvalues
3.220289763758683
5.880492957779885
5.880494420042631
10.686000306400198
10.693916904140842
10.693985649654861
12.317663955277824
12.317841950952001
13.357188228407987
13.421969474652757
13.422017385323521
14.214824068580635

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);