ペアワイズ合計

数値解析において、ペアワイズ加算（カスケード加算とも呼ばれる）は、有限精度 浮動小数点数のシーケンスを加算する手法であり、単純に合計を順番に累積する場合と比較して、累積丸め誤差を大幅に削減します。 ^{[1]通常は丸め誤差がさらに小さい}カハン加算 などの他の手法もありますが、ペアワイズ加算は、計算コストがはるかに低く、ほぼ同じくらい優れています（対数係数のみの違い）。つまり、単純な加算とほぼ同じコスト（およびまったく同じ数の算術演算）になるように実装できます。

特に、n個の数値のシーケンスx _nのペアワイズ加算は、シーケンスを再帰的に2 つの半分に分割し、各半分を合計し、2 つの合計を加算することによって機能します。これは分割統治アルゴリズムです。最悪の場合の丸め誤差は、最大でO (ε log  n )に漸近的に増加します。ここで、 ε はマシン精度です (以下で説明するように、条件数 が固定されていると仮定)。 ^[1] これと比較して、合計を順番に累積する単純な手法 (各x _{i を}1 つずつi  = 1、...、  nに対して加算する) では、丸め誤差は最悪の場合O (ε n ) に増加します。^[1] Kahan 加算では、最悪の場合の誤差はnとは無関係に約O (ε)ですが、数倍の算術演算が必要です。^[1] 丸め誤差がランダムで、特にランダムな符号を持つ場合、それらはランダムウォークを形成し、ペアワイズ加算の誤差の増加は平均にまで減少します。^[2] ${\displaystyle O(\varepsilon {\sqrt {\log n}})}$

非常によく似た再帰的な加算構造は多くの高速フーリエ変換（FFT）アルゴリズムに見られ、それらのFFTの同じように遅い丸め累積の原因となっている。^{[2] [3]}

疑似コードでは、長さn ≥ 0の配列 xのペアワイズ合計アルゴリズムは次のように記述できます。

s = pairwise ( x [1... n ])
       if  n ≤ N     基本ケース: 十分に小さい配列の単純な合計
          s = 0
           for  i = 1 to n 
              s = s + x [ i ]
       else         分割統治: 配列の2つの半分を再帰的に合計
          m = floor ( n / 2)
           s = pairwise ( x [1... m ]) + pairwise ( x [ m +1... n ])
       end if

十分に小さいNに対して、このアルゴリズムは、基本ケースとして単純なループベースの合計に切り替わり、その誤差境界は O(Nε) です。^[4] 合計全体では、最悪の場合の誤差は、与えられた条件数に対して、大きなnに対してO (ε log  n ) として漸近的に増加します(以下を参照)。

この種のアルゴリズムでは（一般的な分割統治アルゴリズム^[5]と同様に）、再帰のオーバーヘッドを償却するために、より大きなベースケースを使用することが望ましい。N = 1 の場合、入力ごと におよそ 1 回の再帰サブルーチン呼び出しがあるが、より一般的には、再帰がちょうどn  =  Nで停止する場合、（およそ） N /2 入力 ごとに 1 回の再帰呼び出しがある。Nを十分に大きくすることで、再帰のオーバーヘッドを無視できる（まさにこの再帰加算のための大きなベースケースの手法は、高性能 FFT 実装で採用されている^[3]）。

Nに関係なく、単純な加算の場合と同じように、合計で正確にn −1 回の加算が実行されるため、再帰のオーバーヘッドが無視できるほど小さい場合、ペアワイズ加算の計算コストは​​単純な加算と本質的に同じになります。

このアイデアのバリエーションとして、各再帰段階で合計をbブロックに分割し、各ブロックを再帰的に合計してから結果を合計するというものがあります。これは提案者によって「スーパーブロック」アルゴリズムと呼ばれています。 ^[6] 上記のペアワイズアルゴリズムは、最後の段階b  =  N を除いて、すべての段階で b = 2に対応します。

Dalton、Wang、Blainey (2014) は、ペアワイズ加算のための反復的な「シフトリデュース」定式化について説明しています。これは、SIMD命令を使用して展開および高速化できます。展開されていないバージョンは次のとおりです。^[7]

double shift_reduce_sum ( double ∗ x , size_t n ) { double stack [ 64 ], v ; size_t p = 0 ; for ( size_t i = 0 ; i < n ; ++ i ) { v = x [ i ]; // シフトfor ( size_t b = 1 ; i & b ; b <<= 1 , −− p ) // 削減v += stack [ p − 1 ]; stack [ p ++ ] = v ; } double sum = 0.0 ; while ( p ) sum += stack [ −− p ]; return sum ; }     
    
     
           
                                     
                
        
      
  
     
   
      
   

i  = 1, ...,  nのn個の値x _{i を}合計するとします。正確な合計は次のようになります。

${\displaystyle S_{n}=\sum _{i=1}^{n}x_{i}}$

(無限の精度で計算されます)。

基本ケースN  = 1のペアワイズ和では、代わりに が得られ、誤差は次のように制限される: ^[1]${\displaystyle S_{n}+E_{n}}$${\displaystyle E_{n}}$

${\displaystyle |E_{n}|\leq {\frac {\varepsilon \log _{2}n}{1-\varepsilon \log _{2}n}}\sum _{i=1}^{n}|x_{i}|}$

ここで、ε は使用される演算の機械精度です（たとえば、標準の倍精度浮動小数点数の場合、ε ≈ 10 ⁻¹⁶）。通常、関心のある量は相対誤差であり、したがって、上側では次の式で制限されます。 ${\displaystyle |E_{n}|/|S_{n}|}$

${\displaystyle {\frac {|E_{n}|}{|S_{n}|}}\leq {\frac {\varepsilon \log _{2}n}{1-\varepsilon \log _{2}n}}\left({\frac {\sum _{i=1}^{n}|x_{i}|}{\left|\sum _{i=1}^{n}x_{i}\right|}}\right).}$

相対誤差境界の式では、分数 (Σ| x _i |/|Σ x _i |) が総和問題の条件数です。基本的に、条件数は、計算方法に関係なく、誤差に対する総和問題の固有の感度を表します。 ^[8]固定精度の固定 アルゴリズム(つまり、任意精度の演算を使用するものや、データに基づいてメモリと時間の要件が変化するアルゴリズムではないもの)によるすべての(後方安定な) 総和方法の相対誤差境界は、この条件数に比例します。 ^[1]悪条件の総和 問題とは、この比率が大きい問題であり、この場合はペアワイズの総和でも相対誤差が大きくなる可能性があります。たとえば、加数x _iが平均がゼロの相関のない乱数である場合、合計はランダムウォークであり、条件数は に比例して増加します。一方、平均がゼロでないランダム入力の場合、条件数は として有限定数に漸近します。入力がすべて非負の場合、条件数は 1 になります。 ${\displaystyle {\sqrt {n}}}$${\displaystyle n\to \infty }$

なお、分母は実際には実質的に 1 です。これは、 n が2 ^{1/ε の}オーダー(倍精度で はおよそ 10 ^{10 15 )}になるまでは 1 よりはるかに小さいためです。${\displaystyle 1-\varepsilon \log _{2}n}$${\displaystyle \varepsilon \log _{2}n}$

比較すると、単純な加算（単純に数字を順番に加算し、各ステップで四捨五入する）の相対誤差境界は、条件数倍になるほど大きくなります。^[1] 実際には、四捨五入誤差は平均がゼロでランダムな符号を持つ可能性の方がはるかに高く、ランダムウォークを形成します。この場合、単純な加算の二乗平均平方根相対誤差は平均して のように大きくなり、ペアワイズ加算の誤差は平均して のように大きくなります。^[2]${\displaystyle O(\varepsilon n)}$${\displaystyle O(\varepsilon {\sqrt {n}})}$${\displaystyle O(\varepsilon {\sqrt {\log n}})}$

ペアワイズ加算はNumPy ^[9]と技術計算言語Julia ^[10]のデフォルトの加算アルゴリズムであり、どちらの場合も単純な加算と同等の速度であることが確認されている（大きな基本ケースを使用しているため）。

その他のソフトウェア実装としては、 C#言語用のHPCsharpライブラリ^[11]やD言語の標準ライブラリsummation ^[12]などがある。