從Neville算法的符號計算多項式

我正在努力依靠計算多項式插值的Neville算法。更多關於它你可以在http://en.wikipedia.org/wiki/Neville%27s_algorithm找到。在某一點計算多項式對我來說不是問題。有很多關於它的來源。我的問題依賴於我不想在某個點計算多項式，我想以這種形式得到多項式：a_0 + a_1x + ... + a_nx^n。我不知道如何開始。你能給我一些提示嗎？從Neville算法的符號計算多項式

來源

2014-10-30 MC2DX

在Numerical Recipes中有這方面的材料，你可能想諮詢一下。然而，他們強調，計算係數是一個微妙的問題，這樣很容易就會失去很多精確度，並且越多的點越多，你越容易失去精度。你真的需要係數嗎？做什麼的？ – dmuir 2014-10-30 17:56:04

我需要係數，因爲我需要計算這個多項式的導數。 – MC2DX 2014-10-30 18:02:20

您可能想閱讀http://en.wikipedia.org/wiki/Runge%27s_phenomenon – dmuir 2014-10-31 12:21:13

Neville算法的變體允許計算一些常數（不是多項式係數），以便與另一個函數一起使用，該函數可以評估任何點處的插值多項式。後一個功能很容易區分。下面的C代碼是我使用的，我相信它工作正常。但是我懷疑你可能會對多項式插值的工作情況感到失望，除非你提供的數據確實來自多項式。如果你可以很容易地在很多點上採樣數據，那麼通過對數據進行最小二乘擬合，找到一個多項式可能會更好（表示爲chebyshev多項式的和）。

// fill C (allocated if null) with params for interpolating polynomial 
// use params with interp_poly_eval 
// !! these are NOT polynomial coefficients. 

double* interp_poly(Int deg, const double* x, const double* y, double* restrict C) 
{ 
double* c = C ?: calloc(deg+1, sizeof *c); 
Int i, j; 
    memcpy(c, y, (deg+1)*sizeof *y); 
    for (i=1; i<=deg; i++) 
    { for (j=deg; j>=i; j--) 
     { c[j] = (c[j]-c[j-1])/(x[j]-x[j-i]); 
     } 
    } 
    return c; 
} 


double interp_poly_eval(Int deg, const double* c, const double* x, double X) 
{ 
double p = c[deg]; 
Int i = deg; 
    while(--i >= 0) 
    { p = c[i] + (X-x[i])*p; 
    } 
    return p; 
} 

// as above but also returns derivative of the polynomial through pdp 
double interp_poly_eval_d(Int deg, const double* c, const double* x, double X, double* pdp) 
{ 
double p = c[deg]; 
double dp = 0.0; 
Int i = deg; 
    while(--i >= 0) 
    { dp = (X-x[i])*dp + p; 
     p = c[i] + (X-x[i])*p; 
    } 
    *pdp = dp; 
    return p; 
}

來源

2014-10-31 12:32:09 dmuir

您是否使用微分商？你可以與我們分享資料嗎？ – MC2DX 2014-10-31 19:28:45

對不起，我不明白你在問什麼。上面的第三個例程計算插值多項式的導數。 – dmuir 2014-11-03 11:25:31

您需要具有可以表示單變量（一個變量）多項式的類/數據結構。那麼很容易計算出實際的多項式用本場比賽的算法本身，即

P[0,0] = y[0] ; constant 
    P[1,1] = y[1] ; constant

，然後遞歸

P[a,b] = P[a,b-1] * x[b]/C + P[a,b-1] * -X/C + P[a+1,b] * -x[a]/C + P[a+1,b] * X/C

其中X是多項式「X（即其係數），而不是僅僅在單個點「（即一階單項式），C是常數X [b] -X [a]。

當您按照慣例在算法中進行遞歸時，這將爲您提供實際的多項式。注意，以上所有算法均爲多項式，即P [a，b-1] * x [b]/C表示多項式P [a，b-1]乘以常數x [b]/C（第b點除以C的x座標）。

如果你想得到一個確切的結果，使用任意精度理性運算包（例如C/C++的GMP）。或者，使用浮點數計算，但最終可能會影響到您的舍入誤差。

來源

2014-12-13 19:37:53

從Neville算法的符號計算多項式

回答

相關問題