matrix.c

Embedded magazine source code. this is the source code in the year 1989

 *		each chain pointer.
 *	    If indices don't match, copy element with lower index into
 *		resultant matrix, and update it's chain pointer.  Leave other
 *		chain as is.
 *	Return resultant matrix.
 *
 ***************************************************************************/
MATRIX *add_matrices( MATRIX *pMatrix1, MATRIX *pMatrix2 )
{
    double dVal;
    USHORT usDimRows;
    USHORT usDimCols;
    USHORT usCol;
    MATRIX *pResult;
    MATRIX_ELEM *pChainM1;
    MATRIX_ELEM *pChainM2;
    BYTE bType;

    if( pMatrix1 == NULL ) {
	pResult = copy_matrix( pMatrix2, ROW );
	return( pResult );
    } /* pMatrix1 == NULL */
    if( pMatrix2 == NULL ) {
	pResult = copy_matrix( pMatrix1, ROW );
	return( pResult );
    } /* pMatrix2 == NULL */

    usDimRows = pMatrix1->bRow;
    usDimCols = pMatrix1->bCol;
    if( ( usDimRows != pMatrix2->bRow ) && ( usDimCols != pMatrix2->bCol ) )
	return( NULL );

    bType = pMatrix1->bType & LOW;
    if( bType != ( pMatrix2->bType & LOW ) )
	change_rep( &pMatrix2 );

    if( ( pResult = dimension( usDimRows, usDimCols, bType ) ) == NULL )
	return( NULL );

    while( 0 < usDimRows-- ) {

	pChainM1 = ( *pMatrix1->pFirst )[ usDimRows ];
	pChainM2 = ( *pMatrix2->pFirst )[ usDimRows ];
	while( ( pChainM1 != NULL ) || ( pChainM2 != NULL ) ) {

	    if( pChainM1->index == pChainM2->index ) {
		dVal = pChainM1->elem;
		dVal += pChainM2->elem;
		usCol = pChainM1->index;
		pChainM1 = pChainM1->pNext;
		pChainM2 = pChainM2->pNext;
	    } else if( pChainM1->index < pChainM2->index ) {
		dVal = pChainM1->elem;
		usCol = pChainM1->index;
		pChainM1 = pChainM1->pNext;
	    } else if( pChainM1->index > pChainM2->index ) {
		dVal = pChainM2->elem;
		usCol = pChainM2->index;
		pChainM2 = pChainM2->pNext;
	    } /* else */
	    if( bType == ROW )
		in_matrix( pResult, usDimRows, usCol, dVal );
	    else
		in_matrix( pResult, usCol, usDimRows, dVal );

	} /* while( ( pChainM1 != NULL ) && ( pChainM2 != NULL ) ) */

    } /* while( 0 < usDimRows-- ) */

    return( pResult );

} /* add_matrices */


/*** sub_matrices() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Subtract the second matrix passed from the first.
 *
 *   INPUT :
 *	result_matrix = sub_matrices( pMatrix1, pMatrix2 )
 *
 *	    pMatrix1 : Pointer to matrix.
 *	    pMatrix2 : Pointer to matrix to be subtracted from matrix 1.
 *
 *   OUTPUT :
 *	MATRIX * : Pointer to resultant matrix = pMatrix1 - pMatrix2.
 *	Error codes are as possible in const_mult_matrix and add_matrices.
 *
 *   OPERATIONS :
 *	Multiply the second matrix by -1.0.
 *	Let add_matrices do the work.
 *	Free the memory taken by intermediate matrix.
 *
 ***************************************************************************/
MATRIX *sub_matrices( MATRIX *pMatrix1, MATRIX *pMatrix2 )
{
    MATRIX *pTemp1;
    MATRIX *pTemp2;

    pTemp1 = const_mult_matrix( pMatrix2, -1.0 );
    pTemp2 = add_matrices( pMatrix1, pTemp1 );
    remove_matrix( pTemp1 );
    return( pTemp2 );

} /* sub_matrices */


