Help us understand the problem. What is going on with this article?

Fortranでinsert/eraseが早いvector作成記

 はじめに

 arrayListは配列より柔軟で実装も軽い便利なデータ構造だが、「末尾以外の要素に対する挿入/削除」の計算量が$O(N)$であり、ボトルネックである。
 それに対してlinkedListでは、このボトルネックは解消されており、末尾以外の要素に対する挿入/削除を$O(1)$で行うことができるが、arrayListで$O(1)$だった「i番目の値の取得」が$O(N)$と悪化している。

 主要なデータ構造に対して、基本的な操作を行った際の計算量を以下に示す。

arrayList linkedList deque(arrayList)
i番目の値を検索する $O(1)$ $O(N)$ $O(1)$
末尾の挿入/削除 $O(1)$ $O(1)$ $O(1)$
先頭への挿入/削除 $O(N)$ $O(1)$ $O(1)$
i番目の値の挿入/削除 $O(N)$ $O(1)$ $O(N)$

 この表から、どのデータ構造も、表に示した4つの操作の少なくとも1つの計算量が$O(N)$だとわかる。
 一方、今回紹介するデータ構造(vectree)に対する同様な表を以下に示す。

vectree
i番目の値を検索する $O(\log(N))$
末尾の挿入/削除 $O(\log(N))$
先頭への挿入/削除 $O(\log(N))$
i番目の値の挿入/削除 $O(\log(N))$

 この表から分かるように、vectreeは表に示した4つの操作に対して$O(\log(N))$で動作する。vectreeは「arrayList,linkedList,dequeのいいとこ取りをして、ちょっとデバフをかけた」ようなデータ構造だといえる。

 競技プログラミング的事情から見ても、$O(N\log N)$の処理は、おそらく$N = 10^5$~$10^6$くらいまで耐えられるため、実用的であろうと考えられる。

vectreeとは

 まず、vectreeという名前だが、これは正式なものではなく、この記事でつけられた名称であり、今回はこの名称を継承している。ここを参考にしたらしい。

 上の記事に図などを用いた丁寧な解説があるため、そちらを参考にしたほうが良いと思うが、vectreeとは、
* 平衡二分木によってデータを保持している。
* i番目の要素は、木を通りがけ順に見ていったときのi番目のノードに保存
* i番目の要素から見たi-1番目の要素のことをpredecessor(前任者)と呼ぶ
* i番目の要素から見たi+1番目の要素のことをsuccesor(後継者)と呼ぶ
* i番目への挿入では、predecessorの右の子ノードに挿入する
* i番目の要素の削除は、predecessorまたはsuccesorをi番目の位置へ移動させる
* i番目の要素を検索するために、各ノードが「自身を根とする部分木の要素数」を保持する
* 「自身を根とする部分木の要素数」は、挿入/削除の度に更新する。
 というものらしい。今回はこのvectreeをFortranで実装することにした。

実装

 平衡二分木にはAA木を用いた。リンク先のコードがFortranに似ていたため、簡易に実装することができた。
 実装は、
* 平衡二分木の構造をvectree_node_mod
* vectree_node_modのラッパーをvectree_mod
 というようにした。平衡二分木の実装は、vectreeの挿入の操作が普通の平衡二分木と異なるため、AA木の実装単体を分離することはできなかった。
 また、procedure名はc++に倣った。

環境

GNU Fortran (Ubuntu 7.5.0-3ubuntu1~18.04) 7.5.0
AtCoderのコードテスト

vectree_node_mod

 各procedureにはvtn_という接頭語をつけた。vectreeを使う時に、このmoduleの内部を気にする必要はない。

module vec_tree_node_mod
    use,intrinsic :: iso_fortran_env
    implicit none
    private
    public:: vtn_insert, vtn_push_back
    public:: vtn_erase
    public:: vtn_to_array
    public:: vtn_print, vtn_print_debug
    public:: vtn_at
    type,public:: vec_tree_node
        integer(int32):: level=1,tree_size=1
        integer(int32):: value
        type(vec_tree_node),pointer:: left_child => null()
        type(vec_tree_node),pointer:: right_child => null()
    end type
contains
    function vtn_has_left_child(n) result(ret)
        type(vec_tree_node),pointer:: n
        logical:: ret

        ret = associated(n%left_child)
    end function


    function vtn_has_right_child(n) result(ret)
        type(vec_tree_node),pointer:: n
        logical:: ret

        ret = associated(n%right_child)
    end function


    function vtn_is_leaf(n) result(ret)
        type(vec_tree_node),pointer:: n
        logical:: ret

        ret = .not. (vtn_has_left_child(n) .or. vtn_has_right_child(n))
    end function


    function vtn_predecessor(n) result(ret)
        !left -> rightend
        type(vec_tree_node),pointer:: n,ret

        ret => n%left_child
        do while(vtn_has_right_child(ret))
            ret => ret%right_child
        end do
    end function


    function vtn_succesor(n) result(ret)
        !right -> leftend
        type(vec_tree_node),pointer:: n,ret

        ret => n%right_child
        do while(vtn_has_left_child(ret))
            ret => ret%left_child
        end do
    end function


    subroutine vtn_size_update(n)
        type(vec_tree_node),pointer:: n
        integer(int32):: old_size, new_size

        old_size = n%tree_size
        new_size = vtn_left_size_of(n)+vtn_right_size_of(n) + 1
        n%tree_size = new_size
    end subroutine


    function vtn_left_size_of(n) result(ret)
        type(vec_tree_node),pointer:: n
        integer(int32):: ret

        ret = 0
        if (vtn_has_left_child(n)) ret = n%left_child%tree_size
    end function


    function vtn_right_size_of(n) result(ret)
        type(vec_tree_node),pointer:: n
        integer(int32):: ret

        ret = 0
        if (vtn_has_right_child(n)) ret = n%right_child%tree_size
    end function


    function vtn_ind_at(n) result(ret)
        type(vec_tree_node),pointer:: n
        integer(int32):: ret

        ret = vtn_left_size_of(n)+1
    end function


    function vtn_skew(n) result(ret)
        type(vec_tree_node),pointer:: n, ret, sub_n

        if (.not. associated(n)) then
            allocate(ret)
            ret => null()
            return
        else if (.not. vtn_has_left_child(n)) then
            ret => n
        else if (n%level == n%left_child%level) then
            sub_n => n%left_child
            n%left_child => sub_n%right_child
            sub_n%right_child => n
            ret => sub_n
        else
            ret => n
        end if
        call vtn_size_update(n)
        call vtn_size_update(ret)
    end function


    function vtn_split(n) result(ret)
        type(vec_tree_node),pointer:: n,ret,sub_n

        if (.not. associated(n)) then
            allocate(ret)
            ret => null()
            return
        else if (.not. vtn_has_right_child(n) .or. .not. vtn_has_right_child(n%right_child)) then
            ret => n
        else if (n%level == n%right_child%right_child%level) then
            sub_n => n%right_child
            n%right_child => sub_n%left_child
            sub_n%left_child => n
            sub_n%level = sub_n%level+1
            ret => sub_n
        else
            ret => n
        end if
        call vtn_size_update(n)
        call vtn_size_update(ret)
    end function


    recursive function vtn_push_back(n,i,val) result(ret)
        type(vec_tree_node),pointer:: n,ret
        integer(int32):: i
        integer(int32):: val

        if (.not. associated(n)) then
            allocate(ret)
            ret = vec_tree_node(value=val)
            return
        else
            n%right_child => vtn_push_back(n%right_child,i,val)
        end if
        n => vtn_skew(n)
        n => vtn_split(n)
        ret => n
    end function


    function vtn_level_of(n) result(ret)
        type(vec_tree_node),pointer:: n
        integer(int32):: ret

        ret = 0
        if (associated(n)) ret = n%level
    end function


    function vtn_decrease_level(n) result(ret)
        type(vec_tree_node),pointer:: n, ret
        integer(int32):: should_be

        should_be = min(vtn_level_of(n%left_child), vtn_level_of(n%right_child)) + 1
        if (should_be < n%level) then
            n%level = should_be
            if (should_be < vtn_level_of(n%right_child)) then
                n%right_child%level = should_be
            end if
        end if
        ret => n
    end function


    recursive function vtn_insert(n,i,val) result(ret)
        type(vec_tree_node),pointer:: n,ret
        integer(int32):: i
        integer(int32):: val

        if (.not. associated(n)) then
            allocate(ret)
            ret = vec_tree_node(value=val,level=1)
            return
        end if
        if (i == vtn_ind_at(n)) then
            n%left_child => vtn_push_back(n%left_child,i,val)
        else if (i > vtn_ind_at(n)) then
            n%right_child => vtn_insert(n%right_child, i-(vtn_ind_at(n)),val)
        else if (i < vtn_ind_at(n)) then
            n%left_child => vtn_insert(n%left_child,i,val)
        end if
        n => vtn_skew(n)
        n => vtn_split(n)
        ret => n
    end function


    recursive function vtn_erase(n,i,ev) result(ret)
        type(vec_tree_node),pointer:: n,ret,nn
        integer(int32):: i
        integer(int32),optional:: ev

        if (.not. associated(n)) then
            allocate(ret)
            ret => null()
            return
        else if (vtn_ind_at(n) > i) then
            n%left_child => vtn_erase(n%left_child, i,ev)
        else if (vtn_ind_at(n) < i) then
            n%right_child => vtn_erase(n%right_child, i-vtn_ind_at(n),ev)
        else
            if (present(ev)) ev = n%value
            if (vtn_is_leaf(n)) then
                allocate(ret)
                ret => null()
                return
            else if (vtn_has_left_child(n)) then
                nn => vtn_predecessor(n)
                n%left_child => vtn_erase(n%left_child,nn%tree_size)
                n%value = nn%value
            else
                nn => vtn_succesor(n)
                n%right_child => vtn_erase(n%right_child,nn%tree_size)
                n%value = nn%value
            end if
        end if
        n => vtn_decrease_level(n)
        n => vtn_skew(n)
        n%right_child => vtn_skew(n%right_child)
        if (vtn_has_left_child(n)) n%right_child%right_child => vtn_skew(n%right_child%right_child)
        n => vtn_split(n)
        n%right_child => vtn_split(n%right_child)
        ret => n
    end function


    recursive subroutine vtn_print_debug(n,i,d)
        type(vec_tree_node),pointer:: n
        integer(int32):: i,j
        character(1):: d

        if (.not. associated(n)) return
        call vtn_print(n%right_child,i+1,"r")
        do j=1,i; write(*,'(a)',advance="no") '  '; end do
        if (d == "r") then
            write(*,'(a)',advance="no") "┌ "
        else if (d == "l") then
            write(*,'(a)',advance="no") "└ "
        else
            write(*,'(a)',advance="no") "  "
        end if
        print'("[ ",3(i0,1x),"]")', n%value, n%level, n%tree_size
        call vtn_print(n%left_child,i+1,"l")
    end subroutine


    recursive subroutine vtn_print(n,i,d)
        type(vec_tree_node),pointer:: n
        integer(int32):: i,j
        character(1):: d

        if (.not. associated(n)) return
        call vtn_print(n%right_child,i+1,"r")
        do j=1,i; write(*,'(a)',advance="no") '  '; end do
        if (d == "r") then
            write(*,'(a)',advance="no") "┌ "
        else if (d == "l") then
            write(*,'(a)',advance="no") "└ "
        else
            write(*,'(a)',advance="no") "  "
        end if
        print'(i0)', n%value
        call vtn_print(n%left_child,i+1,"l")
    end subroutine


    recursive subroutine vtn_to_array(n,ind,ar)
        type(vec_tree_node),pointer:: n
        integer(int32),intent(inout):: ind
        integer(int32),intent(inout):: ar(:)
        integer(int32):: i

        if (.not. associated(n)) return
        i=ind
        call vtn_to_array(n%left_child,i,ar)
        ar(i) = n%value
        i=i+1
        call vtn_to_array(n%right_child,i,ar)
        ind=i
    end subroutine


    recursive function vtn_at(n,ind) result(ret)
        type(vec_tree_node),pointer:: n
        integer(int32):: ind
        integer(int32):: ret

        ret = 0
        if (.not. associated(n)) return
        if (vtn_ind_at(n) > ind) then
            ret = vtn_at(n%left_child, ind)
        else if (vtn_ind_at(n) < ind) then
            ret = vtn_at(n%right_child, ind-vtn_ind_at(n))
        else ! vtn_ind_at(n) == ind
            ret = n%value
        end if
    end function
