24.7. Ranks of discrete operators

from ngsolve import *
from netgen.occ import *

box = Box((0,0,0), (1,1,1))
mesh = Mesh(OCCGeometry(box).GenerateMesh(maxh=0.3))

fesH1 = H1(mesh, order=1)
fesHc = HCurl(mesh, order=0)
fesHd = HDiv(mesh, order=0)
fesL2 = L2(mesh, order=0)

mesh.nv, mesh.nedge, mesh.nface, mesh.GetNE(VOL)

(76, 337, 458, 196)

fesH1.ndof, fesHc.ndof, fesHd.ndof, fesL2.ndof

(76, 337, 458, 196)

ug,vg = fesH1.TnT()
uc,vc = fesHc.TnT()
ud,vd = fesHd.TnT()
ul2,vl2 = fesL2.TnT()

opgrad = ConvertOperator(fesH1, fesHc, grad(ug)).DeleteZeroElements(1e-10)
opcurl = ConvertOperator(fesHc, fesHd, curl(uc)).DeleteZeroElements(1e-10)
opdiv = ConvertOperator(fesHd, fesL2, div(ud)).DeleteZeroElements(1e-10)

compute ranks of the discrete operators:

from scipy.linalg import svd
import numpy as np

sgrad = svd(opgrad.ToDense())[1]    # singular values
print ("rank of grad:", np.sum(sgrad > 1e-10))

scurl = svd(opcurl.ToDense())[1]    # singular values
print ("rank of curl:", np.sum(scurl > 1e-10))

sdiv = svd(opdiv.ToDense())[1]    # singular values
print ("rank of div:", np.sum(sdiv > 1e-10))

rank of grad: 75
rank of curl: 262
rank of div: 196

24.7.1. Exercise:

• interprete these numbers

• create domains with holes and tunnels, and recompute space dimensions and ranks