Mercurial > jwe > octave
diff liboctave/array/fMatrix.cc @ 21120:499b851fbfae
Replace pattern if/err_XXX/else/code with if/err_XXX/ code.
* schur.cc, ov-complex.h, ov-cx-mat.cc, ov-cx-sparse.cc, ov-float.h,
ov-flt-complex.h, ov-flt-cx-mat.cc, ov-flt-re-mat.cc, ov-range.cc,
ov-re-mat.cc, ov-re-sparse.cc, ov-scalar.h, ov-str-mat.cc, ops.h, pt-idx.cc,
Array.cc, CColVector.cc, CMatrix.cc, CRowVector.cc, MSparse.cc, PermMatrix.cc,
Sparse.cc, dColVector.cc, dMatrix.cc, dRowVector.cc, dSparse.cc,
fCColVector.cc, fCMatrix.cc, fCRowVector.cc, fColVector.cc, fMatrix.cc,
fRowVector.cc:
Replace pattern if/err_XXX/else/code with if/err_XXX/ code.
author | Rik <rik@octave.org> |
---|---|
date | Wed, 20 Jan 2016 16:58:59 -0800 |
parents | 358aa7fcbd33 |
children | 228b65504557 |
line wrap: on
line diff
--- a/liboctave/array/fMatrix.cc Wed Jan 20 16:10:23 2016 -0800 +++ b/liboctave/array/fMatrix.cc Wed Jan 20 16:58:59 2016 -0800 @@ -3113,71 +3113,69 @@ if (a_nc != b_nr) err_nonconformant ("operator *", a_nr, a_nc, b_nr, b_nc); + + if (a_nr == 0 || a_nc == 0 || b_nc == 0) + retval = FloatMatrix (a_nr, b_nc, 0.0); + else if (a.data () == b.data () && a_nr == b_nc && tra != trb) + { + octave_idx_type lda = a.rows (); + + retval = FloatMatrix (a_nr, b_nc); + float *c = retval.fortran_vec (); + + const char ctra = get_blas_trans_arg (tra); + F77_XFCN (ssyrk, SSYRK, (F77_CONST_CHAR_ARG2 ("U", 1), + F77_CONST_CHAR_ARG2 (&ctra, 1), + a_nr, a_nc, 1.0, + a.data (), lda, 0.0, c, a_nr + F77_CHAR_ARG_LEN (1) + F77_CHAR_ARG_LEN (1))); + for (int j = 0; j < a_nr; j++) + for (int i = 0; i < j; i++) + retval.xelem (j,i) = retval.xelem (i,j); + + } else { - if (a_nr == 0 || a_nc == 0 || b_nc == 0) - retval = FloatMatrix (a_nr, b_nc, 0.0); - else if (a.data () == b.data () && a_nr == b_nc && tra != trb) + octave_idx_type lda = a.rows (); + octave_idx_type tda = a.cols (); + octave_idx_type ldb = b.rows (); + octave_idx_type tdb = b.cols (); + + retval = FloatMatrix (a_nr, b_nc); + float *c = retval.fortran_vec (); + + if (b_nc == 1) { - octave_idx_type lda = a.rows (); - - retval = FloatMatrix (a_nr, b_nc); - float *c = retval.fortran_vec (); - - const char ctra = get_blas_trans_arg (tra); - F77_XFCN (ssyrk, SSYRK, (F77_CONST_CHAR_ARG2 ("U", 1), - F77_CONST_CHAR_ARG2 (&ctra, 1), - a_nr, a_nc, 1.0, - a.data (), lda, 0.0, c, a_nr - F77_CHAR_ARG_LEN (1) + if (a_nr == 1) + F77_FUNC (xsdot, XSDOT) (a_nc, a.data (), 1, b.data (), 1, *c); + else + { + const char ctra = get_blas_trans_arg (tra); + F77_XFCN (sgemv, SGEMV, (F77_CONST_CHAR_ARG2 (&ctra, 1), + lda, tda, 1.0, a.data (), lda, + b.data (), 1, 0.0, c, 1 + F77_CHAR_ARG_LEN (1))); + } + } + else if (a_nr == 1) + { + const char crevtrb = get_blas_trans_arg (! trb); + F77_XFCN (sgemv, SGEMV, (F77_CONST_CHAR_ARG2 (&crevtrb, 1), + ldb, tdb, 1.0, b.data (), ldb, + a.data (), 1, 0.0, c, 1 F77_CHAR_ARG_LEN (1))); - for (int j = 0; j < a_nr; j++) - for (int i = 0; i < j; i++) - retval.xelem (j,i) = retval.xelem (i,j); - } else { - octave_idx_type lda = a.rows (); - octave_idx_type tda = a.cols (); - octave_idx_type ldb = b.rows (); - octave_idx_type tdb = b.cols (); - - retval = FloatMatrix (a_nr, b_nc); - float *c = retval.fortran_vec (); - - if (b_nc == 1) - { - if (a_nr == 1) - F77_FUNC (xsdot, XSDOT) (a_nc, a.data (), 1, b.data (), 1, *c); - else - { - const char ctra = get_blas_trans_arg (tra); - F77_XFCN (sgemv, SGEMV, (F77_CONST_CHAR_ARG2 (&ctra, 1), - lda, tda, 1.0, a.data (), lda, - b.data (), 1, 0.0, c, 1 - F77_CHAR_ARG_LEN (1))); - } - } - else if (a_nr == 1) - { - const char crevtrb = get_blas_trans_arg (! trb); - F77_XFCN (sgemv, SGEMV, (F77_CONST_CHAR_ARG2 (&crevtrb, 1), - ldb, tdb, 1.0, b.data (), ldb, - a.data (), 1, 0.0, c, 1 - F77_CHAR_ARG_LEN (1))); - } - else - { - const char ctra = get_blas_trans_arg (tra); - const char ctrb = get_blas_trans_arg (trb); - F77_XFCN (sgemm, SGEMM, (F77_CONST_CHAR_ARG2 (&ctra, 1), - F77_CONST_CHAR_ARG2 (&ctrb, 1), - a_nr, b_nc, a_nc, 1.0, a.data (), - lda, b.data (), ldb, 0.0, c, a_nr - F77_CHAR_ARG_LEN (1) - F77_CHAR_ARG_LEN (1))); - } + const char ctra = get_blas_trans_arg (tra); + const char ctrb = get_blas_trans_arg (trb); + F77_XFCN (sgemm, SGEMM, (F77_CONST_CHAR_ARG2 (&ctra, 1), + F77_CONST_CHAR_ARG2 (&ctrb, 1), + a_nr, b_nc, a_nc, 1.0, a.data (), + lda, b.data (), ldb, 0.0, c, a_nr + F77_CHAR_ARG_LEN (1) + F77_CHAR_ARG_LEN (1))); } }