### How to create an array of boundary conditions to later use in an if statement for applying analytic BCs

203
views
0
7 months ago by
Hey everyone I am trying to use Dirichlet BCs on two different boundaries but apply them to a single array

BCs:
def boundary_x0(x, on_boundary):
return near(x[0],x_min) or near(x[0],x_max)

def boundary_x1(x, on_boundary):
return near(x[1],y_min) or near(x[1],y_max)

​

I then apply the BC values to an array bc_array

id = DirichletBC(V, Constant(1.0), boundary_x0)
id = DirichletBC(V, Constant(2.0), boundary_x1)

bc_id = Function(V)
bc_id.vector()[:] = 0.0
id.apply(bc_id.vector())
bc_array = np.array(bc_id.vector()[:])
print('bc', bc_array)
But when I go to print the array it is always the second condition applied:
(fenicsproject) -bash-4.1$python BC_neumann.py bc [ 2. 2. 2. 2. 2. 2. 2. 2.] ​ I need these values (1 or 2) to use in an if statement in a separate c function to calculate the analytic heat flux on the boundary to impose adiabatic conditions. Any help would be greatly appreciated! Best, Terrence Community: FEniCS Project ### 1 Answer 3 7 months ago by You basically reassign what's stored under the name id with these two lines: id = DirichletBC(V, Constant(1.0), boundary_x0) id = DirichletBC(V, Constant(2.0), boundary_x1) The first BC never has a chance to be used. Just give your second BC a different name. MWE: from fenics import * import numpy as np mesh = UnitSquareMesh(4,4) V = FunctionSpace(mesh,FiniteElement('Lagrange',mesh.ufl_cell(),1)) x_min = 0. x_max = 1. y_min = 0. y_max = 1. def boundary_x0(x, on_boundary): return near(x[0],x_min) or near(x[0],x_max) def boundary_x1(x, on_boundary): return near(x[1],y_min) or near(x[1],y_max) id0 = DirichletBC(V, Constant(1.0), boundary_x0) id1 = DirichletBC(V, Constant(2.0), boundary_x1) bc_id = Function(V) id0.apply(bc_id.vector()) id1.apply(bc_id.vector()) bc_array = np.array(bc_id.vector()[:]) print(bc_array) Output: [ 2. 1. 2. 1. 0. 2. 1. 0. 0. 2. 2. 0. 0. 0. 2. 2. 0. 0. 1. 2. 0. 1. 2. 1. 2.] Two notes: 1. The statement bc_id.vector()[:] = 0.0 is redundant, since bc_id = Function(V) already initializes to an all zero field. 2. Maybe refrain from using python reserved keywords such as id as variable names Thank you for the response I copied your code: #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 #define number of plot steps plot_steps=50 #create mesh mesh = BoxMesh(Point(x_min, y_min, z_min), Point(x_max, y_max, z_max), nx, ny, nz) # define function space V = FunctionSpace(mesh, 'P', 1) # Create boundary markers def boundary_bottom(x, on_boundary): return near(x[2],z_min) def boundary_x0(x, on_boundary): return near(x[0],x_min) or near(x[0],x_max) def boundary_x1(x, on_boundary): return near(x[1],y_min) or near(x[1],y_max) id0 = DirichletBC(V, Constant(1.0), boundary_x0) id1 = DirichletBC(V, Constant(2.0), boundary_x1) bc_id = Function(V) id0.apply(bc_id.vector()) id1.apply(bc_id.vector()) bc_array = np.array(bc_id.vector()[:]) print('bc', bc_array) # Initial condition and body forces ​ I still get the same thing: (fenicsproject) -bash-4.1$ python BC_neumann.py
bc [ 2.  2.  2.  2.  2.  2.  2.  2.]

I really don't know what I am doing wrong except that I am in 3D and have a coarser mesh

written 7 months ago by Terrence
1
I can't try it out myself at the moment, but try to increase your element numbers.
All your points lie in both your boundary definitions, since your BoxMesh only consists of 8 points.
written 7 months ago by klunkean