벡터와 행렬

이 단원에서 기술하는 함수들은 일반적인 C 배열에 대한 간단한 벡터와 행렬 인터페이스를 제공합니다. 이러한 배열의 메모리 관리는 후술할 ‘블럭(Block)’이라는 구현체로 이루어집니다. 사용자 정의 함수를 벡터와 행렬을 이용해서 작성하면, 별도의 추가 인자를 구현할 필요 없이 값과 차원을 모두 포함하는 단일 자료형을 인자로 전달할 수 있습니다. 이 구조는 BLAS에서 사용하는 벡터와 행렬과 호환됩니다.

자료형

이 단원에서 기술하는 자료형들은 대응되는 표준 자료형들이 존재합니다. 배정밀도 double 의 정밀도를 가지는 gsl_block , gsl_vector , 그리고 gsl_matrix 자료형이 존재하며, 단정밀도 float 의 경우 gsl_block_float , gsl_vector_float , 그리고 gsl_matrix_float 로 존재합니다.

사용 가능한 전체 자료형은 다음과 같습니다.

GSL 자료형	Type
gsl_block	double
gsl_block_float	float
gsl_block_long_double	long double
gsl_block_int	int
gsl_block_uint	unsigned int
gsl_block_long	long
gsl_block_ulong	unsigned long
gsl_block_short	short
gsl_block_ushort	unsigned short
gsl_block_char	char
gsl_block_uchar	unsigned char
gsl_block_complex	complex double
gsl_block_complex_float	complex float
gsl_block_complex_long_double	complex long double

벡터 gsl_vector 와 행렬 gsl_matrix 에 대해서도 동일한 형태로 존재합니다.

블럭

일관성을 위해 모든 메모리는 gsl_block 구조체를 통해 할당됩니다. 이 구조체는 2가지 변수를 포함합니다. 메모리의 크기를 담는 size 와 실제 데이터가 담긴 메모리를 가르키는 포인터 data 구성됩니다.(*) gsl_block 구조체는 다음과 같이 정의됩니다.

type gsl_block

typedef struct
{
    size_t size;
    double * data;
} gsl_block;

벡터와 행렬들은 이러한 블럭을 자른 조각들로 구현됩니다. 각 조각은 스텝 크기와 인덱스들의 조합, 그리고 시작 설정으로 구성됩니다. 행렬의 경우 열 인덱스의 스텝 크기는 행의 길이를 나타냅니다. 벡터의 스텝 크기는 보폭 stride 로 표기합니다.

이러한 블럭을 메모리에 할당하고, 해제하는 함수들은 gsl_block.h 로 정의되어 있습니다.

블럭 할당

블럭을 메모리에 할당하는 함수들은 malloc 와 free 형태로 이루어 집니다. 여기에 더해, 각자 고유의 오류 검사 기능을 수행합니다. 만약 블럭을 할당하기에 메모리가 충분하지 않다면 각각의 함수들은 GSL 오류 관리자를 호출해 GSL_ENOMEM 를 넘기고 함수는 NULL 포인터를 반환합니다. 따라서 라이브러리의 오류 관리자이 프로그램을 정지하는 형태로 사용한다면 모든 alloc 과정을 검사할 필요가 없습니다.

gsl_block *gsl_block_alloc(size_t n)

n 의 배정밀도 크기의 메모리를 블럭에 할당하고, 해당 블럭을 가르키는 포인터를 반환합니다. 이 블럭이 가진 메모리는 초기화가 이루어지지 않은 상태입니다. gsl_block_calloc() 를 사용하면 모든 메모리가 \(0\) 으로 초기화 된 상태의 블럭을 얻을 수 있습니다.

\(0\) 크기의 블록을 반환하라는 것도 유효합니다.

NULL 포인터가 아닌 정상적인 포인터를 반환합니다. NULL 포인터는 메모리를 할당하지 못했을 때 반환됩니다.

gsl_block *gsl_block_calloc(size_t n)

모든 원소를 \(0\) 으로 초기화한 블록을 가르키는 포인터를 반환합니다.

void gsl_block_free(gsl_block *b)

블럭 포인터 b 가르키는 메모리를 메모리에서 해제합니다.

블럭 읽고 쓰기 (Reading and writing blocks)

이 라이브러리에서는 이진 값이나 규격화 된 문서 파일들을 블럭에 읽고 쓸 수 있는 함수들을 제공합니다.

int gsl_block_fwrite(FILE *stream, const gsl_block *b)

블럭 b 원소 값들을 스트림 stream 바이너리 형태로 작성 합니다. 성공할 경우 \(0\) 값을 반환하고 파일을 작성하는 데 문제가 생길 경우 GSL_EFILED 값을 반환합니다. 데이터가 플랫폼의 기본 바이너리 규격으로 작성되기 때문에 다른 아키텍처에서는 호환이 안될 수 있습니다.

int gsl_block_fread(FILE *stream, gsl_block *b)

열려진 스트림 stream 로 부터 바이너리 규격 값을 읽어 블럭 b 작성합니다. 블럭 b 반드시 사전에 적절한 길이로 할당되어 있어야합니다. 왜냐하면, b 길이를 기반으로 얼마나 많은 데이터를 읽을 지 결정하기 때문입니다. 성공할 경우 \(0\) 값을 반환하고 값을 읽는 데 문제가 생길 경우 GSL_EFILED 값을 반환합니다. 읽는 값들은 동일한 플랫폼의 기본 바이너리 규격으로 작성되었다고 가정합니다.

int gsl_block_fprintf(FILE *stream, const gsl_block *b, const char *format)

블럭 b 원소들을 줄마다 스트림 stream 작성합니다. 이때, 형 결정자 format 로 값들이 작성됩니다. 이 결정자는 부동 소수점을 위한 %g , %e , 그리고 %f 가 존재하며, 정수를 위한 %d 이 있습니다. 성공할 경우 \(0\) 값을 반환하고 파일을 작성하는 데 문제가 생길 경우 GSL_EFILED 값을 반환합니다.

int gsl_block_fscanf(FILE *stream, gsl_block *b)

스트림 stream 로 부터 규격화 된 값을 읽어 블럭 b 에 저장합니다. 블럭 b 반드시 사전에 적절한 길이로 할당되어 있어야합니다. 왜냐하면, 이 함수는 b 길이를 기반으로 얼마나 많은 데이터를 읽을 지 결정하기 때문입니다. 성공할 경우 \(0\) 값을 반환하고 값을 읽는 데 문제가 생길 경우 GSL_EFILED 값을 반환합니다.

블럭 예제

다음 프로그램은 블럭을 할당하는 방법을 보여줍니다.

#include <stdio.h>
#include <gsl/gsl_block.h>

int
main (void)
{
  gsl_block * b = gsl_block_alloc (100);

  printf ("length of block = %zu\n", b->size);
  printf ("block data address = %p\n", b->data);

  gsl_block_free (b);
  return 0;
}

다음은 프로그램의 출력 결과입니다.

length of block = 100
block data address = 0x804b0d8

벡터

벡터는 gsl_vector 구조체로 블럭을 쪼갠 형태로 정의되어있습니다. 각기 다른 벡터가 같은 블록을 가르킬 수도 있습니다. 벡터 슬라이스는 동일한 크기 메모리 영역들의 집합입니다.

gsl_vector 구조체는 \(5\) 개의 성분을 가지고 있습니다. size, stride, 원소들이 저장되어 있는 메모리를 가리키는 포인터 data 벡터가 소유하고 있는 블럭을 나타내는 포인터 block 그리고 소유 상태를 나타내는 owner 가 있습니다. 이 구조체는 간단히 다음과 같이 정의되어 있습니다.

type gsl_vector

typedef struct
{
  size_t size;
  size_t stride;
  double * data;
  gsl_block * block;
  int owner;
} gsl_vector;

