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.22909326185347, 22.63644750413801, 38.759126813443395, 56.95688370738365, 64.68904673363234, 71.61510642142184, 76.89311292906328, 86.36950975162256, 90.92730267869732, 99.90890224792415, 108.07709800844383, 123.53078333222933]
1 : [3.2310126754128734, 6.086234263628024, 6.406899224019449, 12.05453347105542, 12.672421029060487, 13.138599688056795, 14.490032245599606, 17.723681707580635, 18.66258439802927, 19.346120604825614, 23.99916095234862, 25.45185223741474]
2 : [3.220379348186283, 5.886623822602872, 5.905548429293182, 10.853747390545367, 10.901287969870118, 11.033338300507763, 12.490608609795057, 13.71319396803288, 14.184813087378217, 15.348674545541842, 16.52982807616963, 17.964540346949384]
3 : [3.2202534507267315, 5.880768708677501, 5.881908488027696, 10.709133420827508, 10.726052065926849, 10.77865621910847, 12.345960396936647, 12.760676043069353, 13.760281355577183, 14.018326043543617, 14.607718515169996, 16.592214146773827]
4 : [3.2202514669005913, 5.880494279981986, 5.880556917690755, 10.689894953482902, 10.698983267274853, 10.718281105549405, 12.323621310854719, 12.515134574973551, 13.53091563614309, 13.646352490729777, 14.013638403415836, 15.844042097802744]
5 : [3.220251432780629, 5.880477114537668, 5.880480023263932, 10.686851036555511, 10.6947008910243, 10.700890407055072, 12.319109793660692, 12.416240501149808, 13.4094778898671, 13.552567887252483, 13.749216310649249, 15.246452463287355]
6 : [3.220251432184982, 5.880473890511277, 5.880477843318864, 10.686161008730588, 10.69403287675773, 10.695820538589272, 12.31800240065142, 12.365108986961188, 13.375722155017382, 13.488497017447648, 13.613981824241739, 14.813746361903023]
7 : [3.220251432174623, 5.88047367648449, 5.880477751399404, 10.685969287196356, 10.693915072917017, 10.694415397843947, 12.317696219797117, 12.338610207852287, 13.365112608760748, 13.452949540259493, 13.531662728447072, 14.548118344000931]
8 : [3.220251432174446, 5.88047366477277, 5.88047774601528, 10.685917731763386, 10.69389210813749, 10.694043457650482, 12.317604665085016, 12.32609806875343, 13.360878384950315, 13.435960447032443, 13.480566199062721, 14.398036524323059]
9 : [3.220251432174435, 5.880473664135097, 5.880477745705937, 10.68590487291473, 10.693887275314022, 10.693949433864253, 12.317575537838657, 12.320797130742582, 13.35887798663533, 13.428258104349512, 13.451525418295777, 14.315735103180945]
10 : [3.2202514321744427, 5.880473664100341, 5.8804777456885144, 10.685901792193969, 10.693886175010208, 10.693926453072343, 12.317565160998456, 12.318715754648863, 13.357892966814212, 13.424793974406857, 13.436298705898404, 14.270543581394165]
11 : [3.2202514321744404, 5.880473664098444, 5.880477745687534, 10.68590106692505, 10.693885916583648, 10.693920952759077, 12.31755883168241, 12.31793983407066, 13.357414569565583, 13.423228033651379, 13.42874751229699, 14.245626512602954]
12 : [3.2202514321744418, 5.880473664098303, 5.880477745687488, 10.685900897346363, 10.693885856844803, 10.693919648675614, 12.317545562904371, 12.317668127314365, 13.357188433505184, 13.422515875284118, 13.425120096969419, 14.231830210671182]
13 : [3.220251432174439, 5.880473664098292, 5.880477745687481, 10.685900857787226, 10.693885843152513, 10.693919341031883, 12.317512859694952, 12.317598523901948, 13.357083766102429, 13.422189470733239, 13.423410806402112, 14.22417804253281]
14 : [3.220251432174436, 5.880473664098299, 5.88047774568748, 10.685900848563495, 10.69388584001375, 10.693919268631673, 12.317488499745926, 12.317586423609148, 13.35703596440833, 13.422036750983787, 13.422616118830943, 14.219929747685024]
15 : [3.220251432174439, 5.88047366409833, 5.880477745687534, 10.685900846411215, 10.6938858392918, 10.693919251610007, 12.317478501716316, 12.317583475289172, 13.357014276540866, 13.421959835910728, 13.422253878445416, 14.217569542190473]
16 : [3.2202514321744378, 5.880473664098262, 5.880477745687419, 10.685900845908687, 10.693885839124496, 10.693919247606473, 12.317474823695743, 12.317582557534537, 13.357004460574343, 13.421913950688014, 13.4220966672754, 14.216257650183788]
17 : [3.2202514321744427, 5.880473664098334, 5.880477745687393, 10.685900845790444, 10.693885839086462, 10.69391924666297, 12.31747350181891, 12.317582246919407, 13.356999997047856, 13.421883405209485, 13.422032581656998, 14.215519222677496]
18 : [3.2202514321744404, 5.880473664098291, 5.880477745687564, 10.685900845763518, 10.693885839077447, 10.693919246443294, 12.317473040899285, 12.317582140606074, 13.356997985525794, 13.421865761626224, 13.422006263423249, 14.215116842195995]
19 : [3.2202514321744435, 5.880473664098214, 5.880477745687396, 10.685900845757331, 10.693885839075836, 10.693919246393392, 12.317472878697938, 12.317582103544755, 13.356997100164, 13.421856651236084, 13.421996068749117, 14.214894115986747]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.2202514321744435
5.880473664098214
5.880477745687396
10.685900845757331
10.693885839075836
10.693919246393392
12.317472878697938
12.317582103544755
13.356997100164
13.421856651236084
13.421996068749117
14.214894115986747
[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);
[ ]:
[ ]: