fortran result variable Not initiallized

Question

I met a surprsing problem about loacal variable initializing.

I got following function to calculate gammar

function gammar(z) result(gz)
implicit none
real(8),intent(out)::gz
real(8)::z,t,low,up
real(8),parameter::increment=1.0
real(8),parameter::lower_t=0.0,upper_t=10.0
integer(4)::i,n
!gz=0.0
n=(upper_t-lower_t)/increment
do i=1,n
low=lower_t+(i-1)*increment
up=lower_t+(i)*increment
gz=gz+(f(z,low)+f(z,up))*increment/2.0
end do
end function gammar

Then I call this function in main program like
df=9.0
t=0.0
write(*,*) gammar((df+1.0)/2.0)/sqrt(pi*df)/gammar(df/2.0)

I got wrong answer!! 0.126 I found the reason was after gammar((df+1.0)/2.0) was calculated, the local variable gz was not set to 0. Hence ,when calculate gammar(df/2.0), the gz still retained old value 24. Eventually,gammar(df/2.0) got wrong answer 34.. If I add gz=0.0 in the gammar function, this problem was fixed. This is really surprising. Why the local gz was not initiallized to zero when the gammar called every time?

Thanks a lot

Regards Ke

Answer 1

Unless you have a statement to initialize a local variable in a procedure, such as the gz = 0 that you have commented out, those local variables are not initialized when the procedure is invoked. Their values are undefined. They could have a value left from a previous invocation, or some random value.

If you use full warning options of the compiler, it will likely tell you of this problem. gfortran warned of an uninitialized variable at compile time. ifort detected the problem at run time.

Another initialization method is with a declaration. That still doesn't repeat the initialization of additional invocations of the procedure. If you initialize a local variable in a procedure with a declaration, such as integer :: count = 0, that initialization is only done on the first invocation on the procedure. But ... the variable remains defined and on the next invocation will retain that value that it had when the previous invocation exited.

P.S. real(8) is not a portable way to obtain double precision reals. The language standard does not specify specific numeric values for kinds ... compilers are free to whatever value they wish. Most compilers use the number of bytes, but use other numbering methods. It is better to use selected_real_kind or the ISO_FORTRAN_ENV and (for double precision) real64. See quad precision in gfortran

P.P.S. Trying this code with gfortran, that compiler points out another problem with gz:

function gammar(z) result(gz)
                             1
Error: Symbol at (1) is not a DUMMY variable

So remove the intent(out) in the declaration.