Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implement the base of a sphere and test it #42

Merged
merged 1 commit into from
Jan 8, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions meshpy/mesh_creation_functions/__init__.py
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
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
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
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