[アルゴリズム][実装]Fortran で マージソート

概要

• マージソートは，長さ$N$の配列を$O(N\log N)$時間でソートするアルゴリズムです．
• 本稿ではFortranでマージソートを実装します．

アルゴリズム

概要

• マージソートでは，配列を二つに分割していき，それを「マージ」していくことで配列をソートします．
• マージとは２つのソートされた配列を組み合わせて１つのソートされた配列を作る操作のことを言います．
• マージの実装は簡単で，２つの配列の先頭の値を比較し，小さいほうの配列から$\textbf{popleft}$して結果の配列に$\textbf{append}$することを繰り返すだけです．（$\textrm{popleft}$しているので２つの配列の先頭は動的に変わっていきます）計算時間は結果の配列の長さに比例します．
• 下図のように木構造で考えると，分割の回数は$N-1$回で，マージの計算時間はそれぞれの深さでの合計が結局$N$に比例するので，全体の深さがおよそ$\log N$となることを考えると，計算時間は$O(N\log N)$となります．

疑似コード

• マージソートの疑似コードを以下に示します．
• 分割の際に木構造を記憶しておく必要があるので，普通は再起を用いて実装します．
• 標準ライブラリに$\textrm{popleft}$をサポートしているコンテナが存在していない場合を考慮し，ここでは先頭にあたるインデックスを管理しています．

実装

• fortranでの実装は，疑似コードとほぼ同じです．

! ======================================================================================
recursive function sorted(xin) result(xout)
implicit none
    !-------------------------------------------------------------------------
    double precision,intent(in)::xin(:)
    double precision::xout(size(xin))
    double precision,allocatable::x1(:)
    double precision,allocatable::x2(:)
    integer::n
    integer::n1
    integer::n2
    integer::i
    integer::i1
    integer::i2
    !-------------------------------------------------------------------------
    n=size(xin)
    if (n==1) then
        xout(1)=xin(1)
        return
    endif
    n1=n/2
    n2=n-n1
    allocate(x1(n1))
    allocate(x2(n2))
    x1=sorted(xin(1:n1))
    x2=sorted(xin(n1+1:n))
    i1=1
    i2=1
    do i=1,n
        if (x1(i1)<x2(i2)) then
            xout(i)=x1(i1)
            i1=i1+1
            if (i1>n1) then
                xout(i+1:)=x2(i2:)
                exit
            endif
        else
            xout(i)=x2(i2)
            i2=i2+1
            if (i2>n2) then
                xout(i+1:)=x1(i1:)
                exit
            endif
        endif
    enddo
    return
    !-------------------------------------------------------------------------
end function
! ======================================================================================

テスト

• 上記のコードを，gfortranでテストします．
• ここではランダムに生成した配列を並び替えてもらいます．

N=10

• 並び替える配列の長さが10のとき

コード

! ======================================================================================
program test
use sort
implicit none
    !-------------------------------------------------------------------------
    real::y
    integer::i
    double precision::x(10)
    !-------------------------------------------------------------------------
    do i=1,10
        call random_number(y)
        x(i)=y
    enddo
    write(*,"('x=(',10(f5.3,',',1x),')')")x
    write(*,"('sorted(x)=(',10(f5.3,',',1x),')')")sorted(x)
    !-------------------------------------------------------------------------
end program
! ======================================================================================

結果

x=(0.998, 0.567, 0.966, 0.748, 0.367, 0.481, 0.074, 0.005, 0.347, 0.342, 0.218, 0.133, 0.901, 0.387, 0.445, 0.662, 0.016, 0.651, 0.646, 0.323, )
sorted(x)=(0.005, 0.016, 0.074, 0.133, 0.218, 0.323, 0.342, 0.347, 0.367, 0.387, 0.445, 0.481, 0.567, 0.646, 0.651, 0.662, 0.748, 0.901, 0.966, 0.998, )

N=100

• 以下のPythonコードで，直下のフォルダに100個のテストケースを用意します．

