Merge pull request #25 from imcs-compsim/add-spline-interpolation-curve

Add vertex spline interpolation function

cubitpy/geometry_creation_functions.py

@@ -36,6 +36,31 @@
 from . import cupy
 
 
+def create_spline_interpolation_curve(cubit, vertices, *, delete_points=True):
+    """
+    Interpolate multiple vertices with a Cubit internal spline.
+
+    Args
+    ----
+    cubit: Cubit
+        Link to the main cubit object.
+    vertices: list(points in R3)
+        Points to be interpolated by the spline curve
+    delete_points: bool
+        If the created vertices should be kept or should be deleted.
+    """
+
+    # Create the vertices
+    vertices = [cubit.create_vertex(*vertex) for vertex in vertices]
+    vertices_ids = [str(vertex.id()) for vertex in vertices]
+    cubit.cmd(
+        "create curve spline vertex {} {}".format(
+            " ".join(vertices_ids), ("delete" if delete_points else "")
+        )
+    )
+    return cubit.curve(cubit.get_last_id(cupy.geometry.curve))
+
+
 def create_parametric_curve(
     cubit,
     f,
@@ -72,17 +97,10 @@ def create_parametric_curve(
         interval[0] + i * (interval[1] - interval[0]) / float(n_segments)
         for i in range(n_segments + 1)
     ]
-    vertices = [
-        cubit.create_vertex(*f(t, *function_args, **function_kwargs))
-        for t in parameter_points
-    ]
-    vertices_ids = [str(vertex.id()) for vertex in vertices]
-    cubit.cmd(
-        "create curve spline vertex {} {}".format(
-            " ".join(vertices_ids), ("delete" if delete_points else "")
-        )
+    vertices = [f(t, *function_args, **function_kwargs) for t in parameter_points]
+    return create_spline_interpolation_curve(
+        cubit, vertices, delete_points=delete_points
     )
-    return cubit.curve(cubit.get_last_id(cupy.geometry.curve))
 
 
 def create_parametric_surface(

tests/testing.py

@@ -47,7 +47,10 @@
 # CubitPy imports.
 from cubitpy import CubitPy, cupy
 from cubitpy.mesh_creation_functions import create_brick, extrude_mesh_normal_to_surface
-from cubitpy.geometry_creation_functions import create_parametric_surface
+from cubitpy.geometry_creation_functions import (
+    create_spline_interpolation_curve,
+    create_parametric_surface,
+)
 from cubitpy.cubit_utility import get_surface_center, import_fluent_geometry
 
 
@@ -1334,6 +1337,51 @@ def f(u, v, arg, kwarg=-1.0):
     assert np.linalg.norm(connectivity - connectivity_ref) == 0
 
 
+def test_spline_interpolation_curve():
+    """
+    Test the create_spline_interpolation_curve function.
+    """
+
+    cubit = CubitPy()
+
+    x = np.linspace(0, 2 * np.pi, 7)
+    y = np.cos(x)
+    z = np.sin(x)
+    vertices = np.array([x, y, z]).transpose()
+
+    curve = create_spline_interpolation_curve(cubit, vertices)
+    curve.mesh()
+
+    coordinates = [
+        cubit.get_nodal_coordinates(i + 1) for i in range(cubit.get_node_count())
+    ]
+    connectivity = [
+        cubit.get_connectivity("edge", i + 1) for i in range(cubit.get_edge_count())
+    ]
+
+    # fmt: off
+    coordinates_ref = np.array([
+        [0.0, 1.0, 0.0],
+        [6.283185307179586, 1.0, -2.4492935982947064e-16], 
+        [0.6219064247387815, 0.7622034923056742, 0.5808964193893371],
+        [1.2706376409420117, 0.30926608007524203, 0.9532391827102926],
+        [1.8922964421051867, -0.3108980458371118, 0.946952808381383],
+        [2.5151234800888007, -0.8099976142632724, 0.5846200862869367],
+        [3.1415926535897927, -0.9999999999999998, 1.6653345369377348e-16],
+        [3.7680618270907873, -0.8099976142632712, -0.5846200862869384],
+        [4.3908888650744, -0.31089804583711017, -0.9469528083813835],
+        [5.012547666237575, 0.30926608007524364, -0.9532391827102922],
+        [5.661278882440805, 0.7622034923056742, -0.5808964193893369]
+    ])
+
+    connectivity_ref = np.array([[1, 3], [3, 4], [4, 5], [5, 6], [6, 7], [7, 8], [8, 9], [9, 10],
+        [10, 11], [11, 2]])
+    # fmt: on
+
+    assert 0.0 == pytest.approx(np.linalg.norm(coordinates - coordinates_ref), 1e-12)
+    assert np.linalg.norm(connectivity - connectivity_ref) == 0
+
+
 def setup_and_check_import_fluent_geometry(
     fluent_geometry, feature_angle, reference_entitys_number
 ):