size 간단히 벡터의 원소 수를 나타냅니다. 원소에 접근할 수 있는 인덱스의 범주는 \(0\) 에서 size-1 까지 입니다. stride 물리적 메모리에서 한 개의 원소에서 다음 원소로 접근할 때의 보폭 값입니다. 이는 적절한 데이터 형의 단위로 정해집니다. 포인터 data 는 메모리에서 첫번째 벡터 원소를 가리킵니다. 포인터 block 은 벡터 원소들이 위치한 메모리 블럭의 주소를 가리킵니다. 벡터가 해당 블럭을 소유하고 있다면, owner 는 \(1\) 로 설정됩니다. 그리고 벡터가 해제될 때, 블럭도 해제됩니다. 다른 객체가 블럭을 소유하고 있다면, 벡터의 owner 는 \(0\) 으로 설정되고 메모리에서 벡터가 해제되어도, 블럭은 해제되지 않습니다.

벡터를 할당하고 접근하는 함수들은 gsl_vector.h 에 정의되어 있습니다.

벡터 할당 (Vector allocation)

벡터를 메모리에 할당하는 함수들은 malloc , free 와 같은 형식을 따릅니다. 이에 더해 각자 고유의 오류 확인 기능을 수행합니다. 메모리 공간이 충분하지 않아 벡터를 할당할 수 없다면, 함수들은 GSL 오류 관리자를 불러 GSL_ENOMEM 오류 값를 넘기고 NULL 포인터를 반환합니다. 때문에 프로그램을 멈추는 데 오류 관리자를 사용한다면 모든 alloc 과정을 확인할 필요가 없습니다.

gsl_vector *gsl_vector_alloc(size_t n)

길이 \(n\) 의 벡터를 만들고, 벡터 구조체의 시작 지점을 가르키는 포인터를 반환합니다. 벡터의 원소로 새로운 블럭이 할당되고 벡터 구조체의 block 변수에 저장됩니다. 블럭은 벡터에 “소유된” 상태입니다. 만약, 벡터가 해제된다면 블럭도 같이 해제됩니다. 크기가 \(0\) 인 요청도 수행되며 non-null 결과를 반환합니다.

gsl_vector *gsl_vector_calloc(size_t n)

길이 \(n\) 의 벡터를 만들고, 벡터 구조체의 시작 지점을 가르키는 포인터를 반환합니다. 벡터의 모든 원소들은 \(0\) 으로 초기화 되어있습니다.

void gsl_vector_free(gsl_vector *v)

할당된 벡터 v 를 해제합니다. 벡터가 gsl_vector_alloc() 를 통해 만들어졌다면 벡터 안의 블럭도 같이 해제됩니다. 만약, 벡터가 메모리의 다른 객체로부터 만들어졌고 다른 객체가 소유권을 가지고 있다면 블럭은 해제되지 않습니다.

벡터 원소에 접근

포트란 컴파일러들과는 다르게 C 컴파일러들은 벡터와 행렬들에 대해 범위 확인 기능을 지원하지 않습니다 1 . 함수 gsl_vector_get() gsl_vector_set() 간단한 범위 확인 기능을 제공하고, 범위 밖 요소에 접근할 경우 오류를 보고합니다.

벡터와 행렬의 원소에 접근하는 함수들은 gsl_vector.h 에 정의 되어 있으며, extern inline 로 인라인 선언 되어 있습니다. 이러한 구조는 함수 호출 과정에서 오버헤드를 줄여줍니다. 이를 위해서는 프로그램을 컴파일 할 때, 전처리기에서 매크로 HAVE_INLINE 처리 하도록 선언해야합니다.

GSL_RANGE_CHECK_OFF

필요한 경우 GSL_RANGE_CHECK_OFF 정의해 재컴파일하면, 코드의 수정 없이 모든 범위 확인 기능을 무시할 수 있습니다. 인라인 함수를 지원하는 컴파일러에서 이러한 범위 확인 기능을 끄는 것은 gsl_vector_get(v,i) 와 gsl_vector_set(v,i,x) 를 각각 v->data[i*v->stride] 와 v->data[i*v->stride]=x 로 바꾸는 것과 같습니다. 따라서, 범위 확인 기능을 해제한 상황에서 범위 확인 기능이 있는 함수를 사용해도 성능상의 하락은 없습니다.

GSL_C99_INLINE

C99 컴파일러를 사용한다면 extern inline 대신 inline 키워드를 사용해 인라인 함수를 정의합니다. 이 경우 매크로 GSL_C99_INLINE 를 정의해야합니다. GCC의 C99 모드 -std=c99 에서는 자동으로 선택됩니다.

int gsl_check_range

만약 인라인 함수들을 사용하지 않으면, gsl_vector_get() 와 gsl_vector_set() 함수는 라이브러리에 존재하는 컴파일된 버전으로 연결됩니다. 이러한 경우, 범위 확인 기능은 전역 정수 변수 gsl_check_range 에 의해 제어됩니다. 기본적으로 이 옵션은 활성화 되어있으며 gsl_check_range \(0\) 으로 설정해 기능을 끌 수 있습니다. 함수 호출 과정에서 오버헤드의 문제 때문에, 인라인 함수와 비교해서 범위 확인 기능을 끄는 것은 권장되지 않습니다.

double gsl_vector_get(const gsl_vector *v, const size_t i)

벡터 v 의 i 번째 원소를 반환합니다. 만약, i 가 \(0\) 에서 size-1 밖에 있다면, 오류 관리자가 실행되고 \(0\) 이 반환됩니다. HAVE_INLINE 이 정의된 경우, 인라인 함수가 사용됩니다.

void gsl_vector_set(gsl_vector *v, const size_t i, double x)

벡터 v 의 i 번째 원소를 x 로 설정합니다. 만약, i 가 \(0\) 에서 size-1 밖에 있다면, 오류 관리자가 실행됩니다. HAVE_INLINE 이 정의된 경우, 인라인 함수가 사용됩니다.

double *gsl_vector_ptr(gsl_vector *v, size_t i)

const double *gsl_vector_const_ptr(const gsl_vector *v, size_t i)

벡터 v 의 i 번째 원소를 가르키는 포인터를 반환합니다. 만약, i 가 \(0\) 에서 size-1 밖에 있다면, 오류 관리자가 실행되고 NULL 포인터가 반환됩니다. HAVE_INLINE 이 정의된 경우, 인라인 함수가 사용됩니다.

벡터 원소 초기화

void gsl_vector_set_all(gsl_vector *v, double x)

벡터 v 의 모든 원소를 주어진 값 x 으로 초기화합니다.

void gsl_vector_set_zero(gsl_vector *v)

벡터 v 의 모든 원소를 \(0\) 으로 초기화합니다.

int gsl_vector_set_basis(gsl_vector *v, size_t i)

벡터 v 의 원소중 i 번째 원소를 \(1\) 로, 나머지 원소를 \(0\) 으로 초기화해, 기저 벡터를 만듭니다.

벡터 읽고 쓰기

이 라이브러리에서는 이진 값이나 규격화 된 문서 파일들을 벡터에 읽고 쓸 수 있는 함수들을 제공합니다.

int gsl_vector_fwrite(FILE *stream, const gsl_vector *v)

벡터 v 의 원소들을 스트림 stream 에 이진 형식의 내용으로 작성합니다. 작성에 오류가 없으면 \(0\) 값을 반환하고, 작성에서 오류가 발생했을 시 GSL_EFAILED 값을 반환합니다. 해당 값들은 작동하는 시스템의 이진 형식을 따르기 때문에, 다른 시스템끼리의 호환성은 보장되지 않습니다.

int gsl_vector_fread(FILE *stream, gsl_vector *v)

열린 스트림 stream 의 값들을 읽어 주어진 벡터 v 에 이진 형태의 값으로 저장합니다. 벡터 v 는 사전에 할당되어 있어야하며, 읽을 값의 크기에 맞추어 적절한 길이를 가지고 있어야합니다. 왜냐하면, 이 함수가 값을 읽는 길이는 v 의 길이에 의존하기 때문입니다. 작성에 오류가 없으면 \(0\) 값을 반환하고, 작성에서 오류가 발생했을 시 GSL_EFAILED 값을 반환합니다. 읽을 값은 같은 시스템에서 사용하는 이진 형식으로 작성되었다라 가정합니다.

int gsl_vector_fprintf(FILE *stream, const gsl_vector *v, const char *format)

주어진 벡터 v 의 원소들을 정해진 형식 format 로 스트림 stream 에 한 줄씩 기록합니다. format 으로 부동 소수점을 위한 %g , %e , 그리고 %f 가 있고, 정수형은 %d 를 사용할 수 있습니다. 작성에 오류가 없으면 \(0\) 값을 반환하고, 작성에서 오류가 발생했을 시 GSL_EFAILED 값을 반환합니다.

int gsl_vector_fscanf(FILE *stream, gsl_vector *v)

stream 로부터 형식화된 값을 읽어 벡터 v 에 저장합니다. 벡터 v 는 사전에 할당되어 있어야하며, 읽을 값의 크기에 맞추어 적절한 길이를 가지고 있어야합니다. 왜냐하면, 이 함수가 값을 읽는 길이는 v 의 길이에 의존하기 때문입니다. 작성에 오류가 없으면 \(0\) 값을 반환하고, 작성에서 오류가 발생했을 시 GSL_EFAILED 값을 반환합니다.

벡터 view

블럭을 나누어 벡터를 만들었다시피, 벡터를 나누어 벡터 view를 만들 수도 있습니다. 예를 들어, 다른 벡터의 부분 벡터를 벡터 view로 표현할 수 있으며, 벡터의 홀수 번째 원소들과 짝수 번째 원소들을 두 개의 벡터 view로 표현할 수도 있습니다.

type gsl_vector_view

type gsl_vector_const_view

벡터 view는 stack에 저장되는 임시 개체로 벡터 원소들의 부분 집합에 적용할 수 있는 연산입니다. 이러한 벡터 view는 상수 벡터와 비상수 벡터 모두에 적용할 수 있습니다. 각각의 성질 보존을 위해 자료형을 분리해 제공합니다. 일반적인 벡터는 gsl_vector_view , 상수 벡터는 gsl_vector_const_view 자료형을 사용합니다. 두 경우 모두, view의 원소는 view 객체의 vector 성분을 사용해 gsl_vector 에 접근할 수 있습니다. gsl_vector * 와 const gsl_vector * 포인터는 이 성분들에 & 연산을 적용해 얻을 수 있습니다.

이 포인터를 사용할 때는, view가 스코프 안에 있는지 확인하는 것이 매우 중요합니다. 가장 간단한 방법은 포인터를 &view.vector 형태로 사용하고, 다른 변수에 해당 값을 저장하지 않는 것입니다.

gsl_vector_view gsl_vector_subvector(gsl_vector *v, size_t offset, size_t n)

gsl_vector_const_view gsl_vector_const_subvector(const gsl_vector *v, size_t offset, size_t n)

벡터 v 의 부분 벡터에 대한 벡터 view를 반환합니다. 새로운 벡터의 시작 값은 원본 벡터의 시작 지점으로부터 거리 offset 을 이용해 정해집니다. 새 벡터는 n 의 원소를 가지고 있습니다. 수학적으로 새로운 벡터 v' 의 i 번째 원소는 다음과 같이 주어집니다.

v'(i) = v->data[(offset +i)*v -> stride]

i 는 0 에서 n-1 까지 입니다.

반환된 벡터 구조체의 data 포인터는 offset n 계수가 원본 벡터의 끝을 넘을 경우 null 초기화됩니다.

새로운 벡터는 원본 벡터 v 이루는 블럭들의 view일 뿐입니다. b 의 원소를 이루는 블럭들은 새로운 벡터가 소유하고 있지 않습니다. view가 스코프 밖으로 나간다 하더라도, 원본 벡터 v 의 블럭들은 여전히 존재합니다. 원본 벡터가 담긴 메모리는 오직 원본 벡터의 해제로만 해제가 가능합니다. 물론, 원본 벡터도 만약 해당 벡터의 view가 사용되고 있을 때에는 해제가 불가능합니다.

함수 gsl_vector_const_subvector() 는 gsl_vector_subvector() 동일한 기능을 제공하고, const 로 선언된 벡터에 대해서 사용가능합니다.

gsl_vector_view gsl_vector_subvector_with_stride(gsl_vector *v, size_t offset, size_t stride, size_t n)

gsl_vector_const_view gsl_vector_const_subvector_with_stride(const gsl_vector *v, size_t offset, size_t stride, size_t n)

보폭 값 stride 이용해 다른 벡터 v 서 생성된 부분 벡터의 벡터 view를 반환합니다. 이 부분 벡터는 gsl_vector_subvector() 와 같은 방식으로 생성되고 n 의 원소를 가지게 됩니다. 이는 원래 벡터의 한 원소에서 다음 원소까지 stride 간격을 가지도록 추출된 원소들입니다. 수학적으로 새 벡터 v' 의 i 째 원소는 다음과 같이 주어집니다.

v'(i) = v->data[(offset + i*stride)*v->stride]

i 는 \(0\) 부터 n-1 까지의 범위를 가집니다.

유의할 점은 부분 벡터의 view들은 원래 벡터의 원소들에대한 직접 접근을 제공한다는 점입니다.

예를 들어서, 다음의 예제 코드는 길이 n 벡터 b 의 짝수번째 원소들을 \(0\) 으로 만듭니다. 홀수 번째 원소들은 건드리지 않습니다.

gsl_vector_view **v_even = gsl_vector_subvector_with_stride (v, 0, 2, n/2);
gsl_vector_set_zero (&v_even.vector);

벡터 view들은 &view.vector 를 사용해 벡터를 인자로 가지는 모든 함수들에 할당된 다른 벡터처럼 인자로 넣어줄 수 있습니다. 예를 들어 다음의 코드는 벡터 v 홀수 번째 원소들에 대해 BLAS 함수 dnrm2 사용해 노름을 계산합니다.

gsl_vector_view **v_odd = gsl_vector_subvector_with_stride (v, 1, 2, n/2);
double r = gsl_blas_dnrm2 (&v_odd.vector);

함수 gsl_vector_const_subvector_with_stride() 는 gsl_vector_subvector_with_stride() 와 동일한 기능을 제공하고, const 로 정의된 벡터들에서 사용가능합니다.

gsl_vector_view gsl_vector_complex_real(gsl_vector_complex *v)

gsl_vector_const_view gsl_vector_complex_const_real(const gsl_vector_complex *v)

복소수 벡터 v 대해 실수부로 이루어진 벡터 view를 반환합니다.

함수 gsl_vector_complex_const_real() 는 gsl_vector_complex_real() 와 동일한 기능을 제공하고, const 로 정의된 벡터들에서 사용가능합니다.

gsl_vector_view gsl_vector_complex_imag(gsl_vector_complex *v)

gsl_vector_const_view gsl_vector_complex_const_imag(const gsl_vector_complex *v)

복소수 벡터 v 대해 허수부로 이루어진 벡터 view를 반환합니다.

gsl_vector_complex_const_imag() 는 gsl_vector_complex_imag() 와 동일한 기능을 제공하고, const 로 정의된 벡터들에서 사용가능합니다.

gsl_vector_view gsl_vector_view_array(double *base, size_t n)

gsl_vector_const_view gsl_vector_const_view_array(const double *base, size_t n)

배열에 대한 벡터 view를 반환합니다. 새 벡터의 값들은 n 의 원소를 가지는 배열 base 의해 결정됩니다. 수학적으로 새 벡터 v' i 번째 원소는 다음과 같이 주어집니다.

v'(i) = base[i];

i 는 \(0\) 부터 n-1 까지의 범위를 가집니다. v 원소들을 가지는 배열은 새로운 벡터 view에 속해져있지 않습니다. view가 소멸하더라도 배열은 여전히 존재합니다. 해당 배열이 속해있는 메모리는 해당 배열을 가르키는 base 포인터가 소멸할 때 해제됩니다. 물론 반대로 배열이 소멸한 상황에서도 view는 여전히 사용할 수 있습니다.

함수 gsl_vector_const_view_array() 는 gsl_vector_view_array() 와 동일한 기능을 제공하고, const 로 정의된 벡터들에서만 사용가능합니다.

