Missing faces when markering boundaries

I try to mark different parts of boundaries with different markers.I create a unit square mesh with each direction Nx=10,Ny=10 cells.I want to mark the part {x=0,0.3<y<1} as 1 and {x=0,0<y<0.3} as 2.And I use the following code to mark faces

boundaries=[(1,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]>=0.3,x[1]<=1))),
            (2,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]<=0.3,x[1]>=0))),
            (3,lambda x:np.isclose(x[1],0)),
            (4,lambda x:np.isclose(x[0],1)),
            (5,lambda x:np.isclose(x[1],1))]
facet_indices,facet_markers=[],[]
fdim=domain.topology.dim-1
for (marker,locator) in boundaries:
    facets=mesh.locate_entities(domain,fdim,locator)
    print(facets)
    facet_indices.append(facets)
    facet_markers.append(np.full_like(facets,marker))
facets_indices=np.hstack(facet_indices).astype(np.int32)
facets_markers=np.hstack(facet_markers).astype(np.int32)
sorted_facets=np.argsort(facets_indices)
facet_tag=mesh.meshtags(domain,fdim,facets_indices[sorted_facets],facets_markers[sorted_facets])

where I use print to show the indices.
However, I find that when markering marker=1 and marker=2, there has one missing face which doesnt belong to any of them, the following is the result print shows

[243 263 280 294 305 313 318]
[166 194]
[  0   5  14  26  41  59  80 104 131 161]
[  2   7  16  28  43  61  82 106 133 163]
[192 218 241 261 278 292 303 311 316 319]

Is it because the marker functions I define are wrong?
Thanks in advance!

Are you sure you dont want to use smaller equal and bigger equal to get the facets with a vertex at 0.3 and the other facet inside the bounds you specified?

Locate entities checks if all vertices of a given facet satisfies the criteria.

I try to change them to be
boundaries=[(1,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]>0.3,x[1]<=1))), (2,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]<0.3,x[1]>=0))),...]
but the result is the same.
Since I divide each direction with 10, so the third face at left boundary counted from bottom to top has y-coordinate of its two vertices equal to 0.2 and 0.3, so it should be included in the second boundary marker function.And it is really confused why it is not included.
I change the 0.3 to 0.2, the result is right.

Now I use

boundaries=[(1,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]>0.3-1e-10,x[1]<=1))), (2,lambda x:np.logical_and(np.isclose(x[0],0),np.logical_and(x[1]<0.3+1e-10,x[1]>=0))),...]

the result is right,the face is now included with marker =2.Maybe it is due to the round-off error?

You should always use a tolerance when you do smaller equal or greater equal, as numpy arrays usually only have accuracy up to 1e-16, thus things can be 3 ± 1e-16, 0±1e-16 etc

1 Like

Thanks a lot.It is really iimportant to notice.