2-4강 - Sparse Matrix (희소 행렬)

5. Sparse matices

0) 정의

1> sparse matrix : 원소의 다수가 0을 갖는 행렬.

2> 행렬을 2차원 배열로 표현할 때, 저장 공간의 낭비가 있다.

3> 선언 : arr[MAX_ROWS][MAX_COLS]

4> 접근 : arr[i][j] (행렬의 i행 j열에 접근)

 

1) ADT

ADT Spase_Matrix

1> object : <row, column, value> (3원소 쌍의 집합)

- value는 item 집합의 원소이다.

2> function

- Create(max_row, max_col) : maxItems(=max_rowXmax_col)개를 저장할 수 있는 sparse marix return

(최대 max_row개의 행과 max_col개의 열을 저장할 수 있는 sparse matix return)

- Transpose(a) : 3원소 쌍에서 row와 column을 모두 뒤바꿔 만든 matrix

- Add(a, b) : [1] a와 b의 dimension 같은지 확인 -> [2] row와 column이 같은 item들을 합한다.

- Multiply(a, b) : [1] a의 col==b의 row 확인 -> [2] d(i, j) = ∑ a(i, k)·b(k, j)으로 생성된 d가 원소인 matrix

 

2) Sparse Matrix Representation

(create 함수를 효과적으로 구현해보기 위해 명심할 것)

1> 3원소 쌍(row, col, value)의 array를 사용한다.

2> row의 index를 오름차순이 되도록 정렬한다.

3> 특정 row에 대한 3원소 쌍들은 column의 index가 오름차순이 되도록 정렬한다.

4> 연산의 종료를 보장하기 위해 행의 개수, 열의 개수, 0이 아닌 원소의 개수를 알아야 한다.

5> Create 함수

Sparse_Matrix Create(max_row, max_col)
{
	#define MAX_TERMS 101 /*maximum number of terms +1*/
	typedef	struct {
    	int col;
        int row;
		int value;
		} term;
	term a[MAX_TERMS];
}

- MAX_TERM = 최대 원소의 개수 + 1 

- 전체 행의 개수, 전체 열의 개수, 0이 아닌 원소의 개수를 저장할 공간을 하나 더 만들기 위해 1을 더함

6> sparse matrix 형태

 

3) Transposing a Matrix

1> simple algorithm : <i, j, value>를 <j, i, value>로 바꿔야 한다.

2> 문제

- row index 별로 transpose하면

- row가 우선으로 정렬되는 특성상

- transpose된 matrix는 col이 작은 것부터 생성되므로 

- 정확히 (구조체의) 몇 번째 index에 저장해야 할 지 정확히 알 수 없다.

3> 해결

- column index를 기준으로 transpose를 진행한다.

- for all element in column j, place element <i, j, value> in element <j, i, value>

4> Transpose 함수

- [0] 0번째 index에 들어있는 행의 개수, 열의 개수 할당하고 0이 아닌 값의 개수는 교환

- [1] a의 col 개수만큼 loop 만들기 => col 오름차순대로 시행하게 함

- [2] a의 원소 개수만큼 loop 만들기 => a의 모든 원소에 대해 transpose할 수 있게 한다

- [3] 그래서 'a[j].col==i' 일 때마다 transpose => a의 col이 작은 것부터 transpose

void transpose (term a[], term b[])
{/*b is set to the transpose of a */
	int n, i, j, currentb;
	n = a[0].value;			/*total number of elements*/
	b[0].row = a[0].col; 	/*rows in b = columns in a*/
    b[0].col = a[0].row; 	/*columns in b = rows in a*/
	b[0].value = n;
    
    if (n > 0) {
    /* nonzero matrix */
    currentb = 1;
    for (i = 0; i < a[0].col; i++)	/* 0번째 col부터 max_col까지 시도 */
    	for (j = 0; j <= n; j++)	/* a의 원소 개수 n번 만큼 시도*/
        	if (a[j].col == i) {	/* i번째 부터 (오름차순) transpose 연산 시도 */
            	b[currentb].row = a[j].col;
                b[currentb].col = a[j].row;
                b[currentb].value = a[j].value;
                currentb++
            }
    }
}

