quadl crashes on evaluation

[David Bateman]

rmcd wrote:
> 
> Ben Abbott wrote:
>>
>> When you say "crash" what do you get, exactly?
>>
>> I get the following errors repeated many times
>>
>> error: evaluating binary operator `+' near line 169, column 3
>> ...
>>
>> If your function is modified
>>
>> 	f(x) = x .* sin (1 ./ (x+eps)) .* sqrt (abs (1 - x));
>>
>> and an integration tolerance is added
>>
>> 	v = quadl ('f' ,0, 3,. 0001)
>> 	v =  1.9818
>>
>> If any event, please describe what you mean by "crash".
>>
>>
> 
> Sorry for being unclear. By "crash" I mean that when I invoke quadl, octave
> exits with no error message and leaves me at the bash prompt (this is what I
> tried to describe at the start of my message). Presumably this should not
> happen. Under Ubuntu with 3.0.0 and 3.0.1 I just get the error messages, so
> this seems like a cygwin-specific bug. 
> 
> Your suggestions to add eps to the function and a tolerance to quadl work
> even under cygwin, thanks very much.
> 

It shouldn't crash as the issue here is only that the recursion limit of 
Octave is met. I suspect that cygwin doesn't like a 256 deep recursion 
as set by default with Octave and you should try something like

max_recursion_depth (128)
quadl (@(x) x .* sin (1 ./ x) .* sqrt (abs (1 - x)), 0, 3)

and see if it still crashes. If it doesn't, I'd suggest that the cygwin 
maximum recursion depth be reduced.

That being said, quadl is not really adapted to this type of singular 
integral. The QUADPACK based "quad" function is but is a little slow. 
The best for this type of singularity on the edge integrand is probably 
the quadgk function in the development branch of Octave as it explicitly 
extracts these singularities before treating the integrand and also uses 
vectorized calls to the integrand.

You can just take quadgk directly from the development tree and use it 
with Octave 3.0.x as it is not dependent on any other functions and 
neither is it dependent on any 3.1.x features.

Regards
David