Fem: Search elements for electrostatic constraints

This commit is contained in:
marioalexis
2025-05-07 18:47:09 -03:00
committed by Benjamin Nauck
parent 4fe1002baa
commit eb82923869
2 changed files with 225 additions and 0 deletions

View File

@@ -144,11 +144,14 @@ class MeshSetsGetter:
self.get_constraints_sectionprint_faces()
self.get_constraints_transform_nodes()
self.get_constraints_temperature_nodes()
self.get_constraints_electrostatic_nodes()
# constraints sets with constraint data
self.get_constraints_force_nodeloads()
self.get_constraints_pressure_faces()
self.get_constraints_heatflux_faces()
self.get_constraints_electrostatic_faces()
self.get_constraints_electricchargedensity_faces()
setstime = round((time.process_time() - time_start), 3)
FreeCAD.Console.PrintMessage(f"Getting mesh data time: {setstime} seconds.\n")
@@ -251,6 +254,18 @@ class MeshSetsGetter:
print_obj_info(femobj["Object"])
femobj["Nodes"] = meshtools.get_femnodes_by_femobj_with_references(self.femmesh, femobj)
def get_constraints_electrostatic_nodes(self):
if not self.member.cons_electrostatic:
return
# get nodes
for femobj in self.member.cons_electrostatic:
# femobj --> dict, FreeCAD document object is femobj["Object"]
if femobj["Object"].BoundaryCondition == "Dirichlet":
print_obj_info(femobj["Object"])
femobj["Nodes"] = meshtools.get_femnodes_by_femobj_with_references(
self.femmesh, femobj
)
def get_constraints_force_nodeloads(self):
if not self.member.cons_force:
return
@@ -355,6 +370,64 @@ class MeshSetsGetter:
femobj["PressureFaces"] = [(some_string, pressure_faces)]
FreeCAD.Console.PrintLog("{}\n".format(femobj["PressureFaces"]))
def get_constraints_electrostatic_faces(self):
if not self.member.cons_electrostatic:
return
if not self.femnodes_mesh:
self.femnodes_mesh = self.femmesh.Nodes
if not self.femelement_table:
self.femelement_table = meshtools.get_femelement_table(self.femmesh)
if not self.femnodes_ele_table:
self.femnodes_ele_table = meshtools.get_femnodes_ele_table(
self.femnodes_mesh, self.femelement_table
)
for femobj in self.member.cons_electrostatic:
if femobj["Object"].BoundaryCondition == "Neumann":
print_obj_info(femobj["Object"])
charged_faces = meshtools.get_charge_density_obj_faces(
self.femmesh, self.femelement_table, self.femnodes_ele_table, femobj
)
some_string = "{}: face electric flux".format(femobj["Object"].Name)
femobj["ElectricFluxFaces"] = [(some_string, charged_faces)]
FreeCAD.Console.PrintLog("{}\n".format(femobj["ElectricFluxFaces"]))
def get_constraints_electricchargedensity_faces(self):
if not self.member.cons_electricchargedensity:
return
if not self.femnodes_mesh:
self.femnodes_mesh = self.femmesh.Nodes
if not self.femelement_table:
self.femelement_table = meshtools.get_femelement_table(self.femmesh)
if not self.femnodes_ele_table:
self.femnodes_ele_table = meshtools.get_femnodes_ele_table(
self.femnodes_mesh, self.femelement_table
)
for femobj in self.member.cons_electricchargedensity:
if femobj["Object"].Mode in ["Interface", "Total Interface"]:
print_obj_info(femobj["Object"])
charged_faces = meshtools.get_charge_density_obj_faces(
self.femmesh, self.femelement_table, self.femnodes_ele_table, femobj
)
some_string = "{}: face electric charge density".format(femobj["Object"].Name)
femobj["ChargeDensityFaces"] = [(some_string, charged_faces)]
FreeCAD.Console.PrintLog("{}\n".format(femobj["ChargeDensityFaces"]))
elif femobj["Object"].Mode in ["Source", "Total Source"]:
print_obj_info(femobj["Object"])
charged_volumes = meshtools.get_charge_density_obj_elements(
self.femmesh, self.femelement_table, self.femnodes_ele_table, femobj
)
some_string = "{}: Elements with electric charge density".format(
femobj["Object"].Name
)
femobj["ChargeDensityElements"] = (some_string, charged_volumes)
FreeCAD.Console.PrintLog("{}\n".format(femobj["ChargeDensityElements"]))
def get_constraints_contact_faces(self):
if not self.member.cons_contact:
return