end module

vectree_mod

各procedureにはvt_という接頭語をつけた。vectreeを用いるときに使用する関数とその使用方法についてはコードの後に示す。

module vec_tree_mod
    use,intrinsic :: iso_fortran_env
    use vec_tree_node_mod
    implicit none
    private
    public:: vt_print, vt_print_debug
    public:: vt_size
    public:: vt_to_array
    type,public:: vectree
        type(vec_tree_node),pointer,private:: head => null()
        integer(int32),private:: len = 0
    contains
        procedure,public:: insert => vt_insert, erase => vt_erase
        procedure,public:: push_back => vt_push_back, push_front => vt_push_front
        procedure,public:: pop => vt_pop, pop_back => vt_pop_back, pop_front => vt_pop_front
        procedure,public:: at => vt_value_at
    end type
contains
    subroutine vt_insert(vt,i,val)
        class(vectree):: vt
        type(vec_tree_node),pointer:: new
        integer(int32):: i
        integer(int32):: val

        if (i < 1 .or. vt%len+1 < i) then
            print'(a)', "inserting index is too big"
            print'(a)', "index, vectree length"
            print'(*(i0,1x))', i, vt%len
            stop
        else if (.not. associated(vt%head)) then
            allocate(new)
            new = vec_tree_node(value=val)
            vt%head => new
        else if (i == vt%len+1) then
            vt%head => vtn_push_back(vt%head,i,val)
        else
            vt%head => vtn_insert(vt%head,i,val)
        end if
        vt%len=vt%len+1
    end subroutine


    subroutine vt_push_back(vt,val)
        class(vectree):: vt
        integer(int32):: val

        call vt_insert(vt,vt%len+1,val)
    end subroutine


    subroutine vt_push_front(vt,val)
        class(vectree):: vt
        integer(int32):: val

        call vt_insert(vt,1,val)
    end subroutine


    function vt_pop(vt,i) result(erased_value)
        class(vectree):: vt
        integer(int32):: i,erased_value

        if (i < 1 .or. vt%len < i) then
            print'(a)', "remove index is too big"
            print'(a)', "index, vectree length"
            print'(*(i0,1x))', i, vt%len
            stop
        end if
        vt%head => vtn_erase(vt%head,i,erased_value)
        vt%len=vt%len-1
    end function


    subroutine vt_erase(vt,i)
        class(vectree):: vt
        integer(int32):: i,tmp

        tmp = vt%pop(i)
    end subroutine


    function vt_pop_back(vt) result(ret)
        class(vectree):: vt
        integer(int32):: ret

        ret = vt%pop(vt%len)
    end function


    function vt_pop_front(vt) result(ret)
        class(vectree):: vt
        integer(int32):: ret

        ret = vt%pop(1)
    end function


    function vt_value_at(vt,i) result(ret)
        class(vectree):: vt
        integer(int32):: i
        integer(int32):: ret

        ret = vtn_at(vt%head,i)
    end function


    subroutine vt_print(vt)
        type(vectree):: vt

        call vtn_print(vt%head,0,' ')
    end subroutine


    subroutine vt_print_debug(vt)
        type(vectree):: vt

        print*, "==[value, vtn_level_of, tree_size]=="
        call vtn_print_debug(vt%head,0,' ')
    end subroutine


    function vt_to_array(vt) result(ret)
        type(vectree):: vt
        integer(int32):: ret(vt%len),i
        i=1
        call vtn_to_array(vt%head,i,ret)
    end function


    function vt_size(vt) result(ret)
        type(vectree):: vt
        integer(int32):: ret

        ret = vt%len
    end function
