Rotor Unbalance & Critical Speed · Model record

Version: rotor-dynamics-v1.0. Original teaching implementation. All parameter defaults are synthetic, not measured aircraft data.

What is modelled

One disc on a massless flexible shaft with stationary supports. The two transverse directions have identical linear stiffness and viscous damping. Gravity is absorbed into a static equilibrium, and displacement is measured from that equilibrium. This is the isotropic Jeffcott reduction, not a multi-bearing engine rotor.

The supported centre position is (x,y). The centre of mass is displaced from it by the rotating residual eccentricity vector. Initial offset and correction are summed as vectors before computing the forcing. x is positive right; y is positive up. The positive spin angle is anticlockwise when viewed along the shaft from the front.

Input units are kg, kN/m, mm, degrees, rpm and seconds. The model converts to SI internally.

Governing equations

Let e be the magnitude of the residual eccentricity and θ include its initial phase. With prescribed spin speed Ω=dθ/dt,

m \ddot{x} + c \dot{x} + k x = m e (Ω^{2} \cos θ + \dot{Ω} \sin θ)

m \ddot{y} + c \dot{y} + k y = m e (Ω^{2} \sin θ - \dot{Ω} \cos θ)

The tangential terms follow by differentiating the eccentric mass-centre position twice. They are retained during acceleration. They vanish at constant speed.

ω_{n} = \sqrt{\frac{k}{m}} c = 2 ζ m ω_{n}

At constant speed the synchronous orbit radius and phase lag are

\begin{matrix} \frac{X}{e} = \frac{r^{2}}{\sqrt{{(1 - r^{2})}^{2} + {(2 ζ r)}^{2}}} \\ φ = atan2 (2 ζ r, 1 - r^{2}) \\ r = \frac{Ω}{ω_{n}} \end{matrix}

The marker at r=1 is the undamped natural frequency (often called the critical-speed reference in this simple model). For a non-zero unbalance and 0<ζ<1/√2, the maximum of the unbalance-response amplitude occurs at r=1/√(1−2ζ²), slightly above that marker. A transient peak can occur elsewhere. At exactly zero unbalance, the observable phase of the zero response is undefined; the plotted phase is still the mathematical transfer-function phase.

Prescribed drive and time integration

Starting from rest, u=t/t_r and Ω=Ω_f(3u²−2u³) for 0≤t≤t_r. After t_r the speed is fixed at Ω_f. The analytic angle during the ramp is Ω_f t_r(u³−u⁴/2); the analytic acceleration is Ω_f(6u−6u²)/t_r. Both speed and acceleration are continuous at the ramp boundaries.

Classical RK4 integrates x, xdot, y, ydot, cumulative forcing work and cumulative damping dissipation. The maximum step is the smaller of 0.001 s and 1/140 of the smaller of the natural/spin periods. Numerical refinement is checked independently. A saved history has at least 32 samples per shorter natural/spin period for the orbit view. The time plot uses extrema-preserving display reduction. The full integration, not the display reduction, determines the peak.

The transient is not a succession of settled states. Dotted curves on the time plot show the constant-speed radius at the instantaneous Ω only as a comparison. They are not fed into the integrator.

The drive speed is prescribed. Its torque, power supply, motor/controller behaviour and any feedback from lateral vibration to spin speed are outside this model.

Lateral energy account

E = \frac{1}{2} m ({\dot{x}}^{2} + {\dot{y}}^{2}) + \frac{1}{2} k (x^{2} + y^{2})

E (t) - E (0) = \int_{0}^{t} (F_{x} \dot{x} + F_{y} \dot{y}) d t - \int_{0}^{t} c ({\dot{x}}^{2} + {\dot{y}}^{2}) d t

This is an energy account for the two lateral coordinates. It is not the total rotor spin energy or the electrical energy consumed by its drive. The displayed residual is normalised by the largest stored-energy change, forcing-work change or dissipation change during that run. Supplied initial cumulative work and dissipation counters are subtracted before forming these changes. It is checked at saved samples. For the zero solution its normalised value is defined as zero. The peak-displacement search includes the initial state and records its actual initial spin speed.

Visual contract

Orbit: the calculated shaft-centre x,y positions, on a labelled micrometre scale shared by A and B.
Shaft: a schematic sinusoidal outline through fixed end supports and the calculated centre y displacement. It is not a beam finite-element mode shape. Deflection is explicitly magnified relative to a 1 m reference span used only for drawing; the selected stiffness is independent of this illustrative span. The graphic does not predict stresses.
Disc orientation: sampled physical phase. A rapidly rotating graphic can alias; RPM is the speed indication. No judgement of speed should be made from apparent animation direction.
Balance diagram: equivalent eccentricity vectors. The correction endpoint has matching numerical sliders for keyboard access. No real-engine mass placement, tooling, procedure or acceptance threshold is prescribed.
A/B comparison: A retains its complete parameter set; B is live. Each transient has its own ramp duration. Both are integrated to the same end time, long enough for the slower selected run. Their common time axis is absolute time from the start of each run, not a normalised cycle.

CSV export

The Download A/B CSV button reuses the current calculated histories and analytical sweeps. It does not run another integration. Export is disabled while a changed setup is awaiting calculation, so its settings describe the exported runs rather than stale results.

One rectangular UTF-8 CSV contains a record_type column: metadata, setting, summary, steady_final, steady_sweep and transient. The metadata records the model version, export timestamp, comparison availability, interpolation method and limitations. Setting rows retain all eight inputs separately for A and B, with their input units; summary rows include residual eccentricity, natural frequency, damping, numerical step, peak and energy-residual evidence. No saved A is represented explicitly, rather than copied from B.

The steady rows retain the calculated amplitude, phase lag and force against spin ratio and RPM. The integrated histories share one uniform absolute-time grid from zero to the common end time, with at least as many time points as the denser stored history. Stored displacement, velocity, radius, energy, work, dissipation and residual are linearly interpolated; radius and energy are interpolated as stored scalar quantities, not recomputed from interpolated coordinates. The exported forcing components, speed, acceleration and shaft angle are evaluated analytically at each shared time. Column names carry SI units, while setting-row units preserve the displayed input units. The file includes the start and final samples; it does not claim to export every RK4 step. Interpolation can smooth extrema, so the separately exported peak summaries retain the peak found during full integration.

The exported lateral-energy account is not total spin energy or electrical consumption. Phase at zero unbalance is the mathematical transfer-function phase, not an observable motion. The shaft-angle column is the mechanical rotation angle; the residual eccentricity vector has its separate setting-derived initial orientation.

Numerical checks

Run node test-model.js in this directory. It creates verification.json. The tests include zero residual unbalance, resonance amplitude/phase, a high-speed limit, a known synchronous trajectory, finite-difference verification of variable-speed forcing, RK4 refinement, lateral energy closure, crossing event accuracy, independent sweep histories, linear scaling after correction, and invalid input rejection.

These checks establish consistency of the implementation with its equations over selected cases. No experimental validation, rotor certification or aircraft safety assessment is claimed.

Limitations

No anisotropic or non-linear bearings, gyroscopic frequency splitting, internal damping instability, oil whirl, rub, contact, real clearances, blade flexibility, thermal growth, multi-disc modes, crack growth or fatigue. The 2 mm excursion flag is an illustrative scale check, not a real mechanical limit; the model does not turn a large excursion into a diagnosis of damage.

The allowed damping ratio is positive (0.015–0.3), so the interface never presents the undamped resonant singularity as a finite physical prediction. Changing mass or stiffness changes the natural frequency; it does not automatically retune the selected final RPM.