アルゴリズム自体はすでに書いているから、再掲はしない。ただ、並列化するという観点で、もう一度誤差拡散法について、考えてみる。

誤差拡散法では、「現在処理を行っているピクセル」により生じた誤差を、右隣・右下・下・左下のセルに分配する。

その都合上、「現在処理を行っているピクセル」の処理が終了することにより、初めて右隣のピクセルが処理できるようになる。図で表すと下記のようになる。

確定確定確定確定確定

確定確定確定処理中未定

処理可能未定未定未定未定

未定未定未定未定未定

上図の「処理中」と書かれている赤いセルが、「現在処理を行っているピクセル」である。

上図から判るとおり、誤差拡散法は左上から順番に処理を進めなければならない。画像の途中にあるピクセルから突然処理を開始することはできないのだ。

ところが、ここに一点だけ、「処理中」のピクセルとは独立に処理可能なピクセルがある。それは、「処理中」のピクセルがある行の、つぎの行の左隅のピクセルだ。上図で言うと「処理可能」と書かれているピクセルである。

このピクセルには、もうこれ以上左隣のピクセルが存在しないため、直上のピクセルの処理が終われば処理可能となるのである。

並列化

ここでおもむろに、左隅にある「処理可能」なピクセルから右に向かって、可能な限り処理を進めてみることにしてみる。

確定確定確定確定確定確定

確定確定確定確定処理中未定

確定確定処理中未定未定未定

処理可能未定未定未定未定未定

すると、上記のようになる。

先ほどと同様に、左隅のピクセルは「処理可能」となる。また、中程にある「処理中」のピクセルの処理は、これ以上先には進めない。それは、その上にある「処理中」のピクセルが先に進まなければ、その間に挟まれた「未定」のピクセルにかかる誤差が確定しないからだ。

ここで更に、同じ要領で、「処理可能」なピクセルから右に向かって処理を進めてみる。

確定確定確定確定確定確定確定確定

確定確定確定確定確定確定処理中未定

確定確定確定確定処理中未定未定未定

確定確定処理中未定未定未定未定未定

処理可能未定未定未定未定未定未定未定

先ほどと同じである。

でもって、上記の「処理中」のセルは全て並列で処理を行うことが可能である。

まとめて言うなら、ある行の処理を追いかけるようにして、次の行の処理を並列に処理することができる、ということだろうか。

同期

上でも少し書いたが、下側にある行の処理は、上側の行の処理を追い越すことはできない。つまり、「追い越さない」同期処理が必要である。

ところでこの「追い越さない」同期処理とは、どんなものだろうか。もちろん、セマフォとか何とかを駆使すれば、そういった同期制御を実現することも可能である。

だが、無理矢理にそんな処理を作り込んでしまったら、処理速度が遅くなってしまう。それでは並列化の意味が無くなる。だから、もっと自然な、より高速な同期方法を考えなければならない。

ではどうするのか。

こういった場合は、フェーズ化のアプローチが有効である。フェーズ化とは、各フェーズは同期なしで並列に処理を行い、フェーズの最後で全ての処理の終了を待ち合わせるというやり方だ。

同期処理は専門用語でバリアーと言われるらしい。各処理は、バリアーを抜けると一斉に開始される。そして、めいめいに処理を実行し、次のバリアーで、全ての処理（大概はスレッド）の待ち合わせを行う。でもって、全てのスレッドで、当該のフェーズの処理が終了したら、一斉にバリアーを抜けて次のフェーズに進む。

これを、誤差拡散法に適用すると、下記のような順序で実行することになる。

1 2 3 4 5 6 7 8

3 4 5 6 7 8 9 10

5 6 7 8 9 10 11 12

7 8 9 10 11 12 13 14

9 10 11 12 13 14 15 16

数字はフェーズを表している。同じ番号のピクセルは同時に処理を行う。これにより、「追い越さない」同期処理という不自然な事をやらなくても済む。

ブロック化

ところで、上記の処理を1ピクセル単位で実行したのでは、やはり処理が重くなってしまう。同期処理が多すぎるのだ。

だから何とかして、一つのフェーズで処理するピクセルの量を増やしてあげなければならない。1フェーズで多くのピクセルを処理できれば、相対的に同期処理の時間を短くすることができる。

ではどうするのか。

ここも同期処理のセオリー通り、処理をブロック化してあげればいい。

A A A A B B B B

A A A A B B B B

A A A A B B B B

A A A A B B B B

C C C C D D D D

C C C C D D D D

C C C C D D D D

C C C C D D D D

同じアルファベットのセルは、1つのフェーズでまとめて処理を行う。各ブロックの中は、シングルスレッドによる今まで通りの方法で処理を行う。

ブロック間での誤差は、下図のように隣のブロックに引き継ぐ。

X1 X2 X3 X4 A1 A2 A3 A4 B1 B2 B3 B4

X5 X6 X7 X8 A5 A6 A7 A8 B5 B6 B7 B8

X9 X10 X11 X12 A9 A10 A11 A12 B9 B10 B11 B12

X13 X14 X15 X16 A13 A14 A15 A16 B13 B14 B15 B16

Y1 Y2 Y3 Y4 C1 C2 C3 C4 D1 D2 D3 D4

Y5 Y6 Y7 Y8 C5 C6 C7 C8 D5 D6 D7 D8

赤い部分の誤差が、青い部分に伝播される。すなわち、B1には、A4の右隣の誤差が、B5には、A4の右下の誤差と、A8の右隣の誤差が伝播される。他のセルも同じ要領である。

実装

理論的には以上の考察だけで、誤差拡散法の並列化は可能である。ただ、実装するとなるとひどくやっかいである。フェーズ化を実施することにより、明示的な排他制御を行う必要はないのだが、それでも、何となく判りにくいのだ。

その面倒な処理を、簡単に書くと下記のようになる。
for ( i = 0; i < 横方向のブロック数; i++ ) {
  for ( j = 0; j < 縦方向のブロック数; j++ ) {
    int bx = ・・・;  // 処理するブロックを決定
    int by = ・・・;
    // ブロックを処理する
    for ( k = 0; k < ブロックの大きさ; k++ ) {
      for ( l = 0; l < ブロックの大きさ; l++ ) {
        // ブロック内のピクセルを処理する
      }
    }
    // ブロック間の誤差の引き継ぎ処理を行う
  }
}
図で書くと、下記のようなイメージである。

赤いブロックは、kとlのループによって処理される、単一のブロックである。

緑の矢印は、jのループによって各ブロックが処理されてゆく順番を表している。

青い矢印は、iのループにより、緑の矢印全体が一つずつ右にずれて行くことを表している。

このように処理を行う場合、jのループにより処理される各ブロックは、それぞれ並列に処理することが可能となる。すなわち、jのループはOpenMPで分割することが可能になることとなる。

また、ブロック間の誤差の引き継ぎに利用するバッファも、複数スレッドで同一のメモリ領域にアクセスする事はなくなるため、明示的に排他制御を行う必要はない。

後は、誤差の引き継ぎさえうまく実装することができれば、理論上はシングルスレッドで処理したときと、全く同じ結果を得ることができるようになるはずである。まぁ、それがなかなかうまくいかないのだが。