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.213129808414807, 17.799813242788595, 37.21982567795192, 44.01930033081317, 62.182333374489374, 65.5333326216621, 73.98716307679119, 81.0823127532986, 88.06001748837035, 91.25055696413106, 113.61980589119533, 122.21673394069046]
1 : [3.238097539150205, 5.998650630542939, 6.161390072576845, 11.485860551086294, 12.037309680921126, 14.413520500950355, 15.234814500074144, 16.21826936975012, 17.34672695159324, 19.892283927857203, 22.638813219947643, 24.807560939532024]
2 : [3.2205426885755712, 5.8848767089007366, 5.888226095541294, 10.785471632087965, 10.847156783288954, 11.381986377533268, 13.343238319334338, 13.372706505829678, 14.112393284633242, 14.988639564826759, 15.344080016855369, 17.68832957113468]
3 : [3.220295414560143, 5.880713446554866, 5.880839526591378, 10.707466526319546, 10.720069178409286, 10.944773310500405, 12.527973954794785, 13.075036549783984, 13.540697915791776, 14.030695748905012, 14.367494218469172, 15.906741577737597]
4 : [3.220289894102192, 5.880502253291053, 5.880518097321773, 10.695114363079819, 10.699129958357274, 10.78775181139401, 12.354098370710007, 12.87105303761511, 13.445683487776474, 13.649403154676099, 14.213859959968095, 14.60766544363044]
5 : [3.22028976642121, 5.880493375254655, 5.88049596311637, 10.692776367606335, 10.694999509695997, 10.721906317218686, 12.324202578881149, 12.640658944371582, 13.417157954197371, 13.509762215125749, 13.81675927476555, 14.260863114555738]
6 : [3.220289763807026, 5.880492978117834, 5.880494510444941, 10.691097181802727, 10.694161153653752, 10.69864359433527, 12.318883281765356, 12.460746604427628, 13.402394959831765, 13.454795875877398, 13.563293262597744, 14.232709515993362]
7 : [3.2202897637595136, 5.880492958816452, 5.880494425037914, 10.687861627671095, 10.693985793779925, 10.694552792291775, 12.317899571105423, 12.372364426305978, 13.387450842668935, 13.433939412861848, 13.473271157827288, 14.223431358514793]
8 : [3.2202897637586982, 5.880492957834412, 5.880494420312184, 10.686464368100053, 10.693942498545153, 10.694085266048125, 12.31771105532951, 12.337392167536851, 13.373735427609649, 13.426640065256654, 13.441852510509896, 14.219267826288071]
9 : [3.220289763758686, 5.880492957782827, 5.880494420057023, 10.686108698963489, 10.69392537131251, 10.694004893999255, 12.317673740214309, 12.324680164866546, 13.365194188430426, 13.423921247822635, 13.430320729894813, 14.217191023850503]
10 : [3.2202897637586854, 5.880492957780056, 5.880494420043308, 10.68602521652345, 10.693919077435847, 10.69398975720334, 12.31766611341253, 12.320210206598114, 13.360871916843035, 13.422816823173653, 13.425650479973779, 14.216103594779277]
11 : [3.2202897637586836, 5.880492957779918, 5.880494420042636, 10.686006011048574, 10.693917414576472, 10.693986570959073, 12.317664492530897, 12.31865920144805, 13.358849290996655, 13.42234403839467, 13.423640057077026, 14.215520283608797]
12 : [3.2202897637586885, 5.880492957779901, 5.8804944200425755, 10.686001612153287, 10.693917021615105, 10.693985859210748, 12.31766411992412, 12.318123650056194, 13.357930683862795, 13.422136540447829, 13.422748131775291, 14.215203521354521]
13 : [3.220289763758683, 5.880492957779892, 5.880494420042584, 10.686000605645942, 10.69391693107802, 10.693985697546928, 12.317664014901627, 12.31793903883704, 13.357518527917707, 13.42204424113298, 13.422347300982894, 14.21503038440964]
14 : [3.2202897637586854, 5.8804929577799, 5.880494420042546, 10.686000375046454, 10.693916910326474, 10.69398566062122, 12.31766397719333, 12.31787541050535, 13.357334517119517, 13.422002966004527, 13.42216590742287, 14.21493537344909]
15 : [3.2202897637586885, 5.8804929577798974, 5.8804944200423375, 10.686000322167038, 10.69391690556052, 10.69398565217043, 12.317663963002484, 12.317853449086092, 13.357252459997124, 13.421984343339139, 13.422083669638587, 14.214883169540151]
16 : [3.220289763758684, 5.880492957780022, 5.880494420043111, 10.686000310009678, 10.693916904465476, 10.693985650230493, 12.317663957900697, 12.31784585945905, 13.357215865337436, 13.421975907139826, 13.422046282087193, 14.214854421597693]
17 : [3.220289763758687, 5.8804929577799605, 5.880494420042576, 10.686000307208516, 10.693916904213037, 10.693985649783874, 12.317663956136965, 12.317843226610497, 13.357199483944116, 13.421972060633866, 13.422029264205355, 14.214838524589577]
18 : [3.220289763758684, 5.880492957779919, 5.880494420042491, 10.686000306563274, 10.693916904154786, 10.693985649681093, 12.317663955526857, 12.3178423146984, 13.357192170632125, 13.421970310142441, 13.42202150583432, 14.21482976733362]
19 : [3.2202897637586836, 5.88049295777999, 5.880494420042796, 10.686000306415625, 10.693916904141503, 10.693985649657737, 12.317663955317236, 12.31784200132117, 13.357188926902761, 13.421969519700735, 13.422017990412309, 14.21482494224847]
[6]:
print ("Eigenvalues")
for lam in evals:
    print (lam)
Eigenvalues
3.2202897637586836
5.88049295777999
5.880494420042796
10.686000306415625
10.693916904141503
10.693985649657737
12.317663955317236
12.31784200132117
13.357188926902761
13.421969519700735
13.422017990412309
14.21482494224847
[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);
[ ]:

[ ]: