# Correct way to ccall a Fortran subroutine requiring pointers

**URL:** <https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191>\
**Category:** General Usage\
**Tags:** fortran\
**Created:** [May 2, 2023, 10:57am UTC](https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191 "2023-05-02T10:57:36Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![lucifer1004](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lucifer1004/32/22311_2.png) [@lucifer1004](https://discourse.julialang.org/u/lucifer1004)\
**Post date:** [May 2, 2023, 10:57am UTC](https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191/1 "2023-05-02T10:57:36Z")

</div>

Given the following Fortran code:

```fortran
MODULE PTR
   IMPLICIT NONE
   REAL, SAVE, TARGET, ALLOCATABLE :: A(:, :), B(:, :)

CONTAINS

   SUBROUTINE INIT(AP, BP)
      REAL, POINTER :: AP(:, :), BP(:, :)

      ALLOCATE(A(6, 6), B(6, 6))
      AP => A(1:4, 2:5)
      BP => B(2:3, 3:6)

      AP = 3.2
      PRINT *, AP

      BP = 4.5
      PRINT *, BP
   END SUBROUTINE
END MODULE

```

Currently, I can only ccall the subroutine with the following hacky method:

```julia
using Libdl

lib = dlopen("minimal.so")

none = () -> nothing
ap = Ref(@cfunction(none, Cvoid, ()))
bp = Ref(@cfunction(none, Cvoid, ()))

ccall(dlsym(lib, :__ptr_MOD_init), Cvoid, (Ref{Ptr{Cvoid}}, Ref{Ptr{Cvoid}}), ap, bp)

a = unsafe_wrap(Array{Float32,2}, convert(Ptr{Float32}, ap[]), (6, 6))
@show a[1:4, 1:4]

b = unsafe_wrap(Array{Float32,2}, convert(Ptr{Float32}, bp[]), (6, 6))
@show b[1:2, 1:4]

```

Output is:

```julia-auto
   3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005    
   4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000    
a[1:4, 1:4] = Float32[3.2 3.2 3.2 3.2; 3.2 3.2 3.2 3.2; 3.2 3.2 3.2 3.2; 3.2 3.2 3.2 3.2]
b[1:2, 1:4] = Float32[4.5 4.5 4.5 4.5; 4.5 4.5 4.5 4.5]

```

I used `@cfunction` to get a function pointer here and pass them to the Fortran subroutine.

All other attempts failed. For example, if I use `ap = Ref{Ptr{Float32}}(C_NULL)` and change the ccall arguments accordingly,

```julia
ap = Ref{Ptr{Float32}}(C_NULL)
bp = Ref{Ptr{Float32}}(C_NULL)

ccall(dlsym(lib, :__ptr_MOD_init), Cvoid, (Ref{Ptr{Float32}}, Ref{Ptr{Float32}}), ap, bp)

```

The output will be:

```julia-auto
   4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000    
a[1:4, 1:4] = Float32[0.0 1.0f-44 7.0f-45 4.3228914f-38; 0.0 0.0 0.0 0.0; 9.0f-44 0.0 1.1f-44 0.0; 0.0 0.0 0.0 1.0f-45]

signal (11): Segmentation fault
...error messages omitted...

```

Using `ap = zeros(Float32, 6, 6)`, i.e.

```julia
ap = zeros(Float32, 6, 6)
bp = zeros(Float32, 6, 6)

ccall(dlsym(lib, :__ptr_MOD_init), Cvoid, (Ptr{Float32}, Ptr{Float32}), ap, bp)

@show ap
@show bp

```

will not cause a segmentation fault, however, it is obvious that we cannot get the correct data:

```julia-auto
   3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005 3.20000005    
   4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000 4.50000000    
ap = Float32[4.7818537f-38 0.0 1.0f-45 1.0f-45 0.0 0.0; 0.0 1.079f-42 0.0 0.0 0.0 0.0; NaN 6.0f-45 6.0f-45 6.0f-45 0.0 0.0; NaN 0.0 0.0 0.0 0.0 0.0; 6.0f-45 1.0f-45 8.0f-45 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0]
bp = Float32[1.7421598f-38 0.0 1.0f-45 1.0f-45 0.0 0.0; 0.0 1.079f-42 0.0 0.0 0.0 0.0; NaN 6.0f-45 3.0f-45 6.0f-45 0.0 0.0; NaN 0.0 0.0 0.0 0.0 0.0; 6.0f-45 1.0f-45 8.0f-45 0.0 0.0 0.0; 0.0 0.0 0.0 0.0 0.0 0.0]

```

Can anyone explain the internal mechanism? Is there a better way to handle Fortran pointers other than using function pointers?

## Another question

Another question is about the indices. Note that in the first example, `b[1:2, 1:4]` in Julia corresponds to `B(2:3, 3:6)` in Fortran, is it always true that indices will start from 1 at the Julia side?

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [May 2, 2023, 11:16am UTC](https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191/2 "2023-05-02T11:16:42Z")

</div>

I think the answer to your general question is this one: [Calling Fortran variables from Julia - #3 by anowacki](https://discourse.julialang.org/t/calling-fortran-variables-from-julia/27293/3)

> [@lucifer1004](#):
>
> Julia corresponds to `B(2:3, 3:6)` in Fortran, is it always true that indices will start from 1 at the Julia side?

Note that this is not specific of julia, but of the interface. Since you pass a reference, one side does not know about how the indexing was at the other side:

```fortran
program main
    integer :: i, a(3), b
    do i = 1,3
        a(i) = i
    end do
    print *, a
    call test(a(2:3),b)
    print *, b
end

subroutine test(a,b)
    integer :: a(2), b
    b = a(1)
end

```

Will print:

```julia
           1 2 3
           2

```

which means that `b` assumed the value of what was `a(2)` in the main program.

You can, in Julia, use `OffsetArrays` to use arbitrary indexing (but then you need one additional layer of care in passing the array to fortran, since an `OffsetArray` is an array contained in a struct that carries the information about the offsets.

---

<div class="post-metadata">

**Author:** ![lucifer1004](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lucifer1004/32/22311_2.png) [@lucifer1004](https://discourse.julialang.org/u/lucifer1004)\
**Post date:** [May 2, 2023, 11:21am UTC](https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191/3 "2023-05-02T11:21:32Z")

</div>

Thanks! What about the first question? How should we pass Fortran pointers correctly (if the Fortran code cannot be modified)?

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [May 2, 2023, 11:27am UTC](https://discourse.julialang.org/t/correct-way-to-ccall-a-fortran-subroutine-requiring-pointers/98191/4 "2023-05-02T11:27:01Z")

</div>

> [@lucifer1004](#):
>
> How should we pass Fortran pointers correctly (if the Fortran code cannot be modified)?

I would probably write, in Fortran, an interfacing module, following those directives, and call that module from Julia instead of the original one.
