2

I am trying to write a subroutine (for minimisation) that has two arguments:

  • an array x of any length
  • a function f that takes an array of that length and returns a scalar

example module:

module foo

contains

  subroutine solve(x, f)
    real, dimension(:), intent(inout) :: x
    interface
      real pure function f(y)
        import x
        real, dimension(size(x)), intent(in) :: y
      end function
    end interface

    print *, x
    print *, f(x)
  end subroutine

end module

and test program:

use foo

real, dimension(2) :: x = [1.0, 2.0]

call solve(x, g)

contains

  real pure function g(y)
    real, dimension(2), intent(in) :: y

    g = sum(y)
  end function

end

gfortran fails on:

call solve(x, g)
              1
Error: Interface mismatch in dummy procedure 'f' at (1): Shape mismatch in dimension 1 of argument 'y'

If I change size(x) => 2 then it compiles (and runs) fine. It also works fine if I change : => 2. But neither of these solutions gets me what I want.

Any ideas on how I can achieve this?

Vladimir F Героям слава
  • 57,977
  • 4
  • 76
  • 119
js947
  • 23
  • 1
  • 3

4 Answers4

5

How about:

interface
  real pure function f(y)
    real, dimension(:), intent(in) :: y
  end function
end interface

When you pass the argument of solve to the function, the size of the array will automatically be passed. You don't need to make this part of the interface.

M. S. B.
  • 28,968
  • 2
  • 46
  • 73
  • Also doesn't compile I'm afraid, with a similar error: `Error: Interface mismatch in dummy procedure 'f' at (1): Shape mismatch in argument 'y'` – js947 Apr 10 '13 at 13:26
  • 1
    all your dimensions should be `:`, also in the `g(y)` function (otherwise you're mixing explicit and assumed-shape arrays) – steabert Apr 10 '13 at 13:30
  • Problem solved I guess! I was hoping to gain safety by being explicit with with the size of the array in `g(y)`, but maybe that isn't really possible. – js947 Apr 10 '13 at 13:41
2

If you want to gain safety as indicated in your comment of M.S.B's solution, you should use -fcheck=bounds and the compiler will generate run-time checks for assumed and deferred shape arrays. See the gfortran man page for more info on -fcheck. However, you will lose some speed.

Community
  • 1
  • 1
bijancn
  • 462
  • 5
  • 11
1

You have the solution, but for what its worth an explanation... if the dummy argument has an explicit interface (which it does here) then there is a requirement that the characteristics of a procedure passed as an actual argument must match those of the dummy argument, with some exceptions around pureness and elemental intrinsics. The characteristics of a procedure include the characteristics of its dummy arguments, amongst other things.

The characteristics of a dummy argument include its shape, amongst other things. If that shape is not a constant expression - the characteristics include "the exact dependence [of the shape] on the entities in the expression".

The interface block for the dummy argument f declares the array to be of size SIZE(x). x is a host associated assumed shape variable - its size can vary at runtime, so SIZE(x) is not a constant. Hence that expression and the entities in it becomes a characteristic of the dummy argument.

The module procedure g declares the array to be of size 2. That is clearly a constant.

Regardless of the value of the non-constant expression for the size of the dummy argument of f, those array size characteristics (some sort of expression vs a constant) don't match - hence the error.

When you replace SIZE(x) with the constant 2 the characteristics obviously match. When you change the assumed shape x to be a constant size 2 - then SIZE(x) becomes a constant expression of value 2 - because it is a constant expression all that is relevant is its value - hence the characteristics of the two arguments then match. When you change both the dummy argument of f and the dummy argument of g to be assumed shape (:), the characteristics match.

IanH
  • 21,026
  • 2
  • 37
  • 59
0

Here is a demo to show how to pass allocatable array.

Some tips:

  • Use modules to avoid cumbersome interface.
  • Add extra matrix size information when pass array to the actual function. For example f(y, sizeinfo) so that inside your actual function you can declare the size of the input matrix correctly. The allocatable array can be passed to subroutine solve, so the size can be obtained using size(mat) in your subroutine solve.

So a corrected version looks like:

module foo

contains

  subroutine solve(x, f)
    real, dimension(:), intent(inout) :: x
    real,external::f
        integer::sizeinfo

    print *,'x=', x
    sizeinfo = size(x)
    print *, 'f(x)=',f(x,sizeinfo)
  end subroutine

    real function g(y,sizeinfo)
    integer::sizeinfo
    real, dimension(sizeinfo) :: y

    g = sum(y)
  end function
end module

Here is main program:

program main
use foo

real, dimension(2) :: x = (/1.0, 2.0/)

call solve(x, g)

end program

And the result is :

x=   1.000000       2.000000    
f(x)=   3.000000
Community
  • 1
  • 1
Kiyomine
  • 33
  • 5
  • This has the same problem like the other answer. The point of passing a subroutine is that you can pass *any* subroutine. External is old, obsolete and error prone. Using `external` for a function from the same module to save a few lines of code is a bad style, really. For more see my previous comment to http://stackoverflow.com/a/42356839/721644 – Vladimir F Героям слава Feb 21 '17 at 06:41
  • Good point. But, in that case, could you please give an EXAMPLE why external is error prone and what if you have a large amount of subroutines to be called? Would you have to write many and many `interface`s? – Kiyomine Feb 21 '17 at 15:07
  • Why would you? You write one interface and those subroutines all conform to that interface. That is the whole point! You write a library and the user of the library writes those subroutines which are passed. Therefore they cannot be in one module, because someone else writes them (even if both persons happen to be you in the end). – Vladimir F Героям слава Feb 21 '17 at 15:09
  • You are right if there are one or few `callsubs(sub, arg1,arg2)`. Sorry I wasn't clear because what I really wanted to ask is what if have a large amount of `callsubs` calling other `sub` s. Moreover, I still wonder why using module is error-prone. I would like to learn if there is a example showing using interface is indeed a good practice. – Kiyomine Feb 21 '17 at 15:15
  • No, using module is not error prone, using `external` is!!! and as I told you http://stackoverflow.com/a/42356839/721644 you can easily replace the `external` by `procedure(sub)`. And once you did that you will realize that `sub` acts as an interface, that you can just have an `abstract interface` in the module instead of `sub` and you are where you were before. – Vladimir F Героям слава Feb 21 '17 at 15:19
  • But this is not a good place to discuss here because this question is about explaining some error and you are not explaining an error, but showing different approach which I think does not answer the question asked. It is more on-topic where you answered first. – Vladimir F Героям слава Feb 21 '17 at 15:32