Indicator function of a rectangle

General Discussion
1

Hi everyone!
I tried to implement the indicator function of a rectangle Q as follows

// Q = [xq1,xq2] x [yq1,yq2], xq1<xq2, yq1<yq2
func real chiQ(real xq1,real xq2,real yq1,real yq2){
if(x>=xq1 & x<=xq2 & y>=yq1 & y<=yq2) return 1;
else return 0;
}

but I think I am wrong with something.
Indeed, by computing the integral of the indicator function of Q = [0.25,0.75] x [0.25,0.75], I should obtain the area of Q, that is 0.25, but I got results, such as 0.26, 0.22, varying based on the mesh on which I’m computing the integral.
Here the way I compute this integral:
cout<<int2d(Th) (chiQ(0.25,0.75,0.25,0.75))<<endl;

Do you have any idea about why?

2

Remerber , than if the rectangle is not edges of mesh, you have no chance to get the exact value
and the precission is in O(h) .