smatrixn.cpp

basic mathematic classes used for math programming

//
// Sparse Matrix:	Aug 11, 1998
//

#include "mathclass.h"
//################################################################################################
void
smatrixN::allocate(int N)
{
	n = N;

	nnze = 0;

	rows = new entity[n];
	for (int i = 0; i < n; i++)
	{
		rows[i].id = HEAD_ENTITY;
		rows[i].value = MEANINGLESS;
		rows[i].next = NULL;
	}
}

void
smatrixN::deallocate()
{
	// deallocate memories for every non-zero entities
	entity*	curr;
	entity*	next;
	for (int i = 0; i < n; i++)
	{
		next = rows[i].next;
		while (next)
		{
			curr = next;
			next = curr->next;
			delete	curr;
		}
	}
	delete []	rows;

	rows = NULL; n = nnze = 0;
}

m_real
smatrixN::setValue(int row, int col, m_real value, entity** e)
{
	assert(0 <= row && row < n);
	assert(0 <= col && col < n);

	if (row > col)	{ int tmp = row; row = col; col = tmp; }

	// find A(row,col) entity
	for (entity* p = &rows[row]; p; p = p->next)
	{
		if (p->id == col)
		{
			m_real	old = p->value;
			p->value = value;
			if (e)	*e = p;

			return	old;
		}
		else if (!p->next || p->next->id > col)
		{
			nnze++;

			entity*	tmp = new entity;
			tmp->id = col;
			tmp->value = value;
			tmp->next = p->next;
			p->next = tmp;
			if (e)	*e = tmp;

			return	0;
		}
	}

	return	0;
}

m_real
smatrixN::getValue(int row, int col) const
{
	assert(0 <= row && row < n);
	assert(0 <= col && col < n);

	if (row > col)	{ int tmp = row; row = col; col = tmp; }

	// find A(row,col) entity
	for (entity* p = rows[row].next; p; p = p->next)
		if (p->id == col)	return	p->value;

	// A(row,col) entity is zero
	return	0;
}

m_real
smatrixN::add(int row, int col, m_real increment)
{
	entity*	e;
	m_real	old;
	
	old = setValue(row, col, 0, &e);
	e->value = old + increment;

	return	(old + increment);
}

m_real
smatrixN::sub(int row, int col, m_real decrement)
{
	entity*	e;
	m_real	old;

	old = setValue(row, col, 0, &e);
	e->value = old - decrement;

	return	(old - decrement);
}

ostream&
operator<<(ostream& os, smatrixN& A)
{
	os << A.n << " x " << A.n << " ";
	os << A.nnze << " non-zero sparse matrix " << endl;

	for (int i = 0; i < A.n; i++)
	for (entity* p = A.rows[i].next; p; p = p->next)
		os << " [ " << i << " , " << p->id << " ] = " << p->value << endl;
	os << endl;

	return	os;
}

istream&
operator>>(istream& is, smatrixN& A)
{
	static char	buf[256];

	int	nrow, ncol, nnze;

//	is >> nrow >> "x" >> ncol;
	is >> nrow >> buf >> ncol;
//	is >> nnze >> "non-zero" >> "sparse" >> "matrix";
	is >> nnze >> buf >> "sparse" >> buf;

	assert(nrow == ncol);

	int	n = nrow;

	A.deallocate();
	A.allocate(n);

	int		row, col;
	m_real	value;
	for (int i = 0; i < nnze; i++)
	{
		//is >> "[" >> row >> "," >> col >> "]" >> "=" >> value;
		is >> buf >> row >> buf >> col >> buf >> buf >> value;

		if (fabs(value) > EPS)	A.setValue(row, col, value);
	}

	return	is;
}

//################################################################################################

SparseMatrixN::SparseMatrixN(){

	n = 0;
	nnze = 0;
    rows = NULL;


}

void 
SparseMatrixN::setSize(int n){

	deallocate();
	allocate(n);
}

void
SparseMatrixN::allocate(int N)
{
	
	if(rows!=NULL) deallocate();

	n = N;

	nnze = 0;

	rows = new entity[n];
	for (int i = 0; i < n; i++)
	{
		rows[i].id = HEAD_ENTITY;
		rows[i].value = MEANINGLESS;
		rows[i].next = NULL;
	}
}

int     
SparseMatrixN::getColNum(int row){

  int colnum;
  colnum = 0;
  assert(0 <= row && row< n);

  for (entity* p = rows[row].next; p; p = p->next){
		
		colnum++; 

  }

  return colnum;


}

int  
SparseMatrixN::deleteSmallValEntity(int row,m_real valthres){

	assert(0 <= row && row < n);

	int delnum;
	delnum = 0;

	entity * delitem= NULL;
	entity * pre;
	pre = &rows[row];
	entity* p = rows[row].next;
	while(p){

		if(p->value<valthres){
		
			delitem = p;
			pre->next = p->next;
			p = p->next;

			delnum++;

			delete delitem;

		}else{

			pre = p;
			p = p->next;
		}
		
	}

	return delnum;

}
m_real  
SparseMatrixN::deleteEntity(int row, int col){

	assert(0 <= row && row < n);
	assert(0 <= col && col < n);

	entity * delitem= NULL;
	entity * pre;
	pre = &rows[row];
	for (entity* p = rows[row].next; p; p = p->next){
		
		if (p->id == col){

			delitem = p;
			pre->next = p->next;
			nnze--;
			break;
		}else{

			if(p->id<col){

				pre = p;
			}else{

				break;
			}
		}

	}

	if(delitem != NULL){

		m_real val ;
		val = delitem->value; 
		delete delitem;
		return val;

	}else{

		return 0;
	}
}

m_real  
SparseMatrixN::sumRow(int row){

	assert(0 <= row && row < n);

	m_real sumval = 0;

	for (entity* p = rows[row].next; p; p = p->next){
		
		sumval = sumval + p->value; 

	}

	return sumval;

}
int  
SparseMatrixN::multiplyRow(int row, m_real val){

	assert(0 <= row && row < n);

	int colsum=0;

	for (entity* p = rows[row].next; p; p = p->next){
		
		p->value =  p->value*val ; 
        colsum ++;
	}

	return colsum;
}


void
SparseMatrixN::deallocate()
{
	// deallocate memories for every non-zero entities
	entity*	curr;
	entity*	next;
	for (int i = 0; i < n; i++)
	{
		next = rows[i].next;
		while (next)
		{
			curr = next;
			next = curr->next;
			delete	curr;
		}
	}
	delete []	rows;

	rows = NULL; n = nnze = 0;
}

m_real
SparseMatrixN::setValue(int row, int col, m_real value, entity** e)
{
	assert(0 <= row && row < n);
	assert(0 <= col && col < n);

	//if (row > col)	{ int tmp = row; row = col; col = tmp; }

	// find A(row,col) entity
	for (entity* p = &rows[row]; p; p = p->next)
	{
		if (p->id == col)
		{
			m_real	old = p->value;
			p->value = value;
			if (e)	*e = p;

			return	old;
		}
		else if (!p->next || p->next->id > col)
		{
			nnze++;

			entity*	tmp = new entity;
			tmp->id = col;
			tmp->value = value;
			tmp->next = p->next;
			p->next = tmp;
			if (e)	*e = tmp;

			return	0;
		}
	}

	return	0;
}


bool  
SparseMatrixN::isRowEmpty(int row){

	assert(0<=row && row<n);
	entity* p;
	p =  rows[row].next;

	if(p==NULL) return true;
	else return false;
}

m_real
SparseMatrixN::getValue(int row, int col) const
{
	assert(0 <= row && row < n);
	assert(0 <= col && col < n);

//	if (row > col)	{ int tmp = row; row = col; col = tmp; }

	// find A(row,col) entity
	for (entity* p = rows[row].next; p; p = p->next)
		if (p->id == col)	return	p->value;

	// A(row,col) entity is zero
	return	0;
}

m_real
SparseMatrixN::add(int row, int col, m_real increment)
{
	entity*	e;
	m_real	old;
	
	old = setValue(row, col, 0, &e);
	e->value = old + increment;

	return	(old + increment);
}

m_real
SparseMatrixN::sub(int row, int col, m_real decrement)
{
	entity*	e;
	m_real	old;

	old = setValue(row, col, 0, &e);
	e->value = old - decrement;

	return	(old - decrement);
}

ostream&
operator<<(ostream& os, SparseMatrixN& A)
{
	os << A.n << " x " << A.n << " ";
	os << A.nnze << " non-zero sparse matrix " << endl;

	for (int i = 0; i < A.n; i++)
	for (entity* p = A.rows[i].next; p; p = p->next)
		os << " [ " << i << " , " << p->id << " ] = " << p->value << endl;
	os << endl;

	return	os;
}

istream&
operator>>(istream& is, SparseMatrixN& A)
{
	static char	buf[256];

	int	nrow, ncol, nnze;

//	is >> nrow >> "x" >> ncol;
	is >> nrow >> buf >> ncol;
//	is >> nnze >> "non-zero" >> "sparse" >> "matrix";
	is >> nnze >> buf >> "sparse" >> buf;

	assert(nrow == ncol);

	int	n = nrow;

	A.deallocate();
	A.allocate(n);

	int		row, col;
	m_real	value;
	for (int i = 0; i < nnze; i++)
	{
		//is >> "[" >> row >> "," >> col >> "]" >> "=" >> value;
		is >> buf >> row >> buf >> col >> buf >> buf >> value;

		if (fabs(value) > EPS)	A.setValue(row, col, value);
	}

	return	is;
}