高速で安定したx * tanh（log1pexp（x））計算
https://cs.stackexchange.com/questions/125002
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
Russian質問
$$ f（x）= x \ tanh（\ log（1 + e ^ x））$$
関数（Mishの起動）は、正確な損失を大きく失うことなく安定したLOG1PEXPを使用して簡単に実装できます。残念ながら、これは計算的に重いです。
より速い方の数値安定した実装を書くことは可能ですか？
x * std::tanh(std::log1p(std::exp(x)))と同じくらい正確さはいいでしょう。厳密な制約はありませんが、ニューラルネットワークでの使用に合理的に正確であるべきです。
入力の分布は $ [ - \ infty、\ infty] $ からです。それは至る所で働くべきです。
解決
特定のopポイント正確さ仕様のためのmish起動関数の実装だから私はこれを最初に特徴付ける必要がありました。その実装は単精度（float）を使用し、正の半平面で安定して正確です。負のハーフプレーンでは、logfの代わりにlog1pfを使用するため、相対エラーが $ X \ To- \ Infty $ を早く成長させます。正確さの損失は $ - 1 $ で、 $ - 16.635324 $ に該当する
$ \ mathrm {tahn} $ を排除する単純な数学的変換を使用することで、同じ精度と動作を達成することができ、GPUが通常融合増加を提供することを考えると-add（FMA）と高速の逆数、どちらを利用したいと思うでしょう。例示的なCUDAコードは次のように見えます：
__device__ float my_mishf (float x)
{
float r;
float e = expf (x);
r = 1.0f / fmaf (fmaf (-0.5f, e, -1.0f), e, -1.0f);
r = fmaf (r, x, x);
return r;
}
OPが指す基準実装と同様に、これは正の半平面内で優れた精度を有し、負の半平面誤差では $ - 16.635324 $ 実装は誤って $ 0 $ を返します。
これらの精度の問題に対処したいという願望がある場合は、以下の観察を適用できます。十分に小さい $ x $ 、 $ f（x）= x \ esp（x）$ 浮動小数点精度floatの計算の場合、これは $ x <-15 $ です。間隔 $ [ - 15、-1] $ [ - 15、-1] $ については、Rational近似 $ r（x）$を使用できます。 / span> $ f（x）：= r（x）x \ exp（x）$ 。例示的なCUDAコードは次のように見えます：
__device__ float my_mishf (float x)
{
float r;
if (x >= -1.0f) {
float e = expf (x);
r = 1.0f / fmaf (fmaf (-0.5f, e, -1.0f), e, -1.0f);
r = fmaf (r, x, x);
} else {
float eh = expf (0.5f * x);
float p = 1.03628484e-3f; // 0x1.0fa7e6p-10
p = fmaf (p, x, -7.28869531e-3f); // -0x1.ddac04p-8
p = fmaf (p, x, 3.47027816e-2f); // 0x1.1c4902p-5
p = fmaf (p, x, -3.54762226e-1f); // -0x1.6b46cap-2
p = fmaf (p, x, 8.58785570e-1f); // 0x1.b7b2bep-1
p = fmaf (p, x, -1.38065982e+0f); // -0x1.6172ecp+0
p = fmaf (p, x, 5.97694337e-1f); // 0x1.3204fep-1
float q = 1.03527203e-3f; // 0x1.0f63eep-10
q = fmaf (q, x, -7.35638570e-3f); // -0x1.e21bacp-8
q = fmaf (q, x, 3.28683928e-2f); // 0x1.0d4204p-5
q = fmaf (q, x, -3.79927397e-1f); // -0x1.850bb0p-2
q = fmaf (q, x, 6.86127126e-1f); // 0x1.5f4c0ep-1
q = fmaf (q, x, -1.81509292e+0f); // -0x1.d0a9eep+0
q = fmaf (q, x, 1.00000000e+0f); // 0x1.000000p+0
r = (1.0f / q) * p;
if (x < -15.0f) r = 1.0f;
r = r * x * eh * eh;
}
return r;
}
残念ながら、この正確な解決策は性能の著しい低下のコストで達成されます。スムーズに崩壊している間に、滑らかに衰弱している左尾を達成しながら、 $ f（x）\ \ exp（x）\ \ \ exp（x）\ \ \ exp（x）\ \ \ exp（x）\ \ \ exp（x）\ \ \ exp（x）\ \ f. / SPAN>、パフォーマンスの多くを回復します。
__device__ float my_mishf (float x)
{
float r;
float e = expf (x);
if (x >= -6.0625f) {
r = 1.0f / fmaf (fmaf (-0.5f, e, -1.0f), e, -1.0f);
r = fmaf (r, x, x);
} else {
r = fmaf (-0.5f, e, 1.0f);
r = r * x * e;
}
return r;
}
マシン固有の性能向上として、expf()は、Device InterSic __expf()に置き換えられます。
他のヒント
代数的操作（@ ORLPの答えで指摘したように）では、次のことを推論できます。
$$ f（x）= x \ tanh（\ log（1 + E ^ x））\ tag {1} $$ $$= x \ frac {（1 + e ^ x）^ 2 - 1} {（1 + e ^ x）^ 2 + 1}= x \ frac {e ^ { 2x} + 2e ^ x} {e ^ {2x} + 2 ^ x + 2} \タグ{2} $$ $$= x - \ frac {2x} {（1 + e ^ x）^ 2 + 1} \ tag {3} $$
式 $（3）$ は $ x $ が否定的に否定されている精度。式 $（2）$ は、 $ x $ の大きな値には適していません。分子と分母の両方で爆発します。
関数 $（1）$ is hysptotally $ x \ to- \ infty $ としてゼロを押す。 。 $ x $ が大きくなると、式 $（3）$ は壊滅的なキャンセルに苦しみます。 ：2つの大きな用語は、本当に少数を与えるために互いにキャンセルします。式 $（2）$ はこの範囲でより適しています。
これは、 $ - 18 $ にかなりうまく機能します。
関数を近づけてみましょう $ f（x）$ を $ x \に近似しようとしましょう。 - \ infty $ 。
$$ f（x）= x \ frac {e ^ {2x} + 2e ^ x} {e ^ {2x} + 2e ^ x + 2} $$ < / SPAN>
$ e ^ {2x} $ は、 $ e ^ x $ 。 $ E ^ x $ は、 $ 1 $ より小さい桁数です。これら2つの事実を使用して、 $ f（x）$ x に近似することができます：
$ f（x）\ \ \ frac {e ^ x} {e ^ x + 1} \約xe ^ x $
結果：
$ f（x）\ begin {ケース} xe ^ x、＆\ text {$ x \ le -18 $の場合} X \ FRAC {E ^ {2x} + 2E ^ x} {E ^ {2x} + 2 ^ x + 2}＆\ text {$ -18 \ lt x \ \ le -0.6 $の場合} X - \ frac {2x} {（1 + E ^ x）^ 2 + 1}、＆\ text {それ以外の場合} \ end {ケース} $
高速CUDA実装：
__device__ float mish(float x)
{
auto e = __expf(x);
if (x <= -18.0f)
return x * e;
auto n = e * e + 2 * e;
if (x <= -0.6f)
return x * __fdividef(n, n + 2);
return x - 2 * __fdividef(x, n + 2);
}
編集：
さらにより速く正確なバージョン：
$ f（x）\ begin {ケース} X \ FRAC {E ^ {2x} + 2E ^ x} {E ^ {2x} + 2 ^ x + 2}＆\ text {$ x \ le -0.6 $} \\ X - \ frac {2x} {（1 + E ^ x）^ 2 + 1}、＆\ text {それ以外の場合} \ end {ケース} $
__device__ float mish(float x)
{
auto e = __expf(value);
auto n = e * e + 2 * e;
if (value <= -0.6f)
return value * __fdividef(n, n + 2);
return value - 2 * __fdividef(value, n + 2);
}
対数を実行する必要はありません。 $ p= 1以降の\ exp（x）$ を照らした場合 $ f（x）= x \ cdot \ dfrac {p ^ 2-1} {p ^ 2 + 1} $ または代わりに $ f（x）= x - \ dfrac {2x} {p ^ 2 + 1} $ 。
私の印象は、誰かが0から1までのスムーズに行く関数f（x）によってxを掛けたいということです。これを行った基本関数を使用して式を使用して、関数の選択の範囲内で数学的な理由がないことがわかりました。 。
パラメータtを選択した後、 $ p_t（x）= 1/2 +（3/4t）x - x ^ 3 /（4t ^ 3）$ 、 $ p_t（0）= 1/2 $ 、 $ p_t（t）= 1 $ 、 $ p_t（-t）= 0 $ 、および $ p_t '（t）= p_t'（ - t）= 0 $ 。 x <-t、1 if x> + 1の場合はg（x）= 0、および $ p_t（x）$ の場合、-t≤x≤+ tの場合これは、0から1に円滑に変化する関数です。別のパラメータsを選択し、f（x）の代わりにx * g（x - s）を計算します。
t= 3.0、s= -0.3が与えられた関数と非常に合理的に一致し、高度に速い（これは重要と思われる）。もちろん違いです。この関数はいくつかの問題のツールとして使用されるので、元の関数がより良いの数学的な理由を見たいと思うでしょう。
コンテキストここは、コンピュータビジョンとニューラルネットの訓練のための活性化機能です。
チャンスはGPUで実行される予定です。パフォーマンスは一般的な入力の分布に依存しているが、一般的に言えば、GPUコードの分岐を回避することが重要です。ワープの発散はあなたのコードのパフォーマンスを大幅に低下させる可能性があります。たとえば、 cuda toolkitのドキュメント< / a>言う：
注：優先順位は高い：同じ縦語内の異なる実行パスを避けます。 フロー制御命令（if、switch、do、の場合）は、同じWarpのスレッドを発散させることで、命令スループットに大きな影響を与える可能性があります。つまり、異なる実行パスに従うことです。これが起こると、異なる実行パスは別々に実行する必要があります。これにより、このWARPに対して実行された命令の総数が増えます。 ... ほんの数本の指示を含む枝の場合、ワープの発散は一般的に限界的な性能損失をもたらします。たとえば、コンパイラは実際のブランチを回避するために予測を使用することができます。代わりに、すべての命令がスケジュールされていますが、スレッドごとの条件コードまたは述語コントロールは、指示を実行するスレッドを実行します。誤った述語を持つスレッドは結果を書き込まず、またアドレスやオペランドの読み取りを評価しないでください。
2つのブランチフリーの実装
OPの答えは短い枝を持ちますので、分岐予測はいくつかのコンパイラで発生する可能性があります。私が気付いたことのもう1つのことは、呼び出しごとに指数関数を計算することが許容できると思われることです。つまり、指数関数への呼び出しが「高価な」または「遅い」ではないと言って、OPの答えを理解しています。
その場合は、次の簡単なコードをお勧めします。
__device__ float mish(float x)
{
float expx = __expf(x);
return x / (1.0f + 2.0f / (expx * (2.0f + expx)));
}
ブランチ、1つの指数、1つの乗算、および2つの部門はありません。部門はしばしば乗算よりも高価ですので、このコードを試してみました：
__device__ float mish(float x)
{
float expx = __expf(x);
float psi = expx * (2.0f + expx);
return x * (psi / (2.0f + psi));
}
は、分岐はありません、1つの指数関数的、2つの乗算、および1つの分割です。
相対誤差
私はこれら2つの実装とOPの答えのlog10の相対的な精度を計算しました。私は、1/1024の増分で間隔（-100,100）を超えて計算し、次に51以上の値を超える最大値を計算した（視覚的な雑然としているが、それでも正しい印象を与える）。倍精度で最初の実装を計算するのは、参照として十分です。指数関数は1つのULP内で正確であり、一握りの算術演算はあります。残りのビットは、テーブルメーカーのジレンマを非常に低いものにするのに十分な以上のものです。したがって、正しく丸みを帯びた単精度基準値を計算できる可能性が非常に高いです。
緑：最初の実装赤：第二の実装。青：opの実装青と赤はそれらの範囲のほとんどの重なり（約-20の左）。
OPに注意してください。フル精度を維持したい場合は、カットオフを-5より大きく変更します。
性能
これら2つの実装をテストしてどちらが速いかを確認する必要があります。彼らは少なくともOPと同じくらい速くされるべきです、そして、私は彼らが枝が不足しているためにはるかに速くなると思われます。しかし、彼らがあなたにとって十分な速さではないならば、あなたはもっとあなたができることがあります。
重要な質問：
あなたが見た典型的な入力値の分布は何ですか？関数は効果的に計算可能な範囲全体にわたって均一に分布する値です。それともほとんど常に0のままになっているのですか？もしそうなら、どの分散/スプレッドで？
漸近性を向上させることができる。
左側のOPは、CUTOFFを持つx * expxを-18に使用します。このカットオフは、精度の損失なしで約-15.5625に増加することができます。追加の乗算のコストでは、x * expx * (1.0f - 0.5f * expx)と約-4.875のカットオフを使用できます。注：0.5の乗算は、指数から1の減算に最適化することができるので、ここではカウントされていない。
右側には、別の漸近的なものを紹介することができます。 x > 8.75、単にreturn xの場合もう少しコストで、x * (1.0f - 2.0f * __expf(-2.0f * x))の場合はx > 6.0を実行できます。
補間
範囲の中央部（-4.875,6.0）は、補間の表を使用できます。それらの範囲が等間隔に配置されている場合は、1つの分割を使用して直接インデックスをテーブルに計算できます（分岐なし）。そのようなテーブルを計算することはいくらかの努力を取りますが、あなたのニーズによっては価値があるかもしれません：一握りのマルチプリ
ESおよびADD が指数関数よりも安価である可能性があります。そうは言っても、図書館の指数関数の実装者はおそらく多くの時間と労力を費やしたことが彼らの正しいかつ迅速に費やしました。また、「Mish」機能は、範囲の縮小の機会をすべて提示していません。