다차원 배열과 할당 배열 - 열 우선 메모리

이 장에서 배우는 것

앞 장에서는 1차원 배열을 반복문 없이 통째로 계산하는 방법을 익혔다. 이 장에서는 배열의 차원을 늘려 표 모양의 자료를 다룬다. 관측 자료는 대개 "날짜 × 관측 시각"처럼 두 가지 기준으로 정리된다. Fortran 은 이런 표를 2차원 배열로 표현하고, 행렬 계산에 쓰는 내장 함수도 갖추고 있다.

이 장의 핵심은 2차원 배열이 메모리에 놓이는 방식이다. Fortran 은 첫 번째 첨자가 가장 빨리 변하는 순서로 값을 저장한다. 이 사실을 알면 반복문을 어떤 순서로 짜야 하는지, reshape 가 값을 어떤 순서로 채우는지가 한꺼번에 설명된다. 크기를 프로그램이 도는 중에 정하는 할당 배열도 함께 다룬다.

2차원 배열을 선언하고 temp(i, j) 형태로 원소와 행·열 단위 부분을 읽는다.
열 우선 저장 방식과, 그에 맞는 이중 반복문의 순서를 설명한다.
reshape, transpose, matmul 로 배열의 모양을 바꾸고 곱한다.
allocatable 배열을 allocate 로 만들고 deallocate 로 돌려준다.
size 와 shape 로 배열의 크기를 묻는다.

문제 상황

작은 기상 관측소가 사흘 동안 하루 네 번, 06시·12시·18시·24시에 기온을 쟀다. 기록은 12개의 숫자다. 이를 1차원 배열 하나에 담으면 "두 번째 날 18시 값"이 몇 번째 원소인지 매번 계산해야 한다. 날짜별 평균, 시각별 평균, 평균에서 얼마나 벗어났는지를 구하려 해도 어느 구간을 합칠지 일일이 따져야 한다.

자료를 3행 4열 표로 두면 "행은 날짜, 열은 관측 시각"이라는 구조가 코드에 그대로 드러난다. 날짜별 평균은 행마다 더하고, 시각별 평균은 열마다 더하면 된다. 다만 표를 컴퓨터가 어떤 순서로 저장하는지 모르면 읽기와 쓰기는 맞아도 느린 코드를 짜게 된다. 또한 관측 일수를 파일에서 알게 되는 경우에는 배열 크기를 컴파일할 때 정해 둘 수 없다. 이 두 문제를 이 장에서 풀어 본다.

2차원 배열과 열 우선 저장

선언과 접근

2차원 배열은 크기를 쉼표로 구분해 선언한다. real(dp) :: temp(3, 4) 는 3행 4열, 곧 원소 12개짜리 표다. 첫 번째 첨자 i 는 행(row), 두 번째 첨자 j 는 열(column)의 위치다. 이 장에서는 행을 날짜, 열을 관측 시각으로 쓴다. temp(2, 3) 은 2일째 세 번째 시각의 값이다.

앞 장의 부분 배열 문법은 차원마다 적용된다. temp(2, :) 는 2일째 하루 전체(행 하나), temp(:, 3) 은 세 번째 시각의 사흘치(열 하나)다. 콜론은 "그 차원 전부"를 뜻한다.

열 우선(column-major) 저장

메모리는 한 줄로 이어진 칸들이라 표를 그대로 놓을 수 없다. 어느 방향으로 펼칠지 정해야 하는데, Fortran 은 첫 번째 첨자가 가장 빨리 변하는 순서로 펼친다. temp(1,1), temp(2,1), temp(3,1) 이 먼저 놓이고, 첫 열이 끝나야 temp(1,2) 로 넘어간다. 이를 열 우선이라 부른다. 그림 1은 3행 4열 표의 칸마다 메모리 순번을 적은 것이다.

그림 2는 같은 12칸을 메모리 순서대로 한 줄로 늘어놓은 모습이다. 같은 열의 세 칸이 붙어 있고, 같은 행의 칸은 세 칸씩 떨어져 있다.

루프 순서

컴퓨터는 메모리에서 값을 하나씩 가져오지 않고 근처 칸을 묶어서 가져온다. 따라서 이웃한 칸을 연달아 쓰는 반복문이 빠르다. 열 우선 저장에서는 안쪽 반복문이 첫 번째 첨자를 돌아야 메모리를 순서대로 훑는다.

