首页 > 解决方案 > C++中的卡汉算法

问题描述

我已经编写了以下 kahan 算法;n 从 100 万到 1000 万是合理的,但是当 n 为 1 亿时会产生非常大的误差。我对其进行了多次测试,但我无法弄清楚为什么会发生这种情况。

#include <iostream>
#include <vector>
#include <random>
#include <cfloat>
#include <iomanip>
using namespace std;
float KahanAlgorithm(const vector<float> &myarray){
    float sum{0.0f};
    float ac{0.0f};

    for(unsigned int i=0; i<myarray.size();i++){
        float temp{sum+myarray[i]};
        if(sum>=myarray[i]){
            ac=ac+(sum-temp)+myarray[i];
        }
        else{
            ac=ac+(myarray[i]-temp)+sum;
        }
        sum=temp;
    }

    return sum+ac;
}
int main() {

    int n=100000000;
    random_device r;
    default_random_engine g(r());
    uniform_real_distribution<float> d(0.f, nextafter(1.f, DBL_MAX));

    vector<float> a(n);
    vector<double> b(n);

    for(auto i=0;i<n;i++){
        a[i]=d(g);
        b[i]=static_cast<double> (a[i]);

    }

    double exact_sum;
    float kahan_sum;

    exact_sum=accumulate(b.begin(),b.end(),0.0);
    cout<<"exact "<<exact_sum<<endl;
    kahan_sum=KahanAlgorithm(a);
    cout<<" Kahan sum "<<kahan_sum<<endl;
    return 0;
}

总和是:

exact 5.00045e+07
Kahan sum 3.35544e+07

标签: c++randomdoubleroundingprecision

解决方案


我做了一个快速测试,使用了我前段时间写的 Kahan summation 的实现,并将它与你的进行了比较:

#include <cfloat>
#include <iomanip>
#include <iostream>
#include <iterator>
#include <numeric>
#include <random>
#include <vector>

namespace Kahan {
template <class InIt>
typename std::iterator_traits<InIt>::value_type accumulate(InIt begin, InIt end) {
    typedef typename std::iterator_traits<InIt>::value_type real;
    real sum = real();
    real running_error = real();

    for ( ; begin != end; ++begin) {
        real difference = *begin - running_error;
        real temp = sum + difference;
        running_error = (temp - sum) - difference;
        sum = temp;
    }
    return sum;
}
}

using namespace std;
float KahanAlgorithm(const vector<float> &myarray){
    float sum{0.0f};
    float ac{0.0f};

    for(unsigned int i=0; i<myarray.size();i++){
        float temp{sum+myarray[i]};
        if(sum>=myarray[i]){
            ac=ac+(sum-temp)+myarray[i];
        }
        else{
            ac=ac+(myarray[i]-temp)+sum;
        }
        sum=temp;
    }

    return sum+ac;
}
int main() {

    int n=100000000;
    random_device r;
    default_random_engine g(r());
    uniform_real_distribution<float> d(0.f, nextafter(1.f, DBL_MAX));

    vector<float> a(n);
    vector<double> b(n);

    for(auto i=0;i<n;i++){
        a[i]=d(g);
        b[i]=static_cast<double> (a[i]);

    }

    double exact_sum;
    float kahan_sum;

    exact_sum=accumulate(b.begin(),b.end(),0.0);
    cout<<"exact "<<exact_sum<<endl;
    kahan_sum=KahanAlgorithm(a);
    cout<<" Kahan sum "<<kahan_sum<<endl;
    float jerry = Kahan::accumulate(a.begin(), a.end());
    cout << "Jerry's implementation of Kahan: " << jerry << "\n";
    return 0;
}

样本结果:

exact 4.99944e+07
 Kahan sum 3.35544e+07
Jerry's implementation of Kahan: 4.99944e+07

推荐阅读