Vibration Transmissibility Optimization

This was the term project of Mechanical Vibrations at Boğaziçi University, taught by Çetin Yılmaz and completed with Cem Geçgel and Mustafa Çağatay Sipahioğlu. We were given a chain of $n$ disks connected by torsional springs, driven from one end by a rotating base, and asked two questions. First, how much of the base’s oscillation reaches the last disk at each excitation frequency? Second, how should two viscous dampers be placed and tuned so that the worst transmission over the whole frequency range is as small as possible?

The system looks abstract, but that is how mathematical modeling works. Materials in real life vibrate, and most of them can be approximated by the same three primitives: inertias, springs, and dampers, in translational or rotational form. For example, the road shakes a car through its wheels, at frequencies that depend on the car’s speed and the road’s geometry, and the suspension exists so that the passengers barely feel that vibration. That’s also a vibration transmission optimization problem. Reduce the car to inertias, springs, and dampers, and it becomes a linear algebra problem. The optimized parameters then transfer back to the real car as material and geometry choices that match the numerical optimization.

This project practices exactly that on a small, perfectly clean mathematical model of such a system.

Transmissibility

Each disk has inertia $I=100/n$ and each spring has stiffness $k=25n$, with $n$ between 2 and 5. The base oscillates with rotation amplitude $\Phi$ at angular frequency $\omega$, written in complex form as $\varphi(t)=\Phi e^{i\omega t}$, whose real part is the physical motion. Complex exponentials are the natural language of vibration analysis: differentiating one multiplies it by $i\omega$ instead of reshuffling sines into cosines, and a single complex number carries both an amplitude and a phase.

Writing the torque balance for every inertia and collecting coefficients gives the equation of motion in matrix form:

Here $\boldsymbol{\theta}(t)$ is an $n$-dimensional vector that collects the rotation angles of the disks, and $M$, $C$, and $K$ are the $n\times n$ inertia, damping, and stiffness matrices that the torque balance assembles. The right-hand side is the excitation, an $n$-dimensional vector whose components are the external torques on the corresponding disks: the base spring applies the torque $k(\varphi-\theta_1)$ to the first disk and to nothing else, so the known part $k\varphi(t)$ multiplies the first basis vector $\boldsymbol{e}_1=(1,0,\dots,0)^{\mathsf{T}}$ while the $k\theta_1$ part stays inside $K$.

In steady state, everything oscillates at the driving frequency $\omega$, because a linear system cannot invent new frequencies: the equation of motion only adds, scales, and differentiates motions, and none of those operations changes a sinusoid’s frequency, only its amplitude and phase. So $\boldsymbol{\theta}(t)=\boldsymbol{\Theta}e^{i\omega t}$, where the entries of $\boldsymbol{\Theta}$ are complex amplitudes carrying each inertia’s response magnitude and phase. Substituting turns the differential equation into one complex linear system per frequency:

The transmissibility is the amplitude ratio between the last disk and the base, and it comes from a single entry of the inverse of the matrix above, called the receptance matrix, which maps torque amplitudes to rotation amplitudes:

The subscript picks row $n$, column 1 of the inverse: the rotation of the last disk caused by a unit torque on the first. Reading the transmissibility is simple. $T=1$ means the last disk swings exactly as far as the base, larger means amplification, smaller means isolation.

For the bare chain there is no damping, $C=0$, and the frequency response is governed by the natural frequencies, the solutions $\omega_j$ of the eigenvalue problem $K\boldsymbol{v}=\omega^2M\boldsymbol{v}$, whose eigenvectors $\boldsymbol{v}$ are the mode shapes. An $n$-degree-of-freedom system has $n$ of them, and for an end-to-end chain every mode is visible from the base to the last disk, so all $n$ show up as resonances in $T$:

