Mercurial > hg > octave-lyh
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