How to get priorly-unknown array as the output of a function in Fortran

Question

In Python:

def select(x):
    y = []
    for e in x:
        if e!=0:
            y.append(e)
    return y

that works as:

x = [1,0,2,0,0,3]
select(x)
[1,2,3]

to be translated into Fortran:

function select(x,n) result(y)
    implicit none
    integer:: x(n),n,i,j,y(?)
    j = 0
    do i=1,n
        if (x(i)/=0) then
            j = j+1
            y(j) = x(i)
        endif
    enddo
end function

The questions are in Fortran:

1. how to declare y(?)?
2. how to declare predefined values for x
3. how to avoid dimension info n

for 1 if it is defined as y(n) the output will be:

x = (/1,0,2,0,0,3/)
print *,select(x,6)
1,2,3,0,0,0

which is not desired!
!-------------------------------
Comments:
1- All given answers are useful in this post. Specially M.S.B and eryksun's.
2- I tried to adapt the ideas for my problem and compile with F2Py however it was not successful. I had already debugged them using GFortran and all were successful. It might be a bug in F2Py or something that I don't know about using it properly. I will try to cover this issue in another post.

Update: A linked question could be found at here.

Solution 1

I hope a real Fortran programmer comes along, but in the absence of better advice, I would only specify the shape and not the size of x(:), use a temporary array temp(size(x)), and make the output y allocatable. Then after the first pass, allocate(y(j)) and copy the values from the temporary array. But I can't stress enough that I'm not a Fortran programmer, so I can't say if the language has a growable array or if a library exists for the latter.

program test
    implicit none
    integer:: x(10) = (/1,0,2,0,3,0,4,0,5,0/)
    print "(10I2.1)", select(x)

contains

    function select(x) result(y)
        implicit none
        integer, intent(in):: x(:) 
        integer:: i, j, temp(size(x))
        integer, allocatable:: y(:)

        j = 0
        do i = 1, size(x)
            if (x(i) /= 0) then
                j = j + 1
                temp(j) = x(i)
            endif
        enddo

        allocate(y(j))
        y = temp(:j)
    end function select

end program test

---

Edit:

Based on M.S.B.'s answer, here's a revised version of the function that grows temp y with over-allocation. As before it copies the result to y at the end. It turns out i's not necessary to explicitly allocate a new array at the final size. Instead it can be done automatically with assignment.

    function select(x) result(y)
        implicit none
        integer, intent(in):: x(:) 
        integer:: i, j, dsize
        integer, allocatable:: temp(:), y(:)

        dsize = 0; allocate(y(0))

        j = 0
        do i = 1, size(x)
            if (x(i) /= 0) then
                j = j + 1

                if (j >= dsize) then         !grow y using temp
                    dsize = j + j / 8 + 8 
                    allocate(temp(dsize))
                    temp(:size(y)) = y
                    call move_alloc(temp, y) !temp gets deallocated
                endif

                y(j) = x(i)
            endif
        enddo
        y = y(:j)
    end function select

Solution 2

Here is an example of a Fortran function returning a variable length array. This is a feature of Fortran 2003. Also used in the test driver is automatic allocation on assignment, another Fortran 2003 feature.

module my_subs

contains

function select(x) result(y)
    implicit none
    integer, dimension (:), intent (in) :: x
    integer, dimension (:), allocatable :: y
    integer :: i, j

    j = 0
    do i=1, size (x)
        if (x(i)/=0) j = j+1
    enddo

    allocate ( y (1:j) )

    j = 0
    do i=1, size (x)
        if (x(i)/=0) then
            j = j+1
            y(j) = x(i)
        endif
    enddo

    return

end function select

end module my_subs

program test

use my_subs

implicit none
integer, dimension (6) :: array = [ 5, 0, 3, 0, 6, 1 ]
integer, dimension (:), allocatable :: answer

answer = select (array)

write (*, *) size (array), size (answer)
write (*, *) array
write (*, *) answer

stop


end program test

Here is an alternative solution that uses a temporary array to "grow" the output array (function return) as needed. While two passes through the input array are avoided, array copies are required. Another Fortran 2003 feature, move_alloc, reduces the number of copies needed. move_alloc also takes care of the (re)allocation of the output array (here "y") and deallocation of the input array (here "temp"). Perhaps this is more elegant, but it is probably less efficient since multiple copies are used. This version is probably more educational then useful. @eryksun's version uses one pass and one copy, at the expense of making the temporary array full size.

function select(x) result(y)
    implicit none
    integer, dimension (:), intent (in) :: x
    integer, dimension (:), allocatable :: y, temp
    integer :: i, j

    j = 0
    do i=1, size (x)
        if (x(i)/=0) then
            j = j+1
            allocate (temp (1:j))
            if ( allocated (y) ) temp (1:j-1) = y
            call move_alloc (temp, y)
            y(j) = x(i)
        endif
    enddo

    return

end function select

Solution 3

If the example in your question really is what you want to do, you can use the Fortran90 intrinsic `pack':

program pack_example

implicit none

integer, dimension(6) :: x

x = (/ 1,0,2,0,0,3 /)

! you can also use other masks than 'x/=0'
write(*,*) pack(x, x/=0)

end program pack_example

The output of the example program is: 1 2 3

Comments

[Eryk Sun]

Doing 2 full passes is another way to go. I chose a temporary array thinking the work done in the function (not this one, but in general) is probably more expensive than copying the result. Also, aren't allocatable arrays a Fortran 90/95 feature?

[M. S. B.]

allocatable arrays are part of Fortran 90. There have been improvements with Fortran 95 and 2003. In this example, a function return that is an allocatable was not supported prior to Fortran 2003. As was the automatic allocation that takes place in the assignment statement in the test driver, "answer = select (array)".

[Developer]

Thanks for the hint. The example in the question was however just for demonstration of the problem. Anyway your hint is very useful for other simple situations as you provided. Indeed the question refers to how pack does implement returning varying sized array as output depending the condition e.g., 'x/=0' based on x without given size info! Another point is that 'pack' output is a subset of input 'x' and cannot be larger in size, for example. After all it is a good hint.

[Developer]

Thanks for sharing your knowledge with the community. Your answer and M.S.B's were very close to each other and both very helpful. By the way, in part j+j/8+8 was that for size of integer so it require changes for other type of data!

[Developer]

Thanks for your answer. It was however very close to eryksun's answer. This post also was useful to me helping to figure out some improvements in Fortran 2003.

[Eryk Sun]

@Supporter: Maybe that wasn't clear. I added the update as an alternative to M.S.B.'s alternate answer. Instead of growing the array one by one, which I think involves an excess of copying, I had it over-allocate. The current size needed is j items, so I allocate that plus an eighth more (also plus a constant 8 for when it's small). Thus the size of the over-allocation grows in proportion to the current size of the array, which minimizes the number of times it has to be reallocated. At the end I shrink it back down, which hopefully is optimized by the compiler to not require a copy.

[Developer]

To me both answers are helpful and correct. Just I double-checked the time that answers had been sent and then figured out you were the first one (7.19 < 7.26). That is so your answer is chosen now.

[knia]

In the first code block, you could also declare the result as integer :: y(count(x /= 0)). (This may actually be valid Fortran 95, but I'm not sure about the automatic function result.) I'm guessing the program will build an intermediate logical array though, depending on how clever the compiler is.