View File

@@ -224,6 +224,7 @@ def get_bit_pattern_dict(femelement_table, femnodes_ele_table, node_set):
The number in the ele_dict is organized as a bit array.
The corresponding bit is set, if the node of the node_set is contained in the element.
"""
# print("BIT PATTERN", femelement_table, femnodes_ele_table, node_set)
FreeCAD.Console.PrintLog("len femnodes_ele_table: " + str(len(femnodes_ele_table)) + "\n")
FreeCAD.Console.PrintLog("len node_set: " + str(len(node_set)) + "\n")
FreeCAD.Console.PrintLog(f"node_set: {node_set}\n")
@@ -241,6 +242,78 @@ def get_bit_pattern_dict(femelement_table, femnodes_ele_table, node_set):
# ************************************************************************************************
def get_ccxelement_volumes_elements_from_binary_search(bit_pattern_dict):
tet10_mask = {0b1111111111: 1}
tet4_mask = {0b1111: 1}
hex8_mask = {0b11111111: 1}
hex20_mask = {0b11111111111111111111: 1}
pent6_mask = {0b111111: 1}
pent15_mask = {0b111111111111111: 1}
vol_dict = {
4: tet4_mask,
6: pent6_mask,
8: hex8_mask,
10: tet10_mask,
15: pent15_mask,
20: hex20_mask,
}
volumes = []
for ele in bit_pattern_dict:
mask_dict = vol_dict[bit_pattern_dict[ele][0]]
for key in mask_dict:
if (key & bit_pattern_dict[ele][1]) == key:
volumes.append(ele)
# print("VOLUMES:", volumes)
FreeCAD.Console.PrintLog(f"found Volumes: {len(volumes)}\n")
# FreeCAD.Console.PrintMessage("faces: {}\n".format(faces))
return volumes
def get_ccxelement_faces_elements_from_binary_search(bit_pattern_dict):
tria3_mask = {0b111: 1}
tria6_mask = {0b111111: 1}
quad4_mask = {0b1111: 1}
quad8_mask = {0b11111111: 1}
vol_dict = {
3: tria3_mask,
6: tria6_mask,
4: quad4_mask,
8: quad8_mask,
}
faces = []
for ele in bit_pattern_dict:
mask_dict = vol_dict[bit_pattern_dict[ele][0]]
for key in mask_dict:
if (key & bit_pattern_dict[ele][1]) == key:
faces.append(ele)
# print("CARAS:", faces)
FreeCAD.Console.PrintMessage(f"found Edges: {len(faces)}\n")
return faces
def get_ccxelement_edges_from_binary_search(bit_pattern_dict):
tria3_mask = {0b011: 1, 0b110: 2, 0b101: 3}
tria6_mask = {0b001011: 1, 0b010110: 2, 0b100101: 3}
quad4_mask = {0b0011: 1, 0b0110: 2, 0b1100: 3, 0b1001: 4}
quad8_mask = {0b00010011: 1, 0b00100110: 2, 0b01001100: 3, 0b10001001: 4}
vol_dict = {
3: tria3_mask,
6: tria6_mask,
4: quad4_mask,
8: quad8_mask,
}
faces = []
for ele in bit_pattern_dict:
mask_dict = vol_dict[bit_pattern_dict[ele][0]]
for key in mask_dict:
if (key & bit_pattern_dict[ele][1]) == key:
faces.append([ele, mask_dict[key]])
# print("EDGES:", faces)
FreeCAD.Console.PrintMessage(f"found Edges: {len(faces)}\n")
return faces
def get_ccxelement_faces_from_binary_search(bit_pattern_dict):
"""get the CalculiX element face numbers"""
# the forum topic discussion with ulrich1a and others ... Better mesh last instead of mesh first
@@ -268,6 +341,7 @@ def get_ccxelement_faces_from_binary_search(bit_pattern_dict):
for key in mask_dict:
if (key & bit_pattern_dict[ele][1]) == key:
faces.append([ele, mask_dict[key]])
# print("FACES:", faces)
FreeCAD.Console.PrintLog(f"found Faces: {len(faces)}\n")
# FreeCAD.Console.PrintMessage("faces: {}\n".format(faces))
return faces
@@ -1427,6 +1501,84 @@ def get_ref_shape_node_sum_geom_table(node_geom_table):
# ************************************************************************************************
# ***** methods for retrieving element face sets *************************************************
# ***** charged faces ****************************************************************************
def get_charge_density_obj_elements(femmesh, femelement_table, femnodes_ele_table, femobj):
node_set = []
res = []
# TODO get elements from mesh groups
# if femmesh.GroupCount:
# node_set = get_femmesh_groupdata_sets_by_name(femmesh, femobj, "Node")
# # FreeCAD.Console.PrintMessage("node_set_group: {}\n".format(node_set))
# if node_set:
# FreeCAD.Console.PrintLog(
# " Finite element mesh nodes where retrieved "
# "from existent finite element mesh group data.\n"
# )
if not node_set:
FreeCAD.Console.PrintLog(
" Finite element mesh nodes will be retrieved "
"by searching the appropriate nodes in the finite element mesh.\n"
)
for feat, ref in femobj["Object"].References:
for sub_ref in ref:
sub = (feat, (sub_ref,))
node_set = get_femnodes_by_references(femmesh, [sub])
charged_volume_node_set = sorted(set(node_set))
bit_pattern_dict = get_bit_pattern_dict(
femelement_table, femnodes_ele_table, charged_volume_node_set
)
sh = feat.getSubObject(sub_ref)
if sh.ShapeType == "Solid":
charged_elem = get_ccxelement_volumes_elements_from_binary_search(
bit_pattern_dict
)
elif sh.ShapeType == "Face":
charged_elem = get_ccxelement_faces_elements_from_binary_search(
bit_pattern_dict
)
res.append((sub, charged_elem))
return res
def get_charge_density_obj_faces(femmesh, femelement_table, femnodes_ele_table, femobj):
node_set = []
res = []
# TODO get elements from mesh groups
# if femmesh.GroupCount:
# node_set = get_femmesh_groupdata_sets_by_name(femmesh, femobj, "Node")
# # FreeCAD.Console.PrintMessage("node_set_group: {}\n".format(node_set))
# if node_set:
# FreeCAD.Console.PrintLog(
# " Finite element mesh nodes where retrieved "
# "from existent finite element mesh group data.\n"
# )
if not node_set:
FreeCAD.Console.PrintLog(
" Finite element mesh nodes will be retrieved "
"by searching the appropriate nodes in the finite element mesh.\n"
)
for feat, ref in femobj["Object"].References:
for sub_ref in ref:
sub = (feat, (sub_ref,))
node_set = get_femnodes_by_references(femmesh, [sub])
charged_face_node_set = sorted(set(node_set))
bit_pattern_dict = get_bit_pattern_dict(
femelement_table, femnodes_ele_table, charged_face_node_set
)
sh = feat.getSubObject(sub_ref)
if sh.ShapeType == "Face":
charged_faces = get_ccxelement_faces_from_binary_search(bit_pattern_dict)
elif sh.ShapeType == "Edge":
charged_faces = get_ccxelement_edges_from_binary_search(bit_pattern_dict)
res.append((sub, charged_faces))
return res
# ***** pressure faces ***************************************************************************
def get_pressure_obj_faces(femmesh, femelement_table, femnodes_ele_table, femobj):
# see get_ccxelement_faces_from_binary_search for more information