Dear all,
this question is related to boundary layers.
I need to write into a file the hex element id and which side of the hex element is touching a defined boundary layer in Cubit. I defined for this reason a side boundary layer.
I assume, when reading the side numbering chapter in the help, that it should be possible to get the information which of the 6 possilble sides of a HEX elemet does touch a boundary layer. Is this interpretation correct?
How can I write this into an external file with phyton?
If I understand the question correctly this is non-trivial. You need to get the face id associated with a set of nodes in a hex. Here is my solution to your problem. I added a link to the Exodus II documentation. I always have to refer to this documentation when I am working at this level of the mesh model.
This code has been minimally tested. It should at least give you some idea of what needs to be done
#!python
# The face numbering is based on the Exodus element definition which was in turn based on Patran
# See for example, https://gsjaardema.github.io/seacas-docs/html/element_types.html#hex
from collections import deque
# the key is a tuple of nodes in a HEX8 element and the value is the
# face id in the Exodus II definition
face_number_dict = { (0,1,5,4) : 1, \
(1,2,6,5) : 2, \
(2,3,7,6) : 3, \
(0,4,7,3) : 4, \
(0,3,2,1) : 5, \
(4,5,6,7) : 6, \
}
def minimum_index(conn):
"""
Description: Find the index of the given node
Return: the index or -1
"""
min_id = min(conn)
min_index = [i for i, j in enumerate(conn) if j == min_id]
return min_index[0]
def rotate_connectivity(conn):
min_index = minimum_index(conn)
conn_deque = deque(conn)
conn_deque.rotate(-min_index)
return conn_deque
def get_side_from_hex( hex_id, quad_id):
"""
Description: Given the id of a hex and the id of a quad determine
the hex face id associated with that quad. If the quad is not
in the hex, return -1
"""
hex_conn = cubit.get_connectivity("hex", hex_id)
quad_conn = cubit.get_connectivity("quad", quad_id)
hex_rotated = rotate_connectivity( hex_conn )
quad_rotated = rotate_connectivity( quad_conn )
for fn in face_number_dict:
if quad_rotated[ 0 ] == hex_rotated[ fn[0] ] and \
quad_rotated[ 1 ] == hex_rotated[ fn[1] ] and \
quad_rotated[ 2 ] == hex_rotated[ fn[2] ] and \
quad_rotated[ 3 ] == hex_rotated[ fn[3] ]:
return face_number_dict[fn]
return -1
cubit.cmd("reset")
cubit.cmd("brick x 10")
cubit.cmd("mesh vol 1")
cubit.cmd("sideset 2 face in surface 1")
hex_id = 1
dimension = 2
sideset_id = 2
hex_quads = cubit.get_sub_elements( "hex", hex_id, dimension )
sideset_quads = cubit.get_sideset_quads( sideset_id )
quad_in_sideset = [quad in sideset_quads for quad in hex_quads]
if any( quad_in_sideset ):
print( f"hex {hex_id} is in sideset {sideset_id}" )
true_index = [i for i, j in enumerate(quad_in_sideset) if j ][0]
print( f"Face id = {get_side_from_hex( hex_id, hex_quads[true_index])}")
One other thing to note. I made the sideset easier for myself by adding the faces into the sideset instead of the surface. If you add the surface, you will have to use cubit.get_geometric_owner(“quad”, quad_id) to first find out if the quad is in the surface that is in the sideset.
Hello karl,
thank you so far for your support, it is indeed not trivial.
Your recommendation covers my needs, but the posted example does not give the quads or faces on a sideset back of the example (brick).
I’m a bit confused how to solve the problem. Could you please post a code where I get the HEX id and the face number of a hex touching a sideset. The face number dict is already what I’m looking for.
Many thanks in advance.
BR Stefan
If I understand your request, I started with the hex id as a given and you want to find the list of hexes associated with the sideset. One way to trim the set of hex ids is to use sideset_id = 2 hex_in_sideset = cubit.parse_cubit_list("hex", f"in node in face in sideset {sideset_id}")
Is that the other piece you were looking for? If you use surfaces in the sideset change the above
to “surface” instead of “face”.
