Return Array from Function

[bebbo1986]

Hello,
I wrote this programm which tries to use the Function MatrixVektor.

PROGRAM MAIN
REAL, DIMENSION(3) :: v
REAL, DIMENSION(3) :: w
REAL, DIMENSION(3,3) :: A

... filling v and A

w = MatrixVektor(A, v, w)
END PROGRAMM

FUNCTION MatrixVektor(A, v)
REAL, DIMENSION(3) :: MatrixVektor
REAL, DIMENSION(3,3) :: A
REAL, DIMENSION(3) :: v

do i=1,3
MatrixVektor(i) = SUM(A(i,1:3)*v)
enddo
END FUNCTION

But there follows an error: Attempt to call a subroutine as if it were an integer(kind=3) function. In the Line where the Function is called in the main programm. Why? What do I have to change.
Thanks for every help.

Regards,
Martin

 

[gullipe]

Hi, Martin

There are three arguments in the call (A,v,w), but only two in the function (A,v). Drop w from the call.

I am also not sure what the function should do:
The sentence "MatrixVektor(i) = SUM(A(i,1:3)*v)" seems strange to me.

Maybe this is correct:

PROGRAM MAIN
REAL, DIMENSION(3) :: v
REAL, DIMENSION(3) :: w
REAL, DIMENSION(3,3) :: A

... filling v and A

w = MatrixVektor(A,v)
END PROGRAM

FUNCTION MatrixVektor(A,v)
REAL, DIMENSION(3) :: MatrixVektor
REAL, DIMENSION(3,3) :: A
REAL, DIMENSION(3) :: v

do i=1,3
MatrixVektor(i) = a(i,1)*v(1)+a(i,2)*v(2)+a(i,3)*v(3)
enddo
END FUNCTION

 

[bebbo1986]

Thanks for your answer.
I tried as you suggested but the problem stayed the same.
The function should deliver the product of the matrix A with the vector v.
I think it should work with SUM(A(i,1:3)*v) because it's the same only in short.

 

[gullipe]

Hi, Marin

OK, I looked up the SUM function and you are right of course.
To solve the problem, try to use SUBROUTINE instead of FUNCTION:

PROGRAM MAIN
REAL, DIMENSION(3) :: v
REAL, DIMENSION(3) :: w
REAL, DIMENSION(3,3) :: A

... filling v and A

call MatrixVektor(A,v,w)
END PROGRAM

SUBROUTINE MatrixVektor(A,v,w)
REAL, DIMENSION(3) :: w
REAL, DIMENSION(3,3) :: A
REAL, DIMENSION(3) :: v

do i=1,3
w(i) = SUM(A(i,1:3)*v)
enddo
END

 

[bebbo1986]

Your right that helps with the problem. But that was also the reason, why I had (A, v, w) in the text before. Because I tried that already and of course it helps.
But why can't it be done via function?

 

[gullipe]

Are sure that FUNCTION can return more than one variable or more than one dimensional array? In your original program it returns three dimensional array. I have written numerous FUNCTION programs (in Fortran 77 though), but I always declared them as INTEGER or REAL, not INTEGER(3) or REAL(3) for example. Maybe some F90 or F95 experts can answer this.

 

[mikrom]

Here is the function which returns vector

[COLOR=#a020f0]PROGRAM[/color] main
[COLOR=#2e8b57][b]    REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color]([COLOR=#ff00ff]2[/color])   :: v, w
[COLOR=#2e8b57][b]    REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color]([COLOR=#ff00ff]2[/color],[COLOR=#ff00ff]2[/color]) :: A

[COLOR=#2e8b57][b]    REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color]([COLOR=#ff00ff]3[/color])   :: t, u [COLOR=#804040][b]=[/b][/color]([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color]., [COLOR=#ff00ff]1[/color]., [COLOR=#ff00ff]1[/color].[COLOR=#804040][b]/[/b][/color])
[COLOR=#2e8b57][b]    REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color]([COLOR=#ff00ff]3[/color],[COLOR=#ff00ff]3[/color]) :: B [COLOR=#804040][b]=[/b][/color] [COLOR=#008080]reshape[/color](([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].[COLOR=#804040][b]/[/b][/color]),[highlight #ffff00][COLOR=#0000ff]&[/color][/highlight]
                                      [highlight #ffff00][COLOR=#0000ff]&[/color][/highlight] ([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]3[/color], [COLOR=#ff00ff]3[/color][COLOR=#804040][b]/[/b][/color])) [COLOR=#0000ff]! fill columns and reshape[/color]

    [COLOR=#0000ff]!filling v and A[/color]
    A([COLOR=#ff00ff]1[/color],[COLOR=#ff00ff]1[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color]. 
    A([COLOR=#ff00ff]1[/color],[COLOR=#ff00ff]2[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]2[/color].
    A([COLOR=#ff00ff]2[/color],[COLOR=#ff00ff]1[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]3[/color].
    A([COLOR=#ff00ff]2[/color],[COLOR=#ff00ff]2[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]4[/color].
    v([COLOR=#ff00ff]1[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color].
    v([COLOR=#ff00ff]2[/color])[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]2[/color].

    [COLOR=#0000ff]! print *, A ! print columns[/color]
    w [COLOR=#804040][b]=[/b][/color] MatrixVector(A, v)
    [COLOR=#804040][b]print[/b][/color] [COLOR=#804040][b]*[/b][/color], [COLOR=#ff00ff]"u = "[/color], (w(j), j [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], [COLOR=#ff00ff]2[/color])

    [COLOR=#0000ff]! print *, B ! print columns[/color]
    t [COLOR=#804040][b]=[/b][/color] MatrixVector(B, u)
    [COLOR=#804040][b]print[/b][/color] [COLOR=#804040][b]*[/b][/color], [COLOR=#ff00ff]"t = "[/color], (t(j), j [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], [COLOR=#ff00ff]3[/color])

[COLOR=#a020f0]CONTAINS[/color]    
[COLOR=#a020f0]FUNCTION[/color] MatrixVector(M, x)
[COLOR=#2e8b57][b]  REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color](:,:) :: M
[COLOR=#2e8b57][b]  REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color](:) :: x
  [COLOR=#2e8b57][b]INTEGER[/b][/color] ::n,i,j
[COLOR=#2e8b57][b]  REAL[/b][/color], [COLOR=#2e8b57][b]DIMENSION[/b][/color](n) :: MatrixVector
  n[COLOR=#804040][b]=[/b][/color][COLOR=#008080]SIZE[/color](x)
  [COLOR=#0000ff]!print *, n[/color]
  [COLOR=#804040][b]do[/b][/color] i[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n
    MatrixVector(i)[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]0[/color].
    [COLOR=#804040][b]do[/b][/color] j[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n
      [COLOR=#0000ff]! MatrixVector(i) = SUM(A(i,1:n)*x) this gives wrong results ![/color]
      MatrixVector(i) [COLOR=#804040][b]=[/b][/color] MatrixVector(i) [COLOR=#804040][b]+[/b][/color] M(i,j)[COLOR=#804040][b]*[/b][/color]x(j)
    [COLOR=#804040][b]enddo[/b][/color]
  [COLOR=#804040][b]enddo[/b][/color]
[COLOR=#a020f0]END FUNCTION[/color] MatrixVector

[COLOR=#a020f0]END PROGRAM[/color] main

Notice: Using the function SUM() delivered bad result with the second example, so I replaced it with the regular Matrix-Vector multiplication.

Compilation (gfortran with MingW):

$ gfortran vector.f90 -o vector

Running:

$ vector
 u =    5.0000000       11.000000    
 t =    3.0000000       6.0000000       9.0000000

However, it rocks only with this magic word CONTAINS. Without it the program crashes. I don't know why it is so. IMHO, with such secret rules or behaviour Fortran could candidate to be the worst programming language

 

[mikrom]

The complications with the vector function can be treated similar to this

http://www.tek-tips.com/viewthread.cfm?qid=1512160&page=1

You can avoid all such issues and unexpected behaviors of the Fortran compiler introducing a hygienic writing style using Fortran modules.

So, if you don't want to bother use modules - all over!

Here is the example with the usage of module.

[COLOR=#a020f0]module[/color] vector_functions
  [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#a020f0]contains[/color]    
  [COLOR=#a020f0]function[/color] MatrixVector(M, x)
    [COLOR=#0000ff]! Matrix-Vector multiplication /for square Matrices only/[/color]
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color](:,:), [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: M
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color](:), [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: x
    [COLOR=#2e8b57][b]integer[/b][/color] ::n, i, j
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color](n) :: MatrixVector
    n[COLOR=#804040][b]=[/b][/color][COLOR=#008080]SIZE[/color](x)
    [COLOR=#0000ff]!print *, n[/color]
    [COLOR=#804040][b]do[/b][/color] i[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n
      MatrixVector(i)[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]0[/color].
      [COLOR=#804040][b]do[/b][/color] j[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n
        MatrixVector(i) [COLOR=#804040][b]=[/b][/color] MatrixVector(i) [COLOR=#804040][b]+[/b][/color] M(i,j)[COLOR=#804040][b]*[/b][/color]x(j)
      [COLOR=#804040][b]enddo[/b][/color]
    [COLOR=#804040][b]enddo[/b][/color]
  [COLOR=#a020f0]end function[/color] MatrixVector
[COLOR=#a020f0]end module[/color] vector_functions

[COLOR=#a020f0]program[/color] main
    [COLOR=#a020f0]use[/color] vector_functions

    [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]

[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color]([COLOR=#ff00ff]2[/color])   :: v, w
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color]([COLOR=#ff00ff]2[/color],[COLOR=#ff00ff]2[/color]) :: A
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color]([COLOR=#ff00ff]3[/color])   :: u, t 
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]dimension[/b][/color]([COLOR=#ff00ff]3[/color],[COLOR=#ff00ff]3[/color]) :: B
    [COLOR=#2e8b57][b]integer[/b][/color] :: j

    [COLOR=#0000ff]!****************************************************************[/color]
    [COLOR=#0000ff]!*** Example 1: [/color]
    [COLOR=#0000ff]!****************************************************************[/color]
    [COLOR=#0000ff]! construct Vector[/color]
    v[COLOR=#804040][b]=[/b][/color]([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color]., [COLOR=#ff00ff]2[/color].[COLOR=#804040][b]/[/b][/color])
    [COLOR=#0000ff]!print *, v ! print Vector[/color]
    [COLOR=#0000ff]! construct Matrix by columns and reshape[/color]
    A [COLOR=#804040][b]=[/b][/color] [COLOR=#008080]reshape[/color](([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]3[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]4[/color].[COLOR=#804040][b]/[/b][/color]), ([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color], [COLOR=#ff00ff]2[/color][COLOR=#804040][b]/[/b][/color])) 
    [COLOR=#0000ff]!print *, A ! print Matrix by columns[/color]
    w [COLOR=#804040][b]=[/b][/color] MatrixVector(A, v)
    [COLOR=#804040][b]print[/b][/color] [COLOR=#804040][b]*[/b][/color], [COLOR=#ff00ff]"u = "[/color], (w(j), j [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], [COLOR=#ff00ff]2[/color])

    [COLOR=#0000ff]!****************************************************************[/color]
    [COLOR=#0000ff]!*** Example 2: [/color]
    [COLOR=#0000ff]!****************************************************************[/color]
    [COLOR=#0000ff]! construct Vector[/color]
    u [COLOR=#804040][b]=[/b][/color]([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color]., [COLOR=#ff00ff]1[/color]., [COLOR=#ff00ff]1[/color].[COLOR=#804040][b]/[/b][/color])
    [COLOR=#0000ff]!print *, u ! print Vector[/color]
    [COLOR=#0000ff]! construct Matrix by columns and reshape[/color]
    B [COLOR=#804040][b]=[/b][/color] [COLOR=#008080]reshape[/color](([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].,[COLOR=#ff00ff]1[/color].,[COLOR=#ff00ff]2[/color].,[COLOR=#ff00ff]3[/color].[COLOR=#804040][b]/[/b][/color]), ([COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]3[/color], [COLOR=#ff00ff]3[/color][COLOR=#804040][b]/[/b][/color]))
    [COLOR=#0000ff]!print *, B ! print Matrix by columns[/color]
    t [COLOR=#804040][b]=[/b][/color] MatrixVector(B, u)
    [COLOR=#804040][b]print[/b][/color] [COLOR=#804040][b]*[/b][/color], [COLOR=#ff00ff]"t = "[/color], (t(j), j [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color], [COLOR=#ff00ff]3[/color])
[COLOR=#a020f0]end program[/color] main

Compile:

$ gfortran vector.f90 -o vector

Run:

$ vector
 u =    5.0000000       11.000000    
 t =    3.0000000       6.0000000       9.0000000

 

[gullipe]

Hi, mikrom & Martin

I seem to be very lucky that my Fortran 77 compiler does not understand this magic CONTAINS statement!!

 

[mikrom]

Hi gulippe,
Yes Fortran 77 is the Fortran! It's simple and powerful!

The new Fortrans 90/95 have some advantages regarding the new syntax (essential part of which is however supported e.g. by g77 compiler too - as enhacement), but such a tricky issues like this one, working only using magic of CONTAINS with or without module, make a programmer unhappy.