Stability Tutorial

This tutorial covers ship stability calculations, including GZ curves (righting lever curves) and KN curves (cross curves of stability).

Introduction

Ship stability analysis determines how a vessel responds to heeling forces and whether it will return to an upright position.

Key concepts:

  • GZ (Righting Lever): The horizontal distance between G and B when heeled

  • GZ Curve: Plot of GZ vs heel angle - shows stability characteristics

  • Range of Stability: Maximum heel angle with positive righting moment

  • KN Curves: Cross curves assuming KG = 0, used for multiple loading conditions

The Righting Lever (GZ)

When a ship heels, the center of buoyancy (B) moves. The righting lever (GZ) is:

\[\begin{split}GZ = KN - KG \\cdot \\sin(\\phi)\end{split}\]

Where:

  • KN = righting lever assuming KG = 0

  • KG = height of center of gravity above keel

  • φ = heel angle

Setting Up

First, load your hull and create the stability calculator:

from pynavaltoolbox import Hull
from pynavaltoolbox.stability import StabilityCalculator

hull = Hull("path/to/hull.stl")
stability = StabilityCalculator(hull, water_density=1025.0)

Calculating a GZ Curve

The GZ curve is calculated for a specific displacement and center of gravity:

# Define loading condition
displacement = 5000000.0  # kg (5000 tonnes)
cog = (70.0, 0.0, 10.0)  # LCG, TCG, VCG in meters

# Calculate GZ curve from 0° to 90° in 5° increments
gz_curve = stability.calculate_gz_curve(
    displacement_mass=displacement,
    cog=cog,
    heels=list(range(0, 91, 5))
)

# Print results
print("Heel(°)\tGZ(m)\tDraft(m)\tTrim(°)")
print("-" * 50)
for point in gz_curve.points:
    print(f"{point.heel}\t{point.gz:.3f}\t{point.draft:.2f}\t\t{point.trim:.2f}")

Analyzing the GZ Curve

The StabilityCurve provides methods for analysis:

# Key stability characteristics
max_gz = gz_curve.get_max_gz()
angle_max_gz = gz_curve.get_angle_of_max_gz()
range_stability = gz_curve.get_range_of_stability()

print(f"Maximum GZ: {max_gz:.3f} m")
print(f"Angle of Maximum GZ: {angle_max_gz:.1f}°")
print(f"Range of Stability: {range_stability:.1f}°")

# Access raw data for plotting
heels = gz_curve.heels  # numpy array
gz_values = gz_curve.gz_values  # numpy array

Plotting the GZ Curve

Create a publication-quality GZ curve plot:

import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10, 6))

ax.plot(gz_curve.heels, gz_curve.gz_values, 'b-', linewidth=2)
ax.axhline(y=0, color='k', linestyle='-', linewidth=0.5)

# Mark key points
ax.axvline(x=angle_max_gz, color='r', linestyle='--', alpha=0.5)
ax.scatter([angle_max_gz], [max_gz], color='r', s=100, zorder=5)
ax.annotate(f'Max GZ = {max_gz:.3f} m at {angle_max_gz:.0f}°',
            xy=(angle_max_gz, max_gz),
            xytext=(angle_max_gz + 10, max_gz),
            fontsize=10)

ax.set_xlabel('Heel Angle (degrees)', fontsize=12)
ax.set_ylabel('GZ (m)', fontsize=12)
ax.set_title('Righting Lever Curve', fontsize=14)
ax.grid(True, alpha=0.3)
ax.set_xlim(0, 90)

plt.tight_layout()
plt.savefig('gz_curve.png', dpi=300)
plt.show()

Calculating KN Curves

KN curves (cross curves of stability) are calculated assuming KG = 0:

# Calculate KN curve at specific displacement
kn_curve = stability.calculate_kn_curve(
    displacement_mass=5000000.0,
    heels=[0, 5, 10, 15, 20, 30, 45, 60, 75, 90]
)

print("Heel(°)\tKN(m)")
print("-" * 20)
for point in kn_curve.points:
    print(f"{point.heel}\t{point.kn:.3f}")

Using KN Curves for Different KG Values

KN curves allow quick GZ calculation for different loading conditions:

import numpy as np

# Different KG values to evaluate
kg_values = [8.0, 10.0, 12.0]

print("Heel\t", end="")
for kg in kg_values:
    print(f"GZ(KG={kg})\t", end="")
print()

for point in kn_curve.points:
    print(f"{point.heel}°\t", end="")
    for kg in kg_values:
        gz = point.kn - kg * np.sin(np.radians(point.heel))
        print(f"{gz:.3f}\t\t", end="")
    print()

Effect of VCG on Stability

This example shows how moving the center of gravity affects stability:

import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10, 6))

vcg_values = [8.0, 10.0, 12.0, 14.0]
colors = ['green', 'blue', 'orange', 'red']

for vcg, color in zip(vcg_values, colors):
    curve = stability.calculate_gz_curve(
        displacement_mass=5000000.0,
        cog=(70.0, 0.0, vcg),
        heels=list(range(0, 91, 5))
    )
    ax.plot(curve.heels, curve.gz_values, color=color,
            linewidth=2, label=f'VCG = {vcg} m')

ax.axhline(y=0, color='k', linestyle='-', linewidth=0.5)
ax.set_xlabel('Heel Angle (degrees)', fontsize=12)
ax.set_ylabel('GZ (m)', fontsize=12)
ax.set_title('Effect of VCG on Stability', fontsize=14)
ax.legend()
ax.grid(True, alpha=0.3)
ax.set_xlim(0, 90)

plt.tight_layout()
plt.savefig('vcg_comparison.png', dpi=300)
plt.show()

Stability Criteria

Common stability criteria (IMO A.749) to check:

  1. Area under GZ curve from 0° to 30°: ≥ 0.055 m·rad

  2. Area under GZ curve from 0° to 40°: ≥ 0.090 m·rad

  3. Area under GZ curve between 30° and 40°: ≥ 0.030 m·rad

  4. GZ at 30°: ≥ 0.20 m

  5. Angle of maximum GZ: ≥ 25°

  6. Initial GM: ≥ 0.15 m

import numpy as np
from scipy.integrate import trapezoid

# Calculate area under curve (convert to radians for integration)
heels_rad = np.radians(gz_curve.heels)
gz = gz_curve.gz_values

# Find indices for key angles
idx_30 = np.argmin(np.abs(gz_curve.heels - 30))
idx_40 = np.argmin(np.abs(gz_curve.heels - 40))

area_0_30 = trapezoid(gz[:idx_30+1], heels_rad[:idx_30+1])
area_0_40 = trapezoid(gz[:idx_40+1], heels_rad[:idx_40+1])
area_30_40 = trapezoid(gz[idx_30:idx_40+1], heels_rad[idx_30:idx_40+1])

print("=== Stability Criteria Check ===")
print(f"Area 0-30°: {area_0_30:.4f} m·rad (min: 0.055)")
print(f"Area 0-40°: {area_0_40:.4f} m·rad (min: 0.090)")
print(f"Area 30-40°: {area_30_40:.4f} m·rad (min: 0.030)")
print(f"GZ at 30°: {gz[idx_30]:.3f} m (min: 0.20)")
print(f"Angle of max GZ: {gz_curve.get_angle_of_max_gz():.1f}° (min: 25°)")

Next Steps