From f8414bf7003aa3ad8c130df28fa79bdb89ebbdc8 Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 3 Jan 2000 15:02:17 +0000 Subject: [PATCH] Let the FullMatrix allocate no memory by default. Also, release all memory upon reinit(0). git-svn-id: https://svn.dealii.org/trunk@2154 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/full_matrix.h | 18 +- .../lac/include/lac/full_matrix.templates.h | 188 +++++++++++++----- 2 files changed, 151 insertions(+), 55 deletions(-) diff --git a/deal.II/lac/include/lac/full_matrix.h b/deal.II/lac/include/lac/full_matrix.h index 421d5628d1..00ef2217b7 100644 --- a/deal.II/lac/include/lac/full_matrix.h +++ b/deal.II/lac/include/lac/full_matrix.h @@ -73,8 +73,10 @@ class FullMatrix : public Subscriptor * conversion of integers and other types * to a matrix, this constructor is * declared #explicit#. + * + * By default, no memory is allocated. */ - explicit FullMatrix (const unsigned int n = 1); + explicit FullMatrix (const unsigned int n = 0); /** * Constructor. Initialize the matrix as @@ -166,6 +168,8 @@ class FullMatrix : public Subscriptor * Set dimension to $m(B)\times n(B)$ and * allocate memory if necessary. Forget * the previous content of the matrix. + * However, this function does not copy + * the contents of #B#. */ template void reinit (const FullMatrix &B); @@ -547,6 +551,10 @@ class FullMatrix : public Subscriptor void print_formatted (ostream &out, const unsigned int presicion=3) const; + /** + * Exception + */ + DeclException0 (ExcEmptyMatrix); /** * Exception */ @@ -624,14 +632,6 @@ class FullMatrix : public Subscriptor */ unsigned int val_size; - /** - * Initialization. Initialize - * memory for a #FullMatrix# - * of #m# rows and #n# - * columns to zero. - */ - void init (const unsigned int m, const unsigned int n); - /** * Return a read-write reference to the * element #(i,j)#. diff --git a/deal.II/lac/include/lac/full_matrix.templates.h b/deal.II/lac/include/lac/full_matrix.templates.h index 51c9d001a6..d441ee3450 100644 --- a/deal.II/lac/include/lac/full_matrix.templates.h +++ b/deal.II/lac/include/lac/full_matrix.templates.h @@ -12,80 +12,96 @@ template -FullMatrix::FullMatrix (const unsigned int n) +FullMatrix::FullMatrix (const unsigned int n) : + val (0), + dim_range (0), + dim_image (0), + val_size (0) { - init (n,n); + reinit (n,n); }; template -FullMatrix::FullMatrix (const unsigned int m, const unsigned int n) +FullMatrix::FullMatrix (const unsigned int m, + const unsigned int n) : + val (0), + dim_range (0), + dim_image (0), + val_size (0) { - init (m,n); + reinit (m,n); }; template -FullMatrix::FullMatrix (const FullMatrix &m): - Subscriptor() +FullMatrix::FullMatrix (const FullMatrix &m) : + Subscriptor (), + val (0), + dim_range (0), + dim_image (0), + val_size (0) { - init (m.dim_image, m.dim_range); + reinit (m.dim_image, m.dim_range); if (dim_range*dim_image != 0) copy (&m.val[0], &m.val[dim_image*dim_range], &val[0]); }; -template -void -FullMatrix::init (const unsigned int mm, const unsigned int nn) -{ - val_size = nn*mm; - val = new number[val_size]; - dim_range = nn; - dim_image = mm; - clear (); -}; - - template FullMatrix::~FullMatrix () { - delete[] val; + if (val != 0) + delete[] val; }; -template -bool -FullMatrix::all_zero () const -{ - const number *p = &val[0], - *e = &val[n()*m()]; - while (p!=e) - if (*p++ != 0.0) - return false; - - return true; -}; - template void FullMatrix::reinit (const unsigned int mm, const unsigned int nn) { - if (val_size @@ -105,6 +121,22 @@ FullMatrix::reinit (const FullMatrix &B) }; +template +bool +FullMatrix::all_zero () const +{ + Assert (val != 0, ExcEmptyMatrix()); + + const number *p = &val[0], + *e = &val[n()*m()]; + while (p!=e) + if (*p++ != 0.0) + return false; + + return true; +}; + + template template @@ -113,6 +145,8 @@ FullMatrix::vmult (Vector& dst, const Vector& src, const bool adding) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(dst.size() == m(), ExcDimensionMismatch(dst.size(), m())); Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n())); @@ -242,6 +276,8 @@ void FullMatrix::Tvmult (Vector &dst, const Vector &src, const bool adding) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(dst.size() == n(), ExcDimensionMismatch(dst.size(), n())); Assert(src.size() == m(), ExcDimensionMismatch(src.size(), m())); @@ -267,6 +303,8 @@ double FullMatrix::residual (Vector& dst, const Vector& src, const Vector& right) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(dst.size() == m(), ExcDimensionMismatch(dst.size(), m())); Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n())); Assert(right.size() == m(), ExcDimensionMismatch(right.size(), m())); @@ -292,6 +330,8 @@ template void FullMatrix::forward (Vector& dst, const Vector& src) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(dst.size() == m(), ExcDimensionMismatch(dst.size(), m())); Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n())); @@ -313,6 +353,8 @@ template void FullMatrix::backward (Vector& dst, const Vector& src) const { + Assert (val != 0, ExcEmptyMatrix()); + unsigned int j; unsigned int nu = (m() FullMatrix& FullMatrix::operator = (const FullMatrix& m) { - reinit(m); + reinit (m); if (dim_range*dim_image != 0) copy (&m.val[0], &m.val[dim_image*dim_range], &val[0]); @@ -376,6 +418,8 @@ void FullMatrix::add_row (const unsigned int i, const number s, const unsigned int j) { + Assert (val != 0, ExcEmptyMatrix()); + for (unsigned int k=0; k::add_row (const unsigned int i, const number t, const unsigned int k) { + Assert (val != 0, ExcEmptyMatrix()); + const unsigned int size_m = m(); for (unsigned l=0; l void FullMatrix::add_col (const unsigned int i, const number s, const unsigned int j) { + Assert (val != 0, ExcEmptyMatrix()); + for (unsigned int k=0; k::add_col (const unsigned int i, const number s, const unsigned int j, const number t, const unsigned int k) { + Assert (val != 0, ExcEmptyMatrix()); + for (unsigned int l=0; l::add_col (const unsigned int i, const number s, template void FullMatrix::swap_row (const unsigned int i, const unsigned int j) { + Assert (val != 0, ExcEmptyMatrix()); + number s; for (unsigned int k=0; k::swap_row (const unsigned int i, const unsigned int j) template void FullMatrix::swap_col (const unsigned int i, const unsigned int j) { + Assert (val != 0, ExcEmptyMatrix()); + number s; for (unsigned int k=0; k::swap_col (const unsigned int i, const unsigned int j) template void FullMatrix::diagadd (const number src) { + Assert (val != 0, ExcEmptyMatrix()); Assert (m() == n(), ExcDimensionMismatch(m(),n())); + for (unsigned int i=0; i void FullMatrix::mmult (FullMatrix& dst, const FullMatrix& src) const { + Assert (val != 0, ExcEmptyMatrix()); Assert (n() == src.m(), ExcDimensionMismatch(n(), src.m())); + unsigned int i,j,k; number2 s = 1.; dst.reinit(m(), src.n()); @@ -474,6 +532,7 @@ void FullMatrix::mmult (FullMatrix& dst, /*void FullMatrix::mmult (FullMatrix& dst, const FullMatrix& src) const { + Assert (val != 0, ExcEmptyMatrix()); Assert (m() == src.n(), ExcDimensionMismatch(m(), src.n())); unsigned int i,j,k; @@ -496,6 +555,7 @@ template template void FullMatrix::Tmmult (FullMatrix& dst, const FullMatrix& src) const { + Assert (val != 0, ExcEmptyMatrix()); Assert (n() == src.n(), ExcDimensionMismatch(n(), src.n())); unsigned int i,j,k; @@ -515,6 +575,7 @@ void FullMatrix::Tmmult (FullMatrix& dst, const FullMatrix::Tmmult(FullMatrix& dst, const FullMatrix& src) const { + Assert (val != 0, ExcEmptyMatrix()); Assert (m() == src.n(), ExcDimensionMismatch(m(), src.n())); unsigned int i,j,k; @@ -537,6 +598,8 @@ template template double FullMatrix::matrix_norm (const Vector &v) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(m() == v.size(), ExcDimensionMismatch(m(),v.size())); Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size())); @@ -563,8 +626,11 @@ double FullMatrix::matrix_norm (const Vector &v) const template template -double FullMatrix::matrix_scalar_product (const Vector &u, const Vector &v) const +double FullMatrix::matrix_scalar_product (const Vector &u, + const Vector &v) const { + Assert (val != 0, ExcEmptyMatrix()); + Assert(m() == u.size(), ExcDimensionMismatch(m(),v.size())); Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size())); @@ -593,6 +659,8 @@ double FullMatrix::matrix_scalar_product (const Vector &u, cons template number FullMatrix::l1_norm () const { + Assert (val != 0, ExcEmptyMatrix()); + number sum=0, max=0; const unsigned int n_rows = m(), n_cols = n(); @@ -612,6 +680,8 @@ number FullMatrix::l1_norm () const template number FullMatrix::linfty_norm () const { + Assert (val != 0, ExcEmptyMatrix()); + number sum=0, max=0; const unsigned int n_rows = m(), n_cols = n(); @@ -632,6 +702,8 @@ template void FullMatrix::print (ostream& s, int w, int p) const { + Assert (val != 0, ExcEmptyMatrix()); + unsigned int i,j; for (i=0;i void FullMatrix::add (const number s,const FullMatrix& src) { + Assert (val != 0, ExcEmptyMatrix()); + Assert (m() == src.m(), ExcDimensionMismatch(m(), src.m())); Assert (n() == src.n(), ExcDimensionMismatch(n(), src.n())); + if ((n()==3) && (m()==3)) { val[0] += s * src.val[0]; @@ -767,6 +842,8 @@ template void FullMatrix::add_diag (const number s, const FullMatrix& src) { + Assert (val != 0, ExcEmptyMatrix()); + Assert (m() == src.m(), ExcDimensionMismatch(m(), src.m())); Assert (n() == src.n(), ExcDimensionMismatch(n(), src.n())); @@ -888,6 +965,8 @@ template void FullMatrix::Tadd (const number s, const FullMatrix& src) { + Assert (val != 0, ExcEmptyMatrix()); + Assert (m() == n(), ExcNotQuadratic()); Assert (m() == src.m(), ExcDimensionMismatch(m(), src.m())); Assert (n() == src.n(), ExcDimensionMismatch(n(), src.n())); @@ -1012,6 +1091,8 @@ template bool FullMatrix::operator == (const FullMatrix &m) const { + Assert (val != 0, ExcEmptyMatrix()); + bool q = (dim_range==m.dim_range) && (dim_image==m.dim_image); if (!q) return false; @@ -1025,6 +1106,8 @@ template double FullMatrix::determinant () const { + Assert (val != 0, ExcEmptyMatrix()); + Assert (dim_range == dim_image, ExcDimensionMismatch(dim_range, dim_image)); Assert ((dim_range>=1) && (dim_range<=3), ExcNotImplemented(dim_range)); @@ -1053,6 +1136,8 @@ template number FullMatrix::norm2 () const { + Assert (val != 0, ExcEmptyMatrix()); + number s = 0.; for (unsigned int i=0;i::norm2 () const template void FullMatrix::clear () { - fill_n (&val[0], n()*m(), 0); + if (val != 0) + fill_n (&val[0], n()*m(), 0); }; @@ -1073,6 +1159,8 @@ template void FullMatrix::invert (const FullMatrix &M) { + Assert (val != 0, ExcEmptyMatrix()); + Assert (dim_range == dim_image, ExcNotQuadratic()); Assert (dim_range == M.dim_range, ExcDimensionMismatch(dim_range,M.dim_range)); @@ -1230,6 +1318,8 @@ template void FullMatrix::print_formatted (ostream &out, const unsigned int precision) const { + Assert (val != 0, ExcEmptyMatrix()); + out.precision (precision); out.setf (ios::scientific, ios::floatfield); // set output format @@ -1257,7 +1347,9 @@ template void FullMatrix::gauss_jordan() { + Assert (val != 0, ExcEmptyMatrix()); Assert (dim_range == dim_image, ExcNotQuadratic()); + vector p(n()); for (unsigned int i=0; i void FullMatrix::householder(Vector& src) { + Assert (val != 0, ExcEmptyMatrix()); + // m > n, src.n() = m Assert (dim_range <= dim_image, ExcDimensionMismatch(dim_range, dim_image)); Assert (src.size() == dim_image, ExcDimensionMismatch(src.size(), dim_image)); @@ -1365,6 +1459,8 @@ template double FullMatrix::least_squares(Vector& dst, Vector& src) { + Assert (val != 0, ExcEmptyMatrix()); + // m > n, m = src.n, n = dst.n householder(src); -- 2.39.5