Re: [Rd] [R] aligned memory allocation in C
Thanks for your remarks. According to http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/howto-compile.html, I need 16 bit aligned memory when using fill_array64. So I suppose I need 8 bit aligned memory. I will test what you advise me and will come back to R-devel list after. Thanks again Christophe 2008/8/13 Prof Brian Ripley [EMAIL PROTECTED] On Tue, 12 Aug 2008, Jeffrey Horner wrote: Christophe Dutang1 wrote: Hi, I'm currently R porting SF Mersenne Twister algorithm of Matsumoto and Saito. To get the full power of their code, I want to use their fonction fill_array32 which need aligned memory. That is to say I need to use the C function memalign on windows, posix_memalign on linux and classic malloc on Mac OS. In 'writing R extenstion', they recommand to use R_alloc function to allocate memory in C. Does R_alloc return a pointer to aligned memory? if not how can I do this? probably no, because R crashes when I succesively R_alloc and fill_array32 (cf below) on my macbook with R 2.7.1. You can still do this. Just take the address returned from R_alloc and test for alignment. If it's not, then just use an aligned address beyond the one returned. We haven't been told what the desired alignment is (and those functions need to be told). On 32-bit Mac OS X, R_alloc is definitely aligned on 4-byte boundaries (on 64-bit OSes it is usually 8-byte aligned). (But then the question is, which direction beyond the one returned? How does one test for that?) Addresses always go upwards. So if you want 64-byte alignment you need to allocate a block at least 64 bytes longer than required, and go up to the nearest multiple of 64. BTW, this is clearly an R-devel question -- see the posting guide. Jeff Thanks in advance Kind regards Christophe PS : http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/howto-compile.htmlhttp://www.math.sci.hiroshima-u.ac.jp/%7Em-mat/MT/SFMT/howto-compile.htmlprovides an example of memalign. PPS : mac os report [removed] -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/http://www.stats.ox.ac.uk/%7Eripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 [[alternative HTML version deleted]] __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] extending the derivs table/fools rushing in
While you are at it could you add { to the
table so that this works:
# this is ok
f - function(x) x*x
D(body(f), x)
x + x
# but not g which is same as f
# except it has { ... } surrounding its body
g - function(x) { x*x }
D(body(g), x)
Error in D(body(g), x) : Function '`{`' is not in the derivatives table
2008/8/14 Ben Bolker [EMAIL PROTECTED]:
I added plogis to the derivative table in the
development version of R; the patch against yesterday's
R-devel src/deriv/main.c is available at
http://www.zoology.ufl.edu/bolker/deriv_patch.txt .
I pretty much followed the framework of the other symbols;
here was my incantation
- } else if (CAR(expr) == PlogisSymbol) {
- ans = simplify(TimesSymbol,
- PP_S(TimesSymbol,
- PP_S2(ExpSymbol,
- PP_S2(MinusSymbol,CADR(expr))),
- PP_S(PowerSymbol,
-PP_S(PlusSymbol,
- Constant(1.),
- PP_S2(ExpSymbol,
-PP_S2(MinusSymbol,CADR(expr,
-Constant(-2.))),
- PP(D(CADR(expr),var)));
- UNPROTECT(8);
It seems to work:
D(quote(plogis(a)),a)
exp(-a) * (1 + exp(-a))^-2
D(quote(plogis(a+b*x)),x)
exp(-(a + b * x)) * (1 + exp(-(a + b * x)))^-2 * b
Any thoughts? I'm sure there's a cleverer way to do this ...
Ben Bolker
PS I get a minor build error at the end of the compilation
(on Ubuntu 8.10):
make[2]: *** No rule to make target `VR.ts', needed by
`stamp-recommended'. Stop.
make[2]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make[1]: *** [recommended-packages] Error 2
make[1]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make: *** [stamp-recommended] Error 2
... don't know if that is important or not.
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] extending the derivs table/fools rushing in
http://www.zoology.ufl.edu/bolker/deriv_patch2.txt
has this change added as well.
However, I'm not as confident that this is the right thing
to do? Should curly brackets even be appearing in mathematical
expressions?
Ben
Gabor Grothendieck wrote:
While you are at it could you add { to the
table so that this works:
# this is ok
f - function(x) x*x
D(body(f), x)
x + x
# but not g which is same as f
# except it has { ... } surrounding its body
g - function(x) { x*x }
D(body(g), x)
Error in D(body(g), x) : Function '`{`' is not in the derivatives table
2008/8/14 Ben Bolker [EMAIL PROTECTED]:
I added plogis to the derivative table in the
development version of R; the patch against yesterday's
R-devel src/deriv/main.c is available at
http://www.zoology.ufl.edu/bolker/deriv_patch.txt .
I pretty much followed the framework of the other symbols;
here was my incantation
- } else if (CAR(expr) == PlogisSymbol) {
- ans = simplify(TimesSymbol,
- PP_S(TimesSymbol,
- PP_S2(ExpSymbol,
- PP_S2(MinusSymbol,CADR(expr))),
- PP_S(PowerSymbol,
-PP_S(PlusSymbol,
- Constant(1.),
- PP_S2(ExpSymbol,
-PP_S2(MinusSymbol,CADR(expr,
-Constant(-2.))),
- PP(D(CADR(expr),var)));
- UNPROTECT(8);
It seems to work:
D(quote(plogis(a)),a)
exp(-a) * (1 + exp(-a))^-2
D(quote(plogis(a+b*x)),x)
exp(-(a + b * x)) * (1 + exp(-(a + b * x)))^-2 * b
Any thoughts? I'm sure there's a cleverer way to do this ...
Ben Bolker
PS I get a minor build error at the end of the compilation
(on Ubuntu 8.10):
make[2]: *** No rule to make target `VR.ts', needed by
`stamp-recommended'. Stop.
make[2]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make[1]: *** [recommended-packages] Error 2
make[1]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make: *** [stamp-recommended] Error 2
... don't know if that is important or not.
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
signature.asc
Description: OpenPGP digital signature
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] extending the derivs table/fools rushing in
The derivative of plogis is surely dlogis. (And yes, there is a good
reason why we have such a function: take a look at its C code.)
That means we would need an entry for dlogis too, I guess. I am not
convinced that there is a real need for these (and where does this stop?)
What would be much more useful is to make this user-extensible (as Bill
Venables pointed out a decade ago). [pd]norm were added in 2002 to
support MASS, the ability to do all of MASS in R being a goal at the time.
On Thu, 14 Aug 2008, Ben Bolker wrote:
I added plogis to the derivative table in the
development version of R; the patch against yesterday's
R-devel src/deriv/main.c is available at
http://www.zoology.ufl.edu/bolker/deriv_patch.txt .
I pretty much followed the framework of the other symbols;
here was my incantation
- } else if (CAR(expr) == PlogisSymbol) {
- ans = simplify(TimesSymbol,
- PP_S(TimesSymbol,
- PP_S2(ExpSymbol,
- PP_S2(MinusSymbol,CADR(expr))),
- PP_S(PowerSymbol,
-PP_S(PlusSymbol,
- Constant(1.),
- PP_S2(ExpSymbol,
-PP_S2(MinusSymbol,CADR(expr,
-Constant(-2.))),
- PP(D(CADR(expr),var)));
- UNPROTECT(8);
It seems to work:
D(quote(plogis(a)),a)
exp(-a) * (1 + exp(-a))^-2
D(quote(plogis(a+b*x)),x)
exp(-(a + b * x)) * (1 + exp(-(a + b * x)))^-2 * b
Any thoughts? I'm sure there's a cleverer way to do this ...
Ben Bolker
PS I get a minor build error at the end of the compilation
(on Ubuntu 8.10):
make[2]: *** No rule to make target `VR.ts', needed by
`stamp-recommended'. Stop.
make[2]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make[1]: *** [recommended-packages] Error 2
make[1]: Leaving directory
`/usr/local/src/R/R-devel/src/library/Recommended'
make: *** [stamp-recommended] Error 2
... don't know if that is important or not.
A local-to-you problem.
--
Brian D. Ripley, [EMAIL PROTECTED]
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UKFax: +44 1865 272595
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] extending the derivs table/fools rushing in
Prof Brian Ripley wrote: The derivative of plogis is surely dlogis. (And yes, there is a good reason why we have such a function: take a look at its C code.) Doh. That means we would need an entry for dlogis too, I guess. I am not convinced that there is a real need for these (and where does this stop?) What would be much more useful is to make this user-extensible (as Bill Venables pointed out a decade ago). [pd]norm were added in 2002 to support MASS, the ability to do all of MASS in R being a goal at the time. I agree it would be great to make this user-extensible, but it's probably a bit beyond me ... I found a web-reference of Venables saying There is a detailed example towards the end of Ch. 9 of VR on how to extend D() and make.call(), which are the work horses for deriv(), to handle new functions. The new functions handled there are dnorm() and pnorm(), but I() would be even easier, of course. http://www.math.yorku.ca/Who/Faculty/Monette/pub/s-old/0690.html ... but this is from 1997 therefore presumably MASS3? or MASS2? -- I can't find my copy of MASS3 at the moment, and don't own MASS2 ... The reason behind this is that I was trying to write a simple analytic derivative calculator for formulae of the form (e.g.) y ~ dbinom(prob=plogis(a+b*x),size=N) Obviously in this case I could just tell people to write the formula out as y ~ dbinom(prob=1/(1+exp(-(a+b*x))),size=N) ... Ben Bolker signature.asc Description: OpenPGP digital signature __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] extending the derivs table/fools rushing in
You could see if Ryacas package can do what you want: http://ryacas.googlecode.com 2008/8/14 Ben Bolker [EMAIL PROTECTED]: Prof Brian Ripley wrote: The derivative of plogis is surely dlogis. (And yes, there is a good reason why we have such a function: take a look at its C code.) Doh. That means we would need an entry for dlogis too, I guess. I am not convinced that there is a real need for these (and where does this stop?) What would be much more useful is to make this user-extensible (as Bill Venables pointed out a decade ago). [pd]norm were added in 2002 to support MASS, the ability to do all of MASS in R being a goal at the time. I agree it would be great to make this user-extensible, but it's probably a bit beyond me ... I found a web-reference of Venables saying There is a detailed example towards the end of Ch. 9 of VR on how to extend D() and make.call(), which are the work horses for deriv(), to handle new functions. The new functions handled there are dnorm() and pnorm(), but I() would be even easier, of course. http://www.math.yorku.ca/Who/Faculty/Monette/pub/s-old/0690.html ... but this is from 1997 therefore presumably MASS3? or MASS2? -- I can't find my copy of MASS3 at the moment, and don't own MASS2 ... The reason behind this is that I was trying to write a simple analytic derivative calculator for formulae of the form (e.g.) y ~ dbinom(prob=plogis(a+b*x),size=N) Obviously in this case I could just tell people to write the formula out as y ~ dbinom(prob=1/(1+exp(-(a+b*x))),size=N) ... Ben Bolker __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
[Rd] Typo in proc.time.Rd
There's a small typo in proc.time.Rd:
\sQuote{system lime} should be \sQuote{system time}
Thanks,
Stephen
--
Rochester, Minn. USA
__
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] [R] RNG Cycle and Duplication (PR#12538)
Shengqiao Li wrote: Hello all, I am generating large samples of random numbers. The RNG help page says: All the supplied uniform generators return 32-bit integer values that are converted to doubles, so they take at most 2^32 distinct values and long runs will return duplicated values. But I find that the cycles are not the same as the 32-bit integer. This is no bug. Cycle length and the set of possible values are different concepts. And please don't cross-post. Duncan Murdoch My test indicated that the cycles for Knuth's methods were 2^30 while Wichmann-Hill's cycle was larger than 2^32! No numbers were duplicated in 10M numbers generated by runif using Wichmann-Hill. The other three methods had cycle length of 2^32. So, anybody can explain this? And any improvement to the implementation can be made to increase the cycle length like the Wichmann-Hill method? Shengqiao Li Research Associate The Department of Statistics West Virginia University Morgantown, WV 26506-6330 __ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] [R] RNG Cycle and Duplication (PR#12540)
This message is in MIME format. The first part should be readable text, while the remaining parts are likely unreadable without MIME-aware tools. ---559023410-851401618-1218751024=:15885 Content-Type: TEXT/PLAIN; charset=ISO-8859-1; format=flowed Content-Transfer-Encoding: QUOTED-PRINTABLE I didn't describe the problem clearly. It's about the number of distinct=20 values. So just ignore cycle issue. My tests were: RNGkind(kind=3DKnuth-TAOCP); sum(duplicated(runif(1e7))); #return 46552 RNGkind(kind=3DKnuth-TAOCP-2002); sum(duplicated(runif(1e7))); #return 46415 #These collision frequency suggested there were 2^30 distinct values by=20 birthday problem. RNGkind(kind=3DMarsaglia-Multicarry); sum(duplicated(runif(1e7))); #return 11682 RNGkind(kind=3DSuper-Duper); sum(duplicated(runif(1e7))); #return 11542 RNGkind(kind=3DMersenne-Twister); sum(duplicated(runif(1e7))); #return 11656 #These indicated there were 2^32 distinct values, which agrees with the=20 help info. RNGkind(kind=3DWichmann-Hill); sum(duplicated(runif(1e7))); #return 0 #So for this method, there should be more than 2^32 distinct values. You may not get the exact numbers, but they should be close. So how to=20 explain above problem? I need generate a large sample without any ties, it seems to me=20 Wichmann-Hill is only choice right now. =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D= =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D Shengqiao Li The Department of Statistics PO Box 6330 West Virginia University Morgantown, WV 26506-6330 =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D= =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D On Thu, 14 Aug 2008, Peter Dalgaard wrote: Shengqiao Li wrote: Hello all, =20 I am generating large samples of random numbers. The RNG help page says:= =20 All the supplied uniform generators return 32-bit integer values that a= re=20 converted to doubles, so they take at most 2^32 distinct values and long= =20 runs will return duplicated values. But I find that the cycles are not = the=20 same as the 32-bit integer. =20 My test indicated that the cycles for Knuth's methods were 2^30 while=20 Wichmann-Hill's cycle was larger than 2^32! No numbers were duplicated i= n=20 10M numbers generated by runif using Wichmann-Hill. The other three meth= ods=20 had cycle length of 2^32. =20 So, anybody can explain this? And any improvement to the implementation = can=20 be made to increase the cycle length like the Wichmann-Hill method? =20 What test? These are not simple linear congruential generators. Just beca= use=20 you get the same value twice, it doesn't mean that the sequence is repeat= ing.=20 Perhaps you should read the entire help page rather than just the note. --=20 O__ Peter Dalgaard =D8ster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 ---559023410-851401618-1218751024=:15885-- __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] [R] RNG Cycle and Duplication (PR#12540)
Remember Wichmann-Hill is a composite generator. Its composition does take more than 2^32 distinct values. You still haven't identifed a problem here. The note is to warn that runif() does repeat within a cycle, because people wrote code assuming otherwise. It would be a poor use of runif() to rely on the low-order bits, and that's standard advice in the field. For a large sample of uniforms use something like the normal inversion does, e.g. 2^(-30) * (runif(N, 0, 2^30) %% 2^30 + runif(N)) Please do leave R-bugs out of this: we already have 4 entries as a result of your misunderstandings and false claims. On Thu, 14 Aug 2008, [EMAIL PROTECTED] wrote: This message is in MIME format. The first part should be readable text, while the remaining parts are likely unreadable without MIME-aware tools. ---559023410-851401618-1218751024=:15885 Content-Type: TEXT/PLAIN; charset=ISO-8859-1; format=flowed Content-Transfer-Encoding: QUOTED-PRINTABLE I didn't describe the problem clearly. It's about the number of distinct=20 values. So just ignore cycle issue. My tests were: RNGkind(kind=3DKnuth-TAOCP); sum(duplicated(runif(1e7))); #return 46552 RNGkind(kind=3DKnuth-TAOCP-2002); sum(duplicated(runif(1e7))); #return 46415 #These collision frequency suggested there were 2^30 distinct values by=20 birthday problem. RNGkind(kind=3DMarsaglia-Multicarry); sum(duplicated(runif(1e7))); #return 11682 RNGkind(kind=3DSuper-Duper); sum(duplicated(runif(1e7))); #return 11542 RNGkind(kind=3DMersenne-Twister); sum(duplicated(runif(1e7))); #return 11656 #These indicated there were 2^32 distinct values, which agrees with the=20 help info. RNGkind(kind=3DWichmann-Hill); sum(duplicated(runif(1e7))); #return 0 #So for this method, there should be more than 2^32 distinct values. You may not get the exact numbers, but they should be close. So how to=20 explain above problem? I need generate a large sample without any ties, it seems to me=20 Wichmann-Hill is only choice right now. =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D= =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D Shengqiao Li The Department of Statistics PO Box 6330 West Virginia University Morgantown, WV 26506-6330 =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D= =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D On Thu, 14 Aug 2008, Peter Dalgaard wrote: Shengqiao Li wrote: Hello all, =20 I am generating large samples of random numbers. The RNG help page says:= =20 All the supplied uniform generators return 32-bit integer values that a= re=20 converted to doubles, so they take at most 2^32 distinct values and long= =20 runs will return duplicated values. But I find that the cycles are not = the=20 same as the 32-bit integer. =20 My test indicated that the cycles for Knuth's methods were 2^30 while=20 Wichmann-Hill's cycle was larger than 2^32! No numbers were duplicated i= n=20 10M numbers generated by runif using Wichmann-Hill. The other three meth= ods=20 had cycle length of 2^32. =20 So, anybody can explain this? And any improvement to the implementation = can=20 be made to increase the cycle length like the Wichmann-Hill method? =20 What test? These are not simple linear congruential generators. Just beca= use=20 you get the same value twice, it doesn't mean that the sequence is repeat= ing.=20 Perhaps you should read the entire help page rather than just the note. --=20 O__ Peter Dalgaard =D8ster Farimagsgade 5, Entr.B c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K (*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918 ~~ - ([EMAIL PROTECTED]) FAX: (+45) 35327907 ---559023410-851401618-1218751024=:15885-- __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel -- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UKFax: +44 1865 272595 __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
