首页 > 其他 > 详细

SICP读书笔记(一)牛顿法求平方根

时间:2019-03-20 23:26:49      阅读:146      评论:0      收藏:0      [点我收藏+]

LispBox装完后试了一下,感觉还是不错的,于是找到了《计算机程序的构造和解释》,边学基础边练习。

书中使用的是Scheme,我依然使用CCL作为练习,专注思想的学习,同时折腾一下函数式风格。

 

实例:牛顿法求平方根

实数平方根定义:如果y2=x,且y>0,则称y为x的平方根

 

牛顿逐步逼近法求平方根步骤:

1.猜测x的平方根为 y,

2.计算y和x/y的平均值,

3.将步骤2的结果作为猜测值,

4.重复步骤2,3,直到误差允许范围

 

CCL 代码:

(defun squarex (x)
  (* x x))

(defun averagex (x y)
  (/ (+ x y) 2.0))

(defun absx (x)
  (if (< x 0) (- 0 x) x))

(defun good-enough? (x y c)
  (< (absx (- x y)) c))

(defun sqrt-improve (guess x)
  (averagex guess (/ x guess)))

(defun sqrt-iter (guess x c)
  (if (good-enough? (squarex guess) x c)
    guess
    (sqrt-iter (sqrt-improve guess x) x c)))

(defun sqrtx (x &key (c 0.001 c-supplied-p))
  (sqrt-iter x x (if (< c 0) 0.001 c)))

   

平方根逼近过程分析:

假设 x的平方根准确值为 a > 0,猜测值为 b > 0

则 (a - b)≥ 0, x = a2

=> a2 - 2ab + b2 ≥ 0

=> a2 + b2 ≥ 2ab

=> x + b≥ 2ab

=> a ≤ (x + b2) / 2b

=> a ≤ (x/b + b) / 2

由上可知,步骤2处理结果一定大于等于平方根,且仅当a = b时取等号。

当a = b时误差为0,故符合误差要求,直接得结果。

当a ≠ b时,有 a < (x/b + b) / 2

令 bn = (bn-1 + x / bn-1) / 2,b0 > 0

则 a < b1 = (b0 + x / b0) / 2,故从 b1后一定有 a < bn

当 n > 1,

bn =  (bn-1 + x/bn-1 ) / 2

令 dn = bn - bn-1 =  (bn-1 + x/bn-1 ) / 2 - (2 * bn-1 / 2)

=> 2 * dn = bn-1 + x/bn-1 - 2*bn-1 = x/bn-1 - bn-1 

因 x / bn-1 = a2 / bn-1  ,且 a2 ≤ a * bn < bn2

=> x / bn-1 < bn-12 / bn-1 = bn-1

=> x/bn-1 - bn-1 < 0

=> dn < 0

=> bn - bn-1 < 0

=> bn < bn-1

=> a < bn < bn-1

综上所述,当 初始值 不等于 精确值 时,使用 步骤2 迭代,返回结果会不断减小,并且一直大于 精确值。

如此理论平方根的平方应该略大于原值,但执行结果上存在小于原值的情况,可能是 浮点 精度损失导致。


附上:牛顿逼近法

求解某个方程,都可以将其转换为 f(x) = 0 的形式。

在图形上,表述为 函数 f(x) 与 x 轴的交点。

此时,可以任取一点作切线,切线代表了所取点的变化,它与x轴的交点便代表了 f(x) 与 x轴相交的方向。使用交点继续作切线,同理下其,最终会接近真实的交点。

技术分享图片

这个过程的导数形式:xn+1 = xn - f(xn) / f‘(xn)

求平方根 即是 a2=y => f(x) = x2 - y

求 f(x) = 0时的 x值。

f‘(x) = 2 * x

=> xn+1 = xn - (xn2 - y) / (2 * xn)

=> xn+1 = (xn2 + y) / (2 * xn)

=> xn+1 = (xn + y/xn) / 2.

 

lisp自带求平方根函数 sqrt 源码:

(defun sqrt (x &aux a b)
  "Return the square root of NUMBER."
  (cond ((zerop x) x)
        ((complexp x) (* (sqrt (abs x)) (cis (/ (phase x) 2))))          
        ((minusp x) (complex 0 (sqrt (- x))))
        ((floatp x)
         (fsqrt x))
        ((and (integerp x) (eql x (* (setq a (isqrt x)) a))) a)
        ((and (ratiop x)
              (let ((n (numerator x))
                    d)
                (and (eql n (* (setq a (isqrt n)) a))
                     (eql (setq d (denominator x))
                          (* (setq b (isqrt d)) b)))))
         (/ a b))          
        (t
         #+32-bit-target
         (target::with-stack-short-floats ((f1))
           (fsqrt (%short-float x f1)))
         #+64-bit-target
         (fsqrt (%short-float x)))))

  

SICP读书笔记(一)牛顿法求平方根

原文:https://www.cnblogs.com/qq1144054302/p/10568392.html

(0)
(0)
   
举报
评论 一句话评论(0
关于我们 - 联系我们 - 留言反馈 - 联系我们:wmxa8@hotmail.com
© 2014 bubuko.com 版权所有
打开技术之扣,分享程序人生!