Numerical integration/Adaptive Simpson's method: Difference between revisions
Numerical integration/Adaptive Simpson's method (view source)
Revision as of 13:02, 19 May 2023
, 1 year agoAdded Icon.
(→{{header|COBOL}}: Bug fix in the use of delta.) |
(Added Icon.) |
||
Line 940:
Simpson's integration of sine from 0 to 1 = 0.459698
</pre>
=={{header|Icon}}==
<syntaxhighlight lang="icon">
record simpson_result (m, fm, quadval)
procedure main ()
write ("estimate of ∫ sin x dx from 0 to 1: ",
quad_asr (sin, 0.0, 1.0, 1E-9, 100))
end
procedure quad_asr(f, a, b, tol, depth)
local fa, fb, whole
fa := f(a)
fb := f(b)
whole := simpson_rule (f, a, fa, b, fb)
return recursive_simpson(f, a, fa, b, fb, tol, whole.quadval,
whole.m, whole.fm, depth)
end
procedure recursive_simpson (f, a, fa, b, fb, tol, whole,
m, fm, depth)
local left, right, delta, tol2
local leftval, rightval, quadval
left := simpson_rule (f, a, fa, m, fm)
right := simpson_rule (f, m, fm, b, fb)
delta := left.quadval + right.quadval - whole
tol2 := 0.5 * tol
if depth <= 0 | tol2 = tol | abs (delta) <= 15 * tol then
quadval := left.quadval + right.quadval + (delta / 15)
else
{
leftval := recursive_simpson (f, a, fa, m, fm, tol2,
left.quadval, left.m, left.fm,
depth - 1)
rightval := recursive_simpson (f, m, fm, b, fb, tol2,
right.quadval, right.m, right.fm,
depth - 1)
quadval := leftval + rightval
}
return quadval
end
procedure simpson_rule (f, a, fa, b, fb)
local m, fm, quadval
m := 0.5 * (a + b)
fm := f(m)
quadval := ((b - a) / 6) * (fa + (4 * fm) + fb)
return simpson_result (m, fm, quadval)
end
</syntaxhighlight>
{{out}}
<pre>$ icon adaptive_simpson_task_arizona.icn
estimate of ∫ sin x dx from 0 to 1: 0.4596976941</pre>
=={{header|J}}==
| |||