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.229093262919805, 22.636447503123236, 38.75912681267974, 56.956883729055605, 64.68904675330454, 71.61510644513578, 76.89311293584277, 86.3695097586251, 90.92730267240465, 99.90890223997273, 108.07709801482166, 123.53078334403475]
1 : [3.2310126754284525, 6.086234263790956, 6.406899223953482, 12.054533474928775, 12.672421028838158, 13.138599686217587, 14.490032249114888, 17.723681707094162, 18.662584397599993, 19.346120604903856, 23.999160961793955, 25.451852239517347]
2 : [3.2203793481864578, 5.886623822609814, 5.905548429291393, 10.85374739087083, 10.901287970225521, 11.033338300341473, 12.490608609984546, 13.713193967665093, 14.184813087400196, 15.348674546192912, 16.529828077953077, 17.96454034896391]
3 : [3.2202534507267333, 5.8807687086776985, 5.8819084880279, 10.709133420845038, 10.726052066008155, 10.778656219090808, 12.345960396971481, 12.760676043147763, 13.76028135565803, 14.018326043290747, 14.607718516078936, 16.592214150487024]
4 : [3.2202514669005953, 5.8804942799820275, 5.880556917690759, 10.689894953486982, 10.698983267285694, 10.71828110555053, 12.323621310864985, 12.515134575145314, 13.530915636202302, 13.64635249061157, 14.013638404056183, 15.844042102830835]
5 : [3.220251432780631, 5.880477114537635, 5.880480023263936, 10.686851036556758, 10.694700891025724, 10.700890407058617, 12.319109793664268, 12.41624050133359, 13.409477889904897, 13.552567887216538, 13.749216311224064, 15.24645246809157]
6 : [3.2202514321849836, 5.880473890511257, 5.880477843318853, 10.68616100873105, 10.694032876757975, 10.69582053859147, 12.31800240065267, 12.365108987106415, 13.375722155050253, 13.488497017427111, 13.61398182479329, 14.813746365362405]
7 : [3.220251432174619, 5.880473676484494, 5.880477751399419, 10.685969287196523, 10.693915072917076, 10.694415397844855, 12.317696219797535, 12.33861020794073, 13.365112608784575, 13.452949540246664, 13.531662728896505, 14.548118346109739]
8 : [3.2202514321744475, 5.880473664772752, 5.880477746015272, 10.685917731763407, 10.693892108137447, 10.694043457650768, 12.317604665085174, 12.326098068794035, 13.360878384965634, 13.435960447030087, 13.48056619933664, 14.398036525465683]
9 : [3.2202514321744413, 5.880473664135074, 5.880477745705953, 10.68590487291463, 10.693887275313898, 10.693949433864132, 12.317575537839504, 12.320797130806886, 13.358877986668402, 13.428258104898768, 13.451525419547613, 14.315735101372848]
10 : [3.2202514321744378, 5.880473664100336, 5.880477745688493, 10.68590179219366, 10.693886175010414, 10.693926453072134, 12.31756516099857, 12.318715753596004, 13.357892966321543, 13.424793978005066, 13.436298697930154, 14.270543541281574]
11 : [3.22025143217444, 5.880473664098413, 5.880477745687519, 10.68590106692484, 10.693885916583973, 10.693920952758221, 12.317558831671347, 12.317939833165791, 13.3574145692157, 13.4232280369211, 13.428747503427514, 14.245626461663132]
12 : [3.220251432174441, 5.880473664098302, 5.880477745687463, 10.68590089734637, 10.693885856844776, 10.693919648674896, 12.317545562853033, 12.317668127019495, 13.357188433100344, 13.422515876276444, 13.42512010078349, 14.231830196914872]
13 : [3.2202514321744404, 5.880473664098313, 5.880477745687461, 10.685900857788392, 10.69388584315231, 10.69391934103966, 12.31751286214783, 12.317598525626977, 13.357083768451776, 13.422189473084996, 13.423410898306418, 14.224178577189573]
14 : [3.22025143217444, 5.8804736640983055, 5.880477745687458, 10.685900848563529, 10.693885840013628, 10.693919268634227, 12.317488501893822, 12.317586424454962, 13.35703596592148, 13.422036767289326, 13.42261627523339, 14.219930141660981]
15 : [3.2202514321744387, 5.880473664098304, 5.880477745687444, 10.685900846411013, 10.693885839291621, 10.693919251609847, 12.317478504438789, 12.317583476337408, 13.357014272042939, 13.421959825358657, 13.422254160625208, 14.217569954274996]
16 : [3.2202514321744427, 5.880473664098457, 5.8804777456875605, 10.685900845908739, 10.693885839120744, 10.693919247577664, 12.317474763190981, 12.317582538545622, 13.3570044540072, 13.42191368723722, 13.42209573446786, 14.216246171987493]
17 : [3.2202514321744404, 5.880473664098386, 5.880477745687429, 10.685900845791004, 10.693885839084295, 10.693919246650141, 12.317473459006573, 12.317582233054281, 13.357000030590362, 13.421883185208177, 13.422032122181019, 14.215510871298408]
18 : [3.2202514321744427, 5.880473664098422, 5.880477745687494, 10.685900845763829, 10.693885839077295, 10.693919246440636, 12.31747302004677, 12.317582134342974, 13.356998038696517, 13.421865422518145, 13.422006992277918, 14.215106642530591]
19 : [3.2202514321744427, 5.88047366409821, 5.880477745687419, 10.685900845757104, 10.69388583907568, 10.69391924639162, 12.317472868221584, 12.317582100677221, 13.356997066149686, 13.421856259657634, 13.421996310692647, 14.21488204626141]
[6]:
print ("Eigenvalues")
for lam in evals:
    print (lam)
Eigenvalues
3.2202514321744427
5.88047366409821
5.880477745687419
10.685900845757104
10.69388583907568
10.69391924639162
12.317472868221584
12.317582100677221
13.356997066149686
13.421856259657634
13.421996310692647
14.21488204626141
[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);
[ ]:

[ ]: