Modes of free vibration in a fixed-fixed beam I am trying to find the mode shapes of vibration on a fixed-fixed beam.
this can be characterised as a PDE:
$$ EI \frac{\partial^4 v(x,t)}{\partial x^4} + \rho A \frac{\partial^2 v(x,t)}{\partial t^2} = 0 $$
with the boundary conditions
$$v(0,t) = 0  \hspace{2mm} \hspace{2mm} v_x(0,t) = 0 $$
$$v(l,t) = 0 \hspace{2mm}  \hspace{2mm} v_{x}(l,t) = 0 $$
We know that $E$, $I$, $\rho$ and $A$ are constant along the length
I have used the separation of variables method and then substituted in my boundary conditions and reached two different solutions:
$$\text{Re}\left[e^{ItB^2} (\cosh[Bx] - \cos[Bx]) + \frac{\sinh[Bl] + \sin[Bl]}{\cosh[Bl] -\cos[Bl]} (\sin[Bx] - \sinh[Bx])\right]$$
Plot of above showing fixed at l
$$\text{Re}\left[e^{ItB^2} \left((\sinh[Bx] - \sin[Bx]) - \frac{\sinh[Bl] - \sin[Bl]}{\cosh[Bl] 
+ \cos[Bl]} (\cosh[Bx] - \cos[Bx])\right) \right]$$
Plot of above showing no slope at l
where $B$ is a constant that combines the natural frequency and the other constants
Each of these solutions satisfies a different three of the boundary conditions but I cant work out how to satisfy them all at the same time
 A: The procedure to extract the modes of a linear PDE that is solvable by separation of variables is standard.
