Why SubDomain.inside gets called the wrong number or times?

Question

Consider the following example:

from dolfin import *
mesh = UnitSquareMesh(1, 1)
sub_domains = CellFunction("size_t", mesh)
sub_domains.set_all(0)

switched_element_expression = Expression("x[0]")

class SwitchedElementSet(SubDomain):
    def __init__(self, switched_expression):
        SubDomain.__init__(self)
        self.switched_expression = switched_expression
        self.count = 0
    def inside(self, x, on_boundary):
        y = self.switched_expression(x)
        self.count += 1
        print self.count, y, y > 0.2
        return y > 0.2 #with this it gets called only once
        #return True   #with this it gets called six times

switched_element_set = SwitchedElementSet(switched_element_expression)
switched_element_set.mark(sub_domains, 1)

Why inside gets called only once?

Answer 1

You're marking a CellFunction. The UnitSquareMesh(1,1) only has two cells. Neither cell is entirely contained by x[0] > 0.2.

Try replacing return y > 0.2 with return y > -DOLFIN_EPS. This ensures that both cells are contained by the SwitchedElementSet.

You should see each of the vertices of the mesh, as well as the cell midpoints, being checked to confirm that each cell is indeed inside the SwitchedElementSet SubDomain.

Comments

Got it, thank you. It basically assumes that the inside test will be based on a continuous "CG" field of order 1, hence it checks each vertex only once in SubDomain.cpp, regardless of the cell. My original problem came from a discontinuous switched_element_expression field, but I think I will be able to work around it with a C++ function that lops over the cells.

Thank you again.