From a7c96a236e74e53e408c0d728cf0fe2854b24939 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Fri, 2 Feb 2018 12:43:56 +0100 Subject: [PATCH] Replace n_rows() by m(), n_cols() by n() in lapack_full_matrix.cc --- source/lac/lapack_full_matrix.cc | 248 +++++++++++++++---------------- 1 file changed, 124 insertions(+), 124 deletions(-) diff --git a/source/lac/lapack_full_matrix.cc b/source/lac/lapack_full_matrix.cc index b84a417b72..afb761b880 100644 --- a/source/lac/lapack_full_matrix.cc +++ b/source/lac/lapack_full_matrix.cc @@ -99,11 +99,11 @@ template void LAPACKFullMatrix::remove_row_and_column (const size_type row, const size_type col) { - AssertIndexRange (row,this->n_rows()); - AssertIndexRange (col,this->n_cols()); + AssertIndexRange (row,this->m()); + AssertIndexRange (col,this->n()); - const size_type nrows = this->n_rows()-1; - const size_type ncols = this->n_cols()-1; + const size_type nrows = this->m()-1; + const size_type ncols = this->n()-1; TransposeTable copy(std::move(*this)); this->TransposeTable::reinit (nrows, ncols); @@ -136,10 +136,10 @@ template LAPACKFullMatrix & LAPACKFullMatrix::operator = (const FullMatrix &M) { - Assert (this->n_rows() == M.n_rows(), ExcDimensionMismatch(this->n_rows(), M.n_rows())); - Assert (this->n_cols() == M.n(), ExcDimensionMismatch(this->n_cols(), M.n())); - for (size_type i=0; in_rows(); ++i) - for (size_type j=0; jn_cols(); ++j) + Assert (this->m() == M.m(), ExcDimensionMismatch(this->m(), M.m())); + Assert (this->n() == M.n(), ExcDimensionMismatch(this->n(), M.n())); + for (size_type i=0; im(); ++i) + for (size_type j=0; jn(); ++j) (*this)(i,j) = M(i,j); state = LAPACKSupport::matrix; @@ -153,10 +153,10 @@ template LAPACKFullMatrix & LAPACKFullMatrix::operator = (const SparseMatrix &M) { - Assert (this->n_rows() == M.n(), ExcDimensionMismatch(this->n_rows(), M.n())); - Assert (this->n_cols() == M.m(), ExcDimensionMismatch(this->n_cols(), M.m())); - for (size_type i=0; in_rows(); ++i) - for (size_type j=0; jn_cols(); ++j) + Assert (this->m() == M.n(), ExcDimensionMismatch(this->m(), M.n())); + Assert (this->n() == M.m(), ExcDimensionMismatch(this->n(), M.m())); + for (size_type i=0; im(); ++i) + for (size_type j=0; jn(); ++j) (*this)(i,j) = M.el(i,j); state = LAPACKSupport::matrix; @@ -346,7 +346,7 @@ LAPACKFullMatrix::rank1_update(const number a, if (state == LAPACKSupport::matrix) { { - const types::blas_int N = this->n_rows(); + const types::blas_int N = this->m(); const char uplo = LAPACKSupport::U; const types::blas_int lda = N; const types::blas_int incx=1; @@ -354,7 +354,7 @@ LAPACKFullMatrix::rank1_update(const number a, syr(&uplo, &N, &a, v.begin(), &incx, this->values.begin(), &lda); } - const size_type N = this->n_rows(); + const size_type N = this->m(); // FIXME: we should really only update upper or lower triangular parts // of symmetric matrices and make sure the interface is consistent, // for example operator(i,j) gives correct results regardless of storage. @@ -379,8 +379,8 @@ LAPACKFullMatrix::vmult ( const Vector &v, const bool adding) const { - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); const number null = 0.; @@ -394,8 +394,8 @@ LAPACKFullMatrix::vmult ( Assert (adding == false, ExcNotImplemented()); - AssertDimension(v.size(), this->n_cols()); - AssertDimension(w.size(), this->n_rows()); + AssertDimension(v.size(), this->n()); + AssertDimension(w.size(), this->m()); const char diag = 'N'; const char trans = 'N'; @@ -419,8 +419,8 @@ LAPACKFullMatrix::vmult ( case matrix: case inverse_matrix: { - AssertDimension(v.size(), this->n_cols()); - AssertDimension(w.size(), this->n_rows()); + AssertDimension(v.size(), this->n()); + AssertDimension(w.size(), this->m()); gemv("N", &mm, &nn, &alpha, &this->values[0], &mm, v.values.get(), &one, &beta, w.values.get(), &one); break; @@ -428,8 +428,8 @@ LAPACKFullMatrix::vmult ( case svd: { Threads::Mutex::ScopedLock lock (mutex); - AssertDimension(v.size(), this->n_cols()); - AssertDimension(w.size(), this->n_rows()); + AssertDimension(v.size(), this->n()); + AssertDimension(w.size(), this->m()); // Compute V^T v work.resize(std::max(mm,nn)); gemv("N", &nn, &nn, &alpha, &svd_vt->values[0], &nn, v.values.get(), &one, &null, work.data(), &one); @@ -443,8 +443,8 @@ LAPACKFullMatrix::vmult ( case inverse_svd: { Threads::Mutex::ScopedLock lock (mutex); - AssertDimension(w.size(), this->n_cols()); - AssertDimension(v.size(), this->n_rows()); + AssertDimension(w.size(), this->n()); + AssertDimension(v.size(), this->m()); // Compute U^T v work.resize(std::max(mm,nn)); gemv("T", &mm, &mm, &alpha, &svd_u->values[0], &mm, v.values.get(), &one, &null, work.data(), &one); @@ -468,8 +468,8 @@ LAPACKFullMatrix::Tvmult ( const Vector &v, const bool adding) const { - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); const number null = 0.; @@ -483,8 +483,8 @@ LAPACKFullMatrix::Tvmult ( Assert (adding == false, ExcNotImplemented()); - AssertDimension(v.size(), this->n_cols()); - AssertDimension(w.size(), this->n_rows()); + AssertDimension(v.size(), this->n()); + AssertDimension(w.size(), this->m()); const char diag = 'N'; const char trans = 'T'; @@ -509,8 +509,8 @@ LAPACKFullMatrix::Tvmult ( case matrix: case inverse_matrix: { - AssertDimension(w.size(), this->n_cols()); - AssertDimension(v.size(), this->n_rows()); + AssertDimension(w.size(), this->n()); + AssertDimension(v.size(), this->m()); gemv("T", &mm, &nn, &alpha, &this->values[0], &mm, v.values.get(), &one, &beta, w.values.get(), &one); break; @@ -518,8 +518,8 @@ LAPACKFullMatrix::Tvmult ( case svd: { Threads::Mutex::ScopedLock lock (mutex); - AssertDimension(w.size(), this->n_cols()); - AssertDimension(v.size(), this->n_rows()); + AssertDimension(w.size(), this->n()); + AssertDimension(v.size(), this->m()); // Compute U^T v work.resize(std::max(mm,nn)); @@ -534,8 +534,8 @@ LAPACKFullMatrix::Tvmult ( case inverse_svd: { Threads::Mutex::ScopedLock lock (mutex); - AssertDimension(v.size(), this->n_cols()); - AssertDimension(w.size(), this->n_rows()); + AssertDimension(v.size(), this->n()); + AssertDimension(w.size(), this->m()); // Compute V^T v work.resize(std::max(mm,nn)); @@ -580,12 +580,12 @@ LAPACKFullMatrix::mmult(LAPACKFullMatrix &C, Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); Assert(C.state == matrix || C.state == inverse_matrix, ExcState(state)); - Assert (this->n_cols() == B.n_rows(), ExcDimensionMismatch(this->n_cols(), B.n_rows())); - Assert (C.n_cols() == B.n_cols(), ExcDimensionMismatch(C.n_cols(), B.n_cols())); - Assert (C.n_rows() == this->n_rows(), ExcDimensionMismatch(this->n_rows(), C.n_rows())); - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = B.n_cols(); - const types::blas_int kk = this->n_cols(); + Assert (this->n() == B.m(), ExcDimensionMismatch(this->n(), B.m())); + Assert (C.n() == B.n(), ExcDimensionMismatch(C.n(), B.n())); + Assert (C.m() == this->m(), ExcDimensionMismatch(this->m(), C.m())); + const types::blas_int mm = this->m(); + const types::blas_int nn = B.n(); + const types::blas_int kk = this->n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -602,12 +602,12 @@ LAPACKFullMatrix::mmult(FullMatrix &C, { Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); - Assert (this->n_cols() == B.n_rows(), ExcDimensionMismatch(this->n_cols(), B.n_rows())); - Assert (C.n_cols() == B.n_cols(), ExcDimensionMismatch(C.n_cols(), B.n_cols())); - Assert (C.n_rows() == this->n_rows(), ExcDimensionMismatch(this->n_rows(), C.n_rows())); - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = B.n_cols(); - const types::blas_int kk = this->n_cols(); + Assert (this->n() == B.m(), ExcDimensionMismatch(this->n(), B.m())); + Assert (C.n() == B.n(), ExcDimensionMismatch(C.n(), B.n())); + Assert (C.m() == this->m(), ExcDimensionMismatch(this->m(), C.m())); + const types::blas_int mm = this->m(); + const types::blas_int nn = B.n(); + const types::blas_int kk = this->n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -628,12 +628,12 @@ LAPACKFullMatrix::Tmmult(LAPACKFullMatrix &C, Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); Assert(C.state == matrix || C.state == inverse_matrix, ExcState(state)); - Assert (this->n_rows() == B.n_rows(), ExcDimensionMismatch(this->n_rows(), B.n_rows())); - Assert (C.n_cols() == B.n_cols(), ExcDimensionMismatch(C.n_cols(), B.n_cols())); - Assert (C.n_rows() == this->n_cols(), ExcDimensionMismatch(this->n_cols(), C.n_rows())); - const types::blas_int mm = this->n_cols(); - const types::blas_int nn = B.n_cols(); - const types::blas_int kk = B.n_rows(); + Assert (this->m() == B.m(), ExcDimensionMismatch(this->m(), B.m())); + Assert (C.n() == B.n(), ExcDimensionMismatch(C.n(), B.n())); + Assert (C.m() == this->n(), ExcDimensionMismatch(this->n(), C.m())); + const types::blas_int mm = this->n(); + const types::blas_int nn = B.n(); + const types::blas_int kk = B.m(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -650,12 +650,12 @@ LAPACKFullMatrix::Tmmult(FullMatrix &C, { Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); - Assert (this->n_rows() == B.n_rows(), ExcDimensionMismatch(this->n_rows(), B.n_rows())); - Assert (C.n_cols() == B.n_cols(), ExcDimensionMismatch(C.n_cols(), B.n_cols())); - Assert (C.n_rows() == this->n_cols(), ExcDimensionMismatch(this->n_cols(), C.n_rows())); - const types::blas_int mm = this->n_cols(); - const types::blas_int nn = B.n_cols(); - const types::blas_int kk = B.n_rows(); + Assert (this->m() == B.m(), ExcDimensionMismatch(this->m(), B.m())); + Assert (C.n() == B.n(), ExcDimensionMismatch(C.n(), B.n())); + Assert (C.m() == this->n(), ExcDimensionMismatch(this->n(), C.m())); + const types::blas_int mm = this->n(); + const types::blas_int nn = B.n(); + const types::blas_int kk = B.m(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -675,12 +675,12 @@ LAPACKFullMatrix::mTmult(LAPACKFullMatrix &C, Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); Assert(C.state == matrix || C.state == inverse_matrix, ExcState(state)); - Assert (this->n_cols() == B.n_cols(), ExcDimensionMismatch(this->n_cols(), B.n_cols())); - Assert (C.n_cols() == B.n_rows(), ExcDimensionMismatch(C.n_cols(), B.n_rows())); - Assert (C.n_rows() == this->n_rows(), ExcDimensionMismatch(this->n_rows(), C.n_rows())); - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = B.n_rows(); - const types::blas_int kk = B.n_cols(); + Assert (this->n() == B.n(), ExcDimensionMismatch(this->n(), B.n())); + Assert (C.n() == B.m(), ExcDimensionMismatch(C.n(), B.m())); + Assert (C.m() == this->m(), ExcDimensionMismatch(this->m(), C.m())); + const types::blas_int mm = this->m(); + const types::blas_int nn = B.m(); + const types::blas_int kk = B.n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -698,12 +698,12 @@ LAPACKFullMatrix::mTmult(FullMatrix &C, { Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); - Assert (this->n_cols() == B.n_cols(), ExcDimensionMismatch(this->n_cols(), B.n_cols())); - Assert (C.n_cols() == B.n_rows(), ExcDimensionMismatch(C.n_cols(), B.n_rows())); - Assert (C.n_rows() == this->n_rows(), ExcDimensionMismatch(this->n_rows(), C.n_rows())); - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = B.n_rows(); - const types::blas_int kk = B.n_cols(); + Assert (this->n() == B.n(), ExcDimensionMismatch(this->n(), B.n())); + Assert (C.n() == B.m(), ExcDimensionMismatch(C.n(), B.m())); + Assert (C.m() == this->m(), ExcDimensionMismatch(this->m(), C.m())); + const types::blas_int mm = this->m(); + const types::blas_int nn = B.m(); + const types::blas_int kk = B.n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -723,12 +723,12 @@ LAPACKFullMatrix::TmTmult(LAPACKFullMatrix &C, Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); Assert(C.state == matrix || C.state == inverse_matrix, ExcState(state)); - Assert (this->n_rows() == B.n_cols(), ExcDimensionMismatch(this->n_rows(), B.n_cols())); - Assert (C.n_cols() == B.n_rows(), ExcDimensionMismatch(C.n_cols(), B.n_rows())); - Assert (C.n_rows() == this->n_cols(), ExcDimensionMismatch(this->n_cols(), C.n_rows())); - const types::blas_int mm = this->n_cols(); - const types::blas_int nn = B.n_rows(); - const types::blas_int kk = B.n_cols(); + Assert (this->m() == B.n(), ExcDimensionMismatch(this->m(), B.n())); + Assert (C.n() == B.m(), ExcDimensionMismatch(C.n(), B.m())); + Assert (C.m() == this->n(), ExcDimensionMismatch(this->n(), C.m())); + const types::blas_int mm = this->n(); + const types::blas_int nn = B.m(); + const types::blas_int kk = B.n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -745,12 +745,12 @@ LAPACKFullMatrix::TmTmult(FullMatrix &C, { Assert(state == matrix || state == inverse_matrix, ExcState(state)); Assert(B.state == matrix || B.state == inverse_matrix, ExcState(state)); - Assert (this->n_rows() == B.n_cols(), ExcDimensionMismatch(this->n_rows(), B.n_cols())); - Assert (C.n_cols() == B.n_rows(), ExcDimensionMismatch(C.n_cols(), B.n_rows())); - Assert (C.n_rows() == this->n_cols(), ExcDimensionMismatch(this->n_cols(), C.n_rows())); - const types::blas_int mm = this->n_cols(); - const types::blas_int nn = B.n_rows(); - const types::blas_int kk = B.n_cols(); + Assert (this->m() == B.n(), ExcDimensionMismatch(this->m(), B.n())); + Assert (C.n() == B.m(), ExcDimensionMismatch(C.n(), B.m())); + Assert (C.m() == this->n(), ExcDimensionMismatch(this->n(), C.m())); + const types::blas_int mm = this->n(); + const types::blas_int nn = B.m(); + const types::blas_int kk = B.n(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); @@ -768,8 +768,8 @@ LAPACKFullMatrix::compute_lu_factorization() Assert(state == matrix, ExcState(state)); state = LAPACKSupport::unusable; - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); number *const values = &this->values[0]; ipiv.resize(mm); types::blas_int info = 0; @@ -830,8 +830,8 @@ number LAPACKFullMatrix::norm(const char type) const state == LAPACKSupport::inverse_matrix, ExcMessage("norms can be called in matrix state only.")); - const types::blas_int N = this->n_cols(); - const types::blas_int M = this->n_rows(); + const types::blas_int N = this->n(); + const types::blas_int M = this->m(); const number *const values = &this->values[0]; if (property == symmetric) { @@ -861,11 +861,11 @@ number LAPACKFullMatrix::trace() const Assert (state == LAPACKSupport::matrix || state == LAPACKSupport::inverse_matrix, ExcMessage("Trace can be called in matrix state only.")); - Assert (this->n_cols() == this->n_rows(), - ExcDimensionMismatch(this->n_cols(), this->n_rows())); + Assert (this->n() == this->m(), + ExcDimensionMismatch(this->n(), this->m())); number tr = 0; - for (size_type i=0; in_rows(); ++i) + for (size_type i=0; im(); ++i) tr += (*this)(i,i); return tr; @@ -881,8 +881,8 @@ LAPACKFullMatrix::compute_cholesky_factorization() Assert(property == symmetric, ExcProperty(property)); state = LAPACKSupport::unusable; - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); (void) mm; Assert (mm == nn, ExcDimensionMismatch(mm,nn)); @@ -908,7 +908,7 @@ LAPACKFullMatrix::reciprocal_condition_number(const number a_norm) const Assert(state == cholesky, ExcState(state)); number rcond = 0.; - const types::blas_int N = this->n_rows(); + const types::blas_int N = this->m(); const number *values = &this->values[0]; types::blas_int info = 0; const types::blas_int lda = std::max(1,N); @@ -937,7 +937,7 @@ LAPACKFullMatrix::reciprocal_condition_number() const ExcProperty(property)); number rcond = 0.; - const types::blas_int N = this->n_rows(); + const types::blas_int N = this->m(); const number *const values = &this->values[0]; types::blas_int info = 0; const types::blas_int lda = std::max(1,N); @@ -966,8 +966,8 @@ LAPACKFullMatrix::compute_svd() Assert(state == matrix, ExcState(state)); state = LAPACKSupport::unusable; - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); number *const values = &this->values[0]; wr.resize(std::max(mm,nn)); std::fill(wr.begin(), wr.end(), 0.); @@ -1031,8 +1031,8 @@ LAPACKFullMatrix::invert() { Assert(state == matrix || state == lu || state == cholesky, ExcState(state)); - const types::blas_int mm = this->n_rows(); - const types::blas_int nn = this->n_cols(); + const types::blas_int mm = this->m(); + const types::blas_int nn = this->n(); Assert (nn == mm, ExcNotQuadratic()); number *const values = &this->values[0]; @@ -1073,11 +1073,11 @@ void LAPACKFullMatrix::solve(Vector &v, const bool transposed) const { - Assert(this->n_rows() == this->n_cols(), + Assert(this->m() == this->n(), LACExceptions::ExcNotQuadratic()); - AssertDimension(this->n_rows(), v.size()); + AssertDimension(this->m(), v.size()); const char *trans = transposed ? &T : &N; - const types::blas_int nn = this->n_cols(); + const types::blas_int nn = this->n(); const number *const values = &this->values[0]; const types::blas_int n_rhs = 1; types::blas_int info = 0; @@ -1122,13 +1122,13 @@ LAPACKFullMatrix::solve(LAPACKFullMatrix &B, { Assert(B.state == matrix, ExcState(B.state)); - Assert(this->n_rows() == this->n_cols(), + Assert(this->m() == this->n(), LACExceptions::ExcNotQuadratic()); - AssertDimension(this->n_rows(), B.n_rows()); + AssertDimension(this->m(), B.m()); const char *trans = transposed ? &T : &N; - const types::blas_int nn = this->n_cols(); + const types::blas_int nn = this->n(); const number *const values = &this->values[0]; - const types::blas_int n_rhs = B.n_cols(); + const types::blas_int n_rhs = B.n(); types::blas_int info = 0; if (state == lu) @@ -1188,7 +1188,7 @@ template number LAPACKFullMatrix::determinant() const { - Assert(this->n_rows() == this->n_cols(), LACExceptions::ExcNotQuadratic()); + Assert(this->m() == this->n(), LACExceptions::ExcNotQuadratic()); // LAPACK doesn't offer a function to compute a matrix determinant. // This is due to the difficulty in maintaining numerical accuracy, as the @@ -1206,9 +1206,9 @@ LAPACKFullMatrix::determinant() const // http://icl.cs.utk.edu/lapack-forum/viewtopic.php?p=341&#p336 // for further information. Assert(state == lu, ExcState(state)); - Assert(ipiv.size() == this->n_rows(), ExcInternalError()); + Assert(ipiv.size() == this->m(), ExcInternalError()); number det = 1.0; - for (size_type i=0; in_rows(); ++i) + for (size_type i=0; im(); ++i) det *= ( ipiv[i] == types::blas_int(i+1) ? this->el(i,i) : -this->el(i,i) ); return det; } @@ -1220,7 +1220,7 @@ LAPACKFullMatrix::compute_eigenvalues(const bool right, const bool left) { Assert(state == matrix, ExcState(state)); - const types::blas_int nn = this->n_cols(); + const types::blas_int nn = this->n(); wr.resize(nn); wi.resize(nn); if (right) vr.resize(nn*nn); @@ -1284,8 +1284,8 @@ LAPACKFullMatrix::compute_eigenvalues_symmetric(const number lowe FullMatrix &eigenvectors) { Assert(state == matrix, ExcState(state)); - const types::blas_int nn = (this->n_cols() > 0 ? this->n_cols() : 1); - Assert(static_cast(nn) == this->n_rows(), ExcNotQuadratic()); + const types::blas_int nn = (this->n() > 0 ? this->n() : 1); + Assert(static_cast(nn) == this->m(), ExcNotQuadratic()); wr.resize(nn); LAPACKFullMatrix matrix_eigenvectors(nn, nn); @@ -1373,11 +1373,11 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric( const types::blas_int itype) { Assert(state == matrix, ExcState(state)); - const types::blas_int nn = (this->n_cols() > 0 ? this->n_cols() : 1); - Assert(static_cast(nn) == this->n_rows(), ExcNotQuadratic()); - Assert(B.n_rows() == B.n_cols(), ExcNotQuadratic()); - Assert(static_cast(nn) == B.n_cols(), - ExcDimensionMismatch (nn, B.n_cols())); + const types::blas_int nn = (this->n() > 0 ? this->n() : 1); + Assert(static_cast(nn) == this->m(), ExcNotQuadratic()); + Assert(B.m() == B.n(), ExcNotQuadratic()); + Assert(static_cast(nn) == B.n(), + ExcDimensionMismatch (nn, B.n())); wr.resize(nn); LAPACKFullMatrix matrix_eigenvectors(nn, nn); @@ -1461,13 +1461,13 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric ( const types::blas_int itype) { Assert(state == matrix, ExcState(state)); - const types::blas_int nn = this->n_cols(); - Assert(static_cast(nn) == this->n_rows(), ExcNotQuadratic()); - Assert(B.n_rows() == B.n_cols(), ExcNotQuadratic()); - Assert(static_cast(nn) == B.n_cols(), - ExcDimensionMismatch (nn, B.n_cols())); + const types::blas_int nn = this->n(); + Assert(static_cast(nn) == this->m(), ExcNotQuadratic()); + Assert(B.m() == B.n(), ExcNotQuadratic()); + Assert(static_cast(nn) == B.n(), + ExcDimensionMismatch (nn, B.n())); Assert(eigenvectors.size() <= static_cast(nn), - ExcMessage ("eigenvectors.size() > matrix.n_cols()")); + ExcMessage ("eigenvectors.size() > matrix.n()")); wr.resize(nn); wi.resize(nn); //This is set purely for consistency reasons with the @@ -1541,7 +1541,7 @@ LAPACKFullMatrix::print_formatted ( { unsigned int width = width_; - Assert ((!this->empty()) || (this->n_cols()+this->n_rows()==0), + Assert ((!this->empty()) || (this->n()+this->m()==0), ExcInternalError()); Assert (state == LAPACKSupport::matrix || state == LAPACKSupport::inverse_matrix || @@ -1566,10 +1566,10 @@ LAPACKFullMatrix::print_formatted ( width = precision+2; } - for (size_type i=0; in_rows(); ++i) + for (size_type i=0; im(); ++i) { // Cholesky is stored in lower triangular, so just output this part: - const size_type nc = state == LAPACKSupport::cholesky ? i+1 : this->n_cols(); + const size_type nc = state == LAPACKSupport::cholesky ? i+1 : this->n(); for (size_type j=0; j