/*** const_mult_matrix() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Multiply the given matrix by the given constant.
 *
 *   INPUT :
 *	result_matrix = const_mult_matrix( pMatrix, dConstant )
 *
 *   OUTPUT :
 *	MATRIX * : Resultant matrix = dConstant * pMatrix
 *
 *   OPERATIONS :
 *	Dimension the resultant matrix.
 *	Determine the number of chains in the list by the matrix type.
 *	Scan each chain.  For each element found, multiply it by the
 *	    constant, and store it in the resultant matrix.
 *	return resultant matrix.
 *
 ***************************************************************************/
MATRIX *const_mult_matrix( MATRIX *pMatrix, double con )
{
    MATRIX *pResult;
    MATRIX_ELEM *pChain;
    USHORT i;
    USHORT usDimRows;
    USHORT usDimCols;
    BYTE bType;
    BYTE bNumPointers;
    double dTemp;

    usDimRows = pMatrix->bRow;
    usDimCols = pMatrix->bCol;
    bType = pMatrix->bType & LOW;

    if( ( pResult = dimension( usDimRows, usDimCols, bType ) ) == NULL )
	return( NULL );

    switch( bType ) {
	case FULL:
	case ROW :
		    bNumPointers = usDimRows;
		    break;

	case COL :
		    bNumPointers = usDimCols;
		    break;

    } /* switch( type )*/


    for( i = 0; i < bNumPointers; i++ ) {

	 pChain = ( *pMatrix->pFirst )[ i ];
	 while( pChain != NULL ) {

		dTemp = pChain->elem * con;
		if( bType == ROW )
		    in_matrix( pResult, i, pChain->index, dTemp );
		else
		    in_matrix( pResult, pChain->index, i, dTemp );
		pChain = pChain->pNext;

	} /* while */

    } /* for */

    return( pResult );

} /* const_mult_matrix */


/*** mult_matrix() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Multiply the two matrices.  In general, since matrix multiplication
 *   is NOT commutative, the order of multiplaction is important.  It is
 *   pMatrix1 * pMatrix2.
 *
 *   INPUT :
 *	result_matrix = mult_matrices( pMatrix1, pMatrix2 )
 *
 *	    pMatrix1 : Pointer to first matrix.
 *	    pMatrix2 : Pointer to second matrix.
 *
 *   OUTPUT :
 *	MATRIX * : Pointer to resultant matrix = pMatrix1 * pMatrix2.
 *	NULL : Error condition.
 *	       Dimensions not proper for multiplication.
 *	       Out of memory.
 *
 *   OPERATIONS :
 *	Check dimension for proper multimplication.  If pMatrix1 is an
 *	    n x m matrix and pMatrix2 is a k x l matrix, than for multiplcation
 *	    m must equal k and the resultant matrix will have dimensions n x l.
 *	Ensure that matrices are of desired type.  pMatrix1 should be a ROW
 *	    matrix, while pMatrix2 should be a column matrix.  This will allow
 *	    the calling of mult_vector to calculate each element of the result.
 *	Dimension the resultant matrix.  Also dimension dummy matrices to be
 *	    used as a row and column vector.
 *	For each row in pMatrix1:
 *	    multiply current row by each column of pMatrix2.  Store this
 *		value in pResult in location ( row, column ).  This is
 *		occomplished by using the dummy matrices defined above, and
 *		assigning the current row to the row matrix, and the current
 *		column to the column matrix, and calling mult_vector.
 *	Free the memory associated with the tempory matrices.
 *	Return the resultant matrix.
 *
 ***************************************************************************/
