specify intrules for 1d mesh

Main Forum Old-Forum - do not post here
1

Hello,

I am trying to generate the lumped mass matrix by setting integration rules.

For a 2d surface mesh in 3d, the following code works well.

ir = IntegrationRule(points = [(0,0), (1,0), (0,1)], weights = [1/6, 1/6, 1/6] )
ds_lumping = ds(intrules = { TRIG : ir })

However, for a 1d curved mesh in 2d, I guess it should be

ir = IntegrationRule(points = [(0,0), (1,0)], weights = [1/2, 1/2] )
ds_lumping = ds(intrules = { TRIG : ir })

However, it does not work.

Here is my code. I hope that the final result is a diagonal matrix.

import netgen.gui
import netgen.meshing as ngm
from ngsolve import *

p = [[1,0,0],[1,1,0],[0,1,0],[0,0,0]]
nr = len(p)
mymesh= ngm.Mesh(dim=2)
pids = []
for i in range(nr):
    pids.append (mymesh.Add (ngm.MeshPoint(ngm.Pnt(p[i]))))  

idx = mymesh.AddRegion("material", dim=1)
for i in range(nr):
    if i<nr-1:
        mymesh.Add(ngm.Element1D([pids[i],pids[i+1]],index=idx))
    else:
        mymesh.Add(ngm.Element1D([pids[i],pids[0]],index=idx))
    mesh = Mesh(mymesh)

fes = H1(mesh,order=1)
u,v = fes.TnT()

ir = IntegrationRule(points = [(0,0), (1,0)], weights = [1/2,1/2])
ds_lumping = ds(intrules = { TRIG : ir }, definedon=mesh.Boundaries("material"))

a = BilinearForm(fes)
a += u*v*ds_lumping

a.Assemble()
print(a.mat)

Thank you very much for your time in advance.

2

Hi Jiashun,

for 1D Elements (in 2D) you need to specify the integration rule for segments, not for triangles. If you change the code to

ir = IntegrationRule(points = [(0,0), (1,0)], weights = [1/2,1/2])
ds_lumping = ds(intrules = { SEGM : ir }, definedon=mesh.Boundaries("material"))

then you should get a diagonal mass matrix.

Best
Michael

3

It works! Thanks for your explanation.