sparse-matrix.c

Go to the documentation of this file.

/* sparse-matrix.c

*****
Copyright 1994-1997, Brown University
Copyright 1998, 1999, Anthony R. Cassandra

All Rights Reserved

Permission to use, copy, modify, and distribute this software and its
documentation for any purpose other than its incorporation into a
commercial product is hereby granted without fee, provided that the
above copyright notice appear in all copies and that both that
copyright notice and this permission notice appear in supporting
documentation.

ANTHONY CASSANDRA DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR ANY
PARTICULAR PURPOSE.  IN NO EVENT SHALL ANTHONY CASSANDRA BE LIABLE FOR
ANY SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
*****

This module will implementall the routines that pertain to creating and
manipulating a sparse representation of a matrix.

*/

#include <stdio.h>
#include <stdlib.h>
#include <assert.h>

#include "sparse-matrix.h"

/**********************************************************************/
/********************  Routines for row linked lists  *****************/
/**********************************************************************/
I_Matrix_Row_Node newRowNode( int col, REAL_VALUE value ) {
        /*
        Just allocates memory and sets the fields.
        */
        I_Matrix_Row_Node new_node;

        new_node = (I_Matrix_Row_Node ) malloc( sizeof( *new_node ));
        checkAllocatedPointer(new_node);

        new_node->column = col;
        new_node->value = value;
        new_node->next = NULL;

        return( new_node );

}  /* newRowNode */
/**********************************************************************/
void destroyRow( I_Matrix_Row_Node row ) {
        I_Matrix_Row_Node temp;

        while( row != NULL ) {
                temp = row;
                row = row->next;
                free( temp );
        } /* while */
}  /* destroyRow */
/**********************************************************************/
I_Matrix_Row_Node removeRowNode( I_Matrix_Row_Node row, 
                                                                int col, int *count ) {
                                                                        /*
                                                                        Remove the node with column == col, if it exists.  It decrements
                                                                        the count if it is removed.  If there is no entry for this in the
                                                                        list then we simply igore it and return the original list.
                                                                        */
                                                                        I_Matrix_Row_Node temp_node, cur_node, trail_node;

                                                                        /* Case if list is empty */
                                                                        if( row == NULL )
                                                                                return( NULL );

                                                                        /* Case of removing first item in the list */
                                                                        if( row->column == col ) {
                                                                                temp_node = row;
                                                                                row = row->next;
                                                                                free( temp_node );
                                                                                (*count)--;
                                                                                return( row );
                                                                        }  /* if we should remove first node */

                                                                        trail_node = row;
                                                                        cur_node = row->next;
                                                                        while( cur_node != NULL ) {

                                                                                /* Then bingo!  We found the one to remove */
                                                                                if( cur_node->column == col ) {
                                                                                        trail_node->next = cur_node->next;
                                                                                        free( cur_node );
                                                                                        (*count)--;
                                                                                        return( row );
                                                                                }  /* if we found the one to remove */

                                                                                trail_node = cur_node;
                                                                                cur_node = cur_node->next;
                                                                        }  /* while */

                                                                        /* Otherwise we did not find an entry with this column, so just */
                                                                        /* return the list and do nothing. */
                                                                        return( row );

}  /* removeRowNode */
/**********************************************************************/
I_Matrix_Row_Node addEntryToRow( I_Matrix_Row_Node row, 
                                                                int col, REAL_VALUE value,
                                                                int *count, int accumulate ) {
                                                                        /*
                                                                        Will step through the linked list until either the column is found
                                                                        or until we see a strictly smaller column number.  It will either
                                                                        replace the current value (if the column already exists) or insert a
                                                                        new node into the list just before the first column number that is
                                                                        highest.  This way, the list always remains sorted by column
                                                                        number.  The pointer to the int 'count' will be incremented only if
                                                                        we add it to the list, not when we simply replace a value.  *count 
                                                                        should be the length of the list sent in, but this routine does not
                                                                        check to make sure it is.

                                                                        The accumulate flag will specify whether we should replace the 
                                                                        current value (if any) or add the value to the existing value.
                                                                        This can be used to build up probabilities in a matrix piecemeal.
                                                                        */
                                                                        I_Matrix_Row_Node new_node, cur_node, trail_node;
                                                                        trail_node = 0; // avoid warning

                                                                        /* Case if we attempt to add a zero entry.  We need to see if there */
                                                                        /* is already a non-zero entry for it in the list, and if so remove */
                                                                        /* it.  Otherwise we simply ignore this entry request, since we are */
                                                                        /* using a sparse representation. */
                                                                        if( IS_ZERO( value ) ) {
                                                                                if ( accumulate )
                                                                                        return ( row );
                                                                                else
                                                                                        return( removeRowNode( row, col, count ));
                                                                        }

                                                                        /* Case if the list is empty */
                                                                        if( row == NULL ) {
                                                                                new_node = newRowNode( col, value );
                                                                                (*count)++;
                                                                                return( new_node );
                                                                        }  /* if row empty */

                                                                        /* Case if it should be inserted into front of list */
                                                                        if( row->column > col ) {
                                                                                new_node = newRowNode( col, value );
                                                                                new_node->next = row;
                                                                                (*count)++;
                                                                                return( new_node );
                                                                        }  /* if should be inserted in front of list */

                                                                        cur_node = row;
                                                                        while( cur_node != NULL ) {

                                                                                /* Case if we should simply replace (or accumulate) 
                                                                                the current value */
                                                                                if( cur_node->column == col ) {

                                                                                        if( accumulate )
                                                                                                cur_node->value += value;
                                                                                        else
                                                                                                cur_node->value = value;

                                                                                        return( row );
                                                                                } /* if we should just replace value */

                                                                                /* Case if this is the point we should insert it (i.e. between */
                                                                                /* trail_node and cur_node */
                                                                                if( cur_node->column > col ) {
                                                                                        new_node = newRowNode( col, value );
                                                                                        trail_node->next = new_node;
                                                                                        new_node->next = cur_node;
                                                                                        (*count)++;
                                                                                        return( row );
                                                                                }  /* if we have found place to insert it */

                                                                                trail_node = cur_node;
                                                                                cur_node = cur_node->next;
                                                                        }  /* while */

                                                                        /* Case if it should be inserted at the end of the list.  This */
                                                                        /* should be the only case where we exit the bottom of the loop, all */
                                                                        /* other cases should exit the routine from within the loop */

                                                                        new_node = newRowNode( col, value );
                                                                        trail_node->next = new_node;
                                                                        (*count)++;
                                                                        return( row );

}  /* addEntryToRow */
/**********************************************************************/
void displayRow( I_Matrix_Row_Node row ) {

        if( row == NULL )
                printf( "<empty>");

        while( row != NULL ) {
                printf("[%d] %.3f ", row->column, row->value );
                row = row->next;
        }  /* while */

        printf( "\n");
}  /* displayRow */