gsl_vector_view gsl_vector_view_array_with_stride(double *base, size_t stride, size_t n)

gsl_vector_const_view gsl_vector_const_view_array_with_stride(const double *base, size_t stride, size_t n)

보폭 값 stride 이용해 배열 base 에서 생성된 벡터 view를 반환합니다. 이 부분 벡터는 gsl_vector_view_array() 같은 방식으로 생성되고 n 의 원소를 가지게 됩니다. 이는 원래 배열의 한 원소에서 다음 원소까지 stride 간격을 가지도록 추출된 원소들입니다. 수학적으로 새 벡터 v' i 째 원소는 다음과 같이 주어집니다.

v'(i) = base[i*stride];

i 는 \(0\) 부터 n-1 까지 범위를 가집니다.

유의할 점은 해당 벡터 view가 본래 배열에 대해 직접적인 접근 방법을 제공한다는 점입니다. 벡터 view들은 &view.vector 를 사용해 벡터를 인자로 가지는 모든 함수들에 할당된 다른 벡터처럼 인자로 넣어줄 수 있습니다.

함수 gsl_vector_const_view_array_with_stride() 는 gsl_vector_view_array_with_stride() 와 동일한 기능을 제공하고, const 로 정의된 벡터들에서만 사용가능합니다.

벡터 복사

덧셈이나 곱셈 등의 일반적인 벡터 연산들은 BLAS 라이브러리에 구현되어 있습니다( BLAS 지원 을 참고할 수 있습니다). 하지만, BLAS 전체 코드가 필요 없는 몇몇 유용한 기능들이 있을 수 있습니다. 다음의 함수들은 이러한 범주에 드는 기능들을 구현한 함수들입니다.

int gsl_vector_memcpy(gsl_vector *dest, const gsl_vector *src)

벡터 src 원소들을 dest 벡터로 복사합니다. 이 두 벡터들은 같은 길이를 가져야합니다.

int gsl_vector_swap(gsl_vector *v, gsl_vector *w)

두 벡터 v 와 w 의 원소들을 교환합니다. 이 두 벡터들은 같은 길이를 가져야합니다.

원소 교환

다음의 함수들은 벡터의 원소들에 대해 교환, 순열 연산을 취하는 데 이용할 수 있습니다.

int gsl_vector_swap_elements(gsl_vector *v, size_t i, size_t j)

벡터 v 의 i 번째와 j 번째 원소들을 교환합니다.

int gsl_vector_reverse(gsl_vector *v)

벡터 v 의 원소 순서를 역으로 바꿉니다.

벡터 연산

int gsl_vector_add(gsl_vector *a, const gsl_vector *b)

벡터 a 와 벡터 b 의 원소를 합합니다. \(a_i \leftarrow a_i + b_i\) 값이 a 벡터에 저장되며 b 는 변하지 않습니다. 이 두 벡터들은 같은 길이를 가져야합니다.

int gsl_vector_sub(gsl_vector *a, const gsl_vector *b)

벡터 a 와 벡터 b 의 원소를 뺍니다. \(a_i \leftarrow a_i - b_i\) 값이 a 벡터에 저장되며 b 는 변하지 않습니다. 이 두 벡터들은 같은 길이를 가져야합니다.

int gsl_vector_mul(gsl_vector *a, const gsl_vector *b)

벡터 a 와 벡터 b 의 원소를 곱합니다. \(a_i \leftarrow a_i * b_i\) 값이 a 벡터에 저장되며 b 는 변하지 않습니다. 이 두 벡터들은 같은 길이를 가져야합니다.

int gsl_vector_div(gsl_vector *a, const gsl_vector *b)

벡터 a 와 벡터 b 의 원소를 나눕니다. \(a_i \leftarrow a_i / b_i\) 값이 a 벡터에 저장되며 b 는 변하지 않습니다. 이 두 벡터들은 같은 길이를 가져야합니다.

int gsl_vector_scale(gsl_vector *a, const double x)

벡터 a 원소들에 상수 계수 x 를 곱합니다. \(a_i \leftarrow x a_i\) 값이 a 에 저장됩니다.

int gsl_vector_add_constant(gsl_vector *a, const double x)

벡터 a 원소들에 상수 x 를 더합니다. \(a_i \leftarrow a_i + x\) 값이 a 에 저장됩니다.

double gsl_vector_sum(const gsl_vector *a)

주어진 벡터 a 원소들의 합 \(\sum_{i=1}^n a_i\) 을 반환합니다.

int gsl_vector_axpby(const double alpha, const gsl_vector *x, const double beta, gsl_vector *y)

\(y \leftarrow \alpha x + \beta y\) 연산을 수행합니다. 벡터 x 와 y 는 반드시 같은 길이를 가져야합니다.

벡터의 최대, 최소 원소 찾기

다음 연산들은 실수 벡터들에 대해서만 사용 가능합니다.

double gsl_vector_max(const gsl_vector *v)

벡터 v 최댓값을 반환합니다.

double gsl_vector_min(const gsl_vector *v)

벡터 v 최솟값을 반환합니다.

void gsl_vector_minmax(const gsl_vector *v, double *min_out, double *max_out)

벡터 v 최댓값과 최솟값을 반환합니다. 해당 값들은 각각 min_out 과 max_out 포인터가 가르키는 공간에 저장됩니다.

size_t gsl_vector_max_index(const gsl_vector *v)

벡터 v 최댓값인 원소의 위치를 반환합니다. 동일한 최댓값이 여럿 존재하는 경우 가장 작은 값이 반환됩니다.

size_t gsl_vector_min_index(const gsl_vector *v)

벡터 v 최솟값인 원소의 위치를 반환합니다. 동일한 최솟값이 여럿 존재하는 경우 가장 작은 값이 반환됩니다.

void gsl_vector_minmax_index(const gsl_vector *v, size_t *imin, size_t *imax)

벡터 v 최댓값과 최솟값의 위치들을 반환합니다. 해당 값들은 각각 imin 과 imax 포인터가 가르키는 공간에 저장됩니다. 최댓값과 최솟값이 여럿 존재하는 경우 각각 가장 작은 값이 반환됩니다.

벡터 성질

다음 함수들은 실수, 복소수 벡터들에대해 모두 정의되어 있습니다. 복소수 벡터들은 실수, 허수부 모두 각각의 함수들의 조건들을 만족해야 합니다. 예를 들어서 복소수 벡터가 뒤의 나오는 함수 gsl_vector_ispos() 에서 참 값이 \(1\) 이 나왔다면 이는 해당 복소수 벡터의 실수부와 허수부 모두 양수라는 것을 의미합니다. 허수부나 실수부 둘 중 하나라도 음수 값이 있다면 해당 함수는 \(0\) 을 반환할 것입니다. (*)

int gsl_vector_isnull(const gsl_vector *v)

int gsl_vector_ispos(const gsl_vector *v)

int gsl_vector_isneg(const gsl_vector *v)

int gsl_vector_isnonneg(const gsl_vector *v)

벡터 v 의 모든 원소들이 \(0\) 이거나, 양수이거나, 음수이거나, 음수가 아닌 상황들을 판정해 참일 경우 \(1\) 을 반환하고 아니면 \(0\) 을 반환합니다.

int gsl_vector_equal(const gsl_vector *u, const gsl_vector *v)

주어진 두 벡터 u v 에 대해, 원소들을 비교해 이 둘이 같으면 \(1\) 을 아니면 \(0\) 을 반환합니다.

벡터 예제

다음 프로그램은 어떻게 벡터를 할당, 초기화하고 읽는 지를 보여줍니다. gsl_vector_alloc() 와 gsl_vector_set() 그리고 gsl_vector_get() 를 사용합니다.

#include <stdio.h>
#include <gsl/gsl_vector.h>

int
main (void)
{
  int i;
  gsl_vector * v = gsl_vector_alloc (3);

  for (i = 0; i < 3; i++)
    {
      gsl_vector_set (v, i, 1.23 + i);
    }

  for (i = 0; i < 100; i++) /* OUT OF RANGE ERROR */
    {
      printf ("v_%d = %g\n", i, gsl_vector_get (v, i));
    }

  gsl_vector_free (v);
  return 0;
}

