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

" The point is not inside the domain" when using two meshes

0 votes

Hello
I have a problem including two parts. In the first part I have created the mesh in FEniCS and for the second part I imported the mesh from the GMSH.
The goal in the first part is to calculate of variable "c" and use it to find the displacement in the second part. My domain is a rectangle (in the fist part) as shown here:

enter image description here
This is the domain for the second part (the same domain for the first part expanded a little vertically):
enter image description here
The red rectangle is the same domain used for the first part. My complete code for the first part is available here:
https://fenicsproject.org/qa/14053/unable-to-evaluate-function-at-point-error
As you can see I have defined my mesh in FEniCS and moved it upward (only in the y-direction). Here is the geo file of the mesh imported for the second part from gmesh:

Point(1) = {0, 0, 0, 1.0};
Point(2) = {0, 2e-6, 0, 1.0};
Point(3) = {0, 12e-6, 0, 1.0};
Point(4) = {0, 192e-6, 0, 1.0};
Point(5) = {0, 202e-6, 0, 1.0};
Point(6) = {0, 204e-6, 0, 1.0};
Point(7) = {0.03, 0, 0, 1.0};
Point(8) = {0.03, 2e-6, 0, 1.0};
Point(9) = {0.03, 12e-6, 0, 1.0};
Point(10) = {0.03, 192e-6, 0, 1.0};
Point(11) = {0.03, 202e-6, 0, 1.0};
Point(12) = {0.03, 204e-6, 0, 1.0};

Line(1) = {1, 7};
Line(2) = {7, 8};
Line(3) = {8, 2};
Line(4) = {2, 1};
Line(5) = {2, 3};
Line(6) = {3, 9};
Line(7) = {9, 8};
Line(8) = {3, 4};
Line(9) = {4, 10};
Line(10) = {10, 9};
Line(11) = {10, 11};
Line(12) = {11, 5};
Line(13) = {5, 4};
Line(14) = {5, 6};
Line(15) = {6, 12};
Line(16) = {12, 11};

Line Loop(17) = {6, -10, -9, -8};
Plane Surface(18) = {17};

Transfinite Surface{18};
Transfinite Line{10,8}=40;
Transfinite Line{9,6}=2000;

Line Loop(19) = {5, 6, 7, 3};
Plane Surface(20) = {19};


Transfinite Line{5,7}=8;
Transfinite Line{3,6}=2000;

Line Loop(21) = {4, 1, 2, 3};
Plane Surface(22) = {21};


Transfinite Line{2,4}=3;
Transfinite Line{1,3}=2000;

Line Loop(23) = {9, 11, 12, 13};
Plane Surface(24) = {23};


Transfinite Line{11,13}=8;
Transfinite Line{9,12}=2000;

Line Loop(25) = {12, 14, 15, 16};
Plane Surface(26) = {25};


Transfinite Line{14,16}=3;
Transfinite Line{12,15}=2000;

For the second part I have divided the domain to some subdomains:

mesh2 = Mesh("mesh2.xml")

class domain_A(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >=  0. and x[1] <=2.*(10**-6) else False

class domain_B(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >=  2.*(10**-6) and x[1] <= 12.*(10**-6) else False

class domain_C(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >=  12.*(10**-6) and x[1] <= 192.*(10**-6) else False

class domain_D(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >=  192.*(10**-6) and x[1] <= 202.0*(10**-6) else False

class domain_E(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >= 202.0*(10**-6) and x[1] <= 204.*(10**-6) else False

#marking subdomains
domains = CellFunction("size_t", mesh2)
domains.set_all(4)

domain_A = domain_A()
domain_A.mark(domains, 1)
domain_B = domain_B()
domain_B.mark(domains, 2)
domain_C = domain_C()
domain_C.mark(domains, 3)
domain_D = domain_D()
domain_D.mark(domains, 4)
domain_E = domain_E()
domain_E.mark(domains, 5)

#lable domain
dx = Measure('dx', domain=mesh2, subdomain_data=domains)

I am using the variable "c" in my variational form only for the domains B,C and D (covering exactly the domain where this variable was calculated on, in the first part). After running the code I get the error:

*** Reason:  The point is not inside the domain...

This is weird because the subdomains B,C,D defined in the second part cover exactly the domain used for calculation of "c" in the first part. I addition the mesh imported from Gmesh seems to be correct.
I have already noticed the problem is related to the boundary between subdomains D and E (upper edge of the rectangular domain of the first part). When I move this boundary down like:

 class domain_D(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >=  192.*(10**-6) and x[1] <= 200.5*(10**-6) else False

class domain_E(SubDomain):
  def inside(self, x, on_boundary):
      return True if x[0] >= 0.0 and x[0] <= 0.03 and x[1] >= 200.5*(10**-6) and x[1] <= 204.*(10**-6) else False

In this case I do not get the error anymore but this boundary should be at exactly :
202.0(10-6) **not 200.5(10**-6)
In the meantime, I can evaluate the variable "c" in the end of the first part of my code:

print c(0.0,202.0 * (10**(-6)) - DOLFIN_EPS)

And it gives me the value of "c" without any error. I do not understand why this variable is not recognized by FEniCS at that small part of the domain in the second part while it is inside the calculation domain. I have tried almost everything like involving the DOLF_EPS in the definition of the subdomains but it was not helpful.Maybe it is related to the moving my mesh coordinate, defining the subdomains in the second part or imported mesh for the second part. Do you know what I am doing wrong here?

asked Jun 21, 2017 by Leo FEniCS Novice (840 points)
...