二維矩陣的抽象資料結構
二維矩陣是線性代數中基本的組成單位。現在有許多程式語言或函式庫，像是 MATLAB 或 R 等，都內建矩陣運算的功能。實務上，當我們要用 C 或 Fortran 進行向量或矩陣運算時，我們會使用 BLAS 系列函式庫，因為 BLAS 是以矩陣運算聞名於世的 Fortran 實作而成，會比手刻的矩陣來得有效率。因此，本範例程式重點在於學習其原理，而非重造輪子來用。
以下是二維矩陣的抽象資料結構：
註：由於矩陣運算項目較多，本文僅列出一些基本的矩陣操作。
M, N, K are matrices
s is a scalar
sub Col(M): size
sub Row(M): size
sub At(M, i, j): result
sub SetAt(M, i, j, value): void
sub Trans(M): N
sub Add(M, N): K
sub Add(s, M): K
sub Sub(M, N): K
sub Sub(s, M): K
sub Sub(M, s): K
sub Mul(M, N): K
sub Mul(s, M): K
sub Div(M, N): K
sub Div(s, M): K
sub Div(M, s): K
sub Prod(M, N): K以下是和二維矩陣的狀態相關的行為：
Col(M)：矩陣的 (直) 行數Row(M)：矩陣的 (橫) 列數At(M, i, j)：取得矩陣中任意位置的值SetAt(M, i, j, value)：存入矩陣中任意位置的值
以下是二維矩陣的操作：
Trans(M)：將矩陣轉置
根據二維矩陣和純量間的相對位置，可得以下三種情境：
Add(M, N)：矩陣和矩陣相加Add(s, M)：純量和矩陣相加，純量在前Add(M, s)：矩陣和純量相加，純量在後
由加法和乘法具有交換律，在後兩個情境者是相等的；但在減法和除法則不等。
以下是二維矩陣的其他運算：
Prod(M, N)：兩矩陣相乘後得純量
宣告二維矩陣類別
以下是二維矩陣的類別宣告：
class Matrix
row <- m, m > 0
col <- n, n > 0
elements <- a new array which size == m * n在內部儲存矩陣元素時，可用二維陣列或一維陣列來儲存；使用二維陣列時，索引值和內部陣列可以一比一對應；使用一維陣列時，索引值要經過轉換才可對應到內部的陣列。一般來說，使用一維陣列的效率會稍微好一些些，故我們採用一維陣列來儲存矩陣元素。
以下是等效的 C 語言程式碼：
typedef struct matrix_t matrix_t;
struct matrix_t {
size_t col;
size_t row;
double *elements;
};二維矩陣的建構函式
以下是 C 語言版本的矩陣建構函式：
matrix_t * matrix_new(size_t col, size_t row)
{
assert(0 < col);
assert(0 < row);
matrix_t *mtx = malloc(sizeof(matrix_t));
if (!mtx)
return mtx;
mtx->col = col;
mtx->row = row;
mtx->elements = calloc(col * row, sizeof(double));
if (!(mtx->elements)) {
matrix_delete(mtx);
mtx = NULL;
return mtx;
}
return mtx;
}該建構函式會建立維度為 col 和 row 的矩陣，並將其初始化為 0.0。該建構函式是比較傳統的，我們也可以在建立矩陣物件時一併將其初始化，使用方式如下：
matrix_t *mtx = matrix_init(3, 2,
1.0, 2.0, 3.0,
4.0, 5.0, 6.0
);在這個建構函式 matrix_init 中，前兩個參數是矩陣的維度，之後的參數則是實際要填入矩陣的值。這個建構函式的實作如下：
matrix_t * matrix_init( /* 1 */
size_t col, /* 2 */
size_t row, /* 3 */
double value, ...) /* 4 */
{ /* 5 */
assert(0 < col); /* 6 */
assert(0 < row); /* 7 */
matrix_t *mtx = \
malloc(sizeof(matrix_t)); /* 8 */
if (!mtx) /* 9 */
return mtx; /* 10 */
mtx->col = col; /* 11 */
mtx->row = row; /* 12 */
mtx->elements = \
calloc(col * row, sizeof(double)); /* 13 */
if (!(mtx->elements)) { /* 14 */
matrix_delete(mtx); /* 15 */
mtx = NULL; /* 16 */
return mtx; /* 17 */
} /* 18 */
mtx->elements[0] = value; /* 19 */
va_list args; /* 20 */
va_start(args, value); /* 21 */
for (size_t j = 0; j < row; j++) { /* 22 */
for (size_t i = 0; i < col; i++) { /* 23 */
if (i == 0 && j == 0) { /* 24 */
continue; /* 25 */
} /* 26 */
mtx->elements[i+col*j] = \
va_arg(args, double); /* 27 */
} /* 28 */
} /* 29 */
va_end(args); /* 30 */
return mtx; /* 31 */
} /* 32 */第 8 行至第 18 行基本上和傳統的矩陣建構函式無異。但此建構函式在第 19 行至第 30 行間利用不定數量參數將此矩陣物件初始化。
二維矩陣的解構函式
以下是 C 語言版本的矩陣解構函式：
void matrix_delete(void *self)
{
if (!self) {
return;
}
double *elements = \
((matrix_t *) self)->elements;
if (elements) {
free(elements);
}
free(self);
}先釋放內部用來儲存矩陣元素的陣列，再釋放矩陣物件本身即可。
得知二維矩陣的大小
由於我們在建立矩陣物件時，已經儲存矩陣的維度，直接讀取其值即可。以下是等效的 C 語言程式碼：
size_t matrix_col(const matrix_t *self)
{
assert(self);
return self->col;
}
size_t matrix_row(const matrix_t *self)
{
assert(self);
return self->row;
}取得矩陣中任意位置的值
先確認索引值是合法的，再進行下一步動作。由於本範例程式內部是以一維矩陣儲存元素，要將座標進行相對應的轉換。轉換的方式有 column major 和 row major 兩種。轉換方式如下：
(m, n) is the matrix dimension
Column major:
(i, j) => i + m * j
Row major:
(i, j) => i * n + j兩種轉換方式都可以使用，重點是要在整個物件內保持一致。本範例程式採用 column major。
以下是等效的 C 語言實作：
double matrix_at(const matrix_t *self, size_t col, size_t row)
{
assert(col < matrix_col(self));
assert(row < matrix_row(self));
return self->elements[col + row * matrix_col(self)];
}設置二維矩陣中任意位置的值
先確認索引值是合法的，接著將座標進行相對應的轉換即可。由於我們先前採用 column major，這裡也保持一致。以下是等效的 C 語言實作：
void matrix_set_at(matrix_t *self, size_t col, size_t row, double value)
{
assert(col < matrix_col(self));
assert(row < matrix_row(self));
self->elements[col + row * matrix_col(self)] = value;
}檢查兩矩陣是否相等
先確認兩矩陣的維度相等，再逐一確認兩矩陣的每個元素皆相等，即可確認兩矩陣相等。以下是等效的 C 語言實作：
bool matrix_is_equal(const matrix_t *m, const matrix_t *n)
{
if (!(matrix_col(m) == matrix_col(n))) {
return false;
}
if (!(matrix_row(m) == matrix_row(n))) {
return false;
}
for (size_t i = 0; i < matrix_col(m); i++) {
for (size_t j = 0; j < matrix_row(m); j++) {
if (!(fabs(matrix_at(m, i, j) - matrix_at(n, i, j)) < 0.000001)) {
return false;
}
}
}
return true;
}由於本範例程式以浮點數儲存元素，會有一些誤差；本範例程式在兩者誤差小於 10^-6 時，將兩元素視為相等。讀者可視自身的需求修改本段程式碼的定義。
矩陣轉置
乍看程式碼略為複雜，但只是按照其數學定義來操作，追蹤一下程式碼即可了解。以下是等效的 C 語言實作：
matrix_t * matrix_trans(const matrix_t *m)
{
matrix_t *out = matrix_new(matrix_row(m), matrix_col(m));
if (!out) {
return out;
}
for (size_t i = 0; i < matrix_col(m); i++) {
for (size_t j = 0; j < matrix_row(m); j++) {
matrix_set_at(out, j, i, matrix_at(m, i, j));
}
}
return out;
}兩矩陣相加
要先確認兩矩陣的維度相等，再進行下一步動作。在本範例程式中，我們另建新的矩陣，這樣比較耗記憶體，如果想要省記憶體，可以將其相加後的值回存在第一個矩陣中。
在實作矩陣加法時，我們利用 C11 後的泛型巨集敘述來達成多型的特性。以下是其宣告：
#if __STDC_VERSION__ >= 201112L
#define matrix_add(m, n) \
_Generic((m), \
matrix_t *: _Generic((n), \
matrix_t *: matrix_add_mm, \
default: matrix_add_ms), \
default: matrix_add_sm)((m), (n))
#else
matrix_t * matrix_add(const matrix_t *m, const matrix_t *n);
#endif當編譯器未支援 C11 時，我們將 matrix_add 改回一般的矩陣加法，做為退路 (fallback)。
以下是實際的 C 語言程式碼：
#if __STDC_VERSION__ < 201112L
matrix_t * matrix_add(const matrix_t *m, const matrix_t *n)
{
return matrix_add_mm(m, n);
}
#endif
matrix_t * matrix_add_mm(const matrix_t *m, const matrix_t *n)
{
assert(matrix_col(m) == matrix_col(n));
assert(matrix_row(m) == matrix_row(n));
matrix_t *out = matrix_new(matrix_col(m), matrix_row(m));
if (!out) {
return out;
}
for (size_t i = 0; i < matrix_col(m); i++) {
for (size_t j = 0; j < matrix_row(m); j++) {
matrix_set_at(out, i, j, matrix_at(m, i, j) + matrix_at(n, i, j));
}
}
return out;
}矩陣和純量相加
按照其數學定義操作即可。矩陣加法具有交換性，故只需實作一次。以下是等效的 C 語言程式碼：
matrix_t * matrix_add_ms(const matrix_t *m, double s)
{
matrix_t *out = matrix_new(matrix_col(m), matrix_row(m));
if (!out) {
return out;
}
for (size_t i = 0; i < matrix_col(m); i++) {
for (size_t j = 0; j < matrix_row(m); j++) {
matrix_set_at(out, i, j, matrix_at(m, i, j) + s);
}
}
return out;
}
matrix_t * matrix_add_sm(double s, const matrix_t *m)
{
return matrix_add_ms(m, s);
}矩陣內積 (Product)
內積的數學定義並沒有特別困難，但實作時很容易搞錯座標，需仔細追蹤程式碼。以下是等效的 C 語言程式碼：
matrix_t * matrix_prod(const matrix_t *m, const matrix_t *n)
{
assert(matrix_col(m) == matrix_row(n));
matrix_t* out = \
matrix_new(matrix_row(m), matrix_col(n));
if (!out) {
return out;
}
double temp;
for (size_t i = 0; i < matrix_row(m); i++) {
for (size_t j = 0; j < matrix_col(n); j++) {
temp = 0.0;
for (size_t k = 0; k < matrix_col(m); k++) {
temp += \
matrix_at(m, k, i) * matrix_at(n, j, k);
}
matrix_set_at(out, j, i, temp);
}
}
return out;
}演算法上的效率
以下是和矩陣的狀態相關的行為：
Col(M)：O(1)Row(M)：O(1)At(M, i, j)：O(1)SetAt(M, i, j, value)：O(1)
以下是矩陣的基本四則運算：
Add(M, N)：O(mn)Add(s, M)：O(mn)
以下是矩陣操作的效率：
Trans(M)：O(mn)
以下是空間複雜度：
- 空間使用率：
O(mn)
除了少數特殊的演算法可以增進一些些效率外，基本上矩陣的操作和運算的效率和其維度成正比。
一般矩陣的實作中，矩陣中所有的值都會存在物件中，需占用較大的記憶體。有些矩陣大部分的值都為零，這時候就會用稀疏矩陣 (sparse matrix) 來節省記憶體。我們將於後文介紹稀疏矩陣。