end module

実装したprocedureとその使用方法

 以下のコードで、vectreeを操作するためのすべてのprocedureについて紹介する。10個作成した。作りすぎた気がする。

program main
    use,intrinsic :: iso_fortran_env
    use vec_tree_mod
    implicit none
    type(vectree):: vt
    integer(int32):: i,n
    integer(int32),allocatable:: ar(:)

    n=5
    print"(a)", "===insert, push==="
    do i=1,2*n
        v = i*2
        call vt%insert(i,v) ! i番目にvを挿入 O(log N)
        v = i*3
        call vt%push_front(v) ! 先頭にvを挿入 O(log N)
        v = i*4
        call vt%push_back(v) ! 末尾にvを挿入 O(log N)
    end do
    print"(a)", "===vt_print==="
    call vt_print(vt) ! 平衡二分木を可愛らしく表示する O(N)
    print"(a)", "===erase, pop==="
    do i=1,n/2
        call vt%erase(i) ! i番目の要素を削除 O(log N)
        print*, vt%pop(i) ! i番目の要素を削除し、その値が返る O(log N)
        print*, vt%pop_front() ! 先頭の要素を削除し、その値が返る O(log N)
        print*, vt%pop_back() ! 末尾の要素を削除し、その値が返る O(log N)
    end do

    do i=1,n/2
        print*, vt%at(i) ! i番目の要素を参照 O(log N)
    end do

    n = vt_size(vt) ! vtの要素数を返す O(1)
    allocate(ar(n))
    ar(1:n) = vt_to_array(vt) ! vtの要素全体を配列で返す O(N)
    print"(a)", "===to_array==="
    print'(*(i0,1x))', ar(:)
end program main

実行結果

===insert, push===
===vt_print===
            ┌ 40
          ┌ 36
        ┌ 32
          └ 28
      ┌ 24
        └ 20
    ┌ 16
        ┌ 12
      └ 8
        └ 4
  ┌ 2
        ┌ 4
      ┌ 6
        └ 8
    └ 10
          ┌ 12
        ┌ 14
          └ 16
      └ 18
        └ 20
  3
      ┌ 6
    ┌ 9
      └ 12
  └ 15
        ┌ 18
      ┌ 21
        └ 24
    └ 27
      └ 30
===erase, pop===
          27
          24
          40
          15
          21
          36
          12
           9
===to_array===
12 9 6 3 20 18 16 14 12 10 8 6 4 2 4 8 12 16 20 24 28 32

バグなど

ランダムなinsert/eraseのクエリ($10^5$程度)をpythonで作成し、その出力をpythonのlist,fortranの自作vector、vectreeの3つでdiffを取って問題ない事を確認済み。また、AtCoderの「コードテスト」上でも正常に動くことを確認。(自作vectorのinsert/eraseの速度がpythonのリストに負けていて泣きそうになった。)

もしバグをみつけたら教えてください。

まとめ

Fortranでvectreeを作成したよ。

Authns
python fortran atcoder緑 ゲーム作成 分子動力学法
Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away