0

これは私が 1 週間苦労してきた問題であり、無駄な時間を過ごした後にあきらめるために戻ってきました...

次のラゲール多項式の係数を見つけることになっています。

P0(x) = 1

P1(x) = 1 - x

Pn(x) = ((2n - 1 - x) / n) * P(n-1) - ((n - 1) / n) * P(n-2)

何らかの理由で取得した係数が大きすぎるように見えるため、実装にエラーがあると思います。これは、このプログラムが生成する出力です。

a1 = -190.234
a2 = -295.833
a3 = 378.283
a4 = -939.537
a5 = 774.861
a6 = -400.612

コードの説明 (以下):

配列を宣言する部分までコードを少し下にスクロールすると、与えられた x と y が見つかります。

関数 polynomial は、特定の x に対する前記多項式の値で配列を埋めるだけです。再帰関数です。出力値を確認したので、うまくいくと思います。

関数 gauss は、出力配列に対してガウス消去法を実行して係数を見つけます。ここから問題が始まると思います。このコードに間違いがあるのでしょうか、それとも結果を検証する私の方法が悪いのでしょうか? 私はそのようにそれらを確認しようとしています:

-190.234 * 1.5 ^ 5 - 295.833 * 1.5 ^ 4 ... - 400.612 = -3017,817625 =/= 2

コード：

#include "stdafx.h"
#include <conio.h>
#include <iostream>
#include <iomanip>
#include <math.h>

using namespace std;

double polynomial(int i, int j, double **tab)
{
    double n = i;
    double **array = tab;
    double x = array[j][0];

    if (i == 0) {
        return 1;
    } else if (i == 1) {
        return 1 - x;
    } else {
        double minusone = polynomial(i - 1, j, array);
        double minustwo = polynomial(i - 2, j, array);
        double result = (((2.0 * n) - 1 - x) / n) * minusone - ((n - 1.0) / n) * minustwo;
        return result;
    }
}

int gauss(int n, double tab[6][7], double results[7])
{
    double multiplier, divider;

    for (int m = 0; m <= n; m++)
    {
        for (int i = m + 1; i <= n; i++)
        {
            multiplier = tab[i][m]; 
            divider = tab[m][m];

            if (divider == 0) {
                return 1;
            }

            for (int j = m; j <= n; j++)
            {
                if (i == n) {
                    break;
                }

                tab[i][j] = (tab[m][j] * multiplier / divider) - tab[i][j];
            }

            for (int j = m; j <= n; j++) {
                tab[i - 1][j] = tab[i - 1][j] / divider;
            }
        }
    }

    double s = 0;
    results[n - 1] = tab[n - 1][n];
    int y = 0;
    for (int i = n-2; i >= 0; i--)
    {
        s = 0;
        y++;
        for (int x = 0; x < n; x++)
        {
            s = s + (tab[i][n - 1 - x] * results[n-(x + 1)]); 

            if (y == x + 1) { 
                break;
            }
        }
        results[i] = tab[i][n] - s;
    }

}


int _tmain(int argc, _TCHAR* argv[])
{
    int num;
    double **array;

    array = new double*[5];
    for (int i = 0; i <= 5; i++)
    {
        array[i] = new double[2];
    }
                         //i         0      1       2       3       4       5
    array[0][0] = 1.5;  //xi         1.5    2       2.5     3.5     3.8     4.1
    array[0][1] = 2;    //yi         2      5       -1      0.5     3       7
    array[1][0] = 2;
    array[1][1] = 5;
    array[2][0] = 2.5;
    array[2][1] = -1;
    array[3][0] = 3.5;
    array[3][1] = 0.5;
    array[4][0] = 3.8;
    array[4][1] = 3;
    array[5][0] = 4.1;
    array[5][1] = 7;

    double W[6][7]; //n + 1

    for (int i = 0; i <= 5; i++)
    {
        for (int j = 0; j <= 5; j++)
        {
            W[i][j] = polynomial(j, i, array);
        }
        W[i][6] = array[i][1];
    }

    for (int i = 0; i <= 5; i++)
    {
        for (int j = 0; j <= 6; j++)
        {
            cout << W[i][j] << "\t";
        }
        cout << endl;
    }

    double results[6];
    gauss(6, W, results);

    for (int i = 0; i < 6; i++) {
        cout << "a" << i + 1 << " = " << results[i] << endl;
    }

    _getch();
    return 0;
}
4

2 に答える 2

2

再帰的多項式の生成に関するあなたの解釈は、修正が必要か、私には少し巧妙すぎると思います。

指定された P[0][5] = {1,0,0,0,0,...}; P[1][5]={1,-1,0,0,0,...};
次に、P[2] は a*P[0] + convolution(P[1], { c, d }); です。
ここで、a = -((n - 1) / n) c = (2n - 1)/n および d= - 1/n

これは一般化できます: P[n] == a*P[n-2] + conv(P[n-1], { c,d }); すべてのステップで、(c + d*x) による多項式の乗算が行われます。これにより、次数が 1 つ (1 つだけ...) 増加し、スカラー a で乗算された P[n-1] に加算されます。

その場合、補間係数 x の範囲は [0..1] である可能性が最も高くなります。

(たたみ込みとは、多項式の乗算を実装する必要があることを意味します。幸いなことに、これは簡単です...)

              [a,b,c,d]
                * [e,f]
        ------------------
           af,bf,cf,df  +
        ae,be,ce,de, 0  +
--------------------------
 (= coefficients of the final polynomial)
于 2012-11-14T17:35:20.310 に答える
0

の定義はP1(x) = x - 1、記載されているとおりに実装されていません。あなたは1 - x計算にあります。

私はそれ以上見ませんでした。

于 2012-11-14T17:16:23.970 に答える