How can I determine the damping coefficient that minimizes the amplitude of vibrations? This is an extension of Coupled ODEs that model a quad rotor \begin{align} m_{eq}\ddot{y}_1 &= k_{eq}\Bigl[y_1 + \frac{m_{eq}}{m_b}\sum_jy_j\Bigr] + c\Bigl[\dot{y}_1 + \frac{m_{eq}}{m_b}\sum_j\dot{y}_j\Bigr] - m_{eq}gy_1 + F_1(t)\\ m_{eq}\ddot{y}_2 &= k_{eq}\Bigl[y_2 + \frac{m_{eq}}{m_b}\sum_jy_j\Bigr] + c\Bigl[\dot{y}_2 + \frac{m_{eq}}{m_b}\sum_j\dot{y}_j\Bigr] - m_{eq}gy_2 + F_2(t)\\ m_{eq}\ddot{y}_3 &= k_{eq}\Bigl[y_3 + \frac{m_{eq}}{m_b}\sum_jy_j\Bigr] + c\Bigl[\dot{y}_3 + \frac{m_{eq}}{m_b}\sum_j\dot{y}_j\Bigr] - m_{eq}gy_3 + F_3(t)\\ m_{eq}\ddot{y}_4 &= k_{eq}\Bigl[y_4 + \frac{m_{eq}}{m_b}\sum_jy_j\Bigr] + c\Bigl[\dot{y}_4 + \frac{m_{eq}}{m_b}\sum_j\dot{y}_j\Bigr] - m_{eq}gy_4 + F_4(t) \end{align} where $F_i(t) = F_i\cos(\omega t + \phi)$
I tried converting to the state space equation but this doesn't make the system any easier to deal with:
By letting $q_1 = y_1$, $q_2 = \dot{y}_1$, $q_3 = y_2$, $q_4 = \dot{y}_2$, $q_5 = y_3$, $q_6 = \dot{y}_3$, $q_7 = y_4$, and $q_8 = \dot{y}_4$, we obtain the following state equations. \begin{align} \dot{\mathbf{q}} &= \begin{bmatrix} 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\ a & b & C & d & C & d & C & d\\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\ C & d & a & b & C & d & C & d\\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\ C & d & C & d & a & b & C & d\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1\\ C & d & C & d & C & d & a & b \end{bmatrix}\mathbf{q} + \frac{1}{m_{eq}} \begin{bmatrix} F_1(t)\\ F_2(t)\\ F_3(t)\\ F_4(t) \end{bmatrix}^{\intercal} \begin{bmatrix} 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \end{bmatrix}\\ \mathbf{y} &= \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \end{bmatrix}\mathbf{q} \end{align} where $a = \frac{k_{eq}}{m_{eq}} + \frac{k_{eq}}{m_b} - g$, $b = \frac{c}{m_{eq}} + \frac{c}{m_b}$, $C = \frac{k_{eq}}{m_b}$, and $d = \frac{c}{m_b}$.
If I take $\mathcal{L}^{-1}\{(s\mathbb{I} - \mathbf{A})^{-1}\}$ with Mathematica, the solution is ridiculously long.
I have also tried setting up $[m]\ddot{y} + [c]\dot{y} + [k]y = [F]$ but I am not sure how to solve this equation either.