Numerical integration: Difference between revisions
Content added Content deleted
m (→{{header|REXX}}: added whitespace, used a template for the output section.) |
|||
Line 1,161: | Line 1,161: | ||
</lang> |
</lang> |
||
=={{header| |
=={{header|Comal}}== |
||
{{works with|OpenComal on Linux}} |
|||
<lang Comal> |
|||
<lang lisp>(defun left-rectangle (f a b n &aux (d (/ (- b a) n))) |
|||
1000 PRINT "F(X)";" FROM";" TO";" L-Rect";" M-Rect";" R-Rect ";" Trapez";" Simpson" |
|||
(* d (loop for x from a below b by d summing (funcall f x)))) |
|||
1010 fromval:=0 |
|||
1020 toval:=1 |
|||
(defun right-rectangle (f a b n &aux (d (/ (- b a) n))) |
|||
1030 PRINT "X^3 "; |
|||
(* d (loop for x from b above a by d summing (funcall f x)))) |
|||
1040 PRINT USING "#####": fromval; |
|||
1050 PRINT USING "#####": toval; |
|||
(defun midpoint-rectangle (f a b n &aux (d (/ (- b a) n))) |
|||
1060 PRINT USING "###.#########": numint(f1, "L", fromval, toval, 100); |
|||
(* d (loop for x from (+ a (/ d 2)) below b by d summing (funcall f x)))) |
|||
1070 PRINT USING "###.#########": numint(f1, "R", fromval, toval, 100); |
|||
1080 PRINT USING "###.#########": numint(f1, "M", fromval, toval, 100); |
|||
(defun trapezium (f a b n &aux (d (/ (- b a) n))) |
|||
1090 PRINT USING "###.#########": numint(f1, "T", fromval, toval, 100); |
|||
(* (/ d 2) |
|||
1100 PRINT USING "###.#########": numint(f1, "S", fromval, toval, 100) |
|||
(+ (funcall f a) |
|||
1110 // |
|||
(* 2 (loop for x from (+ a d) below b by d summing (funcall f x))) |
|||
1120 fromval:=1 |
|||
1130 toval:=100 |
|||
1140 PRINT "1/X "; |
|||
(defun simpson (f a b n) |
|||
1150 PRINT USING "#####": fromval; |
|||
(loop with h = (/ (- b a) n) |
|||
1160 PRINT USING "#####": toval; |
|||
with sum1 = (funcall f (+ a (/ h 2))) |
|||
1170 PRINT USING "###.#########": numint(f2, "L", fromval, toval, 1000); |
|||
with sum2 = 0 |
|||
1180 PRINT USING "###.#########": numint(f2, "R", fromval, toval, 1000); |
|||
for i from 1 below n |
|||
1190 PRINT USING "###.#########": numint(f2, "M", fromval, toval, 1000); |
|||
do (incf sum1 (funcall f (+ a (* h i) (/ h 2)))) |
|||
1200 PRINT USING "###.#########": numint(f2, "T", fromval, toval, 1000); |
|||
do (incf sum2 (funcall f (+ a (* h i)))) |
|||
1210 PRINT USING "###.#########": numint(f2, "S", fromval, toval, 1000) |
|||
finally (return (* (/ h 6) |
|||
1220 fromval:=0 |
|||
(+ (funcall f a) |
|||
1230 toval:=5000 |
|||
(funcall f b) |
|||
1240 PRINT "X "; |
|||
(* 4 sum1) |
|||
1250 PRINT USING "#####": fromval; |
|||
(* 2 sum2))))))</lang> |
|||
1260 PRINT USING "#####": toval; |
|||
1270 PRINT USING "#########.###": numint(f3, "L", fromval, toval, 5000000); |
|||
1280 PRINT USING "#########.###": numint(f3, "R", fromval, toval, 5000000); |
|||
1290 PRINT USING "#########.###": numint(f3, "M", fromval, toval, 5000000); |
|||
1300 PRINT USING "#########.###": numint(f3, "T", fromval, toval, 5000000); |
|||
1310 PRINT USING "#########.###": numint(f3, "S", fromval, toval, 5000000) |
|||
1320 // |
|||
1330 fromval:=0 |
|||
1340 toval:=6000 |
|||
1350 PRINT "X "; |
|||
1360 PRINT USING "#####": fromval; |
|||
1370 PRINT USING "#####": toval; |
|||
1380 PRINT USING "#########.###": numint(f3, "L", fromval, toval, 6000000); |
|||
1390 PRINT USING "#########.###": numint(f3, "R", fromval, toval, 6000000); |
|||
1400 PRINT USING "#########.###": numint(f3, "M", fromval, toval, 6000000); |
|||
1410 PRINT USING "#########.###": numint(f3, "T", fromval, toval, 6000000); |
|||
1420 PRINT USING "#########.###": numint(f3, "S", fromval, toval, 6000000) |
|||
1430 END |
|||
1440 // |
|||
1450 FUNC numint(FUNC f, type$, lbound, rbound, iters) CLOSED |
|||
1460 delta:=(rbound-lbound)/iters |
|||
1470 integral:=0 |
|||
1480 CASE type$ OF |
|||
1490 WHEN "L", "T", "S" |
|||
1500 actval:=lbound |
|||
1510 WHEN "M" |
|||
1520 actval:=lbound+delta/2 |
|||
1530 WHEN "R" |
|||
1540 actval:=lbound+delta |
|||
1550 OTHERWISE |
|||
1560 actval:=lbound |
|||
1570 ENDCASE |
|||
1580 FOR n:=0 TO iters-1 DO |
|||
1590 CASE type$ OF |
|||
1600 WHEN "L", "M", "R" |
|||
1610 integral:+f(actval+n*delta)*delta |
|||
1620 WHEN "T" |
|||
1630 integral:+delta*(f(actval+n*delta)+f(actval+(n+1)*delta))/2 |
|||
1640 WHEN "S" |
|||
1650 IF n=0 THEN |
|||
1660 sum1:=f(lbound+delta/2) |
|||
1670 sum2:=0 |
|||
1680 ELSE |
|||
1690 sum1:+f(actval+n*delta+delta/2) |
|||
1700 sum2:+f(actval+n*delta) |
|||
1710 ENDIF |
|||
1720 OTHERWISE |
|||
1730 integral:=0 |
|||
1740 ENDCASE |
|||
1750 ENDFOR |
|||
1760 IF type$="S" THEN |
|||
1770 RETURN (delta/6)*(f(lbound)+f(rbound)+4*sum1+2*sum2) |
|||
1780 ELSE |
|||
1790 RETURN integral |
|||
1800 ENDIF |
|||
1810 ENDFUNC |
|||
1820 // |
|||
1830 FUNC f1(x) CLOSED |
|||
1840 RETURN x^3 |
|||
1850 ENDFUNC |
|||
1860 // |
|||
1870 FUNC f2(x) CLOSED |
|||
1880 RETURN 1/x |
|||
1890 ENDFUNC |
|||
1900 // |
|||
1910 FUNC f3(x) CLOSED |
|||
1920 RETURN x |
|||
1930 ENDFUNC |
|||
</lang> |
|||
{{out}} |
|||
<pre> |
|||
F(X) FROM TO L-Rect M-Rect R-Rect Trapez Simpson |
|||
X^3 0 1 0.245025000 0.255025000 0.249987500 0.250025000 0.250000000 |
|||
1/X 1 100 4.654991058 4.556981058 4.604762549 4.605986058 4.605170385 |
|||
X 0 5000 12499997.500 12500002.500 12500000.000 12500000.000 12500000.000 |
|||
X 0 6000 17999997.000 18000003.000 18000000.000 18000000.000 18000000.000 |
|||
</pre> |
|||
=={{header|D}}== |
=={{header|D}}== |