import numpy as np
N=100
T=100
for i in range(T):
    with open(f".\\cases\\{i:02d}","w") as f:
        arr=np.random.uniform(0,1,(N,))
        f.write(" ".join([f"{float(j):.4f}" for j in arr])+"\n")
        arr.sort()
        f.write(" ".join([f"{float(j):.4f}" for j in arr])+"\n")

コード

! ======================================================================================
program test
use sort
implicit none
    !-------------------------------------------------------------------------
    integer::i
    double precision::x(100)
    double precision::y(100)
    character::fn*20
    !-------------------------------------------------------------------------
    do i=0,99
        write(fn,"('.\\cases\\',i2.2)")i
        open(1,file=fn,action="read",status="old")
        read(1,"(100(f6.4,1x))")x
        read(1,"(100(f6.4,1x))")y
        write(*,"('case ',i0,':: ')",advance="no")i+1
        if (all(y==sorted(x))) then
            print*,"Correct Answer"
        else
            print*,"Wrong Answer"
        endif
    enddo
    !-------------------------------------------------------------------------
end program
! ======================================================================================

結果

結果を見る

case 1::  Correct Answer
case 2::  Correct Answer
case 3::  Correct Answer
case 4::  Correct Answer
case 5::  Correct Answer
case 6::  Correct Answer
case 7::  Correct Answer
case 8::  Correct Answer
case 9::  Correct Answer
case 10::  Correct Answer
case 11::  Correct Answer
case 12::  Correct Answer
case 13::  Correct Answer
case 14::  Correct Answer
case 15::  Correct Answer
case 16::  Correct Answer
case 17::  Correct Answer
case 18::  Correct Answer
case 19::  Correct Answer
case 20::  Correct Answer
case 21::  Correct Answer
case 22::  Correct Answer
case 23::  Correct Answer
case 24::  Correct Answer
case 25::  Correct Answer
case 26::  Correct Answer
case 27::  Correct Answer
case 28::  Correct Answer
case 29::  Correct Answer
case 30::  Correct Answer
case 31::  Correct Answer
case 32::  Correct Answer
case 33::  Correct Answer
case 34::  Correct Answer
case 35::  Correct Answer
case 36::  Correct Answer
case 37::  Correct Answer
case 38::  Correct Answer
case 39::  Correct Answer
case 40::  Correct Answer
case 41::  Correct Answer
case 42::  Correct Answer
case 43::  Correct Answer
case 44::  Correct Answer
case 45::  Correct Answer
case 46::  Correct Answer
case 47::  Correct Answer
case 48::  Correct Answer
case 49::  Correct Answer
case 50::  Correct Answer
case 51::  Correct Answer
case 52::  Correct Answer
case 53::  Correct Answer
case 54::  Correct Answer
case 55::  Correct Answer
case 56::  Correct Answer
case 57::  Correct Answer
case 58::  Correct Answer
case 59::  Correct Answer
case 60::  Correct Answer
case 61::  Correct Answer
case 62::  Correct Answer
case 63::  Correct Answer
case 64::  Correct Answer
case 65::  Correct Answer
case 66::  Correct Answer
case 67::  Correct Answer
case 68::  Correct Answer
case 69::  Correct Answer
case 70::  Correct Answer
case 71::  Correct Answer
case 72::  Correct Answer
case 73::  Correct Answer
case 74::  Correct Answer
case 75::  Correct Answer
case 76::  Correct Answer
case 77::  Correct Answer
case 78::  Correct Answer
case 79::  Correct Answer
case 80::  Correct Answer
case 81::  Correct Answer
case 82::  Correct Answer
case 83::  Correct Answer
case 84::  Correct Answer
case 85::  Correct Answer
case 86::  Correct Answer
case 87::  Correct Answer
case 88::  Correct Answer
case 89::  Correct Answer
case 90::  Correct Answer
case 91::  Correct Answer
case 92::  Correct Answer
case 93::  Correct Answer
case 94::  Correct Answer
case 95::  Correct Answer
case 96::  Correct Answer
case 97::  Correct Answer
case 98::  Correct Answer
case 99::  Correct Answer
case 100::  Correct Answer

• 100個すべてのテストケースで正答しています．