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.229093264088796, 22.63644750021968, 38.7591268074329, 56.95688373189905, 64.68904674067471, 71.61510645567323, 76.89311294098684, 86.36950973809988, 90.92730267815365, 99.90890223837845, 108.07709800706414, 123.53078334439867]
1 : [3.2310126754290005, 6.086234263752179, 6.406899224043229, 12.054533474924122, 12.672421028411293, 13.138599686815688, 14.490032247184844, 17.72368170574987, 18.662584397407336, 19.346120605404124, 23.999160963409903, 25.451852239301758]
2 : [3.220379348186456, 5.886623822607963, 5.905548429307376, 10.853747390985935, 10.901287970134389, 11.03333830047591, 12.490608609927285, 13.71319396637087, 14.184813087489335, 15.348674547736547, 16.529828081490365, 17.96454034960235]
3 : [3.2202534507267337, 5.880768708677619, 5.881908488028921, 10.709133420843477, 10.72605206603876, 10.778656219123098, 12.345960396965067, 12.760676042627274, 13.760281355821052, 14.018326044808092, 14.607718515929712, 16.592214150627004]
4 : [3.2202514669005953, 5.88049427998201, 5.880556917690821, 10.689894953484133, 10.698983267293832, 10.71828110555591, 12.32362131086209, 12.515134574883007, 13.530915636662087, 13.646352490647383, 14.013638403923258, 15.844042102884288]
5 : [3.220251432780632, 5.8804771145376336, 5.880480023263901, 10.686851036555376, 10.694700891027196, 10.700890407059877, 12.31910979366284, 12.416240501191043, 13.409477890036053, 13.552567887229943, 13.749216311176202, 15.246452468104627]
6 : [3.2202514321849836, 5.880473890511293, 5.880477843318865, 10.68616100873049, 10.694032876758227, 10.695820538591883, 12.318002400652107, 12.36510898703465, 13.375722155096774, 13.488497017424129, 13.613981824788235, 14.813746365350404]
7 : [3.2202514321746194, 5.880473676484502, 5.8804777513993995, 10.685969287196327, 10.693915072917113, 10.694415397844987, 12.317696219797337, 12.338610207908276, 13.36511260880564, 13.452949540240063, 13.531662728896436, 14.548118346083903]
8 : [3.2202514321744475, 5.880473664772748, 5.880477746015268, 10.685917731763391, 10.693892108137481, 10.694043457650748, 12.317604665085106, 12.326098068790301, 13.360878384976116, 13.435960447024712, 13.480566199468338, 14.398036525649388]
9 : [3.220251432174443, 5.8804736641351, 5.880477745705931, 10.68590487291453, 10.693887275312727, 10.693949433866182, 12.317575537840744, 12.32079713143904, 13.358877986662575, 13.428258103729869, 13.451525421854582, 14.315735113904907]
10 : [3.220251432174444, 5.8804736641003155, 5.880477745688506, 10.68590179219364, 10.69388617500978, 10.693926453073297, 12.31756516099858, 12.318715754394306, 13.357892966264878, 13.424793976977238, 13.436298704118734, 14.270543561283421]
11 : [3.2202514321744444, 5.880473664098402, 5.880477745687525, 10.685901066924883, 10.693885916583579, 10.69392095275884, 12.31755883167368, 12.317939833650952, 13.357414569189489, 13.423228036207517, 13.428747507155766, 14.24562647937341]
12 : [3.2202514321744444, 5.880473664098303, 5.8804777456874655, 10.68590089734635, 10.693885856844744, 10.693919648675251, 12.317545562872866, 12.317668127206618, 13.357188433262388, 13.422515876862162, 13.42512010294933, 14.231830202025922]
13 : [3.2202514321744404, 5.8804736640983215, 5.880477745687496, 10.685900857786553, 10.693885843152772, 10.693919341024397, 12.317512856712282, 12.317598521558184, 13.357083766211474, 13.422189470464152, 13.42341073472244, 14.224177628050795]
14 : [3.220251432174439, 5.880473664098299, 5.880477745687474, 10.685900848563225, 10.693885840013854, 10.693919268627607, 12.317488499158648, 12.31758642359948, 13.357035963784853, 13.42203675463739, 13.422616307483743, 14.219929801303488]
15 : [3.220251432174444, 5.880473664098363, 5.880477745687188, 10.685900846411, 10.693885839291843, 10.693919251608317, 12.317478503467783, 12.317583476033713, 13.357014275646097, 13.421959823729392, 13.42225422305614, 14.217570177744546]
16 : [3.2202514321744427, 5.880473664098332, 5.880477745687485, 10.685900845908515, 10.693885839125034, 10.693919247603196, 12.317474822373725, 12.317582557769454, 13.357004464909275, 13.421913878740591, 13.422096867113751, 14.21625660286332]
17 : [3.220251432174442, 5.88047366409822, 5.880477745687368, 10.685900845790712, 10.693885839085802, 10.693919246655923, 12.317473496384858, 12.31758224570242, 13.357000016663953, 13.421883426343364, 13.422033241757564, 14.215521447095176]
18 : [3.2202514321744427, 5.880473664098514, 5.880477745687513, 10.685900845763632, 10.693885839077847, 10.69391924644308, 12.317473039000927, 12.317582140358553, 13.35699801278686, 13.421865743187334, 13.422007669692468, 14.215116544165959]
19 : [3.220251432174444, 5.880473664098146, 5.88047774568748, 10.685900845757219, 10.693885839075763, 10.693919246393488, 12.317472877742524, 12.31758210342281, 13.356997120071576, 13.421856610628058, 13.42199689925418, 14.214892713222143]
[6]:
print ("Eigenvalues")
for lam in evals:
print (lam)
Eigenvalues
3.220251432174444
5.880473664098146
5.88047774568748
10.685900845757219
10.693885839075763
10.693919246393488
12.317472877742524
12.31758210342281
13.356997120071576
13.421856610628058
13.42199689925418
14.214892713222143
[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);
[ ]:
[ ]: