1

平面上の 2 点が与えられたx, y場合:

x, f(x)
1, 3
2, 5

Lagrange と find を使用してそれらを補間できf(1.5)ます4。少し考えてみると、方程式の係数を見つける方法を見つけることができました。

void l1Coefficients(const vector<double> &x, const vector<double> &y) {

    double a0 = y[0]/(x[0]-x[1]);
    double a1 = y[1]/(x[1]-x[0]);

    double b0 = (-x[1]*y[0])/(x[0]-x[1]);
    double b1 = (-x[0]*y[1])/(x[1]-x[0]);

    double a = a0 + a1;
    double b = b0 + b1;

    cout << "P1(x) = " << a << "x +" << b << endl;
}

それは私に与えます P1(x) = 2x +1。

もう少し考えて、それを拡張して2nd方程式を注文することができました。したがって、ポイントを考えると：

1, 1
2, 4
3, 9

P2(x) = 1x^2 +0x +0次の式を見つけました。

void l2Coefficients(const vector<double> &x, const vector<double> &y) {

    double a0 =              y[0] / ((x[0]-x[1])*(x[0]-x[2]));
    double a1 =              y[1] / ((x[1]-x[0])*(x[1]-x[2]));
    double a2 =              y[2] / ((x[2]-x[0])*(x[2]-x[1]));

    double b0 = -(x[1]+x[2])*y[0] / ((x[0]-x[1])*(x[0]-x[2]));
    double b1 = -(x[0]+x[2])*y[1] / ((x[1]-x[0])*(x[1]-x[2]));
    double b2 = -(x[0]+x[1])*y[2] / ((x[2]-x[0])*(x[2]-x[1]));

    double c0 =  (x[1]*x[2])*y[0] / ((x[0]-x[1])*(x[0]-x[2]));
    double c1 =  (x[0]*x[2])*y[1] / ((x[1]-x[0])*(x[1]-x[2]));
    double c2 =  (x[0]*x[1])*y[2] / ((x[2]-x[0])*(x[2]-x[1]));

    double a = a0 + a1 + a2;
    double b = b0 + b1 + b2;
    double c = c0 + c1 + c2;

    cout << "P2(x) = " << a << "x^2 +" << b << "x +" << c << endl;
}

一生懸命働いて、実際に4次までの次数の方程式の係数を見つけることができました。

n次数方程式の係数を見つけるには? どこ

Pn(x) = c_2x^2 + c_1x^1 + c_0x^0 + ...
4

2 に答える 2

1

単純な線形代数の問題です。

x k -> f(x k )という形式の N 個のサンプルのセットがあり、関数 f(x) の一般的な形式は次のとおりです。

f(x) = c 0 x 0 + c 1 x 1 + ... + c N-1 x N-1

係数 c 0 ... c N-1を見つけたいとします。それを達成するために、次の形式の N 個の方程式のシステムを構築します。

c 0 x k 0 + c 1 x k 1 + ... + c N-1 x k N-1 = f(x k )

ここで、kはサンプル番号です。x kと f(x k ) は変数ではなく定数であるため、線形連立方程式が得られます。

線形代数で表現すると、次を解決する必要があります。

Ac = b

ここで、A はxの累乗のヴァンデルモンド行列で、bは f(x k ) 値のベクトルです。

このような系を解くには、 などの線形代数ライブラリが必要ですEigen。コード例については、こちらを参照してください。

このようなアプローチで問題が発生する可能性がある唯一のことは、線形方程式系が過小決定されることです。これは、N サンプルが N-1 未満の次数の多項式で適合できる場合に発生します。このような場合でも、次のようなムーア-ペンローズ疑似逆行列を使用してこのシステムを解くことができます。

c = pinv(A)*b

残念ながら、実装Eigenはありませんが、pinv()Singular Value Decomposition (SVD) に関して自分でコーディングするのは非常に簡単です。

于 2017-04-23T11:02:45.650 に答える