다음은 프로그램의 출력 결과입니다. 반복문의 마지막 단계에서 함수는 벡터 v 의 길이 이상을 읽으려 시도하고 해당 시도는 gsl_vector_get() 함수의 길이 확인 코드에 의해 감지되어 오류를 반환합니다.

$./a.out
v_0 = 1.23
v_1 = 2.23
v_2 = 3.23
gsl: vector_source.12: ERROR: index out of range
Default GSL error handler invoked.
Aborted (core dumped)

다음 프로그램은 벡터를 어떻게 파일에 기록하는지를 보여줍니다.

#include <stdio.h>
#include <gsl/gsl_vector.h>

int
main (void)
{
  int i;
  gsl_vector * v = gsl_vector_alloc (100);

  for (i = 0; i < 100; i++)
    {
      gsl_vector_set (v, i, 1.23 + i);
    }

  {
     FILE * f = fopen ("test.dat", "w");
     gsl_vector_fprintf (f, v, "%.5g");
     fclose (f);
  }

  gsl_vector_free (v);
  return 0;
}

프로그램을 실행하고 나면 test.dat 파일이 생성됩니다. 이 파일은 v 원소들을 값으로 가지고 있습니다. 해당 값들은 %.5g 형식으로 저장되어져 있습니다. gsl_vector_fscanf(f,v) 를 사용해 다음과 같이 읽을 수 있습니다.

#include <stdio.h>
#include <gsl/gsl_vector.h>

int
main (void)
{
  int i;
  gsl_vector * v = gsl_vector_alloc (10);

  {
     FILE * f = fopen ("test.dat", "r");
     gsl_vector_fscanf (f, v);
     fclose (f);
  }

  for (i = 0; i < 10; i++)
    {
      printf ("%g\n", gsl_vector_get(v, i));
    }

  gsl_vector_free (v);
  return 0;
}

행렬 (Matrices)

참고

번역중

행렬은 gsl_matrix 구조체로 정의됩니다. 이 구조체는 일반화된 블록들의 분할체 입니다. 벡터와 같이 이 구조체는 메모리의 특정 원소 집합을 나타냅니다. 그러나 벡터와는 달리 2개의 인덱스를 사용합니다.

type gsl_matrix

gsl_matrix 구조체는 여섯개의 변수를 가지고 있습니다. 행렬의 두 차원, 물리적 차원, 행렬 원소들이 저장된 메모리를 가르키는 포인터 data , 행렬이 소유한 블럭 block 그리고 소유권을 나타내는 owner 입니다. 물리적 차원은 메모리의 배치를 나타내며 행렬의 차원과는 다를 수 있습니다. 이는 부분 행렬의 사용을 위함합니다. gsl_matrix 구조체는 다음과 같이 간단하게 정의되어 있습니다.

typedef struct
{
  size_t size1;
  size_t size2;
  size_t tda;
  double * data;
  gsl_block * block;
  int owner;
} gsl_matrix;

행렬은 행-우선 순서(row-major order)로 되어 있습니다. 이 의미는 한 행의 원소들은 메모리 내에 연속적으로 배치되어 있음을 의미합니다. 이 방법은 2차원 배열에 대한 표준 “C 언어 순서”입니다. 포트란의 경우 열-우선 순서(column-major order)로 되어 있습니다. 행의 갯수는 size1 로 표기됩니다. 사용가능한 행 인덱스는 \(0\) 에서 size1-1 까지 입니다. 같은 방식으로 열의 갯수는 size2 표기됩니다. 열 인덱스는 \(0\) 에서 size2-1 지 입니다. 물리적인 행 차원은 tda 로 표기됩니다. 이는 trailing dimension 이라고도 불리며, 메모리 내 행렬의 배치에서 행의 크기를 나타냅니다.

예를 들어서 다음 행렬은 size1 가 \(3\) , size2 가 \(4\) , 그리고 tda 기 \(8\) 인 행렬입니다. 메모리 내의 행렬의 물리적 배치는 좌상단에서 부터 시작해 좌에서 우로 각 행들을 따라 값들이 배치됩니다.

00 01 02 03 XX XX XX XX
10 11 12 13 XX XX XX XX
20 21 22 23 XX XX XX XX

메모리 내에서 사용되지 않는 영역은 ” XX “로 표기됩니다. 포인터 data 메모리 내 행렬의 첫번째 원소의 위치를 가르킵니다. 포인터 block 행렬의 원소들이 저장된 메모리 블럭의 위치를 나타냅니다. 만약 행렬이 해당 블럭을 소유하고 있다면, owner 는 \(1\) 로 설정됩니다. 이 경우 블럭은 행렬이 해제되면 같이 해제됩니다. 만약, 행렬이 다른 객체가 소유한 블럭의 조각이라면 owner 는 \(0\) 으로 설정되고 해당 행렬이 가르키는 블럭은 해제되지 않습니다.

행렬을 할당하고 접근하는 함수들은 gsl_matrix.h 헤더 파일에 정의되어 있습니다.

행렬 할당

메모리에 행렬을 할당하는 함수들은 malloc 과 free 형식을 따릅니다. 그 자체로 오류 검증 기능을 가지고 있습니다. 행렬을 할당하기에 적절한 메모리가 없는 경우 GSL 오류 관리자를 호출해 GSL_ENOMEM 오류 값를 전달하고 Null 포인터를 반환합니다. 따라서, 프로그램을 정지시키기기 위해 라이브러리의 오류 관리자를 사용한다면 모든 alloc 과정을 검사할 필요가 없습니다.

gsl_matrix *gsl_matrix_alloc(size_t n1, size_t n2)

This function creates a matrix of size n1 rows by n2 columns, returning a pointer to a newly initialized matrix struct. A new block is allocated for the elements of the matrix, and stored in the block component of the matrix struct. The block is “owned” by the matrix, and will be deallocated when the matrix is deallocated. Requesting zero for n1 or n2 is valid and returns a non-null result.

gsl_matrix *gsl_matrix_calloc(size_t n1, size_t n2)

This function allocates memory for a matrix of size n1 rows by n2 columns and initializes all the elements of the matrix to zero.

void gsl_matrix_free(gsl_matrix *m)

This function frees a previously allocated matrix m. If the matrix was created using :fun`gsl_matrix_alloc` then the block underlying the matrix will also be deallocated. If the matrix has been created from another object then the memory is still owned by that object and will not be deallocated.

행렬 원소 접근

The functions for accessing the elements of a matrix use the same range checking system as vectors. You can turn off range checking by recompiling your program with the preprocessor definition GSL_RANGE_CHECK_OFF .

The elements of the matrix are stored in “C-order”, where the second index moves continuously through memory. More precisely, the element accessed by the function gsl_matrix_get(m,i,j) and gsl_matrix_set(m,i,j,x) is:

m->data[i * m->tda + j]

where tda is the physical row-length of the matrix.

double gsl_matrix_get(const gsl_matrix *m, const size_t i, const size_t j)

This function returns the \((i,j)\)-th element of a matrix m. If i or j lie outside the allowed range of 0 to n1 - 1 and 0 to n2 - 1 then the error handler is invoked and 0 is returned. An inline version of this function is used when HAVE_INLINE is defined.

void gsl_matrix_set(gsl_matrix *m, const size_t i, const size_t j, double x)

This function sets the value of the \((i,j)\)-th element of a matrix m to x. If i or j lies outside the allowed range of 0 to n1 - 1 and 0 to n2 - 1 then the error handler is invoked. An inline version of this function is used when HAVE_INLINE is defined.

double *gsl_matrix_ptr(gsl_matrix *m, size_t i, size_t j)

const double *gsl_matrix_const_ptr(const gsl_matrix *m, size_t i, size_t j)

These functions return a pointer to the \((i,j)\)-th element of a matrix m. If i or j lie outside the allowed range of 0 to n1 - 1 and 0 to n2 - 1 then the error handler is invoked and a null pointer is returned. Inline versions of these functions are used when HAVE_INLINE is defined.

행렬 원소 초기화