MATRIX *mult_matrices( MATRIX *pMatrix1, MATRIX *pMatrix2 )
{
    MATRIX *pResult;
    MATRIX *pRowMat;
    MATRIX *pColMat;
    USHORT usRows;
    USHORT usCols;
    USHORT usCurCol;
    double dTemp;
    BYTE bType;

    if( pMatrix1->bCol != pMatrix2->bRow )
	return( NULL );

    bType = pMatrix1->bType & LOW;
    if( bType != ROW )
	change_rep( &pMatrix1 );

    bType = pMatrix2->bType & LOW;
    if( bType != COL )
	change_rep( &pMatrix2 );

    usRows = pMatrix1->bRow;
    usCols = pMatrix2->bCol;
    if( (pResult = dimension( usRows, usCols, ROW ) ) == NULL )
	return( NULL );
    if( (pRowMat = dimension( 1, usCols, ROW ) ) == NULL )
	return( NULL );
    if( (pColMat = dimension( usRows, 1, COL ) ) == NULL )
	return( NULL );

    while( usRows-- > 0 ) {

	usCurCol = usCols;
	while( usCurCol-- > 0 ) {

	    ( *pRowMat->pFirst )[ 0 ] = ( *pMatrix1->pFirst )[ usRows ];
	    ( *pColMat->pFirst )[ 0 ] = ( *pMatrix2->pFirst )[ usCurCol ];
	    dTemp = mult_vector( pRowMat, pColMat );
	    in_matrix( pResult, usRows, usCurCol, dTemp );

	} /* while( usCurCol-- > 0 ) */

    } /* while( usRows-- > 0 ) */

    free( pRowMat );
    free( pColMat );

    return( pResult );

} /* mult_matrices */


/*** mult_vector() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Multiply the given row vector by the given column vector.
 *
 *   INPUT :
 *	dResult = mult_vector( pRowVect, pColVect )
 *
 *	    pRowVect : Pointer to a ROW type matrix with only one row.
 *	    pColVect : Pointer to a COL type matrix with only one column.
 *
 *   OUTPUT :
 *	dResult : Result of vector multiplication.  This is often called
 *		  a dot product, or inner product.
 *	Error conditions should be reported using errorno global.
 *
 *   OPERATIONS :
 *	Varify matrix types, and that matrix dimensions match.
 *	Initialize result to zero.
 *	Scan the row matrix, for each element, if there is a corresponding
 *	    element in the column vector, multiply them and add to result.
 *	    If there isn't a corresponding element, proceed to the next
 *	    row element.  This is because corresponding column element must
 *	    be zero, and 0 = 0 * any_number.
 *	Return result.
 *
 ***************************************************************************/
double mult_vector( MATRIX *pRowVect, MATRIX *pColVect )
{
    MATRIX_ELEM *pNextRow;
    MATRIX_ELEM *pNextCol;
    double dResult = 0.0;
    USHORT usRow;
    USHORT usCol;
    BYTE bType;

    bType = pRowVect->bType;
    if( ( bType & LOW ) != ROW ) {
	printf( "Error multiplying vectors\n" );
	return( 0.0 );
    }
    bType = pColVect->bType;
    if( ( bType & LOW ) != COL ) {
	printf( "Error multiplying vectors\n" );
	return( 0.0 );
    }

    usRow = pRowVect->bRow;
    if( usRow != pColVect->bCol ) {

	printf( "Error multiplying vectors\n" );
	return( 0.0 );

    } /* if( usRow != pColVect->bCol ) */

    pNextRow = ( *pRowVect->pFirst )[0];
    usRow = pNextRow->index;
    pNextCol = ( *pColVect->pFirst )[0];
    usCol = pNextCol->index;
    while( pNextRow != NULL ) {

	if( usRow == usCol ) {

	    dResult += pNextCol->elem * pNextRow->elem;
	    pNextCol = pNextCol->pNext;
	    usCol = pNextCol->index;
	    pNextRow = pNextRow->pNext;
	    usRow = pNextRow->index;

	} else if( usRow < usCol ) {

	    pNextRow = pNextRow->pNext;
	    usRow = pNextRow->index;

	} else if( usRow > usCol ) {

	    pNextCol = pNextCol->pNext;
	    usCol = pNextCol->index;

	} /* else if( usRow > usCol ) */

    } /* while( pNextRow != NULL ) */

    return( dResult );

} /* mult_vector */


/*** change_rep() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Change the given matrix representation from COL to ROW or vice versa.
 *
 *   INPUT :
 *	ret_condition = change_rep( &pMatrix )
 *
 *	    pMatrix : Pointer to a pointer to a matrix.
 *
 *   OUTPUT :
 *	0 : Representation was changed.
 *	1 : Error, Out of memory.
 *
 *   OPERATIONS :
 *	From the matrices current type, determine the new type, and the
 *	    number of chains in the current matrix.
 *	Dimension the resultant matrix, with the same dimension as the old
 *	    matrix, but with the new type.
 *	For each chain in the old matrix :
 *	    Scan the chain.  For each element, make a new element and insert
 *		it into the new matrix.
 *	Free memory occupied by the old matrix.
 *	Return the pResult matrix.
 *
 ***************************************************************************/
