Re: [Rd] Replace isnan and lgamma in Fortran subroutine in R package
On Tue, 2014-09-23 at 07:43 +0200, Berend Hasselman wrote: > On 23-09-2014, at 00:33, Wang, Zhu wrote: > > > Hello, > > > > I submitted a package which used Fortran functions isnan and lgamma. > > However, I was told that: > > > > isnan and lgamma are not Fortran 95 functions. > > > > I was asked to write 'cross-platform portable code' and so should not be > > writing GNU extensions to Fortran. > > > > See http://cran.r-project.org/web/checks/check_results_mpath.html, which > > will shortly show installation failures under Solaris. > > > > I will appreciate advice on how to replace these two functions to avoid > > failure on some platforms. > > > > I don’t know about lgamma. > > Instead of isnan you could use Lapack’s logical function disnan to test for > NaN (it’s in lapack 3.4.2; I don’t know about earlier versions). > > Another way would be to write a C function using functions provided by R to > test for NaN. > That function should be declared with F77_NAME so that it can be called from > a Fortran routine. > > I haven’t tried this so I don’t know if this would be foolproof or is the > best way to do it.. > > Berend > As described by Berend (See also section 6.6 of Writing R Extensions) you can define Fortran-callable wrappers around ISNAN and lgammafn, both of which are provided with the C interface to R (lgammafn is declared in Rmath.h). Then declare these wrappers as external subroutines in your Fortran code and call them in place of isnan and lgamma. Martyn __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Replace isnan and lgamma in Fortran subroutine in R package
Try this patch. Martyn On Mon, 2014-09-22 at 22:33 +, Wang, Zhu wrote: > Hello, > > I submitted a package which used Fortran functions isnan and lgamma. However, > I was told that: > > isnan and lgamma are not Fortran 95 functions. > > I was asked to write 'cross-platform portable code' and so should not be > writing GNU extensions to Fortran. > > See http://cran.r-project.org/web/checks/check_results_mpath.html, which will > shortly show installation failures under Solaris. > > I will appreciate advice on how to replace these two functions to avoid > failure on some platforms. > > Thanks in advance, > > Zhu Wang > > > > [[alternative HTML version deleted]] > > __ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel --- This message and its attachments are strictly confidential. If you are not the intended recipient of this message, please immediately notify the sender and delete it. Since its integrity cannot be guaranteed, its content cannot involve the sender's responsibility. Any misuse, any disclosure or publication of its content, either whole or partial, is prohibited, exception made of formally approved use --- diff -uNr mpath/src/cfuns.c mpath-fixed/src/cfuns.c --- mpath/src/cfuns.c 1970-01-01 01:00:00.0 +0100 +++ mpath-fixed/src/cfuns.c 2014-09-23 09:56:24.303347209 +0200 @@ -0,0 +1,5 @@ +#include +#include + +int F77_SUB(cisnan)(double *x) { return ISNAN(*x); } +double F77_SUB(rlgamma)(double *x) { return lgammafn(*x); } diff -uNr mpath/src/midloop.f mpath-fixed/src/midloop.f --- mpath/src/midloop.f 2014-09-22 19:32:25.0 +0200 +++ mpath-fixed/src/midloop.f 2014-09-23 09:56:55.455497773 +0200 @@ -255,6 +255,8 @@ integer n, family, i double precision dev, y(n), mu(n), theta, weights(n),tmp + integer cisnan + external cisnan dev = 0.0D0 do 10 i=1, n @@ -282,7 +284,7 @@ dev=dev+2*(weights(i)*(y(i)*dlog(max(1.0D0,y(i))/mu(i))- + (y(i)+theta)*dlog((y(i)+theta)/(mu(i) + theta endif - if(isnan(dev)) then + if(cisnan(dev).NE.0) then call intpr("i=", -1, i, 1) call dblepr("y(i)=", -1, y(i), 1) call dblepr("mu(i)=", -1, mu(i), 1) @@ -300,6 +302,8 @@ integer n, family, i, y0 double precision ll, y(n), mu(n), theta, w(n) + double precision rlgamma + external rlgamma ll = 0 do 10 i=1,n @@ -309,8 +313,8 @@ else y0=0 endif - ll=ll + w(i) * (lgamma(theta+y(i))-lgamma(theta)- - + lgamma(y(i)+1)+ + ll=ll + w(i) * (rlgamma(theta+y(i))-rlgamma(theta)- + + rlgamma(y(i)+1)+ + theta*log(theta) + y(i)*log(mu(i)+y0) - (theta + y(i))* + log(theta + mu(i))) else if(family .EQ. 1)then !gaussian @@ -320,7 +324,7 @@ ll=ll + w(i)*(y(i)*log(mu(i)/(1-mu(i)))+log(1-mu(i))) endif else if(family .EQ. 3) then !poisson -ll=ll + w(i)*(-mu(i) + y(i)*log(mu(i))-lgamma(y(i)+1)) +ll=ll + w(i)*(-mu(i) + y(i)*log(mu(i))-rlgamma(y(i)+1)) endif 10 continue return __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Replace isnan and lgamma in Fortran subroutine in R package
That patch works flawlessly. Thanks Martyn. Zhu -Original Message- From: Martyn Plummer [mailto:plumm...@iarc.fr] Sent: Tuesday, September 23, 2014 8:39 AM To: Wang, Zhu Cc: r-devel@r-project.org Subject: Re: [Rd] Replace isnan and lgamma in Fortran subroutine in R package Try this patch. Martyn On Mon, 2014-09-22 at 22:33 +, Wang, Zhu wrote: > Hello, > > I submitted a package which used Fortran functions isnan and lgamma. However, > I was told that: > > isnan and lgamma are not Fortran 95 functions. > > I was asked to write 'cross-platform portable code' and so should not be > writing GNU extensions to Fortran. > > See http://cran.r-project.org/web/checks/check_results_mpath.html, which will > shortly show installation failures under Solaris. > > I will appreciate advice on how to replace these two functions to avoid > failure on some platforms. > > Thanks in advance, > > Zhu Wang > > > > [[alternative HTML version deleted]] > > __ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel --- This message and its attachments are strictly confidenti...{{dropped:8}} __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Bug in new behaviour for all.equal and environments?
> Duncan Murdoch > on Sun, 21 Sep 2014 15:57:00 -0400 writes: > On 21/09/2014, 1:38 PM, Rui Barradas wrote: >> Hello, >> >> In R 3.1.1 on Windows 7 it's ok. >> >> > all.equal(baseenv(), baseenv()) [1] TRUE > >> sessionInfo() R version 3.1.1 (2014-07-10) Platform: >> x86_64-w64-mingw32/x64 (64-bit) > I'm not sure if that's really "ok", since it returns TRUE > for any pair of environments, e.g. > all.equal(baseenv(), emptyenv()) > The current R-devel behaviour is temporary and certain to > change. > Duncan Murdoch Indeed. A bit more than an hour ago, I have committed several changes to R-devel's behavior in dealing with objects such as (rp <- getClass("refClass")@prototype) (ner <- new("envRefClass")) both of which would not even correctly *print* in current R. They now do in R-devel, also work with str() there, and all.equal(baseenv(), baseenv()) and many more cases of all.equal() with environments now work there as well. I still expect that creative useRs will be able to construct cases where the new all.equal() methods will fail, possibly even end in an infinite loop (which is caught and leads to an error). Please use R-devel svn rev >= 8 for such experiments. Thanking you for testing ..., Martin Maechler ETH Zurich >> [] >> Rui Barradas >> >> Em 21-09-2014 18:06, Kevin Ushey escreveu: >>> Hi R-devel, >>> >>> The following code: >>> >>> all.equal(baseenv(), baseenv()) >>> >>> gives the error when run in a clean R session with >>> latest R-devel (r66650): >>> >>> kevin:~$ R --vanilla --slave -e "all.equal(baseenv(), >>> baseenv())" Error in all.equal.envRefClass(target[[i]], >>> current[[i]], check.attributes = check.attributes, : >>> attempt to apply non-function Calls: all.equal >>> ... all.equal.list -> all.equal -> all.equal.envRefClass >>> Execution halted >>> >>> Although I don't know if it's helpful -- it appears that >>> packages that include some S4 machinery will effect the >>> outcome of this error, e.g. if we load 'Biobase' first: >>> >>> kevin:~$ R --vanilla --slave -e >>> "suppressPackageStartupMessages(library(Biobase)); >>> all.equal(baseenv(), baseenv())" Error in >>> target$getClass() : object '.refClassDef' not found >>> Calls: all.equal ... all.equal.list -> all.equal -> >>> all.equal.envRefClass -> Execution halted >>> >>> We were bumping into an error with this in Rcpp -- we >>> used all.equal (through RUnit) to confirm that >>> baseenv(), and a function returning the base >>> environment, would return TRUE. >>> >>> For completeness: >>> >>> kevin:~$ R --vanilla --slave -e "sessionInfo()" R Under >>> development (unstable) (2014-09-20 r66650) Platform: >>> x86_64-apple-darwin13.3.0 (64-bit) >>> >>> locale: [1] >>> en_CA.UTF-8/en_CA.UTF-8/en_CA.UTF-8/C/en_CA.UTF-8/en_CA.UTF-8 >>> >>> attached base packages: [1] stats graphics grDevices >>> utils datasets methods base >>> >>> Thanks, Kevin __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
[Rd] Patch for R to fix some buffer overruns and add a missing PROTECT().
This patch is against current svn and contains three classes of fix: - Ensure the result is properly terminated after calls to strncpy() - Replace calls of sprintf() with snprintf() - Added a PROTECT() call in do_while which could cause memory errors if evaluating the condition results in a warning. Thanks, Karl __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Patch for R to fix some buffer overruns and add a missing PROTECT().
On 23/09/2014 3:20 PM, Karl Millar wrote: This patch is against current svn and contains three classes of fix: - Ensure the result is properly terminated after calls to strncpy() - Replace calls of sprintf() with snprintf() - Added a PROTECT() call in do_while which could cause memory errors if evaluating the condition results in a warning. Nothing was attached. Generally fixes like this are best sent to bugs.r-project.org, and they receive highest priority if accompanied by code demonstrating why they are needed, i.e. crashes or incorrect results in current R. Those will likely be incorporated as regression tests. Duncan Murdoch __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Patch for R to fix some buffer overruns and add a missing PROTECT().
Bug submitted. Thanks. On Tue, Sep 23, 2014 at 12:42 PM, Duncan Murdoch wrote: > On 23/09/2014 3:20 PM, Karl Millar wrote: >> >> This patch is against current svn and contains three classes of fix: >> - Ensure the result is properly terminated after calls to strncpy() >> - Replace calls of sprintf() with snprintf() >> - Added a PROTECT() call in do_while which could cause memory >> errors if evaluating the condition results in a warning. > > > Nothing was attached. > > Generally fixes like this are best sent to bugs.r-project.org, and they > receive highest priority if accompanied by code demonstrating why they are > needed, i.e. crashes or incorrect results in current R. Those will likely > be incorporated as regression tests. > > Duncan Murdoch __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
Re: [Rd] Pointer ownership with GECreateDevDesc/GEDestroyDevDesc
This seems like the best solution (free and NULL the dev pointer in dev_Close(), BEFORE GEdestroyDevDesc() has a go). I think this should be safe because dev_Close() and GEdestroyDevDesc() are only called in one place (so your free() should always happen before GEdestroyDevDesc() and your free() shouldn't happen in any other, potentially inappropriate, context). Paul On 09/23/14 03:15, Ritch Melton wrote: I perused the source for RStudio and noticed that they had a similar problem. I didn't think that the GEDevDesc structure was exported, but apparently this is a worthwhile workaround to the defect. void GD_Close(pDevDesc dev) { ...elided... // explicitly free and then null out the dev pointer of the GEDevDesc // This is to avoid incompatabilities between the heap we are compiled with // and the heap R is compiled with (we observed this to a problem with // 64-bit R) std::free(s_pGEDevDesc->dev); s_pGEDevDesc->dev = NULL; // set GDDevDesc to NULL so we don't reference it again s_pGEDevDesc = NULL; } } I will On Fri, Sep 19, 2014 at 10:37 AM, Ritch Melton wrote: According to the "R Internals" document, for a custom device, I should create a pDevDesc structure that gets passed to GECreateDevDesc. ..elided... pDevDesc dev; /* Allocate and initialize the device driver data */ if (!(dev = (pDevDesc) calloc(1, sizeof(DevDesc return 0; /* or error() */ /* set up device driver or free ’dev’ and error() */ gdd = GEcreateDevDesc(dev); GEaddDevice2(gdd, "dev_name"); ...elided... which indicates to me that the calling code owns the pDevDesc structure and should be responsible for freeing it. However, in GEdestroyDevDesc, that particular code calls free(dd->dev), freeing the pointer. In my case, I have a device implemented and created in a language that is using a different allocator. The pDevDesc I pass in is not owned by R and causes a crash on free with "dev.off()". I thought I could mitigate the issue, by making a call to the allocator R is using via an exported R function, such as R_alloc, but I did not see such a function available to me. I would like to know if I should be creating the pDevDesc via an imported alloc routine somehow? or is there a way to prevent this code from freeing the structure I pass in. Blue Skies, Ritch [[alternative HTML version deleted]] __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel -- Dr Paul Murrell Department of Statistics The University of Auckland Private Bag 92019 Auckland New Zealand 64 9 3737599 x85392 p...@stat.auckland.ac.nz http://www.stat.auckland.ac.nz/~paul/ __ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel