Numerical integration: Difference between revisions
Content added Content deleted
(Nimrod -> Nim) |
|||
Line 875: | Line 875: | ||
double s = integrate(f, 0.0, 1.0, 10, simpson());</lang> |
double s = integrate(f, 0.0, 1.0, 10, simpson());</lang> |
||
=={{header|Chapel}}== |
|||
<lang chapel> |
|||
proc f1(x:real):real { |
|||
return x**3; |
|||
} |
|||
proc f2(x:real):real { |
|||
return 1/x; |
|||
} |
|||
proc f3(x:real):real { |
|||
return x; |
|||
} |
|||
proc leftRectangleIntegration(a: real, b: real, N: int, f): real{ |
|||
var h: real = (b - a)/N; |
|||
var sum: real = 0.0; |
|||
var x_n: real; |
|||
for n in 0..N-1 { |
|||
x_n = a + n * h; |
|||
sum = sum + f(x_n); |
|||
} |
|||
return h * sum; |
|||
} |
|||
proc rightRectangleIntegration(a: real, b: real, N: int, f): real{ |
|||
var h: real = (b - a)/N; |
|||
var sum: real = 0.0; |
|||
var x_n: real; |
|||
for n in 0..N-1 { |
|||
x_n = a + (n + 1) * h; |
|||
sum = sum + f(x_n); |
|||
} |
|||
return h * sum; |
|||
} |
|||
proc midpointRectangleIntegration(a: real, b: real, N: int, f): real{ |
|||
var h: real = (b - a)/N; |
|||
var sum: real = 0.0; |
|||
var x_n: real; |
|||
for n in 0..N-1 { |
|||
x_n = a + (n + 0.5) * h; |
|||
sum = sum + f(x_n); |
|||
} |
|||
return h * sum; |
|||
} |
|||
proc trapezoidIntegration(a: real(64), b: real(64), N: int(64), f): real{ |
|||
var h: real(64) = (b - a)/N; |
|||
var sum: real(64) = f(a) + f(b); |
|||
var x_n: real(64); |
|||
for n in 1..N-1 { |
|||
x_n = a + n * h; |
|||
sum = sum + 2.0 * f(x_n); |
|||
} |
|||
return (h/2.0) * sum; |
|||
} |
|||
proc simpsonsIntegration(a: real(64), b: real(64), N: int(64), f): real{ |
|||
var h: real(64) = (b - a)/N; |
|||
var sum: real(64) = f(a) + f(b); |
|||
var x_n: real(64); |
|||
for n in 1..N-1 by 2 { |
|||
x_n = a + n * h; |
|||
sum = sum + 4.0 * f(x_n); |
|||
} |
|||
for n in 2..N-2 by 2 { |
|||
x_n = a + n * h; |
|||
sum = sum + 2.0 * f(x_n); |
|||
} |
|||
return (h/3.0) * sum; |
|||
} |
|||
var exact:real; |
|||
var calculated:real; |
|||
writeln("f(x) = x**3 with 100 steps from 0 to 1"); |
|||
exact = 0.25; |
|||
calculated = leftRectangleIntegration(a = 0.0, b = 1.0, N = 100, f = f1); |
|||
writeln("leftRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = rightRectangleIntegration(a = 0.0, b = 1.0, N = 100, f = f1); |
|||
writeln("rightRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = midpointRectangleIntegration(a = 0.0, b = 1.0, N = 100, f = f1); |
|||
writeln("midpointRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = trapezoidIntegration(a = 0.0, b = 1.0, N = 100, f = f1); |
|||
writeln("trapezoidIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = simpsonsIntegration(a = 0.0, b = 1.0, N = 100, f = f1); |
|||
writeln("simpsonsIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
writeln(); |
|||
writeln("f(x) = 1/x with 1000 steps from 1 to 100"); |
|||
exact = 4.605170; |
|||
calculated = leftRectangleIntegration(a = 1.0, b = 100.0, N = 1000, f = f2); |
|||
writeln("leftRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = rightRectangleIntegration(a = 1.0, b = 100.0, N = 1000, f = f2); |
|||
writeln("rightRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = midpointRectangleIntegration(a = 1.0, b = 100.0, N = 1000, f = f2); |
|||
writeln("midpointRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = trapezoidIntegration(a = 1.0, b = 100.0, N = 1000, f = f2); |
|||
writeln("trapezoidIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = simpsonsIntegration(a = 1.0, b = 100.0, N = 1000, f = f2); |
|||
writeln("simpsonsIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
writeln(); |
|||
writeln("f(x) = x with 5000000 steps from 0 to 5000"); |
|||
exact = 12500000; |
|||
calculated = leftRectangleIntegration(a = 0.0, b = 5000.0, N = 5000000, f = f3); |
|||
writeln("leftRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = rightRectangleIntegration(a = 0.0, b = 5000.0, N = 5000000, f = f3); |
|||
writeln("rightRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = midpointRectangleIntegration(a = 0.0, b = 5000.0, N = 5000000, f = f3); |
|||
writeln("midpointRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = trapezoidIntegration(a = 0.0, b = 5000.0, N = 5000000, f = f3); |
|||
writeln("trapezoidIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = simpsonsIntegration(a = 0.0, b = 5000.0, N = 5000000, f = f3); |
|||
writeln("simpsonsIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
writeln(); |
|||
writeln("f(x) = x with 6000000 steps from 0 to 6000"); |
|||
exact = 18000000; |
|||
calculated = leftRectangleIntegration(a = 0.0, b = 6000.0, N = 6000000, f = f3); |
|||
writeln("leftRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = rightRectangleIntegration(a = 0.0, b = 6000.0, N = 6000000, f = f3); |
|||
writeln("rightRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = midpointRectangleIntegration(a = 0.0, b = 6000.0, N = 6000000, f = f3); |
|||
writeln("midpointRectangleIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = trapezoidIntegration(a = 0.0, b = 6000.0, N = 6000000, f = f3); |
|||
writeln("trapezoidIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
calculated = simpsonsIntegration(a = 0.0, b = 6000.0, N = 6000000, f = f3); |
|||
writeln("simpsonsIntegration: calculated = ", calculated, "; exact = ", exact, "; difference = ", abs(calculated - exact)); |
|||
writeln(); |
|||
</lang> |
|||
output |
|||
<lang> |
|||
f(x) = x**3 with 100 steps from 0 to 1 |
|||
leftRectangleIntegration: calculated = 0.245025; exact = 0.25; difference = 0.004975 |
|||
rightRectangleIntegration: calculated = 0.255025; exact = 0.25; difference = 0.005025 |
|||
midpointRectangleIntegration: calculated = 0.249988; exact = 0.25; difference = 1.25e-05 |
|||
trapezoidIntegration: calculated = 0.250025; exact = 0.25; difference = 2.5e-05 |
|||
simpsonsIntegration: calculated = 0.25; exact = 0.25; difference = 5.55112e-17 |
|||
f(x) = 1/x with 1000 steps from 1 to 100 |
|||
leftRectangleIntegration: calculated = 4.65499; exact = 4.60517; difference = 0.0498211 |
|||
rightRectangleIntegration: calculated = 4.55698; exact = 4.60517; difference = 0.0481889 |
|||
midpointRectangleIntegration: calculated = 4.60476; exact = 4.60517; difference = 0.000407451 |
|||
trapezoidIntegration: calculated = 4.60599; exact = 4.60517; difference = 0.000816058 |
|||
simpsonsIntegration: calculated = 4.60517; exact = 4.60517; difference = 3.31627e-06 |
|||
f(x) = x with 5000000 steps from 0 to 5000 |
|||
leftRectangleIntegration: calculated = 1.25e+07; exact = 1.25e+07; difference = 2.5 |
|||
rightRectangleIntegration: calculated = 1.25e+07; exact = 1.25e+07; difference = 2.5 |
|||
midpointRectangleIntegration: calculated = 1.25e+07; exact = 1.25e+07; difference = 0.0 |
|||
trapezoidIntegration: calculated = 1.25e+07; exact = 1.25e+07; difference = 1.86265e-09 |
|||
simpsonsIntegration: calculated = 1.25e+07; exact = 1.25e+07; difference = 3.72529e-09 |
|||
f(x) = x with 6000000 steps from 0 to 6000 |
|||
leftRectangleIntegration: calculated = 1.8e+07; exact = 1.8e+07; difference = 3.0 |
|||
rightRectangleIntegration: calculated = 1.8e+07; exact = 1.8e+07; difference = 3.0 |
|||
midpointRectangleIntegration: calculated = 1.8e+07; exact = 1.8e+07; difference = 7.45058e-09 |
|||
trapezoidIntegration: calculated = 1.8e+07; exact = 1.8e+07; difference = 3.72529e-09 |
|||
simpsonsIntegration: calculated = 1.8e+07; exact = 1.8e+07; difference = 0.0 |
|||
</lang> |
|||
=={{header|CoffeeScript}}== |
=={{header|CoffeeScript}}== |
||
{{trans|python}} |
{{trans|python}} |