void gsl_matrix_set_all(gsl_matrix *m, double x)

This function sets all the elements of the matrix m to the value x.

void gsl_matrix_set_zero(gsl_matrix *m)

This function sets all the elements of the matrix m to zero.

void gsl_matrix_set_identity(gsl_matrix *m)

This function sets the elements of the matrix m to the corresponding elements of the identity matrix, \(m(i,j) = \delta(i,j)\), i.e. a unit diagonal with all off-diagonal elements zero. This applies to both square and rectangular matrices.

행렬 읽고 쓰기

The library provides functions for reading and writing matrices to a file as binary data or formatted text.

int gsl_matrix_fwrite(FILE *stream, const gsl_matrix *m)

This function writes the elements of the matrix m to the stream stream in binary format. The return value is 0 for success and GSL_EFAILED if there was a problem writing to the file. Since the data is written in the native binary format it may not be portable between different architectures.

int gsl_matrix_fread(FILE *stream, gsl_matrix *m)

This function reads into the matrix m from the open stream stream in binary format. The matrix m must be preallocated with the correct dimensions since the function uses the size of m to determine how many bytes to read. The return value is 0 for success and GSL_EFAILED if there was a problem reading from the file. The data is assumed to have been written in the native binary format on the same architecture.

int gsl_matrix_fprintf(FILE *stream, const gsl_matrix *m, const char *format)

This function writes the elements of the matrix m line-by-line to the stream stream using the format specifier format, which should be one of the %g, %e or %f formats for floating point numbers and %d for integers. The function returns 0 for success and GSL_EFAILED if there was a problem writing to the file.

int gsl_matrix_fscanf(FILE *stream, gsl_matrix *m)

This function reads formatted data from the stream stream into the matrix m. The matrix m must be preallocated with the correct dimensions since the function uses the size of m to determine how many numbers to read. The function returns 0 for success and GSL_EFAILED if there was a problem reading from the file.

행렬 view

Writing {: .label .label-red}

type gsl_matrix_view

type gsl_matrix_const_view

A matrix view is a temporary object, stored on the stack, which can be used to operate on a subset of matrix elements. Matrix views can be defined for both constant and non-constant matrices using separate types that preserve constness. A matrix view has the type gsl_matrix_view and a constant matrix view has the type gsl_matrix_const_view . In both cases the elements of the view can by accessed using the matrix component of the view object. A pointer gsl_matrix * or const gsl_matrix * can be obtained by taking the address of the matrix component with the & operator. In addition to matrix views it is also possible to create vector views of a matrix, such as row or column views.

gsl_matrix_view gsl_matrix_submatrix(gsl_matrix *m, size_t k1, size_t k2, size_t n1, size_t n2)

gsl_matrix_const_view gsl_matrix_const_submatrix(const gsl_matrix *m, size_t k1, size_t k2, size_t n1, size_t n2)

These functions return a matrix view of a submatrix of the matrix m. The upper-left element of the submatrix is the element ( k1, k2 ) of the original matrix. The submatrix has n1 rows and n2 columns. The physical number of columns in memory given by tda is unchanged. Mathematically, the \((i,j)\) -th element of the new matrix is given by:

m'(i,j) = m->data[(k1*m->tda + k2) + i*m->tda + j]

where the index i runs from 0 to n1 - 1 and the index j runs from 0 to n2 - 1.

The data pointer of the returned matrix struct is set to null if the combined parameters (i, j, n1, n2, tda) overrun the ends of the original matrix.

The new matrix view is only a view of the block underlying the existing matrix, m. The block containing the elements of m is not owned by the new matrix view. When the view goes out of scope the original matrix m and its block will continue to exist. The original memory can only be deallocated by freeing the original matrix. Of course, the original matrix should not be deallocated while the view is still in use.

The function :fun`gsl_matrix_const_submatrix` is equivalent to :fun`gsl_matrix_submatrix` but can be used for matrices which are declared const.

gsl_matrix_view gsl_matrix_view_array(double *base, size_t n1, size_t n2)

gsl_matrix_const_view gsl_matrix_const_view_array(const double *base, size_t n1, size_t n2)

These functions return a matrix view of the array base . The matrix has n1 rows and n2 columns. The physical number of columns in memory is also given by n2 . Mathematically, the \((i,j)\) -th element of the new matrix is given by:

m'(i,j) = base[i*n2 + j]

where the index i runs from 0 to n1 - 1 and the index j runs from 0 to n2 - 1 .

The new matrix is only a view of the array base . When the view goes out of scope the original array base will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use.

The function :fun`gsl_matrix_const_view_array` is equivalent to :fun`gsl_matrix_view_array` but can be used for matrices which are declared const .

gsl_matrix_view gsl_matrix_view_array_with_tda(double *base, size_t n1, size_t n2, size_t tda)

gsl_matrix_const_view gsl_matrix_const_view_array_with_tda(const double *base, size_t n1, size_t n2, size_t tda)

These functions return a matrix view of the array base with a physical number of columns tda which may differ from the corresponding dimension of the matrix. The matrix has n1 rows and n2 columns, and the physical number of columns in memory is given by tda . Mathematically, the \((i,j)\)-th element of the new matrix is given by:

m'(i,j) = base[i*tda + j]

where the index i runs from 0 to n1 - 1 and the index j runs from 0 to n2 - 1 .

The new matrix is only a view of the array base. When the view goes out of scope the original array base will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use.

The function :fun`gsl_matrix_const_view_array_with_tda` is equivalent to :fun`gsl_matrix_view_array_with_tda` but can be used for matrices which are declared const .

gsl_matrix_view gsl_matrix_view_vector(gsl_vector *v, size_t n1, size_t n2)

gsl_matrix_const_view gsl_matrix_const_view_vector(const gsl_vector *v, size_t n1, size_t n2)

These functions return a matrix view of the vector v. The matrix has n1 rows and n2 columns. The vector must have unit stride. The physical number of columns in memory is also given by n2. Mathematically, the \((i,j)\)-th element of the new matrix is given by:

m'(i,j) = v->data[i*n2 + j]

where the index i runs from 0 to n1 - 1 and the index j runs from 0 to n2 - 1 .

The new matrix is only a view of the vector v. When the view goes out of scope the original vector v will continue to exist. The original memory can only be deallocated by freeing the original vector. Of course, the original vector should not be deallocated while the view is still in use.

The function :fun`gsl_matrix_const_view_vector` is equivalent to :fun`gsl_matrix_view_vector` but can be used for matrices which are declared const .

gsl_matrix_view gsl_matrix_view_vector_with_tda(gsl_vector *v, size_t n1, size_t n2, size_t tda)

gsl_matrix_const_view gsl_matrix_const_view_vector_with_tda(const gsl_vector *v, size_t n1, size_t n2, size_t tda)

These functions return a matrix view of the vector v with a physical number of columns tda which may differ from the corresponding matrix dimension. The vector must have unit stride. The matrix has n1 rows and n2 columns, and the physical number of columns in memory is given by tda. Mathematically, the \((i,j)\) -th element of the new matrix is given by:

m'(i,j) = v->data[i*tda + j]

where the index i runs from 0 to n1 - 1 and the index j runs from 0 to n2 - 1.

The new matrix is only a view of the vector v. When the view goes out of scope the original vector v will continue to exist. The original memory can only be deallocated by freeing the original vector. Of course, the original vector should not be deallocated while the view is still in use.

The function :fun`gsl_matrix_const_view_vector_with_tda` is equivalent to :fun`gsl_matrix_view_vector_with_tda` but can be used for matrices which are declared const.

행, 열 view 생성

Writing {: .label .label-red}

In general there are two ways to access an object, by reference or by copying. The functions described in this section create vector views which allow access to a row or column of a matrix by reference. Modifying elements of the view is equivalent to modifying the matrix, since both the vector view and the matrix point to the same memory block.

gsl_vector_view gsl_matrix_row(gsl_matrix *m, size_t i)

gsl_vector_const_view gsl_matrix_const_row(const gsl_matrix *m, size_t i)

These functions return a vector view of the i-th row of the matrix m. The data pointer of the new vector is set to null if i is out of range.

