diff liboctave/fCMatrix.cc @ 9661:afcf852256d2

optimize / and '\ for triangular matrices
author Jaroslav Hajek <highegg@gmail.com>
date Wed, 23 Sep 2009 10:00:16 +0200
parents 3429c956de6f
children 7e5b4de5fbfe
line wrap: on
line diff
--- a/liboctave/fCMatrix.cc
+++ b/liboctave/fCMatrix.cc
@@ -1845,7 +1845,7 @@
 FloatComplexMatrix::utsolve (MatrixType &mattype, const FloatComplexMatrix& b, 
 			octave_idx_type& info, float& rcon, 
 			solve_singularity_handler sing_handler,
-			bool calc_cond) const
+			bool calc_cond, blas_trans_type transt) const
 {
   FloatComplexMatrix retval;
 
@@ -1921,7 +1921,7 @@
 		  FloatComplex *result = retval.fortran_vec ();
 
 		  char uplo = 'U';
-		  char trans = 'N';
+		  char trans = get_blas_char (transt);
 		  char dia = 'N';
 
 		  F77_XFCN (ctrtrs, CTRTRS, (F77_CONST_CHAR_ARG2 (&uplo, 1), 
@@ -1946,7 +1946,7 @@
 FloatComplexMatrix::ltsolve (MatrixType &mattype, const FloatComplexMatrix& b, 
 			octave_idx_type& info, float& rcon, 
 			solve_singularity_handler sing_handler,
-			bool calc_cond) const
+			bool calc_cond, blas_trans_type transt) const
 {
   FloatComplexMatrix retval;
 
@@ -2022,7 +2022,7 @@
 		  FloatComplex *result = retval.fortran_vec ();
 
 		  char uplo = 'L';
-		  char trans = 'N';
+		  char trans = get_blas_char (transt);
 		  char dia = 'N';
 
 		  F77_XFCN (ctrtrs, CTRTRS, (F77_CONST_CHAR_ARG2 (&uplo, 1), 
@@ -2253,10 +2253,10 @@
 FloatComplexMatrix
 FloatComplexMatrix::solve (MatrixType &typ, const FloatMatrix& b, octave_idx_type& info, 
 		      float& rcon, solve_singularity_handler sing_handler,
-		      bool singular_fallback) const
+		      bool singular_fallback, blas_trans_type transt) const
 {
   FloatComplexMatrix tmp (b);
-  return solve (typ, tmp, info, rcon, sing_handler, singular_fallback);
+  return solve (typ, tmp, info, rcon, sing_handler, singular_fallback, transt);
 }
 
 FloatComplexMatrix
@@ -2286,7 +2286,7 @@
 FloatComplexMatrix::solve (MatrixType &mattype, const FloatComplexMatrix& b, 
 		      octave_idx_type& info, float& rcon,
 		      solve_singularity_handler sing_handler,
-		      bool singular_fallback) const
+		      bool singular_fallback, blas_trans_type transt) const
 {
   FloatComplexMatrix retval;
   int typ = mattype.type ();
@@ -2296,9 +2296,13 @@
 
   // Only calculate the condition number for LU/Cholesky
   if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
-    retval = utsolve (mattype, b, info, rcon, sing_handler, false);
+    retval = utsolve (mattype, b, info, rcon, sing_handler, false, transt);
   else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
-    retval = ltsolve (mattype, b, info, rcon, sing_handler, false);
+    retval = ltsolve (mattype, b, info, rcon, sing_handler, false, transt);
+  else if (transt == blas_trans)
+    return transpose ().solve (mattype, b, info, rcon, sing_handler, singular_fallback);
+  else if (transt == blas_conj_trans)
+    retval = hermitian ().solve (mattype, b, info, rcon, sing_handler, singular_fallback);
   else if (typ == MatrixType::Full || typ == MatrixType::Hermitian)
     retval = fsolve (mattype, b, info, rcon, sing_handler, true);
   else if (typ != MatrixType::Rectangular)
@@ -2343,9 +2347,9 @@
 FloatComplexColumnVector
 FloatComplexMatrix::solve (MatrixType &typ, const FloatColumnVector& b, 
 		      octave_idx_type& info, float& rcon,
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
-  return solve (typ, FloatComplexColumnVector (b), info, rcon, sing_handler);
+  return solve (typ, FloatComplexColumnVector (b), info, rcon, sing_handler, transt);
 }
 
 FloatComplexColumnVector
@@ -2374,11 +2378,11 @@
 FloatComplexColumnVector
 FloatComplexMatrix::solve (MatrixType &typ, const FloatComplexColumnVector& b,
 		      octave_idx_type& info, float& rcon,
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
 
   FloatComplexMatrix tmp (b);
-  return solve (typ, tmp, info, rcon, sing_handler).column(static_cast<octave_idx_type> (0));
+  return solve (typ, tmp, info, rcon, sing_handler, transt).column(static_cast<octave_idx_type> (0));
 }
 
 FloatComplexMatrix
@@ -2404,10 +2408,10 @@
 
 FloatComplexMatrix
 FloatComplexMatrix::solve (const FloatMatrix& b, octave_idx_type& info, float& rcon,
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
   FloatComplexMatrix tmp (b);
-  return solve (tmp, info, rcon, sing_handler);
+  return solve (tmp, info, rcon, sing_handler, transt);
 }
 
 FloatComplexMatrix
@@ -2433,10 +2437,10 @@
 
 FloatComplexMatrix
 FloatComplexMatrix::solve (const FloatComplexMatrix& b, octave_idx_type& info, float& rcon,
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
   MatrixType mattype (*this);
-  return solve (mattype, b, info, rcon, sing_handler);
+  return solve (mattype, b, info, rcon, sing_handler, true, transt);
 }
 
 FloatComplexColumnVector
@@ -2464,9 +2468,9 @@
 FloatComplexColumnVector
 FloatComplexMatrix::solve (const FloatColumnVector& b, octave_idx_type& info, 
 		      float& rcon, 
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
-  return solve (FloatComplexColumnVector (b), info, rcon, sing_handler);
+  return solve (FloatComplexColumnVector (b), info, rcon, sing_handler, transt);
 }
 
 FloatComplexColumnVector
@@ -2494,10 +2498,10 @@
 FloatComplexColumnVector
 FloatComplexMatrix::solve (const FloatComplexColumnVector& b, octave_idx_type& info,
 		      float& rcon,
-		      solve_singularity_handler sing_handler) const
+		      solve_singularity_handler sing_handler, blas_trans_type transt) const
 {
   MatrixType mattype (*this);
-  return solve (mattype, b, info, rcon, sing_handler);
+  return solve (mattype, b, info, rcon, sing_handler, transt);
 }
 
 FloatComplexMatrix