10

線形相補性問題を解決するためのProjectedGauss-SeidelアルゴリズムのC＃実装を探しています。これまでのところ、BulletライブラリでC ++で記述されたものを見つけましたが、残念ながら高度に最適化されています（したがって、C＃に変換するのは難しいでしょう）。

同様の質問で、.NETの数値ライブラリを調べることを提案しました。それらはすべて、連立一次方程式を解くためのアルゴリズムのみを含んでいます。

編集：私が見つけたとしても、それは完全ではないように見えるので、質問はまだ開いています。

4

2 に答える 2

11

プロジェクションなしで Gaus Seidel を実装しました。射影されたガウス ザイデルの場合、解を下限と上限内に射影 (またはクランプ) する必要があります。

public static double[] Solve (double[,] matrix, double[] right,
                              double relaxation, int iterations, 
                              double[] lo, double[] hi)
{
    // Validation omitted
    var x = right;
    double delta;

    // Gauss-Seidel with Successive OverRelaxation Solver
    for (int k = 0; k < iterations; ++k) {
        for (int i = 0; i < right.Length; ++i) {
            delta = 0.0f;

            for (int j = 0; j < i; ++j)
                delta += matrix [i, j] * x [j];
            for (int j = i + 1; j < right.Length; ++j)
                delta += matrix [i, j] * x [j];

            delta = (right [i] - delta) / matrix [i, i];
            x [i] += relaxation * (delta - x [i]);
    // Project the solution within the lower and higher limits
            if (x[i]<lo[i])
                x[i]=lo[i];
            if (x[i]>hi[i])
                x[i]=hi[i];
        }
    }
    return x;
}

軽微な修正です。Bullet Physics Library から A 行列と b ベクトルを抽出し、射影された Gauss Seidel を使用してそれを解く方法を示す要点は次のとおりです: https://gist.github.com/erwincoumans/6666160

于 2013-09-23T03:32:39.737 に答える
10

1 週間の検索の後、私はついにこの出版物を見つけました (ロシア語で、Kenny Erleben の作品に基づいています)。そこでは、予測された Gauss-Seidel アルゴリズムが説明されており、SORと終了条件で拡張されています。この基本的な C# 実装に使用した C++ の例をすべて示します。

public static double[] Solve (double[,] matrix, double[] right,
                              double relaxation, int iterations)
{
    // Validation omitted
    var x = right;
    double delta;

    // Gauss-Seidel with Successive OverRelaxation Solver
    for (int k = 0; k < iterations; ++k) {
        for (int i = 0; i < right.Length; ++i) {
            delta = 0.0f;

            for (int j = 0; j < i; ++j)
                delta += matrix [i, j] * x [j];
            for (int j = i + 1; j < right.Length; ++j)
                delta += matrix [i, j] * x [j];

            delta = (right [i] - delta) / matrix [i, i];
            x [i] += relaxation * (delta - x [i]);
        }
    }

    return x;
}
于 2012-08-02T09:53:51.290 に答える