int change_rep( MATRIX **mat )
{
    MATRIX *pResult;
    MATRIX *pMat;
    MATRIX_ELEM *pBeg;
    MATRIX_ELEM *pOld;
    MATRIX_ELEM *pNew;
    MATRIX_ELEM *pTemp;
    BYTE bRow;
    BYTE bCol;
    BYTE bType;
    BYTE bNewType;
    BYTE bLinks;
    USHORT i;

    pMat = *mat;
    bRow = pMat->bRow;
    bCol = pMat->bCol;
    bType = pMat->bType & LOW;

    if( bType == ROW ) {
	bNewType = pMat->bType & HIGH;
	bNewType += COL;
	bLinks = bRow;
    } else {
	bNewType = pMat->bType & HIGH;
	bNewType += ROW;
	bLinks = bCol;
    }
    if( ( pResult = dimension( bRow, bCol, bNewType ) ) == NULL )
	return( 1 );

    for( i = 0; i < bLinks; i++ ) {
	pBeg = ( *pMat->pFirst )[ i ];
	while( pBeg != NULL ) {

	    pNew = make_elem( i, pBeg->elem );
	    pTemp = ( *pResult->pFirst )[ pBeg->index ];
	    if( pTemp == NULL )
		( *pResult->pFirst )[ pBeg->index ] = pNew;
	    else {
		pOld = find_location( pTemp, i );
		pOld->pNext = pNew;
	    }
	    pBeg = pBeg->pNext;

	} /* while( pBeg != NULL ) */

    } /* for( i = 0; i < usRow; i++ ) */

    remove_matrix( pMat );
    *mat = pResult;
    return( 0 );

} /* change_rep */


/*** remove_matrix() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Remove ALL storage associated with a given matrix.
 *
 *   INPUT :
 *	remove_matrix( pMatrix )
 *
 *	    pMatrix : Pointer to matrix to be removed from system.
 *
 *   OUTPUT :
 *	None.
 *
 *   OPERATIONS :
 *	Use the matrix type to determine the number of chains in a matrix.
 *	Use the pFirst array to get the starting location of each chain, and
 *	    call remove_chain to remove it.
 *	Free the memory associated with the matrix body itself.
 *
 ***************************************************************************/
void remove_matrix( MATRIX *mat )
{
    MATRIX_ELEM *pNext;
    USHORT i;
    BYTE bLength;
    BYTE bType;

    bType = mat->bType & LOW;
    switch( bType ) {
	case ROW :
		    bLength = mat->bRow;
		    break;
	case COL :
		    bLength = mat->bCol;
		    break;
	case FULL :
		    bLength = mat->bCol + mat->bRow;
		    break;
    } /* switch( type ) */

    for( i = 0; i < bLength; i++ ) {

	pNext = ( *mat->pFirst )[i];
	remove_chain( pNext );

    } /* for( i = 0; i < usLength; i++ ) */

    free( mat );

} /* remove_matrix */


/*** remove_chain() ***********************************************************
 *
 *   Mike Foley
 *   January 21, 1989
 *
 *   FUNCTIONAL SUMMARY :
 *	Free ALL memory associated with the given chain.
 *
 *   INPUT :
 *	remove_chain( pChain )
 *
 *	    pChain : Pointer to chain to be freed.
 *
 *   OUTPUT :
 *	None.
 *
 *   OPERATIONS :
 *	Until the chain is depleated :
 *	    Save a pointer to the next element in the chain, for when the
 *		current element is freed, it is lost.
 *	    Free the current element in the chain.
 *	    Update the chain pointer to the stored next value.
 *
 ***************************************************************************/
void remove_chain( MATRIX_ELEM *pChain )
{
    MATRIX_ELEM *pTemp;

    while( pChain != NULL ) {

	pTemp = pChain->pNext;
	free( pChain );
	pChain = pTemp;