At low frequency the chain simply follows the base and $T\to1$. At each natural frequency the response is unbounded, since nothing dissipates energy. Past the last resonance the chain can no longer keep up with the base and genuinely isolates. If the operating range contains resonances, the worst case over the range is what matters, and for the bare chain the worst case is infinite.

Houdaille dampers

The assigned remedy is the Houdaille damper: a free inertia $I_a$ sealed inside the host disk, coupled to it only by a viscous fluid with damping constant $c_a$. Each damper adds one degree of freedom, the free inertia’s own rotation: $I_a$ joins the diagonal of $M$, the block $c_a\bigl[\begin{smallmatrix}1&-1\\-1&1\end{smallmatrix}\bigr]$ joins $C$, coupling absorber and host through their relative velocity, and $K$ is untouched, so the full system has $n+2$ degrees of freedom. Both extremes of $c_a$ are useless. With $c_a=0$ the inertia spins freely and does nothing. With $c_a\to\infty$ it rigidly locks to its host, adding inertia that shifts the resonances slightly but dissipating nothing, so the peaks stay unbounded. Every useful setting is in between, where the inertia lags its host and fluid friction converts the relative motion into heat. The best setting is a genuine optimization problem.

Two such dampers share a fixed inertia budget, $I_{a1}+I_{a2}=\mu$ with $\mu$ between 0.1 and 0.3. The objective is minimax, over the frequency range prescribed by the assignment:

where $\omega_{\max}$ is the largest natural frequency. In the first stage the dampers sit on disks 1 and $n$ with $I_{a1}=I_{a2}=\mu/2$ and only $c_{a1},c_{a2}$ are tuned; in the final stage the inertia split and the two host disks are design variables too.

The optimization

The peaks are extremely narrow, so evaluating $T$ on a dense frequency grid at every optimizer step is both expensive and unreliable. What is needed is the exact location of every peak, and that means finding the resonance frequencies of the damped system. The classical route is closed: two dampers touching only two disks make the damping non-proportional, meaning the undamped mode shapes do not diagonalize $C$, and the modal analysis that works for the bare chain breaks down. The resonance condition itself still generalizes: peaks live where the dynamic stiffness $K-\omega^2M+i\omega C$ comes closest to singular. The systematic way to find them is the state-space transformation.

Stacking the angles and the angular velocities into one state vector turns the $n+2$ second-order equations into $2(n+2)$ first-order ones:

where $\mathbf{1}$ and $\mathbf{0}$ are the identity and zero matrices. Unlike the bare chain’s eigenvalue problem, $A$ is not symmetric, and its eigenvalues come in complex pairs $\lambda=-\sigma\pm i\omega_d$: the decay rate $\sigma$ says how fast that mode dies out, and $\omega_d$ is its damped natural frequency, where a transmissibility peak lives. This eigenvalue problem is solved for every candidate design the optimizer tries, and $T$ is evaluated at the resulting $\omega_d$ values: a handful of exact evaluations per step instead of thousands of grid points. Getting to this formulation was a small research problem of its own; I asked about it on Mathematics Stack Exchange at the time and ended up answering it myself: math.stackexchange.com/questions/4464736

The outer problem, minimizing the largest of the $n$ peaks, went to MATLAB’s minimax solver fminimax. In the free-placement stage the inertia budget entered as a linear equality constraint, all $\binom{n}{2}$ host pairs were enumerated, and the best design was kept.

The result

For the representative case, $n=5$ and $\mu=0.3$, tuning equal absorbers fixed to disks 1 and 5 caps every peak at a finite height. Releasing the placement and the split beats it: the optimizer zeroes the first absorber entirely and spends the whole budget on the last disk. The interpretation is physical: the middle disks only pass energy along, and dissipation is most effective where the motion being minimized actually occurs. Splitting the budget wastes half of it on a disk nobody asked about.

The complete derivation, all tested $n$ and $\mu$ cases, and the MATLAB listing are in the report: PDF

Source code: github.com/sinaatalay/MultipleTorsionalViscousDamperDesign