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 = 44384

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.2
  warnings.warn(f"A NumPy version >={np_minversion} and <{np_maxversion}"

0 : [5.921595123472257, 16.354970608244265, 19.605603006681456, 42.58472025783241, 47.46753828791688, 53.51420708010899, 67.47212884647678, 73.88926300011428, 82.3553172805869, 91.5356762617824, 101.12543397816916, 108.83008313968564]
1 : [3.2324570677857767, 6.017167188489162, 6.219423865566229, 11.719693778492864, 12.177834973217287, 12.681820756451707, 14.477850948322892, 15.596821107717835, 16.850696279850073, 17.2909634030003, 20.775093837442515, 24.10381246738526]
2 : [3.220424486542232, 5.884453035364446, 5.895380849604207, 10.805685750624702, 10.948315781132681, 11.003449902891706, 12.687851534241654, 13.149348652834448, 13.947130501938767, 14.165023198626022, 15.461066417762387, 18.936725139123396]
3 : [3.2202857963162868, 5.880622718645475, 5.881215482934005, 10.703348609960335, 10.751617499951246, 10.770399529855702, 12.441884158191746, 12.578846777704573, 13.57585743009314, 13.652071364164213, 14.353411562514411, 16.70122547292047]
4 : [3.2202839009361597, 5.880495053922039, 5.880522638350744, 10.688469332223855, 10.708276467471567, 10.714592150223202, 12.36355261938914, 12.409317987009931, 13.44283323351695, 13.50927095849069, 13.962225999867286, 15.353134485132866]
5 : [3.2202838706937658, 5.880486677010456, 5.8804906271290305, 10.686373706201401, 10.697364433358455, 10.699377121660776, 12.334461398700928, 12.3516752437036, 13.39858661986391, 13.452338700257691, 13.730574951609224, 14.724791911742741]
6 : [3.220283870202488, 5.880484832710461, 5.880490431252413, 10.686060068944412, 10.694716372916917, 10.695289566610578, 12.323561396325562, 12.330081936185064, 13.377256157257936, 13.43262202993094, 13.579314969035071, 14.453539725640526]
7 : [3.2202838701945957, 5.880484735958815, 5.880490423365975, 10.686008883820245, 10.694102643883191, 10.694266634976774, 12.319694706522874, 12.322036881067788, 13.366371150703559, 13.426058619405325, 13.496332840503298, 14.334088341352379]
8 : [3.220283870194476, 5.880484730924844, 5.880490423014433, 10.685999809329148, 10.693953956054113, 10.694032711798025, 12.318347615609646, 12.319158839328802, 13.361210677722958, 13.423642238601444, 13.455798099479795, 14.27684905060523]
9 : [3.2202838701944696, 5.880484730662485, 5.880490422997885, 10.685998084852145, 10.693915270703894, 10.69398366916357, 12.317869754879496, 12.318171918659553, 13.35885607322059, 13.422668881407537, 13.437079419761012, 14.247753580076342]
10 : [3.220283870194468, 5.880484730648793, 5.880490422997049, 10.685997738750748, 10.693905985117295, 10.693973001228917, 12.317678925906096, 12.317860021078241, 13.357800345505858, 13.422254543379097, 13.428645256011587, 14.232476427321721]
11 : [3.220283870194465, 5.880484730648054, 5.880490422997009, 10.685997666378451, 10.693903874831467, 10.693970583363068, 12.3175990911151, 12.317767809225293, 13.35733016755795, 13.422072584789936, 13.42488636536481, 14.224316938023254]
12 : [3.220283870194469, 5.880484730648058, 5.880490422997027, 10.685997650808497, 10.693903402353394, 10.693970028463722, 12.317569384942354, 12.317738636475191, 13.357121483676053, 13.421991026419231, 13.42321982953708, 14.219921948383178]
13 : [3.2202838701944687, 5.880484730648056, 5.880490422997017, 10.685997647391048, 10.693903296793275, 10.69396990068894, 12.317559010371685, 12.317728834552064, 13.357028932951536, 13.421953300556401, 13.422483119581504, 14.217543183255321]
14 : [3.2202838701944665, 5.88048473064803, 5.8804904229970045, 10.685997646630994, 10.693903273198794, 10.69396987123413, 12.317555459505215, 12.317725462390923, 13.356987910765485, 13.42193365118518, 13.422159993054358, 14.21625258334222]
15 : [3.220283870194473, 5.880484730648062, 5.880490422997067, 10.685997646460411, 10.693903267915188, 10.693969864438365, 12.317554250638917, 12.317724292221477, 13.356969720657116, 13.42191842378299, 13.422023396471154, 14.215551363116958]
16 : [3.2202838701944705, 5.880484730648085, 5.880490422997044, 10.685997646422127, 10.69390326672944, 10.693969862870839, 12.317553838951106, 12.317723884787483, 13.356961657690812, 13.421901139677642, 13.421973201609816, 14.215169689726736]
17 : [3.220283870194468, 5.880484730647119, 5.880490422997236, 10.685997646412423, 10.693903266462925, 10.693969862506277, 12.31755369883548, 12.317723741859023, 13.356958069422976, 13.421886908827524, 13.421957371827602, 14.214961824952377]
18 : [3.220283870194471, 5.880484730648109, 5.88049042299703, 10.685997646411437, 10.693903266401941, 10.693969862423314, 12.317553649601159, 12.317723692087744, 13.356956467079325, 13.421878479046136, 13.421951777311849, 14.21484408246443]
19 : [3.2202838701944705, 5.8804847306480275, 5.880490422997007, 10.685997646410607, 10.693903266388594, 10.693969862402811, 12.317553633238667, 12.317723675309914, 13.356955710387597, 13.421874671440822, 13.421949822111344, 14.214776309092017]

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

Eigenvalues
3.2202838701944705
5.8804847306480275
5.880490422997007
10.685997646410607
10.693903266388594
10.693969862402811
12.317553633238667
12.317723675309914
13.356955710387597
13.421874671440822
13.421949822111344
14.214776309092017

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