First separate the variables setting $\omega^4=\frac{\rho A}{EI}$ and set $u=X(x)T(t)$, where $X,T$ can be shown to satisfy the equations:
$$X''''-(\lambda\omega)^4X=0\\T''+\lambda^4T=0$$
where $\lambda$ is an, as of now, undetermined real constant. We choose to parametrize the separated solutions this way, since we would like to obtain a wave-like solution to the equations. Solving these equations we find that:
$$X(x)=A\cos(\lambda\omega x)+B\sin(\lambda\omega x)+C\cosh(\lambda\omega x)+D\sinh(\lambda\omega x)\\T(t)=A'\cos(\lambda^2t)+B'\sin(\lambda^2t)$$
However $X(x)$ must satisfy some constraints:
$$X(0)=X(l)=X'(0)=X'(l)=0$$
which must be simultaneously satisfied for them to be a solution to the fixed fixed beam problem (note the close analogy to the modes of oscillation of a string fixed on both ends). From this we obstain a linear system of equations for the parameters $A,B,C,D$:
$$M\begin{pmatrix}A\\B\\C\\D\end{pmatrix}=\begin{pmatrix}1&0&1&0\\\cos(\lambda\omega l)&\sin(\lambda\omega l)&\cosh(\lambda\omega l)&\sinh(\lambda\omega l)\\0&1&0&1\\-\sin(\lambda\omega l)&\cos(\lambda\omega l)&\sinh(\lambda\omega l)&\cosh(\lambda\omega l)\end{pmatrix}\begin{pmatrix}A\\B\\C\\D\end{pmatrix}=0$$
From basic linear algebra, we know that for this system to have a solution a necessary and sufficient condition is:
$$\det M=\cos (\lambda\omega l)\cosh(\lambda\omega l)+\sin (\lambda\omega l)\sinh (\lambda\omega l)-1=0$$ 
and therefore the eigenvalues can be written $\lambda_n=\frac{x_n}{\omega l}$, where $x_n$ represents the solutions to the equation
$$\cos x\cosh x+\sin x \sinh x=1$$
There doesn't seem to be a analytical solution for this equation, but it is easy to show that there is an infinite number of solutions to it as follows. Set $f(x)=\cos(x)\cosh(x)+\sin(x)\sinh(x)-1$. Since this function is even, we will restrict our attention to positive eigenvalues. Then $f'(x)=2\cos x\sinh x$. It is easy to see that this posesses minima at $2n\pi-\pi/2,n\geq1, n\in \mathbb{N} $ and maxima at $2n\pi+\pi/2, n\geq0$. Since 
$$f(2n\pi-\pi/2)=-\sinh(2n\pi-\pi/2)-1<0~~,~~f(2n\pi+\pi/2)=\sinh(2n\pi+\pi/2)-1>0, f(0)=0$$ 
we find that it possesses exactly one root in each one of those intervals. Specifically, we see that:
$$x_1\in(\pi/2,3\pi/2), x_{n+1}\in(2n\pi-\pi/2, 2n\pi+\pi/2), n>0$$
A: If you use the separation of variables, such that:
$$v(x,t)=V(x)\exp(i\omega t),$$
it is only necessary to write the space equation
$$\frac{d^4}{dx^4}V(x)+k^4V(x)=0$$
where $k=\left(\omega^2 \frac{\rho A}{EI}\right)^{1/4}$.
the solution can be written as
$$V(x)=A\cos(kx)+B\sin(kx)+C\cosh(kx)+D\sinh(kx)$$
and applying your boundary conditions:
$$V(0)=0, V'(0)=0, V(l)=0, V'(l)=0$$
The four equations can be written in matrix format:
$$\begin{pmatrix}1 & 0 & 1 & 0\\
\cos{\left( kl\right) } & \sin{\left( kl\right) } & \cosh{\left( kl\right) } & \operatorname{sinh}\left( kl\right) \\
0 & 1 & 0 & 1\\
-\sin{\left( kl\right) } & \cos{\left( kl\right) } & \operatorname{sinh}\left( kl\right)  & \cosh{\left( kl\right) }\end{pmatrix}\begin{pmatrix} A\\B\\C\\D\end{pmatrix}=0$$
For nontrival solution, the determinant of the matrix in this system of equations reduces to
$$\cos(kl)\cosh(kl)-1=0$$
which needs to be solved numerically. For numerical computation of the solution is better to rewrite this equation as
$$\cos(kl)=\frac{1}{\cosh(kl)}$$
so you can also plot a graph of the functions on the both sides of equation and check where they are equal. The first root, for instance, occurs when $kl\approx 4.73004$. It is possible to check that for large $kl$, the roots approximate to $\cos(kl)=0$.
Now, to calculate the mode shapes, it is necessary to use three of the four equations and write the mode shape as a function of one of the constants ($A,B,C$ or $D$).
If we choose to write the constants in terms of $A$, it is found that
$$B=-\frac{\cosh{\left( kl\right) }-\cos{\left( L k\right) }}{\operatorname{sinh}\left( kl\right) -\sin{\left( kl\right) }}A $$
$$C=-A$$
$$D=\frac{\cosh{\left( kl\right) }-\cos{\left( kl\right) }}{\operatorname{sinh}\left(  kl\right) -\sin{\left( kl\right) }}A $$
The expression for $B,C$ and $D$ are than replaced in the solution $V(x)$
$$V(x)=\left[\cos(kx)-\cosh(kx) + \frac{\cosh(kl)-\cos(kl)}{\sinh(kl)-\sin(kl)}\left(\sinh(kx)-\sin(kx)\right)\right]A$$
since $V(x)$ is a mode shape, it can be scaled by making $A=1$, to give:
$$\Psi(x)=\cos(kx)-\cosh(kx) + \frac{\cosh(kl)-\cos(kl)}{\sinh(kl)-\sin(kl)}\left(\sinh(kx)-\sin(kx)\right)$$
which is not unique.
The first mode shape is shown:
A: Too long for a late comment.
Starting from @DinosaurEgg's answer, if we look for the positive zero's of
$$f(x)=\cos (x) \cosh (x)+\sin (x) \sinh (x)-1$$
$$f'(x)=2 \cos (x) \sinh (x)$$ we have minima and maxima for $x=(2n+1)\frac \pi 2$.
Building the $[2,2]$ Padé approximant round these points, we just need to solve a quadratic equation in $\big[x-(2n+1)\frac \pi 2\big]$ to have a rather good estimate. For the first solutions
$$\left(
\begin{array}{ccc}
n & \text{approximation} & \text{solution} \\
 0 & 2.19426 & 2.18354 \\
 1 & 5.52114 & 5.50353 \\
 2 & 8.65614 & 8.63913 \\
 3 & 11.7980 & 11.7810
\end{array}
\right)$$ and
$$x_{n+1} \sim x_n+\pi$$ One iteration of Newton method is more than sufficient.
