目次
1. はじめに
C言語はシステム開発や組み込みプログラミングで広く使われる言語であり、行列計算もその一部として重要な役割を果たします。本記事では、C言語を用いた行列積の実装方法を解説し、基本から最適化、高速化まで幅広くカバーします。
本記事の対象読者：
- C言語を学習中の初心者
- 行列計算の実装方法を知りたい人
- 行列積の最適化や高速化に興味がある人
この記事を読むことで、C言語で行列積を実装する基本的な方法や、効率的に計算するための手法を習得できます。まずは、行列積の基礎から学んでいきましょう。
2. 行列積の基礎
2.1 行列とは？行列積の基本をわかりやすく解説
行列（Matrix）は、数値が縦横に整列した構造を持つデータの集合です。一般的に、行と列を持つ2次元の配列として表現されます。
例えば、以下のような2×3の行列があります：
A =
[ 1 2 3 ]
[ 4 5 6 ]行列積は、2つの行列を掛け合わせた結果を求める演算であり、次の条件を満たす必要があります：
- 左側の行列 A の列数 と 右側の行列 B の行数 が等しいこと。
例えば、次のように2×3行列と3×2行列を掛けると、2×2の行列が得られます。
A =
[ 1 2 3 ]
[ 4 5 6 ]
B =
[ 7 8 ]
[ 9 10 ]
[11 12 ]このときの行列積 C = A × B は以下のように計算されます。
C =
[ 58 64 ]
[139 154 ]2.2 行列積の計算方法｜具体例と数式で理解する
行列積の計算は、各行の要素と列の要素を掛け合わせ、合計することで求められます。
一般に、行列 A（サイズ m × n）と行列 B（サイズ n × p）を掛ける場合、得られる行列 C のサイズは m × p となります。
行列積の一般的な計算式：
C[i][j] = Σ(A[i][k] * B[k][j]) (k=0 から n-1)ここで、
C[i][j]は行列Cのi行j列の要素A[i][k]は行列Aのi行k列の要素B[k][j]は行列Bのk行j列の要素
2.3 行列積の性質（結合法則、分配法則など）
行列積には以下のような性質があります。
- 結合法則（Associativity）
(A × B) × C = A × (B × C)ただし、行列の掛け算は交換可能ではありません。
- 分配法則（Distributive Law）
A × (B + C) = A × B + A × C- 単位行列（Identity Matrix）
単位行列Iは、すべての対角要素が 1、その他の要素が 0 の行列で、任意の行列Aに対して次の関係が成り立ちます。
A × I = I × A = A行列積の性質を理解することで、計算の最適化や応用の幅が広がります。
3. C言語で行列積を実装する方法｜サンプルコード付き
3.1 C言語で行列を2次元配列で表現する方法
固定サイズの2次元配列を使う
行列を 固定サイズの2次元配列 で表現する方法はシンプルでわかりやすく、初心者向けです。
#include <stdio.h>
#define ROWS 2
#define COLS 3
int main() {
int matrix[ROWS][COLS] = {
{1, 2, 3},
{4, 5, 6}
};
// 行列の表示
for (int i = 0; i < ROWS; i++) {
for (int j = 0; j < COLS; j++) {
printf("%d ", matrix[i][j]);
}
printf("\n");
}
return 0;
}動的メモリ確保を使う
行列のサイズを柔軟に扱いたい場合は、 malloc() を使って動的にメモリを確保 します。
#include <stdio.h>
#include <stdlib.h>
int main() {
int rows = 2, cols = 3;
int **matrix;
// メモリ確保
matrix = (int **)malloc(rows * sizeof(int *));
for (int i = 0; i < rows; i++) {
matrix[i] = (int *)malloc(cols * sizeof(int));
}
// データの入力
int value = 1;
for (int i = 0; i < rows; i++) {
for (int j = 0; j < cols; j++) {
matrix[i][j] = value++;
}
}
// 行列の表示
for (int i = 0; i < rows; i++) {
for (int j = 0; j < cols; j++) {
printf("%d ", matrix[i][j]);
}
printf("\n");
}
// メモリ解放
for (int i = 0; i < rows; i++) {
free(matrix[i]);
}
free(matrix);
return 0;
}3.2 行列積をC言語で実装｜基本コードを解説
ここでは、2つの行列の積を計算する基本的なプログラム を紹介します。
#include <stdio.h>
#define N 3 // 行列のサイズ
// 行列の積を計算する関数
void matrixMultiplication(int A[N][N], int B[N][N], int C[N][N]) {
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
C[i][j] = 0; // 初期化
for (int k = 0; k < N; k++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
}
// 行列を表示する関数
void printMatrix(int matrix[N][N]) {
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
printf("%d ", matrix[i][j]);
}
printf("\n");
}
}
int main() {
int A[N][N] = {
{1, 2, 3},
{4, 5, 6},
{7, 8, 9}
};
int B[N][N] = {
{9, 8, 7},
{6, 5, 4},
{3, 2, 1}
};
int C[N][N]; // 結果行列
// 行列積の計算
matrixMultiplication(A, B, C);
// 結果の表示
printf("行列 A:\n");
printMatrix(A);
printf("行列 B:\n");
printMatrix(B);
printf("行列積 C:\n");
printMatrix(C);
return 0;
}このコードでは、matrixMultiplication() 関数で3重ループを使用し、各要素ごとの積と和を計算しています。
4. 行列積の最適化手法
4.1 ループ最適化で行列積を高速化
基本的な行列積のコードでは、3重ループ を使って計算を行います。しかし、このループの順序を工夫することで キャッシュ効率を向上させ、処理速度を速くする ことが可能です。
メモリアクセスとキャッシュ効率
CPUのキャッシュは、メモリアクセスの局所性が高いと効率的に動作します。行列積の計算では、次のようなループ順序の変更が効果的です。
#include <stdio.h>
#define N 3 // 行列のサイズ
// ループ順序を変更してキャッシュ効率を改善
void optimizedMatrixMultiplication(int A[N][N], int B[N][N], int C[N][N]) {
for (int i = 0; i < N; i++) {
for (int k = 0; k < N; k++) { // k を外側にする
for (int j = 0; j < N; j++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
}
int main() {
int A[N][N] = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};
int B[N][N] = {{9, 8, 7}, {6, 5, 4}, {3, 2, 1}};
int C[N][N] = {0};
optimizedMatrixMultiplication(A, B, C);
// 結果表示
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
printf("%d ", C[i][j]);
}
printf("\n");
}
return 0;
}4.2 ストラッセン法で行列積を最適化
ストラッセン法（Strassen Algorithm）は、行列積の計算量を O(N³) → O(N^{2.81}) に削減できるアルゴリズムです。
ストラッセン法のC言語実装
#include <stdio.h>
#include <stdlib.h>
// 2×2 行列に対するストラッセン法
void strassen(int A[2][2], int B[2][2], int C[2][2]) {
int M1 = (A[0][0] + A[1][1]) * (B[0][0] + B[1][1]);
int M2 = (A[1][0] + A[1][1]) * B[0][0];
int M3 = A[0][0] * (B[0][1] - B[1][1]);
int M4 = A[1][1] * (B[1][0] - B[0][0]);
int M5 = (A[0][0] + A[0][1]) * B[1][1];
int M6 = (A[1][0] - A[0][0]) * (B[0][0] + B[0][1]);
int M7 = (A[0][1] - A[1][1]) * (B[1][0] + B[1][1]);
C[0][0] = M1 + M4 - M5 + M7;
C[0][1] = M3 + M5;
C[1][0] = M2 + M4;
C[1][1] = M1 - M2 + M3 + M6;
}
int main() {
int A[2][2] = {{1, 2}, {3, 4}};
int B[2][2] = {{5, 6}, {7, 8}};
int C[2][2];
strassen(A, B, C);
// 結果表示
for (int i = 0; i < 2; i++) {
for (int j = 0; j < 2; j++) {
printf("%d ", C[i][j]);
}
printf("\n");
}
return 0;
}4.3 OpenMPで並列処理を行う
行列積の計算は 各要素が独立して計算できる ため、並列処理 による高速化が可能です。C言語では OpenMP を使うことで簡単に並列処理を実装できます。
#include <stdio.h>
#include <omp.h>
#define N 3 // 行列サイズ
// OpenMPを使用した並列行列積
void parallelMatrixMultiplication(int A[N][N], int B[N][N], int C[N][N]) {
#pragma omp parallel for collapse(2)
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
C[i][j] = 0;
for (int k = 0; k < N; k++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
}
int main() {
int A[N][N] = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};
int B[N][N] = {{9, 8, 7}, {6, 5, 4}, {3, 2, 1}};
int C[N][N] = {0};
parallelMatrixMultiplication(A, B, C);
// 結果表示
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
printf("%d ", C[i][j]);
}
printf("\n");
}
return 0;
}
5. C言語の行列積で発生しやすいエラーと対処法
5.1 メモリリークを防ぐ方法
メモリリークとは？
動的メモリ確保 (malloc) を行った後に free で解放しない と、不要なメモリがプログラム終了後も確保されたままとなり、メモリリークが発生します。これにより、メモリ使用量が増え続け、プログラムの動作が不安定になる可能性があります。
誤ったコード例
#include <stdio.h>
#include <stdlib.h>
int main() {
int rows = 3, cols = 3;
int **matrix;
// メモリ確保
matrix = (int **)malloc(rows * sizeof(int *));
for (int i = 0; i < rows; i++) {
matrix[i] = (int *)malloc(cols * sizeof(int));
}
return 0; // freeを忘れているため、メモリリークが発生
}対処法
確保したメモリは、必ず解放する必要があります。
#include <stdio.h>
#include <stdlib.h>
int main() {
int rows = 3, cols = 3;
int **matrix;
// メモリ確保
matrix = (int **)malloc(rows * sizeof(int *));
for (int i = 0; i < rows; i++) {
matrix[i] = (int *)malloc(cols * sizeof(int));
}
// メモリ解放
for (int i = 0; i < rows; i++) {
free(matrix[i]);
}
free(matrix);
return 0;
}5.2 配列外参照エラーの原因と解決策
配列外参照とは？
配列の範囲外のメモリにアクセスすると、予期しない動作が発生し、場合によっては セグメンテーションフォルト（Segmentation Fault） というエラーになります。
誤ったコード例
#include <stdio.h>
#define N 3
int main() {
int matrix[N][N] = {
{1, 2, 3},
{4, 5, 6},
{7, 8, 9}
};
// 範囲外アクセス (N = 3 のとき、matrix[3][0] は存在しない)
printf("%d\n", matrix[3][0]);
return 0;
}対処法
#include <stdio.h>
#define N 3
int main() {
int matrix[N][N] = {
{1, 2, 3},
{4, 5, 6},
{7, 8, 9}
};
for (int i = 0; i < N; i++) {
for (int j = 0; j < N; j++) {
printf("%d ", matrix[i][j]);
}
printf("\n");
}
return 0;
}5.3 デバッグツールGDBを使ったエラー解析
エラーの発生箇所を特定するために、C言語のデバッグツール GDB（GNU Debugger） を使うと便利です。
GDBの基本的な使い方
- コンパイル時にデバッグ情報を含める
gcc -g program.c -o program- GDBを起動
gdb ./program- プログラムを実行
run- エラー発生時のスタックトレースを確認
backtraceGDBを使った実例
(gdb) run
Starting program: ./program
Program received signal SIGSEGV, Segmentation fault.
0x0000555555555142 in main () at program.c:10
10 printf("%d\n", matrix[3][0]); // ここでエラー発生6. 行列積の応用例
6.1 グラフィックスにおける変換行列
行列を使った2D・3D変換
コンピュータグラフィックス（CG）では、オブジェクトの移動（平行移動）、回転、スケーリング（拡大・縮小） などを行列を使って計算します。
2Dの変換行列
2次元の座標変換は、次のような 3×3の行列 を使います。
[x'] [ a b tx ] [x]
[y'] = [ c d ty ] * [y]
[ 1] [ 0 0 1 ] [1]C言語で2D変換行列を適用
#include <stdio.h>
#include <math.h>
#define PI 3.14159265
// 2D座標を回転する関数
void rotate2D(float x, float y, float angle, float *x_out, float *y_out) {
float rad = angle * PI / 180.0; // 角度をラジアンに変換
float rotationMatrix[2][2] = {
{cos(rad), -sin(rad)},
{sin(rad), cos(rad)}
};
// 行列積を計算
*x_out = rotationMatrix[0][0] * x + rotationMatrix[0][1] * y;
*y_out = rotationMatrix[1][0] * x + rotationMatrix[1][1] * y;
}
int main() {
float x = 1.0, y = 0.0;
float angle = 90.0;
float x_new, y_new;
rotate2D(x, y, angle, &x_new, &y_new);
printf("回転後の座標: (%.2f, %.2f)\n", x_new, y_new);
return 0;
}6.2 機械学習におけるニューラルネットワークの重み計算
ニューラルネットワークでは、入力データ X に対して、重み W を掛けて次の層の入力を計算します。
Y = X × WC言語で単純なニューラルネットワークの重み計算
#include <stdio.h>
// 行列積を計算する関数
void matrixMultiply(float X[2][3], float W[3][2], float Y[2][2]) {
for (int i = 0; i < 2; i++) {
for (int j = 0; j < 2; j++) {
Y[i][j] = 0;
for (int k = 0; k < 3; k++) {
Y[i][j] += X[i][k] * W[k][j];
}
}
}
}
int main() {
float X[2][3] = {
{1, 2, 3},
{4, 5, 6}
};
float W[3][2] = {
{0.1, 0.2},
{0.3, 0.4},
{0.5, 0.6}
};
float Y[2][2];
matrixMultiply(X, W, Y);
// 結果表示
for (int i = 0; i < 2; i++) {
for (int j = 0; j < 2; j++) {
printf("%.2f ", Y[i][j]);
}
printf("\n");
}
return 0;
}6.3 物理シミュレーションにおける状態遷移
状態遷移行列
例えば、2Dの物体の位置を (x, y)、速度を (vx, vy) とすると、時間 t の経過後の新しい状態 S' は、以下の行列積で表せます。
S' = T × S[x'] [ 1 0 t 0 ] [ x ]
[y'] = [ 0 1 0 t ] * [ y ]
[vx'] [ 0 0 1 0 ] [vx ]
[vy'] [ 0 0 0 1 ] [vy ]C言語での実装例
#include <stdio.h>
// 状態遷移を計算
void stateTransition(float T[4][4], float S[4], float S_new[4]) {
for (int i = 0; i < 4; i++) {
S_new[i] = 0;
for (int j = 0; j < 4; j++) {
S_new[i] += T[i][j] * S[j];
}
}
}
int main() {
float T[4][4] = {
{1, 0, 1, 0},
{0, 1, 0, 1},
{0, 0, 1, 0},
{0, 0, 0, 1}
};
float S[4] = {0, 0, 1, 1}; // 初期位置(0,0) 速度(1,1)
float S_new[4];
stateTransition(T, S, S_new);
printf("新しい状態: %.2f %.2f %.2f %.2f\n", S_new[0], S_new[1], S_new[2], S_new[3]);
return 0;
}7. FAQ（よくある質問）
7.1 行列のサイズが大きくなると計算速度が遅くなるのはなぜ？
計算量が O(N³) で増加するため
基本的な行列積の計算は、3重ループを用いるため、計算量が O(N³) になります。
行列サイズが2倍になると、計算回数は 8倍 になるため、サイズが大きくなるほど処理時間が急増します。
解決策
- 最適化手法を活用
- ループの順序を変更（キャッシュ効率向上）
- ストラッセン法を使用（計算量を O(N^{2.81}) に削減）
- OpenMP などの並列化技術を活用
- 専用の数値計算ライブラリを使用
- BLAS（Basic Linear Algebra Subprograms）
- Intel MKL（Math Kernel Library）
- OpenBLAS（オープンソースの高速線形代数ライブラリ）
7.2 ストラッセン法はどのような場合に有効？
大規模な行列の積を計算する場合
ストラッセン法は、行列のサイズが大きいほど高速化の効果が大きくなる アルゴリズムです。
ただし、次のような制約があります。
メリット
- 計算量が O(N^{2.81}) になり、大規模な行列積を高速に計算可能
- 一般的な O(N³) の方法と比べて、特に N が大きい場合に有利
デメリット
- 小さな行列ではオーバーヘッドが大きく、むしろ遅くなる
- 再帰呼び出しのため、メモリ使用量が増加する可能性
7.3 行列積の計算で注意すべき数値誤差は？
浮動小数点演算による誤差の蓄積
C言語では、浮動小数点（float や double）を用いた計算で 丸め誤差が発生 することがあります。
行列積では、多くの乗算と加算を行うため、この誤差が蓄積しやすいです。
数値誤差を減らす方法
double型を使用
float（32bit）よりdouble（64bit）の方が精度が高い- 例:
floatの精度は約7桁、doubleは約15桁
- 高精度演算ライブラリを活用
- GNU MP（GMP） や Intel MKL を利用して、高精度な演算を行う
7.4 C言語以外の言語での行列積の実装は？
Python（NumPyを使用）
Pythonでは NumPy を使うことで、非常に簡単に行列積を計算できます。
import numpy as np
A = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
B = np.array([[9, 8, 7], [6, 5, 4], [3, 2, 1]])
C = np.dot(A, B) # 行列積
print(C)Java（標準ライブラリを使用）
Javaでも 2次元配列を用いた実装 が可能。
public class MatrixMultiplication {
public static void main(String[] args) {
int[][] A = { {1, 2, 3}, {4, 5, 6}, {7, 8, 9} };
int[][] B = { {9, 8, 7}, {6, 5, 4}, {3, 2, 1} };
int[][] C = new int[3][3];
for (int i = 0; i < 3; i++) {
for (int j = 0; j < 3; j++) {
for (int k = 0; k < 3; k++) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
for (int i = 0; i < 3; i++) {
for (int j = 0; j < 3; j++) {
System.out.print(C[i][j] + " ");
}
System.out.println();
}
}
}
