Operators as Matrices: Towards Dirac Notation

 In order to make the jump from the wave function representation to Dirac Notation more obvious and natural let us consider the simplest problem usually discussed at great length in quantum mechanics introductory text book, namely a particle in an infinite box. We shall show that we can write the position and momentum operators as matrices.
To start off, we have a particle in an box with infinite potential at some finite boundary i.e the particle can't escape or tunnel and we are assured that at the boundaries the particle's wave function becomes zero. After going through the calculation(we shall not be done here) we arrive at the fact that the wave functions for the different energy states are given by \(\psi_n=  \sqrt(\frac{2}{L}) \sin (\frac{n\pi x}{L}) \) where L is the length of the box and the energies are  \(E_n =  \frac{n^2 \pi^2 \hbar ^2}{2mL^2} \). Now the Hamiltonian is now merely the kinetic energy operator : \( \frac{p^2}{2m} \) that means that,
\begin{align}
\frac{p^2}{2m} \Psi(t) &= \sum a_n \frac{p^2}{2m} \psi_n \\
                                      &= \sum_n a_n E_n \psi_n \hspace{10mm} \text{eq.1}
\end{align}

Remember that we could represent the wave function as a column vector, the elements of which are the coefficients:

\begin{equation}
\begin{pmatrix} a_1 \\ a_2 \\ a_3 \\ \vdots \end{pmatrix}
\end{equation}

Looking at eq.1 we see that we can represent the kinetic energy operator as
\begin{equation}
\begin{pmatrix}
 E_1 &  0 &  0 & 0 & \ldots \\
 0  &  E_2 & 0 & 0 &  \ldots \\
 0  &   0  &  E_3  & 0 & \ldots \\
\vdots & \vdots & \vdots \\
\end{pmatrix}
\end{equation}

What about the position operator, \( \hat{x}\)?
\begin{align}
 \hat{x}\Psi &= \sum_n a_n \hat{x}\psi_n \\
  \hat{x} \psi_n &= \sum _m x_{nm} \psi_ m  \hspace{10mm} \text{eq.2}
\end{align}
We now use the orthogonality of the \( \psi_m \) to find the matrix elements \( x_{mn} \):
\begin{equation}
 \int \psi_m ^{*} \hat{x} \psi_ n = x_{mn} \hspace{10mm} \text{eq.3}
\end{equation}
We multipled both sides by \(\psi^* \)  in eq.2 and use the fact that  \(\int \psi_n* \psi_m = \delta_{mn} \)
Doing the integral using the eigenfunctions of the particle in  a box mentioned in the introductory second paragraph we get,  \(x_{mn} =  \frac{L}{\pi^2}\frac{4mn}{(m^2-n^2)^2}((-1)^{(m-n) }-1)\). For m=n we get the diagonal elements to be L/2. Remember that m represents the row number and n represents the column number. So when we say elements where m=n we mean the diagonal elements. We then go to our equation make m=n and then find the number to go into the diagonal spots in our matrix. If you want to find the number in the first row second column, make m=1 an n=2 in the equation. We therefore have that our matrix representing the position operator to be:
\begin{equation}
\begin{pmatrix}
\frac{L}{2} & \frac{-16L}{9 \pi^2} & 0 \ldots \\
\frac{-16L}{9 \pi^2}& \frac{L}{2} & \frac{-48L}{25 \pi^2}&\ldots \\
0 & \frac{-48L}{25 \pi^2} & \frac{L}{2}& \ldots \\
\vdots & \vdots &  \vdots  & \ldots
\end{pmatrix}
\end{equation}

Notice in eq.3 we sandwiched the operator between the eigenfunction and its complex conjugate in order to find the matrix element.
The same can be done with the momentum operator \hat{p}or represented concretely as \( \frac{-i\hbar d}{dx} \). We have :
\begin{align}
 \hat{p}\Psi &= \sum_n a_n \hat{p} \psi_n \\
   \hat{p}\psi_n &= \sum_m  p_{nm} \psi_ m \\
\end{align}
Using the orthogonality of the \( \psi_n \)  we get another integral to do namely, \(\int \psi_m ^{*} \hat{p} \psi_ n \). This time the integrals gives us \( p_{mn} = \frac{\hbar}{iL} \frac{2mn}{m^2-n^2} (1 - (-1)^{m-n})\). So the matrix representing our momentum operator is:
\begin{equation}
\begin{pmatrix}
 0 & \frac{8i\hbar}{3L} & 0 \ldots \\
\frac{-8i\hbar}{3L}& 0 & \frac{24 i\hbar}{5L}&\ldots \\
0 & \frac{-24 i \hbar}{5L} & 0 & \ldots \\
\vdots & \vdots &  \vdots  & \ldots
\end{pmatrix}
\end{equation}

So as we can see  our matrices are really infinite and they act on our infinite dimensional column.
We have now laid the ground for Dirac notation. We have seen that in this simple problem there is a representation which can be arrived at where we think of ourselves as living in some sort of vector space where our operators(matrices) and our vectors can be infinite dimensional.
<script type="text/javascript">MathJax.Hub.Queue(["Typeset",MathJax.Hub]);</script>
 

Transition from wave functions to Dirac notation

We begin as one would expect with Schrodinger's equation in one dimension; which is
\( -\frac{\hbar}{2m}\frac{d^2 \psi }{dx^2} + V \psi = \frac{i \hbar \partial \psi}{\partial t} \) . If we solve the solve this using separation of variable we get a term that is only time dependent and the other that is only spatially dependent. The steps to this abound everywhere so I shall not go into it (It would not be in the spirit of this blog).  But we could introduce what is called the Hamiltonian operator and rewrite the time dependent Schrodinger's equation as follows:
\begin{equation}
H\psi = i \hbar  \frac{\partial \psi}{\partial t}
\end{equation}
where H =  \(-\frac{-\hbar}{2m}\frac{ d^2 }{dx^2} + V  \).
The spatially dependent part would look like this; (If we write in terms of our newly defined Hamilton, H)
\begin{equation}
H\psi = E \psi
\end{equation}
Notice that E(energy) is just a number and this looks like some eigenvalue problem which in fact it is. So if we have a system with different energy levels , these would correspond to different eigenfunctions.  So in that case our eigenvalue problem would look like this:
\begin{equation}
H\psi_n = E_n \psi_n
\end{equation}
The next statement is very important and will be stated without proof. The \(\psi_n \) live in a vector space that physicist refer to as the Hilbert Space, not only do they live in it but also can form a basis set for it (after of course they have been normalized). In other words if we have some solution to the schrodinger's equation, we can express it as a linear combination of these  \(\psi_n \). That is to say that:
\begin{equation}
 \Phi  = \sum_n c_n \psi_n