The function :fun`gsl_matrix_const_row` is equivalent to :fun`gsl_matrix_row` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_column(gsl_matrix *m, size_t j)

gsl_vector_const_view gsl_matrix_const_column(const gsl_matrix *m, size_t j)

These functions return a vector view of the j-th column of the matrix m. The data pointer of the new vector is set to null if j is out of range.

The function :fun`gsl_matrix_const_column` is equivalent to :fun`gsl_matrix_column` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_subrow(gsl_matrix *m, size_t i, size_t offset, size_t n)

gsl_vector_const_view gsl_matrix_const_subrow(const gsl_matrix *m, size_t i, size_t offset, size_t n)

These functions return a vector view of the i-th row of the matrix m beginning at offset elements past the first column and containing n elements. The data pointer of the new vector is set to null if i, offset, or n are out of range.

The function :fun`gsl_matrix_const_subrow` is equivalent to :fun`gsl_matrix_subrow` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_subcolumn(gsl_matrix *m, size_t j, size_t offset, size_t n)

gsl_vector_const_view gsl_matrix_const_subcolumn(const gsl_matrix *m, size_t j, size_t offset, size_t n)

These functions return a vector view of the j-th column of the matrix m beginning at offset elements past the first row and containing n elements. The data pointer of the new vector is set to null if j, offset, or n are out of range.

The function :fun`gsl_matrix_const_subcolumn` is equivalent to :fun`gsl_matrix_subcolumn` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_diagonal(gsl_matrix *m)

gsl_vector_const_view gsl_matrix_const_diagonal(const gsl_matrix *m)

These functions return a vector view of the diagonal of the matrix m. The matrix m is not required to be square. For a rectangular matrix the length of the diagonal is the same as the smaller dimension of the matrix.

The function :fun`gsl_matrix_const_diagonal` is equivalent to :fun`gsl_matrix_diagonal` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_subdiagonal(gsl_matrix *m, size_t k)

gsl_vector_const_view gsl_matrix_const_subdiagonal(const gsl_matrix *m, size_t k)

These functions return a vector view of the k-th subdiagonal of the matrix m. The matrix m is not required to be square. The diagonal of the matrix corresponds to \(k = 0\).

The function :fun`gsl_matrix_const_subdiagonal` is equivalent to :fun`gsl_matrix_subdiagonal` but can be used for matrices which are declared const.

gsl_vector_view gsl_matrix_superdiagonal(gsl_matrix *m, size_t k)

gsl_vector_const_view gsl_matrix_const_superdiagonal(const gsl_matrix *m, size_t k)

These functions return a vector view of the k-th superdiagonal of the matrix m . The matrix m is not required to be square. The diagonal of the matrix corresponds to \(k = 0\) .

The function :fun`gsl_matrix_const_superdiagonal` is equivalent to :fun`gsl_matrix_superdiagonal` but can be used for matrices which are declared const .

행렬 복사

Writing {: .label .label-red}

int gsl_matrix_memcpy(gsl_matrix *dest, const gsl_matrix *src)

This function copies the elements of the matrix src into the matrix dest. The two matrices must have the same size.

int gsl_matrix_swap(gsl_matrix *m1, gsl_matrix *m2)

This function exchanges the elements of the matrices m1 and :da

행, 열 복사

Writing {: .label .label-red}

The functions described in this section copy a row or column of a matrix into a vector. This allows the elements of the vector and the matrix to be modified independently. Note that if the matrix and the vector point to overlapping regions of memory then the result will be undefined. The same effect can be achieved with more generality using :fun`gsl_vector_memcpy` with vector views of rows and columns.

int gsl_matrix_get_row(gsl_vector *v, const gsl_matrix *m, size_t i)

This function copies the elements of the i-th row of the matrix m into the vector v. The length of the vector must be the same as the length of the row.

int gsl_matrix_get_col(gsl_vector *v, const gsl_matrix *m, size_t j)

This function copies the elements of the j-th column of the matrix m into the vector v. The length of the vector must be the same as the length of the column.

int gsl_matrix_set_row(gsl_matrix *m, size_t i, const gsl_vector *v)

This function copies the elements of the vector v into the i-th row of the matrix m. The length of the vector must be the same as the length of the row.

int gsl_matrix_set_col(gsl_matrix *m, size_t j, const gsl_vector *v)

This function copies the elements of the vector v into the j-th column of the matrix m. The length of the vector must be the same as the length of the column.

행과 열 교환하기

Writing {: .label .label-red}

The following functions can be used to exchange the rows and columns of a matrix.

int gsl_matrix_swap_rows(gsl_matrix *m, size_t i, size_t j)

This function exchanges the i-th and j-th rows of the matrix m in-place.

int gsl_matrix_swap_columns(gsl_matrix *m, size_t i, size_t j)

This function exchanges the i-th and j-th columns of the matrix m in-place.

int gsl_matrix_swap_rowcol(gsl_matrix *m, size_t i, size_t j)

This function exchanges the i-th row and j-th column of the matrix m in-place. The matrix must be square for this operation to be possible.

int gsl_matrix_transpose_memcpy(gsl_matrix *dest, const gsl_matrix *src)

This function makes the matrix dest the transpose of the matrix src by copying the elements of src into dest. This function works for all matrices provided that the dimensions of the matrix dest match the transposed dimensions of the matrix src.

int gsl_matrix_transpose(gsl_matrix *m)

This function replaces the matrix m by its transpose by copying the elements of the matrix in-place. The matrix must be square for this operation to be possible.

int gsl_matrix_complex_conjtrans_memcpy(gsl_matrix *dest, const gsl_matrix *src)

This function makes the matrix dest the conjugate transpose of the matrix src by copying the complex conjugate elements of src into dest. This function works for all complex matrices provided that the dimensions of the matrix dest match the transposed dimensions of the matrix src.

행렬 연산

Writing {: .label .label-red}

The following operations are defined for real and complex matrices.

int gsl_matrix_add(gsl_matrix *a, const gsl_matrix *b)

This function adds the elements of matrix b to the elements of matrix a. The result \(a(i,j) \leftarrow a(i,j) + b(i,j)\) is stored in a and b remains unchanged. The two matrices must have the same dimensions.

int gsl_matrix_sub(gsl_matrix *a, const gsl_matrix *b)

This function subtracts the elements of matrix b from the elements of matrix a. The result \(a(i,j) \leftarrow a(i,j) - b(i,j)\) is stored in a and b remains unchanged. The two matrices must have the same dimensions.

int gsl_matrix_mul_elements(gsl_matrix *a, const gsl_matrix *b)

This function multiplies the elements of matrix a by the elements of matrix b. The result \(a(i,j) \leftarrow a(i,j) * b(i,j)\) is stored in a and b remains unchanged. The two matrices must have the same dimensions.

int gsl_matrix_div_elements(gsl_matrix *a, const gsl_matrix *b)

This function divides the elements of matrix a by the elements of matrix b. The result \(a(i,j) \leftarrow a(i,j) / b(i,j)\) is stored in a and b remains unchanged. The two matrices must have the same dimensions.

int gsl_matrix_scale(gsl_matrix *a, const double x)

This function multiplies the elements of matrix a by the constant factor x. The result \(a(i,j) \leftarrow x a(i,j)\) is stored in a.

int gsl_matrix_scale_columns(gsl_matrix *A, const gsl_vector *x)

This function scales the columns of the \(M\)-by-\(N\) matrix A by the elements of the vector x, of length \(N\). The \(j\) -th column of A is multiplied by \(x_j\). This is equivalent to forming

\[A \rightarrow A X\]

where \(X = \textrm{diag}(x)\).

int gsl_matrix_scale_rows(gsl_matrix *A, const gsl_vector *x)

This function scales the rows of the \(M\)-by-\(N\) matrix A by the elements of the vector x, of length \(M\). The \(i\) -th row of A is multiplied by \(x_i\). This is equivalent to forming

\[A \rightarrow X A\]

where \(X = \textrm{diag}(x)\).

int gsl_matrix_add_constant(gsl_matrix *a, const double x)

