From 7b8894561fa5e185bf38c1771b51981be073790a Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 9 Mar 1998 18:32:33 +0000 Subject: [PATCH] One more reinit for dSMatrix git-svn-id: https://svn.dealii.org/trunk@51 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/dsmatrix.h | 93 ++++++++++++++++-------------- deal.II/lac/source/dsmatrix.cc | 92 +++++++++++++++++------------ 2 files changed, 104 insertions(+), 81 deletions(-) diff --git a/deal.II/lac/include/lac/dsmatrix.h b/deal.II/lac/include/lac/dsmatrix.h index 4588e46af1..60db825da9 100644 --- a/deal.II/lac/include/lac/dsmatrix.h +++ b/deal.II/lac/include/lac/dsmatrix.h @@ -110,54 +110,61 @@ CLASS */ class dSMatrix { - dSMatrixStruct& cols; - double* val; - int max_len; - public: - ////////// - void reinit(); - ////////// - dSMatrix(dSMatrixStruct& c) - : cols(c), val(0), max_len(0) - { - reinit(); - } - ////////// - ~dSMatrix() - { - delete[] val; - } + dSMatrixStruct * cols; + double* val; + int max_len; + public: + // + void reinit(); + // + void reinit (dSMatrixStruct &); + // + dSMatrix(dSMatrixStruct& c) + : cols(&c), val(0), max_len(0) + { + reinit(); + } + // + ~dSMatrix() + { + delete[] val; + } - ////////// - int m() const { return cols.rows; } - ////////// - int n() const { return cols.cols; } + // + int m() const { return cols->rows; } + // + int n() const { return cols->cols; } - ////////// - void set(int i,int j,double value) { val[cols(i,j)] = value; } - ////////// - void add(int i,int j,double value) { val[cols(i,j)]+= value; } + // + void set(int i,int j,double value) { val[cols->operator()(i,j)] = value; } + // + void add(int i,int j,double value) { val[cols->operator()(i,j)]+= value; } - ////////// - void vmult (dVector& dst,const dVector& src); - ////////// - void Tvmult(dVector& dst,const dVector& src); + // + void vmult (dVector& dst,const dVector& src); + // + void Tvmult(dVector& dst,const dVector& src); - ////////// - double residual (dVector& dst,const dVector& x,const dVector& b); - ////////// - void Jacobi_precond(dVector& dst,const dVector& src,double om = 1.); - ////////// - void SSOR_precond(dVector& dst,const dVector& src,double om = 1.); - ////////// - void SOR_precond(dVector& dst,const dVector& src,double om = 1.); - ////////// - void SSOR(dVector& dst,double om = 1.); - ////////// - void SOR(dVector& dst,double om = 1.); - ////////// - void precondition(dVector& dst,const dVector& src) { dst=src; } + // + double residual (dVector& dst,const dVector& x,const dVector& b); + // + void Jacobi_precond(dVector& dst,const dVector& src,double om = 1.); + // + void SSOR_precond(dVector& dst,const dVector& src,double om = 1.); + // + void SOR_precond(dVector& dst,const dVector& src,double om = 1.); + // + void SSOR(dVector& dst,double om = 1.); + // + void SOR(dVector& dst,double om = 1.); + // + void precondition(dVector& dst,const dVector& src) { dst=src; } void print (ostream &) const; + + /** + * Exception + */ + DeclException0 (ExcNotCompressed); }; #endif diff --git a/deal.II/lac/source/dsmatrix.cc b/deal.II/lac/source/dsmatrix.cc index b05430401c..d2e344f3f3 100644 --- a/deal.II/lac/source/dsmatrix.cc +++ b/deal.II/lac/source/dsmatrix.cc @@ -255,18 +255,6 @@ dSMatrixStruct::add_matrix(const iVector& rows, const iVector& cols) add(rows(i), cols(j)); } -void -dSMatrix::reinit() -{ - if(max_len=0 ; i--) val[i] = 0; -} void dSMatrixStruct::print_gnuplot (ostream &out) const @@ -296,6 +284,34 @@ dSMatrixStruct::bandwidth () const return b; } + + +void +dSMatrix::reinit() +{ + Assert (cols->compressed, ExcNotCompressed()); + + 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; +} + + + +void +dSMatrix::reinit (dSMatrixStruct &sparsity) { + cols = &sparsity; + reinit (); +}; + + + + void dSMatrix::vmult(dVector& dst,const dVector& src) { @@ -305,9 +321,9 @@ dSMatrix::vmult(dVector& dst,const dVector& src) for (int i=0;irowstart[i]; j < cols->rowstart[i+1] ;j++) { - int p = cols.colnums[j]; + int p = cols->colnums[j]; s += val[j] * src(p); } dst(i) = s; @@ -326,9 +342,9 @@ dSMatrix::Tvmult(dVector& dst,const dVector& src) for (i=0;irowstart[i]; jrowstart[i+1] ;j++) { - int p = cols.colnums[j]; + int p = cols->colnums[j]; dst(p) += val[j] * src(i); } } @@ -345,9 +361,9 @@ dSMatrix::residual(dVector& dst,const dVector& u,const dVector& b) for (int i=0;irowstart[i]; jrowstart[i+1] ;j++) { - int p = cols.colnums[j]; + int p = cols->colnums[j]; s -= val[j] * u(p); } dst(i) = s; @@ -363,7 +379,7 @@ dSMatrix::Jacobi_precond(dVector& dst,const dVector& src,double om) for (int i=0;irowstart[i]]; } } @@ -376,23 +392,23 @@ dSMatrix::SSOR_precond(dVector& dst,const dVector& src,double om) for (i=0;irowstart[i]; jrowstart[i+1] ;j++) { - p = cols.colnums[j]; + p = cols->colnums[j]; if (prowstart[i]]; } - for (i=0;irowstart[i]]; for (i=n-1;i>=0;i--) { - for (j=cols.rowstart[i];jrowstart[i];jrowstart[i+1];j++) { - p = cols.colnums[j]; + p = cols->colnums[j]; if (p>i) dst(i) -= om* val[j] * dst(p); } - dst(i) /= val[cols.rowstart[i]]; + dst(i) /= val[cols->rowstart[i]]; } } @@ -412,13 +428,13 @@ dSMatrix::SOR(dVector& dst, double om) for (int i=0;irowstart[i]; jrowstart[i+1] ;j++) { - int p = cols.colnums[j]; + int p = cols->colnums[j]; if (prowstart[i]]; } } @@ -432,36 +448,36 @@ dSMatrix::SSOR(dVector& dst, double om) for (i=0;irowstart[i]; jrowstart[i+1] ;j++) { - p = cols.colnums[j]; + p = cols->colnums[j]; if (p>=0) { if (i>j) s += val[j] * dst(p); } } dst(i) -= s * om; - dst(i) /= val[cols.rowstart[i]]; + dst(i) /= val[cols->rowstart[i]]; } for (i=n-1;i>=0;i--) { s = 0.; - for (j=cols.rowstart[i]; jrowstart[i]; jrowstart[i+1] ;j++) { - p = cols.colnums[j]; + p = cols->colnums[j]; if (p>=0) { if (irowstart[i]]; } } void dSMatrix::print (ostream &out) const { - for (int i=0; irows; ++i) + for (int j=cols->rowstart[i]; jrowstart[i+1]; ++j) + out << "(" << i << "," << cols->colnums[j] << ") " << val[j] << endl; }; -- 2.39.5