\end{equation}
where \( \Phi(t) \) is a solution to Schrodinger's equation.
It is here that the student may think his or herself, " well that's neat, what is the next topic?"  Pose for a while and contemplate this, we have a infinite dimensional vector space, spanned by these \( \psi_n \). A little but deep thought could occur to one: why represent the vector by writing the tedious sum expressed in the above equation? Could we not just use the \(c_n \)  place them in a column vector (which in this case would be infinite dimensional) and then our operators would be matrices rather than differential operators. After all matrices act upon column vectors.  This is the key insight to Dirac notation. Our quantum state will no longer be written by some function but instead it will be some column vector, with the elements in it the \(c_n\).
In other words or symbols:
\begin{equation}
 \Phi =  \begin{pmatrix}  c_1\\ c_2 \\ c_3 \\ \vdots \end{pmatrix}
\end{equation}

NOTE: If we are in an infinite dimensional space we have to worry about convergence but at the most we live in a textbook world everything is well behaved. The next step is to see the transition from differential operators to matrices explicitly.

A late realization.(Dirac Notation)

Since I would like to have some sort of structure to the posts on this blog. I realized that I jumped into Density matrix theory without ever talking about Dirac Notation which was heavily used. Since I would like this blog to be completely self contained to the curious reader, I shall stop , back up and slowly introduced Dirac notation and may be in the process introduce Quantum Mechanics.

Density Matrix Theory: An Introduction III

One of the main reason from introducing the density matrix is that it allows to describe statistical mixture as opposed to merely pure states. We therefore have to make some comments as to how we generally extract out information about a system from quantum mechanics. We know from introductory quantum mechanics that the wave function has all the information we could want or that quantum mechanics can provide.
Let's for a moment talk more generally about what is knowable in quantum mechanics. We have all heard about Heisenberg's uncertainty principle, namely that we can't know position and momentum exactly simultaneously. One often hears this in public discussions, but a more precise statement and in fact a more correct statement is that we can't know the position in a certain direction and the momentum in that same direction simultaneously to arbitrary accuracy. The mathematical version of this statement would be  that \([ \hat{q_i},\hat{p}_j ] = i\hbar \delta_{ij}\) . Here \(q_i\) and \(p_i\) are generalized co-ordinates and generalized momentum respectively. So the more familiar version would be \( [\hat{x}, \hat{p_x}] = i\hbar \).
The underlying reason is that the position operator in the \(\hat{x}\) does not commute with the momentum operator in the same position. This implies that I can indeed now the position exactly in the \(\hat{y}\) and the momentum for example in the \(\hat{x}\) .
Now to jump to more general statements, we may proceed as follows. If we have a set T=\( ( Q_1,Q_2,Q_3 ...Q_n) \)  consisting of operators any two of which commute with each other  and a correponding set U = \((q_1, q_2,...q_n)\)  of eigenvalues for each operator ( we assume that this sets are the largest possible size they can be so the ket  \(| q_1, q_2,q_3,...q_n \rangle \) describes our system), then these two sets represent the "maximum knowledge" of the system that we can have. These states of "maximum knowledge"  are what we called the pure states earlier.
One must keep in mind that these two sets may not be unique and in fact  are rarely  reproducible in experiments. Instead what we usually have is that we know with a certain probability  \( W_n \) that our system is in a pure state \( \psi_n \). We therefore deal with statistical mixtures more than we do with pure states.

NOTE: Do not confuse a superposition of states with a statistical mixture. With a superposition there is a phase relationship between the quantum states, we can thus write down a pure state wave function describing a quantum superposition of these two states. On the other hand, statistical states  \( \psi_n \) have no phase relationship. To put it another way, if I have a state and I can't tell in principle whether it is spin up or spin down until I measure then I have a superposition but if I have a number of states that are either spin up or spin down or superposition of the two (so there is probability of me picking up a state with either spin or spin down or superposition of the two) then I have a mixed state.

In quantum mechanics we connect with the "real world" by calculating expectation values  of some Operator. So for example, we continuously prepare a state in some specific configuration and keep measuring its energy. At the end we shall have a probability of having the state with energy \(E_n\)  in the corresponding eigenfunction  \( \psi_n\). We can calculate the expectation value for a pure state (in braket notation) like this:
\begin{equation}
\langle O \rangle = \langle \psi| O | \psi \rangle
\end{equation}
and for mixed states
\begin{equation}
\langle O \rangle = \sum _n \langle \psi_n| O | \psi_n \rangle
\end{equation}

Density Matrix
We can now describe our density matrix of our statistical mixture in the following manner:
\begin{equation}
 \rho = \sum_n W_n |\psi_n \rangle \langle \psi_n|
\end{equation}
Here \( W_n \) are the statistical weights and\( |psi_n \rangle\) are the independently prepared states. These  independently prepared states are not necessarily orthonormal  so we could re-write in terms of states that are so that :
\begin{equation}
\psi_n = \sum_m a_m^n \phi_m
\end{equation}

We are now in a position to re-write our density matrix in the following manner:
\begin{equation}
 \rho = \sum_{nmk} W_n a_m^{n} a_k^{n*} |\phi_m \rangle \langle \phi_k|
\end{equation}

Thus we are able to put matrix in density matrix by find how to get its matrix elements. Using the orthogonality condition of the \( \phi_n \) state we observe that:
\begin{equation}
\rho_{ij} = \langle \phi_i |\rho | \phi_n \rangle = \sum _n W_n a_j ^{n}a_j^{n*}
\end{equation}

Introducing these orthonormal states allows us to arrive at a simple expression or the expectation value for an operator of a mixed state. We proceed as follows:
\begin{align}
 \langle O \rangle &= \sum _n \langle \psi_n| O | \psi_n \rangle \\
                         &= \sum_{mk} \sum_{n} W_n a_m^{n} a_k^{n*} \langle \phi_k |O | \phi_m \rangle \\
                           & = \sum_{n}\langle \phi_m|\rho | \phi_k \rangle  \langle \phi_k |O | \phi_m \rangle  \\
                            &= tr (\rho O)
\end{align}

Look at the power of the density matrix. If I know the density matrix (remember something I can define for pure or statistical mixtures), I can arrive at the very thing I need to connect with  the "real " world namely the expectation of the operator or observable in question. This something that I can't easily get if I stick with the wave function and have a statistical mixture in my lab.

Density Matrix Theory: An introduction II

We have now introduced the polarization vector and how to calculate its components. To see the relationship between it and the density we shall move as follows:
Remember that if we have  a pure state then we can assign a state vector to our system but if we have a mixture then we instead have to consider a statistical mixture. We can thus write the density matrix for a statistical mixture in the following manner:\chi
\begin{equation} \label{eq.1}
\rho = \sum_k W_k |\chi \rangle  \langle \chi | \hspace{10mm}  \text{eq.1}
\end{equation}
where \(W_k \) is the statistical weight. From this I hope one can see that a statistical mixture is a generalization of a pure state. Let's multiply the above equation by the pauli matrix \( \sigma_i \) to get:

\begin{equation} \label{eq.2}
 \rho \sigma_i =  \sum_k W_k |\chi \rangle  \langle \chi |  \sigma_i   \hspace{10mm}  \text{eq.2}
\end{equation}
We shall now take the trace of each term i.e the \(tr (\rho \sigma_i \) to get:
\begin{equation} \label{eq.3}
tr ( \rho \sigma_i ) = \sum_k W_k \langle \chi |  \sigma_i | \chi \rangle = P_i   \hspace{10mm}  \text{eq.3}
 \end{equation}

Now, the pauli matrices in combination with the identity matrix form a basis for a 4 dimensional space. Thus if we stick with two level systems, then we can write our density matrix as:
\begin{equation}
 \rho = a_o I + \sum_i a_i \sigma_i
\end{equation}
 We can the use the following easily derivable facts to find the \( a_o , a_i  \)constants:
 1) pauli spin algebra- \(\sigma_i \sigma_j = i \sum_i \epsilon_ {ijk} \sigma_k  + \delta_{ij}I\),
 2) the \(tr(\sigma_i) =0 \)
 3)\( tr (\rho ^2 )= \frac{(1+ P^2)}{2}\)

Using  2) and 3) in combination tells us that \(a_o \) is actually 1/2. Multiply our density matrix again by  \(sigma_i \)  and taking the trace and then using 1) and 2) gives us that
\begin{eqnarray}
 tr ( \rho \sigma_ i) &=& 2 \sum_i b_i \delta_{ij} \\
                              &=& 2 b_i
\end{eqnarray}
Now referring back to eq.3  we see that  \( b_i  = \frac{P_i}{2} \)
But now let's take a close look at eq.1 in conjuction with fact 3). In general eq.1 will be less than 1 since the statistical weights themselves are less than one. Which means that in general  \( \rho^2 \) for  a statistical mixture will be less than one while if we have a pure state then \(  \rho^2\) is equal to one. Now we finally have a mathematical way of distinguishing between a density matrix for a pure state and a density matrix for a statistical mixture.

In our discussion of the density matrix we considered states \( | \chi \rangle \) . These states were not necessarily orthonormal. But we could re -write the density matrix in terms of orthornormal states. That is to say we can first re-write the \( | \chi \rangle \) like this:
\begin{equation}
 |\chi  \rangle  =  \sum _n a_n | \psi \rangle
\end{equation}

So if we then re-write the density matrix (which is in general represent a mixture) we have that:
\begin{equation}
 \rho = \sum_k W_k a_k a_k^{*} | \psi \rangle \langle \psi |
\end{equation}

This last equation gives us a hint of a much more general formulation for the density matrix with the eventual goal of considering open quantum systems and the density matrix that may describe them.

Density matrix theory: An Introduction, Preliminary Considerations

A very powerful and useful instrument that will be discussed is the density matrix. Consider the outer product of a ket and a bra. i.e \( |\psi \rangle \langle \psi | \). This might look like a needless complication but as we shall see this entity tells us quite a bit of information.
Before we start with the wonders we shall start with the Stern Gerlach experiment. A beam of electron is sent into a magnetic field pointing in the \(\hat{z} \) direction. What was expected was that electron would be bent or deflected in all directions instead what was seen was that either the electrons were deflected up or down. This was the first indication of what physicists now call spin. I shall not go on about this experiment since it is one that is rehashed in almost every modern physics textbook or quantum mechanics text book. Instead what shall is consider the following scenario:

Supposing a beam of light is let go from a source and it is able to move through the Stern Gerlach instrunmentation by turning the instrument in some direction, then we shall call it a pure state. Note that we have not specified a direction, we have merely said a direction or orientation can be found in which the beam entirely goes through.
This is important because if  particles  in the  beam  are in a superposition of spin up and spin down, as long as all the particles are in this state then this beam is considered to be in a pure state. If this is the case we can then describe the system by a single state vector.

If one the other had we have beam of particles in one pure state and another beam of particles in another pure state, then this ensemble is considered to be in a mixed state. In order for this to occur the beams have to be prepared independently and by this we mean that no phase relation exists between these two beams.

Now we move to introduce the polarization vector which shall be defined as follows:
$$ P_i = \langle \sigma_i \rangle =\langle \chi |\sigma_i | \chi \rangle  $$
If we are talking about a two level system then the sigma are the pauli spin matrices. Before we continue, let's get an inkling as to why it is called the polarization vector. If we consider the general two dimensional state vector \( | \chi \rangle = \cos \theta + e^{\delta}\sin \theta \)  and then calculate \(P_i\), one gets the following column vector \(  (\cos \delta \sin \theta, \sin \delta \sin \theta, \cos \theta)^{T} \) This should remind one of spherical co-ordinates. So we picture the bloch sphere the polarization is telling us the "composition" of our state vector. Hence how it is polarized, in very much the analogous way we might ask how the electric field for example is polarized. In fact there is an analogy to be made between filter EM waves and spins the Stern-Gerlach experiment.

The definition above should just remind one of a pure beam, well can we generalize it for mixtures and as a point in fact it can, namely:

$$ P_i = \sum _a W_a \langle \chi_a |\sigma_i | \chi _a\rangle $$

Here \( W_a \) are simply the statistical weights namely the proportion of particles in state \( \chi_a \)

The next step is to connect this with density matrix.

The Schrodinger Equation and The Form of the Solution

One might wonder Schrodinger's equation is in the form it is in. It is called a wave equation but we are not using the wave equation. Well, I can't claim to understand what I am about to expostulate( in other words give an intuitive understanding) but rather take the cope out way and argue with a combination of mathematical and experimental reasoning.

Our story begins with an electron going through the double slit and hitting a screen somewhere behind it. The surprising result was that an interference pattern was observed as a series of these electrons passed through one at a time. The obvious explanation is that the electron has a wave property so we quickly write down our trusty sine function( or cosine if you are inclined that way). Now the next step is to get the interference pattern through our results. So we do the naive thing and write down two sine functions differing by a phase factor and add them together . One must recall that we want an interference pattern that has a maximum at the center and decreases "sinusoidally" in amplitude as we go away from it from either side.  Unfortunately when we go with our naive guess and square our wave function(addition of the sines) since E.M tells us that the intensity is proportional to wave function squared our interference is sinusoidal but does not decrease as we go away from the center. thus the interference pattern is still illusive.

In fact we realize something worse. Supposing two waves are approaching each other and let's stick to one dimension(we are modelling two electrons approaching each other). We expect a combination of destructive and constructive interference . So we re-write the sine functions in this manner  i.e \( \sin (kx - \omega t) + \sin (kx + \omega t) \); the result is this \(2 \sin kx \cos \omega t \). But this is a terrifying result, at multiples of t = \( \frac{ \pi}{2 \omega} \) the particles disappear everywhere. This sort of tom-foolery is allowed for the Cheshire cat but nor for our dear electrons. So we want waves but adding sines or cosines together does not work!

Finally we remember something about Euler's identity and it being related to waves. So we decide to represent our electrons by this wave function \( e^{(ikx - i\omega t)} \).  We add one coming from left \( e^{(ikx - i\omega t)} \) and one coming from the right  \( e^{(-ikx - i\omega t)} \) and get \( 2 e^{(-i\omega t) } \cos \omega t\). Ah! the function is still sinusoidal and never disappears , but with a weird complex amplitude. Remember our goal is to describe the interference patterns on the screen, but these patterns have real amplitudes. Before we move on we take note that we either take the wave function \( e^{(ikx - i\omega t)} \) or  \( e^{(-ikx + i\omega t)} \) but not both. Reason is if we take both and add them then we are back to cosine and sine functions which have all the wrong properties.  So we shall take the former and assume the amplitude is A.

Back to the puzzling imaginary amplitude. This is no cause for concern because remember the intensity is the square of the wavefunction, which means taking the modulus of our wave function.  First let's add two wave functions that will describe our electron (since we are looking for an interference pattern):\(A_1 e^{(ikx_1 - i\omega t)}+A_2e^{(ikx_2 - i\omega t)} \)
Now let's look at the modulus

\begin{align}
&\left(A_1 e^{(ikx_1 - i\omega t)}+A_2e^{(ikx_2 - i\omega t)}\right) \left(A_1 e^{(-ikx_1 + i\omega t)}+A_2e^{(-ikx_2 + i\omega t)}\right)\\
&\left(A_1 e^{ikx_1}+A_2e^{ikx_2 }\right) \left(A_1e^{-ikx_1}+A_2e^{-ikx_2 }\right)\\
&\text{a bit of algebra....}\\
&A_1^2 + A_2^2 + 2A_1A_2 \cos k(x_2 - x_1)
\end{align}

The formula above is quite astonishing, if we close slit 2 in our experiment and assume our wave function is the right one while holding firm to the stipulation that its modulus is the thing that matters(avoiding complex amplitudes) then we expect \(A_1^2\) namely just a mound of electrons in front of slit 1 and similarly \(A_2^2\) if we close slit 1 and leave slit 2 open. In both of these situations the cosine term is not there. Once we open both slits the interference term appears and we get our interference pattern we were looking for.

Notice the classical guess from E.M gives a cosine term squared for the intensity, which of course does not explain the places on the screen where there are no electrons hitting.

Posts on Quantum Mechanics

I recently re-looked at my Quantum Mechanics text-book by Ohanian, and I was rather impressed by his style of presentation of the material. There will be a new tab with my examination and notes on the textbook. It offers a nice complementary look to Griffiths introduction. His method avoids solving partial differential equations and instead opts for the method of operators throughout. The first post will be examining the wave picture and the matrix picture and showing their equivalence. 

The Summa : Euler's Rigid Body Equations

We finally arrive at the top, namely the Euler equations for rigid bodies. I mentioned the inertia tensor in previous section, but I have decided to eschew from it since there are a multitude of good resources on it. What matters to us at this point is that it is a matrix and not only that but a symmetric matrix. Given that this is the case, we can always diagonalize it so that it simply has entries on the diagonal. These entries we shall call \( I_x, I_y, I_z \). These are the so called principal axes. Now if we jump into a co-ordinate system which coincide with these axes then the angular momentum L = \( \left( I_x\omega_x, I_y \omega_y , I_z \omega_z \right) \). The angular moment as observed in the inertial reference is $$ \dot{L_{i.f}} = \dot{L_{r.f}}+ \Omega \times L = \Gamma $$ where i.f is "inertial reference frame" , r.f is "rotating reference frame" and gamma is the external torque. If we plug in what L is we arrive at the blessed equations
\begin{align}
\Gamma_x & =  \dot{L_x} - (I_y - I_z) \omega_z \omega_y\\
\Gamma_y &=   \dot{L_y}- (I_z - I_x) \omega_y \omega_z\\
\Gamma_z &=  \dot{L_z}-(I_x - I_y) \omega_x \omega_ y \\
\end{align}

Notice that the time derivate does not affect I, that is the purpose of using these principal axes. Now supposing external torque is zero, then the angular moment L should be constant. Here is the proof:
1) If external torque is zero, angular momentum is constant

Multiply with ith equation above by  \( I_i\omega_i\) and add all three equations
\begin{align}
I_x^2 \omega_x \dot{\omega_x} + I_y^2  \omega_y \dot{\omega_y}+ I_z^2  \omega_z \dot{\omega_z} &=  \left(I_x(I_y - I_z)  + (I_z - I_x)I_y + (I_x - I_y)I_z \right) \omega_z \omega_y\omega_x \\
I_x^2 \omega_x \dot{\omega_x} + I_y^2  \omega_y \dot{\omega_y}+ I_z^2  \omega_z \dot{\omega_z} &= 0\\
L_x\dot{L_x}+ L_y\dot{L_y}+L_z\dot{L_z} &= 0\\
 \frac{d}{dt}\frac{\left( L \cdot L\right)}{2} &= 0\\
\frac{d}{dt}\frac{L^2}{2} &= 0\\
\dot{L} &= 0 \hspace{20 mm} QED
\end{align}

2) If we have a lamina rotating freely(no torques) about a point O of the lamina. The component of \( \omega \) in the plane is constant.

Again since we have no torque we start from Euler's equations with no torque, but we use a slight of hand which should not look like one since we have seen it before, namely use the perpendicular axis theorem. This gets rid of the \( I_x \) and \(I_y\)  but we are getting ahead of ourselves:
\begin{align}
I_x \dot{\omega_x} &= (I_y - I_z) \omega_y \omega_z\\
I_y \dot{\omega_y} &= (I_z - I_x) \omega_z \omega_x\\
I_z \dot{\omega_z} &= (I_x - I_y) \omega_y \omega_x\\
\end{align}
Use perpendicular axis theorem

\begin{align}
\frac{I_z}{2} \dot{\omega_x} &= (\frac{I_z}{2} - I_z) \omega_y \omega_z\\
\frac{I_z}{2} \dot{\omega_y} &= (I_z - \frac{I_z}{2}) \omega_z \omega_x\\
I_z \dot{\omega_z}& = (\frac{I_z - I_z}{2}) \omega_y \omega_x\\
\end{align}
The new set of equations
\begin{align}
\dot{\omega_x}& = -  \omega_y \omega_z\\
 \dot{\omega_y} &= \omega_z \omega_x\\
\end{align}
Multiply the top equation above by \(\omega_x\) and the bottom by \( \omega_y\) and then add the two to get \(\omega_x \dot{\omega_x}+ \omega_y\dot{\omega_y} = 0 \). But this is merely \(\frac{d}{dt}\left( \omega_x^2 + \omega_y^2\right) \)= 0. Voila

A few solved problems

Since this blog is meant to concentrate on these not emphasized in textbooks or offer a slightly different way of looking at old things, we now come to place where I think a great deal of material has been written on the topic and most of it quite good. Now, we derive Euler's rigid body equations in all their glory, I found a few interesting problems, which relate to rotational motion.

1) You are to imagine a rod pivoted to a horizontal beam but free to swing in a plane perpendicular to the plane containing the beam. The rod has zero width, and a density not necessarily constant. Its moment of inertia about an axis parallel to the beam is I and the center o f mass is located distance b from the pivot.  Find the equation of motion for small swings, its period and the length of a simple pendulum that would have the same period.

This is called a compound pendulum( contrast it with the simple pendulum where all the mass is concentrated at a point). Let \(\Gamma\) be the torque and L is the angular momentum.
\begin{align}
 \dot{L}&= \Gamma \\
  I \dot{\omega}& = - bmg\sin \theta \\
  \ddot{\theta} &= -bg \theta m \\
 \ddot{\theta} &= \frac{-bg\theta m}{I}\\
\end{align}

It should b:e obvious that the \(\omega\) or angular frequency is\( \sqrt{\frac{mbg}{I}}\) with period \(2 \pi\sqrt{\frac{I}{mbg}}\). Now the angular frequency of a simple pendulum is \( \sqrt{\frac{g}{l}}\). Setting the expressions for the different angular frequencies equal gives that the length of the simple pendulum must be \(\frac{I}{mb}\).

2) Imagine the same rod as above but this time it is hit with a horizontal force a distance c from the pivot with impulse F\(\Delta \)T = \(\xi\). What is the angular momentum delivered to the pendulum, pivot just after it has been hit? What is special distance \(b_o\) such that the impulse delivered to the pivot is zero?

Mathematically, I do not think this is a challenging problem. It is the conceptual physics involved that makes it quite challenging. So here is my solution:

We can treat the pendulum as though all its mass is centered at the center of mass. Now since it is not moving before we hit it, the initial velocity is zero. But this implies that the impulse is equal to the moment just after we hit it.(Just consider newton's second law for discrete times). Thus the angular momentum after it has been hit is simply \(b \xi\). Now the impulse delivered to the pivot will be m\(v_p\) where \(v_p\) is the velocity at the pivot. The velocity at the pivot is the velocity of the center of mass \( v_{cm} - v_{h}\) where \(v_h\) is the tangential velocity at the point  it is hit. Therefore  impulse is : \( mv_{cm} -m c\omega\) = \( \xi - \frac{mbc\xi}{I}\). Remember L = \(I\omega\),  c is the distance from the pivot to where it is hit and b is the distance from pivot to center of mass. If we want to find that special distance \(b_o\)  where impulse at the pivot is zero just set the expression for the impulse equal to zero and you will get that \( b_o = \frac{I}{mb}\)

Moments of Inertia and conclusion (Part II)

We continue our journey with moments of inertia by calculating the I for a solid cone about two different axes i.e through the axis that goes through the center of mass and the axis that goes perpendicular to its pointy edge.

5) Axis through the center of mass
\begin{align}
I_z &= \rho \left( \int r^2\, dV \right)\\
  &=  \rho \left( \int r^3 \, d\phi dr dz \right)\\
  &=  \rho  2\pi \left(\int_0^{\frac{zR}{L}} r^3 \, dr \,\int_0^ L\, dz\right)\\
 &= \frac{\pi \rho R^4 L }{10} = \frac{3MR^2}{10}
\end{align}

Axis that is perpendicular to the point edge.
Now this interesting, and it worth slowing down the pace and considering. This problem can be done two ways: One is to realize that if we put in the cone upside(not necessary) then this axis is \(I_{yy}\) in the moment of inertia tensor; now calculate away the integrals. But supposing you know nothing about the inertial tensor, is there hope? Indeed, but it requires a little foray into the perpendicular axis theorem.

Let's suppose you have a two dimensional figure(located in the xy place). Then the perpendicular axis theorem says \(I_z\) = \(I_x + I_y\). The proof is quite simple and goes as follows:
\begin{align}
 I_z &= \left(\int r^2 dm \right) = \left( \int x^2 + y^2 dm\right)\\
       &= \left(\int x^2 dm\right) + \left( \int y^2 dm\right) = I_y + I_x \hspace{20mm}QED\\
\end{align}
Its most common application is for example finding \(I_x\) or \(I_y\) of a disk. Since we know \(I_z\)  because it is easy to calculate, then finding  \(I_x\) or  \(I_y\)  is trivial. Again you do not need integrals, the theorem saves you.
 This theorem  suggest another way of calculating our desired results. Here is the strategy:
1) There is symmetry in the x- y plane so \(I_x\) and \(I_y\) are the same and we calculated \(I_z\) in the problem above(axis through center of mass)
2) Starting from the perpendicular axis theorem we can add \(I_z\) to both sides to get that \(2I = I_x+I_y+I_z\). Since we know \(I_x = I_y\) then \( 2I_y + I_z = 2I\). Now all we need to do is calculate 2I. This is a generalization of the perpendicular axis theorem.

\begin{align}
 2I &= 2\rho \left( \int (r^2 \cos^2 \phi + r^2 \sin^2\phi + z^2 ) r dr d\phi dz \right)\\
      & = 4\pi \rho \left(\int_0^L dz \int_0^{\frac{zR}{L}}dr(r^3 + z^2r) \right)\\
      & = 4\pi \rho \left( \frac{R^4L}{20}+ \frac{L^3R^2}{10}\right)\\
    & = 12M\left( \frac{R^2}{20} + \frac{L^2}{10} \right)\\
& \text{plugging into previous results for } I_z \text{ we get } \frac{3M}{5}\left(  \frac{R^2}{4} + L^2\right)
\end{align}

 This theorem is very useful if you are in three dimensions and have symmetry in two or more axes. For example take the moment of inertia for a spherical shell, \( I_z, I_y \) and \(I_x\) are all the same so we have that \( 3I_z = 2 I\)  and I is easy to calculate. "I"here is the moment of inertia about an axis in the \( 4^{th}\) dimension!!! Poooffff ( Your mind being blown away)

Finding center of mass and moments of Inertia

The ultimate goal is still to arrive at the Euler equations for rigid body dynamics. I often find that equations are neatly derived in books and in class but they often prove to be almost useless to the student because students are not comfortable in the actual calculation of the quantities in the equations. We so far have introduced rotational motion and gone through it in some detail. Thus the idea of jumping into the rotational reference frame in order to use these equations should not be foreign. The next step is to learn how to find center of mass, and moments of inertia of various shapes. The goal is to try and make the derivations as straight forward as possible, which means using brute force methods. We shall leave the trickery and elegance for pretentious physics professors and authors. I am a simple minded man and I hope you are too.

Our first goal is to be comfortable finding the center of mass, then we shall move to finding the moment of inertia. Along this journey we shall meet interesting concepts like the parallel axis theorem and the less known perpendicular axis theorem.
First let's begin simply, we shall find the center of mass for a 1D line, cylinder, sphere and a cone. The equation is easily recalled namely: \(\int_{V}\vec{r}dm\). The important thing to remember is to take \(\vec{r}\) seriously. This is the mistake I have often made in the past, because \(\vec{r}\) is not in cartesian co-ordinates, before we take the integral this will have to be changed to Cartesian co-ordinates. This point becomes very important when we go to 3D dimensions.
1) 1D line.
  \begin{align}
   &\frac{\rho}{M}\left( \int_0^L x dx \right) \\
 &  \rho \frac{x^2}{2}|_{0}^L = \rho\left(\frac{L^2}{2}\right) = L/2\\
&\rho \text{ here is }\frac{M}{L}
 \end{align}

2) Cylinder
  \begin{align}
 &\frac{\rho}{M}\left( \int x \, dV \, \hat{x} + \int y\, dV\, \hat{y}+ \int z \, dV \, \hat{z} \right)\\
&\frac{\rho}{M}\left(\int r\,\cos\phi\, r d\phi dr dz\, \hat{x} + \int r\,\sin\phi\, rd\phi dr dz\,\hat{y}+ \int_0^L \int_0^R\, z\,r\,dr d\phi dz\, \hat{z} \right)
\end{align}
The first two integrals  will  go to zero since the integral of \(\sin \phi\) and \(\cos \phi \) over 2\(\pi\) radians is zero. So all we need  is to concentrate on the last integral.
\begin{align}
&\frac{\rho}{M}\left( \frac{R^2}{2} 2\left(\pi\right) \frac{L^2}{2}\right)\\
&\text{put in the expression for the density of a cylinder  to get} \frac{L}{2}
\end{align}
With the cylindrical example it is now clear what people mean when they say,"by symmetry, the center of mass is on this or that axis". It is important that one knows what occurs when options are eliminated by reasons of symmetry. This is why I am doing all the integrals, one has to justify one's intuitions however right they may be.
3) Spherical Shell
\begin{align}
&\frac{\rho}{M}\left( \int R\,\cos \phi \,R^2\, (\ sin\theta)^2\, d\phi d\theta\, \hat{x} + \int \, R\, \sin \phi \, (\sin\theta)^2 \,  R^2 d\theta d\phi \hat{y} + \int  R \cos \, R^2 \, \sin \theta\, d\theta d\phi\, \hat{z}\right)
\end{align}
I have already plugged in what x, y and z are in spherical co-ordinates to get the above expression
Again, the first two integrals go to zero and doing the third integral also gives you zero, i.e the center of mass is at the origin.
4)Solid Cone
\begin{align}
\frac{\rho}{M}\left( \int \, x \, r d\phi dr dz \, \hat{x}+ \int y \, r d\phi dr dz \, \hat{y}+ \int z\, r\, drd\phi dz \, \hat{z}\right)
\end{align}
Again the first two integrals disappear (could have argued this from the beginning by using symmetry arguments).  Now r depends on z the integrals over r and z can't be done completely independently, we will need to write r in terms of z. This will be done by considering the limits. If the cone is balance on its pointed edge then r and z are related in the following manner \(r = \frac{zR}{L}\) . Thus the last integral is the following $$ \frac{\rho}{M}2\pi \left(\int_0^L z\, dz \int_0 ^{\frac{zR}{L}}r\, dr \right)$$ The integral gives \( \frac{\rho \pi R^2  L^2}{4M}\) and plugging for the density gives \(\frac{3L}{4}\)

All other shapes usually found in  text books can be attacked from these examples by playing around with the limits.

Moments of Inertia

Again it is quite easy to remember the general equation which is \(\int r^2 \, dm \) but as with most thing s the devil is in the details. The above form is slightly misleading if you have had vector calculus. One is led to be believe that r is from the point dm to the origin but this is wrong. The r in the equation is from the dm to axis of rotation. Do not get the impression that calculation for moment of inertia is like the calculation for center of mass except the r is squared that could not be further from the truth. This point will be illustrated when we calculate the moment of inertia for a solid sphere.

Secondly, it matters what our axis of rotation is, you have to know that before you begin writing down integrals. This point will be illustrated by deriving the parallel axis theorem.

Let \(r_i\) be the distance from some origin to a \(m_i\) and \(r_i^{'}\) be the distance from the center of mass to that same mass and lastly let \(d\) be the distance from the center of mass to the fore mentioned origin. Then the moment of inertia about this origin is clearly, \(\sum r_i^2 m_i\) which is equal to  \(\sum (r_i^{'}+d)^2 m_i\) since r = r' + d.
\begin{align}
     \sum r_i^2 m_i &= \sum  r{'}^2 m_i + \sum d^2 m_i + 2d \cdot \left( \sum r{'}m_i\right) \\
                            &= \sum r{'}^2 m_i  +d^2M +  2d \cdot \left( \sum r' m_i \right)
\end{align}
Look at the last term on the right in parentheses. It is weighted term defined in order that it be zero. Look at the definition of center of mass and  put all the terms on one side  and you will see what I mean. We thus arrive at the parallel axis theorem. Whenever you find yourself trying to calculate the moment of inertia about some weird axis. STOP! do not just write down integrals. It might be easier to do the integrals about the axis that goes through the center of mass and then use the parallel axis theorem.

1) I of line about the end and about center of mass
\begin{align}
I &= \int r^2 dm  = \rho \int r^2 dV \\
  &= \rho \int_0^L x^2 dx  = \rho\frac{L^3}{3}\\
  &= \frac{ML^2}{3}
\end{align}
Now suppose we want to find the moment of inertia about the center of mass, do we re-do the integrals. NO! we apply the parallel axis theorem: \(I = I_{cm}+ I_d\)  we just calculated I and the end is separated from the center of mass by length d. Plug and chug to get \(\frac{ML^2}{12}\).

2) I of cylinder about axis
\begin{align}
I &= \int r^2 dm  = \rho \int r^2 dV \\
  &= \rho \int r^3 \, dr d\phi dz\\
  &= \frac{\rho \pi L R^4}{2}\\
 &= \frac{MR^2}{2}
\end{align}

3) I of solid sphere about axis (any really)
  \begin{align}
  I &= \rho \int (r \sin \theta)^2 \, r^2 \sin \theta \, dr d\phi d\theta \\
     &= \rho \int (\sin \theta)^3 r^4 \, d\theta d\phi \\
     & = \frac{\rho R^5 8\pi}{15}\\
    & = \frac{2MR^2}{5}
  \end{align}
Notice the \(r \sin\theta\)  in the first integral which would not have been there if r was the usual r in spherical co-ordinates and also I have not broken up r into its components x, y and z as in the center of mass case.

4) I of a Spherical shell
   There is a rigorous way of assuming it has a finite but appreciable thickness , finding the moment of inertia then taking the limit  of the thickness dropping to zero. I will not do that instead I shall I assume it has a finite mass in an infinitesimally thin space. Hence,
\begin{align}
  I &= \rho \int r^2 dA\\
  &= \rho \int (R \sin \theta)^2\, R^2 \sin \theta \,d\phi d\theta \\
  &= \rho \int (\sin \theta)^3 R^4 \, d\theta d\phi\\
 I &= \frac{2MR^2}{3}\\
\text { remember} \rho \text{ is }\frac{M}{4\pi R^2}
\end{align} 

Non-inertial reference frames (example calculations)

In this section we shall do the usual example of free-fall motion but the distinction will be the use of spherical co-ordinates to ease the use of mental geometric gymnastics and also a classical introduction to Larmor precession will be introduced. Because the earth is for the most part a sphere, spherical co-ordinates should be used with only very minor modifications. Gravity points in the -\(\hat{r}\) direction the \(\hat{\theta}\) is going from north pole to south pole and \(\hat{\phi}\) is going along a latitude. Ultimately we want to describe the motion as a function of latitude so the simple transformation of \(\frac{\pi}{2}- \theta = \lambda\) will be introduced ; \(\lambda\) is the latitude defined as being zero at the equator and goes to \(45^{\circ}\) in either direction.

First a note on unit vectors:
The \(\hat{r}\) is the one commonly introduced in text books and is reproduced below along with \(\hat{\theta}\) and \(\hat{\phi}\)

\begin{equation}
\hat{r}= \cos \phi \sin \theta \hat{x}+ \sin\phi \sin \theta \hat{y} + \cos \theta \hat{z}\\
\hat{\theta}= \cos \phi \cos \theta \hat{x}+ \sin\phi \cos \theta \hat{y} - \sin \theta \hat{z}\\
\hat{\phi}= -\sin \theta \hat{x}+  \cos \theta \hat{y}
\end{equation}

From these three equations one can easily find the inverse transformations that give \(\hat{x}\), \(\hat{y}\) and  \(\hat{z}\) in terms of the spherical unit vectors in this manner. To find what  \(\hat{x}\) is simply collect the terms next to the  \(\hat{x}\) in the above equations and group them together.
So  \(\hat{x} = \cos \phi \sin \theta \hat{r} + \cos \phi \cos \theta  \hat{\theta} -\sin \theta \hat{\phi}\). Same applies for the rest.
Also, if you know what the \(\hat{r }\) unit vector looks like in Cartesian unit vectors then \( \hat{\theta} \)is not very far away. First note that it is perpendicular to the \(\hat{r}\) unit vector and therefore while we had to project the \(\hat{r}\) onto the x-y axis by taking the \(\sin \theta \) in order to project the \(\hat{\theta}\) we use \(\cos \theta\). Since \(\hat{\theta}\) is going in the -\(\hat{z}\) direction we take the  -\(\sin \theta\) .

So why is all this important? Well, our \(\Omega\) or angular velocity will be in the \(\hat{z}\) direction and we will need to put it into spherical unit vectors. But from the above discussion, rather than playing with right triangles and twisting and contorting your head in unseemly manners. Simply collect the \(\hat{z}\) in the \(\hat{r}\) and  \(\hat{\theta}\) equations. So we have the following:
\begin{equation}
 \mathbf{\Omega \hat{z}} = \Omega \cos \theta \hat{r} - \Omega \sin \theta \hat{\theta}
\end{equation}

Imagine if you were on a test, you just saved yourself five minutes of agony. Now since \(\hat{\theta}\) moves in the -\(\hat{z}\)  we will need to re-write \(\mathbf{\Omega}\) in the following manner:

\begin{equation}
 \mathbf{\Omega \hat{z}} = \Omega \cos \theta \hat{r} + \Omega \sin \theta \hat{\theta}
\end{equation}

Also because we will put the origin at the surface of the earth and not the center we shall use the effective gravitational acceleration that already includes the effects of the centrifugal force. Thus our expression for equation of motion will have no expression involving position and we will not have to bite our lips for this time-saving step.

Finally we come to the general equation of motion:
\begin{equation}
 F_{rot} =  mg_{eff}+ F_{cor}\\
 F_{rot} = mg_{eff} + v_{rot} \times 2\Omega
\end{equation}

In the above expressions the subscript "rot" stands for rotational reference frame. Next we take the determinant to get the expression below:
\begin{equation}
v \times \mathbf{\Omega} = \hat{r} \left(v_{\phi}\Omega \sin \theta \right) + \hat{\theta}\left(v_{\phi}\Omega \cos \theta \right) + \hat{\phi}\left( v_{r}\Omega \sin \theta - v_{\theta}\Omega \cos \theta \right)
\end{equation}

We know have our equations of motion:
\begin{equation}
\ddot{r} = v_{\phi} \Omega \sin \theta - g_{eff} \\
\ddot{\theta} = v_{\phi}\Omega \cos \theta \\
\ddot{\phi} = \left( v_{r} \Omega \sin \theta - v_{\phi} \Omega \cos \theta \right)
\end{equation}

The procedure of getting the solution is by using a bit of perturbation theory. First solve the equations by assuming \(\Omega\) is zero to get the following equations:
\begin{equation}
\theta = 0\\
\phi = 0 \\
 r = h - \frac{g t^2}{2}\\
v_{\theta} = 0\\
v_{\phi} = 0\\
v_{r} = -gt
\end{equation}

