Re: Wrong Integral result for a Piecewise function

*To*: mathgroup at smc.vnet.net*Subject*: [mg58561] Re: Wrong Integral result for a Piecewise function*From*: Maxim <ab_def at prontomail.com>*Date*: Thu, 7 Jul 2005 05:35:51 -0400 (EDT)*References*: <dag1bg$5j4$1@smc.vnet.net>*Sender*: owner-wri-mathgroup at wolfram.com

On Wed, 6 Jul 2005 07:30:56 +0000 (UTC), Dean Nairn <dnairn at udel.edu> wrote: > This integral gives gives the wrong result for the interval [2,3] > > h[x_] := Integrate[Boole[x - 1 < 2 y + 2 z < x], {y, 0, 1}, {z, 0, 1}] > > Plot[h[x],{x,0,5}] > and > Plot[Evaluate[h[x]],{x,0,5}] > > give different plots, The curve should be smooth and bell shaped from > 0 to 5, so the first looks correct. The second has a jump discontinuity > at 2 and 3. Also > > h[5/2] > and > h[x]/.x-> 5/2 > > give different answers > > This is using some new features in Mathematica 5.1. Same result on a > Mac (10.4) and SunOS (5.9). > > Breaking into a difference of two integrals gives the correct answer: > Integrate[Boole[ 2 y + 2 z < x], {y, 0, 1}, {z, 0, 1}] - > Integrate[Boole[ 2 y + 2 z < x-1], {y, 0, 1}, {z, 0, 1}] > > Finally > Integrate[h[x], {x, 0, 5}] and NIntegrate[h[x], {x, 0, 5}] > both give the wrong answer, it should be 1. The triple integral is > correct > > Integrate[Boole[x - 1 < 2 y + 2z < x], {y, 0, 1}, {z, 0, 1}, {x, 0, 5}] > > Any suggestions on integrating over regions with linear constraints? > Versions 5.1 has powerful new piecewise functions > Have a look at this thread from a few days ago: http://forums.wolfram.com/mathgroup/archive/2005/Jul/msg00042.html . The same two recipes apply here: either use PiecewiseIntegrate In[1]:= PiecewiseIntegrate[Boole[x - 1 < 2*y + 2*z < x], {y, 0, 1}, {z, 0, 1}] Out[1]= If[x == 3, 3/8, 0] + If[0 < x < 1, x^2/8, 0] + If[1 <= x <= 2, (1/8)*(-1 + 2*x), 0] + If[2 < x < 3, (1/8)*(-9 + 10*x - 2*x^2), 0] + If[Inequality[3, Less, x, LessEqual, 4], (1/8)*(9 - 2*x), 0] + If[4 < x < 5, (1/8)*(-5 + x)^2, 0] or integrate separately for ranges of the parameter values where the answer can be given as a single expression without conditions: In[2]:= Piecewise[ {Assuming[#, Integrate[Boole[x - 1 < 2*y + 2*z < x], {y, 0, 1}, {z, 0, 1}]], #}& /@ (Flatten[#, 1]&)@ Table[{i < x < i + 1, x == i + 1}, {i, 0, 4}] ] // PiecewiseExpand Out[2]= Piecewise[{{1/8, x == 1 || x == 4}, {3/8, x == 2 || x == 3}, {(1/8)*(9 - 2*x), 3 < x < 4}, {(1/8)*(-5 + x)^2, 4 < x < 5}, {x^2/8, 0 < x < 1}, {(1/8)*(-1 + 2*x), 1 < x < 2}, {(1/8)*(-9 + 10*x - 2*x^2), 2 < x < 3}}] This integral also can be derived from the general convolution formula. If xi_1, ..., xi_n are independent and uniformly distributed on [0, 1], then the CDF of their sum is given as F[n_, x_] = (1/n!)*Sum[(-1)^k*Binomial[n, k]*(x - k)^n*UnitStep[x - k], {k, 0, n}] Then the integral in question equals F[2, x/2] - F[2, (x - 1)/2]. Also note that if h[x] is called inside Plot, it still uses Integrate but with an approximate integrand (after substituting a numerical value for x), which is not the same as using NIntegrate. Besides, NIntegrate in versions 5.0 and 5.1 pre-evaluates the integrand, so the numerical integration proceeds with the wrong symbolic value of h[x]. Maxim Rytin m.r at inbox.ru