From hwborchers at googlemail.com  Fri Jan 27 13:23:45 2012
From: hwborchers at googlemail.com (Hans W Borchers)
Date: Fri, 27 Jan 2012 13:23:45 +0100
Subject: [Rd] Numerical instability in new R Windows development version
Message-ID: <CAML4n3MNRsp_kTpPHarDh7_kLCVwGZJkweOxkoo6S8OMwENfjg@mail.gmail.com>

I have a question concerning the new Windows toolchain for R >= 2.14.2.
When trying out my package 'pracma' on the win-builder development version
it will stop with the following error message:

  > f3 <- function(x, y) sqrt((1 - (x^2 + y^2)) * (x^2 + y^2 <= 1))
  > dblquad(f3, -1, 1, -1, 1)     #   2.094395124 , i.e. 2/3*pi , err = 2e-8
  Warning in sqrt((1 - (x^2 + y^2)) * (x^2 + y^2 <= 1)) : NaNs produced
  Warning in sqrt((1 - (x^2 + y^2)) * (x^2 + y^2 <= 1)) : NaNs produced
  Error in integrate(function(y) f(x, y), ya, yb, subdivisions = subdivs,  :
    non-finite function value
  Calls: dblquad ...
         <Anonymous> -> f -> do.call -> mapply -> <Anonymous> -> integrate
  Execution halted
  ** running examples for arch 'x64' ... ERROR
  Running examples in 'pracma-Ex.R' failed

This probably means that the following expression got negative for some
values x, y:

  (1 - (x^2 + y^2)) * (x^2 + y^2 <= 1)

It appears to be an often used trick in numerical analysis. One advantage is
that a function using it is immediately vectorized while an expression such
as, e.g., "max(0, 1 - (x^2 + y^2))" is not.

The example runs fine on Debian Linux and Mac OS X 32-/64-bit architectures.
In my understanding the approach is correct and, as said above, often used in
numerical applications.

Can someone explain to me why this fails for the Windows 64-bit compiler and
what I should use instead. Thanks.

Hans Werner Borchers
ABB Corporate Research



[Rd] Numerical instability in new R Windows development version
Hans W Borchers