Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Debug Fortran code from within R

I have a long piece of inherited FORTRAN77 code that I call from R using .Fortran(). The Fortran code contains a set of subroutines, and works when embedded in a Fortran program and subsequently compiled and run from the command line. However, when I call it from R, it crashes R the second time I call the function.

As the Fortran code uses a lot of indices and array dimensions stored as a variable, I reckon something goes wrong there. At some point, the Fortran code is looking somewhere in the memory where it's not supposed to be. So I need to step through the Fortran code and check whether all that came from R is what I think it is, and the code does what I think it does.

If it would be an R function, I'd have a choice of using debug(), adding browser() statements and printing out any value I'd like to see at one point in the code. But the Fortran code doesn't allow me any of these things afaik. If I understood it right, Fortran's screen output is not captured by R.

So has anybody an idea how exactly I can check the type and value of the arguments that R passes to the Fortran code. If you can explain how I can subsequently debug that code when called from R, that would be splendid.

Here's an example to illustrate what I mean.

C An example program
C
      PROGRAM EXAMPLE
      INTEGER N
      PARAMETER (N=10)
      REAL X0, X(N),MEAN

C
      X0 = 14
      DO 10 I = 1,10
         FI = FLOAT(I)
         X(I) = X0 + FI
   10 CONTINUE
      CALL MYSUB(X0,MEAN)
      END
C
C Mysub the subroutine
C
      SUBROUTINE MYSUB(X,N,MEAN)
      INTEGER N
      REAL X(N), MEAN
      MEAN = 0
      DO 20 I = 1,N
         MEAN = MEAN + X(I)
   20 CONTINUE
      MEAN = MEAN / N
      RETURN
      END

Say I want to call the subroutine mysub from R, and I want to make sure that I get X and N correctly. I use following function :

mysub <- function(x){
    if(!is.vector(x) | is.numeric(x)) stop("X has to be a numeric vector")
    n <- length(x)
    res <- .Fortran('mysub',X=as.single(x), N=as.integer(n), MEAN=single(1))
    return(res$MEAN)
}
like image 365
Joris Meys Avatar asked Oct 06 '22 08:10

Joris Meys


1 Answers

In your example program you call MYSUB with just 2 arguments instead of the 3 found in MYSUB's definition: (X, N, MEAN)

This probably has nothing to do with your question, but since you are asking about debugging FORTRAN and arguments, I thought I should point it out. FORTRAN subroutines are independent compilations. There is nothing to prevent you from passing incorrect number of arguments (no compiler error or linkage hint) and this can cause problems.

like image 60
Steve Rawlins Avatar answered Oct 13 '22 00:10

Steve Rawlins