Department of Scientific Computing    Institute for Numerical Simulation    University of Bonn

Documentation

Download

Programming References

Bug Reports / Suggestions

FAQ

Authors

           

---

Matrix2D.cc

# include "Matrix2D.hpp"

//extern int debugRefine ;

Matrix2D::Matrix2D() {
  size[0]=size[1]=0 ; 
  a=NULL ;
}

Matrix2D::Matrix2D(const int s0,const int s1) {
  Init(s0,s1) ;
}

void Matrix2D::Init(const int s0,const int s1) {
  size[0]=s0;size[1]=s1 ; 
  if (s0*s1<=0) return ;
  a=(double **)malloc(sizeof(double *)*s0); assert(a);
  for (int i=0; i<s0; i++) { a[i]=(double *)calloc(s1,sizeof(double)); assert(a[i]) ; }
}

Matrix2D::~Matrix2D()
{ 
  if (a) for (int i=0; i<size[0]; i++) free(a[i]) ;
  free(a) ;
  a=NULL  ;
  size[0]=size[1]=0 ;
}

void Matrix2D::Print() const
{
  std::cout <<"size: "<<size[0]<<' '<<size[1]<<' '<<a<<'\n' ;
  for (int i=0; i<size[0]; i++) {
    for (int j=0; j<size[1]; j++) std::cout <<a[i][j]<<' '  ;
    std::cout << '\n' ;
  } 
}

void Matrix2D::Copy(Matrix2D *B) 
{
  assert(a==NULL) ;
  Init(B->size[0],B->size[1]) ;
  
  for (int i=0; i<size[0]; i++)
    for (int j=0; j<size[1]; j++) a[i][j]=B->a[i][j] ;
}

void Matrix2D::Mirror(Matrix2D *B) 
{
  assert(a==NULL) ;
  Init(B->size[0],B->size[1]) ;
  
  for (int i=0; i<size[0]; i++)
    for (int j=0; j<size[1]; j++) a[size[0]-1-i][size[1]-1-j]=B->a[i][j] ;
}


void Matrix2D::Read(FILE *fid)
{
  fscanf(fid,"%d %d",&size[0],&size[1]) ;
  if (size[0]*size[1]<=0) { 
    size[0]=size[1]=0 ;
  }
  Init(size[0],size[1]) ;
  
  for (int i=0;i<size[0];i++)
    for (int j=0;j<size[1];j++) fscanf(fid,"%le",&a[i][j]) ;
} 

/********************************************************************
  Matrix Multiplications, which allow for calculations in place,    *
  i.e. V and W may be the same                                      */

void Matrix2D::MatVec(double *const V,double *const W) const
{
  double s,t[100] ;
  int i ; 
  
  assert(size[1]<100) ;
  for (i=0;i<size[1];i++) t[i]=V[i] ;
  
  for (i=0;i<size[0];i++)
    {s=0; for (int j=0;j<size[1];j++) s+=a[i][j]*t[j] ; W[i]=s ;}
}


void Matrix2D::MatTVec(double *const V,double *const W) const 
{
  double s,t[100] ;
  int i ;
  
  assert(size[1]<100) ;
  for (i=0;i<size[0];i++) t[i]=V[i] ;
  
  for (i=0;i<size[1];i++)
    {s=0; for (int j=0;j<size[0];j++) s+=a[j][i]*t[j] ; W[i]=s ;}
}

void Matrix2D::APlusBMatVec(const double al,const double bl,double *const V,double *const W) const 
{
  double s,t[100] ;
  int i ;
  
  assert(size[1]<100) ;
  for (i=0;i<size[1];i++) t[i]=V[i] ; 
  
  for (i=0;i<size[0];i++)
    {s=0.0; for (int j=0;j<size[1];j++) s+=a[i][j]*t[j] ; W[i] =al*W[i]+bl*s ;}
}

void Matrix2D::APlusBMatTVec(const double al,const double bl,double *const V,double *const W) const 
{
  double s,t[100] ;
  int i ;
  
  assert(size[1]<100) ;
  for (i=0;i<size[0];i++) t[i]=V[i] ;
  
  for (i=0;i<size[1];i++)
    {s=0.0; for (int j=0;j<size[0];j++) s+=a[j][i]*t[j] ; W[i] =al*W[i]+bl*s ;}
}

void Matrix2D::AdaptiveAPlusBMatVec(const double al,const double bl,double *const V,double *const W,int j0,int j1) const 
{
  double s,t[100] ;
  int i,j ;
  
  assert(size[1]<100) ;
  for (i=0;i<size[1];i++) t[i]=V[i] ; 
  
  for(i=0;i<size[0];i++) W[i]*=al ;
  
  for (j=j0;j<=min(j1,size[1]-1);j++) {
    s=t[j]*bl ;
    for (i=0;i<size[0];i++) W[i] += a[i][j]*s ;
  }
}

void Matrix2D::AdaptiveAPlusBMatTVec(const double al,const double bl,double *const V,double *const W,int j0,int j1) const 
{
  double s,t[100] ;
  int i,j ;
  
  assert(size[0]<100) ;
  for (i=0;i<size[0];i++) t[i]=V[i] ; 
  
  for(i=0;i<size[1];i++) W[i]*=al ;
  
  for (j=j0;j<=min(j1,size[0]-1);j++) {
    s=t[j]*bl ;
    for (i=0;i<size[1];i++) W[i] += a[j][i]*s ;
  }
}

---