/**********************************************************************/
/********************  Routines for intermediate matrix    ************/
/**********************************************************************/
I_Matrix newIMatrix( int num_rows ) {
        I_Matrix i_matrix;

        i_matrix = (I_Matrix) malloc( sizeof( *i_matrix));
        checkAllocatedPointer(i_matrix);

        i_matrix->num_rows = num_rows;

        i_matrix->row = (I_Matrix_Row_Node *)
                calloc( num_rows, sizeof( *(i_matrix->row) ));

        i_matrix->row_length = (int *) calloc( num_rows, sizeof( int ));

        return( i_matrix );
}  /* newIMatrix */
/**********************************************************************/
void destroyIMatrix( I_Matrix i_matrix ) {
        int i;

        free( i_matrix->row_length );

        for( i = 0; i < i_matrix->num_rows; i++ )
                destroyRow( i_matrix->row[i] );
        free( i_matrix->row );

        free( i_matrix );

}  /* destroyIMatrix */
/**********************************************************************/
int addEntryToIMatrix( I_Matrix i_matrix, int row, 
                                          int col, REAL_VALUE value ) {

                                                  assert(( i_matrix != NULL) 
                                                          && (row >=0) && ( row < i_matrix->num_rows ));

                                                  i_matrix->row[row] = addEntryToRow( i_matrix->row[row], col, value, 
                                                          &(i_matrix->row_length[row]), 0 );

                                                  return( 1 );
}  /* addEntryToIMatrix */
/**********************************************************************/
int accumulateEntryInIMatrix( I_Matrix i_matrix, int row, 
                                                         int col, REAL_VALUE value ) {
                                                                 /*
                                                                 This routine is the same as addEntryToIMatrix() except it will call 
                                                                 the addEntryToRow() routine with the accumulate flag set.  Thus, if
                                                                 there is a value already at this row and column, the new value will 
                                                                 be the sum of the old and the new value.
                                                                 */
                                                                 assert(( i_matrix != NULL) 
                                                                         && (row >=0) && ( row < i_matrix->num_rows ));

                                                                 i_matrix->row[row] = addEntryToRow( i_matrix->row[row], col, value, 
                                                                         &(i_matrix->row_length[row]), 1 );

                                                                 return( 1 );
}  /* addEntryToIMatrix */
/**********************************************************************/
int countEntriesInIMatrix( I_Matrix i_matrix ) {
        int i;
        int total = 0;

        for( i = 0; i < i_matrix->num_rows; i++ )
                total += i_matrix->row_length[i];

        return( total );
}  /* countEntriesInIMatrix */
/**********************************************************************/
REAL_VALUE sumIMatrixRowValues( I_Matrix i_matrix, int row ) {
        REAL_VALUE sum = 0.0;
        I_Matrix_Row_Node list;

        list = i_matrix->row[row];

        while( list != NULL ) {
                sum += list->value;
                list = list->next;
        }  /* while */
        return( sum );
}  /* sumIMatrixRowValues */
/**********************************************************************/
void displayIMatrix( I_Matrix i_matrix ) {
        int i;

        for( i = 0; i < i_matrix->num_rows; i++ ) {
                printf( "(len=%d, sum =%.1f)Row=%d: ", i_matrix->row_length[i], 
                        sumIMatrixRowValues( i_matrix, i ), i );
                displayRow( i_matrix->row[i] );

        }  /* for i */
        
} /* displayIMatrix */

