The Rotation Library

The rotation library provides comprehensive tools for 3D rotations in continuum mechanics. The Rotation class inherits from scipy.spatial.transform.Rotation, so all scipy features (batch operations, mean(), Slerp, RotationSpline, etc.) are available directly, while simcoon adds continuum-mechanics methods for Voigt-notation tensors.

Architecture

The rotation system is organized in three layers:

C++ core (simcoon::Rotation in include/simcoon/Simulation/Maths/rotation.hpp):

The full rotation class with quaternion internals, all factory methods, conversions, algebra (composition, inverse, SLERP), and the Voigt mechanics operations. This is what C++ code uses directly.

pybind11 bridge (_CppRotation):

A minimal binding that exposes only what the Python layer needs from C++: a single factory method (from_quat) and the nine mechanics operations (apply_stress, apply_stiffness, etc.). Everything else — creation, conversion, composition, interpolation — is handled by scipy on the Python side.

Python class (simcoon.Rotation):

Inherits from scipy.spatial.transform.Rotation. Users get the full scipy API for free (factory methods, conversions, batch operations, mean(), Slerp, RotationSpline). The mechanics methods delegate to C++ by converting the quaternion on each call:

r.apply_stiffness(L)
    |
    v
Python:  self._to_cpp()  -->  _CppRotation.from_quat(self.as_quat())
    |                              ^ scipy provides the quaternion
    v
C++:     simcoon::Rotation::from_quat(q)  -->  .apply_stiffness(L)
    |
    v
Python:  returns 6x6 numpy array

The quaternion transfer (4 doubles) is negligible. This design avoids duplicating scipy’s well-tested rotation algebra while keeping the performance-critical Voigt operations in C++.

The Rotation Class

The Rotation class extends scipy.spatial.transform.Rotation with methods for rotating stress, strain, stiffness, and compliance tensors. It uses unit quaternions internally for numerical stability and efficient composition.

Creating Rotations

All scipy factory methods are available and return a simcoon.Rotation instance:

import simcoon as smc
import numpy as np

# Identity rotation (no rotation)
r = smc.Rotation.identity()

# From Euler angles (scipy convention: uppercase = intrinsic)
r = smc.Rotation.from_euler('ZXZ', [psi, theta, phi])

# From Euler angles in degrees, extrinsic
r = smc.Rotation.from_euler('xyz', [roll, pitch, yaw], degrees=True)

# From rotation matrix
R = np.array([[1, 0, 0], [0, 0, -1], [0, 1, 0]])  # 90° around x
r = smc.Rotation.from_matrix(R)

# From quaternion [qx, qy, qz, qw] (scalar-last)
q = np.array([0, 0, np.sin(np.pi/4), np.cos(np.pi/4)])  # 90° around z
r = smc.Rotation.from_quat(q)

# From rotation vector (axis × angle)
rotvec = np.array([0, 0, np.pi/2])  # 90° around z
r = smc.Rotation.from_rotvec(rotvec)

# From axis and angle (simcoon-specific)
r = smc.Rotation.from_axis_angle(np.pi/2, 3)  # 90° around z (axis 3)

# Random rotation (uniform distribution)
r = smc.Rotation.random()

Converting Between Representations

r = smc.Rotation.from_euler('ZXZ', [0.5, 0.3, 0.7])

# To rotation matrix
R = r.as_matrix()  # 3×3 numpy array

# To quaternion
q = r.as_quat()  # [qx, qy, qz, qw]

# To Euler angles (scipy convention)
angles = r.as_euler('ZXZ')  # [psi, theta, phi]
angles_deg = r.as_euler('ZXZ', degrees=True)

# To rotation vector
rotvec = r.as_rotvec()
rotvec_deg = r.as_rotvec(degrees=True)

# To Voigt rotation matrices (simcoon-specific)
QS = r.as_voigt_stress_rotation()  # 6×6 for stress
QE = r.as_voigt_strain_rotation()  # 6×6 for strain

Applying Rotations

To 3D vectors:

r = smc.Rotation.from_axis_angle(np.pi/2, 3)

v = np.array([1.0, 0.0, 0.0])
v_rot = r.apply(v)  # [0, 1, 0]

# Inverse rotation
v_back = r.apply(v_rot, inverse=True)  # [1, 0, 0]

To 3×3 tensors:

