From 7760c5cce1fc11ce70c0382fc03bdbdfba8c55c1 Mon Sep 17 00:00:00 2001 From: Bernd Hahnebach Date: Mon, 1 Aug 2016 21:56:57 +0100 Subject: [PATCH] FEM: constraint plane rotation: add implementation for solver CalculiX --- src/Mod/Fem/FemInputWriter.py | 15 ++++++++- src/Mod/Fem/FemInputWriterCcx.py | 56 ++++++++++++++++++++++++++++++-- src/Mod/Fem/FemInputWriterZ88.py | 4 +-- src/Mod/Fem/FemMeshTools.py | 44 +++++++++++++++++++++++++ src/Mod/Fem/FemTools.py | 8 +++++ src/Mod/Fem/FemToolsCcx.py | 2 +- src/Mod/Fem/FemToolsZ88.py | 2 +- 7 files changed, 124 insertions(+), 7 deletions(-) diff --git a/src/Mod/Fem/FemInputWriter.py b/src/Mod/Fem/FemInputWriter.py index 39e6eb373..242843bf1 100644 --- a/src/Mod/Fem/FemInputWriter.py +++ b/src/Mod/Fem/FemInputWriter.py @@ -45,7 +45,7 @@ class FemInputWriter(): analysis_obj, solver_obj, mesh_obj, mat_obj, fixed_obj, displacement_obj, - contact_obj, + contact_obj, planerotation_obj, selfweight_obj, force_obj, pressure_obj, beamsection_obj, shellthickness_obj, analysis_type, eigenmode_parameters, @@ -58,6 +58,7 @@ class FemInputWriter(): self.fixed_objects = fixed_obj self.displacement_objects = displacement_obj self.contact_objects = contact_obj + self.planerotation_objects = planerotation_obj self.selfweight_objects = selfweight_obj self.force_objects = force_obj self.pressure_objects = pressure_obj @@ -80,16 +81,28 @@ class FemInputWriter(): self.femmesh = self.mesh_object.FemMesh self.femnodes_mesh = {} self.femelement_table = {} + self.constraint_conflict_nodes = [] def get_constraints_fixed_nodes(self): # get nodes for femobj in self.fixed_objects: # femobj --> dict, FreeCAD document object is femobj['Object'] femobj['Nodes'] = FemMeshTools.get_femnodes_by_references(self.femmesh, femobj['Object'].References) + # add nodes to constraint_conflict_nodes, needed by constraint plane rotation + for node in femobj['Nodes']: + self.constraint_conflict_nodes.append(node) def get_constraints_displacement_nodes(self): # get nodes for femobj in self.displacement_objects: # femobj --> dict, FreeCAD document object is femobj['Object'] femobj['Nodes'] = FemMeshTools.get_femnodes_by_references(self.femmesh, femobj['Object'].References) + # add nodes to constraint_conflict_nodes, needed by constraint plane rotation + for node in femobj['Nodes']: + self.constraint_conflict_nodes.append(node) + + def get_constraints_planerotation_nodes(self): + # get nodes + for femobj in self.planerotation_objects: # femobj --> dict, FreeCAD document object is femobj['Object'] + femobj['Nodes'] = FemMeshTools.get_femnodes_by_references(self.femmesh, femobj['Object'].References) def get_constraints_force_nodeloads(self): # check shape type of reference shape diff --git a/src/Mod/Fem/FemInputWriterCcx.py b/src/Mod/Fem/FemInputWriterCcx.py index 360c11a78..a6f4702b4 100644 --- a/src/Mod/Fem/FemInputWriterCcx.py +++ b/src/Mod/Fem/FemInputWriterCcx.py @@ -40,7 +40,7 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): analysis_obj, solver_obj, mesh_obj, mat_obj, fixed_obj, displacement_obj, - contact_obj, + contact_obj, planerotation_obj, selfweight_obj, force_obj, pressure_obj, beamsection_obj, shellthickness_obj, analysis_type=None, eigenmode_parameters=None, @@ -50,7 +50,7 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): analysis_obj, solver_obj, mesh_obj, mat_obj, fixed_obj, displacement_obj, - contact_obj, + contact_obj, planerotation_obj, selfweight_obj, force_obj, pressure_obj, beamsection_obj, shellthickness_obj, analysis_type, eigenmode_parameters, @@ -71,10 +71,14 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): self.write_node_sets_constraints_fixed(inpfile) if self.displacement_objects: self.write_node_sets_constraints_displacement(inpfile) + if self.planerotation_objects: + self.write_node_sets_constraints_planerotation(inpfile) if self.contact_objects: self.write_surfaces_contraints_contact(inpfile) self.write_materials(inpfile) self.write_femelementsets(inpfile) + if self.planerotation_objects: + self.write_constraints_planerotation(inpfile) if self.contact_objects: self.write_constraints_contact(inpfile) self.write_step_begin(inpfile) @@ -93,6 +97,7 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): self.write_analysis_frequency(inpfile) self.write_outputs_types(inpfile) self.write_step_end(inpfile) + self.write_footer(inpfile) inpfile.close() return self.file_name @@ -158,6 +163,43 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): for n in femobj['Nodes']: f.write(str(n) + ',\n') + def write_node_sets_constraints_planerotation(self, f): + # get nodes + self.get_constraints_planerotation_nodes() + # write nodes to file + if not self.femnodes_mesh: + self.femnodes_mesh = self.femmesh.Nodes + f.write('\n***********************************************************\n') + f.write('** Node set for plane rotation constraint\n') + f.write('** written by {} function\n'.format(sys._getframe().f_code.co_name)) + # info about self.constraint_conflict_nodes: + # is used to check if MPC and constraint fixed and constraint displacement share same nodes, + # because MPC's and constriants fixed an constraints displacement can't share same nodes. + # thus call write_node_sets_constraints_planerotation has to be after constraint fixed and constraint displacement + for femobj in self.planerotation_objects: # femobj --> dict, FreeCAD document object is femobj['Object'] + l_nodes = femobj['Nodes'] + fric_obj = femobj['Object'] + f.write('*NSET,NSET=' + fric_obj.Name + '\n') + # Code to extract nodes and coordinates on the PlaneRotation support face + nodes_coords = [] + for node in l_nodes: + nodes_coords.append((node, self.femnodes_mesh[node].x, self.femnodes_mesh[node].y, self.femnodes_mesh[node].z)) + node_planerotation = FemMeshTools.get_three_non_colinear_nodes(nodes_coords) + for i in range(len(l_nodes)): + if l_nodes[i] not in node_planerotation: + node_planerotation.append(l_nodes[i]) + MPC_nodes = [] + for i in range(len(node_planerotation)): + cnt = 0 + for j in range(len(self.constraint_conflict_nodes)): + if node_planerotation[i] == self.constraint_conflict_nodes[j]: + cnt = cnt + 1 + if cnt == 0: + MPC = node_planerotation[i] + MPC_nodes.append(MPC) + for i in range(len(MPC_nodes)): + f.write(str(MPC_nodes[i]) + ',\n') + def write_surfaces_contraints_contact(self, f): # get surface nodes and write them to file f.write('\n***********************************************************\n') @@ -338,6 +380,16 @@ class FemInputWriterCcx(FemInputWriter.FemInputWriter): stick = (slope / 10.0) f.write(str(friction) + ', ' + str(stick) + ' \n') + def write_constraints_planerotation(self, f): + f.write('\n***********************************************************\n') + f.write('** PlaneRotation Constraints\n') + f.write('** written by {} function\n'.format(sys._getframe().f_code.co_name)) + for femobj in self.planerotation_objects: # femobj --> dict, FreeCAD document object is femobj['Object'] + fric_obj_name = femobj['Object'].Name + f.write('*MPC\n') + f.write('PLANE,' + fric_obj_name + '\n') + f.write('\n') + def write_constraints_selfweight(self, f): f.write('\n***********************************************************\n') f.write('** Self weight\n') diff --git a/src/Mod/Fem/FemInputWriterZ88.py b/src/Mod/Fem/FemInputWriterZ88.py index 9b30e7c73..60a7abfaf 100644 --- a/src/Mod/Fem/FemInputWriterZ88.py +++ b/src/Mod/Fem/FemInputWriterZ88.py @@ -37,7 +37,7 @@ class FemInputWriterZ88(FemInputWriter.FemInputWriter): analysis_obj, solver_obj, mesh_obj, mat_obj, fixed_obj, displacement_obj, - contact_obj, + contact_obj, planerotation_obj, selfweight_obj, force_obj, pressure_obj, beamsection_obj, shellthickness_obj, analysis_type=None, eigenmode_parameters=None, @@ -47,7 +47,7 @@ class FemInputWriterZ88(FemInputWriter.FemInputWriter): analysis_obj, solver_obj, mesh_obj, mat_obj, fixed_obj, displacement_obj, - contact_obj, + contact_obj, planerotation_obj, selfweight_obj, force_obj, pressure_obj, beamsection_obj, shellthickness_obj, analysis_type, eigenmode_parameters, diff --git a/src/Mod/Fem/FemMeshTools.py b/src/Mod/Fem/FemMeshTools.py index 60429ab31..91e64867d 100644 --- a/src/Mod/Fem/FemMeshTools.py +++ b/src/Mod/Fem/FemMeshTools.py @@ -578,6 +578,50 @@ def is_edge_femmesh(femmesh): return True +def get_three_non_colinear_nodes(nodes_coords): + # Code to obtain three non-colinear nodes on the PlaneRotation support face + # nodes_coords --> [(nodenumber, x, y, z), (nodenumber, x, y, z), ...] + if not nodes_coords: + print(len(nodes_coords)) + print('Error: No nodes in nodes_coords') + return [] + dum_max = [1, 2, 3, 4, 5, 6, 7, 8, 0] + for i in range(len(nodes_coords)): + for j in range(len(nodes_coords) - 1 - i): + x_1 = nodes_coords[j][1] + x_2 = nodes_coords[j + 1][1] + y_1 = nodes_coords[j][2] + y_2 = nodes_coords[j + 1][2] + z_1 = nodes_coords[j][3] + z_2 = nodes_coords[j + 1][3] + node_1 = nodes_coords[j][0] + node_2 = nodes_coords[j + 1][0] + distance = ((x_1 - x_2) ** 2 + (y_1 - y_2) ** 2 + (z_1 - z_2) ** 2) ** 0.5 + if distance > dum_max[8]: + dum_max = [node_1, x_1, y_1, z_1, node_2, x_2, y_2, z_2, distance] + node_dis = [1, 0] + for i in range(len(nodes_coords)): + x_1 = dum_max[1] + x_2 = dum_max[5] + x_3 = nodes_coords[i][1] + y_1 = dum_max[2] + y_2 = dum_max[6] + y_3 = nodes_coords[i][2] + z_1 = dum_max[3] + z_2 = dum_max[7] + z_3 = nodes_coords[i][3] + node_3 = int(nodes_coords[j][0]) + distance_1 = ((x_1 - x_3) ** 2 + (y_1 - y_3) ** 2 + (z_1 - z_3) ** 2) ** 0.5 + distance_2 = ((x_3 - x_2) ** 2 + (y_3 - y_2) ** 2 + (z_3 - z_2) ** 2) ** 0.5 + tot = distance_1 + distance_2 + if tot > node_dis[1]: + node_dis = [node_3, tot] + node_1 = int(dum_max[0]) + node_2 = int(dum_max[4]) + print([node_1, node_2, node_3]) + return [node_1, node_2, node_3] + + def make_femmesh(mesh_data): ''' makes an FreeCAD FEM Mesh object from FEM Mesh data ''' diff --git a/src/Mod/Fem/FemTools.py b/src/Mod/Fem/FemTools.py index 8c01bbbde..e6c90967b 100644 --- a/src/Mod/Fem/FemTools.py +++ b/src/Mod/Fem/FemTools.py @@ -190,6 +190,10 @@ class FemTools(QtCore.QRunnable, QtCore.QObject): # set of displacements for the analysis. Updated with update_objects # Individual displacement_constraints are Proxy.Type "FemConstraintDisplacement" self.displacement_constraints = [] + ## @var planerotation_constraints + # set of plane rotation constraints from the analysis. Updated with update_objects + # Individual constraints are "Fem::ConstraintPlaneRotation" type + self.planerotation_constraints = [] ## @var contact_constraints # set of contact constraints from the analysis. Updated with update_objects # Individual constraints are "Fem::ConstraintContact" type @@ -240,6 +244,10 @@ class FemTools(QtCore.QRunnable, QtCore.QObject): displacement_constraint_dict = {} displacement_constraint_dict['Object'] = m self.displacement_constraints.append(displacement_constraint_dict) + elif m.isDerivedFrom("Fem::ConstraintPlaneRotation"): + planerotation_constraint_dict = {} + planerotation_constraint_dict['Object'] = m + self.planerotation_constraints.append(planerotation_constraint_dict) elif m.isDerivedFrom("Fem::ConstraintContact"): contact_constraint_dict = {} contact_constraint_dict['Object'] = m diff --git a/src/Mod/Fem/FemToolsCcx.py b/src/Mod/Fem/FemToolsCcx.py index f61b49379..83ea7f8ad 100644 --- a/src/Mod/Fem/FemToolsCcx.py +++ b/src/Mod/Fem/FemToolsCcx.py @@ -92,7 +92,7 @@ class FemToolsCcx(FemTools.FemTools): self.analysis, self.solver, self.mesh, self.materials, self.fixed_constraints, self.displacement_constraints, - self.contact_constraints, + self.contact_constraints, self.planerotation_constraints, self.selfweight_constraints, self.force_constraints, self.pressure_constraints, self.beam_sections, self.shell_thicknesses, self.analysis_type, self.eigenmode_parameters, diff --git a/src/Mod/Fem/FemToolsZ88.py b/src/Mod/Fem/FemToolsZ88.py index 460a18428..083a28bb6 100644 --- a/src/Mod/Fem/FemToolsZ88.py +++ b/src/Mod/Fem/FemToolsZ88.py @@ -84,7 +84,7 @@ class FemToolsZ88(FemTools.FemTools): self.analysis, self.solver, self.mesh, self.materials, self.fixed_constraints, self.displacement_constraints, - self.contact_constraints, + self.contact_constraints, self.planerotation_constraints, self.selfweight_constraints, self.force_constraints, self.pressure_constraints, self.beam_sections, self.shell_thicknesses, self.analysis_type, None,