Integration over a 2D region

Dear NGSolve community,

I am trying to integrate the function 1 over a small region within a 2D mesh.
Unfortunately, I get different values for finer or coarser meshes.
Should I do it in a different manner? I understand that the number of triangles inside the region might vary refining and, very much likely, they won’t be conforming to the region. However, it is clear that the area, which is the result of my computation, doesn’t change.

import ngsolve as ng
from netgen.geom2d import SplineGeometry
geo = SplineGeometry()
geo.AddRectangle((-1,-1),(1,1),bc="rectangle")
meshA = ng.Mesh(geo.GenerateMesh(maxh=0.2))

geo = SplineGeometry()
geo.AddRectangle((-1,-1),(1,1),bc="rectangle")
meshB = ng.Mesh(geo.GenerateMesh(maxh=0.1))

ng.Integrate(ng.IfPos(ng.x- -0.9393939393939394 , 1., 0.) * 
             ng.IfPos( -0.8787878787878788 -ng.x, 1., 0.) * 
             ng.IfPos(ng.y- -0.8181818181818181 , 1., 0.) * 
             ng.IfPos( -0.7575757575757576 -ng.y, 1., 0.), meshA)

ng.Integrate(ng.IfPos(ng.x- -0.9393939393939394 , 1., 0.) * 
             ng.IfPos( -0.8787878787878788 -ng.x, 1., 0.) * 
             ng.IfPos(ng.y- -0.8181818181818181 , 1., 0.) * 
             ng.IfPos( -0.7575757575757576 -ng.y, 1., 0.), meshB)

Thanks very much!
Francesco

You are integrating a characteristic function which is pretty unsmooth. Numerical integration will be inexact and depend on the integration points. It’s no surprise.

Try to integrate a step function in 1D with one or two elements on paper. You will observe the same behavior.

If you want to respect the geometry of your characteristic function, try the multilevelset function cut integration of ngsxfem. That will give you the same integral value numerically on different meshes (as long as the geometry can be described exactly).

Best,
Christoph

Dear Christoph,

Is this what you meant? I am sorry but I don’t know this add-on library. I tried to make it work looking at the examples. Unfortunately now I get a value of 0.

x0 = -0.9393939393939394
x1 = -0.8787878787878788
y0 = -0.8181818181818181
y1 = -0.7575757575757576

region = ng.IfPos(ng.x- x0 , 1., 0.) * ng.IfPos( x1 -ng.x, 1., 0.) *\
         ng.IfPos(ng.y- y0 , 1., 0.) * ng.IfPos( y1 -ng.y, 1., 0.)

import xfem
order=3
V = ng.H1(meshA,order=order,autoupdate=True)
region_approx = ng.GridFunction(V,autoupdate=True)
xfem.InterpolateToP1(region, region_approx)
ng.Integrate(ng.CoefficientFunction(1.) * xfem.dCut(levelset=region_approx, domain_type=xfem.NEG, order=order), mesh=meshA)

Best,
Francesco

level_sets = (ng.x- x0, x1 -ng.x, ng.y- y0, y1 -ng.y)
nr_ls = len(level_sets)
level_sets_p1 = tuple(ng.GridFunction(ng.H1(meshA, order=1)) for i in range(nr_ls))

for i, lset_p1 in enumerate(level_sets_p1):
    xfem.InterpolateToP1(level_sets[i], lset_p1)
  
ng.Integrate(ng.CoefficientFunction(1) * xfem.dCut(level_sets_p1, (xfem.POS, xfem.POS, xfem.POS, xfem.POS), order=1), mesh=meshA)

This is the way that works for me. I’m sorry, probably this is what you meant in the first place.
Thank you,
Francesco

exactly, that’s what I meant.

Dear Christoph,
I am having some difficulty installing xfem locally. I have a local installation of ngsolve that works, done through pip. I am doing the same for xfem. Unfortunately, even though I have no error, I can’t import the module in python. I get the following error

  File "<stdin>", line 1, in <module>
  File "/home/mantegaz/venv/lib/python3.10/site-packages/xfem/__init__.py", line 23, in <module>
    from xfem.ngsxfem_py import *
ImportError: libsolve.so: cannot open shared object file: No such file or directory

This is consistent with the following

ldd /home/mantegaz/venv/lib/python3.10/site-packages/xfem/ngsxfem_py.so
        linux-vdso.so.1 (0x00007ffd9a3ef000)
        libsolve.so => not found
        libngcomp.so => not found
        libngfem.so => not found
        libngla.so => not found
        libngbla.so => not found
        libngstd.so => not found
        libnglib.so => not found
        libngcore.so => not found

What am I missing?
I have looked for all these libraries and I am sure that they do exist, except libsolve.so . Evidently the installation through pip couldn’t find them.

Bests,
Francesco

Fixed! The recommended version of ngsolve was not the one actually needed.
It looks like you need 6.2.2302 for ngsolve, instead of 6.2.2301 as recommended in ngs_check.py.

Francesco