bugGNU Octave - Bugs: bug #58636, Integral fails to call quadgk for...

 
 

bug #58636: Integral fails to call quadgk for complex integration

Submitted by:  Nicholas Jankowski <nrjank>
Submitted on:  Sat 20 Jun 2020 08:19:36 PM UTC  
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Unexpected Error
Status:  None Assigned to:  None
Originator Name:  Nicholas Jankowski Open/Closed:  Open
Release:  5.2.0 Operating System:  Any

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

Sun 21 Jun 2020 03:50:41 AM UTC, comment #5: 

ok. this version of integral has the limit checks, the feval(a,b) check, and a try/catch if those checks pass but the function still returns a complex answer. It currently works for the existing quadcc.  I added some complex input tests and it passes.

pending patch will tweak quadcc docstring to note it doesn't handle complex inputs or outputs, and once I can test that error messages split in quadcc, I can do that too and update integral to catch that particular error message.

doc fix in quad.txi and that should do it.

(file #49335)

Nicholas Jankowski <nrjank>
Sun 21 Jun 2020 01:40:05 AM UTC, comment #4: 

looking at quadcc, there's really only one error message that needs to be caught.  there are two that involve errors for complex numbers, but one is for complex values in the singularities list, and integral doesn't include that detail so never passes it to quadcc. The other one is:

if (fvals.length () != 1 || ! fvals(0).is_real_matrix ())
  error ("quadcc: integrand F must return a single, real-valued vector");

at three places in the code.  I don't know if I want this to catch that and run quadgk for both. looking at the code I'm not sure what sort of input would trigger that that first length error. I tried a couple inputs but couldn't get it to pop up. Maybe if I figure out how to trigger that I could see if it needs to be separated out.

Separating them is easy enough into:

if (fvals.length () != 1)
  error ("quadcc: integrand F must return a single vector");

if (! fvals(0).is_real_matrix ())
  error ("quadcc: integrand F must return a real-valued vector");

So if it doesn't cause issues integral could try/catch based on the first one.  I also noticed there are no BISTs for that error, so I'd want to add one for each after the split, if I can find a fn for the latter.

will play with the try/catch in integral now.

Nicholas Jankowski <nrjank>
Sat 20 Jun 2020 09:35:32 PM UTC, comment #3: 

need to see if there's a way to catch complex functions. something link

integral (exp (i*x), 0, 1) evaluates to real values at the real endpoints, so it'll sneak past the check. Also, not sure about whether there are other not-so-obvious complex functions. even something like

>> fx = @(x) complex(x)
fx =

@(x) complex (x)

>> quadgk(fx,0,1)
ans =  0.50000
>> quadcc(fx,0,1)
error: quadcc: integrand F must return a single, real-valued vector

needs to be caught.

Nicholas Jankowski <nrjank>
Sat 20 Jun 2020 09:15:51 PM UTC, comment #2: 

and this version checks the value of the function at the limits to catch things like
integral (sqrt(x), -1,1)

(file #49333)

Nicholas Jankowski <nrjank>
Sat 20 Jun 2020 08:54:28 PM UTC, comment #1: 

attached version of integral now works for complex limits using an iscomplex? check at the beginning. if there's an easy way to test the integrand, should be easy to add to that check and the whole function will work.  need to add some complex tests as well then this would be fixed.

(file #49332)

Nicholas Jankowski <nrjank>
Sat 20 Jun 2020 08:19:36 PM UTC, original submission:  

As reported by Brett on the mailing list [1]:

Integral is a wrapper for quadcc, quadgk, and quadv, and defaults to quadcc where possible. quadgk is the only quadrature method built for complex integration and limits, including complex path integration, and should be called.  There is currently no check for complex functions or limits, so if waypoints are not specified the script chooses quadcc which issues an error, while both quadgk and quadv would have successfully computed the result.

Complex function case:

>> fx = @(x) (x+i.*3).^2;
>> cmplx_int_val = integral(fx,-1,2)
error: quadcc: integrand F must return a single, real-valued vector
error: called from
    integral at line 114 column 7
>> cmplx_int_val = quadcc(fx,-1,2)
error: quadcc: integrand F must return a single, real-valued vector
>> cmplx_int_val = quadgk(fx,-1,2)
cmplx_int_val = -24.0000 +  9.0000i
>> cmplx_int_val = quadv(fx,-1,2)
cmplx_int_val = -24 +  9i

Complex limits case:

>> fx = @(x) (x).^2;
>> cmplx_int_val = integral(fx,-i,i)
error: quadcc: lower limit of integration (A) must be a real scalar
error: called from
    integral at line 114 column 7
>> cmplx_int_val = quadcc(fx,-i,i)
error: quadcc: lower limit of integration (A) must be a real scalar
>> cmplx_int_val = quadgk(fx,-i,i)
cmplx_int_val =  0.00000 - 0.66667i
>> cmplx_int_val = quadv(fx,-i,i)
cmplx_int_val = -0.00000 - 0.66667i

the isreal check on the limits should be easy enough.  the check on the function may be a bit more complicated, and I'm not sure the best way to check that.

it would be inefficient, but lacking something better could just if (catch specific complex error message from quadcc) --> kick it to quadgk. 

better suggestions welcome

(also, need to fix the docs for any change including fixing the preamble to the one function integration section that doesn't mention integral calls quadcc)

[1] https://octave.1599824.n4.nabble.com/Why-is-there-no-function-for-complex-integrals-tp4697602.html

Nicholas Jankowski <nrjank>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #49335:  integral.m added by nrjank (11KiB - text/plain)
file #49333:  integral.m added by nrjank (11KiB - text/plain - works for cmplex limits and functions where f(a) or f(b) are complex)
file #49332:  integral.m added by nrjank (10KiB - text/plain - integral function now working for complex limits, still not for cmplx integrand)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by nrjank (Submitted the item)
  • -email is unavailable- added by nrjank (reported on mailing list)
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

     

     

    Follow 4 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2020-06-21 nrjank Attached File- => Added integral.m, #49335
    2020-06-20 nrjank Attached File- => Added integral.m, #49333
    2020-06-20 nrjank Attached File- => Added integral.m, #49332
    2020-06-20 nrjank Carbon-Copy- => Added green.brett.r@gmail.com

    Back to the top


    Powered by Savane 3.5