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

0 : [5.921595127472971, 16.354970609601462, 19.605603000278688, 42.58472025357642, 47.4675382921202, 53.51420708192481, 67.47212885821452, 73.88926300148832, 82.3553172556532, 91.53567627124974, 101.12543400516255, 108.83008317845227]
1 : [3.232457067800617, 6.017167188445006, 6.219423865699035, 11.71969377985415, 12.177834973417813, 12.681820756102205, 14.477850948676139, 15.596821107487202, 16.850696283117973, 17.290963401950982, 20.775093836284224, 24.103812463429616]
2 : [3.220424486542368, 5.884453035363423, 5.895380849615198, 10.805685750899539, 10.94831578092535, 11.00344990290829, 12.68785153393187, 13.149348654046499, 13.94713049997284, 14.165023197709555, 15.46106641814932, 18.936725138466013]
3 : [3.2202857963162868, 5.880622718645461, 5.881215482934627, 10.703348610002125, 10.75161749986115, 10.770399529913274, 12.441884157962702, 12.578846778088776, 13.575857429614196, 13.65207136345768, 14.353411562568903, 16.70122547565227]
4 : [3.220283900936159, 5.880495053922065, 5.88052263835076, 10.688469332231065, 10.708276467451062, 10.714592150244147, 12.363552619290276, 12.409317987187931, 13.44283323335727, 13.509270958221487, 13.962225999862392, 15.353134488509168]
5 : [3.220283870693764, 5.880486677010455, 5.88049062712904, 10.686373706202739, 10.697364433354574, 10.699377121668544, 12.33446139866717, 12.351675243796294, 13.398586619781945, 13.452338700183684, 13.730574951812375, 14.724791913872028]
6 : [3.22028387020248, 5.8804848327104775, 5.88049043125242, 10.686060068944673, 10.69471637291609, 10.69528956661312, 12.323561396315283, 12.330081936228705, 13.377256157219298, 13.432622029910245, 13.57931496928683, 14.45353972673245]
7 : [3.2202838701945966, 5.880484735958813, 5.880490423365995, 10.68600888382032, 10.694102643883063, 10.694266634977517, 12.319694706519774, 12.322036881085822, 13.366371150686746, 13.426058619398201, 13.496332840674674, 14.33408834191198]
8 : [3.2202838701944714, 5.880484730924877, 5.880490423014426, 10.685999809329171, 10.693953956054024, 10.69403271179812, 12.318347615610174, 12.319158839332061, 13.361210677717057, 13.423642238598887, 13.455798099553965, 14.276849050847884]
9 : [3.220283870194469, 5.880484730662487, 5.8804904229978625, 10.685998084852157, 10.69391527070272, 10.693983669163355, 12.31786975486964, 12.318171918634357, 13.35885607323577, 13.422668881397883, 13.437079419103755, 14.24775357874087]
10 : [3.2202838701944674, 5.880484730648779, 5.880490422997055, 10.685997738750745, 10.693905985118082, 10.693973001228947, 12.317678925929073, 12.31786002109632, 13.357800345512146, 13.422254543374812, 13.428645256880394, 14.232476430696767]
11 : [3.22028387019447, 5.880484730648088, 5.880490422997009, 10.685997666378482, 10.693903874830465, 10.693970583363784, 12.317599091053154, 12.31776780925397, 13.357330167652105, 13.42207258490099, 13.424886359116194, 14.2243169341841]
12 : [3.220283870194467, 5.88048473064805, 5.880490422997022, 10.685997650808671, 10.693903402361505, 10.693970028465328, 12.317569385530337, 12.317738636539527, 13.357121486344832, 13.421991026833636, 13.423219895399743, 14.219922036431779]
13 : [3.220283870194469, 5.880484730648021, 5.880490422997001, 10.685997647391126, 10.693903296796817, 10.693969900689169, 12.317559010755657, 12.317728834697473, 13.357028935039125, 13.421953300324692, 13.422483178671015, 14.21754326185706]
14 : [3.2202838701944674, 5.880484730648008, 5.880490422997021, 10.685997646631032, 10.693903273193964, 10.693969871233449, 12.31755545861548, 12.317725462497584, 13.356987908624923, 13.42193365006137, 13.422159789139032, 14.216252301451698]
15 : [3.220283870194473, 5.880484730647985, 5.880490422997159, 10.685997646460413, 10.6939032679125, 10.693969864436689, 12.317554249777071, 12.317724292325678, 13.356969712758232, 13.42191840140823, 13.422023197253099, 14.215550825355113]
16 : [3.2202838701944683, 5.880484730647817, 5.880490422996893, 10.685997646421963, 10.693903266728379, 10.69396986286691, 12.317553838568264, 12.31772388354239, 13.356961644545835, 13.421901122327737, 13.421973137667504, 14.215168803828073]
17 : [3.2202838701944674, 5.880484730648136, 5.880490422996978, 10.685997646413334, 10.693903266461836, 10.69396986250417, 12.317553697828167, 12.317723741918531, 13.356958037518046, 13.42188689600829, 13.421957484188265, 14.21495847492692]
18 : [3.2202838701944683, 5.8804847306482895, 5.880490422997081, 10.685997646411323, 10.693903266402438, 10.6939698624219, 12.317553650415034, 12.31772369260791, 13.356956458067256, 13.42187900123067, 13.421952137454493, 14.214845160580287]
19 : [3.2202838701944696, 5.88048473064765, 5.880490422996957, 10.685997646410906, 10.693903266389034, 10.693969862403591, 12.317553634570947, 12.317723675536895, 13.356955762607456, 13.42187526129148, 13.421950055617765, 14.214784930642123]

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

Eigenvalues
3.2202838701944696
5.88048473064765
5.880490422996957
10.685997646410906
10.693903266389034
10.693969862403591
12.317553634570947
12.317723675536895
13.356955762607456
13.42187526129148
13.421950055617765
14.214784930642123

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