I added a new function for finding rotation matrices.

daniobot/daniobot.py

@@ -1,3 +1,5 @@
+import math
+
 import pterasoftware as ps
 
 daniobot = ps.geometry.Airplane(
@@ -12,6 +14,7 @@
             y_le=0.0,
             z_le=0.0,
             symmetric=True,
+            unit_chordwise_vector=[1 / math.sqrt(2), 0, 1 / math.sqrt(2)],
             num_chordwise_panels=16,
             chordwise_spacing="uniform",
             wing_cross_sections=[
@@ -54,16 +57,16 @@
 
 tail_mid_wing_cross_section_movement = ps.movement.WingCrossSectionMovement(
     base_wing_cross_section=daniobot.wings[0].wing_cross_sections[1],
-    pitching_amplitude=15.0,
-    pitching_period=1 / 8,
-    pitching_spacing="sine",
+    # pitching_amplitude=15.0,
+    # pitching_period=1 / 8,
+    # pitching_spacing="sine",
 )
 
 tail_tip_wing_cross_section_movement = ps.movement.WingCrossSectionMovement(
     base_wing_cross_section=daniobot.wings[0].wing_cross_sections[2],
-    pitching_amplitude=15.0,
-    pitching_period=1 / 8,
-    pitching_spacing="sine",
+    # pitching_amplitude=15.0,
+    # pitching_period=1 / 8,
+    # pitching_spacing="sine",
 )
 
 tail_movement = ps.movement.WingMovement(

pterasoftware/functions.py

@@ -14,6 +14,11 @@
     angle_axis_rotation_matrix: This function is used to find the rotation matrix for
     a given axis and angle.
 
+    axis_axis_rotation_matrix: This function is used to find the rotation matrix that
+    would rotate a vector in the direction of one axis to the direction of another
+    axis. Both input axes will be normalized, so this rotation matrix will preserve
+    length.
+
     numba_centroid_of_quadrilateral: This function is used to find the centroid of a
     quadrilateral. It has been optimized for JIT compilation using Numba.
 
@@ -50,6 +55,7 @@
 """
 
 import logging
+from audioop import cross
 
 import numpy as np
 from numba import njit
@@ -151,6 +157,49 @@ def angle_axis_rotation_matrix(angle, axis, axis_already_normalized=False):
     return rotation_matrix
 
 
+def axis_to_axis_rotation_matrix(start_axis, end_axis):
+    """This function is used to find the rotation matrix that would rotate a vector
+    in the direction of one axis to the direction of another axis. Both input axes
+    will be normalized, so this rotation matrix will preserve length.
+
+    Citation:
+        Adapted from:         ttps://math.stackexchange.com/q/476311
+        Author:               Jur van den Berg
+        Date of Retrieval:    12/05/2024
+
+    :param start_axis: (3,) array of floats.
+        This is the initial state of the axis.
+    :param end_axis: (3,) array of floats.
+        This is the final state of the axis.
+    :return (3, 3) array of floats
+        This is the rotation matrix that
+    """
+    # Normalize both axes.
+    start_axis = start_axis / np.linalg.norm(start_axis)
+    end_axis = end_axis / np.linalg.norm(end_axis)
+
+    cross_product = np.cross(start_axis, end_axis)
+
+    # Find the cosine of the angle between the two axes in radians.
+    cosine = np.dot(start_axis, end_axis)
+
+    # Find the skew-symmetric matrix version of the cross product.
+    cross_product_skew_sym = np.cross(np.eye(3), cross_product)
+
+    if cosine != -1:
+        # If the vectors aren't anti-parallel, apply the formula, and return the
+        # rotation matrix.
+        return (
+            np.eye(3)
+            + cross_product_skew_sym
+            + (cross_product_skew_sym @ cross_product_skew_sym) / (1 + cosine)
+        )
+    else:
+        # If the vectors are anti-parallel, return the negative Identity matrix,
+        # which is the rotation matrix.
+        return -np.eye(3)
+
+
 @njit(cache=True, fastmath=False)
 def numba_centroid_of_quadrilateral(
     front_left_vertex, front_right_vertex, back_left_vertex, back_right_vertex

pterasoftware/geometry.py

@@ -285,7 +285,15 @@ def __init__(
 
         # Check that the root wing cross section's leading edge isn't offset from the
         # wing's leading edge.
-        if np.any(self.wing_cross_sections[0].leading_edge):
+        if np.any(
+            np.array(
+                [
+                    self.wing_cross_sections[0].x_le,
+                    self.wing_cross_sections[0].y_le,
+                    self.wing_cross_sections[0].z_le,
+                ]
+            )
+        ):
             raise Exception(
                 "The root wing cross section's leading edge must not be offset from "
                 "the wing's leading edge."
@@ -330,9 +338,6 @@ def __init__(
         # Calculate the number of panels on this wing.
         self.num_panels = self.num_spanwise_panels * self.num_chordwise_panels
 
-        for wing_cross_section in self.wing_cross_sections:
-            wing_cross_section.wing_unit_chordwise_vector = self.unit_chordwise_vector
-
         # Initialize the panels attribute. Then mesh the wing, which will populate
         # this attribute.
         self.panels = None
@@ -648,7 +653,7 @@ def unit_up_vector(self):
         """This method defines a property for the wing cross section's unit up vector.
 
         :return: (3,) array of floats
-            This is the unit vector for the wing cross section's chordwise direction.
+            This is the unit vector for the wing cross section's up direction.
             The units are meters.
         """
         return np.cross(self.unit_chordwise_vector, self.unit_normal_vector)

pterasoftware/movement.py

@@ -511,6 +511,9 @@ def __init__(
         z_le_amplitude=0.0,
         z_le_period=0.0,
         z_le_spacing="sine",
+        # wing_pitching_amplitude=0.0,
+        # wing_pitching_period=0.0,
+        # wing_pitching_spacing="sine",
     ):
         """This is the initialization method.
 
@@ -722,7 +725,7 @@ def generate_wings(self, num_steps=10, delta_time=0.1):
 
                 if base_y_le != inner_y_les[0]:
                     # Find the base sweep angle of this wing cross section compared
-                    # to the inner wing cross section at the first time step.
+                    # to its inner wing cross section at the first time step.
                     base_wing_section_sweep = (
                         np.arctan(
                             (base_z_le - inner_z_les[0]) / (base_y_le - inner_y_les[0])
@@ -732,7 +735,7 @@ def generate_wings(self, num_steps=10, delta_time=0.1):
                     )
 
                     # Find the base heave angle of this wing cross section compared
-                    # to the inner wing cross section at the first time step.
+                    # to its inner wing cross section at the first time step.
                     base_wing_section_heave = (
                         np.arctan(
                             (base_x_le - inner_x_les[0]) / (base_y_le - inner_y_les[0])