Subscribe to GameDev.net Direct to receive the latest updates and exclusive content.
Guest, the last post of this topic is over 60 days old and at this point you may not reply in this topic. If you wish to continue this conversation start a new topic.
Posted 05 March 2004 - 05:38 AM
Posted 05 March 2004 - 09:43 AM
Posted 06 March 2004 - 01:12 AM
Posted 06 March 2004 - 01:53 AM
(1/det(A)) * TRANSFORM( A )
Posted 06 March 2004 - 06:44 PM
Posted 06 March 2004 - 06:54 PM
Posted 07 March 2004 - 12:01 AM
Posted 07 March 2004 - 12:02 AM
quote:
Original post by optics
darookie: The matrix is of abitrary size( i.e. 4x4 ) but the values are not.
Posted 07 March 2004 - 07:33 AM
quote:
Original post by optics
then solve for x_{1} to x_{k}. I am having a hard time trying to figure out what b_{1}, b_{2}, b_{3}, b_{4} is suppose to be set to. Do we pull these values out of thin air or is there something I overlooked?
| e11 e12 e13 e14 1 0 0 0 |
| e21 e22 e23 e24 0 1 0 0 |
| e31 e32 e33 e34 0 0 1 0 |
| e41 e42 e43 e44 0 0 0 1 |
| 1 0 0 0 e''11 e''12 e''13 e''14 |
| 0 1 0 0 e''21 e''22 e''23 e''24 |
| 0 0 1 0 e''31 e''32 e''33 e''34 |
| 0 0 0 1 e''41 e''42 e''43 e''44 |
Posted 07 March 2004 - 08:58 AM
quote:
Original post by optics
then solve for x<sub>1</sub> to x<sub>k</sub>. I am having a hard time trying to figure out what b<sub>1</sub>, b<sub>2</sub>, b<sub>3</sub>, b<sub>4</sub> is suppose to be set to. Do we pull these values out of thin air or is there something I overlooked?quote:
The problem is to solve the system Ax = b
If you calculate the inverse of A, then you can calculate inv(A)b to find x.
You can factor A into upper and lower triangular matrices, to get LUx = b. Then you solve the problem Lc = b, then solve Ux = c. Because the matrices are triangular, the systems can be easily solved with substitution.
You can also attach b as an extra column to A, which is called an augmented matrix, and then use elimination to turn A into an upper triangular matrix, applying the same transformations to b that you do to the rows of A. Then you can solve for x by using the triangular matrix, the new values for b, and substitution.
Posted 07 March 2004 - 11:20 AM
Posted 07 March 2004 - 11:27 AM
quote:
Original post by optics
1) Could you show me the mathematical solution to obtain I(A^{-1})?
quote:
2)If I post my code would you look over and let me know if I am going in the right direction or am just way off base?
Posted 07 March 2004 - 11:52 AM
| e11 e12 e13 e14 1 0 0 0 |
| e21 e22 e23 e24 0 1 0 0 |
| e31 e32 e33 e34 0 0 1 0 |
| e41 e42 e43 e44 0 0 0 1 |
| 1 0 0 0 e''11 e''12 e''13 e''14 |
| 0 1 0 0 e''21 e''22 e''23 e''24 |
| 0 0 1 0 e''31 e''32 e''33 e''34 |
| 0 0 0 1 e''41 e''42 e''43 e''44 |
class CMatrix4x4
{
public:
CMatrix4x4(void); // default constructor
CMatrix4x4(double c11, double c12, double c13, double c14,
double c21, double c22, double c23, double c24,
double c31, double c32, double c33, double c34,
double c41, double c42, double c43, double c44); // constructor
~CMatrix4x4(); // destructor
double Determinant(); // determinant of the matrix
CMatrix4x4& Inverse(); // inverse of the matrix
private:
double e11, e12, e13, e14,
e21, e22, e23, e24,
e31, e32, e33, e34,
e41, e42, e43, e44;
};
double CMatrix4x4::Determinant()
{
return e11*e22*e33*e44 - e11*e22*e34*e43 + e11*e23*e34*e42 - e11*e23*e32*e44
+ e11*e24*e32*e43 - e11*e24*e33*e42 - e12*e23*e34*e41 + e12*e23*e31*e44
- e12*e24*e31*e43 + e12*e24*e33*e41 - e12*e21*e33*e44 + e12*e21*e34*e43
+ e13*e24*e31*e42 - e13*e24*e32*e41 + e13*e21*e32*e44 - e13*e21*e34*e42
+ e13*e22*e34*e41 - e13*e22*e31*e44 - e14*e21*e32*e43 + e14*e21*e33*e42
- e14*e22*e33*e41 + e14*e22*e31*e43 - e14*e23*e31*e42 + e14*e23*e32*e41;
}
CMatrix4x4& CMatrix4x4::Inverse()
{
double s1, s2, s3, s4;
double negX;
if(!this->Determinant()) // if determinant is 0.0 then matrix has no inverse
{ return *this; }
else // otherwise there is an inverse and we need to find it
{
/************************** Step 1 of Gaussian Elimination *******************************/
if( e21 == 1 ) // swap row 2 with row 1 so that e11 == 1
{
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e21; e12 = e22; e13 = e23; e14 = e24;
e21 = s1; e22 = s2; e23 = s3; e24 = s4;
}
else if( e31 == 1 ) // swap row 3 with row 1 so that e11 == 1
{
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e31; e12 = e32; e13 = e33; e14 = e34;
e31 = s1; e32 = s2; e33 = s3; e34 = s4;
}
else if( e41 == 1 ) // swap row 4 with row 1 so that e11 == 1
{ // do Gaussian step 1
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e41; e12 = e42; e13 = e43; e14 = e44;
e41 = s1; e42 = s2; e43 = s3; e44 = s4;
}
// if column 1 does not have a 1 in it then we will do math to make e11 == 1
// sort of hackish but we need to get the code working then we can fix it
//row 2 math
negX = -e21;
e21 = e11/e11*negX + e21;
e22 = e12/e11*negX + e22;
e23 = e13/e11*negX + e23;
e24 = e14/e11*negX + e24;
//row 3 math
negX = -e31;
e31 = e11/e11*negX + e31;
e32 = e12/e11*negX + e32;
e33 = e13/e11*negX + e33;
e34 = e14/e11*negX + e34;
//row 4 math
negX = -e41;
e41 = e11/e11*negX + e41;
e42 = e12/e11*negX + e42;
e43 = e13/e11*negX + e43;
e44 = e14/e11*negX + e44;
/********************** Step 2 of Gaussian Elimination *************************************/
//row 3 math
negX = -e32;
e32 = e22/e22*negX + e32;
e33 = e23/e22*negX + e33;
e34 = e24/e22*negX + e34;
//row 4 math
negX = -e42;
e42 = e22/e22*negX + e42;
e43 = e23/e22*negX + e43;
e44 = e24/e22*negX + e44;
/***************************** Step 3 of Gaussian Elimination *******************************/
//row 4 math
negX = -e43;
e43 = e33/e33*negX + e43;
e44 = e34/e33*negX + e44;
}
return *this;
}
Posted 07 March 2004 - 12:18 PM
to here:
| e11 e12 e13 e14 1 0 0 0 |
| e21 e22 e23 e24 0 1 0 0 |
| e31 e32 e33 e34 0 0 1 0 |
| e41 e42 e43 e44 0 0 0 1 |
| 1 0 0 0 e'11 e'12 e'13 e'14 |
| 0 1 0 0 e'21 e'22 e'23 e'24 |
| 0 0 1 0 e'31 e'32 e'33 e'34 |
| 0 0 0 1 e'41 e'42 e'43 e'44 |
| 1 e12/e11 e13/e11 e14/e11 1/e11 0 0 0 |
| e21 e22 e23 e24 0 1 0 0 |
| e31 e32 e33 e34 0 0 1 0 |
| e41 e42 e43 e44 0 0 0 1 |
| 1 e12/e11 e13/e11 e14/e11 1/e11 0 0 0 |
| 0 e*22 e*23 e*24 -e21/e11 1 0 0 |
| 0 e*32 e*33 e*34 -e31/e11 0 1 0 |
| 0 e*42 e*43 e*44 -e31/e11 0 0 1 |
Posted 08 March 2004 - 07:39 AM
class CMatrix4x4
{
public:
CMatrix4x4(void); // default constructor
CMatrix4x4(double c11, double c12, double c13, double c14,
double c21, double c22, double c23, double c24,
double c31, double c32, double c33, double c34,
double c41, double c42, double c43, double c44); // constructor
~CMatrix4x4(); // destructor
Identity(); // identity matrix
double Determinant(); // determinant of the matrix
CMatrix4x4& Inverse(); // inverse of the matrix
private:
double e11, e12, e13, e14,
e21, e22, e23, e24,
e31, e32, e33, e34,
e41, e42, e43, e44;
};
CMatrix4x4::Identity()
{
e11 = e22 = e33 = e44 = 1.0f;
e12 = e13 = e14 = e21 = e23 = e24 =
e31 = e32 = e34 = e41 = e42 = e43 = 0.0f;
}
double CMatrix4x4::Determinant()
{
return e11*e22*e33*e44 - e11*e22*e34*e43 + e11*e23*e34*e42 - e11*e23*e32*e44
+ e11*e24*e32*e43 - e11*e24*e33*e42 - e12*e23*e34*e41 + e12*e23*e31*e44
- e12*e24*e31*e43 + e12*e24*e33*e41 - e12*e21*e33*e44 + e12*e21*e34*e43
+ e13*e24*e31*e42 - e13*e24*e32*e41 + e13*e21*e32*e44 - e13*e21*e34*e42
+ e13*e22*e34*e41 - e13*e22*e31*e44 - e14*e21*e32*e43 + e14*e21*e33*e42
- e14*e22*e33*e41 + e14*e22*e31*e43 - e14*e23*e31*e42 + e14*e23*e32*e41;
}
CMatrix4x4& CMatrix4x4::Inverse()
{
double s1, s2, s3, s4;
double negX, div;
CMatrix4x4 m4x4Inv;
m4x4Inv.Identity();
if(!this->Determinant()) // if determinant is 0.0 then matrix has no inverse
{ return *this; }
else // otherwise there is an inverse and we need to find it
{
/************************** Step 1 of Gaussian Elimination *******************************/
if( e21 == 1 ) // swap row 2 with row 1 so that e11 == 1
{
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e21; e12 = e22; e13 = e23; e14 = e24;
e21 = s1; e22 = s2; e23 = s3; e24 = s4;
}
else if( e31 == 1 ) // swap row 3 with row 1 so that e11 == 1
{
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e31; e12 = e32; e13 = e33; e14 = e34;
e31 = s1; e32 = s2; e33 = s3; e34 = s4;
}
else if( e41 == 1 ) // swap row 4 with row 1 so that e11 == 1
{
s1 = e11; s2 = e12; s3 = e13; s4 = e14;
e11 = e41; e12 = e42; e13 = e43; e14 = e44;
e41 = s1; e42 = s2; e43 = s3; e44 = s4;
}
//if column 1 does not have a 1 in it then we will do math to make e11 == 1
//sort of hackish but we need to get the code working then we can fix it
// row 1 math
div = e11;
e11 = e11/div;
e12 = e12/div;
e13 = e13/div;
e14 = e14/div;
m4x4Inv.e11 = 1/div;
// row 2 math
negX = -e21;
e21 = e11*negX + e21;
e22 = e12*negX + e22;
e23 = e13*negX + e23;
e24 = e14*negX + e24;
m4x4Inv.e21 = negX/div;
// row 3 math
negX = -e31;
e31 = negX + e31;
e32 = e12*negX + e32;
e33 = e13*negX + e33;
e34 = e14*negX + e34;
m4x4Inv.e31 = negX/div;
// row 4 math
negX = -e41;
e41 = negX + e41;
e42 = e12*negX + e42;
e43 = e13*negX + e43;
e44 = e14*negX + e44;
m4x4Inv.e41 = negX/div;
/********************* Step 2 of Gaussian Elimination *************************************/
div = e22;
// row 2 math
e22 = e22/div;
e23 = e23/div;
e24 = e24/div;
m4x4Inv.e22 = 1/div;
// row 1 math
negX = -e12;
e12 = negX + e12;
e13 = e23*negX + e13;
e14 = e24*negX + e14;
m4x4Inv.e12 = negX/div;
// row 3 math
negX = -e32;
e32 = negX + e32;
e33 = e23*negX + e33;
e34 = e24*negX + e34;
m4x4Inv.e32 = negX/div;
// row 4 math
negX = -e42;
e42 = negX + e42;
e43 = e23*negX + e43;
e44 = e24*negX + e44;
m4x4Inv.e42 = negX/div;
/***************************** Step 3 of Gaussian Elimination *******************************/
// row 3 math
div = e33;
e33 = e33/div;
e34 = e34/div;
m4x4Inv.e33 = 1/div;
//row 1 math
negX = -e13;
e13 = negX + e13;
e14 = e34*negX + e14;
m4x4Inv.e13 = negX/div;
// row 2 math
negX = -e23;
e23 = negX + e23;
e24 = e34*negX + e24;
m4x4Inv.e23 = negX/div;
// row 4 math
negX = -e43;
e43 = negX + e43;
e44 = e34*negX + e44;
m4x4Inv.e43 = negX/div;
/***************************** Step 4 of Gaussian Elimination *******************************/
// row 4 math
div = e44;
e44 = e44/div;
m4x4Inv.e44 = 1/div;
// row 1 math
negX = -e14;
e14 = negX + e14;
m4x4Inv.e14 = negX/div;
// row 2 math
negX = -e24;
e24 = negX + e24;
m4x4Inv.e24 = negX/div;
// row 3 math
negX = -e34;
e34 = negX + e34;
m4x4Inv.e34 = negX/div;
// assign m4x4Inv to the matrix element by element
e11 = m4x4Inv.e11; e12 = m4x4Inv.e12; e13 = m4x4Inv.e13; e14 = m4x4Inv.e14;
e21 = m4x4Inv.e21; e22 = m4x4Inv.e22; e23 = m4x4Inv.e23; e24 = m4x4Inv.e24;
e31 = m4x4Inv.e31; e32 = m4x4Inv.e32; e33 = m4x4Inv.e33; e34 = m4x4Inv.e34;
e41 = m4x4Inv.e41; e42 = m4x4Inv.e42; e43 = m4x4Inv.e43; e44 = m4x4Inv.e44;
}
return *this;
}
Posted 08 March 2004 - 09:12 AM
quote:
Original post by optics
On the first part, I thank you for that but I also want to ask what you meant by naive method.
quote:
Second part: you have to find out of the matrix has a determinant otherwise you will be doing calculations for nothing.
Posted 08 March 2004 - 11:31 AM
/** Step 2 of Gaussian Elimination **/
// row 3 math
// div = e22
negX = -e32;
e32 = negX + e32;
e33 = e23*negX + e33;
e34 = e24*negX + e34;
m4x4Inv.e31 = negX/div*m4x4Inv.e31;
m4x4Inv.e32 = negX/div;
Guest, the last post of this topic is over 60 days old and at this point you may not reply in this topic. If you wish to continue this conversation start a new topic.
GameDev.net™, the GameDev.net logo, and GDNet™ are trademarks of GameDev.net, LLC.