From: Wolfgang Bangerth Date: Wed, 3 Jun 1998 00:36:10 +0000 (+0000) Subject: Implement a scheme to free unused memory for dVector and dSMatrix. The memory is... X-Git-Tag: v8.0.0~22881 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=39c2152718ad2d3b9fb06a00c0c89b01357bffe7;p=dealii.git Implement a scheme to free unused memory for dVector and dSMatrix. The memory is not freed (as before) if the object is made smaller, but is now freed altogether if the object is assigned the size zero. This functionality is needed if you really want to free all memory, e.g. to let an object sleep for a while but don't want to use those unhandy pointers for half a dozen vectors and matrices. -- Hopefully the new scheme doesn't break too many other things; the previous state of the library is a bit hacked together (my fault also, I admit...). git-svn-id: https://svn.dealii.org/trunk@373 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/dsmatrix.h b/deal.II/lac/include/lac/dsmatrix.h index 064cf26a48..9291ff0833 100644 --- a/deal.II/lac/include/lac/dsmatrix.h +++ b/deal.II/lac/include/lac/dsmatrix.h @@ -28,13 +28,21 @@ template class DoFHandler; CLASS dSMatrixStruct - @author Original version by Roland Becker, Guido Kanschat, - Franz-Theo Suttmeier; lots of enhancements, reorganisation - and documentation by Wolfgang Bangerth + @author Original version by Roland Becker, Guido Kanschat, Franz-Theo Suttmeier; lots of enhancements, reorganisation and documentation by Wolfgang Bangerth */ class dSMatrixStruct { public: + /** + * Initialize the matrix empty, i.e. with + * no memory allocated. This is useful if + * you want such objects as member + * variables in other classes. You can make + * the structure usable by calling the + * #reinit# function. + */ + dSMatrixStruct (); + /** * Initialize a rectangular matrix with * #m# rows and #n# columns, @@ -50,7 +58,8 @@ class dSMatrixStruct * #n# with at most #max_per_row# * nonzero entries per row. */ - dSMatrixStruct (const unsigned int n, const unsigned int max_per_row); + dSMatrixStruct (const unsigned int n, + const unsigned int max_per_row); /** * Destructor. @@ -63,6 +72,14 @@ class dSMatrixStruct * #m# rows and #n# columns, * with at most #max_per_row# * nonzero entries per row. + * + * If #m*n==0# all memory is freed, + * resulting in a total reinitialization + * of the object. If it is nonzero, new + * memory is only allocated if the new + * size extends the old one. This is done + * to save time and to avoid fragmentation + * of the heap. */ void reinit (const unsigned int m, const unsigned int n, @@ -86,6 +103,15 @@ class dSMatrixStruct */ void compress (); + /** + * Return whether the object is empty. It + * is empty if no memory is allocated, + * which is the same as that both + * dimensions are zero. + */ + bool empty () const; + + /** * Return the index of the matrix * element with row number #i# and @@ -289,6 +315,10 @@ class dSMatrixStruct /** * Exception */ + DeclException0 (ExcEmptyObject); + /** + * Exception + */ DeclException0 (ExcInternalError); private: @@ -351,7 +381,7 @@ class dSMatrix * long as #reinit# is not called with a * new sparsity structure. */ - dSMatrix (dSMatrixStruct& c); + dSMatrix (const dSMatrixStruct &sparsity); /** * Destructor. Free all memory, but do not @@ -361,15 +391,31 @@ class dSMatrix ~dSMatrix (); - // + /** + * Reinitialize the object but keep to + * the sparsity pattern previously used. + * This may be necessary if you #reinit#'d + * the sparsity structure and want to + * update the size of the matrix. + * + * Note that memory is only reallocated if + * the new size exceeds the old size. If + * that is not the case, the allocated + * memory is not reduced. However, if the + * sparsity structure is empty (i.e. the + * dimensions are zero), then all memory + * is freed. + */ void reinit (); - // - void reinit (dSMatrixStruct &); /** - * Release all memory and return to a state - * just like after having called the - * default constructor. + * Reinitialize the sparse matrix with the + * given sparsity pattern. The latter tells + * the matrix how many nonzero elements + * there need to be reserved. + * + * Regarding memory allocation, the same + * applies as said above. * * You have to make sure that the lifetime * of the sparsity structure is at least @@ -377,6 +423,15 @@ class dSMatrix * long as #reinit# is not called with a * new sparsity structure. */ + void reinit (const dSMatrixStruct &sparsity); + + /** + * Release all memory and return to a state + * just like after having called the + * default constructor. It also forgets the + * sparsity pattern it was previously tied + * to. + */ void clear (); /** diff --git a/deal.II/lac/include/lac/dvector.h b/deal.II/lac/include/lac/dvector.h index e306f96df4..a1286f3201 100644 --- a/deal.II/lac/include/lac/dvector.h +++ b/deal.II/lac/include/lac/dvector.h @@ -108,6 +108,9 @@ class dVector : public VectorBase /** * Set all entries to zero. Equivalent to * #v = 0#, but more obvious and faster. + * Note that this function does not change + * the size of the vector, unlike the + * STL's #vector<>::clear# function. */ void clear (); @@ -160,10 +163,16 @@ class dVector : public VectorBase /** * Change the dimension of the vector to * #N#. The reserved memory for this vector - * remains unchanged, however, to make + * remains unchanged if possible, to make * things faster, but this may waste some * memory, so take this in the back of your * head. + * However, if #N==0# all memory is freed, + * i.e. if you want to resize the vector + * and release the memory not needed, you + * have to first call #reinit(0)# and then + * #reinit(N)#. This cited behaviour is + * analogous to that of the STL containers. * * On #fast==false#, the vector is filled by * zeros. @@ -175,7 +184,9 @@ class dVector : public VectorBase * vector #V#. The same applies as for * the other #reinit# function. * - * The elements of #V# are not copied. + * The elements of #V# are not copied, i.e. + * this function is the same as calling + * #reinit (V.size(), fast)#. */ void reinit (const dVector& V, const bool fast=false); @@ -457,6 +468,10 @@ class dVector : public VectorBase * Exception */ DeclException0 (ExcOutOfMemory); + /** + * Exception + */ + DeclException0 (ExcEmptyVector); }; diff --git a/deal.II/lac/source/dsmatrix.cc b/deal.II/lac/source/dsmatrix.cc index 95c0525077..f1109d6bf1 100644 --- a/deal.II/lac/source/dsmatrix.cc +++ b/deal.II/lac/source/dsmatrix.cc @@ -11,6 +11,17 @@ +dSMatrixStruct::dSMatrixStruct () : + max_dim(0), + max_vec_len(0), + rowstart(0), + colnums(0) +{ + reinit (0,0,0); +}; + + + dSMatrixStruct::dSMatrixStruct (const unsigned int m, const unsigned int n, const unsigned int max_per_row) : max_dim(0), @@ -37,10 +48,8 @@ dSMatrixStruct::dSMatrixStruct (const unsigned int n, dSMatrixStruct::~dSMatrixStruct () { - if (rowstart != 0) - delete[] rowstart; - if (colnums != 0) - delete[] colnums; + if (rowstart != 0) delete[] rowstart; + if (colnums != 0) delete[] colnums; } @@ -50,28 +59,37 @@ void dSMatrixStruct::reinit (const unsigned int m, const unsigned int n, const unsigned int max_per_row) { - Assert (m>0, ExcInvalidNumber(m)); - Assert (n>0, ExcInvalidNumber(n)); - Assert (max_per_row>0, ExcInvalidNumber(max_per_row)); + Assert ((max_per_row>0) || ((m==0) && (n==0)), ExcInvalidNumber(max_per_row)); rows = m; cols = n; vec_len = m * max_per_row; max_row_len = max_per_row; + if (m*n == 0) + { + if (rowstart) delete[] rowstart; + if (colnums) delete[] colnums; + rowstart = 0; + colnums = 0; + max_vec_len = vec_len = max_dim = rows = cols = 0; + compressed = false; + return; + }; + if (rows > max_dim) { if (rowstart) delete[] rowstart; max_dim = rows; rowstart = new unsigned int[max_dim+1]; - } - + }; + if (vec_len > max_vec_len) { if (colnums) delete[] colnums; max_vec_len = vec_len; colnums = new int[max_vec_len]; - } - + }; + for (unsigned int i=0; i<=rows; i++) rowstart[i] = i * max_per_row; fill_n (&colnums[0], vec_len, -1); @@ -87,6 +105,8 @@ dSMatrixStruct::reinit (const unsigned int m, const unsigned int n, void dSMatrixStruct::compress () { + Assert ((rowstart!=0) && (colnums!=0), ExcEmptyObject()); + if (compressed) return; unsigned int next_free_entry = 0, next_row_start = 0, @@ -160,9 +180,38 @@ dSMatrixStruct::compress () +bool +dSMatrixStruct::empty () const { + // let's try to be on the safe side of + // life by using multiple possibilities in + // the check for emptiness... (sorry for + // this kludge -- emptying matrices and + // freeing memory was not present in the + // original implementation and I don't + // know at how many places I missed + // something in adding it, so I try to + // be cautious. wb) + if ((rowstart==0) || (rows==0) || (cols==0)) + { + Assert (rowstart==0, ExcInternalError()); + Assert (rows==0, ExcInternalError()); + Assert (cols==0, ExcInternalError()); + Assert (colnums==0, ExcInternalError()); + Assert (vec_len==0, ExcInternalError()); + Assert (max_vec_len==0, ExcInternalError()); + Assert (vec_len==0, ExcInternalError()); + + return true; + }; + return false; +}; + + + int dSMatrixStruct::operator () (const unsigned int i, const unsigned int j) const { + Assert ((rowstart!=0) && (colnums!=0), ExcEmptyObject()); Assert (i=0) @@ -278,6 +333,7 @@ dSMatrixStruct::print_gnuplot (ostream &out) const unsigned int dSMatrixStruct::bandwidth () const { + Assert ((rowstart!=0) && (colnums!=0), ExcEmptyObject()); unsigned int b=0; for (unsigned int i=0; icompressed, ExcNotCompressed()); - + Assert (cols->compressed || cols->empty(), ExcNotCompressed()); + + if (cols->empty()) + { + if (val) delete[] val; + val = 0; + max_len = 0; + return; + }; + if (max_lenvec_len) - { - if (val) delete[] val; - val = new double[cols->vec_len]; - max_len = cols->vec_len; - } -// memset(val, 0, sizeof(*val) * cols->vec_len); - for (int i = cols->vec_len-1 ; i>=0 ; i--) val[i] = 0; + { + if (val) delete[] val; + val = new double[cols->vec_len]; + max_len = cols->vec_len; + }; + + if (val) + fill_n (&val[0], cols->vec_len, 0); } void -dSMatrix::reinit (dSMatrixStruct &sparsity) { +dSMatrix::reinit (const dSMatrixStruct &sparsity) { cols = &sparsity; reinit (); }; @@ -377,6 +443,7 @@ dSMatrix::n_nonzero_elements () const { dSMatrix & dSMatrix::copy_from (const dSMatrix &matrix) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert (cols == matrix.cols, ExcDifferentSparsityPatterns()); for (unsigned int i = 0 ; ivec_len; ++i) @@ -390,6 +457,7 @@ dSMatrix::copy_from (const dSMatrix &matrix) { void dSMatrix::add_scaled (const double factor, const dSMatrix &matrix) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert (cols == matrix.cols, ExcDifferentSparsityPatterns()); for (unsigned int i = 0 ; ivec_len; ++i) @@ -402,6 +470,7 @@ void dSMatrix::vmult (dVector& dst, const dVector& src) const { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert(m() == dst.size(), ExcDimensionsDontMatch(m(),dst.size())); Assert(n() == src.size(), ExcDimensionsDontMatch(n(),src.size())); @@ -422,6 +491,7 @@ dSMatrix::vmult (dVector& dst, const dVector& src) const void dSMatrix::Tvmult (dVector& dst, const dVector& src) const { + Assert (val != 0, ExcMatrixNotInitialized()); Assert (cols != 0, ExcMatrixNotInitialized()); Assert(n() == dst.size(), ExcDimensionsDontMatch(n(),dst.size())); Assert(m() == src.size(), ExcDimensionsDontMatch(m(),src.size())); @@ -444,6 +514,7 @@ double dSMatrix::matrix_norm (const dVector& v) const { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert(m() == v.size(), ExcDimensionsDontMatch(m(),v.size())); Assert(n() == v.size(), ExcDimensionsDontMatch(n(),v.size())); @@ -470,6 +541,7 @@ double dSMatrix::residual (dVector& dst, const dVector& u, const dVector& b) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert(m() == dst.size(), ExcDimensionsDontMatch(m(),dst.size())); Assert(m() == b.size(), ExcDimensionsDontMatch(m(),b.size())); Assert(n() == u.size(), ExcDimensionsDontMatch(n(),u.size())); @@ -494,6 +566,7 @@ void dSMatrix::Jacobi_precond (dVector& dst, const dVector& src, const double om) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); const unsigned int n = src.size(); @@ -507,6 +580,7 @@ void dSMatrix::SSOR_precond (dVector& dst, const dVector& src, const double om) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); int p; const unsigned int n = src.size(); @@ -541,6 +615,7 @@ void dSMatrix::SOR_precond (dVector& dst, const dVector& src, const double om) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); dst = src; SOR(dst,om); } @@ -549,6 +624,7 @@ void dSMatrix::SOR (dVector& dst, const double om) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); Assert(n() == m(), ExcDimensionsDontMatch(n(),m())); Assert(m() == dst.size(), ExcDimensionsDontMatch(m(),dst.size())); @@ -569,6 +645,7 @@ void dSMatrix::SSOR (dVector& dst, const double om) { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); int p; const unsigned int n = dst.size(); @@ -616,6 +693,7 @@ const dSMatrixStruct & dSMatrix::get_sparsity_pattern () const { void dSMatrix::print (ostream &out) const { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); for (unsigned int i=0; irows; ++i) for (unsigned int j=cols->rowstart[i]; jrowstart[i+1]; ++j) @@ -626,6 +704,7 @@ void dSMatrix::print (ostream &out) const { void dSMatrix::print_formatted (ostream &out, const unsigned int precision) const { Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (val != 0, ExcMatrixNotInitialized()); out.precision (precision); out.setf (ios::scientific, ios::floatfield); // set output format @@ -642,3 +721,4 @@ void dSMatrix::print_formatted (ostream &out, const unsigned int precision) cons out.setf (0, ios::floatfield); // reset output format }; + diff --git a/deal.II/lac/source/dvector.cc b/deal.II/lac/source/dvector.cc index 4b67a950ed..8543c19d4f 100644 --- a/deal.II/lac/source/dvector.cc +++ b/deal.II/lac/source/dvector.cc @@ -9,7 +9,7 @@ #include -inline double sqr (const double x) { +static inline double sqr (const double x) { return x*x; }; @@ -24,17 +24,10 @@ dVector::dVector () : dVector::dVector (const unsigned int n) : dim(n), - maxdim(n), + maxdim(0), val(0) { - Assert (n>0, ExcInvalidNumber(n)); - - if (n) - { - val = new double[maxdim]; - Assert (val != 0, ExcOutOfMemory()); - clear (); - } + reinit (n, false); } @@ -54,38 +47,33 @@ dVector::dVector (const dVector& v) : -void dVector::reinit (const unsigned int n, const bool fast) -{ - Assert (n>0, ExcInvalidNumber(n)); - +void dVector::reinit (const unsigned int n, const bool fast) { + if (n==0) + { + if (val) delete[] val; + val = 0; + maxdim = dim = 0; + return; + }; + if (n>maxdim) - { - if (val) delete[] val; - val = new double[n]; - Assert (val != 0, ExcOutOfMemory()); - maxdim = n; - } + { + if (val) delete[] val; + val = new double[n]; + Assert (val != 0, ExcOutOfMemory()); + maxdim = n; + }; dim = n; - if (!fast) + if (fast == false) clear (); } -void dVector::reinit (const dVector& v, const bool fast) -{ - const unsigned int n = v.size(); - if (n>maxdim) - { - if (val) delete[] val; - val = new double[n]; - Assert (val != 0, ExcOutOfMemory()); - maxdim = n; - } - dim = n; - if (!fast) - clear (); -} +void dVector::reinit (const dVector& v, const bool fast) { + reinit (v.size(), fast); +}; + @@ -97,12 +85,15 @@ dVector::~dVector () void dVector::clear () { - fill (begin(), end(), 0.); + if (dim>0) + fill (begin(), end(), 0.); } bool dVector::all_zero () const { + Assert (dim!=0, ExcEmptyVector()); + const_iterator p = begin(), e = end(); while (p!=e) @@ -115,6 +106,8 @@ bool dVector::all_zero () const { double dVector::operator * (const dVector& v) const { + Assert (dim!=0, ExcEmptyVector()); + if (&v == this) return norm_sqr(); @@ -149,6 +142,8 @@ double dVector::operator * (const dVector& v) const double dVector::norm_sqr () const { + Assert (dim!=0, ExcEmptyVector()); + double sum0 = 0, sum1 = 0, sum2 = 0, @@ -177,6 +172,8 @@ double dVector::norm_sqr () const double dVector::mean_value () const { + Assert (dim!=0, ExcEmptyVector()); + double sum0 = 0, sum1 = 0, sum2 = 0, @@ -205,6 +202,8 @@ double dVector::mean_value () const double dVector::l1_norm () const { + Assert (dim!=0, ExcEmptyVector()); + double sum0 = 0, sum1 = 0, sum2 = 0, @@ -239,6 +238,8 @@ double dVector::l2_norm () const double dVector::linfty_norm () const { + Assert (dim!=0, ExcEmptyVector()); + double max0=0., max1=0., max2=0., @@ -265,6 +266,8 @@ double dVector::linfty_norm () const { dVector& dVector::operator += (const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); + add (v); return *this; } @@ -273,7 +276,9 @@ dVector& dVector::operator += (const dVector& v) dVector& dVector::operator -= (const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); + iterator i_ptr = begin(), i_end = end(); const_iterator v_ptr = v.begin(); @@ -287,6 +292,8 @@ dVector& dVector::operator -= (const dVector& v) void dVector::add (const double v) { + Assert (dim!=0, ExcEmptyVector()); + iterator i_ptr = begin(), i_end = end(); while (i_ptr!=i_end) @@ -297,7 +304,9 @@ void dVector::add (const double v) void dVector::add (const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); + iterator i_ptr = begin(), i_end = end(); const_iterator v_ptr = v.begin(); @@ -309,7 +318,9 @@ void dVector::add (const dVector& v) void dVector::add (const double a, const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); + iterator i_ptr = begin(), i_end = end(); const_iterator v_ptr = v.begin(); @@ -322,6 +333,7 @@ void dVector::add (const double a, const dVector& v) void dVector::add (const double a, const dVector& v, const double b, const dVector& w) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); Assert (dim == w.dim, ExcDimensionsDontMatch(dim, w.dim)); iterator i_ptr = begin(), @@ -336,6 +348,7 @@ void dVector::add (const double a, const dVector& v, void dVector::sadd (const double x, const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); iterator i_ptr = begin(), i_end = end(); @@ -348,6 +361,7 @@ void dVector::sadd (const double x, const dVector& v) void dVector::sadd (const double x, const double a, const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); iterator i_ptr = begin(), i_end = end(); @@ -361,6 +375,7 @@ void dVector::sadd (const double x, const double a, const dVector& v) void dVector::sadd (const double x, const double a, const dVector& v, const double b, const dVector& w) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); Assert (dim == w.dim, ExcDimensionsDontMatch(dim, w.dim)); iterator i_ptr = begin(), @@ -377,6 +392,7 @@ void dVector::sadd (const double x, const double a, const dVector& v, const double b, const dVector& w, const double c, const dVector& y) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); Assert (dim == w.dim, ExcDimensionsDontMatch(dim, w.dim)); Assert (dim == y.dim, ExcDimensionsDontMatch(dim, y.dim)); @@ -394,6 +410,8 @@ void dVector::sadd (const double x, const double a, void dVector::scale (const double factor) { + Assert (dim!=0, ExcEmptyVector()); + iterator ptr=begin(), eptr=end(); while (ptr!=eptr) *ptr++ *= factor; @@ -404,6 +422,7 @@ void dVector::scale (const double factor) void dVector::equ (const double a, const dVector& u, const double b, const dVector& v) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == u.dim, ExcDimensionsDontMatch(dim, u.dim)); Assert (dim == v.dim, ExcDimensionsDontMatch(dim, v.dim)); iterator i_ptr = begin(), @@ -418,6 +437,7 @@ void dVector::equ (const double a, const dVector& u, void dVector::equ (const double a, const dVector& u) { + Assert (dim!=0, ExcEmptyVector()); Assert (dim == u.dim, ExcDimensionsDontMatch(dim, u.dim)); iterator i_ptr = begin(), i_end = end(); @@ -429,6 +449,7 @@ void dVector::equ (const double a, const dVector& u) void dVector::ratio (const dVector &a, const dVector &b) { + Assert (dim!=0, ExcEmptyVector()); Assert (a.dim == b.dim, ExcDimensionsDontMatch (a.dim, b.dim)); // no need to reinit with zeros, since @@ -446,6 +467,7 @@ void dVector::ratio (const dVector &a, const dVector &b) { dVector& dVector::operator = (const double s) { + Assert (dim!=0, ExcEmptyVector()); fill (begin(), end(), s); return *this; } @@ -457,7 +479,8 @@ dVector& dVector::operator = (const dVector& v) if (v.dim != dim) reinit (v.dim, true); - copy (v.begin(), v.end(), begin()); + if (dim!=0) + copy (v.begin(), v.end(), begin()); return *this; } @@ -466,6 +489,7 @@ dVector& dVector::operator = (const dVector& v) void dVector::cadd (const unsigned int i, const VectorBase& V, const double s, const unsigned int j) { + Assert (dim!=0, ExcEmptyVector()); const dVector& v = (const dVector&) V; Assert (i