46 lines
1.2 KiB
Nim
46 lines
1.2 KiB
Nim
type Function = proc(x: float): float
|
|
type Rule = proc(f: Function; x, h: float): float
|
|
|
|
proc leftRect(f: Function; x, h: float): float =
|
|
f(x)
|
|
|
|
proc midRect(f: Function; x, h: float): float =
|
|
f(x + h/2.0)
|
|
|
|
proc rightRect(f: Function; x, h: float): float =
|
|
f(x + h)
|
|
|
|
proc trapezium(f: Function; x, h: float): float =
|
|
(f(x) + f(x+h)) / 2.0
|
|
|
|
proc simpson(f: Function, x, h: float): float =
|
|
(f(x) + 4.0*f(x+h/2.0) + f(x+h)) / 6.0
|
|
|
|
proc cube(x: float): float =
|
|
x * x *x
|
|
|
|
proc reciprocal(x: float): float =
|
|
1.0 / x
|
|
|
|
proc identity(x: float): float =
|
|
x
|
|
|
|
proc integrate(f: Function; a, b: float; steps: int; meth: Rule): float =
|
|
let h = (b-a) / float(steps)
|
|
for i in 0 .. <steps:
|
|
result += meth(f, a+float(i)*h, h)
|
|
result = h * result
|
|
|
|
for fName, a, b, steps, fun in items(
|
|
[("cube", 0, 1, 100, cube),
|
|
("reciprocal", 1, 100, 1000, reciprocal),
|
|
("identity", 0, 5000, 5_000_000, identity),
|
|
("identity", 0, 6000, 6_000_000, identity)]):
|
|
|
|
for rName, rule in items({"leftRect": leftRect, "midRect": midRect,
|
|
"rightRect": rightRect, "trapezium": trapezium, "simpson": simpson}):
|
|
|
|
echo fName, " integrated using ", rName
|
|
echo " from ", a, " to ", b, " (", steps, " steps) = ",
|
|
integrate(fun, float(a), float(b), steps, rule)
|