Re: Integration of function containing if-test

[Torquil Macdonald Sørensen]

Søren Hauberg wrote:
> tor, 05 02 2009 kl. 18:17 +0100, skrev Torquil Macdonald Sørensen:
> > Hi!
> >
> > I am having some problems with integration of the following function 
> > (this is really a simplification of my function, but the problem I 
> > illustrate is the same).
> >
> > I have a function nu(x) that is defined using an if-test. It has a 
> > special value (0) at x = 0, and is defined in terms of an ordinary 
> > function f(x) everywhere else. Therefore I have:
> >
> > function y = nu(x)
> >         if (x == 0)
> >                 y = 0;
> >         else
> >                 y = f(x);
> >         endif
> > endfunction
> >
> > When I pass it a 0 it returns a zero. But when I pass it a vector [0 1] 
> > it doesn't work, because the if-test is apparently not executed on each 
> > element of the vector individually, but on the vector x as a whole. 
> > Since x = [0 1] is not equal to 0, the first part of the if-test is 
> > never entered, and therefore f(x) is evaluated at both 0 and 1, thereby 
> > returning the wrong result. In my case it returns [ NaN 1], since f(0) = 
> > Nan and f(1) = 1.
>
> I'm not quite I understand your problem, but:
>
> 1) Can you just do
>
>    function y = nu(x)
>      if any (x == 0)
>        y = 0;
>      else
>        y = f(x);
>      endif
>    endfunction
>
> ? Notice the 'any'.

Thanks Søren, actually this method didn't work for me as is stands, 
because I need it to return a vector even if one of the input values are 
0. This returns "0" of at least one element of x is 0.

> 2) Since you only have a finite number of points where 'x' is zero, then
> I would have thought you could just ignore them when doing integration.
> Can't you just integrate 'f (x)' directly?

Yes, I decided to do this and it works fine. I had looked at the help 
for quadl, not quad, so I didn't think of the possibility to specify a 
point where the function is "singular", as described in "help quad". 
Actually the function I want to integrate is well-defined at x = 0, but 
its defining expression in octave is given  as x*x*g(x) where g(0) = 
NaN. And since 0*NaN is also NaN, I need to treat it as a special case.

Thanks again
Torquil Sørensen