do j = 1, nslots        ! 바깥: 두 번째 첨자
   do i = 1, ndays      ! 안쪽: 첫 번째 첨자
      anomaly(i, j) = temp(i, j) - daily_mean(i)
   end do
end do

3×4 표에서는 어느 순서든 차이를 느낄 수 없다. 격자가 수만 칸 이상으로 커지면 순서가 실행 시간에 영향을 준다. 결과 값은 어느 순서든 같다는 점도 알아 두자. 순서가 영향을 주는 것은 속도뿐이다. 다른 언어에서 행 우선 습관을 들였다면 이 점을 특히 바꿔야 한다.

size 와 shape 로 크기 묻기

배열의 크기는 직접 외워 두지 않고 내장 함수로 묻는다. 코드 곳곳에 숫자 3과 4를 적어 두면 표 크기가 바뀔 때 모두 고쳐야 하기 때문이다.

3행 4열 배열 temp 에 대한 크기 질의 결과
호출	값	의미
`size(temp)`	12	전체 원소 수
`size(temp, 1)`	3	1번 차원(행) 길이
`size(temp, 2)`	4	2번 차원(열) 길이
`shape(temp)`	[3, 4]	차원별 길이를 담은 정수 배열

모양을 바꾸고 곱하는 내장 함수

reshape

reshape(원본, 새_모양) 은 원소 순서는 그대로 두고 모양만 바꾼 배열을 돌려준다. 이때 새 배열은 열 우선 순서로 채워진다. 원본 값 12개를 [3, 4] 모양으로 바꾸면 처음 3개가 첫 열이 된다. 행 단위로 값을 적어 놓은 표를 옮길 때는 이 순서를 뒤집어야 하는데, 이는 뒤의 "자주 틀리는 것"에서 다룬다.

transpose

transpose(a) 는 행과 열을 맞바꾼 배열을 돌려준다. 3행 4열 표는 4행 3열이 되고 b(j, i) 가 a(i, j) 와 같은 값이다. "행은 시각, 열은 날짜"인 표가 필요할 때 쓴다.

matmul

matmul(a, b) 는 행렬 곱이다. 결과의 (i, j) 원소는 a 의 i번째 행과 b 의 j번째 열을 원소끼리 곱해 더한 값이다. 그래서 a 의 열 수와 b 의 행 수가 같아야 한다. 3행 4열 표에 길이 4짜리 가중치 벡터를 곱하면 길이 3짜리 벡터가 나온다. 날짜마다 시각별 가중치를 곱해 더한 값, 곧 가중 평균을 한 번에 구하는 셈이다. 여기서 * 는 같은 위치끼리 곱하는 연산일 뿐 행렬 곱이 아니다.

할당 배열

관측 일수가 파일마다 다르다면 배열 크기를 컴파일할 때 정할 수 없다. 이럴 때 할당 배열(allocatable array)을 쓴다. 선언에서는 차원 수만 정하고 각 차원은 콜론으로 비워 둔다.

real(dp), allocatable :: anomaly(:, :)

이 시점의 anomaly 는 메모리가 없는 빈 상태다. 크기가 정해졌을 때 allocate 로 메모리를 잡는다.

allocate(anomaly(size(temp, 1), size(temp, 2)), stat=ierr)
if (ierr /= 0) stop 'allocation failed'

stat= 를 붙이면 메모리가 모자랄 때 프로그램이 곧바로 멈추지 않고 0이 아닌 값이 변수에 들어온다. 쓰고 난 배열은 deallocate 로 돌려준다. 할당 배열은 프로그램 단위가 끝나면 자동으로 해제되지만, 큰 배열을 중간에 정리하고 싶을 때 직접 호출한다. 지금 할당되어 있는지는 allocated(배열) 이 논리값으로 알려 준다.

대입으로도 할당된다. 왼쪽이 할당 배열이면 tt = transpose(temp) 처럼 모양이 맞는 값을 대입할 때 알아서 알맞은 크기로 할당된다. 이 장의 코드에서도 tt 와 weighted 는 이 방식을 쓴다. 한편 a(:) = ... 처럼 첨자를 적은 대입은 이미 할당된 배열을 전제로 한다는 점에 주의하자.

완성 코드

사흘 × 네 시각 기온 표로 날짜별 평균, 시각별 평균, 평균과의 차를 계산하고, 전치와 가중 평균까지 구한다. 한 파일로 되어 있다.

main.f90

program station_grid
  implicit none
  integer, parameter :: dp = selected_real_kind(15, 307)
  integer, parameter :: ndays = 3, nslots = 4
  real(dp) :: temp(ndays, nslots)
  real(dp) :: weights(nslots)
  real(dp), allocatable :: daily_mean(:), slot_mean(:), weighted(:)
  real(dp), allocatable :: anomaly(:, :), tt(:, :)
  integer :: order(ndays, nslots)
  integer :: i, j, k, ierr

  ! 값은 열 우선 순서로 나열한다: 같은 시각의 사흘치가 먼저 온다
  temp = reshape([11.0_dp, 12.5_dp, 10.0_dp, &
                  18.0_dp, 20.5_dp, 16.0_dp, &
                  16.0_dp, 19.0_dp, 14.0_dp, &
                  12.0_dp, 13.0_dp, 11.0_dp], [ndays, nslots])

  print '(a)', '[1] temp(day, slot)'
  write (*, '(5x,4a8)') '06h', '12h', '18h', '24h'
  do i = 1, ndays
    write (*, '("day ",i1,4f8.2)') i, temp(i, :)
  end do

  print '(a)', '[2] size and shape'
  write (*, '("size=",i0," rows=",i0," cols=",i0)') &
    size(temp), size(temp, 1), size(temp, 2)
  write (*, '("shape=",2i3)') shape(temp)

  print '(a)', '[3] column-major order'
  k = 0
  fill_col: do j = 1, nslots
    fill_row: do i = 1, ndays
      k = k + 1
      order(i, j) = k
    end do fill_row
  end do fill_col
  write (*, '("fill order:",12i3)') order
  write (*, '("temp flat:",12f6.1)') temp
  do i = 1, ndays
    write (*, '("day ",i1,4i3)') i, order(i, :)
  end do

  print '(a)', '[4] means'
  allocate (daily_mean(ndays), slot_mean(nslots), stat=ierr)
  if (ierr /= 0) stop 'allocation failed'
  daily_mean = sum(temp, dim=2) / real(nslots, dp)
  slot_mean = sum(temp, dim=1) / real(ndays, dp)
  write (*, '("daily mean:",3f8.2)') daily_mean
  write (*, '("slot mean: ",4f8.2)') slot_mean

  print '(a)', '[5] anomaly = temp - daily mean'
  write (*, '("allocated(anomaly) before allocate: ",l1)') allocated(anomaly)
  allocate (anomaly(size(temp, 1), size(temp, 2)), stat=ierr)
  if (ierr /= 0) stop 'allocation failed'
  write (*, '("allocated(anomaly) after allocate: ",l1)') allocated(anomaly)
  an_col: do j = 1, nslots
    an_row: do i = 1, ndays
      anomaly(i, j) = temp(i, j) - daily_mean(i)
    end do an_row
  end do an_col
  do i = 1, ndays
    write (*, '("day ",i1,4f8.2)') i, anomaly(i, :)
  end do

  print '(a)', '[6] transpose'
  tt = transpose(temp)
  write (*, '("shape(tt)=",2i3)') shape(tt)
  do j = 1, nslots
    write (*, '("slot ",i1,3f8.2)') j, tt(j, :)
  end do

  print '(a)', '[7] matmul'
  weights = [0.1_dp, 0.4_dp, 0.4_dp, 0.1_dp]
  weighted = matmul(temp, weights)
  write (*, '("weighted: ",3f8.2)') weighted
  write (*, '("shape(matmul(tt,temp))=",2i3)') shape(matmul(tt, temp))

  print '(a)', '[8] deallocate'
  deallocate (anomaly, daily_mean, slot_mean, weighted, tt)
  write (*, '("allocated(anomaly) after deallocate: ",l1)') allocated(anomaly)
end program station_grid

줄별 해설

선언부. dp 는 배정밀도 실수의 kind 값이다. ndays 와 nslots 는 이름 붙은 상수로 두어 숫자 3과 4가 코드에 흩어지지 않게 했다. temp 는 고정 크기 2차원 배열이고, 나머지 real(dp) 배열 다섯 개는 allocatable 이라 아직 메모리가 없다. 정수 배열 order 는 열 우선 순서를 눈으로 확인하는 용도다.

