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 = 42562
[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 : [6.229093264016037, 22.636447487152928, 38.759126814004915, 56.95688373128489, 64.68904675251993, 71.61510642964112, 76.89311294878544, 86.36950974785135, 90.92730267845576, 99.90890224154136, 108.07709801268318, 123.53078334386488]
1 : [3.231012675425504, 6.086234263668395, 6.406899224308958, 12.054533474436528, 12.672421029552494, 13.138599688742014, 14.490032246800967, 17.7236817076432, 18.66258439869991, 19.34612060416855, 23.999160955795507, 25.451852238753077]
2 : [3.2203793481864347, 5.886623822608349, 5.905548429302, 10.853747390797615, 10.90128797011855, 11.03333830068371, 12.490608609935009, 13.713193967572682, 14.184813087352996, 15.34867454484223, 16.52982807319949, 17.96454034795094]
3 : [3.2202534507267306, 5.880768708677551, 5.881908488027979, 10.709133420833274, 10.726052065992498, 10.778656219085692, 12.345960396962447, 12.760676042846962, 13.76028135552251, 14.018326042065077, 14.607718515461615, 16.59221414926929]
4 : [3.2202514669005975, 5.880494279981993, 5.880556917690756, 10.689894953486856, 10.698983267279418, 10.718281105536073, 12.323621310863306, 12.515134574925504, 13.530915635830729, 13.646352490654957, 14.013638403754763, 15.844042101366297]
5 : [3.2202514327806324, 5.880477114537638, 5.880480023263909, 10.686851036555908, 10.694700891024304, 10.700890407052356, 12.319109793663774, 12.416240501185367, 13.409477889811722, 13.552567887243343, 13.749216310982947, 15.246452466761616]
6 : [3.2202514321849836, 5.880473890511258, 5.880477843318854, 10.686161008730483, 10.69403287675768, 10.69582053858934, 12.318002400652404, 12.365108987016681, 13.375722155016895, 13.488497017446313, 13.613981824594303, 14.813746364445949]
7 : [3.220251432174621, 5.88047367648449, 5.8804777513994155, 10.685969287196304, 10.693915072917, 10.694415397844162, 12.317696219797451, 12.33861020789386, 13.365112608765903, 13.452949540258942, 13.531662728749435, 14.548118345565937]
8 : [3.2202514321744444, 5.880473664772746, 5.8804777460152975, 10.685917731763363, 10.693892108137433, 10.694043457650599, 12.317604665085128, 12.326098068770873, 13.360878384952601, 13.43596044702771, 13.480566199145386, 14.3980365250577]
9 : [3.220251432174435, 5.880473664135121, 5.8804777457059405, 10.685904872914808, 10.693887275314248, 10.693949433864507, 12.31757553783861, 12.320797130603848, 13.358877986605034, 13.428258103833539, 13.451525414951107, 14.315735101361371]
10 : [3.2202514321744427, 5.880473664100316, 5.880477745688507, 10.685901792193821, 10.693886175010512, 10.693926453072415, 12.317565160998418, 12.31871575379432, 13.357892966495195, 13.424793976652712, 13.436298699102318, 14.270543555603394]
11 : [3.2202514321744387, 5.880473664098439, 5.880477745687515, 10.685901066924925, 10.693885916583975, 10.693920952758551, 12.31755883167538, 12.31793983342744, 13.357414569334887, 13.423228035798676, 13.428747507415528, 14.24562648260975]
12 : [3.220251432174439, 5.8804736640983215, 5.880477745687476, 10.685900897346398, 10.693885856844844, 10.693919648675092, 12.31754556287718, 12.317668127164769, 13.357188433283259, 13.422515876156938, 13.425120106257358, 14.23183021438993]
13 : [3.220251432174441, 5.880473664098298, 5.880477745687468, 10.685900857789003, 10.69388584315233, 10.69391934104211, 12.317512863811661, 12.317598526707068, 13.357083767333409, 13.422189471915251, 13.42341102133127, 14.224178983902583]
14 : [3.220251432174438, 5.880473664098296, 5.880477745687473, 10.685900848564195, 10.693885840013051, 10.693919268634566, 12.317488503271386, 12.317586424792355, 13.357035966561412, 13.422036750537055, 13.422616280128718, 14.219930662060026]
15 : [3.2202514321744453, 5.880473664098317, 5.8804777456874495, 10.68590084641151, 10.693885839291504, 10.6939192516096, 12.317478504562649, 12.317583476185536, 13.35701427880973, 13.421959832783699, 13.422254166968807, 14.217570646410673]
16 : [3.2202514321744413, 5.880473664098246, 5.880477745687427, 10.685900845908526, 10.693885839125343, 10.69391924760622, 12.317474828730191, 12.317582559762524, 13.357004462300361, 13.421913902905228, 13.422096978427177, 14.21625782860729]
17 : [3.2202514321744418, 5.88047366409832, 5.8804777456875446, 10.685900845791156, 10.693885839086736, 10.693919246665889, 12.31747351392447, 12.317582250132384, 13.35700002684088, 13.4218837610775, 13.422033304751844, 14.21552892105699]
18 : [3.220251432174442, 5.880473664098314, 5.880477745687456, 10.685900845763177, 10.693885839077913, 10.693919246442308, 12.317473041641328, 12.31758214158059, 13.35699792208318, 13.42186555648969, 13.422007716926473, 14.215112614482354]
19 : [3.220251432174442, 5.880473664098385, 5.880477745687463, 10.685900845757075, 10.693885839075776, 10.693919246392827, 12.31747287834129, 12.31758210389758, 13.356997049608681, 13.421856451617586, 13.421996929632373, 14.214888171460952]
[6]:
print ("Eigenvalues")
for lam in evals:
    print (lam)
Eigenvalues
3.220251432174442
5.880473664098385
5.880477745687463
10.685900845757075
10.693885839075776
10.693919246392827
12.31747287834129
12.31758210389758
13.356997049608681
13.421856451617586
13.421996929632373
14.214888171460952
[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);
[ ]:

[ ]: