This post is the second part of a series about creating a matrix class in C#. In this post I will add matrix addition, subtraction, multiplication, and scalar multiplication. Matrix creation, cloning, transposing, creating an identity matrix, and converting a matrix to and from an array were previously covered in the first part of this series. In a future post I plan to add matrix equation solution through LU decomposition with pivoting (LUP), and matrix inversion.
I want to be able to take two matrices and add them. For those unfamiliar with this, there are reviews on Wikipedia Matrix Addition and Wolfram MathWorld Matrix Addition . For this I added the code below.
public static MatrixNumeric Add(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumColumns);
System.Diagnostics.Debug.Assert(mLeft.NumRows == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] + mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator +(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Add(mLeft, mRight);
}
The add method takes two matrices, mLeft and mRight and adds them together. It first checks to insure that the number of rows and columns in the matrix are the same for each matrix, otherwise matrix addition doesn’t really make sense. It then creates a new matrix with the same number of rows and columns as the matrices passed in. It then iterates over the rows and columns of the matrices, adds the elements from each matrix, and puts the results into the new matrix. The new matrix is then returned. I also added the addition operator overload to make matrix addition easier in code. This overload simply calls the Add method and passes in the two matrices.
Another thing that I want to be able to do is to take two matrices and subtract them. This is rather similar to matrix addition, except the values are subtracted instead of added. For this I added the code below.
public static MatrixNumeric Subtract(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumColumns);
System.Diagnostics.Debug.Assert(mLeft.NumRows == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] - mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator -(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Subtract(mLeft, mRight);
}
This code is quite similar to that for matrix addition, except that instead of adding the element values together, they are subtracted. As with matrix addition I added the subtraction operator overload to make matrix subtraction easier. The overload simple calls the Subtract method and passes in the two matrices.
The next thing that I want to be able to do is to take two matrices and multiply them. More information about matrix multiplication can be found on Wikipedia Matrix Multiplication and Wolfram MathWorld Matrix Multiplication For this I added the code below.
public static MatrixNumeric Multiply(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mRight.NumColumns; i++)
{
for (int j = 0; j < mLeft.NumRows; j++)
{
double dValue = 0.0;
for (int k = 0; k < mRight.NumRows; k++)
{
dValue += mLeft[j, k] * mRight[k, i];
}
mRet[j, i] = dValue;
}
}
return mRet;
}
public static MatrixNumeric operator *(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Multiply(mLeft, mRight);
}
The static Multiply function accepts two matrices and checks to ensure that the number of columns in the left matrix are equal to the number of columns in the right matrix. It then creates a new matrix with the number of rows equal to the number of rows in the left matrix and the number of columns equal to the number of columns in the right matrix. The multiplication starts with the first column in the left matrix and the first row in the right matrix. It creates and zeroes an accumulator dValue. It then multiplies the top most element in the first column of the right matrix by the left most element in the first row in the left matrix and adds the result to dValue. It then continues moving down the first column in the right matrix and across the first row in the left matrix multiplying the elements and adding it to dValue. After it reaches the end of the column and row, it stores the value of dValue into the newly created matrix. It then proceeds to repeat the process with the next row of the left matrix and the first column of the right matrix. After it passes the last row in the left matrix, it then moves to the next column in the right matrix and repeats. This process is continued until it reaches the last row in the left matrix and the last column in the right matrix, at which point the matrix multiplication is complete. The function then returns the new matrix that contains the multiplication result. As with matrix addition and subtraction, I added the multiplication operator overload to make matrix multiplication easier. The overload simple calls the Multiply method and passes in the two matrices.
The next thing that I want to be able to do is to multiply a matrix by a scalar. For this I added the code below.
public static MatrixNumeric Multiply(double dLeft, MatrixNumeric mRight)
{
MatrixNumeric mRet = new MatrixNumeric(mRight.NumRows, mRight.NumColumns);
for (int i = 0; i < mRight.NumRows; i++)
{
for (int j = 0; j < mRight.NumColumns; j++)
{
mRet[i, j] = dLeft * mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator *(double dLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Multiply(dLeft, mRight);
}
This is a fairly simple function. It takes a double which represents the scalar and a matrix which is multiplied by the scalar. It creates a new matrix which has the same number of rows and columns as the matrix passed in. It then iterates over the rows and columns and sets the element value of the new matrix to the element value of the original matrix times the scalar. As with the above cases, I added the multiplication operator overload to make this matrix multiplication easier. The overload simply calls the Multiply method and passes in the scalar and matrix.
The next thing that I want to be able to do is to perform the above operation when the matrix and scalar orders are reveresed. For this I added the code below.
public static MatrixNumeric Multiply(MatrixNumeric mLeft, double dRight)
{
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mLeft.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] * dRight;
}
}
return mRet;
}
public static MatrixNumeric operator *(MatrixNumeric mLeft, double dRight)
{
return MatrixNumeric.Multiply(mLeft, dRight);
}
This code operates very similar to the above section, with the only difference being that in this case the matrix is passed in as the first parameter and the scalar as the second. I also created an overloaded operator to make it easier to multiply a matrix by a scalar.
The next thing that I want to be able to do is to take a matrix and divide it by a scalar. For this I added the code below.
public static MatrixNumeric Divide(MatrixNumeric mLeft, double dRight)
{
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mLeft.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] / dRight;
}
}
return mRet;
}
public static MatrixNumeric operator /(MatrixNumeric mLeft, double dRight)
{
return MatrixNumeric.Divide(mLeft, dRight);
}
This code works quite similar to the above multiplication of a matrix with a scalar, except this time the matrix elements are divided by the scalar instead of multiplied by them. As with the above cases, I also overrode the division operator to make division of a matrix by a scalar easier to do in code.
That wraps up this post in the series. In a future post I plan to add matrix equation solution using LU with pivoting decomposition and matrix inversion. The full code of this article plus the previous one in this series is shown below.
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
namespace NumericalMethods
{
public class MatrixNumeric
{
private int m_nNumRows = 3;
private int m_nNumColumns = 3;
private double[,] m_dValues;
public MatrixNumeric()
{
m_dValues = new double[m_nNumRows, m_nNumColumns];
}
public MatrixNumeric(int nNumRows, int nNumColumns)
{
m_nNumRows = nNumRows;
m_nNumColumns = nNumColumns;
m_dValues = new double[m_nNumRows, m_nNumColumns];
}
public double this[int nRow, int nColumn]
{
get { return m_dValues[nRow, nColumn]; }
set { m_dValues[nRow, nColumn] = value; }
}
public int NumRows
{
get { return m_nNumRows; }
}
public int NumColumns
{
get { return m_nNumColumns; }
}
public MatrixNumeric Clone()
{
MatrixNumeric mRet = new MatrixNumeric(m_nNumRows, m_nNumColumns);
for (int i = 0; i < m_nNumRows; i++)
{
for (int j = 0; j < m_nNumColumns; j++)
{
mRet[i, j] = this[i, j];
}
}
return mRet;
}
public static MatrixNumeric Identity(int nSize)
{
MatrixNumeric mRet = new MatrixNumeric(nSize, nSize);
for (int i = 0; i < nSize; i++)
{
for (int j = 0; j < nSize; j++)
{
mRet[i, j] = (i == j) ? 1.0 : 0.0;
}
}
return mRet;
}
public MatrixNumeric Transpose()
{
MatrixNumeric mRet = new MatrixNumeric(m_nNumColumns, m_nNumRows);
for (int i = 0; i < m_nNumRows; i++)
{
for (int j = 0; j < m_nNumColumns; j++)
{
mRet[j, i] = this[i, j];
}
}
return mRet;
}
public static MatrixNumeric FromArray(double[] dLeft)
{
int nLength = dLeft.Length;
MatrixNumeric mRet = new MatrixNumeric(nLength, 1);
for (int i = 0; i < nLength; i++)
{
mRet[i, 0] = dLeft[i];
}
return mRet;
}
public static implicit operator MatrixNumeric(double[] dLeft)
{
return FromArray(dLeft);
}
public static double[] ToArray(MatrixNumeric mLeft)
{
System.Diagnostics.Debug.Assert((mLeft.NumColumns == 1 && mLeft.NumRows >= 1) || (mLeft.NumRows == 1 && mLeft.NumColumns >= 1));
double[] dRet = null;
if (mLeft.NumColumns > 1)
{
int nNumElements = mLeft.NumColumns;
dRet = new double[nNumElements];
for (int i = 0; i < nNumElements; i++)
{
dRet[i] = mLeft[0, i];
}
}
else
{
int nNumElements = mLeft.NumRows;
dRet = new double[nNumElements];
for (int i = 0; i < nNumElements; i++)
{
dRet[i] = mLeft[i, 0];
}
}
return dRet;
}
public static implicit operator double[](MatrixNumeric mLeft)
{
return ToArray(mLeft);
}
public static MatrixNumeric FromDoubleArray(double[,] dLeft)
{
int nLength0 = dLeft.GetLength(0);
int nLength1 = dLeft.GetLength(1);
MatrixNumeric mRet = new MatrixNumeric(nLength0, nLength1);
for (int i = 0; i < nLength0; i++)
{
for (int j = 0; j < nLength1; j++)
{
mRet[i, j] = dLeft[i, j];
}
}
return mRet;
}
public static implicit operator MatrixNumeric(double[,] dLeft)
{
return FromDoubleArray(dLeft);
}
public static double[,] ToDoubleArray(MatrixNumeric mLeft)
{
double[,] dRet = new double[mLeft.NumRows, mLeft.NumColumns];
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
dRet[i, j] = mLeft[i, j];
}
}
return dRet;
}
public static implicit operator double[,](MatrixNumeric mLeft)
{
return ToDoubleArray(mLeft);
}
public static MatrixNumeric Add(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumColumns);
System.Diagnostics.Debug.Assert(mLeft.NumRows == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] + mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator +(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Add(mLeft, mRight);
}
public static MatrixNumeric Subtract(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumColumns);
System.Diagnostics.Debug.Assert(mLeft.NumRows == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] - mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator -(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Subtract(mLeft, mRight);
}
public static MatrixNumeric Multiply(MatrixNumeric mLeft, MatrixNumeric mRight)
{
System.Diagnostics.Debug.Assert(mLeft.NumColumns == mRight.NumRows);
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mRight.NumColumns);
for (int i = 0; i < mRight.NumColumns; i++)
{
for (int j = 0; j < mLeft.NumRows; j++)
{
double dValue = 0.0;
for (int k = 0; k < mRight.NumRows; k++)
{
dValue += mLeft[j, k] * mRight[k, i];
}
mRet[j, i] = dValue;
}
}
return mRet;
}
public static MatrixNumeric operator *(MatrixNumeric mLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Multiply(mLeft, mRight);
}
public static MatrixNumeric Multiply(double dLeft, MatrixNumeric mRight)
{
MatrixNumeric mRet = new MatrixNumeric(mRight.NumRows, mRight.NumColumns);
for (int i = 0; i < mRight.NumRows; i++)
{
for (int j = 0; j < mRight.NumColumns; j++)
{
mRet[i, j] = dLeft * mRight[i, j];
}
}
return mRet;
}
public static MatrixNumeric operator *(double dLeft, MatrixNumeric mRight)
{
return MatrixNumeric.Multiply(dLeft, mRight);
}
public static MatrixNumeric Multiply(MatrixNumeric mLeft, double dRight)
{
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mLeft.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] * dRight;
}
}
return mRet;
}
public static MatrixNumeric operator *(MatrixNumeric mLeft, double dRight)
{
return MatrixNumeric.Multiply(mLeft, dRight);
}
public static MatrixNumeric Divide(MatrixNumeric mLeft, double dRight)
{
MatrixNumeric mRet = new MatrixNumeric(mLeft.NumRows, mLeft.NumColumns);
for (int i = 0; i < mLeft.NumRows; i++)
{
for (int j = 0; j < mLeft.NumColumns; j++)
{
mRet[i, j] = mLeft[i, j] / dRight;
}
}
return mRet;
}
public static MatrixNumeric operator /(MatrixNumeric mLeft, double dRight)
{
return MatrixNumeric.Divide(mLeft, dRight);
}
}
}