/**********************************************************************/
/********************  Routines for sparse matrix    ******************/
/**********************************************************************/
Matrix newMatrix( int num_rows, int num_non_zero ) {
        Matrix matrix;

        matrix = (Matrix) malloc( sizeof( *matrix ));
        checkAllocatedPointer(matrix);

        matrix->num_rows = num_rows;
        matrix->num_non_zero = num_non_zero;

        matrix->mat_val = (REAL_VALUE *) 
                calloc( num_non_zero, sizeof( REAL_VALUE ));
        matrix->col = (int *) 
                calloc( num_non_zero, sizeof( int ));
        matrix->row_start = (int *) 
                calloc( num_rows, sizeof( int ));
        matrix->row_length = (int *) 
                calloc( num_rows, sizeof( int ));

        return( matrix );
}  /* newMatrix */
/**********************************************************************/
void destroyMatrix( Matrix matrix ) {

        if(matrix != NULL)
        {
                if(matrix->row_length != NULL)
                        free( matrix->row_length );
                
                if(matrix->row_start != NULL)
                        free( matrix->row_start );
                
                if(matrix->col != NULL)
                        free( matrix->col );

                if(matrix->mat_val != NULL)
                        free( matrix->mat_val );

                free( matrix );
        }
}  /* destroyMatrix */
/**********************************************************************/
Matrix transformIMatrix( I_Matrix i_matrix ) {
        /*
        Will convert a matrix object in intermediate form into a sparse
        representation.  It does not free the memory of the intermediate
        representation.
        */
        Matrix matrix;
        int row;
        I_Matrix_Row_Node cur_node;
        int index = 0;  /* holds the position in the mat_val and col arrays */

        /* Allocate the appropriate amount of memory */
        matrix = newMatrix( i_matrix->num_rows, 
                countEntriesInIMatrix( i_matrix ));

        /* Now go through and set the values */
        for( row = 0; row < i_matrix->num_rows; row++ ) {

                matrix->row_start[row] = index;
                matrix->row_length[row] = i_matrix->row_length[row];

                cur_node = i_matrix->row[row];
                while( cur_node != NULL ) {

                        matrix->col[index] = cur_node->column;
                        matrix->mat_val[index] = cur_node->value;
                        index++;

                        cur_node = cur_node->next;
                }  /* while cur_node */

        }  /* for row */

        assert( index == matrix->num_non_zero );

        return( matrix );
}  /* transformIMatrix */
/**********************************************************************/
REAL_VALUE sumRowValues( Matrix matrix, int row ) {
        REAL_VALUE sum = 0.0;
        int col;

        for( col = matrix->row_start[row];
                col < matrix->row_start[row] + matrix->row_length[row];
                col++ )
                sum += matrix->mat_val[col];

        return( sum );
}  /* sumRowValues */
/**********************************************************************/
void displayMatrix( Matrix matrix ) {
        int i, j;

        for( i = 0; i < matrix->num_rows; i++ ) {
                printf( "(len=%d, sum=%.1f)Row=%d: ", matrix->row_length[i], 
                        sumRowValues( matrix, i ), i );


                if( matrix->row_length[i] == 0 )
                        printf( "<empty>");

                for( j = matrix->row_start[i];
                        j < matrix->row_start[i] + matrix->row_length[i];
                        j++ )
                        printf("[%d] %.3f ", matrix->col[j], matrix->mat_val[j] );

                printf( "\n");

        }  /* for i */

} /* displayMatrix */
/**********************************************************************/
REAL_VALUE getEntryMatrix( Matrix matrix, int row, int col ) {
        /*
        Returns the value for a particular entry in the matrix.  It looks in
        the row for the specified column, and if it is not found it will
        return zero, since after all this is a sparse representation.
        */
        int j;

        for( j = matrix->row_start[row];
                j < matrix->row_start[row] + matrix->row_length[row];
                j++ )

                if( matrix->col[j] == col )
                        return( matrix->mat_val[j] );

        return( 0.0 );
}  /* getEntryMatrix */
/**********************************************************************/