G3D::Matrix::Impl Class Reference

#include <Matrix.h>

Public Member Functions
	~Impl ()
Public Member Functions inherited from G3D::ReferenceCountedObject
virtual	~ReferenceCountedObject ()

Static Public Member Functions
static void *	operator new (size_t size)
static void	operator delete (void *p)

Private Member Functions
void	setSize (int newRows, int newCols)
	Impl ()
	Impl (const Matrix3 &M)
	Impl (const Matrix4 &M)
	Impl (int r, int c)
Impl &	operator= (const Impl &m)
	Impl (const Impl &B)
void	setZero ()
void	set (int r, int c, T v)
const T &	get (int r, int c) const
void	mul (const Impl &B, Impl &out) const
void	add (const Impl &B, Impl &out) const
void	add (T B, Impl &out) const
void	sub (const Impl &B, Impl &out) const
void	sub (T B, Impl &out) const
void	lsub (T B, Impl &out) const
void	arrayMul (const Impl &B, Impl &out) const
void	mul (T B, Impl &out) const
void	arrayDiv (const Impl &B, Impl &out) const
void	div (T B, Impl &out) const
void	negate (Impl &out) const
void	inverseViaAdjoint (Impl &out) const
void	inverseInPlaceGaussJordan ()
void	adjoint (Impl &out) const
void	cofactor (Impl &out) const
T	cofactor (int r, int c) const
void	transpose (Impl &out) const
T	determinant () const
T	determinant (int r, int c) const
void	arrayLog (Impl &out) const
void	arrayExp (Impl &out) const
void	arraySqrt (Impl &out) const
void	arrayCos (Impl &out) const
void	arraySin (Impl &out) const
void	swapRows (int r0, int r1)
void	swapAndNegateCols (int c0, int c1)
void	mulRow (int r, const T &v)
void	abs (Impl &out) const
void	withoutRowAndCol (int excludeRow, int excludeCol, Impl &out) const
bool	anyNonZero () const
bool	allNonZero () const
void	setRow (int r, const T *vals)
void	setCol (int c, const T *vals)

Private Attributes
T **	elt
T *	data
int	R
int	C
int	dataSize

Friends
class	Matrix

Detailed Description

Used internally by Matrix.

Does not throw exceptions– assumes the caller has taken care of argument checking.

Constructor & Destructor Documentation

G3D::Matrix::Impl::~Impl

(

)

 {
 //delete[] elt;
 System::free(elt);
 System::alignedFree(data);
}

Here is the call graph for this function:

G3D::Matrix::Impl::Impl ( )

inlineprivate

: elt(NULL), data(NULL), R(0), C(0), dataSize(0) {}

G3D::Matrix::Impl::Impl ( const Matrix3 &  M )

private

 : elt(NULL), data(NULL), R(0), C(0), dataSize(0){
 setSize(3, 3);
 for (int r = 0; r < 3; ++r) {
 for (int c = 0; c < 3; ++c) {
 set(r, c, M[r][c]);
 }
 }

}

Here is the call graph for this function:

G3D::Matrix::Impl::Impl ( const Matrix4 &  M )

private

 : elt(NULL), data(NULL), R(0), C(0), dataSize(0) {
 setSize(4, 4);
 for (int r = 0; r < 4; ++r) {
 for (int c = 0; c < 4; ++c) {
 set(r, c, M[r][c]);
 }
 }
}

Here is the call graph for this function:

G3D::Matrix::Impl::Impl ( int  r, int  c  )

inlineprivate

 : elt(NULL), data(NULL), R(0), C(0), dataSize(0) {
 setSize(r, c);
 }

Here is the call graph for this function:

G3D::Matrix::Impl::Impl ( const Impl &  B )

inlineprivate

 : elt(NULL), data(NULL), R(0), C(0), dataSize(0) {
 // Use the assignment operator
 *this = B;
 }

Member Function Documentation

void G3D::Matrix::Impl::abs ( Impl &  out ) const

private

