R help - Mar 2009 - Bug? FORTTRAN help

I was feeling masochistic the other day and we have been having some wierd
memory problems so I started digging into the source for L-BFGS-B. In the
lbgfsb.c file I see the following code:

/*     Cholesky factorization of (2,2) block of wn. */
    F77_CALL(dpofa)(&wn[*col + 1 + (*col + 1) * wn_dim1], &m2, col,
info);
    if (*info != 0) {
	*info = -2;
	return;
    }

If I am not mistaken this says that there is a m2 * col matrix that starts at
'col + 1 + (col + 1) * wn_dm1. Where wn_dm1 is 2 * m. My first question is
to verify that statement.
Say I am trying to optimize the "banana function" as given in the
documentation. In that case n = 2 and the default m = 5. So m2 is 10 and wn_dim1
is 20 and the dimension of wn is 100 (this is all by deduction. So if col is 5
then the offset into the array is 55 and there is not room in the vector for a
10 x 5 array. I am worried that the optimizer will silently write info memory
that it shouldn't but more than likely it is something that I don't
understand. So please vefify my first statement.

Thank you.

Kevin

G'day Kevin,

On Thu, 26 Mar 2009 13:42:20 -0700
<rkevinburton at charter.net> wrote:
> I was feeling masochistic the other day [...]
Welcome to the club. :)
> and we have been having some wierd memory problems so I started
> digging into the source for L-BFGS-B. In the lbgfsb.c file I see the
> following code:
> 
> /*     Cholesky factorization of (2,2) block of wn. */
>     F77_CALL(dpofa)(&wn[*col + 1 + (*col + 1) * wn_dim1], &m2, col,
> info); if (*info != 0) {
> 	*info = -2;
> 	return;
>     }
> 
> If I am not mistaken this says that there is a m2 * col matrix that
> starts at 'col + 1 + (col + 1) * wn_dm1. Where wn_dm1 is 2 * m.
I think your interpretation is not quite correct.  Note that it makes
only a sense to calculate a Cholesky factorization of a square matrix.
The interface of dpofa (in the linpack library) is available at:
     http://www.netlib.org/linpack/dpofa.f

Thus, the call above says, calculate the Cholesky factorization of a
col * col matrix whose (1,1) element is stored at &wn[*col+1+(*col+1)]
and that matrix is stored within a matrix which was allocated such that
it has m2 rows.

Or, in other words, calculate the Cholesky factorization of a col * col
matrix whose (1,1) element is stored at &wn[*col+1+(*col+1)] and to
move from the (1,1) element to the (1,2) element you have to move to
the memory location m2*sizeof(double) ahead/behind of (1,1).

Fortran uses a column major form to store arrays, i.e. element (1,1) is
followed by element (2,1), (3,1) and so forth.  To know where to find
element (1,2) of the matrix, you have to tell Fortran with how many
rows the "big matrix" that holds "your matrix" was
allocated.
> I am worried that the optimizer will silently write info memory that
> it shouldn't [...]
If you are worried about such issues, you should read chapter 4 of
"Writing R extensions", in particular Section 4.3 on gctorture and
valgrind.  Then run R on a platform that supports valgrind.  It is very
useful to catch problems such as accessing or writing into memory that
you should not access or write to.

HTH.

Cheers,
	
	Berwin

=========================== Full address ============================Berwin A
Turlach                            Tel.: +65 6516 4416 (secr)
Dept of Statistics and Applied Probability        +65 6516 6650 (self)
Faculty of Science                          FAX : +65 6872 3919       
National University of Singapore     
6 Science Drive 2, Blk S16, Level 7          e-mail: statba at nus.edu.sg
Singapore 117546                    http://www.stat.nus.edu.sg/~statba