LoginSignup
5
2

More than 3 years have passed since last update.

JuliaからFortranの構造体を引数にもつsubroutineを呼び出す

Last updated at Posted at 2019-07-09

「JuliaからFortranのsubroutineを呼び出す」
https://qiita.com/cometscome_phys/items/d098b00ab97f44916783
では、FortranのsubroutineをJuliaから呼び出していました。
しかし、Fortranで構造体(type)を使っている場合にどのようにすればいいのかについては述べていませんでした。
この記事では、JuliaのStructをFortranのsubroutineの引数として突っ込む方法について述べます。

最後に、自分がわからなかった点について書いています。structの中の配列サイズを任意にする方法がわかりませんでした。

バージョン

Julia 1.1.0

Fortranコード

まず、呼び出したいFortranコードとして、

testfunc.f90
module test

  contains
    subroutine testpoint(p) 
        implicit none
        type point         
            real(8)::x,y,z
            real(8)::d(10)
            real(8)::d2(3,3)
        end type point

        type(point)::p 
        real(8)::t

        t = p%x
        p%x = p%y
        p%y = t
        p%z = -2*p%z
        p%d(1) = t 

        p%d2(1,3) = 100

        write(*,*) "p = ",p
        write(*,*) "d = ",p%d
        write(*,*) "d2 = ",p%d2
        return
    end subroutine testpoint


end module

というものを考えましょう。このコードを、

gfortran testfunc.f90 -o test.so -shared -fPIC 

でコンパイルします。
そして、

nm test.so 

をして、

0000000000000af0 T ___test_MOD_testpoint
                 U __gfortran_st_write
                 U __gfortran_st_write_done
                 U __gfortran_transfer_array_write
                 U __gfortran_transfer_character_write
                 U __gfortran_transfer_real_write
                 U dyld_stub_binder

というものを見て、___test_MOD_testpointという関数があることを確認します。
なお、以前の記事にあるようなbind(c, name = 'mat')のようなものをつけようとすると、

testfunc.f90:4:26:

    4 |     subroutine testpoint(p) bind(c, name = 'mat')
      |                          1
Error: Variable 'p' at (1) is a dummy argument to the BIND(C) procedure 'testpoint' but is not C interoperable because derived type 'point' is not C interoperable

というエラーでできませんでしたので、今回はnmを使って関数を探しておくことにします。

[19/07/10 追記]
@cure_honey さんにやり方を教えてもらいました。記事の最後に追記します。

Fortranを呼び出すJuliaコード

次に、上のコードを呼び出すJuliaコードを

mutable struct Point
    x::Float64
    y::Float64
    z::Float64
    d::NTuple{10,Float64}
    d2::NTuple{9,Float64}
end

A = rand(3,3)
println("A = ", A)
foo = Point(1,2,3,ntuple(i -> 2*i, 10),tuple(A...))

ccall((:__test_MOD_testpoint,"./test.so"),Nothing, 
    (Ref{Point},), #b
    foo)
println(foo)

としました。
Juliaでmutable structで作っておくと、Fortran側の構造体として引き渡された後に値を変更することができます。

なお、出力は、

A = [0.964356 0.451978 0.0943697; 0.304721 0.417815 0.110606; 0.235112 0.916069 0.454091]
 p =    2.0000000000000000        1.0000000000000000       -6.0000000000000000        1.0000000000000000        4.0000000000000000        6.0000000000000000        8.0000000000000000        10.000000000000000        12.000000000000000        14.000000000000000        16.000000000000000        18.000000000000000        20.000000000000000       0.96435610564740015       0.30472062004468348       0.23511159444040275       0.45197848999368229       0.41781472623058180       0.91606942024800531        100.00000000000000       0.11060585340863494       0.45409130957036381     
 d =    1.0000000000000000        4.0000000000000000        6.0000000000000000        8.0000000000000000        10.000000000000000        12.000000000000000        14.000000000000000        16.000000000000000        18.000000000000000        20.000000000000000     
 d2 =   0.96435610564740015       0.30472062004468348       0.23511159444040275       0.45197848999368229       0.41781472623058180       0.91606942024800531        100.00000000000000       0.11060585340863494       0.45409130957036381     
Point(2.0, 1.0, -6.0, (1.0, 4.0, 6.0, 8.0, 10.0, 12.0, 14.0, 16.0, 18.0, 20.0), (0.9643561056474002, 0.3047206200446835, 0.23511159444040275, 0.4519784899936823, 0.4178147262305818, 0.9160694202480053, 100.0, 0.11060585340863494, 0.4540913095703638))

となります。ちゃんとFortran側で書き換えていることがわかります。

配列の扱い(未解決点)

配列に関しては、mutable structで作った時にサイズを指定しなければダメでした。そのためには、Ntupleでなければならないようです。
(参考:https://docs.julialang.org/en/v1/manual/calling-c-and-fortran-code/index.html)
Array{Float64,1}ではできませんでした。
structの中の配列サイズを任意にする方法がわかりませんでした。
マクロを使えば、できるのでしょうか。

追記

Fortranでbindでエラーが出た件ですが、typeの部分にもつけることで解決できるそうです (@cure_honey さんありがとうございます)。

以下のコードがFortranとJuliaのコードです。

module test

  contains
    subroutine testpoint(p) bind(c, name = 'testpoint')
        implicit none
        type, bind(c) :: point      
            real(8)::x,y,z
            real(8)::d(10)
            real(8)::d2(3,3)
        end type point

        type(point)::p 
        real(8)::t

        t = p%x
        p%x = p%y
        p%y = t
        p%z = -2*p%z
        p%d(1) = t 

        p%d2(1,3) = 100

        write(*,*) "p = ",p
        write(*,*) "d = ",p%d
        write(*,*) "d2 = ",p%d2
        return
    end subroutine testpoint

end module
mutable struct Point
    x::Float64
    y::Float64
    z::Float64
    d::NTuple{10,Float64}
    d2::NTuple{9,Float64}
end

A = rand(3,3)
println("A = ", A)
foo = Point(1,2,3,ntuple(i -> 2*i, 10),tuple(A...))

ccall((:testpoint,"./test.so"),Nothing, 
    (Ref{Point},), #b
    foo)
println(foo)


5
2
2

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
5
2