[R] x*x*x*... vs x^n
davidr@rhotrading.com
davidr at rhotrading.com
Wed Jun 29 17:23:54 CEST 2005
Looking at the code for gsl_pow_int, I see they do use that method.
David L. Reiner
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch [mailto:r-help-
> bounces at stat.math.ethz.ch] On Behalf Of David Reiner
> <davidr at rhotrading.com>
> Sent: Wednesday, June 29, 2005 9:50 AM
> To: r-help
> Subject: Re: [R] x*x*x*... vs x^n
>
> 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
>
> ______________________________________________
> 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