Basically
To get an understanding about rotation Matrixes we start in the 2 dimensional room:
To rotate a vector p1(x,y) in the 2 dimensional room by the angle φ, its x and y coordinate have to be moved from p1 to p2.
That means:
With the addition theorem cos(α + β) = cos(α) * cos(β) - sin(α) * sin(β) And the fact that
we can say
and that can be reduced to
The same can be done for y2:
With the addition theorem sin(α + β) = sin(α) * cos(β) + cos(α) * sin(β)
We finally get
These two formulas we can as well write in Matrix form
and get the same result
The Matrix
Is the rotation Matrix in the 2 dimensional room
3 dimensional room
Now let’s move forward to the 3 dimensional room. Here a rotation Matrix looks like
If we multiply a vector (x1, y1, z1) by this Matrix, we get a vector
The y coordinate is not touched and the x and z coordinate look the same as the x and y coordinate in the 2 dimensional room further above. That means this rotation Matrix rotates the (x1, y1, z1) by the angle φ in the x-z plane.
To rotate in all 3 dimensions, 3 different rotation Matrixes would be required.
For the x-y plane
For the x-z plane
For the y-z plane.
These 3 Matrixes could be put into one single 3 dimensional rotation Matrix. But that wouldn’t be too clear any more.
N dimensional room
For the y-z plane.
Now, in the N*N room the situation is basically the same. A rotation Matrix rotates a vector in one plane. Only the number of dimensions increases.
A N*N rotation Matrix U(p,q, φ) is a Matrix with the elements:
uij = 1 for I <> p, q
upp = uqq = cos(φ)
upq = sin(φ)
uqp = -sin(φ)
uij = 0 for all other elements In a 5x5 Matrix that looks like:
Where p and q are the 2 dimensions of the plane the rotation takes place.
Givens Transformation
To transform a Matrix into certain form it usually is advantageous to use a simple transformation Matrix. Such a simple Matrix is the rotation Matrix U(p,q, φ) as mentioned above.
This rotation Matrix can be used to eliminate elements of a Matrix as is done in the Gaussian elimination algorithm. If we regard the elements a11 as y and a12 as x of a Matrix A as a vector and want to eliminate a12, we might have to rotate this vector in a way the a12 coordinate becomes 0. This can be achieved if the new angle of this vector becomes 90° or 270°
If we rotate to 90° and let a11 and a12 be > 0, we get
Multiplying the Matrix A by a rotation Matrix with this angle eliminates a12 and if we perform the same calculation with.
a13, a14,… a23, a24,…, a34,.. we finally get a right triangle Matrix. And this is the Givens transformation with the condition that a12, a13, a14,… a23, a24,…, a34,..are all >0.
Now in the real Givens transformation the formulation is a bit different. For a11 and a12 it looks like:
The difference between these formulations and the twos further above is that here the y value of the rotated vector becomes negative if a11 is negative and with the twos further above it will always become positive. The x value always becomes 0. If we want to use the rotation Matrix just to eliminate elements in a Matrix equation, we can leave that and use the formulas further above. That does not make a big difference. Only the sign can switch in certain lines of the Matrix equation.
Solving a Matrix equation by a Givens transformations
To solve a Matrix equation by applying Givens transformations is basically the same as using the Gaussian elimination algorithm. We have to eliminate the elements below the main diagonal of the Matrix first to be able to solve the equations afterwards. For this we have to build a rotation Matrix U(p,q, φ) for each element we want to eliminate and multiply first the Matrix and the the solution vector by U(p,q, φ):
and build the rotation matrix as
That looks like
public void CalcRMatrix(int p, int q)
{
int i, j;
double radius;
radius = Math.Sqrt(m[p, p] * m[p, p] + m[q, p] * m[q, p]);
for (i = 0 ; i < maxOrder; i++)
{
for (j = 0; j < maxOrder; j++)
{
r[i,j] = 0.0;
}
r[i,i] = 1.0;
}
if (radius != 0.0)
{
if (m[p, p] == 0.0)
{
r[p, p] = 0.0;
r[p, q] = 1.0;
}
else
{
r[p, p] = m[p, p] / radius;
r[p, q] = m[q, p] / radius;
}
r[q,q] = r[p,p];
r[q,p] = -r[p,q];
}
}
with
double[,] m = double[MaxOrder, MaxOrder]; // equation Matrix
double[] r = new double[MaxOrder]; // rotation Matrix
Now p basically represents the main diagonal of the Matrix and q the elements below. With one index running down the main diagonal and a second running down from the firs index to the bottom of the Matrix we can build the right triangle Matrix
public void RunGivens(int maxOrder)
{
int i, k;
for (i = 0; i < maxOrder; i++)
{
for (k = i + 1; k < maxOrder; k++)
{
if (m[k, i] != 0.0)
{
CalcRMatrix(i, k);
MultMatrix(r, i, k);
}
}
}
}
Now for the multiplication we have the special case that only the elements of the p and q columns are affected by the rotation. Therefore we can just keep all the other elements and save a lot of calculation effort. For this multiplication I add the function
MultMatrix(double[,] matrix2, int
p, int q)
to the TMatrix class. It multiplies the own Matrix equation by matrix2 considering p and q. I made a small speed optimisation in this multiplication: As there are only the elements of row p and q affected by this multiplication, I only manipulate these two rows. All the other rows can be left as they are. That reduces the calculation effort by the factor n – 2 and means that this multiplication may not be used for a standard matrix multiplication. The function uses a temporary Matrix al and the temporary vector yl to buffer the data. The result remains in the own matrix as I want to go on working with this matrix.
public void MultMatrix(double[,] matrix, int p, int q)
{
int i, j, k;
double[,] al = new double[maxOrder, maxOrder];
double[] yl = new double[maxOrder];
for (i = 0; i < maxOrder; i++)
{
if ((i == p) || (i == q)) // only calculate row p and q
{
yl[i] = 0.0;
for (j = 0; j < maxOrder; j++)
{
yl[i] = yl[i] + (y[j] * matrix[i, j]);
al[i, j] = 0.0;
for (k = 0; k < maxOrder; k++)
{
al[i, j] = al[i, j] + (m[k, j] * matrix[i, k]);
}
}
}
else
{
yl[i] = y[i];
for (j = 0; j < maxOrder; j++)
{
al[i, j] = m[i, j];
}
}
}
for (i = 0; i < maxOrder; i++)
{
y[i] = yl[i];
for (j = 0; j < maxOrder; j++)
{
m[i, j] = al[i, j];
}
}
}
The function copies the multiplied Matrix into its own Matrix c and the solution vector into d. From here the Matrix equation can be solved now by a small function
public double[] SolveRMatrix()
{
int k, l;
double[] x = new double[maxOrder];
for (k = maxOrder - 1; k >= 0; k--)
{
for (l = maxOrder - 1; l > k; l--)
{
y[k] = y[k] - x[l] * m[k, l];
}
if (m[k, k] != 0)
x[k] = y[k] / m[k, k];
else
x[k] = 0;
}
return x;
}
Solving a Matrix equation like this yields a bit more accuracy and is a bit faster than using the Gaussian elimination algorithm. The higher accuracy is the result of the smaller values in the main Matrix after building the right triangle Matrix.
A sample Matrix equation like
Has, solved by Givens transformations, the result
The same Matrix equation solved by Gauss results in
o.k. it’s not that much. But it’s just a 4*4 Matrix
I implemented an online Matrix equation solver with the Givens rotation matrix for elimination. This Matrix equation solver can be found here