CODEWITHSUNDEEP
fortranoperator-overloadingfortran90dynamic-arrays
jorispilot
jorispilot

Reputation: 483

Fortran's interface operator behavior on allocatable arrays

I have a derived type that contains a set allocatable arrays and I'm trying to overload some operators. As my arrays can become very large I don't want Fortran to do implicit copy of my arrays, but I can't figure out what is it actually done. Here is a minimal program that illustrate my point:

module type_vector_field
  implicit none


  type vector_field
     double precision, dimension(:), allocatable :: f
  end type vector_field

  interface operator(+)
     module procedure vector_field_add
  end interface


contains


  function vector_field_add(vf1, vf2) result(vf3)
    type(vector_field), intent(in) :: vf1, vf2
    type(vector_field) :: vf3

    print*, "vf* addresses in add:  ", loc(vf1), loc(vf2), loc(vf3)
    vf3%f = vf1%f + vf2%f

  end function vector_field_add


  subroutine vector_field_add2(vf1, vf2, vf3)
    type(vector_field), intent(in) :: vf1, vf2
    type(vector_field), intent(out) :: vf3

    print*, "vf* addresses in add2: ", loc(vf1), loc(vf2), loc(vf3)
    vf3%f = vf1%f + vf2%f

  end subroutine vector_field_add2


end module type_vector_field

The addition of two vector_field can be done either by the + operator or by the vector_field_add2 subroutine. I tested this with the following sample program:

program main
  use type_vector_field
  implicit none

  integer :: n = 6283
  type(vector_field) :: vf1, vf2, vf3

  allocate(vf1%f(n), vf2%f(n))
  vf1%f = 3.1415
  vf2%f = 3.1415

  !! first version

  allocate(vf3%f(n))
  print*, "vf* addresses in main: ", loc(vf1), loc(vf2), loc(vf3)

  vf3 = vf1 + vf2
  print*, "vf* addresses in main: ", loc(vf1), loc(vf2), loc(vf3)
  print*, vf3%f(n)

  deallocate(vf3%f)

  !! second version

  allocate(vf3%f(n))
  print*, "vf* addresses in main: ", loc(vf1), loc(vf2), loc(vf3)

  call vector_field_add2(vf1, vf2, vf3)
  print*, "vf* addresses in main: ", loc(vf1), loc(vf2), loc(vf3)
  print*, vf3%f(n)

end program main

And this is what it outputs when compiled with gfortran-4.8.2:

$ gfortran -g -Wall test.f90 -o test && test
 vf* addresses in main:           6303968          6304032          6304096
 vf* addresses in add:            6303968          6304032  140733663003232
 vf* addresses in main:           6303968          6304032          6304096
   6.2829999923706055     
 vf* addresses in main:           6303968          6304032          6304096
 vf* addresses in add2:           6303968          6304032          6304096
 vf* addresses in main:           6303968          6304032          6304096
   6.2829999923706055

It seems that the subroutine works directly with the vf3 instance defined in the main scope, while the operator create an other vf3 instance in its scope and copy back the result. When I print the arrays addresses instead of the derived type addresses, I get:

 vf*%f addresses in main:          39440240         39490512         39540784
 vf*%f addresses in add:           39440240         39490512                0
 vf*%f addresses in main:          39440240         39490512         39591056
   6.2829999923706055     
 vf*%f addresses in main:          39440240         39490512         39540784
 vf*%f addresses in add2:          39440240         39490512                0
 vf*%f addresses in main:          39440240         39490512         39540784
   6.2829999923706055

Hence the operator seems to copy the address of its vf3%f in the vf3%f instance of the main scope, while the subroutine works directly with the vf3%f of the main scope. In both cases, I don't understand where points the vf3%f of both the subroutine and the operator.

My questions:

Upvotes: 1

Views: 317

Answers (1)

user4490638
user4490638

Reputation:

The problem you encounter has, as far as I can see, nothing to do with operator overloading.

More or less by accident (I presume) you had INTENT(OUT) in your subroutine vector_field_add2, and you don't have it in vector_field_add.

One often-overlooked effect of INTENT(OUT) is that all allocatable arrays get deallocated if the dummy arguments have the allocatable attribute upon entry into vector_field_add2. The same goes for types having components with the allocatable attribute.

If you have

   real, dimension(:), allocatable: my_array

   allocate(my_array(some:thing))

and then do

 call sub1(my_array)

with

  sub1(a)
    real, dimension(:), intent(out) :: a

nothing surprising happens.

If you do instead

   call sub2(a)

where sub2 is

   subroutine sub2(a)
   real, dimension(:), allocatable, intent(out) :: a
   end subroutine sub2

the call to sub2 is equivalent to a call to deallocate. If you leave out the intent(out) or replace it, e.g. with intent(inout), nothing happens.

The same thing happened to you, only you had the allocatable component inside your type, so it was even less visible.

Otherwise, there is no difference if you use an operator or if you call a subroutine directly. It only matters for clarity of the code. The compiler translates the assignment statement with the overloaded operator to the subroutine call.

If you want to see how it is done, translate your program using the -fdump-tree-original option. This will generate a file with the main name like your source file, with an added .003.original, where you can see a C-like representation of your Fortran code.

Upvotes: 1

Related Questions