Re: chisqsim broken on alpha - solved + next bug in fft found (PR#277)

About this list Date view Thread view Subject view Author view Other groups

Subject: Re: chisqsim broken on alpha - solved + next bug in fft found (PR#277)
albrecht.gebhardt@uni-klu.ac.at
Date: Wed 15 Sep 1999 - 23:33:09 EST


Message-Id: <199909151333.PAA17503@pubhealth.ku.dk>

On Wed, 15 Sep 1999 ripley@stats.ox.ac.uk wrote:

> On Tue, 14 Sep 1999, Albrecht Gebhardt wrote:
>
> [...]
>
> > chisqsim (called from chisq.test() crashes. First I found that the crash
> > is at
> > chisqsim.c:167
> > 167 fact[i] = x;
>
> [...]
>
> chisqsum.c has its interface in terms of long *, called by .C that supplies
> int *. I think those are different on an Alpha. Could you please
> edit long to int in that file and try again.
>
> Brian
>

Yes this helps.
The patch is:
--- ./src/appl/chisqsim.c.chisqsim-patch Wed Sep 15 09:17:56 1999
+++ ./src/appl/chisqsim.c Wed Sep 15 14:30:56 1999
@@ -11,10 +11,15 @@
 #include "S.h"
 
 static void
+#ifdef __osf__
+rcont2(int *nrow, int *ncol, int *nrowt, int *ncolt, int *ntotal,
+ double *fact, int *jwork, int *matrix)
+#else
 rcont2(long *nrow, long *ncol, long *nrowt, long *ncolt, long *ntotal,
        double *fact, long *jwork, long *matrix)
+#endif
 {
- long nlmp, j, l, m, ia, ib, ic, jc, id, ie, ii, nrowtl, iap, idp,
+ int nlmp, j, l, m, ia, ib, ic, jc, id, ie, ii, nrowtl, iap, idp,
        igp, ihp, iip, nll, nlm, nrowm, ncolm, lsm, lsp;
     double x, y, dummy, sumprb;
 
@@ -151,12 +156,18 @@
    */
 
 void
+#ifdef __osf__
+chisqsim(int *nrow, int *ncol, int *nrowt, int *ncolt, int *n,
+ int *b, double *expected, int *observed, double *fact,
+ int *jwork, double *results)
+#else
 chisqsim(long *nrow, long *ncol, long *nrowt, long *ncolt, long *n,
         long *b, double *expected, long *observed, double *fact,
         long *jwork, double *results)
+#endif
 {
     /* Local variables */
- long i, j, iter;
+ int i, j, iter;
     double chi, e, o, x;
 
     /* Calculate log-factorials */

And I get to the next bug:

> data(faithful)
> x <- xx <- faithful$eruptions
> x[i.out <- sample(length(x), 10)] <- NA
  I modifed this to
  x[i.out <- c(1,10,20,50,100,200)] <- NA
> doR <- density(x, bw=0.15, na.rm = TRUE) # works
> doN <- density(x, bw=0.15, na.rm = FALSE)
segfault at approx.c:

    for(i=0 ; i<*nxy ; i++)
        if(ISNA(x[i]) || ISNA(y[i]))
            error("approx(): attempted to interpolate NA values\n");

Both x and y are NULL pointers.
The error is somewhere before, comparing with R on SuSE 5.3:
somewhere in density():

debug: kords[(n + 2):(2 * n)] <- -kords[n:2]
linux and osf both get:
kords
[1] 0.00000000 0.01094819 0.02189638 0.03284457 0.04379277 ...
[1021] -0.04379277 -0.03284457 -0.02189638 -0.01094819

then in the following if-statement we go to
debug: kords <- convolve(y, kords)[1:n]
with linux result
kords
  [1] 2.523589e-13 4.199167e-13 6.950994e-13 1.144613e-12 ...
[511] 5.032797e-13 3.017881e-13
and alpha:
kords
  [1] NA NA NA NA NA NA NA NA ...
[501] NA NA NA NA NA NA NA NA NA NA NA NA
... now on the alpha x and y are set to numeric(0) and approx is called
with this arguments, bang.

Now in convolve(x,y) on entry x=y and y=kords on linux and osf are
identical.
x:
   [1] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 ...
 [151] 4.866728e-03 6.204044e-03 3.258272e-03 3.552390e-03 4.218750e-03
[1021] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 ...
y:
   [1] 2.659615e+00 2.652540e+00 2.631429e+00 2.596616e+00 ...
[1021] 2.548649e+00 2.596616e+00 2.631429e+00 2.652540e+00

now the error is in line
debug: x <- fft(fft(x) * (if (conj) Conj(fft(y)) else fft(y)), inv = TRUE)
of convolve() which returns only NAs on the alpha.

... so tomorrow I will have a look into fft.c

Albrecht

-------------------------------------------------------------------------------
Albrecht Gebhardt email : albrecht.gebhardt@uni-klu.ac.at
Institut fuer Mathematik Tel. : (++43 463) 2700/837
Universitaet Klagenfurt Fax : (++43 463) 2700/834
Villacher Str. 161
A-9020 Klagenfurt, Austria
-------------------------------------------------------------------------------

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._


About this list Date view Thread view Subject view Author view Other groups

This archive was generated by hypermail 2b25 : Tue 04 Jan 2000 - 14:16:08 EST