# Rotate a tensor: R · T · R^T
T = np.eye(3)
T_rot = r.apply_tensor(T)

To Voigt notation tensors:

# Stress vector [σ11, σ22, σ33, σ12, σ13, σ23]
sigma = np.array([100.0, 50.0, 25.0, 10.0, 5.0, 2.0])
sigma_rot = r.apply_stress(sigma)

# Strain vector [ε11, ε22, ε33, 2ε12, 2ε13, 2ε23]
epsilon = np.array([0.01, -0.005, -0.005, 0.002, 0.001, 0.0])
epsilon_rot = r.apply_strain(epsilon)

# Stiffness matrix (6×6)
L = smc.L_iso([210e9, 0.3], "Enu")
L_rot = r.apply_stiffness(L)

# Compliance matrix (6×6)
M = smc.M_iso([210e9, 0.3], "Enu")
M_rot = r.apply_compliance(M)

# Strain concentration tensor (6×6): A' = QE * A * QS^T
A_global = r.apply_strain_concentration(A_local)

# Stress concentration tensor (6×6): B' = QS * B * QE^T
B_global = r.apply_stress_concentration(B_local)

Composing Rotations

r1 = smc.Rotation.from_axis_angle(np.pi/4, 3)  # 45° around z
r2 = smc.Rotation.from_axis_angle(np.pi/6, 1)  # 30° around x

# Compose: apply r1 first, then r2
r_combined = r2 * r1

# This is equivalent to:
v = np.array([1.0, 0.0, 0.0])
v1 = r2.apply(r1.apply(v))
v2 = r_combined.apply(v)
# v1 == v2

# Inverse
r_inv = r1.inv()
r_identity = r1 * r_inv  # Should be identity

Interpolating Rotations

SLERP (Spherical Linear Interpolation) is available via scipy’s Slerp class:

from scipy.spatial.transform import Slerp

r_start = smc.Rotation.identity()
r_end = smc.Rotation.from_euler('ZXZ', [np.pi/2, np.pi/4, 0])

key_rots = smc.Rotation.concatenate([r_start, r_end])
slerp = Slerp([0, 1], key_rots)

# Interpolate along path
for t in np.linspace(0, 1, 11):
    r_t = slerp(t)
    print(f"t={t:.1f}: angle={np.degrees(r_t.magnitude()):.1f}°")

# Or use the convenience method
r_mid = r_start.slerp_to(r_end, 0.5)

Utility Methods

r = smc.Rotation.from_euler('ZXZ', [0.5, 0.3, 0.7])

# Get rotation magnitude (angle)
angle = r.magnitude()  # radians
angle_deg = np.degrees(r.magnitude())

# Check if identity
is_id = r.is_identity()

# Check equality (accounts for quaternion sign ambiguity)
r2 = smc.Rotation.from_euler('ZXZ', [0.5, 0.3, 0.7])
are_equal = r.equals(r2, tol=1e-12)

Scipy Integration

Since simcoon.Rotation inherits from scipy.spatial.transform.Rotation, all scipy features work directly:

from scipy.spatial.transform import Rotation as R

# simcoon.Rotation IS a scipy Rotation
r = smc.Rotation.from_axis_angle(np.pi/4, 3)
isinstance(r, R)  # True

# All scipy methods work directly
r.as_mrp()           # Modified Rodrigues parameters
r.as_davenport(...)  # Davenport angles

# Batch operations
rots = smc.Rotation.random(100)
rots.mean()

# RotationSpline
from scipy.spatial.transform import RotationSpline
times = [0, 1, 2, 3]
key_rots = smc.Rotation.random(4)
spline = RotationSpline(times, key_rots)
r_interp = spline(1.5)

Passing plain scipy Rotations to simcoon:

If you already have a scipy.spatial.transform.Rotation (e.g. from another library), you can upgrade it to a simcoon.Rotation with from_scipy():

from scipy.spatial.transform import Rotation as R

# Created somewhere else, e.g. by a third-party library
scipy_rot = R.from_euler('z', 45, degrees=True)

# Upgrade to simcoon.Rotation to unlock mechanics methods
r = smc.Rotation.from_scipy(scipy_rot)
L_rot = r.apply_stiffness(L)   # now available

Active vs Passive Rotations

