How to use Boundarymesh(mesh,'interior' or 'local')

6 months ago by
Hey everyone. I am working with a 1x1x1 cube. I am currently using boundarymesh(mesh,'exterior') to calculate the xyz components of the normal,centroid,each node, and the cell area of each 2-D cell on the boundary. The problem is I only need this information for the sides at x[0]=0,1 and x[1]=0,1. I read in the documentation that there are options for 'interior' and 'local' but I am unsure how to implement these and if this is what I want. Here is a short working example:
from dolfin import *
import time
import numpy as np
import neumann

#define domain
x_min = y_min = z_min = 0.0
x_max = y_max = 1.0
z_max = 1.98/5.0

#define number of elements
nx = ny =1
nz =1

#create mesh
mesh = BoxMesh(Point(x_min, y_min, z_min), Point(x_max, y_max, z_max), nx, ny, nz)
coordinates = mesh.coordinates()
# define function space
V = FunctionSpace(mesh, 'P', 1)

u = TrialFunction(V)
v = TestFunction(V)
u0 = Function(V)

bmesh = BoundaryMesh(mesh, 'exterior')
dof_mesh = V.tabulate_dof_coordinates().reshape(V.dim(),mesh.geometry().dim())
h = []
norm_x0 = []
norm_x1 = []
norm_x2 = []
coord_x0 = []
coord_x1 = []
coord_x2 = []
area_array = []
msg = 'Cell %d in bmesh(%d in mesh) has normal %r'
for bci, cell in enumerate(cells(bmesh)):
    # calculate normal vectors
    x = cell.cell_normal().x()
    y = cell.cell_normal().y()
    z = cell.cell_normal().z()
    # calculate areas
    area = cell.volume()
    # calc centroid coordinates
    cent_x0 = cell.midpoint().x()
    cent_x1 = cell.midpoint().y()
    cent_x2 = cell.midpoint().z()
    #get coordinates of vertex of bmesh

h = np.array(h)
coord_x0 = np.array(coord_x0)
coord_x1 = np.array(coord_x1)
coord_x2 = np.array(coord_x2)
norm_x0 = np.array(norm_x0)
norm_x1 = np.array(norm_x1)
norm_x2 = np.array(norm_x2)
area_array = np.array(area_array)
#then a bunch of stuff that is not relevant 

Thanks in advance,

Community: FEniCS Project

1 Answer

6 months ago by
Hi Terrence,

you can find (if) and which cell center exists at (or near) x=0.1 and y=0.1 with a simple if statement. Then you can get the cell index in order to use it later. Try something like this ( I haven't tried it):

for cell in cells(bmesh):
    if cell.midpoint().x() = 0.1 and cell.midpoint().y() = 0.1:
        c = cell.index()
        x = cell.cell_normal.x()
        y = cell.cell_normal.y()
        z = cell.cell_normal.z()
        area = cell.volume()
        h = cell.get_vertex_coordinates()

If you are not sure if the center of the cell lies on (0.1, 0.1, 0.0) then you should use as criterion the fact that in order a cell to lie near that point, then at least one vertex x-coordinate has to be greater than 0.1 and at least one vertex x-coordinate has to be smaller than 0.1. The same holds for y-coordinate.

for cell in cells(bmesh):
    for vertex in vertices(cell):
         x = vertex.x()
         y = vertex.y()
         if x > 0.1 and y < 0.1:
         c = cell.index()
         if x < 0.1 and y > 0.1:
         c = cell.index()
         print("Next vertices combination")
Please login to add an answer/comment or follow this question.

Similar posts:
Search »