diff --git a/meshpy/mesh_creation_functions/__init__.py b/meshpy/mesh_creation_functions/__init__.py index 33cde52c..eebdfb2c 100644 --- a/meshpy/mesh_creation_functions/__init__.py +++ b/meshpy/mesh_creation_functions/__init__.py @@ -69,6 +69,7 @@ create_nurbs_hollow_cylinder_segment_2d, create_nurbs_flat_plate_2d, create_nurbs_brick, + create_nurbs_sphere_surface, ) # Define the items that will be exported by default. @@ -99,4 +100,5 @@ "create_nurbs_hollow_cylinder_segment_2d", "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 3fb1786a..370291e6 100644 --- a/meshpy/mesh_creation_functions/nurbs_geometries.py +++ b/meshpy/mesh_creation_functions/nurbs_geometries.py @@ -218,6 +218,128 @@ def create_nurbs_flat_plate_2d(width, length, *, n_ele_u=1, n_ele_v=1): return surf +def create_nurbs_sphere_surface(radius, n_ele_u=1, n_ele_v=1): + """ + Generates a patch of a sphere as a NURBS surface. + This function constructs a segment of a spherical + surface using Non-Uniform Rational B-Splines (NURBS) + based on the specified radius and the number of elements + in the parametric u and v directions. + + Args + --- + radius: double + radius of the sphere + 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 + """ + + # 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 + + dummy = 6.0 * ( + 5.0 / 12.0 + + 0.5 * np.sqrt(2.0 / 3.0) + - 0.25 / np.sqrt(3.0) + - 0.5 * np.sqrt(2.0 / 3.0) * np.sqrt(3.0) / 2.0 + ) + + ctrlpts = [ + [ + 2.0 * radius / np.sqrt(6.0) * -np.sin(1.0 / 4.0 * np.pi), + -radius / np.sqrt(3.0), + 2.0 * radius / np.sqrt(6.0) * np.cos(1.0 / 4.0 * np.pi), + ], + [ + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.cos(1.0 / 4.0 * np.pi) + + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * -np.sin(1.0 / 4.0 * np.pi), + -np.sqrt(3.0) / 2 * radius, + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.cos(1.0 / 4.0 * np.pi) + + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.sin(1.0 / 4.0 * np.pi), + ], + [ + 2.0 * radius / np.sqrt(6.0) * np.cos(1.0 / 4.0 * np.pi), + -radius / np.sqrt(3.0), + 2.0 * radius / np.sqrt(6.0) * np.sin(1.0 / 4.0 * np.pi), + ], + [ + radius * np.sqrt(6.0) / 2.0 * -np.sin(1.0 / 4.0 * np.pi), + 0.0, + radius * np.sqrt(6.0) / 2.0 * np.cos(1.0 / 4.0 * np.pi), + ], + [ + radius * dummy * np.sqrt(2.0) / 2.0 * np.cos(1.0 / 4.0 * np.pi) + + radius * dummy * np.sqrt(2.0) / 2.0 * -np.sin(1.0 / 4.0 * np.pi), + 0.0, + radius * dummy * np.sqrt(2.0) / 2.0 * np.cos(1.0 / 4.0 * np.pi) + + radius * dummy * np.sqrt(2.0) / 2.0 * np.sin(1.0 / 4.0 * np.pi), + ], + [ + radius * np.sqrt(6.0) / 2.0 * np.cos(1.0 / 4.0 * np.pi), + 0.0, + radius * np.sqrt(6.0) / 2.0 * np.sin(1.0 / 4.0 * np.pi), + ], + [ + 2.0 * radius / np.sqrt(6.0) * -np.sin(1.0 / 4.0 * np.pi), + 2.0 * radius / np.sqrt(6.0) * np.cos(1.0 / 4.0 * np.pi), + radius / np.sqrt(3.0), + ], + [ + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.cos(1.0 / 4.0 * np.pi) + + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * -np.sin(1.0 / 4.0 * np.pi), + np.sqrt(3.0) / 2 * radius, + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.cos(1.0 / 4.0 * np.pi) + + radius * np.sqrt(3.0) / (np.sqrt(2.0) * 2.0) * np.sin(1.0 / 4.0 * np.pi), + ], + [ + 2.0 * radius / np.sqrt(6.0) * np.cos(1.0 / 4.0 * np.pi), + radius / np.sqrt(3.0), + 2.0 * radius / np.sqrt(6.0) * np.sin(1.0 / 4.0 * np.pi), + ], + ] + + weights = [ + 1.0, + 2.0 / np.sqrt(6.0), + 1.0, + 2.0 / np.sqrt(6.0), + 2.0 / 3.0, + 2.0 / np.sqrt(6.0), + 1.0, + 2.0 / np.sqrt(6.0), + 1.0, + ] + + t_ctrlptsw = compat.combine_ctrlpts_weights(ctrlpts, weights) + n_ctrlptsw = compat.flip_ctrlpts_u(t_ctrlptsw, p_size_u, p_size_v) + + surf.ctrlpts_size_u = p_size_u + surf.ctrlpts_size_v = p_size_v + surf.ctrlptsw = n_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_brick(width, length, height, *, n_ele_u=1, n_ele_v=1, n_ele_w=1): """ Creates a patch of a 3 dimensional brick. diff --git a/tests/reference-files/test_nurbs_sphere_surface_reference.dat b/tests/reference-files/test_nurbs_sphere_surface_reference.dat new file mode 100644 index 00000000..5a243f0a --- /dev/null +++ b/tests/reference-files/test_nurbs_sphere_surface_reference.dat @@ -0,0 +1,110 @@ +// ----------------------------------------------------------------------------- +// This input file was created with MeshPy. +// Copyright (c) 2018-2023 +// 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 -1.0 NUE 0.0 DENS 0.0 +-----------------------------------------------------------STRUCTURE KNOTVECTORS +BEGIN NURBSPATCH +ID 1 +NUMKNOTS 8 +DEGREE 2 +TYPE Interpolated +0.0 +0.0 +0.0 +0.3333333333333333 +0.6666666666666666 +1.0 +1.0 +1.0 +NUMKNOTS 7 +DEGREE 2 +TYPE Interpolated +0.0 +0.0 +0.0 +0.5 +1.0 +1.0 +1.0 +END NURBSPATCH +-------------------------------------------------------------DNODE-NODE TOPOLOGY +NODE 1 DNODE 1 +NODE 5 DNODE 2 +NODE 16 DNODE 3 +NODE 20 DNODE 4 +-------------------------------------------------------------DLINE-NODE TOPOLOGY +NODE 1 DLINE 1 +NODE 2 DLINE 1 +NODE 3 DLINE 1 +NODE 4 DLINE 1 +NODE 5 DLINE 1 +NODE 5 DLINE 2 +NODE 10 DLINE 2 +NODE 15 DLINE 2 +NODE 20 DLINE 2 +NODE 16 DLINE 3 +NODE 17 DLINE 3 +NODE 18 DLINE 3 +NODE 19 DLINE 3 +NODE 20 DLINE 3 +NODE 1 DLINE 4 +NODE 6 DLINE 4 +NODE 11 DLINE 4 +NODE 16 DLINE 4 +-------------------------------------------------------------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 COORDS +CP 1 COORD -0.57735026919 -0.57735026919 0.57735026919 1.0 +CP 2 COORD -0.409977610553 -0.661036598508 0.661036598508 0.9388321936425754 +CP 3 COORD 0 -0.723160822194 0.723160822194 0.8980536560709591 +CP 4 COORD 0.409977610553 -0.661036598508 0.661036598508 0.9388321936425754 +CP 5 COORD 0.57735026919 -0.57735026919 0.57735026919 1.0 +CP 6 COORD -0.707106781187 -0.317837245196 0.707106781187 0.9082482904638631 +CP 7 COORD -0.50211797591 -0.363907427874 0.896007962233 0.8526927349083075 +CP 8 COORD 0 -0.398107450235 1.03623802516 0.8156556978712706 +CP 9 COORD 0.50211797591 -0.363907427874 0.896007962233 0.8526927349083075 +CP 10 COORD 0.707106781187 -0.317837245196 0.707106781187 0.9082482904638631 +CP 11 COORD -0.707106781187 0.317837245196 0.707106781187 0.9082482904638631 +CP 12 COORD -0.50211797591 0.363907427874 0.896007962233 0.8526927349083075 +CP 13 COORD 0 0.398107450235 1.03623802516 0.8156556978712706 +CP 14 COORD 0.50211797591 0.363907427874 0.896007962233 0.8526927349083075 +CP 15 COORD 0.707106781187 0.317837245196 0.707106781187 0.9082482904638631 +CP 16 COORD -0.57735026919 0.57735026919 0.57735026919 1.0 +CP 17 COORD -0.409977610553 0.661036598508 0.661036598508 0.9388321936425754 +CP 18 COORD 0 0.723160822194 0.723160822194 0.8980536560709591 +CP 19 COORD 0.409977610553 0.661036598508 0.661036598508 0.9388321936425754 +CP 20 COORD 0.57735026919 0.57735026919 0.57735026919 1.0 +--------------------------------------------------------------STRUCTURE ELEMENTS +1 WALLNURBS NURBS9 1 2 3 6 7 8 11 12 13 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +2 WALLNURBS NURBS9 2 3 4 7 8 9 12 13 14 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +3 WALLNURBS NURBS9 3 4 5 8 9 10 13 14 15 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +4 WALLNURBS NURBS9 6 7 8 11 12 13 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 12 13 14 17 18 19 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +6 WALLNURBS NURBS9 8 9 10 13 14 15 18 19 20 MAT 1 KINEM linear EAS none THICK 1.0 STRESS_STRAIN plane_strain GP 3 3 +------------------------------------------------------------------FLUID ELEMENTS +-----------------------------------------------------------------------------END \ No newline at end of file diff --git a/tests/testing_nurbs.py b/tests/testing_nurbs.py index 51262a26..27f0941e 100644 --- a/tests/testing_nurbs.py +++ b/tests/testing_nurbs.py @@ -50,6 +50,7 @@ create_nurbs_hollow_cylinder_segment_2d, create_nurbs_flat_plate_2d, create_nurbs_brick, + create_nurbs_sphere_surface, ) # Testing imports @@ -258,6 +259,34 @@ def test_nurbs_couple_nurbs_meshes(self): # Compare with the reference file compare_test_result(self, input_file.get_string(header=False)) + def test_nurbs_sphere_surface(self): + """Test the creating of the base patch of the surface of a sphere.""" + + # Create input file + input_file = InputFile() + + # Create the base of a sphere + surf_obj = create_nurbs_sphere_surface(1, n_ele_u=3, n_ele_v=2) + + # Create first patch set + mat = MaterialStVenantKirchhoff() + + 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)) + if __name__ == "__main__": # Execution part of script.