计算高阶函数的黎曼积分公式

Riemann integral formula to compute high-order function

SICP在1.3.1章节介绍了黎曼积分公式

(define (integral f a b dx)
  (define (add-dx x) (+ x dx))
  (* (sum f (+ a (/ dx 2.0)) add-dx b)
     dx))

将其应用于特定案例

#+name: case-1.3.1-integral.scm
#+BEGIN_SRC scheme :session sicp
(define pi 3.141592653589793)

(define (integral2 f a b dx)
  (define (add-dx x) (+ x dx))
  (* (sum (f b)
          (+ a (/ dx 2.0))
          (lambda (x) (+ x dx))
          b)
     dx))

(define (f b)
  (lambda (x) (/ 1 (sqrt
                    (- (sin x)
                       (sin b))))))


(* (integral2 f 0 (/ pi 6) 0.00001)
   (sqrt (/ 40
            (* 3 9.8))))
#+END_SRC

#+RESULTS: case-1.3.1-integral.scm
: 0.0-1.777598336021436i

得到完美答案 1.777598336021436

然后翻译成elisp

从小做起:

#+name: case-1.3.1-integral.el
#+begin_src emacs-lisp :session sicp :lexical t
(defun integral (f a b dx)
  (* (sum f
     (+ a (/ dx 2.0))
     (lambda (x) (+ x dx))
     b)
     dx))

(defun sum(term a next b)
  (if (> a b)
      0
    (+ (funcall term a)
       (sum term (funcall next a) next b))))

(integral #'cube 0 1 0.01)
#+end_src

#+RESULTS: case-1.3.1-integral.el
: 0.24998750000000042

它有效,因此用它来解决前面的问题

#+begin_src emacs-lisp :session sicp :lexical t 
(defvar pi 3.141592653589793)

(defun integral (f a b dx)
  (* (sum f
     (+ a (/ dx 2.0))
     (lambda (x) (+ x dx))
     b)
     dx))

(defun f (b)
  (lambda (x) (/ 1 (sqrt
                    (- (sin x)
                       (sin b))))))

(defun integral2 (f a b dx)
  (* (sum (funcall f b)
          (+ a (/ dx 2.0))
          (lambda (x) (+ x dx))
          b)
      dx))

(integral2 #'f 0 (/ pi 6) 0.01)

#+end_src

但它 return 一个毫无意义的结果

ELISP> (integral2 #'f 0 (/ pi 6) 0.01)
-0.0e+NaN

有什么问题?

你在使用Scheme时得到的答案是一个复数,调用sqrt的结果(你确定Scheme代码一开始是正确的吗?你应该仔细检查一下):

0.0-1.777598336021436i

不幸的是,Elisp 复数,这就是我们在那里得到 NaN 的原因。但这不是真正的问题;您应该调查 为什么 在 Scheme 代码中得到复杂的结果,积分不应该 return 复杂的值!