Mercurial > hg > octave-lyh

changeset 15437:c9954a15bc03 stable

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Fix quadcc when there are NaNs in the integrand (bug #37414) * quadcc.cc: Use fx[nans[i]] = octave_NaN instead of fx[i]. Add %!test for integrand with NaNs in it.
author Rik <rik@octave.org>
date Sun, 23 Sep 2012 10:36:25 -0700
parents fd5c0159b588
children e39a51e0d54b d174210ce1ec
files src/DLD-FUNCTIONS/quadcc.cc
diffstat 1 files changed, 17 insertions(+), 4 deletions(-) [+]
line wrap: on
line diff
--- a/src/DLD-FUNCTIONS/quadcc.cc
+++ b/src/DLD-FUNCTIONS/quadcc.cc
@@ -1715,7 +1715,7 @@
       Vinvfx (iv->fx, &(iv->c[idx[2]]), 2);
       Vinvfx (iv->fx, &(iv->c[0]), 0);
       for (i = 0; i < nnans; i++)
-        iv->fx[i] = octave_NaN;
+        iv->fx[nans[i]] = octave_NaN;
       iv->a = iivals[j];
       iv->b = iivals[j + 1];
       iv->depth = 3;
@@ -1844,7 +1844,7 @@
             {
               downdate (&(iv->c[idx[d]]), n[d], d, nans, nnans);
               for (i = 0; i < nnans; i++)
-                iv->fx[i] = octave_NaN;
+                iv->fx[nans[i]] = octave_NaN;
             }
 
           /* Compute the error estimate. */
@@ -1989,7 +1989,7 @@
             {
               downdate (ivl->c, n[0], 0, nans, nnans);
               for (i = 0; i < nnans; i++)
-                ivl->fx[i] = octave_NaN;
+                ivl->fx[nans[i]] = octave_NaN;
             }
           for (i = 0; i <= n[d]; i++)
             {
@@ -2085,7 +2085,7 @@
             {
               downdate (ivr->c, n[0], 0, nans, nnans);
               for (i = 0; i < nnans; i++)
-                ivr->fx[i] = octave_NaN;
+                ivr->fx[nans[i]] = octave_NaN;
             }
           for (i = 0; i <= n[d]; i++)
             {
@@ -2251,6 +2251,19 @@
 %!assert (quadcc (@(x) exp(-x .^ 2), -Inf, Inf), sqrt(pi), 1e-6)
 %!assert (quadcc (@(x) exp(-x .^ 2), -Inf, 0), sqrt(pi)/2, 1e-6)
 
+## Test function with NaNs in interval 
+%!function y = __nansin (x)
+%!  nan_locs = [-3*pi/4, -pi/4, 0, pi/3, pi/2, pi];
+%!  y = sin (x);
+%!  idx = min (abs (bsxfun (@minus, x(:), nan_locs)), [], 2); 
+%!  y(idx < 1e-10) = NaN;
+%!endfunction 
+
+%!test
+%! [q, err, npoints] = quadcc ("__nansin", -pi, pi); 
+%! assert (q, 0, eps);
+%! assert (err, 0, 15*eps);
+
 %% Test input validation
 %!error (quadcc ())
 %!error (quadcc (@sin))