Using levelset to integrate a fespace variable in a subregion of the mesh

Hi there!
I was wondering if it is possible to integrate a FESpace variable over a subregion of the original mesh.
Let’s suppose Th3D in my example is a cylinder with a height of 0.0 to 100 cm. I want to integrate the Pow FESpace variable over just the axial region where z≥10.0z \geq 10.0z≥10.0 and z≤20.0z \leq 20.0z≤20.0. I tried the following approach, but it only generates a vector composed entirely of zeros.

Thank you in advance!

Christian

ofstream fffout("power_density_results/PowerDensity_levelset_testing.txt");
real[int] PowerDensityRegionsLevelsetTest(64);

// Defintion of the function that has to be used in levelset
func flevelset = (z>= 10.0 && z<=20.0);
//

for (int i=0; i<64; i++)
{
  real xc =vecCenters(i,0);
  real yc =vecCenters(i,1);
  int thisregion = Th3D(xc,yc,10.0).region;
  PowerDensityRegionsLevelsetTest[i] = int3d(Th3D,levelset=flevelset)(PowerDensityFespaceVariable)/int3d(Th3D,levelset=flevelset)(1.0);
  fffout.scientific<<PowerDensityRegionsLevelsetTest[i]<<endl;

}

Hello,
according to the FreeFem documentation, when you write
int3d(Th3D,levelset=phi) the integral is performed where phi<0.
Thus you must change your definition of flevelset
Otherwise you can do also as in Vectorial PDE, conditions, stiffness matrix - #4 by fb77