reshape 대입. 대괄호 안의 12개 값은 3개씩 한 줄이다. 한 줄이 한 열, 곧 한 관측 시각의 사흘치다. 기온이 날짜 순서로 나열되므로 첫 줄은 06시의 1일·2일·3일 값이다. 줄 끝의 & 는 다음 줄로 이어진다는 표시다.

[1] 표 출력. (5x,4a8) 은 공백 5칸 뒤에 폭 8의 문자열 4개를 쓴다. 문자열이 폭보다 짧으면 오른쪽에 맞춰지므로 아래 숫자와 열이 정렬된다. 행 하나를 temp(i, :) 로 넘기고 4f8.2 로 소수 둘째 자리까지 쓴다.

[2] 크기 묻기. i0 는 자릿수에 맞춰 최소 폭으로 정수를 쓴다. shape(temp) 는 길이 2의 정수 배열이므로 2i3 으로 한 번에 쓴다.

[3] 열 우선 확인. 안쪽 반복문이 i 를 돌며 1부터 번호를 매긴다. 번호가 열 단위로 이어져 들어간다. 이름 붙은 반복문(fill_col, fill_row)은 한 프로그램 안에서 이름이 겹치지 않아야 하므로 뒤의 반복문과 다른 이름을 썼다. write 에 배열 전체를 주면 메모리 순서대로 풀려 나오므로 order 는 1부터 12까지 차례로 출력된다. temp 도 같은 순서로 출력해 reshape 에 적은 값의 순서와 비교할 수 있다.

[4] 평균. sum(temp, dim=2) 는 두 번째 차원을 따라 더해서 길이 3(날짜별) 배열을 만든다. dim=1 은 첫 번째 차원을 따라 더해 길이 4(시각별) 배열을 만든다. 합을 칸 수로 나눌 때는 real(…, dp) 로 실수로 바꿔 두었다. 두 평균 배열은 앞에서 allocate 로 만든 뒤 대입한다. stat=ierr 로 할당 실패를 검사한다.

[5] 평균과의 차. allocated 의 결과가 allocate 를 기준으로 F 에서 T 로 바뀌는 것을 출력한다. 크기는 size(temp, 1) 과 size(temp, 2) 로 정했으므로 temp 의 크기가 바뀌어도 이 줄은 고칠 필요가 없다. 이중 반복문은 안쪽이 i 이므로 열 우선 순서에 맞는다. 하루 평균 daily_mean(i) 는 행 번호로 고른다.

[6] 전치. tt 는 할당하지 않고 transpose 결과를 대입했다. 대입하면서 4×3 모양으로 할당된다. shape(tt) 가 3, 4 가 아닌 4, 3 으로 나오는지 확인하자.

[7] 행렬 곱. 가중치 합이 1이 되도록 낮 시간대에 무게를 두었다. matmul(temp, weights) 는 3×4 와 길이 4 벡터의 곱이라 길이 3 배열이 되어 weighted 에 할당된다. 3행째 줄은 결과 모양만 확인한다. 4×3 인 tt 와 3×4 인 temp 를 곱하면 4×4 가 된다.

[8] 해제. deallocate 에 배열을 여러 개 나열해 한 번에 해제하고, allocated 가 F 가 된 것을 확인한다. 해제된 배열을 다시 쓰려면 allocate 를 새로 해야 한다.

실행 결과

$ gfortran -std=f2018 -Wall main.f90 -o station_grid
$ ./station_grid
[1] temp(day, slot)
          06h     12h     18h     24h
