Fortran Program Fails Depending On A Write Statement Before Subroutine Call

It's been a number of years since I've worked with Fortran, so maybe I'm missing a fundamental issue, but here it goes. I'm not even sure how to properly describe this issue, so I

Solution 1:

In the subroutine CURVE_DERIVS_ALG1, the dummy argument CK seems not initialized, so could you check its initial value? If I set it to 0.0d0 before entering the loop, the code seems to work fine, but I am not sure if this initial value is OK. (Please also note that if INTENT(OUT) is given, all the elements must be defined.)

SUBROUTINE CURVE_DERIVS_ALG1(N, P, UK, CPW, U, D, CK)
    ...
    DOUBLE PRECISION, INTENT(OUT) :: CK(0:D, 0:3)
    ...    
    CALL FIND_SPAN(N, P, U, UK, SPAN)
    CALL DERS_BASIS_FUNS(SPAN, U, P, DU, UK, M, NDERS)

    CK(:,:) = 0.0d0          !! <--- here

    DO K = 0, DU
        DO J = 0, P
            CK(K, :) = CK(K, :) + NDERS(K, J) * CPW(SPAN - P + J, :)
            ...

Another potential issue is

IF (U .EQ. UK(N + 1)) THEN

which compares two floating-point numbers. Although this condition seems not met in this program, it is probably safer to rewrite this as, e.g.

IF ( abs( U - UK(N + 1) ) < 1.0d-10 ) THEN    !! threshold depends on your need...

---

EDIT: To detect the above error of CK automatically, it may be useful to compile as

gfortran -finit-real=snan -ffpe-trap=invalid evaluate.f90 main.f90

which gives (with gfortran4.8 on Linux x86_64)

Program received signal 8 (SIGFPE): Floating-point exception.

Backtrace for this error:
#0  0x00000039becac5f4 inwait () from /lib64/libc.so.6#1  0x00000039c501400d in ?? () from /usr/lib64/libgfortran.so.3#2  0x00000039c501582e in ?? () from /usr/lib64/libgfortran.so.3#3  0x00000039c50146ca in ?? () from /usr/lib64/libgfortran.so.3#4  <signal handler called>#5  0x0000000000401787 in __evaluate_MOD_curve_derivs_alg1 ()    <--- here#6  0x0000000000400fce in __evaluate_MOD_rat_curve_derivs ()#7  0x0000000000402b26 in MAIN__ ()#8  0x0000000000402cbb in main ()