Chapter Eight of "Maxima By Example"



On July 28, 2009, Mike Rotenberry wrote:
-------------------------------------



> I have been reading through chapter eight of your excellent "Maxima By
> Example" tutorial. In section 8.2.5 it could be inferred that
> quad_qags cannot be used to calculate double integrals, but this has
> not been my experience. Using your example 2 in section 8.2.6:
>
> Maxima 5.18.1 http://maxima.sourceforge.net
> Using Lisp SBCL 1.0.18.debian

> (%i1) 
> quad_qags(lambda([x],quad_qags(lambda([y],exp(x-y^2)),y,1,2+x)[1]),x,0,1);
> (%o1)          [.2384683615524653, 2.647530656401834e-15, 21, 0]
>
>
> I have an interest in determining knot energies, and double integrals
> arise naturally in these calculations. My tests show that for double
> integrals quad_qags usually provides a better estimate than romberg.
>
> Regards,
> James M Rotenberry
------------------------------------
Hi Mike,

Thanks for the example of use of quad_qags for a double integral in which
the inner integral has variable limits.

I will modify Ch. 8 to include your example.

Here is a comparison for the first example in that section which shows
 quad_qags is a little more accurate for this example.
-----------------------------------------------------------------------
(%i1) fpprintprec:8$
(%i2) g : x*y/(x+y)$
(%i3) tval : bfloat( integrate( integrate(g,y,0,x/2),x,1,3) ),fpprec:20;
Is  x  positive, negative, or zero?

p;
(%o3)                            8.1930239b-1
(%i4) load(romberg)$
(%i5) [rombergtol,rombergabs,rombergit,rombergmin];
(%o5)                        [1.0E-4, 0.0, 11, 0]
(%i6) rombergtol:1e-15$
(%i7) romberg( romberg(g,y,0,x/2),x,1,3);
(%o7)                              0.819302

romberg absolute error:
(%i8) abs(% - tval),fpprec:20;
(%o8)                            1.7298445b-16

romberg relative error:
(%i9) %/tval,fpprec:20;
(%o9)                            2.1113627b-16

(%i10) quad_qags(lambda([x], 
first(quad_qags(lambda([y],g),y,0,x/2))),x,1,3);

(%o10) quad_qags(lambda([x], first(quad_qags(lambda([y], g), y, 0, x/2))), 
x, 1,
                                 3, epsrel = 1.0E-8, epsabs = 0.0, limit = 
200)

(%i11) quad_qags(lambda([x], 
first(quad_qags(lambda([y],x*y/(x+y)),y,0,x/2))),
                                                                      x,1,3);
(%o11)                 [0.819302, 9.09608385E-15, 21, 0]

quad_qags absolute error:
(%i12) abs(%[1] - tval), fpprec:20;
(%o12)                           6.1962154b-17

quad_qags relative error:
(%i13) %/tval,fpprec:20;
(%o13)                           7.5627942b-17
---------------------------------------------------
Evidently, we cannot use the assigned expression "g" as we can
in romberg and integrate, but we can use a Maxima function instead:
-----------------------------------------

(%i14) define  ( f (x, y),  g)$

(%i16) quad_qags(lambda([x], 
first(quad_qags(lambda([y],f(x,y)),y,0,x/2))),x,1,3);

(%o16)                 [0.819302,  9.09608385E-15,  21,  0]
---------------------------------------------------------

Best Wishes,
Ted