void G3D::Matrix::Impl::add ( const Impl &  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::add ( T  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::adjoint ( Impl &  out ) const

private

 { 
 cofactor(out);
 // Transpose is safe to perform in place
 out.transpose(out);
}

Here is the call graph for this function:

bool G3D::Matrix::Impl::allNonZero ( ) const

private

 {
 int N = R * C;
 for (int i = 0; i < N; ++i) {
 if (data[i] == 0.0) {
 return false;
 }
 }
 return true;
}

bool G3D::Matrix::Impl::anyNonZero ( ) const

private

 {
 int N = R * C;
 for (int i = 0; i < N; ++i) {
 if (data[i] != 0.0) {
 return true;
 }
 }
 return false;
}

void G3D::Matrix::Impl::arrayCos ( Impl &  out ) const

private

void G3D::Matrix::Impl::arrayDiv ( const Impl &  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::arrayExp ( Impl &  out ) const

private

void G3D::Matrix::Impl::arrayLog ( Impl &  out ) const

private

void G3D::Matrix::Impl::arrayMul ( const Impl &  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::arraySin ( Impl &  out ) const

private

void G3D::Matrix::Impl::arraySqrt ( Impl &  out ) const

private

void G3D::Matrix::Impl::cofactor ( Impl &  out ) const

private

Matrix of all cofactors

 {
 debugAssert(&out != this);
 for(int r = 0; r < R; ++r) {
 for(int c = 0; c < C; ++c) {
 out.set(r, c, cofactor(r, c));
 }
 } 
}

Here is the call graph for this function:

Matrix::T G3D::Matrix::Impl::cofactor ( int  r, int  c  ) const

private

Cofactor [r][c] is defined as C[r][c] = -1 ^(r+c) * det(A[r][c]), where A[r][c] is the (R-1)x(C-1) matrix formed by removing row r and column c from the original matrix.

 {
 // Strang p. 217
 float s = isEven(r + c) ? 1.0f : -1.0f;

 return s * determinant(r, c);
}

Here is the call graph for this function:

Matrix::T G3D::Matrix::Impl::determinant ( ) const

private

 {

 debugAssert(R == C);

 // Compute using cofactors
 switch(R) {
 case 0:
 return 0;

 case 1:
 // Determinant of a 1x1 is the element
 return elt[0][0];

 case 2:
 // Determinant of a 2x2 is ad-bc
 return elt[0][0] * elt[1][1] - elt[0][1] * elt[1][0];

 case 3:
 {
 // Determinant of an nxn matrix is the dot product of the first
 // row with the first row of cofactors. The base cases of this
 // method get called a lot, so we spell out the implementation
 // for the 3x3 case.

 float cofactor00 = elt[1][1] * elt[2][2] - elt[1][2] * elt[2][1];
 float cofactor10 = elt[1][2] * elt[2][0] - elt[1][0] * elt[2][2];
 float cofactor20 = elt[1][0] * elt[2][1] - elt[1][1] * elt[2][0];
 
 return Matrix::T(
 elt[0][0] * cofactor00 +
 elt[0][1] * cofactor10 +
 elt[0][2] * cofactor20);
 }
 
 default:
 {
 // Determinant of an n x n matrix is the dot product of the first
 // row with the first row of cofactors
 T det = 0;

 for (int c = 0; c < C; ++c) {
 det += elt[0][c] * cofactor(0, c);
 }

 return det;
 }
 }
}

Here is the caller graph for this function:

Matrix::T G3D::Matrix::Impl::determinant ( int  r, int  c  ) const

private

Determinant computed without the given row and column

 {
 debugAssert(R > 0);
 debugAssert(C > 0);
 Impl A(R - 1, C - 1);
 withoutRowAndCol(nr, nc, A);
 return A.determinant();
}

Here is the call graph for this function:

void G3D::Matrix::Impl::div ( T  B, Impl &  out  ) const

private

Ok if out == this or out == B

Here is the caller graph for this function:

const T& G3D::Matrix::Impl::get ( int  r, int  c  ) const

inlineprivate

 {
 debugAssert(r < R);
 debugAssert(c < C);
 return elt[r][c];
 }

Here is the caller graph for this function:

void G3D::Matrix::Impl::inverseInPlaceGaussJordan ( )

private

Use Gaussian elimination with pivots to solve for the inverse destructively in place.

 {
 debugAssertM(R == C, 
 format(
 "Cannot perform Gauss-Jordan inverse on a non-square matrix."
 " (Argument was %dx%d)",
 R, C));

 // Exchange to float elements
# define SWAP(x, y) {float temp = x; x = y; y = temp;}

 // The integer arrays pivot, rowIndex, and colIndex are
 // used for bookkeeping on the pivoting
 static Array<int> colIndex, rowIndex, pivot;

 int col = 0, row = 0;

 colIndex.resize(R);
 rowIndex.resize(R);
 pivot.resize(R);

 static const int NO_PIVOT = -1;

 // Initialize the pivot array to default values.
 for (int i = 0; i < R; ++i) {
 pivot[i] = NO_PIVOT;
 }

 // This is the main loop over the columns to be reduced
 // Loop over the columns.
 for (int c = 0; c < R; ++c) {

 // Find the largest element and use that as a pivot
 float largestMagnitude = 0.0;

 // This is the outer loop of the search for a pivot element
 for (int r = 0; r < R; ++r) {

 // Unless we've already found the pivot, keep going
 if (pivot[r] != 0) {

 // Find the largest pivot
 for (int k = 0; k < R; ++k) {
 if (pivot[k] == NO_PIVOT) {
 const float mag = fabs(elt[r][k]);

 if (mag >= largestMagnitude) {
 largestMagnitude = mag;
 row = r; col = k;
 }
 }
 }
 }
 }

 pivot[col] += 1;

 // Interchange columns so that the pivot element is on the diagonal (we'll have to undo this
 // at the end)
 if (row != col) {
 for (int k = 0; k < R; ++k) {
 SWAP(elt[row][k], elt[col][k])
 }
 }

 // The pivot is now at [row, col]
 rowIndex[c] = row; 
 colIndex[c] = col;
 
 double piv = elt[col][col];

 debugAssertM(piv != 0.0, "Matrix is singular");

 // Divide everything by the pivot (avoid computing the division
 // multiple times).
 const double pivotInverse = 1.0 / piv;
 elt[col][col] = 1.0;

 for (int k = 0; k < R; ++k) {
 elt[col][k] *= Matrix::T(pivotInverse);
 }

 // Reduce all rows
 for (int r = 0; r < R; ++r) {
 // Skip over the pivot row
 if (r != col) {

 double oldValue = elt[r][col];
 elt[r][col] = 0.0;

 for (int k = 0; k < R; ++k) {
 elt[r][k] -= Matrix::T(elt[col][k] * oldValue);
 }
 }
 }
 }


 // Put the columns back in the correct locations
 for (int i = R - 1; i >= 0; --i) {
 if (rowIndex[i] != colIndex[i]) {
 for (int k = 0; k < R; ++k) {
 SWAP(elt[k][rowIndex[i]], elt[k][colIndex[i]]);
 }
 }
 } 
 
# undef SWAP
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G3D::Matrix::Impl::inverseViaAdjoint ( Impl &  out ) const

private

Slow way of computing an inverse; for reference

 {
 debugAssert(&out != this);

 // Inverse = adjoint / determinant

 adjoint(out);

 // Don't call the determinant method when we already have an
 // adjoint matrix; there's a faster way of computing it: the dot
 // product of the first row and the adjoint's first col.
 double det = 0.0;
 for (int r = R - 1; r >= 0; --r) {
 det += elt[0][r] * out.elt[r][0];
 }

 out.div(Matrix::T(det), out);
}

Here is the call graph for this function:

void G3D::Matrix::Impl::lsub ( Matrix::T  B, Impl &  out  ) const

private

B - this

 {
 const Impl& A = *this;

 debugAssert(A.C == out.C);
 debugAssert(A.R == out.R);

 for (int i = R * C - 1; i >= 0; --i) {
 out.data[i] = B - A.data[i];
 }
}

void G3D::Matrix::Impl::mul ( const Impl &  B, Impl &  out  ) const

private

Multiplies this by B and puts the result in out.

 {
 const Impl& A = *this;

 debugAssertM(
 (this != &out) && (&B != &out),
 "Output argument to mul cannot be the same as an input argument.");

 debugAssert(A.C == B.R);
 debugAssert(A.R == out.R);
 debugAssert(B.C == out.C);

 for (int r = 0; r < out.R; ++r) {
 for (int c = 0; c < out.C; ++c) {
 T sum = 0.0;
 for (int i = 0; i < A.C; ++i) {
 sum += A.get(r, i) * B.get(i, c);
 }
 out.set(r, c, sum);
 }
 }
}

Here is the call graph for this function:

void G3D::Matrix::Impl::mul ( T  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::mulRow ( int  r, const T &  v  )

private

 {
 T* row = elt[r];

 for (int c = 0; c < C; ++c) {
 row[c] *= v;
 }
}

Here is the call graph for this function:

void G3D::Matrix::Impl::negate ( Impl &  out ) const

private

static void G3D::Matrix::Impl::operator delete ( void *  p )

inlinestatic

 {
 System::free(p);
 }

Here is the call graph for this function:

static void* G3D::Matrix::Impl::operator new ( size_t  size )

inlinestatic

 {
 return System::malloc(size);
 }

Here is the call graph for this function:

Matrix::Impl & G3D::Matrix::Impl::operator= ( const Impl &  m )

private

 {
 setSize(m.R, m.C);
 System::memcpy(data, m.data, R * C * sizeof(T));
 ++Matrix::debugNumCopyOps;
 return *this;
}

Here is the call graph for this function:

void G3D::Matrix::Impl::set ( int  r, int  c, T  v  )

inlineprivate

 {
 debugAssert(r < R);
 debugAssert(c < C);
 elt[r][c] = v;
 }

Here is the caller graph for this function:

void G3D::Matrix::Impl::setCol ( int  c, const T *  vals  )

private

 {
 for (int r = 0; r < R; ++r) {
 elt[r][c] = vals[r];
 }
}

void G3D::Matrix::Impl::setRow ( int  r, const T *  vals  )

private

 {
 debugAssert(r >= 0);
 System::memcpy(elt[r], vals, sizeof(T) * C);
}

Here is the call graph for this function:

void G3D::Matrix::Impl::setSize ( int  newRows, int  newCols  )

private

If R*C is much larger or smaller than the current, deletes all previous data and resets to random data. Otherwise re-uses existing memory and just resets R, C, and the row pointers.

 {
 if ((R == newRows) && (C == newCols)) {
 // Nothing to do
 return;
 }

 int newSize = newRows * newCols;

 R = newRows; C = newCols;

 // Only allocate if we need more space
 // or the size difference is ridiculous
 if ((newSize > dataSize) || (newSize < dataSize / 4)) {
 System::alignedFree(data);
 data = (float*)System::alignedMalloc(R * C * sizeof(T), 16);
 ++Matrix::debugNumAllocOps;
 dataSize = newSize;
 }

 // Construct the row pointers
 //delete[] elt;
 System::free(elt);
 elt = (T**)System::malloc(R * sizeof(T*));// new T*[R];

 for (int r = 0; r < R; ++ r) {
 elt[r] = data + r * C;
 }
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G3D::Matrix::Impl::setZero ( )

private

 {
 System::memset(data, 0, R * C * sizeof(T));
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G3D::Matrix::Impl::sub ( const Impl &  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::sub ( T  B, Impl &  out  ) const

private

Ok if out == this or out == B

void G3D::Matrix::Impl::swapAndNegateCols ( int  c0, int  c1  )

private

 {
 for (int r = 0; r < R; ++r) {
 T* row = elt[r];

 const T temp = -row[c0];
 row[c0] = -row[c1];
 row[c1] = temp;
 }
}

Here is the call graph for this function:

void G3D::Matrix::Impl::swapRows ( int  r0, int  r1  )

private

 {
 T* R0 = elt[r0];
 T* R1 = elt[r1];

 for (int c = 0; c < C; ++c) {
 T temp = R0[c];
 R0[c] = R1[c];
 R1[c] = temp;
 }
}

void G3D::Matrix::Impl::transpose ( Impl &  out ) const

private

Ok if out == this or out == B

 {
 debugAssert(out.R == C);
 debugAssert(out.C == R);

 if (&out == this) {
 // Square matrix in place
 for (int r = 0; r < R; ++r) {
 for (int c = r + 1; c < C; ++c) {
 T temp = get(r, c);
 out.set(r, c, get(c, r));
 out.set(c, r, temp);
 }
 }
 } else {
 for (int r = 0; r < R; ++r) {
 for (int c = 0; c < C; ++c) {
 out.set(c, r, get(r, c)); 
 }
 }
 }
}

Here is the call graph for this function:

Here is the caller graph for this function:

void G3D::Matrix::Impl::withoutRowAndCol ( int  excludeRow, int  excludeCol, Impl &  out  ) const

private

Makes a (R-1)x(C-1) copy of this matrix

 {
 debugAssert(out.R == R - 1);
 debugAssert(out.C == C - 1);

 for (int r = 0; r < out.R; ++r) {
 for (int c = 0; c < out.C; ++c) {
 out.elt[r][c] = elt[r + ((r >= excludeRow) ? 1 : 0)][c + ((c >= excludeCol) ? 1 : 0)];
 }
 }
}

Friends And Related Function Documentation

friend class Matrix

friend

Member Data Documentation

int G3D::Matrix::Impl::C

private

The number of columns

T* G3D::Matrix::Impl::data

private

Row major data for the entire matrix.

int G3D::Matrix::Impl::dataSize

private

T** G3D::Matrix::Impl::elt

private

elt[r][c] = the element. Pointers into data.

int G3D::Matrix::Impl::R

private

The number of rows

---

The documentation for this class was generated from the following files:

• dep/g3dlite/include/G3D/Matrix.h
• dep/g3dlite/source/Matrix.cpp