2

ここに記載されているものを作成する必要があるこのコードを書きました。

ウィキペディアからの勾配の共役

しかし、いくつかの反復の後、変数denomAlphaはゼロになるため、 NAN on を取得しalphaます。では、私のアルゴリズムの何が問題なのですか?

import Jama.Matrix;

public class ConjugateGrad {

    private static final int MAX_IT = 20;
    private static final int MAX_SIZE = 50;

    public static void main(String[] args)  {
        Matrix A = Matrix.random(MAX_SIZE, MAX_SIZE);
        Matrix b = Matrix.random(MAX_SIZE, 1);

        double[][] d = new double[MAX_SIZE][1];
        for(int ii=0;ii<MAX_SIZE;ii++) {
            d[ii][0] =0;
        }

        Matrix x = Matrix.constructWithCopy(d);
        Matrix r = b.minus(A.times(x));
        Matrix p = r;
        Matrix rTrasp_r = r.transpose().times(p);

        for (int i = 0; i < MAX_IT; i++) {
            Matrix denomAlpha = p.transpose().times(A.times(p));
            double numeratorAlpha = rTrasp_r.getArray()[0][0];
            double Alpha = numeratorAlpha / denomAlpha.getArray()[0][0];
            x = x.plus(p.times(Alpha)); 
            r = r.minus(A.times(p));  
            Matrix rNew = r.transpose().times(r); 
            if (Math.sqrt(rNew.getArray()[0][0]) <1.0e-6) {
                break;
            }
            double Beta = rNew.getArray()[0][0] / rTrasp_r.getArray()[0][0];
            p = r.plus(p.times(Beta));
            rTrasp_r = rNew;           
        }
    }
}

これらのパラメータと同じです：

double[][] matrixA = {{4,1},{1,3}};
Matrix A = Matrix.constructWithCopy(matrixA);    
double[][] vectorb = {{1},{2}};
Matrix b = Matrix.constructWithCopy(vectorb);
double[][] d = {{2},{1}};
Matrix x = Matrix.constructWithCopy(d);

アルゴリズムの最初のステップでは問題はありませんが、2 番目のステップではそうではありません...

r: -8.0, -3.0
Alpha: 0.22054380664652568
Beta:   12.67123287671233
x:  0.2356495468277946, 0.33836858006042303,

第二段階 ：

Alpha: 0.0337280177221555
Beta:  159.11259655226627
x:  -2.2726985108925097, -0.47156587291133856,

OK、エラーが 1 つ見つかりました:

r = r.minus(A.times(p).times(Alpha));

今では動作します：

r: -8.0, -3.0, 
Alpha: 0.22054380664652568
rNew: 0.6403099643121183, 
Beta: 0.008771369374138607
p:  -0.3511377223647101, 0.7229306048685207, 
x:  0.2356495468277946, 0.33836858006042303,
4

1 に答える 1

0

ハックの答えで申し訳ありませんが...ウィキペディアの記事の数値例のパラメーターを使用し、各ステップで行列をターミナルに出力すると、矛盾が見つかる可能性があります。

于 2012-08-19T17:22:06.943 に答える