File: dblint.dem

package info (click to toggle)
maxima 5.6-17
  • links: PTS
  • area: main
  • in suites: woody
  • size: 30,572 kB
  • ctags: 47,715
  • sloc: ansic: 154,079; lisp: 147,553; asm: 45,843; tcl: 16,744; sh: 11,057; makefile: 7,198; perl: 1,842; sed: 334; fortran: 24; awk: 5
file content (55 lines) | stat: -rw-r--r-- 1,741 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
/* Example of use of double integrator (DBLINT) */

/* F must be a function of two variables, R and S must be functions of
one variable, and A and B must be floating point numbers (all of the 
functions must be translated or compiled prior to using DBLINT) */

/* Get the fasl file containing DBLINT(F,R,S,A,B) */ 

batchload("dblint"); 

/* To get the double integral of exp(x-y^2) over the region bounded by
y=1 and y=2+x and x=0 to x=1 */

/* Define the integrand as a function of two variables */

F(X,Y):=(MODE_DECLARE([X,Y],FLOAT),EXP(X-Y^2)); 

/* Define the lower and upper limits on the inner (y in this case)
integral as a function of the outer variable (x in this case) */

R(X):=(MODE_DECLARE(X,FLOAT),1.0); 
S(X):=(MODE_DECLARE(X,FLOAT),2.0+X);

/* Now translate these functions for the sake of efficiency */

TRANSLATE(F,R,S); 

/* Call the DBLINT function with quoted arguments for function names, and
floating point values for the endpoints of the outer (x) integration
*/ 

DBLINT_ANS:DBLINT('F,'R,'S,0.0,1.0); 

/* Compare with the exact answer */

INTY:RISCH(EXP(X-Y^2),Y); 
XINT:EV(INTY,Y:2+X)-EV(INTY,Y:1);
DINT:RISCH(XINT,X); 
ANS:EV(DINT,X:1.0)-EV(DINT,X:0.0),NUMER; 

/* Relative error check here */
(ANS-DBLINT_ANS)/ANS;
/* Quite reasonable */
/* Of course, the DBLINT routine will still work even when we choose an
area over which the closed-form integral fails to be expressible in
terms of standard transcendental functions */ 

S1(X):=(MODE_DECLARE(X,FLOAT),2.0+X^(3/2));
TRANSLATE(S1);
DBLINT('F,'R,'S1,0.0,1.0); 

/* To see what happens when this integral is attempted in closed-form,
you may BATCH(DBLINT,DEMO1,DSK,SHARE1), and that will also show sample use
of FLOATDEFUNK and QUANC8 */
"end of dblint.dem"$