Next assume that the angular velocity is non-zero but very small so that our zero order approximation are fairly accurate enough and plug the results in the equations of motions. The only equation that survives is \(\ddot{\phi} = -gt \Omega \sin \theta\). This gives the deviation in the \hat{\phi} direction. To get how this deviation varies according to latitudes simply put in the transformation mentioned at the beginning of this section.

The next example introduces the phenomenon of Larmor precession. Imagine a negative electric charge q is moving in an elliptical orbit about a positive charge q. What are the equations of motion in the rotating reference frame?
\begin{equation}
F_{rot} = 2m\left(v_{rot} \times \Omega \right)  + q\left(v_{rot} + \Omega \times r \right) \times \vec{B} + \left(\Omega \times r \right) \times \Omega   + \frac{q^2}{4\pi \epsilon r^2} \hat{r} \\
F_{rot} = 2m\left(v_{rot} \times \Omega \right)  + qv_{rot} \times \vec{B} + q\left(\Omega \times r \right) \times \vec{B}+ \left(\Omega \times r \right) \times \Omega   + \frac{q^2}{4\pi \epsilon r^2}  \hat{r}\\
F_{rot} = v_{rot} \times \left(2m \Omega - q\vec{B}\right) - q\left(\Omega \times r \right) \times \vec{B}+\left(\Omega \times r \right) \times \Omega   + \frac{q^2}{4\pi \epsilon r^2} \hat{r}
\end{equation}
Now in order not to have linear velocity in the rotating reference frame the first term in the last equation has to disappear and therefore \Omega has to be \(\frac{qB}{2m}\). This is the Larmor frequency at which the precession occurs. The particle has no translational velocity. We are then left with the following expression:
\begin{equation}
F_{rot} = -\frac{q^2}{2m}\left(\vec{B} \times r \right) \times \vec{B} + \frac{q^2}{2m}\left(\vec{B}\times r \right) \times \vec{B} - \frac{q^2}{4\pi \epsilon r^2} \hat{r}
\end{equation}

Non-inertial frames: Matrix formulation

In the Lagrangian formulation  we came across this equation \( \left(\frac{de^{'}_i}{dt}\right)_i = \left(\frac{d\Omega}{dt}\right)_{ij}e^{'}_j  \) and it was mentioned that this matrix anti-symmetric. First we may give a heuristic argument as to why we may expect it to be. Notice that it arrives we take the derivative of a unit vector. Now, in what way can a unit vector change? It certainly can't change its length therefore the matrix in the equation can't have real eigenvalues. Remember a diagonalizable matrix acting upon a vector returns that same vector multiplied by a constant. As we have mentioned before this constant cannot be real for our matrix, it certainly can not be zero  so the only other option is that this eigenvalue is imaginary. But we also know that anti-symmetric matrices have imaginary eigenvalues hence we may expect our matrix to be anti-symmetric. This argument is not a proof because I am sure there are normal matrices with imaginary eigenvalues. But since we are in SO(3) we have some hope. Here then is a more rigorous argument to show that it is anti-symmetric. The argument goes as follows:
\begin{eqnarray}
 0 =& \frac{d}{dt}\left(e_i \cdot e_j \right)\\
 0 = &\frac{de_i}{dt}e_j + \frac{de_j}{dt}e_i\\
 0 =&\left(\frac{d\Omega}{dt}\right)_{ij}e_j e_j + \left(\frac{d\Omega}{dt}\right)_{ji}e_i e_i\\
 0 =  &\left(\frac{d\Omega}{dt}\right)_{ij} + \left(\frac{d\Omega}{dt}\right)_{ji}
\end{eqnarray}
 Looking at the last step we see that taking the transpose of the matrix and adding it to the matrix gives us the zero matrix. Hence we conclude it is an anti-symmetric matrix. QED
Since we start out by knowing that the unit vector cannot change its length in time we know that this change is actually related to a rotation. But at this point you should be disturbed, the equation we started with in the first paragraph has one specific direction namely, \( e_j\). This is where knowing a lot about matrices comes in. It turns out that all cross products can be represented as anti-symmetric matrices. The  equation we started off with is actually a  cross product! Re-writing \( \left(\frac{d\Omega}{dt}\right)\) as \(\omega \times \) where \(\omega = \frac{d\Omega}{dt}\) we have our connection to cross products.
Now that we know we have anti-symmetric matrix on our hand, let's reformulate newton's second law in non-inertial frames using matrices.
Lets suppose we have a vector  r, then we know that there is a standard matrix such that the vector's representation in another reference frame is given by Mr' where M is the standard matrix and r' is the vector in terms of another set of basis vectors. Since we are looking at orthogonal matrices then \(M^{-1}= M^{T}\) . That means that \(r' = M^{T}r\).Now look at the first derivative of r.
\begin{eqnarray}
 \dot{r} =& \dot{M}r' + M\dot{r^{'}}\\
  =&\dot{M}M^{T}r + M\dot{r^{'}}
\end{eqnarray}

Before we move any further lets take note that \(MM^{T} = I\)  and let's take the derivative of this expression to get \(\dot{M}M^{T} + M\dot{M^{T}} = 0\). It is important to notice that the matrix \(\dot{M}M^{T}\) is anti-symmetric based on the preceding equation. Taking its transpose and adding the two together gives us 0. Marvels of all marvels this is the same matrix we discussed earlier under a different guise. If you still have not witnessed the miracle, hold my statement with your grains of salt and read on. Let's define \(\dot{M}M^{T}\) to be the operator O. this will make our work later on easier.
Let's take the second derivative:
\begin{eqnarray}
\ddot{r} =& \dot{O}r + O\dot{r} + \dot{M}\dot{r}' + M \ddot{r}'\\
=& \dot{O}r + O\dot{r} - O^{T}\dot{r}+ O^{T}\dot{M}r' + M\ddot{r}'\\
=& \dot{O}r + 2 O\dot{r} + O^{T}\dot{M}r' + O\ddot{r}'\\
=& \dot{O}r + 2 O\dot{r} + O^{T}\dot{M}M^{T}r + M\ddot{r}'\\
=& \dot{O}r + 2 O\dot{r} - OOr + M\ddot{r}'\\
\end{eqnarray}
For the second step I went back to the first derivative of r to substitute what \(\dot{r}'\) was.

Finally, let's take a closer look at the last step we get when we calculate the second derivative. My claim was that O which is \(\dot{M}M^{T}\) is really \(\omega \times \) which is also the matrix \(\left(\frac{d\Omega}{dt}\right)_{ij}\). Let's substitute into the last step for the second derivative the expression \(\omega \times \) for O. We get the following:
\begin{equation}
\ddot{r} = \dot{\omega} \times r + 2 \omega \times \dot{r} - \omega \times \omega \times r + M\ddot{r}'
\end{equation}

What does all this mean? Well, mostly that matrices are unbelievably  awesome !