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.229093264873581, 22.636447497960596, 38.75912682011573, 56.95688373252868, 64.68904674457876, 71.61510642577963, 76.89311296711642, 86.36950975632872, 90.92730267885543, 99.90890224785282, 108.0770980138028, 123.53078334097053]
1 : [3.231012675430815, 6.086234263751503, 6.406899224299318, 12.054533474069045, 12.67242103010273, 13.138599688515992, 14.490032248663178, 17.72368170754059, 18.662584400170466, 19.34612060372826, 23.99916095311159, 25.451852239638114]
2 : [3.2203793481864817, 5.886623822609479, 5.905548429299362, 10.853747390711698, 10.901287970203162, 11.033338300766479, 12.49060860997046, 13.713193967479356, 14.184813087421976, 15.348674545471054, 16.529828073177566, 17.964540347054847]
3 : [3.220253450726732, 5.880768708677545, 5.881908488027794, 10.709133420831666, 10.726052065996917, 10.778656219125445, 12.345960396967433, 12.760676042782727, 13.76028135557727, 14.018326042410783, 14.607718515540343, 16.5922141479089]
4 : [3.2202514669005917, 5.880494279981969, 5.8805569176907335, 10.68989495348774, 10.698983267281399, 10.71828110554932, 12.323621310865128, 12.515134574876944, 13.530915635872402, 13.646352490695305, 14.01363840387991, 15.84404209968021]
5 : [3.220251432780628, 5.880477114537641, 5.880480023263928, 10.686851036556705, 10.694700891025006, 10.700890407055608, 12.31910979366469, 12.416240501140667, 13.409477889811503, 13.55256788724539, 13.749216311056772, 15.246452465162703]
6 : [3.2202514321849844, 5.880473890511284, 5.88047784331887, 10.686161008730831, 10.694032876757856, 10.695820538589883, 12.318002400652832, 12.365108986981532, 13.375722155012593, 13.48849701744168, 13.613981824595317, 14.813746363271134]
7 : [3.220251432174622, 5.880473676484498, 5.880477751399405, 10.685969287196382, 10.693915072917065, 10.694415397844237, 12.3176962197976, 12.338610207872089, 13.365112608763097, 13.452949540256789, 13.531662728716293, 14.548118344835192]
8 : [3.220251432174442, 5.880473664772775, 5.880477746015278, 10.685917731763396, 10.693892108137483, 10.69404345765055, 12.317604665085202, 12.32609806877186, 13.360878384950926, 13.435960447037093, 13.480566199359023, 14.398036524965532]
9 : [3.2202514321744404, 5.8804736641351365, 5.880477745705936, 10.685904872914037, 10.693887275313262, 10.693949433863503, 12.317575537838549, 12.320797130574396, 13.35887798662936, 13.428258105375482, 13.451525417077468, 14.315735094027982]
10 : [3.220251432174442, 5.880473664100282, 5.880477745688511, 10.685901792193503, 10.693886175010238, 10.693926453071917, 12.31756516099856, 12.31871575398835, 13.357892966248826, 13.424793978612122, 13.436298706131357, 14.27054355151956]
11 : [3.2202514321744435, 5.880473664098413, 5.880477745687532, 10.685901066924892, 10.693885916583804, 10.693920952758857, 12.3175588316756, 12.317939833690394, 13.357414569195097, 13.423228037278756, 13.428747514453748, 14.245626485687902]
12 : [3.2202514321744404, 5.880473664098331, 5.880477745687489, 10.685900897346357, 10.693885856844478, 10.693919648675267, 12.317545562886854, 12.317668127256633, 13.357188433147362, 13.422515876210108, 13.425120100781639, 14.231830206012276]
13 : [3.220251432174438, 5.880473664098298, 5.88047774568747, 10.68590085778918, 10.693885843152671, 10.693919341038288, 12.317512862414628, 12.31759852561525, 13.357083769938585, 13.422189472994024, 13.423410986469026, 14.224178858958897]
14 : [3.2202514321744404, 5.8804736640982975, 5.880477745687464, 10.685900848563987, 10.693885840013364, 10.693919268633383, 12.31748850217273, 12.317586424529747, 13.357035967824656, 13.422036759989606, 13.42261626439888, 14.219930321360685]
15 : [3.220251432174443, 5.880473664098306, 5.880477745687474, 10.685900846411329, 10.69388583929156, 10.69391925160967, 12.317478503442619, 12.31758347590467, 13.35701427905396, 13.421959838675496, 13.422254067423328, 14.217570076837106]
16 : [3.220251432174442, 5.880473664098296, 5.880477745687508, 10.6859008459087, 10.693885839125182, 10.693919247607116, 12.31747482767613, 12.317582559263894, 13.357004470357978, 13.421913935773118, 13.422096846655153, 14.216257792386632]
17 : [3.2202514321744418, 5.880473664098327, 5.880477745687428, 10.685900845789417, 10.6938858390868, 10.693919246649498, 12.317473454588328, 12.317582235144716, 13.35700000326222, 13.421882201470778, 13.422030273838304, 14.21548561242283]
18 : [3.2202514321744427, 5.880473664098299, 5.880477745687465, 10.68590084576292, 10.693885839077968, 10.693919246439673, 12.317473015102557, 12.317582134427582, 13.356997997498462, 13.421864545359407, 13.422005787796177, 14.215085116023557]
19 : [3.2202514321744404, 5.880473664098209, 5.880477745687459, 10.685900845757136, 10.693885839075929, 10.693919246392024, 12.317472865234624, 12.317582100366804, 13.356997082502032, 13.421855675240014, 13.421995684151984, 14.214866118578447]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.2202514321744404
5.880473664098209
5.880477745687459
10.685900845757136
10.693885839075929
10.693919246392024
12.317472865234624
12.317582100366804
13.356997082502032
13.421855675240014
13.421995684151984
14.214866118578447
[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);
[ ]:
[ ]: