Changeset 5675 in orxonox.OLD for trunk/src/lib/math

Timestamp:

Nov 21, 2005, 9:59:09 PM (19 years ago)

Author:

bensch

Message:

orxonox/trunk: optimisations

Location:

trunk/src/lib/math

Files:

• matrix.cc (modified) (2 diffs)
• matrix.h (modified) (1 diff)

trunk/src/lib/math/matrix.cc

@@ -19 +19 @@
 
 
-Vector Matrix::getEigenValues() const
+int Matrix::getEigenValues(Vector& eigenValues) const
 {
-  Vector eigVl;
-
+  int retVal = -1;
   float a = 0;
   float b = 0;
@@ -61 +60 @@
     float p = sqrt((b/2.0)*(b/2.0) - Q);
 
-    eigVl.x = c[2]/3.0 + 2 * pow(p, 1/3.0) * cos(psi/3.0);
-    eigVl.y = c[2]/3.0 - pow(p, 1/3.0) * (cos(psi/3.0)
+    eigenValues.x = c[2]/3.0 + 2 * pow(p, 1/3.0) * cos(psi/3.0);
+    eigenValues.y = c[2]/3.0 - pow(p, 1/3.0) * (cos(psi/3.0)
         + sqrt(3.0) * sin(psi/3.0));
-    eigVl.z = c[2]/3.0 - pow(p, 1/3.0) * (cos(psi/3.0)
+    eigenValues.z = c[2]/3.0 - pow(p, 1/3.0) * (cos(psi/3.0)
         - sqrt(3.0) * sin(psi/3.0));
-
+    retVal = 3;
   }
   // 2 Distinct Roots
   else if (Q == 0)
   {
-    eigVl.x = c[2]/3.0 + pow(b/2.0, 1.0/3.0);
-    eigVl.y = c[2]/3.0 + pow(b/2.0, 1.0/3.0);
-    eigVl.z = c[2]/3.0 + 2* pow(b/2.0, 1.0/3.0);
+    eigenValues.x = eigenValues.y = c[2]/3.0 + pow(b/2.0, 1.0/3.0);
+    eigenValues.z = c[2]/3.0 + 2* pow(b/2.0, 1.0/3.0);
+    retVal = 2;
   }
   // 1 Root (not calculating anything.)
   else if (Q > 0)
   {
-    printf("This Matrix is a multiple of the Identity matrix (lambda * I3))\n");
-    eigVl.x = eigVl.y = eigVl.z = 1;
+    eigenValues.x = eigenValues.y = eigenValues.z = 1;
+    retVal = 1;
   }
-  return eigVl;
+  return retVal;
 }
 
-void Matrix::getEigenVectors(Vector& a, Vector& b, Vector& c) const
+void Matrix::getEigenVectors(Vector& eigVc1, Vector& eigVc2, Vector& eigVc3) const
 {
-  Vector eigVl = this->getEigenValues();
+  Vector eigenValues;
+  int eigenValuesCount = this->getEigenValues(eigenValues);
 
-  float eigVal[3] = { eigVl.x, eigVl.y, eigVl.z };
+  /* eigenvec creation */
+  eigVc1.x = -1/this->m13*(this->m33 - eigenValues.x) +
+      (this->m32*(-this->m31*this->m32 + this->m12*this->m33 - this->m12*eigenValues.x)) /
+      this->m13*(-this->m13*this->m22 - this->m12*this->m23 + this->m13*eigenValues.x);
 
-  Vector eigVc[3];
-  /* eigenvec test */
-  for(int i = 0; i < 2; i++)
-  {
-    eigVc[i].x = -1/this->m13*(this->m33 - eigVal[i]) + (this->m32*(-this->m31*this->m32 + this->m12*this->m33 - this->m12*eigVal[i])) /
-        this->m13*(-this->m13*this->m22 - this->m12*this->m23 + this->m13*eigVal[i]);
+  eigVc1.y = -( -this->m13*this->m23 + this->m12*this->m33 - this->m12*eigenValues.x) /
+      (-this->m31*this->m22 + this->m12*this->m23 + this->m13*eigenValues.x);
 
-    eigVc[i].y = -( -this->m13*this->m23 + this->m12*this->m33 - this->m12*eigVal[i]) /
-        (-this->m31*this->m22 + this->m12*this->m23 + this->m13*eigVal[i]);
+  eigVc1.z = 1.0f;
 
-    eigVc[i].z = 1.0f;
+  eigVc2.x = -1/this->m13*(this->m33 - eigenValues.y) +
+      (this->m32*(-this->m31*this->m32 + this->m12*this->m33 - this->m12*eigenValues.y)) /
+      this->m13*(-this->m13*this->m22 - this->m12*this->m23 + this->m13*eigenValues.y);
 
-    eigVc[i] /= eigVc[i].len();
-  }
-  eigVc[2] = eigVc[0].cross(eigVc[1]);
+  eigVc2.y = -( -this->m13*this->m23 + this->m12*this->m33 - this->m12*eigenValues.y) /
+      (-this->m31*this->m22 + this->m12*this->m23 + this->m13*eigenValues.y);
 
-  a = eigVc[0];
-  b = eigVc[1];
-  c = eigVc[2];
+  eigVc2.z = 1.0f;
+
+  eigVc3 = eigVc1.cross(eigVc2);
+
+  eigVc1.normalize();
+  eigVc2.normalize();
+  eigVc3.normalize();
 }
 

trunk/src/lib/math/matrix.h

@@ -58 +58 @@
     }
 
-    Vector getEigenValues() const;
+    int getEigenValues(Vector& eigenVectors) const;
     void getEigenVectors(Vector& a, Vector& b, Vector& c) const;