以下是我的 SICP exercise 1.29 的代码。该练习要求我们实现
使用高阶过程 sum
的辛普森法则。应该是更多
比原来的integral
程序更准确。但我不知道为什么不是
我的代码中的情况:
(define (simpson-integral f a b n)
(define h (/ (- b a) n))
(define (next x) (+ x (* 2 h)))
(* (/ h 3) (+ (f a)
(* 4 (sum f (+ a h) next (- b h)))
(* 2 (sum f (+ a (* 2 h)) next (- b (* 2 h))))
(f b))))
我的代码的一些解释:如
h/3 * (y_{0} + 4*y_{1} + 2*y_{2} + 4*y_{3} + 2*y_{4} + ... + 2*y_{n-2} + 4*y_{n-1} + y_{n})
等于
h/3 * (y_{0}
+ 4 * (y_{1} + y_{3} + ... + y_{n-1})
+ 2 * (y_{2} + y_{4} + ... + y_{n-2})
+ y_{n})
我只是使用 sum
来计算 y_{1} + y_{3} + ... + y_{n-1}
和 y_{2 } +
y_{4} + ... + y_{n-2}
。
完整代码在这里:
#lang racket
(define (cube x) (* x x x))
(define (sum term a next b)
(if (> a b)
0
(+ (term a)
(sum term (next a) next b))))
(define (integral f a b dx)
(define (add-dx x) (+ x dx))
(* (sum f (+ a (/ dx 2.0)) add-dx b)
dx))
(define (simpson-integral f a b n)
(define h (/ (- b a) n))
(define (next x) (+ x (* 2 h)))
(* (/ h 3) (+ (f a)
(* 4 (sum f (+ a h) next (- b h)))
(* 2 (sum f (+ a (* 2 h)) next (- b (* 2 h))))
(f b))))
一些测试(准确值应该是0.25):
> (integral cube 0 1 0.01)
0.24998750000000042
> (integral cube 0 1 0.001)
0.249999875000001
> (simpson-integral cube 0 1.0 100)
0.23078806666666699
> (simpson-integral cube 0 1.0 1000)
0.24800798800666748
> (simpson-integral cube 0 1.0 10000)
0.2499999999999509
最佳答案
在您的解决方案中,x 值的计算方式如下:
h = (b-a)/n
x1 = a+1
x3 = x1 +2*h
x5 = x3 +2*h
...
这意味着舍入误差会慢慢累积。
当 (b-a)/n
无法表示为浮点时,就会发生这种情况。
如果我们将 xi
计算为 a+ (i*(b-a))/n
,您将获得更准确的结果。
此解决方案变体使用上述方法来计算 xi
。
(define (simpson-integral3 f a b n)
(define h (/ (- b a) n))
(define (next i) (+ i 2))
(define (f* i) (f (+ a (/ (* i (- b a)) n))))
(* (/ h 3)
(+ (f a)
(* 4 (sum f* 1 next n))
(* 2 (sum f* 2 next (- n 1)))
(f b))))
关于algorithm - 辛普森规则的实现(SICP 练习 1.29),我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/19706893/