Jaki jest wzór na interpolację kwadratową?
Z dwóch wartości mamy liniową, zatem kwadratowa wymaga trzech wartości.
no, może trochę inaczej - chodzi mi o takie coś:
mamy trzy punkty:
y1, y2 i y3, dla x1, x2 i x3
i teraz chcę wyliczyć x0, znaczy punkt w którym y = y0.
x0 = ?
W czym problem? Rozwiązujesz układ trzech równań liniowych z trzema niewiadomymi (współczynnikami funkcji kwadratowej).
bogdans napisał(a):
W czym problem? Rozwiązujesz układ trzech równań liniowych z trzema niewiadomymi (współczynnikami funkcji kwadratowej).
Pewnie tak można.
Wyznaczamy wielomian: p(x) = ax2 + bx + c;
czyli te trzy stałe: a,b,c;
a potem rozwiązujemy równanie kwadratowe: y(x1) = y1;
Tylko że to jest chyba jakieś prymitywne - są pewnie na to lepsze metody..
zwłaszcza że mogę potem użyć lepszej interpolacji: 3 - kubicznej, 4, 5...
a wtedy to już tak łatwo nie rozwiążę tego równania: p(x) = y1, dla wielomianu np. stopnia 7.
bogdans napisał(a):
Zajrzyj tu https://en.wikipedia.org/wiki/Polynomial_interpolation.
Coś w pobliżu.
Nie ma tam o mim przypadku, znaczy o interpolacji odwrotnej.
Ja muszę wyliczyć x dla znanego y, a standardowa interpolacja pozwala wyliczać y dla danego x.
W zasadzie można to wprost odwrócić, tz. zamieniamy x z y, a wtedy interpolujemy funkcję odwrotną:
y(x) = x => x = f(y), gdzie f jest funkcją odwrotną y.
Obawiam się jednak że to nie zadziała w moim przypadku.
Przykładowo. mając dane:
y x
2 1
1 2
1 3
i w przypadku y(x) jest to prosta parabola, gdzie minimum leży pomiędzy x = 2 i 3,
natomiast funkcja odwrotna x(y) byłaby niejednoznaczna.
a takich raczej nie da rady interpolować.
byłby to schemat typu;
x y (po odwróceniu)
1 2
2 2 i 3 naraz
3 2 i 3
vpiotr napisał(a):
https://4programmers.net/Forum/Newbie/173006-programowe_rozwiazywanie_ukladow_rownan_kwadratowych
raczej prymitywne to jest.
dawno temu robiłem ogólną metodę do rozwiązywania wielomianów.
Mogę wam podać gotowca - metoda Laguera (na zespolonych):
double TSolvW::laguer(double *wsp, int n) // wsp - współczynniki wielomianu rzędu n
{
poly012(&x, wsp, n);
// ny/(y' ± sqrt(H)); H=(n-1)[(n-1)y'^2 - nyy'']
complexd ny(y), h(yp);
ny *= n;
h.sqr(); h *= n-1; // h = y'^2 * (n-1)
yb *= ny; // nyy''
h -= yb; h *= n-1;
h.sqrt();
d = yp;
double c = d.re*h.re + d.im*h.im; // cos(a,b) > 0 -> |a+b| > |a-b|
if( c > 0 ) d += h; else d -= h;
if( d.Zero() ) // gdy: y' = y'' = 0
{
if( y.Zero() ) ;
else d.re = -1;
}
else d = ny /= d;
x -= d;
return d.nt();
}
Zawsze jest dokładnie n rozwiązań wielomianu stopnia n.
Przykładowe wyniki:
x2 + x + 1 = 0, czyli współczynniki: 1, 1, 1
rozwiązania:
-0.5 + 0.866025403784439i
-0.5 - 0.866025403784439i
teraz bierzemy coś z rzeczywistych, np.:
x^2 - x - 1 = 0, współczynniki: 1, -1, -1
rozw.:
-0.618033988749895
+0.618033988749895
weźmy coś wyższego, np. wielomian stopnia 7;
x7 + x6 + x5 + x4 + x3 + x2 + x + 1 = 0,, czyli współczynniki: 1, 1, 1 ...
rozwiązania:
0 + i
0 - i
-0.707106781186548 - 0.707106781186548i
-0.707106781186548 + 0.707106781186548i
-1
0.707106781186548 - 0.707106781186548i
0.707106781186548 + 0.707106781186548i
siedem sztuk;
ta moja funkcja co zwraca?
double TSolvW::laguer(double *wsp, int n)
nie pamiętam... jakaś tam tylko pomocniczą rzecz..
to jest tylko jeden krok algorytmu Laguera...
cały kod jest trochę dłuższy - deflacja, i inne rzeczy.
Mogę podać cały kod jak chcesz...
class TSolvW // rozywiązywasz wielomianów
{
public: TSolvW(int n=128) { wsp_a = 0; reset(n); }
~TSolvW() { delete []wsp_a; }
void reset(int n);
int solvone();
int solv();
double laguer(double *wsp, int n);
double newtonW();
int deflate(int m);
int nA, nb, nz, maxW;
int limitI;
double beps;
double *wsp_a; // wsp. początkowe wielomianu
double *wsp_b; // wsp. po deflacji
complexd *zera; // znalezione zera
complexd x,y,yp,yb;
complexd d; // ostatni przyrost: x' = x + d
};
To pozostaje interpolacja kwadratowa (wyznaczenie współczynników wielomianu) a potem rozwiązanie równania drugiego stopnia. Nadal nie widzę w czym tkwi problem.
bogdans napisał(a):
To pozostaje interpolacja kwadratowa (wyznaczenie współczynników wielomianu) a potem rozwiązanie równania drugiego stopnia. Nadal nie widzę w czym tkwi problem.
Jest to chyba zbyt kosztowne rozwiązanie, ponieważ mi nie potrzeba żadnego wielomianu, a tylko jedną wartość.
Na pewno istnieje prostszy sposób...
Zadowalałoby mnie coś takiego:
robimy najpierw interpolację liniową,
a dopiero potem ją korygujemy odpowiednio, uzyskując kwadratową.
A może coś z krzywych beziera?
Liniowa interpolacja wygląda tam tak:
P1(A,B) = At + B(1-t);
a kwadratowa jest tu znowu kombinacją linową - dwóch liniowych:
A, B, C - zadane punkty, a wtedy:
P2(A,B,C) = P1(A,B)t + P1(B,C)(1-t);
co można sobie podstawić i wyliczyć wprost, i chyba takie coś wyjdzie:
P2(A,B,C) = (1-t2)A + 2tB + t2*C; t = <0.1>
Ale teraz to chyba mam dwie (albo i trzy) funkcje kwadratowe: jedna dla x, druga dla y i dla z.