day 1   11.00   18.00   16.00   12.00
day 2   12.50   20.50   19.00   13.00
day 3   10.00   16.00   14.00   11.00
[2] size and shape
size=12 rows=3 cols=4
shape=  3  4
[3] column-major order
fill order:  1  2  3  4  5  6  7  8  9 10 11 12
temp flat:  11.0  12.5  10.0  18.0  20.5  16.0  16.0  19.0  14.0  12.0  13.0  11.0
day 1  1  4  7 10
day 2  2  5  8 11
day 3  3  6  9 12
[4] means
daily mean:   14.25   16.25   12.75
slot mean:    11.17   18.17   16.33   12.00
[5] anomaly = temp - daily mean
allocated(anomaly) before allocate: F
allocated(anomaly) after allocate: T
day 1   -3.25    3.75    1.75   -2.25
day 2   -3.75    4.25    2.75   -3.25
day 3   -2.75    3.25    1.25   -1.75
[6] transpose
shape(tt)=  4  3
slot 1   11.00   12.50   10.00
slot 2   18.00   20.50   16.00
slot 3   16.00   19.00   14.00
slot 4   12.00   13.00   11.00
[7] matmul
weighted:    15.90   18.35   14.10
shape(matmul(tt,temp))=  4  4
[8] deallocate
allocated(anomaly) after deallocate: F

[3] 의 마지막 표는 행 방향으로는 순번이 3씩 뛰고, 열 방향으로는 1씩 늘어남을 보여 준다. 그림 1과 같은 배치다.

실무에서 자주 틀리는 것

1. 안쪽 반복문이 두 번째 첨자를 돈다

다른 언어에서 배운 습관대로 행을 바깥, 열을 안쪽에 두면 값은 맞아도 메모리를 건너뛰며 읽는다.

! 틀린 순서: 안쪽이 j 라서 메모리를 ndays 칸씩 건너뛴다
do i = 1, ndays
  do j = 1, nslots
    anomaly(i, j) = temp(i, j) - daily_mean(i)
  end do
end do

! 고친 순서: 안쪽이 i 라서 메모리를 차례로 읽는다
do j = 1, nslots
  do i = 1, ndays
    anomaly(i, j) = temp(i, j) - daily_mean(i)
  end do
end do

작은 표에서는 차이가 드러나지 않아 발견이 늦다. 격자 크기를 키운 뒤에야 느리다는 것을 알게 되는 경우가 많다.

2. reshape 가 행 순서로 채운다고 생각한다

종이에 가로로 적은 표를 그대로 옮기면 행과 열이 뒤바뀐다.

! 틀린 코드: 날짜별로 4개씩 적었지만 열 우선으로 채워진다
temp = reshape([11.0_dp, 18.0_dp, 16.0_dp, 12.0_dp, &
                12.5_dp, 20.5_dp, 19.0_dp, 13.0_dp, &
                10.0_dp, 16.0_dp, 14.0_dp, 11.0_dp], [ndays, nslots])

값 12개가 3개씩 끊겨 열이 되므로 temp(1, 1), temp(2, 1), temp(3, 1) 에 1일의 06시·12시·18시 값이 들어간다. 고치는 방법은 두 가지다. 한 가지는 order=[2, 1] 로 두 번째 차원이 먼저 변하도록 지정하는 것이다.

temp = reshape([11.0_dp, 18.0_dp, 16.0_dp, 12.0_dp, &
                12.5_dp, 20.5_dp, 19.0_dp, 13.0_dp, &
                10.0_dp, 16.0_dp, 14.0_dp, 11.0_dp], &
               [ndays, nslots], order=[2, 1])

다른 한 가지는 4×3 모양으로 만든 뒤 transpose 하는 것이다.

3. 할당하지 않은 배열을 쓰거나 두 번 할당한다

real(dp), allocatable :: daily_mean(:)
daily_mean(1) = 14.25_dp            ! 틀림: 메모리가 아직 없다

allocate (daily_mean(ndays))
allocate (daily_mean(ndays))        ! 틀림: 이미 할당되어 있다

첫 대입은 정의되지 않은 동작이다. 오류 없이 지나가거나 이상한 값을 만들거나 멈출 수 있다. 두 번째 allocate 는 실행 중에 오류로 멈춘다. 쓰기 전에 allocate 하고, 반복 호출될 수 있는 곳에서는 상태를 먼저 확인한다.

if (allocated(daily_mean)) deallocate (daily_mean)
allocate (daily_mean(ndays))
daily_mean(1) = 14.25_dp

4. * 를 행렬 곱으로 착각한다

real(dp) :: a(3, 3), b(3, 3), c(3, 3)
c = a * b            ! 틀림: 같은 위치끼리의 곱이다

c = matmul(a, b)     ! 행렬 곱

* 는 두 배열 모양이 같으면 컴파일되므로 오류 메시지 없이 다른 값이 나온다. 행렬 곱이 필요한 자리인지 원소별 곱이 필요한 자리인지 식을 쓸 때마다 구분하자.

한눈에 보기

이 장에서 쓴 도구와 기억할 점
도구	쓰임	예	주의
2차원 선언	행·열 표	`real(dp) :: t(3, 4)`	첫 첨자가 행, 둘째가 열
열 우선 저장	메모리 배치	`t(1,1), t(2,1), t(3,1), t(1,2)…`	안쪽 반복문이 첫 첨자
`reshape`	모양 바꾸기	`reshape(v, [3, 4])`	열 우선으로 채움, 행 순서는 `order=[2,1]`
`transpose`	행·열 맞바꿈	`tt = transpose(t)`	3×4 가 4×3 이 됨
`matmul`	행렬 곱	`matmul(t, w)`	앞의 열 수 = 뒤의 행 수
`allocatable`	크기를 나중에 정함	`real(dp), allocatable :: a(:,:)`	차원은 콜론으로 비움
`allocate` / `deallocate`	메모리 잡기와 돌려주기	`allocate(a(n, m), stat=ierr)`	두 번 할당하면 오류
`size` / `shape`	크기 묻기	`size(t, 1)`, `shape(t)`	`shape` 는 배열을 돌려줌

첨자와 부분 배열이 가리키는 범위 (3행 4열 temp)
표기	가리키는 것	길이	메모리에서
`temp(2, 3)`	원소 하나	1	한 칸
`temp(:, 3)`	세 번째 열	3	붙어 있음
`temp(2, :)`	2일째 행	4	세 칸씩 떨어져 있음
`temp`	표 전체	12	열 순서로 이어짐

세부 규칙이 궁금하면 GNU Fortran 공식 문서에서 reshape, matmul, allocate 항목을 확인할 수 있다.

연습 문제

a = reshape([1, 2, 3, 4, 5, 6], [2, 3]) 으로 만든 정수 배열에서 a(1, 3) 과 a(2, 2) 의 값은 무엇인가?
위의 a 에 대해 b = transpose(a) 의 shape 와 b(3, 1) 의 값을 구하라.
transpose 를 쓰지 않고 이중 반복문으로 전치를 구현하라. 결과를 담을 배열은 할당 배열로 선언하고, size 로 크기를 정하라. 안쪽 반복문은 어느 배열의 메모리 순서에 맞추게 되는지도 설명하라.
a 가 2×3, x 가 길이 3 인 실수 배열일 때 matmul(a, x) 와 matmul(x, a) 중 오류가 나는 것은 어느 쪽이며 그 이유는 무엇인가?

정답과 해설

a(1, 3) 은 5, a(2, 2) 는 4 다. 값은 열 단위로 채워지므로 첫 열이 1, 2, 둘째 열이 3, 4, 셋째 열이 5, 6 이다. 셋째 열의 1행이 5, 둘째 열의 2행이 4 다.
shape(b) 는 [3, 2] 이고 b(3, 1) 은 a(1, 3) 과 같은 5 다. 전치된 배열의 b(j, i) 는 원본의 a(i, j) 이기 때문이다.
예시 답은 다음과 같다.
```
real, allocatable :: b(:, :)
integer :: i, j

allocate (b(size(a, 2), size(a, 1)))
do j = 1, size(a, 2)
  do i = 1, size(a, 1)
    b(j, i) = a(i, j)
  end do
end do
```
안쪽 i 는 a 의 첫 번째 첨자라 a 는 메모리 순서대로 읽는다. 반면 b 는 두 번째 첨자가 i 라서 건너뛰며 쓴다. 전치에서는 읽기와 쓰기 양쪽을 동시에 순서대로 만들 수 없으므로 어느 한쪽을 고르게 된다. 이런 이유로 내장 transpose 를 쓰는 것이 간단하다.
matmul(x, a) 가 오류다. 앞 인수가 길이 3 벡터이므로 뒤 인수의 행 수가 3 이어야 하는데 a 의 행 수는 2 다. 크기가 상수이면 컴파일 때, 아니면 실행 때 크기 불일치로 보고된다. matmul(a, x) 는 a 의 열 수 3 과 x 의 길이 3 이 맞아 길이 2 벡터를 돌려준다.