- Time complexity = O(column*elements) <2개의 for loop 곱한다.>

- If elements = rows·columns, (모든 원소가 0이 아닌 경우)

     Time complexity = O(column*elements)=O(rows*columns^2)

 

5> fast transpose 예시

- row_terms[MAX_COL] : col index 별로 몇 개의 nonzero element가 있는지 array 형태로 기록

- starting_pos[MAX_COL] : col index 별로 transpose 후에 b의 몇 번째 index부터 각 element를 기록할지

- [0] 적절한 값 할당

- [1] row_terms array를 일단 col 존재하는 것 0으로 초기화

- [2] row_terms array에 index 별로 row_terms 값 저장

- [3] starting_pos 초기화 : 0번째 index에 1 할당

- [4] starting_pos 값 할당 : 이 방식으로 자리 확보 => starting_pos[i]=starting_pos[i-1]+row_terms[i-1]

- [5] tranpose

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

#define MAX_COLS 101 
typedef struct{
	int col;
	int row;
	int value;
}term;

void fast_transpose(term a[], term b[])
{ /* the transpose of a is placed in b */
	int row_terms[MAX_COL]; // col index 별로 몇 개의 nonzero element가 있는지
	int starting_pos[MAX_COL]; // col index 별로 transpose 후 b에 저장된 위치를 알려주는 array
	int i, j;
	int num_cols = a[0].col;
	int num_terms = a[0].value;
	b[0].row = num_cols;
	b[0].col = a[0].row;
	b[0].value = num_terms;

	if (num_terms > 0) { /* nonzero matrix */
		for (i = 0; i < num_cols; i++) row_terms[i] = 0; // row_terms 초기화
		for (i = 1; i <= num_terms; i++) row_terms[a[i].col]++; // index 별로 row_terms 값 저장
		starting_pos[0] = 1;

		for (i = 1; i < num_cols; i++)
			starting_pos[i] = starting_pos[i - 1] + row_terms[i - 1]; // col별 원소 개수 만큼 여유 공간 확보
		for (i = 1; i <= num_terms; i++) {
			j = starting_pos[a[i].col]++; // col index별로 중복을 피하기 위해 j에 대입하고 1 증가시킨다.
			// transpose
			b[j].row = a[i].col;
			b[j].col = a[i].row;
			b[j].value = a[i].value;
		}
	}
}

- Time complexity = O(columns + elements) <columns, elements 만큼의 loop가 여러개 있어서>

- If elements = rows·columns, (모든 원소가 0이 아닌 경우)

     Time complexity = O(columns + elements) = O(rows·columns)

 

6> dense representation 형태의 transpose

for (int j = 0; j < col; j++)
	for (int i = 0; i < row; i++)
    	b[j][i] = a[i][j];

- Time complexity = O(rows·columns)

 

4) Matrix Multiplication

0> 정의

1> dense represenatation 형태의 matrix multiplication

int a[rows_a][cols_a];
int b[rows_b][cols_b];

int d[rows_a][cols_b; ]

for (int i = 0; i < rows_a; i++) {
	for (int j = 0; j < cols_b; j++) {
		sum = 0;
		for (int k = 0; k < cols_a; k++) {
			sum += a[i][k] * b[k][j];
		}
		d[i][j] = sum;
	}
}

- Time complexity = O(rows_a*cols_a*cols_b)

- 2개의 sparse matrix의 곱의 결과는 dense할 수 있다.

 

2> sparse matrix의 multiplication 원리

- Matrix A, B, D : term array a, b, d 내부에 3원소 쌍 (추가로 B_T도 term array new_b에 저장)

- row : B의 column에 곱하고 있는 A의 row

- row_begin : 현재 행에서 첫 번째 element가 term array a에서 위치

- column : A의 row가 곱하고 있는 B의 column

- totald : D matrix에서(term array d에서) 원소의 개수

- i, j : A의 row, B의 column을 연속적으로 검사하기 위한 pointer들

 

3> sparse marix의 multiplication

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

