LoginSignup
2
0
@PannacottaFGK(パンナコッタ・F・G・カズンド)

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

Last updated at Posted at 2024-06-22

概要

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

アルゴリズム

概要

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

Whiteboard 3.png

疑似コード

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

sodapdf-converted (1).jpg

実装

  • 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個すべてのテストケースで正答しています．
2
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
2
0