From 56cf72359be982bbeffded34fc3ccd7a893cb198 Mon Sep 17 00:00:00 2001 From: Gabriela Loera Date: Thu, 19 Dec 2024 17:47:05 +0100 Subject: [PATCH] Implement a nurbs cylindrical shell sector --- meshpy/mesh_creation_functions/__init__.py | 2 + .../nurbs_geometries.py | 101 +++++++++++ .../test_nurbs_cylindrical_shell_sector.dat | 163 ++++++++++++++++++ tests/testing_nurbs.py | 33 ++++ 4 files changed, 299 insertions(+) create mode 100644 tests/reference-files/test_nurbs_cylindrical_shell_sector.dat diff --git a/meshpy/mesh_creation_functions/__init__.py b/meshpy/mesh_creation_functions/__init__.py index d85baeb2..337a13ab 100644 --- a/meshpy/mesh_creation_functions/__init__.py +++ b/meshpy/mesh_creation_functions/__init__.py @@ -68,6 +68,7 @@ # NURBS geometry functions from .nurbs_geometries import ( create_nurbs_brick, + create_nurbs_cylindrical_shell_sector, create_nurbs_flat_plate_2d, create_nurbs_hemisphere_surface, create_nurbs_hollow_cylinder_segment_2d, @@ -102,6 +103,7 @@ "add_geomdl_nurbs_to_mesh", # NURBS geometry functions "create_nurbs_hollow_cylinder_segment_2d", + "create_nurbs_cylindrical_shell_sector", "create_nurbs_flat_plate_2d", "create_nurbs_brick", "create_nurbs_sphere_surface", diff --git a/meshpy/mesh_creation_functions/nurbs_geometries.py b/meshpy/mesh_creation_functions/nurbs_geometries.py index a5f27d77..b9d87725 100644 --- a/meshpy/mesh_creation_functions/nurbs_geometries.py +++ b/meshpy/mesh_creation_functions/nurbs_geometries.py @@ -149,6 +149,107 @@ def create_nurbs_hollow_cylinder_segment_2d( return surf +def create_nurbs_cylindrical_shell_sector( + radius, angle, length, *, n_ele_u=1, n_ele_v=1 +): + """Creates a patch of a surface of a 3-dimensional sector of a cylindrical + shell. The center of such cylindrical shell sector is located in [0, 0, 0] + + Args + ---- + radius: double + cylindrical shell radius + angle: double + angle of the cylindrical shell (radians) + length: double + length of the cylindrical shell + n_ele_u: int + number of elements in the parametric u-direction + n_ele_v: int + number of elements in the parametric v-direction + + + Return + ---- + surf: geomdl object + geomdl object that contains the surface information + """ + + # Check the validity of the input values: + if (angle >= np.pi) or (angle < 0): + raise ValueError( + "The following algorithm for creating a cylindrical shell sector is only valid for 0 < angle <= pi." + ) + + # Create a NURBS surface instance + surf = NURBS.Surface() + + # Set degrees + surf.degree_u = 2 + surf.degree_v = 2 + + # Control points and set them to the surface + p_size_u = 3 + p_size_v = 3 + + # Obtaining the control points + cp_1 = [-radius * np.sin(angle / 2), -length / 2, -radius * np.cos(angle / 2)] + cp_3 = [radius * np.sin(angle / 2), -length / 2, -radius * np.cos(angle / 2)] + + # Calculating position of middle points. This is done by + # obtaining the tangents on the points cp_1 and cp_3 and + # calculating the intersection point between the tangents. + m_1 = -np.tan(-angle / 2) + m_3 = -np.tan(angle / 2) + b_1 = cp_1[2] + m_1 * cp_1[0] + b_3 = cp_3[2] + m_3 * cp_3[0] + + inter_point_x = (b_3 - b_1) / (m_3 - m_1) + inter_point_z = m_1 * inter_point_x + b_1 + + # The intersection point is assigned to the middle points cp_2, cp_5 and cp_8 + cp_2 = [inter_point_x, -length / 2, inter_point_z] + + cp_4 = [-radius * np.sin(angle / 2), 0.0, -radius * np.cos(angle / 2)] + cp_5 = [inter_point_x, 0.0, inter_point_z] + cp_6 = [radius * np.sin(angle / 2), 0.0, -radius * np.cos(angle / 2)] + + cp_7 = [ + -radius * np.sin(angle / 2), + length / 2, + -radius * np.cos(angle / 2), + ] + cp_8 = [inter_point_x, length / 2, inter_point_z] + cp_9 = [radius * np.sin(angle / 2), length / 2, -radius * np.cos(angle / 2)] + + ctrlpts = [cp_1, cp_4, cp_7, cp_2, cp_5, cp_8, cp_3, cp_6, cp_9] + + weights = [ + 1.0, + 1.0, + 1.0, + np.cos(angle / 2), + np.cos(angle / 2), + np.cos(angle / 2), + 1.0, + 1.0, + 1.0, + ] + + t_ctrlptsw = compat.combine_ctrlpts_weights(ctrlpts, weights) + + surf.ctrlpts_size_u = p_size_u + surf.ctrlpts_size_v = p_size_v + surf.ctrlptsw = t_ctrlptsw + + surf.knotvector_u = [0.0, 0.0, 0.0, 1.0, 1.0, 1.0] + surf.knotvector_v = [0.0, 0.0, 0.0, 1.0, 1.0, 1.0] + + do_uniform_knot_refinement_surface(surf, n_ele_u, n_ele_v) + + return surf + + def create_nurbs_flat_plate_2d(width, length, *, n_ele_u=1, n_ele_v=1): """Creates a patch of a 2 dimensional flat plate. diff --git a/tests/reference-files/test_nurbs_cylindrical_shell_sector.dat b/tests/reference-files/test_nurbs_cylindrical_shell_sector.dat new file mode 100644 index 00000000..9e061469 --- /dev/null +++ b/tests/reference-files/test_nurbs_cylindrical_shell_sector.dat @@ -0,0 +1,163 @@ +// ----------------------------------------------------------------------------- +// This input file was created with MeshPy. +// Copyright (c) 2018-2024 +// Ivo Steinbrecher +// Institute for Mathematics and Computer-Based Simulation +// Universitaet der Bundeswehr Muenchen +// https://www.unibw.de/imcs-en +// ----------------------------------------------------------------------------- +-----------------------------------------------------------------------MATERIALS +MAT 1 MAT_Struct_StVenantKirchhoff YOUNG 50 NUE 0.19 DENS 5.3e-07 +-----------------------------------------------------------STRUCTURE KNOTVECTORS +NURBS_DIMENSION 2 +BEGIN NURBSPATCH +ID 1 +NUMKNOTS 9 +DEGREE 2 +TYPE Interpolated +0.0 +0.0 +0.0 +0.25 +0.5 +0.75 +1.0 +1.0 +1.0 +NUMKNOTS 8 +DEGREE 2 +TYPE Interpolated +0.0 +0.0 +0.0 +0.3333333333333333 +0.6666666666666666 +1.0 +1.0 +1.0 +END NURBSPATCH +-------------------------------------------------------------DNODE-NODE TOPOLOGY +NODE 30 DNODE 1 +NODE 6 DNODE 2 +NODE 25 DNODE 3 +NODE 1 DNODE 4 +-------------------------------------------------------------DLINE-NODE TOPOLOGY +NODE 6 DLINE 1 +NODE 12 DLINE 1 +NODE 18 DLINE 1 +NODE 24 DLINE 1 +NODE 30 DLINE 1 +NODE 5 DLINE 2 +NODE 11 DLINE 2 +NODE 17 DLINE 2 +NODE 23 DLINE 2 +NODE 29 DLINE 2 +NODE 1 DLINE 3 +NODE 7 DLINE 3 +NODE 13 DLINE 3 +NODE 19 DLINE 3 +NODE 25 DLINE 3 +NODE 2 DLINE 4 +NODE 8 DLINE 4 +NODE 14 DLINE 4 +NODE 20 DLINE 4 +NODE 26 DLINE 4 +NODE 25 DLINE 5 +NODE 26 DLINE 5 +NODE 27 DLINE 5 +NODE 28 DLINE 5 +NODE 29 DLINE 5 +NODE 30 DLINE 5 +NODE 19 DLINE 6 +NODE 20 DLINE 6 +NODE 21 DLINE 6 +NODE 22 DLINE 6 +NODE 23 DLINE 6 +NODE 24 DLINE 6 +NODE 1 DLINE 7 +NODE 2 DLINE 7 +NODE 3 DLINE 7 +NODE 4 DLINE 7 +NODE 5 DLINE 7 +NODE 6 DLINE 7 +NODE 7 DLINE 8 +NODE 8 DLINE 8 +NODE 9 DLINE 8 +NODE 10 DLINE 8 +NODE 11 DLINE 8 +NODE 12 DLINE 8 +-------------------------------------------------------------DSURF-NODE TOPOLOGY +NODE 1 DSURFACE 1 +NODE 2 DSURFACE 1 +NODE 3 DSURFACE 1 +NODE 4 DSURFACE 1 +NODE 5 DSURFACE 1 +NODE 6 DSURFACE 1 +NODE 7 DSURFACE 1 +NODE 8 DSURFACE 1 +NODE 9 DSURFACE 1 +NODE 10 DSURFACE 1 +NODE 11 DSURFACE 1 +NODE 12 DSURFACE 1 +NODE 13 DSURFACE 1 +NODE 14 DSURFACE 1 +NODE 15 DSURFACE 1 +NODE 16 DSURFACE 1 +NODE 17 DSURFACE 1 +NODE 18 DSURFACE 1 +NODE 19 DSURFACE 1 +NODE 20 DSURFACE 1 +NODE 21 DSURFACE 1 +NODE 22 DSURFACE 1 +NODE 23 DSURFACE 1 +NODE 24 DSURFACE 1 +NODE 25 DSURFACE 1 +NODE 26 DSURFACE 1 +NODE 27 DSURFACE 1 +NODE 28 DSURFACE 1 +NODE 29 DSURFACE 1 +NODE 30 DSURFACE 1 +---------------------------------------------------------------------NODE COORDS +CP 1 COORD -0.5 -0.5 -0.866025403784 1.0 +CP 2 COORD -0.38799538113 -0.5 -0.930691300639 0.9665063509461097 +CP 3 COORD -0.133974596216 -0.5 -1 0.9330127018922194 +CP 4 COORD 0.133974596216 -0.5 -1 0.9330127018922194 +CP 5 COORD 0.38799538113 -0.5 -0.930691300639 0.9665063509461097 +CP 6 COORD 0.5 -0.5 -0.866025403784 1.0 +CP 7 COORD -0.5 -0.333333333333 -0.866025403784 1.0 +CP 8 COORD -0.38799538113 -0.333333333333 -0.930691300639 0.9665063509461098 +CP 9 COORD -0.133974596216 -0.333333333333 -1 0.9330127018922194 +CP 10 COORD 0.133974596216 -0.333333333333 -1 0.9330127018922194 +CP 11 COORD 0.38799538113 -0.333333333333 -0.930691300639 0.9665063509461098 +CP 12 COORD 0.5 -0.333333333333 -0.866025403784 1.0 +CP 13 COORD -0.5 0 -0.866025403784 1.0 +CP 14 COORD -0.38799538113 0 -0.930691300639 0.9665063509461098 +CP 15 COORD -0.133974596216 0 -1 0.9330127018922194 +CP 16 COORD 0.133974596216 0 -1 0.9330127018922194 +CP 17 COORD 0.38799538113 0 -0.930691300639 0.9665063509461098 +CP 18 COORD 0.5 0 -0.866025403784 1.0 +CP 19 COORD -0.5 0.333333333333 -0.866025403784 1.0 +CP 20 COORD -0.38799538113 0.333333333333 -0.930691300639 0.9665063509461097 +CP 21 COORD -0.133974596216 0.333333333333 -1 0.9330127018922194 +CP 22 COORD 0.133974596216 0.333333333333 -1 0.9330127018922194 +CP 23 COORD 0.38799538113 0.333333333333 -0.930691300639 0.9665063509461097 +CP 24 COORD 0.5 0.333333333333 -0.866025403784 1.0 +CP 25 COORD -0.5 0.5 -0.866025403784 1.0 +CP 26 COORD -0.38799538113 0.5 -0.930691300639 0.9665063509461097 +CP 27 COORD -0.133974596216 0.5 -1 0.9330127018922194 +CP 28 COORD 0.133974596216 0.5 -1 0.9330127018922194 +CP 29 COORD 0.38799538113 0.5 -0.930691300639 0.9665063509461097 +CP 30 COORD 0.5 0.5 -0.866025403784 1.0 +--------------------------------------------------------------STRUCTURE ELEMENTS +1 WALLNURBS NURBS9 1 2 3 7 8 9 13 14 15 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +2 WALLNURBS NURBS9 2 3 4 8 9 10 14 15 16 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +3 WALLNURBS NURBS9 3 4 5 9 10 11 15 16 17 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +4 WALLNURBS NURBS9 4 5 6 10 11 12 16 17 18 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +5 WALLNURBS NURBS9 7 8 9 13 14 15 19 20 21 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +6 WALLNURBS NURBS9 8 9 10 14 15 16 20 21 22 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +7 WALLNURBS NURBS9 9 10 11 15 16 17 21 22 23 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +8 WALLNURBS NURBS9 10 11 12 16 17 18 22 23 24 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +9 WALLNURBS NURBS9 13 14 15 19 20 21 25 26 27 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +10 WALLNURBS NURBS9 14 15 16 20 21 22 26 27 28 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +11 WALLNURBS NURBS9 15 16 17 21 22 23 27 28 29 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +12 WALLNURBS NURBS9 16 17 18 22 23 24 28 29 30 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 diff --git a/tests/testing_nurbs.py b/tests/testing_nurbs.py index 6a0b9748..76dde2a1 100644 --- a/tests/testing_nurbs.py +++ b/tests/testing_nurbs.py @@ -37,6 +37,7 @@ from meshpy import ( InputFile, + InputSection, MaterialString, MaterialStVenantKirchhoff, Rotation, @@ -44,6 +45,7 @@ from meshpy.mesh_creation_functions import ( add_geomdl_nurbs_to_mesh, create_nurbs_brick, + create_nurbs_cylindrical_shell_sector, create_nurbs_flat_plate_2d, create_nurbs_hemisphere_surface, create_nurbs_hollow_cylinder_segment_2d, @@ -86,6 +88,37 @@ def test_nurbs_hollow_cylinder_segment_2d(self): # Compare with the reference file compare_test_result(self, input_file.get_string(header=False)) + def test_nurbs_cylindrical_shell_sector(self): + """Test the creation of a 3-dimensional cylindrical shell sector.""" + + # Create the surface of a quarter of a hollow cylinder + surf_obj = create_nurbs_cylindrical_shell_sector( + 1, np.pi / 3, 1, n_ele_u=4, n_ele_v=3 + ) + + # Create input file + input_file = InputFile() + + # Add material + mat = MaterialStVenantKirchhoff(youngs_modulus=50, nu=0.19, density=5.3e-7) + + # Create patch set + element_description = ( + "KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3" + ) + + patch_set = add_geomdl_nurbs_to_mesh( + input_file, + surf_obj, + material=mat, + element_description=element_description, + ) + + input_file.add(patch_set) + + # Compare with the reference file + compare_test_result(self, input_file.get_string(header=False)) + def test_nurbs_flat_plate_2d(self): """Test the creation of a two dimensional flat plate."""