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

171
views
0
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)

face_numbers=mesh.num_faces()
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()
norm_x0.append(x)
norm_x1.append(y)
norm_x2.append(z)
# calculate areas
area = cell.volume()
area_array.append(area)
# calc centroid coordinates
cent_x0 = cell.midpoint().x()
cent_x1 = cell.midpoint().y()
cent_x2 = cell.midpoint().z()
coord_x0.append(cent_x0)
coord_x1.append(cent_x1)
coord_x2.append(cent_x2)
#get coordinates of vertex of bmesh
h.append(cell.get_vertex_coordinates())

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
​

Terrece

Community: FEniCS Project

0
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()
else:
print("Next vertices combination")