The active parameter on Voigt methods (apply_stress, apply_strain, apply_stiffness, apply_compliance, as_voigt_stress_rotation, as_voigt_strain_rotation) controls the rotation convention:

• active=True (default): Alibi rotation — rotates the physical object while the coordinate system stays fixed.

• active=False: Alias rotation — rotates the coordinate system while the object stays fixed. This is equivalent to the inverse active rotation.

r = smc.Rotation.from_axis_angle(np.pi/4, 3)

# Active: rotate the stress tensor itself
sigma_active = r.apply_stress(sigma, active=True)

# Passive: express the same stress in a rotated coordinate system
sigma_passive = r.apply_stress(sigma, active=False)

Input Validation

All mechanics methods that accept NumPy arrays validate the input shape and raise ValueError with a descriptive message if the dimensions are wrong:

r = smc.Rotation.from_axis_angle(np.pi/4, 3)

r.apply_stress(np.array([1.0, 2.0, 3.0]))
# ValueError: sigma must have 6 elements, got 3

r.apply_stiffness(np.eye(3))
# ValueError: L must have shape (6, 6), got (3, 3)

Expected input shapes:

• apply: 1D array, 3 elements (or Nx3 array for batch)

• apply_tensor: 2D array, shape (3, 3)

• apply_stress, apply_strain: 1D array, 6 elements

• apply_stiffness, apply_compliance: 2D array, shape (6, 6)

Note

Quaternions q and -q represent the same rotation. The equals() method accounts for this antipodal equivalence.

Euler Angle Conventions

Simcoon uses scipy’s Euler angle conventions:

• Uppercase letters ('ZXZ', 'XYZ', etc.) → intrinsic rotations (rotating frame)

• Lowercase letters ('zxz', 'xyz', etc.) → extrinsic rotations (fixed frame)

Proper Euler angles (axis sequence where first = last):

• ZXZ, ZYZ, XYX, XZX, YXY, YZY

Tait-Bryan angles (all axes different):

• XYZ (roll-pitch-yaw), XZY, YXZ, YZX, ZXY, ZYX

# Intrinsic ZXZ
r1 = smc.Rotation.from_euler('ZXZ', [psi, theta, phi])

# Extrinsic xyz (same physical rotation as intrinsic ZYX)
r2 = smc.Rotation.from_euler('xyz', [roll, pitch, yaw])

Practical Examples

Example 1: Rotating Orthotropic Material Properties

import simcoon as smc
import numpy as np

# Define orthotropic material in local coordinates
# [E1, E2, E3, nu12, nu13, nu23, G12, G13, G23]
props = [150e9, 10e9, 10e9, 0.3, 0.3, 0.45, 5e9, 5e9, 3.5e9]
L_local = smc.L_ortho(props, "EnuG")

# Fiber orientation: 45° rotation around z-axis
r = smc.Rotation.from_axis_angle(np.pi/4, 3)

# Rotate stiffness to global coordinates
L_global = r.apply_stiffness(L_local)

print("Local stiffness L11:", L_local[0, 0] / 1e9, "GPa")
print("Global stiffness L11:", L_global[0, 0] / 1e9, "GPa")

Example 2: Stress Transformation in a Rotated Element

import simcoon as smc
import numpy as np

# Stress in global coordinates (Voigt notation)
sigma_global = np.array([100.0, 50.0, 0.0, 25.0, 0.0, 0.0])  # MPa

# Element orientation (Euler angles)
psi, theta, phi = np.radians([30, 45, 60])

# Create rotation (scipy convention: uppercase = intrinsic)
r = smc.Rotation.from_euler('ZXZ', [psi, theta, phi])

# Transform stress to local element coordinates
sigma_local = r.inv().apply_stress(sigma_global)

print("Global stress:", sigma_global)
print("Local stress:", sigma_local)

Example 3: Averaging Orientations

import simcoon as smc
import numpy as np

# Generate random orientations
rotations = smc.Rotation.random(10)

# Use scipy's built-in mean
r_mean = rotations.mean()

print("Mean rotation angle:", np.degrees(r_mean.magnitude()), "degrees")

See Also

• The Constitutive Library - Constitutive matrices (L_iso, L_ortho, etc.)

• Transfer helpers (Voigt <-> matrix) - Voigt notation conversions

• The Kinematics Library - Finite strain kinematics

• The Rotation Module - Full C++ API reference