C で行列を扱う: 動的なメモリ確保

計算工学及演習でメモリの動的確保とかポインタとか扱ってるけど、皆そこそこ苦労してるっぽいんで参考になればなー、と課題 1 の内容を動的なメモリ確保を使って書き直してみた。 参考になるかどうかはわからないが

というか C でクラスベースっぽくプログラムを書くサンプルになっちゃったかもしれない。 Java とか Perl とか JavaScript とかのオブジェクト指向言語 (Perl は微妙？) を日ごろ使ってる人が C を書くとこうなります的な。 C はあんまり得意じゃないんで 「ここはこうした方がいいんじゃない？」 という意見があればお願いします。

mail: webmaster[あっと]r-definition.com ([あっと] の部分は半角のアットマーク @ に変えてください)

コンパイル方法

"kadai1.5.c"、"Matrix.c"、"Matrix.h" を同じディレクトリに置いて、ターミナル上で

% gcc kadai1.5.c Matrix.c -lm

でコンパイルできます。 Windows XP + Cygwin + gcc で確認しています。 その他環境では適当にコマンド等を変更してください。

source code

kadai1.5.c のソースコード
/**
 * 電気電子計算工学及演習 課題 1 改変
 * 動的なメモリ領域確保
 */

#include <stdio.h>

#include "Matrix.h"

int main(void) {
    /*--------------------*
     * 変数宣言 
     *--------------------*/
    int n = 5; /* 行列のサイズ */
    Matrix matF = NULL;
    Matrix matA = NULL;
    Matrix matG = NULL;
    
    /*--------------------*
     * 処理 
     *--------------------*/
    matF = Matrix_createPascalMatrix(n);
    matA = Matrix_add(matF, matF);
    matG = Matrix_multiply(matF, matF);
    
    printf("P = ");
    Matrix_print(matF, "%-#17.14g ");
    printf("A = P+P = ");
    Matrix_print(matA, "%-#17.14g ");
    printf("G = P*P = ");
    Matrix_print(matG, "%-#17.14g ");
    
    printf( "infty norm of A = %#20.16g\n", Matrix_evalInftyNorm(matA) );
    
    /*--------------------*
     * メモリ領域開放
     *--------------------*/
    delete_Matrix(&matF);
    delete_Matrix(&matA);
    delete_Matrix(&matG);
    
    return 0;
}
Matrix.h のソースコード
/**
 * 電気電子計算工学及演習
 * 擬似クラス Matrix (header file)
 */

#ifndef MATRIX_H__
#define MATRIX_H__

#include <stdio.h>
#include <math.h>

/* ====================
 *   構造体
 * ==================== */
typedef struct {
    int _numRow; /* 行数 */
    int _numCol; /* 列数 */
    double** elem;
} _Matrix;
typedef _Matrix* Matrix;

/* ====================
 *   関数宣言
 * ==================== */

/**
 * n×m 行列を生成する. 行列の各要素の値は未定義とする. 
 * メモリ確保失敗時、NULL を返す. 
 * @param n 生成する行列の行数
 * @param m 生成する行列の列数
 * @return 生成した行列
 */
Matrix Matrix_createNullMatrix( int n, int m );

/**
 * n×n の Pascal Matrix を生成する.
 * メモリ確保失敗時、NULL を返す. 
 * @param n 行列のサイズ
 * @return 生成した行列
 */
Matrix Matrix_createPascalMatrix( int n );

/**
 * 行列の使用メモリを開放する.
 * @param matPt メモリを開放する行列へのポインタ
 */
void delete_Matrix( Matrix* matPt );


/**
 * 行列の内容を表示する.
 * @param mat 表示する対象の行列
 * @param format 各要素の値を表示する際のフォーマット
 */
void Matrix_print( Matrix mat, char* format );

/**
 * 行列同士の和算を行う.
 * @param mat1 演算対象の行列
 * @param mat2 演算対象の行列
 * @return 演算結果として生成される行列
 */
Matrix Matrix_add( Matrix mat1, Matrix mat2 );

/**
 * 行列同士の積算を行う.
 * @param mat1 演算対象の行列. 掛けられる側
 * @param mat2 演算対象の行列. 掛ける側
 * @return 演算結果として生成される行列
 */
Matrix Matrix_multiply( Matrix mat1, Matrix mat2 );

/**
 * 行列の無限大ノルムを求める.
 * @param mat 演算対象の行列
 * @return 求めた無限大ノルム
 */
double Matrix_evalInftyNorm( Matrix mat );

/**
 * 行列の 1-ノルムを求める
 * @param mat1 演算対象の行列
 * @return 求めた 1-ノルム
 */
double Matrix_eval1Norm( Matrix mat );

#endif
Matrix.c のソースコード
/**
 * 電気電子計算工学及演習
 * 擬似クラス Matrix
 */

#include "Matrix.h"

/* ==============================
 *   Private methods 宣言
 * ============================== */

/**
 * n×m 行列用のメモリ領域を確保する. 
 * 失敗時は NULL を返す. 
 */
Matrix _new( int n, int m );

/* ==============================
 *   Public methods 定義
 * ============================== */
Matrix Matrix_createNullMatrix( int aNumRow, int aNumCol ) {
    return _new( aNumRow, aNumCol );
}