This function adds the constant value x to the elements of the matrix a. The result \(a(i,j) \leftarrow a(i,j) + x\) is stored in a.

행렬 원소의 최대, 최소 찾기

Writing {: .label .label-red}

The following operations are only defined for real matrices.

double gsl_matrix_max(const gsl_matrix *m)

This function returns the maximum value in the matrix m.

double gsl_matrix_min(const gsl_matrix *m)

This function returns the minimum value in the matrix m.

void gsl_matrix_minmax(const gsl_matrix *m, double *min_out, double *max_out)

This function returns the minimum and maximum values in the matrix m, storing them in min_out and max_out.

void gsl_matrix_max_index(const gsl_matrix *m, size_t *imax, size_t *jmax)

This function returns the indices of the maximum value in the matrix m, storing them in imax and jmax. When there are several equal maximum elements then the first element found is returned, searching in row-major order.

void gsl_matrix_min_index(const gsl_matrix *m, size_t *imin, size_t *jmin)

This function returns the indices of the minimum value in the matrix m, storing them in imin and jmin. When there are several equal minimum elements then the first element found is returned, searching in row-major order.

void gsl_matrix_minmax_index(const gsl_matrix *m, size_t *imin, size_t *jmin, size_t *imax, size_t *jmax)

This function returns the indices of the minimum and maximum values in the matrix m, storing them in (imin, jmin) and (imax, jmax). When there are several equal minimum or maximum elements then the first elements found are returned, searching in row-major order.

행렬의 성질

다음 함수들은 실수와 복소수 행렬 모두 사용할 수 있습니다. 복소수 행렬의 경우 각 함수에서 검사하는 성질이 실수부와 허수부가 모두 일치해야 해당 행렬이 함수의 검사 조건을 만족한다 간주합니다.

int gsl_matrix_isnull(const gsl_matrix *m)

int gsl_matrix_ispos(const gsl_matrix *m)

int gsl_matrix_isneg(const gsl_matrix *m)

int gsl_matrix_isnonneg(const gsl_matrix *m)

각각 행렬의 모든 원소들이 \(0\) 이거나, 양수이거나, 음수이거나, 음수가 아니면 \(1\) 을 반환하고 아니면 \(0\) 을 반환합니다. 행렬이 positive-definite임을 검사하려면 숄레스키 분해(Cholesky decomposition)를 사용해야 합니다.

int gsl_matrix_equal(const gsl_matrix *a, const gsl_matrix *b)

행렬 a b 원소를 비교해서 같은 행렬이면 \(1\) 을 아니면 \(0\) 을 반환합니다.

double gsl_matrix_norm1(const gsl_matrix *A)

\(m\) - \(n\) 행렬 A 대해, \(1\) -노름값을 반환합니다. 이는 다음과 같이 정의됩니다.

\[||A||_1 = \text{max}_{1\leq j \leq n} \sum_{i=1}^m |A_{ij}|\]

행렬 예제

아래의 프로그램은 함수 gsl_matrix_alloc() 와 gsl_mtrix_set() 그리고 gsl_matrix_get() 를 이용해 행렬의 할당, 초기화, 그리고 값을 읽는 방법을 보여줍니다.

#include <stdio.h>
#include <gsl/gsl_matrix.h>

int
main (void)
{
  int i, j;
  gsl_matrix * m = gsl_matrix_alloc (10, 3);

  for (i = 0; i < 10; i++)
    for (j = 0; j < 3; j++)
      gsl_matrix_set (m, i, j, 0.23 + 100*i + j);

  for (i = 0; i < 100; i++)  /* OUT OF RANGE ERROR */
    for (j = 0; j < 3; j++)
      printf ("m(%d,%d) = %g\n", i, j,
              gsl_matrix_get (m, i, j));

  gsl_matrix_free (m);

  return 0;
}

다음은 프로그램의 출력 결과입니다. 마지막 반복문은 행렬 m 범위를 넘어 읽으려 시도하고, gsl_matrix_get() 의 범위 확인 코드가 이를 감지했습니다.

$./a.out
 m(0,0) = 0.23
 m(0,1) = 1.23
 m(0,2) = 2.23
 m(1,0) = 100.23
 m(1,1) = 101.23
 m(1,2) = 102.23
 ...
 m(9,2) = 902.23
 gsl: matrix_source.13: ERROR: first index out of range
 Default GSL error handler invoked.
 Aborted (core dumped)

다음 프로그램은 행렬을 어떻게 파일로 저장하는지 보여줍니다.

#include <stdio.h> #include <gsl/gsl_matrix.h>

int main (void) {

int i, j, k = 0; gsl_matrix * m = gsl_matrix_alloc (100, 100); gsl_matrix * a = gsl_matrix_alloc (100, 100);

for (i = 0; i < 100; i++)

for (j = 0; j < 100; j++)

gsl_matrix_set (m, i, j, 0.23 + i + j);

{

FILE * f = fopen (“test.dat”, “wb”); gsl_matrix_fwrite (f, m); fclose (f);

}

{

FILE * f = fopen (“test.dat”, “rb”); gsl_matrix_fread (f, a); fclose (f);

}

for (i = 0; i < 100; i++)

for (j = 0; j < 100; j++)

{

double mij = gsl_matrix_get (m, i, j); double aij = gsl_matrix_get (a, i, j); if (mij != aij) k++;

}

gsl_matrix_free (m); gsl_matrix_free (a);

printf (“differences = %d (should be zero)n”, k); return (k > 0);

}

프로그램의 실행이 끝나면, m 행렬의 원소들을 담고 있는 파일 test.dat 파일이 생성됩니다. 이 파일은 바이너리 형식으로 작성됩니다. 함수 gsl_matrix_fread() 사용해 파일을 읽어 본래 행렬과 같은 행렬을 얻을 수 있습니다.

다음 프로그램은 벡터 view를 사용하는 법을 기술합니다. 이 프로그램은 행렬의 열 노름을 계산합니다.

#include <math.h>
#include <stdio.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_blas.h>

int
main (void)
{
  size_t i,j;

  gsl_matrix *m = gsl_matrix_alloc (10, 10);

  for (i = 0; i < 10; i++)
    for (j = 0; j < 10; j++)
      gsl_matrix_set (m, i, j, sin (i) + cos (j));

  for (j = 0; j < 10; j++)
    {
      gsl_vector_view column = gsl_matrix_column (m, j);
      double d;

      d = gsl_blas_dnrm2 (&column.vector);

      printf ("matrix column %zu, norm = %g\n", j, d);
    }

  gsl_matrix_free (m);

  return 0;
}

결과는 다음과 같습니다.

matrix column 0, norm = 4.31461
matrix column 1, norm = 3.1205
matrix column 2, norm = 2.19316
matrix column 3, norm = 3.26114
matrix column 4, norm = 2.53416
matrix column 5, norm = 2.57281
matrix column 6, norm = 4.20469
matrix column 7, norm = 3.65202
matrix column 8, norm = 2.08524
matrix column 9, norm = 3.07313

GNU octave를 이용해 이 결과를 검증해볼 수 있습니다.

$octave
GNU Octave, version 2.0.16.92
octave> m = sin(0:9)' * ones(1,10)
               + ones(10,1) * cos(0:9);
octave> sqrt(sum(m.^2))
ans =
  4.3146  3.1205  2.1932  3.2611  2.5342  2.5728
  4.2047  3.6520  2.0852  3.0731

참고 문헌과 추가 자료

GSL의 블럭, 벡터 그리고 행렬 객체들은 C++의 valarraay 모델을 따릅니다. 이 모델의 자세한 기술은 다음을 참고할 수 있습니다.

• B. Stroustrup, The C++ Programming Language (3rd Ed), Section 22.4 Vector Arithmetic. Addison-Wesley 1997, ISBN 0-201-88954-4.

1

GNC C 컴파일러에서는 확장기능을 통해 범주 확인 기능을 지원합니다. 하지만, GCC의 기본 설치 목록에 포함되어 있지는 않습니다. 메모리 접근은 GCC의 gcc -fmudflap 라는 메모리 보호 옵션과 Valgrind로 확인할 수 있습니다.