double* val;
int max_len;
public:
- //
- void reinit();
- //
- void reinit (dSMatrixStruct &);
+
+ /**
+ * Constructor; initializes the matrix to be
+ * empty, without any structure, i.e. the
+ * matrix is not usable at all. This constructor
+ * is therefore only useful for matrices which
+ * are members of a class. You have to initialize
+ * the matrix before usage with
+ * #reinit(dSMatrixStruct)#.
+ */
+ dSMatrix ();
+
//
dSMatrix(dSMatrixStruct& c)
: cols(&c), val(0), max_len(0)
delete[] val;
}
+
+ //
+ void reinit();
+ //
+ void reinit (dSMatrixStruct &);
+
//
int m() const { return cols->rows; }
//
* Exception
*/
DeclException0 (ExcNotCompressed);
+ /**
+ * Exception
+ */
+ DeclException0 (ExcMatrixNotInitialized);
+ /**
+ * Exception
+ */
+ DeclException2 (ExcDimensionsDontMatch,
+ int, int,
+ << "The dimensions " << arg1 << " and " << arg2
+ << " do not match properly.");
};
#endif
}
-/*----------------- from sort.h -------------------------*/
+/*----------------- end: from sort.h -------------------------*/
+
+
+/*-------------------------------------------------------------------------*/
+
+
+dSMatrix::dSMatrix () :
+ cols(0),
+ val(0),
+ max_len(0) {};
+
+
void
dSMatrix::reinit()
{
+ Assert (cols != 0, ExcMatrixNotInitialized());
Assert (cols->compressed, ExcNotCompressed());
if(max_len<cols->vec_len)
void
dSMatrix::vmult(dVector& dst,const dVector& src)
{
- // Assert(m() = dst.n())
- // Assert(n() = src.n())
+ Assert (cols != 0, ExcMatrixNotInitialized());
+ Assert(m() == dst.n(), ExcDimensionsDontMatch(m(),dst.n()));
+ Assert(n() == src.n(), ExcDimensionsDontMatch(n(),src.n()));
for (int i=0;i<m();i++)
{
void
dSMatrix::Tvmult(dVector& dst,const dVector& src)
{
- // Assert(n() = dst.n())
- // Assert(m() = src.n())
+ Assert (cols != 0, ExcMatrixNotInitialized());
+ Assert(n() == dst.n(), ExcDimensionsDontMatch(n(),dst.n()));
+ Assert(m() == src.n(), ExcDimensionsDontMatch(m(),src.n()));
int i;
double
dSMatrix::residual(dVector& dst,const dVector& u,const dVector& b)
{
- // Assert(m() = dst.n())
- // Assert(n() = src.n())
+ Assert (cols != 0, ExcMatrixNotInitialized());
+ Assert(m() == dst.n(), ExcDimensionsDontMatch(m(),dst.n()));
+ Assert(m() == b.n(), ExcDimensionsDontMatch(m(),b.n()));
+ Assert(n() == u.n(), ExcDimensionsDontMatch(n(),u.n()));
double s,norm=0.;
void
dSMatrix::Jacobi_precond(dVector& dst,const dVector& src,double om)
{
+ Assert (cols != 0, ExcMatrixNotInitialized());
+
int n = src.n();
for (int i=0;i<n;++i)
void
dSMatrix::SSOR_precond(dVector& dst,const dVector& src,double om)
{
+ Assert (cols != 0, ExcMatrixNotInitialized());
+
int p,n = src.n();
int i,j;
void
dSMatrix::SOR_precond(dVector& dst,const dVector& src,double om)
{
+ Assert (cols != 0, ExcMatrixNotInitialized());
dst = src;
SOR(dst,om);
}
void
dSMatrix::SOR(dVector& dst, double om)
{
- // Assert(m()==n())
- // Assert(m()==dst.n())
+ Assert (cols != 0, ExcMatrixNotInitialized());
+ Assert(n() == m(), ExcDimensionsDontMatch(n(),m()));
+ Assert(m() == dst.n(), ExcDimensionsDontMatch(m(),dst.n()));
for (int i=0;i<m();i++)
{
void
dSMatrix::SSOR(dVector& dst, double om)
{
+ Assert (cols != 0, ExcMatrixNotInitialized());
+
int p,n = dst.n();
int i,j;
double s;
void dSMatrix::print (ostream &out) const {
+ Assert (cols != 0, ExcMatrixNotInitialized());
+
for (int i=0; i<cols->rows; ++i)
for (int j=cols->rowstart[i]; j<cols->rowstart[i+1]; ++j)
out << "(" << i << "," << cols->colnums[j] << ") " << val[j] << endl;