modifying large R objects in place

I have a C function, which performs a transformation
of a large integer matrix. The matrix may be of size 1.6GB,
so I can have only one copy in RAM and have to modify it
in place. This is possible using .Call and works fine. For
debugging, I need two copies of a smaller matrix and modify only
one of them. This may also be done, for example, by
  A <- some integer matrix
  B <- A + as.integer(0)
  .Call("transform", A)
Then, B is still the original and A is the transformed one.
Up to now, I do not have any real problem with this, but there
are things, which could help me. Namely the following ones:

1. Is there a general way to duplicate an R object on the
   level of R language? I would use it instead of
      B <- A + as.integer(0)
   Looking for the keywords "duplicate" and "copy" in help pages
   and R-exts did not find what I need.

2. Is there a way how a C function can verify that its
   argument does not share data with other objects?
   I included a test of NAMED(A) in my C function, but I do
   not know, whether it may be used for the purpose which I need.
   I did not find a way how to generate an R object with
   NAMED equal 0. Is this possible on the level of R language?
   If yes, I will require that the argument passed to my
   function is generated in this way and the function
   will refuse to modify it if it is not so.
   The documentation suggests to call "duplicate", if the
   data are shared, but I cannot afford this due to memory.
   So, I can only stop if NAMED is not 0.

I appreciate, if anybody could give me advice on the above things.

Petr Savicky.

For the most part, doing anything to an R object result in it's
duplication. You generally have to do a lot of work to NOT copy an R
object.

  a <- matrix(as.integer(1),nrow=14100,ncol=14100) # 774m
  a[1,1] <- 0 # 3.0g
  gc() # 1.5g

On Wed, Sep 26, 2007 at 10:52:28AM -0700, Byron Ellis wrote:
  

For the most part, doing anything to an R object result in it's
duplication. You generally have to do a lot of work to NOT copy an R
object.
    

Thank you for your response. Unfortunately, you are right. For example,
the allocated memory determined by top command on Linux may change during
a session as follows:
  a <- matrix(as.integer(1),nrow=14100,ncol=14100) # 774m
  a[1,1] <- 0 # 3.0g
  gc() # 1.5g

In the current applicatin, I modify the matrix only using my own C code
and only read it on R level. So, the above is not a big problem for me
(at least not now).

However, there is a related thing, which could be a bug. The following
code determines the value of NAMED field in SEXP header of an object:

  SEXP getnamed(SEXP a)
  {
      SEXP out;
      PROTECT(out = allocVector(INTSXP, 1));
      INTEGER(out)[0] = NAMED(a);
      UNPROTECT(1);
      return(out);
  }

Now, consider the following session

  u <- matrix(as.integer(1),nrow=5,ncol=3) + as.integer(0)
  .Call("getnamed",u) # 1 (OK)

  length(u)
  .Call("getnamed",u) # 1 (OK)

  dim(u)
  .Call("getnamed",u) # 1 (OK)

  nrow(u)
  .Call("getnamed",u) # 2 (why?)

  u <- matrix(as.integer(1),nrow=5,ncol=3) + as.integer(0)
  .Call("getnamed",u) # 1 (OK)
  ncol(u)
  .Call("getnamed",u) # 2 (so, ncol does the same)

Is this a bug?
  

Petr Savicky.

On Wed, Sep 26, 2007 at 10:52:28AM -0700, Byron Ellis wrote:

For the most part, doing anything to an R object result in it's
duplication. You generally have to do a lot of work to NOT copy an R
object.

Thank you for your response. Unfortunately, you are right. For example,
the allocated memory determined by top command on Linux may change during
a session as follows:
 a <- matrix(as.integer(1),nrow=14100,ncol=14100) # 774m
 a[1,1] <- 0 # 3.0g
 gc() # 1.5g

In the current applicatin, I modify the matrix only using my own C code
and only read it on R level. So, the above is not a big problem for me
(at least not now).

However, there is a related thing, which could be a bug. The following
code determines the value of NAMED field in SEXP header of an object:

 SEXP getnamed(SEXP a)
 {
     SEXP out;
     PROTECT(out = allocVector(INTSXP, 1));
     INTEGER(out)[0] = NAMED(a);
     UNPROTECT(1);
     return(out);
 }

Now, consider the following session

 u <- matrix(as.integer(1),nrow=5,ncol=3) + as.integer(0)
 .Call("getnamed",u) # 1 (OK)

 length(u)
 .Call("getnamed",u) # 1 (OK)

 dim(u)
 .Call("getnamed",u) # 1 (OK)

 nrow(u)
 .Call("getnamed",u) # 2 (why?)

 u <- matrix(as.integer(1),nrow=5,ncol=3) + as.integer(0)
 .Call("getnamed",u) # 1 (OK)
 ncol(u)
 .Call("getnamed",u) # 2 (so, ncol does the same)

Is this a bug?

Petr Savicky.

2) I expected NAMED on 'a' to be incremented by nrow(a): here is my 
understanding.

When you called nrow(a) you created another reference to 'a' in the 
evaluation frame of nrow.  (At a finer level you first created a promise 
to 'a' and then dim(x) evaluated that promise, which did SET_NAMED(<SEXP>) 
= 2.)  So NAMED(a) was correctly bumped to 2, and it is never reduced.

More generally, any argument to a closure that actually gets used will 
get NAMED set to 2.

Thank you very much for all the explanations. In particular for pointing
out that nrow is not a .Primitive unlike dim, which is the
reason for the difference in their behavior. (I rised the question
of possible bug due to this difference, not just being unsatisfied
with nrow). Also, thanks for:

On Thu, Sep 27, 2007 at 05:59:05PM +0100, Prof Brian Ripley wrote:
[...]
  

2) I expected NAMED on 'a' to be incremented by nrow(a): here is my 
understanding.

When you called nrow(a) you created another reference to 'a' in the 
evaluation frame of nrow.  (At a finer level you first created a promise 
to 'a' and then dim(x) evaluated that promise, which did SET_NAMED(<SEXP>) 
= 2.)  So NAMED(a) was correctly bumped to 2, and it is never reduced.

More generally, any argument to a closure that actually gets used will 
get NAMED set to 2.
    

[...]

This explains a lot.

I appreciate also the patch to matrix by Henrik Bengtsson, which saved
me time formulating a further question just about this.

I do not know, whether there is a reason to keep nrow, ncol not .Primitive,
but if there is such, the problem may be solved by rewriting
them as follows:

nrow <- function(...) dim(...)[1]
ncol <- function(...) dim(...)[2]

At least in my environment, the new versions preserved NAMED == 1.
  

It has a side effect that this unifies the error messages generated
by too many arguments to nrow(x) and dim(x). Currently
  a <- matrix(1:6,nrow=2)
  nrow(a,a) # Error in nrow(a, a) : unused argument(s) (1:6)
  dim(a,a) # Error: 2 arguments passed to 'dim' which requires 1

May be, also other solutions exist.

Petr Savicky.

nrow <- function(...) dim(...)[1]
ncol <- function(...) dim(...)[2]

At least in my environment, the new versions preserved NAMED == 1.
 

Yes, but changing the formal arguments is a bit messy, is it not?

Presumably, nrow <- function(x) eval.parent(substitute(dim(x)[1])) works 
too, but if the gain is important enough to warrant that sort of 
programming, you might as well make nrow a .Primitive.

Longer-term, I still have some hope for better reference counting, but 
the semantics of environments make it really ugly -- an environment can 
contain an object that contains the environment, a simple example being 

f <- function()
   g <- function() 0
f()

At the end of f(), we should decide whether to destroy f's evaluation 
environment. In the present example, what we need to be able to see is 
that this would remove all refences to g and that the reference from g 
to f can therefore be ignored.  Complete logic for sorting this out is 
basically equivalent to a new garbage collector, and one can suspect 
that applying the logic upon every function return is going to be 
terribly inefficient. However, partial heuristics might apply.

