保持精度的小trick:Kahan's summation Formula

由于最近用GPU编程,涉及到了float数组,就不得不涉及精度问题。在 CPU 上进行计算时,我们使用 double(即 64 bits 浮点数)来累进计算过程,而在 GPU 上则只能用 float(32 bits 浮点数)。在累加大量数字的时候,由于累加结果很快会变大,因此后面的数字很容易被舍去过多的位数。对于双精度如C中double以及Fortran中real(kind = 8),一般运算的精度足以保持,但是单精度数组,在大量操作后极易出现“大数吃小数”等不稳定现象。在不能使用更高精度数组的前提下,可以用一个小技巧来保持精度:Kahan求和。

function KahanSum(input)
    var sum = 0.0
    var c = 0.0             //A running compensation for lost low-order bits.
    for i = 1 to input.length do
        y = input[i] - c    //So far, so good: c is zero.
        t = sum + y         //Alas, sum is big, y small, so low-order digits of y are lost.
        c = (t - sum) - y   //(t - sum) recovers the high-order part of y; subtracting y recovers -(low part of y)
        sum = t             //Algebraically, c should always be zero. Beware eagerly optimising compilers!
        //Next time around, the lost low part will be added to y in a fresh attempt.
    return sum

见下面一段Fortran代码:

program main
implicit none
    integer, parameter:: N = 1000000
    integer i
    real(kind = 4), parameter:: ELEMENT = 0.001
    real(kind = 4) s, eps, y, t

write (*, "('Theoretical value: ', F10.5)") N*ELEMENT

s = 0.0
    do i = 1, N
        s = s+ELEMENT
    enddo
    write (*, "('Naive method value: ', F10.5)") s

s = 0.0
    eps = 0.0
    do i = 1, N
        y = ELEMENT-eps
        t = s+y
        eps = (t-s)-y
        s = t
    enddo
    write (*, "('Kahan method value: ',F10.5)") s

stop
end

运行结果为:

Theoretical value: 1000.00006
Naive method value:  991.14154
Kahan method value: 1000.00006

对于N个0.001,普通方法累加到991左右就已经丢失精度了。可以看到用“Kahan method”能够得到近乎于理论的精度数值。分析一下他的原理。我们发现,如果没有精度损失,eps永远为0,y就是ELEMENT=0.001。一旦在 i 到了某个数值出现了大数吃小数 的情形时,不妨激进的设小数部分全部被截断,则如s = 991.0000时,由于eps之前为0,则y=0.0010.之后t=s+y,得到的就是“吃掉”的结果,如991.0000,绝对误差达0.001.此时:eps=(t-s)-y=(991.0000-991.0000)-0.0010=-0.001,可见eps起了保存“损失位”的作用。此时s=t=991.0000.下个循环:y = 0.001--0.001=0.002,t = s+y=991.0000,eps=-0.002,如此反复,这样足够多循环后,eps足可以复现大的校正值,从而保证结果的高精度。当eps足够大时候,(t-s)-y=0,从而使eps重新为0,继续起保存损失的作用。

