finding specific indices with pointer array

2020-04-20 09:01发布

问题:

I am relatively new to Fortran and break my head about one thing for hours now:

I want to write a subroutine for finding the indexes for specific elements in a real 1D array (given to the routine as input).

I have generated an array with 100 random reals, called arr, and now want to determine the indexes of those elements which are greater than a real value min, which is also passed to subroutine.

Plus, in the end I would like to have a pointer I'd allocate in the end, which I was said would be better than using an array indices containing the indexes once found.

I just didn't find how to solve that, I had following approach:

SUBROUTINE COMP(arr, min)
   real, intent(in)                 :: arr(:)
   real, intent(in)                 :: min
   integer, pointer, dimension(:)   :: Indices
   integer                          :: i, j

   ! now here I need a loop which assigns to each element of the pointer 
   ! array the Indices one after another, i don't know how many indices there
   ! are to be pointed at 
   ! And I dont know how to manage that the Indices are pointed at one after another,
   ! like Indices(1) => 4
   !      Indices(2) => 7
   !      Indices(3) => 32
   !      Indices(4) => 69
   !      ...
   ! instead of
   !      Indices(4) => 4
   !      Indices(7) => 7
   !      Indices(32) => 32
   !      Indices(69) => 69
   !      ...


   DO i = 1, size(arr)
      IF (arr(i) > min) THEN
         ???
      ENDIF
   ENDDO
allocate(Indices)
END SUBROUTINE COMP

回答1:

I think you're on the right track, but you're ignoring some intrinsic Fortran functions, specifically count, and you aren't returning anything!

subroutine comp(arr, min)
   real, intent(in) :: arr(:)
   real, intent(in) :: min
! local variables
   integer, allocatable :: indices(:)
   integer :: i,j, indx


! count counts the number of TRUE elements in the array 
   indx = count(arr > min)
   allocate(indices(indx))

! the variable j here is the counter to increment the index of indices
   j=1
   do i=1,size(arr)
      if(arr(i) > min) then
         indices(j) = i
         j = j+1
      endif
   enddo
end subroutine comp

Then you can use the array indices as

do i=1,size(indices)
   var = arr(indices(i))
enddo

Note that since you are not returning indices, you will lose all the information found once you return--unless you plan on using it in the subroutine, then you're okay. If you're not using it there, you could define it as a global variable in a module and the other subroutines should see it.



回答2:

If succinctness (rather than performance) floats your boat... consider:

FUNCTION find_indexes_for_specific_elements_in_a_real_1D_array(array, min)  &
    RESULT(indices)
  REAL, INTENT(IN) :: array(:)
  REAL, INTENT(IN) :: min
  INTEGER, ALLOCATABLE :: indices(:)
  INTEGER :: i
  indices = PACK([(i,i=1,SIZE(array))], array >= min)
END FUNCTION find_indexes_for_specific_elements_in_a_real_1D_array

[Requires F2003. Procedures that have assumed shape arguments and functions with allocatable results need to have an explicit interface accessible where they are referenced, so all well behaved Fortran programmers put them in a module.]



回答3:

A simple way to get the indices of a rank 1 array arr for elements greater than value min is

indices = PACK([(i, i=LBOUND(arr,1), UBOUND(arr,1))], arr.gt.min)

where indices is allocatable, dimension(:). If your compiler doesn't support automatic (re-)allocation than an ALLOCATE(indices(COUNT(arr.gt.min)) would be needed before that point (with a DEALLOCATE before that if indices is already allocated).

As explanation: the [(i, i=...)] creates an array with the numbers of the indices of the other array, and the PACK intrinsic selects those corresponding to the mask condition.

Note that if you are doing index calculations in a subroutine you have to be careful:

subroutine COMP(arr, min, indices)
  real, intent(in) :: arr(:)
  real, intent(in) :: min
  integer, allocatable, intent(out) :: indices(:)

  !...
end subroutine

arr in the subroutine will have lower bound 1 regardless of the bounds of the actual argument (the array passed) (which could be, say VALS(10:109). You will also then want to pass the lower bound to the subroutine, or address that later.

[Automatic allocation is not an F90 feature, but in F90 one also has to think about allocatable subroutine arguments