Matrix Matrix_createPascalMatrix( int aSize ) {
    /* 変数宣言 */
    int i, j;
    Matrix mat;
    
    mat = _new( aSize, aSize );
    if( mat == NULL ) {
        return NULL;
    }
    /* 値を代入 */
    for( i = 0; i < mat->_numRow; i++ ) {
        for( j = 0; j < mat->_numCol; j++ ) {
            if( i == 0 || j == 0 ) {
                mat->elem[i][j] = 1;
            } else {
                mat->elem[i][j] = mat->elem[i-1][j] + mat->elem[i][j-1];
            }
        }
    }
    
    return mat;
}

void delete_Matrix(Matrix* aMatPt) {
    if( *aMatPt != NULL ) {
        if( (*aMatPt)->elem != NULL && (*aMatPt)->elem[0] != NULL ) {
            free((*aMatPt)->elem[0]);
            (*aMatPt)->elem[0] = NULL;
        }
        if( (*aMatPt)->elem != NULL ) {
            free((*aMatPt)->elem);
            (*aMatPt)->elem = NULL;
        }
    }
    free(*aMatPt);
    *aMatPt = NULL;
}

void Matrix_print( Matrix aMat, char* aFormat ) {
    int i, j;
    
    printf("[\n");
    for( i = 0; i < aMat->_numRow; i++ ) {
        for( j = 0; j < aMat->_numCol; j++ ) {
            printf( aFormat, aMat->elem[i][j] );
        }
        printf("\n");
    }
    printf("]\n");
}

Matrix Matrix_add( Matrix in1, Matrix in2 ) {
    int i, j;
    Matrix result;
    
    /* in1 と in2 の行列のサイズが等しいことを確認 */
    if( in1->_numRow != in2->_numRow || in1->_numCol != in2->_numCol ) {
        puts("[Matrix] EXCEPTION: 引数により与えられた 2 つの行列のサイズが異なります...");
        return NULL;
    }
    
    result = _new( in1->_numRow, in1->_numCol );
    
    for( i = 0; i < in1->_numRow; i++ ) {
        for( j = 0; j < in2->_numCol; j++ ) {
            result->elem[i][j] = in1->elem[i][j] + in2->elem[i][j];
        }
    }
    
    return result;
}

Matrix Matrix_multiply( Matrix in1, Matrix in2 ) {
    int i, j, k;
    Matrix result;
    
    /* in1 の列数と in2 の行数が等しいことを確認 */
    if( in1->_numCol != in2->_numRow ) {
        puts("[Matrix] EXCEPTION: 積計算において、1 つ目の行列の列数と 2 つ目の行列の行数が一致しません...");
        return NULL;
    }
    
    result = _new( in1->_numRow, in2->_numCol );
    
    for( i = 0; i < result->_numRow; i++ ) {
        for( j = 0; j < result->_numCol; j++ ) {
            result->elem[i][j] = 0;
            for( k = 0; k < result->_numRow; k++ ) {
                result->elem[i][j] += in1->elem[i][k] * in2->elem[k][j];
            }
        }
    }
    
    return result;
}

double Matrix_evalInftyNorm( Matrix aMat ) {
    int i, j, n;
    double max = 0.;
    double temp;
    
    for( i = 0; i < aMat->_numRow; i++ ) {
        temp = 0;
        for( j = 0; j < aMat->_numCol; j++ ) {
            /* 絶対値和をとる */
            if( aMat->elem[i][j] >= 0. ) {
                temp += aMat->elem[i][j];
            } else {
                temp -= aMat->elem[i][j];
            }
        }
        if( max < temp ) {
            max = temp;
        }
    }
    return max;
}

double Matrix_eval1Norm_1( Matrix aMat ) {
    int i, j;
    double max = 0;
    double temp;
    
    for( j = 0; j < aMat->_numCol; j++ ) {
        temp = 0;
        for( i = 0; i < aMat->_numRow; i++ ) {
            /* 絶対値和をとる */
            if( aMat->elem[i][j] >= 0 ) {
                temp += aMat->elem[i][j];
            } else {
                temp -= aMat->elem[i][j];
            }
        }
        if( max < temp ) {
            max = temp;
        }
    }
    return max;
}


/* ==============================
 *   Private methods 定義
 * ============================== */

Matrix _new( int aNumRow, int aNumCol ) {
    int i;
    Matrix mat;
    mat = (Matrix) malloc( sizeof(_Matrix) );
    if( mat == NULL ) {
        puts("[Matrix] メモリ不足...");
        return NULL;
    }
    mat->elem = (double**) malloc( aNumRow * sizeof(double*) );
    if( mat->elem == NULL ) {
        puts("[Matrix] メモリ不足...");
        free(mat);
        return NULL;
    }
    mat->elem[0] = (double*) malloc( aNumRow * aNumCol * sizeof(double) );
    if( mat->elem[0] == NULL ) {
        puts("[Matrix] メモリ不足...");
        free(mat->elem);
        free(mat);
        return NULL;
    }
    for( i = 1; i < aNumRow; i++ ) {
        mat->elem[i] = mat->elem[i-1] + aNumCol;
    }
    mat->_numRow = aNumRow;
    mat->_numCol = aNumCol;
    return mat;
}