This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

AutoSubDomain fail to get marked

+1 vote

Hi, I've got a quite simple problem with boundary description. I have the following simplified code

from dolfin import *
from mshr import *

# Define Domain

diam    = 5.0
mainBox = Rectangle(Point(0,0),Point(diam,diam))

geo = mainBox

# Prepare domain labels

inletRadius     = 0.5
INmid       = Point(1.25,diam)
OUTmid      = Point(3.75,0.0)

def dist(mid,x):
    dx = mid[0] - x[0]
    dy = mid[1] - x[1]
    return sqrt(dx*dx + dy*dy)

def inlet(x):   
    if ( near(x[1],diam) and  dist(INmid,x) < inletRadius+DOLFIN_EPS ): 
        return True
    else:
        return False

def outlet(x):
    if ( near(x[1],0.0) and  dist(OUTmid,x) < inletRadius+DOLFIN_EPS ): 
        return True
    else:
        return False


def noslip(x):
    if ( not outlet(x) and not inlet(x) ) \
    and (   near(x[0],0.0) or near(x[0],diam)\
        or near(x[1],0.0) or near(x[1],diam) ):
        return True
    else:
        return False

# Build mesh

mesh3d = generate_mesh(geo, 10)

# Boundary description
bndry   = FacetFunction("size_t", mesh3d)

noV     = AutoSubDomain(lambda x,bndry: bndry and noslip(x) )
inletV      = AutoSubDomain(lambda x,bndry: bndry and inlet(x) )
outletV     = AutoSubDomain(lambda x,bndry: bndry and outlet(x) )

noV     .mark(bndry,1)
inletV      .mark(bndry,2)
outletV     .mark(bndry,3)

plot(bndry,interactive=True)

The funny thing is that 'noV' part of boundary is plotted/marked properly, i.e. functions 'inlet(x)' and 'outlet(x)' are evaluated correctly, BUT, the 'inletV' and 'outletV' parts are not marked at all.

I'm totally lost, since the functions inlet(x) and outlet(x) obviously works for 'noV' but not for inletV nor for outletV. I hope it is something really trivial that I'm missing (like typo) - but I'm looking at this for too long now and can't find the problem.

Thanks for your help.

asked Jun 20, 2015 by Ianx86 FEniCS User (1,050 points)
edited Jun 21, 2015 by Ianx86

1 Answer

0 votes

Ok, finally I realize the correct solution. It suffices to use the following definition of inletV and outletV

inletV      = AutoSubDomain(lambda x,bndry: inlet(x) )
outletV     = AutoSubDomain(lambda x,bndry: outlet(x) )

Nevertheless the reason why this works is still hidden to me.

answered Jun 25, 2015 by Ianx86 FEniCS User (1,050 points)
edited Jun 26, 2015 by Ianx86
...