#define MAX_TERMS 101 /*maximum number of terms 1*/
typedef struct{
	int col;
	int row;
	int value;
}term;

void storesum(term d[], int *totald, int row, int column, int *sum)
{ /* if *sum != 0, then it along with its row and column position is stored as the *totald+1 entry in d */
	if(*sum)
		if (*totald < MAX_TERMS) {
			d[++*totald].row = row;
			d[*totald].col = column;
			d[*totald].value = *sum;
			*sum = 0;
		}
		else {
			fprintf(stderr, "Numbers of terms exceeds %d" , MAX_TERMS);
			exit(1);
		}
}

// a의 원소 개수 = 1(행렬 정보 저장) + nonzero element 개수 + 1(여유분)
void mmult(term a[], term b[], term d[])
/* multiply two sparse matrices */
{
	// 0. 변수 설정
	int i, j, column;
	int rows_a = a[0].row, cols_a = a[0].col, totala = a[0].value; // a의 행렬 정보 변수화
	int cols_b = b[0].col, totalb = b[0].value; // b의 행렬 정보 변수화
	if (cols_a != b[0].row) {
		fprintf(stderr, "Incompatible matrices n");
		exit(1);
	} // 1> b의 row 개수 확인 2> 행렬곱 가능 여부 확인 
	term new_b[MAX_TERMS];
	fast_transpose(b, new_b);
	// 1. 초기화
	int totald = 0; // d의 nonzero element 개수 초기화
	int row = a[1].row; // row 초기화: a의 원소 중 첫 번째 원소의 row로 (B의 column에 곱하고 있는 A의 row)
	int row_begin = 1; // row_begin 초기화: row 행의 1번째 원소가 nonzero라 가정하고 1로 초기화
	int sum = 0; // A의 특정 행과 B의 특정 열의 행렬곱 결과가 저장될 장소 (일단 0으로 초기화)
	/* set boundary condition */
	a[totala + 1].row = rows_a;
	new_b[totalb + 1].row = cols_b; 
	new_b[totalb + 1].col = -1;

	for (i = 1; i <= totala; ) { // a의 nonzero element 개수 만큼의 loop
		column = new_b[1].row; // column 초기화 : A의 row가 곱하고 있는 B의 column
		for (j = 1; j <= totalb + 1; ) { // new_b의 nonzero element 개수 만큼의 loop
			/* multiply row of a by column of b */
			// D의 원소 하나 생성할 때마다 if나 else if로 간다.
			if (a[i].row != row) { // 현재 a의 row가 내가 연산하려고 하는 row가 아니라면
				storesum(d, totald, row, column, &sum); 
				// 1> totald의 값 1 증가 => d의 nonzero element개수 counting
				// 2> d의 원소 중 증가된 totald번째에 이전에 진행한 row, column을 저장
				// 3> value에는 행렬곱 연산한 sum을 저장하고 0으로 초기화
				i = row_begin; // a의 첫번째 원소(a[1])부터 다시 시작
				for (; new_b[j].row == column; j++); 
                // 이전에 연산한 d의 column과 이제 연산할 b의 column이 다를 때까지 j 증가
				column = new_b[j].row; // 이제 연산할 b의 column으로 d의 column도 초기화
			}
			else if (new_b[j].row != column) { // 현재 b의 col이 내가 연산하려고 하는 col이 아니라면 
				storesum(d, totald, row, column, &sum);
				i = row_begin;
				column = new_b[j].row;
			}
			else switch (COMPARE(a[i].col, new_b[j].col)) { // a의 열(i)과 b의 행(j)을 비교 (자리가 맞아야 곱함)
			case -1: /* go to next term in a */
				i++; break;
			case 0: /* add terms, go to next term in a and b */
				sum += (a[i++].value * new_b[j++].value);
				break;
			case 1: /* go to next term in b */
				j++;
			}
		} /* end of for j <= totalb+1 */
		for (; a[i].row == row; i++);
		row_begin = i;
		row = a[i].row;
	} /* end of for i <= totala */
	d[0].row = rows_a;
	d[0].col = cols_b;
	d[0].value = totald;
}

- Time complexity = O(cols_b*totala+rows_a*totalb)