Kahan's summation Formula相关推荐

  1. Kahan's Summation Formula原理—它是如何处理大数吃小数的

    Kahan's Summation Formula原理-它是如何避免大数吃小数的 Kahan求和公式原理: 首先,这个算法就是用来求和的,求a1+a2+a3+...为什么不直接相加呢,而要用Kahan ...

  2. 深入浅出CUDA编程

    标签: cuda编程threadfloatconflictexpress 2010-12-10 13:29 44960人阅读 评论(7) 收藏 举报 CUDA 是 NVIDIA 的 GPGPU 模型, ...

  3. 深入浅出谈CUDA(二)

    前面介绍的计算平方和的程序,似乎没有什么实用价值.所以我们的第二个 CUDA 程序,要做一个确实有(某些)实用价值的程序,也就是进行矩阵乘法.而且,这次我们会使用浮点数. 虽然矩阵乘法有点老套,不过因 ...

  4. bookmarks_2021_9_28

    书签栏 通讯 s7-1200与s7-200smart通讯-工业支持中心-西门子中国 IO_device S7-1200PROFINET通信 ET 200SP 安装视频 - ID: 95886218 - ...

  5. CUDA: 矩阵乘法优化

    矩阵乘法是有实用价值的程序,我们会使用浮点数. 虽然矩阵乘法有点老套,不过因为它相当简单,而且也可以用来介绍一些有关 CUDA 的有趣性质. 矩阵乘法 为了单纯起见,我们这里以方形的矩阵为例子.基本上 ...

  6. CUDA 深入浅出谈[转]

    CUDA 深入浅出谈           "CUDA 是 NVIDIA 的 GPGPU 模型,它使用 C 语言为基础,可以直接以大多数人熟悉的 C 语言,写出在显示芯片上执行的程序,而不需要 ...

  7. CUDA 深入浅出谈

    CUDA 深入浅出谈           "CUDA 是 NVIDIA 的 GPGPU 模型,它使用 C 语言为基础,可以直接以大多数人熟悉的 C 语言,写出在显示芯片上执行的程序,而不需要 ...

  8. CUDA编程深入浅出,案列讲解

    CUDA 是 NVIDIA 的 GPGPU 模型,它使用 C 语言为基础,可以直接以大多数人熟悉的 C 语言,写出在显示芯片上执行的程序,而不需要去学习特定的显示芯片的指令或是特殊的结构." ...

  9. 矩阵乘法——CUDA 优化记录

    CUDA 是 NVIDIA 的 GPGPU 模型,它使用 C 语言为基础,可以直接以大多数人熟悉的 C 语言,写出在显示芯片上执行的程序,而不需要去学习特定的显示芯片的指令或是特殊的结构." ...

最新文章

  1. 删除数组中指定元素_如何删除PHP数组元素键值并重新排序
  2. html+引导,html – 引导点的CSS样式
  3. 服务器中毒后老板差点把我开除了。。。
  4. Docker系列教程01-Centos7安装新版Docker教程(10步)
  5. java 查询线程池_[代码全屏查看]-我的 Java 线程池测试类
  6. CMU Database Systems - Sorting,Aggregation,Join
  7. NETGEAR拒绝连接请求_详解 Tomcat 的连接数与线程池
  8. 3.8 - Using the Print Function
  9. 基于自动图像分割算法和扩展数据集深度学习的经济作物病害识别
  10. scala 类中的对象是类_Scala类和对象– Singleton对象,伴侣类
  11. 结对作业_core组
  12. (转)OpenStack Kilo 版本中 Neutron 的新变化
  13. anaconda 卸载_Windows安装Anaconda使用教程
  14. 在CAD中容易混淆的概念
  15. 韩顺平java30天Utils包下的工具类
  16. Mac下Chrome添加.crx浏览器插件
  17. 泛化误差,偏差方差分解
  18. tas5717php手册,TAS5715 具有扬声器均衡、双频带 DRC 和 DC 保护的 25W 立体声 I2S 音频放大器...
  19. 北京航空航天大学计算机学院保研,北京航空航天大学计算机学院(专业学位)计算机技术保研...
  20. 中国数字乳房断层合成(DBT)设备市场趋势报告、技术动态创新及市场预测

热门文章

  1. 金匮要略重点整理 笔记
  2. 韩国顶级舞台剧《爱上街舞少年的芭蕾少女》掀起街舞狂潮
  3. [生存志] 第5节 第一篇 以史为鉴 罗振宇的思维逻辑
  4. 233_S32DS共性特征学习
  5. 毕设模块之一:最新版 微博网络爬虫(可登录版)
  6. C++万年历课程设计
  7. 四川企立方:拼多多标题要怎么组成
  8. Android模拟Windows10,windows10模拟器
  9. Axure RP8 下载、安装、破解、汉化一条龙服务
  10. Unity Shader 麻将平面阴影高光