[R] x*x*x*... vs x^n

davidr@rhotrading.com davidr at rhotrading.com
Wed Jun 29 16:49:31 CEST 2005


In general, the "Russian peasant algorithm", which requires only O(log
n) multiplications, is very good. Section 4.6.3 of Knuth's The Art of
Computer Programming. Volume 2: Seminumerical Algorithms has an in depth
discussion.
I have had to use this in the past, when computers were slower and
compilers were not so clever. It is also better when x is not just a
real number, say complex or matrix (as has been mentioned.)
In many cases though, one needs many powers sequentially, and then it
may be more efficient to just multiply the previous power by x and use
the power, etc. (unless you have a parallel computer.)
So 

pows <- x^(1:1000) 
# use pows in calculations

could be sped up by employing a faster algorithm, but probably a loop
will be faster:

pows <- 1
for(i in 1:1000) {
  pows <- pows * x
  # use this power
}

David L. Reiner, Ph.D.
Rho Trading
 

> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch [mailto:r-help-
> bounces at stat.math.ethz.ch] On Behalf Of Robin Hankin
> Sent: Wednesday, June 29, 2005 9:13 AM
> To: Duncan Murdoch
> Cc: r-help; Robin Hankin
> Subject: Re: [R] x*x*x*... vs x^n
> 
> 
> On Jun 29, 2005, at 02:47 pm, Duncan Murdoch wrote:
> 
> > On 6/29/2005 9:31 AM, Robin Hankin wrote:
> >
> >> Hi  Duncan
> >>
> >>
> >> library(gsl)
> >> system.time(ignore <- pow_int(a,8))
> >> [1] 1.07 1.11 3.08 0.00 0.00
> >>
> >> <why the slow execution time?>
> >>
> >
> > Shouldn't you ask the gsl maintainer that? :-)
> >
> 
> well  I did ask myself, but  in this case the gsl maintainer
> told me to ask the gsl co-author, who
> is no doubt much better informed in these matters ;-)
> 
> >>
> >> (Of course, I'm not suggesting that other programming tasks be
> >> suspended!  All I'm pointing
> >> out is that there may exist a user to whom fast integer powers are
> >> very very important)
> >>
> >
> > Then that user should submit the patch, not you.  But whoever does
it
> > should include an argument to convince an R core member that the
> > change
> > is worth looking at, and do it well enough that the patch is
accepted.
> > Changes to the way R does arithmetic affect everyone, so they had
> > better
> > be done right, and checking them takes time.
> >
> 
> yes, that's a fair point.
> But including a native R command pow.int(), say, wouldn't  affect
> anyone, would it?
> One could even use the (tested) GSL code, as it is GPL'ed.
> 
> This would just be a new function that users could use at their
> discretion, and no
> existing code would break.
> 
> I assume that such a function would  not suffer whatever performance
> disadvantage
> that the GSL package approach had, so it may well be quite a
> significant improvement
> over the method used by R_pow_di() in arithmetic.c especially for
> large n.
> 
> 
> best wishes
> 
> rksh
> 
> > Duncan Murdoch
> >
> > ______________________________________________
> > R-help at stat.math.ethz.ch mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide! http://www.R-project.org/posting-
> > guide.html
> >
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-
> guide.html




More information about the R-help mailing list