On Fri, Sep 28, 2007 at 12:39:30AM +0200, Peter Dalgaard wrote:
[...]

nrow <- function(...) dim(...)[1]
ncol <- function(...) dim(...)[2]

At least in my environment, the new versions preserved NAMED == 1.

Yes, but changing the formal arguments is a bit messy, is it not?

Specifically for nrow, ncol, I think not much, since almost nobody needs
to know (or even knows) that the name of the formal argument is "x".

However, there is another argument against the ... solution: it solves
the problem only in the simplest cases like nrow, ncol, but is not
usable in other, like colSums, rowSums. These functions also increase
NAMED of its argument, although their output does not contain any
reference to the original content of their arguments.

I think that a systematic solution of this problem may be helpful.
However, making these functions Internal or Primitive would
not be good in my opinion. It is advantageous that these functions
contain an R level part, which
makes the basic decisions before a call to .Internal.
If nothing else, this serves as a sort of documentation.

For my purposes, I replaced calls to "colSums" and "matrix" by the
corresponding calls to .Internal in my script. The result is that
now I can complete several runs of my calculation in a cycle instead
of restarting R after each of the runs.

This leads me to a question. Some of the tests, which I did, suggest
that gc() may not free all the memory, even if I remove all data
objects by rm() before calling gc(). Is this possible or I must have
missed something?

A possible solution to the unwanted increase of NAMED due to temporary
calculations could be to give the user the possibility
to store NAMED attribute of an object before a call to a function
and restore it after the call. To use this, the user should be
confident that no new reference to the object persists after the
function is completed.

Presumably, nrow <- function(x) eval.parent(substitute(dim(x)[1])) works
too, but if the gain is important enough to warrant that sort of
programming, you might as well make nrow a .Primitive.

You are right. This indeed works.

Longer-term, I still have some hope for better reference counting, but
the semantics of environments make it really ugly -- an environment can
contain an object that contains the environment, a simple example being

f <- function()
   g <- function() 0
f()

At the end of f(), we should decide whether to destroy f's evaluation
environment. In the present example, what we need to be able to see is
that this would remove all refences to g and that the reference from g
to f can therefore be ignored.  Complete logic for sorting this out is
basically equivalent to a new garbage collector, and one can suspect
that applying the logic upon every function return is going to be
terribly inefficient. However, partial heuristics might apply.

I have to say that I do not understand the example very much.
What is the input and output of f? Is g inside only defined or
also used?

Let me ask the following question. I assume that gc() scans the whole
memory and determines for each part of data, whether a reference
to it still exists or not. In my understanding, this is equivalent to
determine, whether NAMED of it may be dropped to zero or not.
Structures, for which this succeeds are then removed. Am I right?
If yes, is it possible during gc() to determine also cases,
when NAMED may be dropped from 2 to 1? How much would this increase
the complexity of gc()?

Thank you in advance for your kind reply.

Petr Savicky.

On Fri, Sep 28, 2007 at 12:39:30AM +0200, Peter Dalgaard wrote:

Longer-term, I still have some hope for better reference counting, but 
the semantics of environments make it really ugly -- an environment can 
contain an object that contains the environment, a simple example being 

f <- function()
   g <- function() 0
f()

At the end of f(), we should decide whether to destroy f's evaluation 
environment. In the present example, what we need to be able to see is 
that this would remove all refences to g and that the reference from g 
to f can therefore be ignored.  Complete logic for sorting this out is 
basically equivalent to a new garbage collector, and one can suspect 
that applying the logic upon every function return is going to be 
terribly inefficient. However, partial heuristics might apply.

I have to say that I do not understand the example very much.
What is the input and output of f? Is g inside only defined or
also used?

On Fri, 28 Sep 2007, Petr Savicky wrote:

This leads me to a question. Some of the tests, which I did, suggest
that gc() may not free all the memory, even if I remove all data
objects by rm() before calling gc(). Is this possible or I must have
missed something?

Not impossible but very unlikely givent he use gc gets. There are a
few internal tables that are grown but not shrunk at the moment but
that should not usually cause much total growth.  If you are ooking at
system memopry use then that is a malloc issue -- there was a thread
about this a month or so ago.

gc()

for(s in search()) lsf.str(s)
gc()

A possible solution to the unwanted increase of NAMED due to temporary
calculations could be to give the user the possibility
to store NAMED attribute of an object before a call to a function
and restore it after the call. To use this, the user should be
confident that no new reference to the object persists after the
function is completed.

This would be too dangerous for general use. Some more structured
approach may be possible. A related issue is that user-defined
assignment functions always see a NAMED of 2 and hence cannot modify
in place. We've been trying to come up with a reasonable solution to
this, so far without success but I'm moderately hopeful.

On 9/28/2007 7:45 AM, Petr Savicky wrote:
  

On Fri, Sep 28, 2007 at 12:39:30AM +0200, Peter Dalgaard wrote:
    

  ...
  

Longer-term, I still have some hope for better reference counting, but 
the semantics of environments make it really ugly -- an environment can 
contain an object that contains the environment, a simple example being 

f <- function()
   g <- function() 0
f()

At the end of f(), we should decide whether to destroy f's evaluation 
environment. In the present example, what we need to be able to see is 
that this would remove all refences to g and that the reference from g 
to f can therefore be ignored.  Complete logic for sorting this out is 
basically equivalent to a new garbage collector, and one can suspect 
that applying the logic upon every function return is going to be 
terribly inefficient. However, partial heuristics might apply.
      

I have to say that I do not understand the example very much.
What is the input and output of f? Is g inside only defined or
also used?
    

f has no input; it's output is the function g, whose environment is the 
evaluation environment of f.  g is never used, but it is returned as the 
value of f.  Thus we have the loop:

g refers to the environment.
the environment contains g.

Even though the result of f() was never saved, two things (the 
environment and g) got created and each would have non-zero reference 
count.

In a more complicated situation you might want to save the result of the 
function and then modify it.  But because of the loop above, you would 
always think there's another reference to the object, so every in-place 
modification would require a copy first.

  

Duncan Murdoch

[...] A related issue is that user-defined
assignment functions always see a NAMED of 2 and hence cannot modify
in place. We've been trying to come up with a reasonable solution to
this, so far without success but I'm moderately hopeful.

Longer-term, I still have some hope for better reference counting, but 
the semantics of environments make it really ugly -- an environment can 
contain an object that contains the environment, a simple example being 

f <- function()
   g <- function() 0
f()

f has no input; it's output is the function g, whose environment is the 
evaluation environment of f.  g is never used, but it is returned as the 
value of f.  Thus we have the loop:

g refers to the environment.
the environment contains g.

Even though the result of f() was never saved, two things (the 
environment and g) got created and each would have non-zero reference 
count.

If yes, is it possible during gc() to determine also cases,
when NAMED may be dropped from 2 to 1? How much would this increase
the complexity of gc()?

Probably not impossible but would be a fair bit of work with probably
not much gain as the NAMED values would still be high until the next
gc of the appropriate level, which will probably be a fair time as an
object being modified is likely to be older, but the interval in which
there would be a benefit is short.

On Fri, 28 Sep 2007, Luke Tierney wrote:

approach may be possible. A related issue is that user-defined
assignment functions always see a NAMED of 2 and hence cannot modify
in place. We've been trying to come up with a reasonable solution to
this, so far without success but I'm moderately hopeful.

I am not persuaded that the difference between NAMED=1/2 makes much 
difference in general use of R, and I recall Ross saying that he no longer 
believed that this was a worthwhile optimization.  It's not just 
'user-defined' replacement functions, but also all the system-defined 
closures (including all methods for the generic replacement functions 
which are primitive) that are unable to benefit from it.