Skip to content

Commit

Permalink
Merge pull request #42 from imcs-compsim/implement-base-sphere
Browse files Browse the repository at this point in the history 
Implement the base of a sphere and test it
  • Loading branch information
@isteinbrecher
isteinbrecher authored Jan 8, 2024
2 parents 3173716 + 3aeea4f commit 5935cf7
Show file tree
Hide file tree
Showing 4 changed files with 263 additions and 0 deletions.
2 changes: 2 additions & 0 deletions meshpy/mesh_creation_functions/__init__.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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.
Expand Down Expand Up @@ -99,4 +100,5 @@
"create_nurbs_hollow_cylinder_segment_2d",
"create_nurbs_flat_plate_2d",
"create_nurbs_brick",
"create_nurbs_sphere_surface",
]
122 changes: 122 additions & 0 deletions meshpy/mesh_creation_functions/nurbs_geometries.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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.
Expand Down
110 changes: 110 additions & 0 deletions tests/reference-files/test_nurbs_sphere_surface_reference.dat
View file
Original file line number Diff line number Diff line change
@@ -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
29 changes: 29 additions & 0 deletions tests/testing_nurbs.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,7 @@
create_nurbs_hollow_cylinder_segment_2d,
create_nurbs_flat_plate_2d,
create_nurbs_brick,
create_nurbs_sphere_surface,
)

# Testing imports
Expand Down Expand Up @@ -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.
Expand Down

0 comments on commit 5935cf7

Please sign in to comment.