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)
0 : [8.213129817349921, 17.7998132400286, 37.21982565076547, 44.01930032728756, 62.182333369408234, 65.53333262079437, 73.9871630782097, 81.08231276873838, 88.06001749342057, 91.25055696920221, 113.61980588433377, 122.21673389691519]
1 : [3.238097539168596, 5.998650630503755, 6.161390072110305, 11.48586055147025, 12.0373096809818, 14.413520500604227, 15.234814501132997, 16.218269369657687, 17.346726953342955, 19.89228392651031, 22.638813210092565, 24.80756094219778]
2 : [3.2205426885757795, 5.884876708897553, 5.888226095524722, 10.785471632196765, 10.847156783330668, 11.381986377559613, 13.34323831889798, 13.372706506679213, 14.112393285860007, 14.98863956309097, 15.344080015841682, 17.688329571272813]
3 : [3.2202954145601526, 5.880713446554425, 5.880839526590908, 10.707466526328226, 10.720069178418957, 10.94477331075434, 12.5279739547512, 13.075036550608978, 13.540697915830394, 14.03069574816355, 14.367494217892018, 15.906741582814744]
4 : [3.220289894102191, 5.880502253291033, 5.880518097321799, 10.695114363080107, 10.699129958358903, 10.78775181162691, 12.354098370696558, 12.87105303876328, 13.445683487765075, 13.649403154306295, 14.21385995940532, 14.607665448576801]
5 : [3.220289766421207, 5.8804933752546535, 5.880495963116381, 10.692776367606495, 10.694999509696284, 10.72190631735076, 12.324202578878726, 12.640658945580281, 13.417157954184695, 13.509762214957393, 13.816759276193379, 14.260863115194164]
6 : [3.2202897638070183, 5.880492978117793, 5.880494510444926, 10.691097181810724, 10.694161153653855, 10.698643594374873, 12.318883281765203, 12.46074660512757, 13.40239495982758, 13.454795875816066, 13.563293263180418, 14.232709516175236]
7 : [3.2202897637595163, 5.880492958816448, 5.88049442503793, 10.687861627679942, 10.693985793780026, 10.694552792296216, 12.317899571105517, 12.372364426600308, 13.387450842686718, 13.433939412846312, 13.473271158032608, 14.223431358594558]
8 : [3.2202897637587, 5.88049295783439, 5.880494420312164, 10.686464368102458, 10.693942498545262, 10.694085266048775, 12.317711055329575, 12.33739216763752, 13.37373542762483, 13.426640065252839, 13.441852510574916, 14.219267826323446]
9 : [3.220289763758686, 5.880492957782836, 5.880494420057038, 10.68610869895639, 10.693925371312453, 10.694004893998049, 12.317673740213685, 12.324680164655264, 13.36519418831342, 13.42392124781682, 13.430320729423793, 14.217191023859078]
10 : [3.2202897637586836, 5.8804929577800475, 5.880494420043343, 10.68602521651851, 10.69391907743561, 10.693989757202504, 12.317666113412292, 12.320210206283717, 13.360871916540283, 13.422816823158326, 13.425650479457058, 14.216103594739844]
11 : [3.2202897637586814, 5.880492957779912, 5.88049442004261, 10.686006011029626, 10.693917414576521, 10.693986570956575, 12.317664492531842, 12.318659201330982, 13.358849291401986, 13.422344037685738, 13.423640051707773, 14.215520283899515]
12 : [3.2202897637586827, 5.880492957779909, 5.880494420042564, 10.686001612158089, 10.693917021615775, 10.693985859212555, 12.317664119926183, 12.318123653399686, 13.357930692680727, 13.42213653988282, 13.422748129735327, 14.215203524593415]
13 : [3.2202897637586836, 5.88049295777991, 5.880494420042561, 10.686000605583322, 10.693916931075414, 10.693985697539018, 12.317664014902157, 12.317939028826782, 13.357518511143082, 13.422044243388502, 13.422347251100593, 14.215030393598846]
14 : [3.220289763758682, 5.8804929577799125, 5.880494420042574, 10.686000375002935, 10.693916910323367, 10.693985660615963, 12.31766397719323, 12.317875402301393, 13.357334515991584, 13.422002981963155, 13.422165851307666, 14.21493541085034]
15 : [3.2202897637586814, 5.8804929577799445, 5.880494420042837, 10.68600032215073, 10.693916905558906, 10.69398565216849, 12.317663963000582, 12.317853445810808, 13.35725246267589, 13.42198435678698, 13.42208362458691, 14.214883178536455]
16 : [3.220289763758686, 5.880492957779887, 5.880494420042669, 10.686000309998978, 10.693916904464647, 10.693985650228994, 12.317663957906527, 12.317845852968825, 13.357215844589632, 13.421975911896313, 13.422046230788036, 14.214854429125403]
17 : [3.220289763758686, 5.8804929577799125, 5.880494420042538, 10.686000307205848, 10.693916904212312, 10.693985649783444, 12.317663956132094, 12.317843225571751, 13.357199489055166, 13.421972039816293, 13.422029225948219, 14.214838502371103]
18 : [3.2202897637586854, 5.88049295777986, 5.880494420042537, 10.68600030655121, 10.693916904154227, 10.693985649678723, 12.317663955503614, 12.317842287981303, 13.357191865758454, 13.421970282250724, 13.422021305219504, 14.214829336487933]
19 : [3.2202897637586836, 5.880492957779674, 5.880494420042459, 10.686000306410595, 10.693916904141066, 10.693985649656394, 12.31766395529581, 12.317841975319016, 13.357188482748915, 13.42196944557237, 13.42201786466157, 14.214824054237692]
[6]:
print ("Eigenvalues")
for lam in evals:
    print (lam)
Eigenvalues
3.2202897637586836
5.880492957779674
5.880494420042459
10.686000306410595
10.693916904141066
10.693985649656394
12.31766395529581
12.317841975319016
13.357188482748915
13.42196944557237
13.42201786466157
14.214824054237692
[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);
[ ]:

[ ]: