Re: st: Question concerning pointers in Mata

[Matthew J Baker <matthew.baker@hunter.cuny.edu>]

Dear Listers --

Thanks to Alan for the very insightful answer! By way of a quick 
follow up, the reason why I asked the question is because 
sometimes one wishes to do a series of complicated 
computations, and then collect results systematically in a 
pointer. In these instances, it is easy to get tripped up on the 
aforementioned problem. 

Anyways, given the explanation, my thought is that a good 
approach is to put the "complicated" computations in a stand-
alone routine and reference this with the pointer, along the lines 
of the following example (the computations aren't complicated in 
the example, but hopefully the idea is clear):

// Begin
clear all
mata
real matrix foo(real scalar k)
{
real matrix Z
X=J(k,k,1)
return(X)
}
Z=J(10,1,NULL)
for (i=1;i<=10;i++) {
Z[i]=&foo(i)
}
*Z[1]
*Z[10]
// End example

This seems to me to be a systematic way of being sure that the 
possibility of a case of mistaken pointer identity is avoided, but 
there may be alternative means of dealing with the issue. If so, 
I'd like to hear about them!

Matt Baker

Dr. Matthew J. Baker
Department of Economics
Hunter College and the Graduate Center, CUNY

---- Original message ----
>Date: Thu, 24 Mar 2011 02:12:08 -0500
>From: owner-statalist@hsphsun2.harvard.edu (on behalf of 
Alan Riley <ariley@stata.com>)
>Subject: Re: st: Question concerning pointers in Mata  
>To: statalist@hsphsun2.harvard.edu
>
>Matthew J. Baker (matthew.baker@hunter.cuny.edu) asked a 
question
>about pointers in Mata:
>
>> I'm curious why the following two snippets of mata code 
produce
>> different results. The observation is that if one defines a
>> temporary matrix and then "points" to this matrix, the result 
is
>> different than if one just simply "points" to the definition.
>> Accordingly, I'm a bit curious about how exactly the pointer 
and
>> the pointed is storing results. Here is an example:
>> 
>> /* Beginning of example */
>> 
>>    mata
>>    Y=J(10,1,NULL)
>>    for (k=1;k<=10;k++) {
>>       Y[k]=&J(k,k,1)
>>    }
>>    Z=J(10,1,NULL)
>>    for (k=1;k<=10;k++) {
>>       Temp=J(k,k,1)
>>       Z[k]=&Temp
>>    }
>>    /* Compare contents of Y[1] and Z[1]  */
>>    *Y[1]
>>    *Z[1]
>>    
>>    /* End of example */
>
>Matthew defines two vectors, Y and Z.  Each contains 10 
pointers.
>He also repeatedly defines a matrix named Temp:
>
>   Y[1] contains a pointer to the result of J(1,1,1)
>   Y[2] contains a pointer to the result of J(2,2,1)
>   ...
>   Y[10] contains a pointer to the result of J(10,10,1)
>
>   Z[1] contains a pointer to Temp.  At the time Z[1] was set to
>      be a pointer to Temp, Temp contained the result of 
J(1,1,1)
>   Z[2] contains a pointer to Temp.  At the time Z[2] was set to
>      be a pointer to Temp, Temp contained the result of 
J(2,2,1)
>   ...
>   Z[10] contains a pointer to Temp.  At the time Z[10] was set 
to
>      be a pointer to Temp, Temp contained the result of 
J(10,10,1)
>
>When Matthew displays what Y[1] is pointing to, he sees what
>J(1,1,1) returned, which is a 1x1 matrix containing the value 
'1'.
>
>When Matthew displays what Z[1] is pointing to, however, he 
sees
>what J(10,10,1) returned, which is a 10x10 matrix of 1s.
>
>I suspect that Matthew expected *Y[1] through *Y[10] to 
contain
>exactly what they do, but that he is surprised by what *Z[1]
>through *Z[10] contain.  *Y[1] through *Y[1] contain the result
>of J(1,1,1) through J(10,10,1).  *Z[1] through *Z[10] all 
contain
>J(10,10,1)!
>
>The loop
>
>   for (k=1;k<=10;k++) {
>      Y[k]=&J(k,k,1)
>   }
>
>and the loop
>
>   for (k=1;k<=10;k++) {
>      Temp=J(k,k,1)
>      Z[k]=&Temp
>   }
>
>certainly look similar.  Why aren't *Y[1] and *Z[1] the same?  
And,
>why are *Z[1] through *Z[10] all J(10,10,1)?
>
>In the first loop, Matthew coded
>
>   Y[k] = &J(k,k,1)
>
>For each value of k, Mata created an unnamed temporary 
matrix
>containing the result of J(k,k,1).  The '&' in front of J()
>told Mata to take the memory address of that unnamed 
temporary
>matrix and to store that memory address in Y[k].  Each of the
>10 times through the first loop, a brand new unnamed 
temporary
>matrix was created, and the address of that brand new matrix
>was placed in Y[k].
>
>I said that the matrices were temporary, but because Matthew
>obtained a pointer (memory address) for each one and stored
>those pointers in vector Y, those matrices will stick around
>until Y is dropped or memory is cleared.  Each of them is
>stored separately, and thus *Y[1] is the result of J(1,1,1),
>*Y[2] is the result of J(2,2,1), and so on.
>
>In the second loop, Matthew coded
>
>   Temp=J(k,k,1)
>   Z[k]=&Temp
>
>For each value of k, Mata stored the result of J(k,k,1) in a
>matrix named Temp.  Then, the '&' in front of Temp told Mata
>to take the memory address of Temp and to store that 
memory
>address in Z[k].
>
>Note that Temp didn't go away in each iteration of the loop;
>it simply got redefined each time.  It didn't go anywhere, and
>thus, its memory address didn't change.  So, Z[1] contains the
>memory address of Temp, Z[2] contains the memory address 
of
>Temp, ..., and Z[10] contains the memory address of Temp.
>Z[1], Z[2], ..., Z[10] all contain exactly the same pointer:
>the memory address of Temp!
>
>What did Temp have in it at the end of the loop?  J(10,10,1)
>Because Z[1] through Z[10] all contain exactly the same 
pointer
>(the memory address of Temp), whether we look at *Z[1] or 
*Z[10]
>doesn't matter.  They will all show you what Temp has in it now
>after the loop:  J(10,10,1)
>
>
>Just in case all the loops, vectors, and so forth above
>have confused the issue at all, let me illustrate the
>same thing with two simpler examples.
>
>First, I'm going to define a pointer to the result of (2+2),
>then I'm going to define a pointer to the result of (3+3), and
>then I'm going to see what is stored in the memory addresses
>pointed to by each of those pointers:
>
>   . mata
>   ------------------------------ mata (type end to exit) -----
>   : PtrToResOne = &(2+2)
>
>   : PtrToResTwo = &(3+3)
>
>   : *PtrToResOne
>     4
>
>   : *PtrToResTwo
>     6
>
>   : end
>
>
>That should be exactly what you expected to see.
>
>Now, rather than directly obtaining pointers to the results of
>(2+2) and (3+3), lets first store the results in X and get
>pointers to X instead:
>
>   . mata
>   ------------------------------ mata (type end to exit) -----
>   : X = (2+2)
>
>   : PtrToX = &X
>
>   : X = (3+3)
>
>   : AnotherPtrToX = &X
>
>   : X
>     6
>
>   : *PtrToX
>     6
>
>   : *AnotherPtrToX
>     6
>
>   : end
>
>
>Both PtrToX and AnotherPtrToX contain pointers to X.  So, it
>shouldn't surprise us that no matter what we store in X in the
>future, both PtrToX and AnotherPtrToX will point to whatever it
>is that we store.  That's the whole point of pointers.
>
>This is exactly what happened with Z[1] through Z[10] in 
Matthew's
>second loop.
>
>
>--Alan
>ariley@stata.com
>
>*
>*   For searches and help try:
>*   http://www.stata.com/help.cgi?search
>*   http://www.stata.com/support/statalist/faq
>*   http://www.ats.ucla.edu/stat/stata/
*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/