interpolacja kubiczna po eipsie

0

Problem:

mam obliczyć pozycję peryhelium orbity, tz. r(t_min) = (x,y,z),
mając dane tylko dwa punkty: r0, v0, oraz r1 i v1, i wiem że minimum jest pomiędzy tymi punktami.

|r| = r_min.

Zatem biorę interpolację kubiczną, bo taka akurat tu podchodzi:
mamy tylko dwa punkty i dwie pochodne, czyli prędkości, w tych samych punktach.

y(t=0), y'(0), oraz y(1), y'(1); gdzie y jak i v to wektory: x,y,z;

Zatem obliczam ten wielomian 3-go stopnia, i tam wyjdzie coś takiego:

y(t) = at<sup>3 + bt</sup>2 + ct + d
gdzie: a, b, c i d należy wyznaczyć, co jest dość proste, np. od razu można zgadnąć, że:
y(0) = d, oraz: y'(0) = c = v(0);

Pozostają tylko a i b do wyznaczenia, co też można sobie łatwo wyliczyć:
y(1) = a+b+c+d
y'(1) = 3a + 2b + d;
...

OK. Mamy ten wielomian 3-stopnia... i cóż z tego?

Gdzie jest to minimum:
|r(t_min)| = r_min, jak wyznaczyć: t = t_min, żeby wyliczyć pozycję: r(t_min) = ?

0

Wyznacz sobie ekstrema tej funkcji po prostu i wybierz minimum. Ekstrema są oczywiście tam gdzie pochodna się zeruje.

0
Shalom napisał(a):

Wyznacz sobie ekstrema tej funkcji po prostu i wybierz minimum. Ekstrema są oczywiście tam gdzie pochodna się zeruje.

Wiem że tak to się robi, ale chyba dla funkcji skalarnych jedynie.

Ale tu masz funkcję wektorową!

Te współczynniki: a, b, c i d są wektorami typu: a = (ax,ay,az).

Zatem tu są w zasadzie trzy funkcje kubiczne, a należy wyliczyć: r(t)=(x,y,z)(t) dla minimum typu: |r|.

0

cd.
zatem podchodząc z marszu otrzymamy tu problem ekstremum warunkowego dla funkcji:
|r|, przy warunkach r = r(t) = trzy równania kubiczne;

można wyliczyć minimum |r|^2, zamiast |r|:

|r|<sup>2 = |at</sup>3 + bt<sup>2 + ct + d|</sup>2 = a_xt<sup>3 + b_x t</sup>2 + c_x t + d_x)<sup>2 + ... = p t</sup>6 + ...

zatem pochodna z tego będzie wielomianem aż 5-go stopnia, którego nie da rady w ogóle rozwiązać analitycznie!

Problem jest raczej prosty, więc coś tu nie gra - musi być inna metoda.

1

Zatem biorę interpolację kubiczną, bo taka akurat tu podchodzi:

Dlaczego? Orbita jest stożkową (elipsą, hiperbolą lub parabolą), jest więc krzywą drugiego stopnia. Na dodatek, orbita jest płaska, więc przy odpowiednim wyborze układu współrzędnych zagadnienie jest dwuwymiarowe.

0
  1. v i r to odpowiednio prędkość (powinna być liniowa) i promień - czyli pozycja 2D (o ile ta notacja nie jest wyssana z palca)
  2. czegoś tu brakuje: masz 2 punkty i wzór wskazujący na pozycję "y" w czasie - nie ma łącznika. Nie wszystko zanotowałeś czy robisz sobie jaja? ("można sobie łatwo wyliczyć")
  3. źle policzyłeś y'(1)

Moja wiedza w zakresie astronomii jest niewielka, ale pewnie to się przyda:
https://pl.wikipedia.org/wiki/Prawa_Keplera

0
bogdans napisał(a):

Zatem biorę interpolację kubiczną, bo taka akurat tu podchodzi:

Dlaczego? Orbita jest stożkową (elipsą, hiperbolą lub parabolą), jest więc krzywą drugiego stopnia. Na dodatek, orbita jest płaska, więc przy odpowiednim wyborze układu współrzędnych zagadnienie jest dwuwymiarowe.

Niezły żart.

To że krzywe stożkowe są wielomianami II-go rzędu wcale nie oznacza, że to są krzywe II-go stopnia w ramach teorii wielomianów.

Przykładowo, weź sobie okrąg, czyli takie coś:
r(t) = 1

a teraz popróbuj sobie interpolować to parabolą... haha!

0
vpiotr napisał(a):
  1. v i r to odpowiednio prędkość (powinna być liniowa) i promień - czyli pozycja 2D (o ile ta notacja nie jest wyssana z palca)
  2. czegoś tu brakuje: masz 2 punkty i wzór wskazujący na pozycję "y" w czasie - nie ma łącznika. Nie wszystko zanotowałeś czy robisz sobie jaja? ("można sobie łatwo wyliczyć")
  3. źle policzyłeś y'(1)

Moja wiedza w zakresie astronomii jest niewielka, ale pewnie to się przyda:
http://pl.wikipedia.org/wiki/Prawa_Keplera

Problem jest taki:
jest Słońce, dookoła którego krąży Ziemia - zgodnie z prawem Newtona.

Zatem biorę sobie i wyliczam kolejne punkty po tej orbicie... no i chcę przy okazji wyliczyć to minimum, zwane peryhelium w astrologii.

Jednak orbita nie jest wcale elipsą, lecz czymś zbliżonym do elipsy, z uwagi na inna planet - Wenus, Mars, itd.

Tym samym nie mogę tu użyć wzoru na elipsę...
no, ale zawsze będzie te jakieś minimum, który chcę właśnie wyznaczyć, i z dwóch punktów bo to jest łatwe!

Obliczam tylko: |r|' w kolejnych punktach, i gdy otrzymam w dwóch kolejnych punktach sytuację typu:
|r0|' < 0, oraz |r1|' > 0 no to pomiędzy nimi musiał być punkt r, w którym |r(t_pomiędzy)| = 0,
czyli tam było właśnie minimum odległości |r|!

|r|' = v_r - prędkość radialna, czyli przyrosty odległości, które liczymy tak: |r|.

0

Dobra. Może taki sposób byłby tu dostatecznie dobry:

zamiast interpolować te wektory, znaczy: r(t) = (x,y,z)
interpolujemy samą odległość, czyli to y = |r|.

wtedy otrzymam też kubiczną, bo mam zarówno: y(t) jak i y' = |r|' = v_r, w obu punktach.

Wyliczę sobie to minimum odległości, wprost z zera pochodnej, czyli z rów. kwadratowego, co jest łatwe.
v_r(t_min) = 3at2 + 2bt + c = 0 -> t_min = ...

Tylko że mi potrzeba cały: r - pozycję, nie samo |r| !
Jak to wyliczyć?

Chyba wracam z powrotem do kubicznej wektorowej:
r(t) = ... i tam podstawiam wyliczony czas: t_min?

Jakieś pojebane to jest.

0

bogdans pisze: "O czym Ty piszesz? Okrąg będę interpolował okręgiem, który jest krzywą drugiego stopnia. Skąd pomysł paraboli?

Masz nasrane.

Okrąg w funkcji jednej zmmiennej - czyli t - czas to takie coś:

x = Rcos(wt);
y = R
sin(wt);
gdzie w = const.

a jest wiadome od kilku wieków, że:
sinx = x - x<sup>3/6 + x</sup>5/120 + ...
i kosinus podobnie...

czyli nie jest to wielomian 2 rzędu, ani też trzeciego, lecz nieskończonego, niestety!

Okrąg w wersji typu: r = r(t), byłoby to po prostu: r(t) = R = const;

No, a po przewaleniu do szkolnego układu Kartezjusza x-y masz takie coś:

y = sqrt(R2-x2);
jak widać nadal nie przypomina to wielomianu II-stopnia... i z oczywistych przyczyn - okrąg nie jest parabolą!

1

Przestań się kompromitować, wbrew swojemu mniemaniu nie masz pojęcia o matematyce i zupełnie brakuje Ci kultury. Równanie okręgu wygląda tak: y2 + x2 = R2, jak widać jest to wielomian 2 stopnia.
Dla rozrywki masz cztery równania, jedno opisuje hiperbolę, drugie parabolę, trzecie elipsę nie będącą okręgiem, czwarte okrąg. Spróbuj powiązać równania z krzywymi.
x<sup>2+y</sup>2-x\cdot y+x+2\cdot y-1=0
x<sup>2+y</sup>2-0,9\cdot x\cdot y+x+2\cdot y-1=0
x<sup>2+y</sup>2-1,1\cdot x\cdot y+x+2\cdot y-1=0
x<sup>2+y</sup>2-2\cdot x-y+6=0
Lubisz współrzędne biegunowe, równanie orbity (w polu centralnym) we współrzędnych biegunowych wygląda tak;r = \frac{p}{1+e\cdot\cos\phi}. Parametr p nie wpływa na rodzaj krzywej, parametr e (mimośród) wpływa. Dla pewnych e otrzymamy elipsę, dla innych parabolę lub hiperbolę. Druga zagadka, dla jakich wartości mimośrodu orbita jest parabolą?
O astronomii też chyba wiesz nie wiele, skoro nie wiesz, że orbita w układzie słonecznym (w każdym układzie centralnym) może być parabolą lub hiperbolą. I dla wielu ciał niebieskich jest taką właśnie krzywą.

0

Naprawdę bardzo odkrywczy jesteś - niemalże jak moja.. myszka domowa! :)

0

To musi być frustrujące, mieszkać z myszką, której się nie dorasta do pięt - intelektualnie oczywiście.

0

No to pokaż w końcu jak rozwiązałeś to minimum na krzywej w 3D.. z dwóch punktów - krzywa Hermita itd.

Twoje wywody są prześwietne... ale bezużyteczne.

0

Nie rozwiązywałem tego problemu, wyrażałem tylko bardzo poważne wątpliwości co do poprawności Twojego podejścia.
Mogę spróbować, ale potrzebuję wyjaśnienia kilku wątpliwości. Najważniejsza jest taka: interpolujemy na podstawie dwóch tylko punktów => dokładność jest bardzo mała, jednocześnie mamy uwzględniać wpływ innych planet, co wymaga dużej dokładności. Możesz podać przykładowe dane wejściowe (r_1,r_2,\vec {v_1}, \vec {v_2}) , chciałbym zobaczyć czy wektory \vec {v_1},\vec {v_2} leżą w płaszczyźnie wyznaczonej przez r_1, r_2, a jeśli nie leżą, to jak duże są ich składowe prostopadłe do tej płaszczyzny.

0
bogdans napisał(a):

Nie rozwiązywałem tego problemu, wyrażałem tylko bardzo poważne wątpliwości co do poprawności Twojego podejścia.
Mogę spróbować, ale potrzebuję wyjaśnienia kilku wątpliwości. Najważniejsza jest taka: interpolujemy na podstawie dwóch tylko punktów => dokładność jest bardzo mała, jednocześnie mamy uwzględniać wpływ innych planet, co wymaga dużej dokładności. Możesz podać przykładowe dane wejściowe (r_1,r_2,\vec {v_1}, \vec {v_2}) , chciałbym zobaczyć czy wektory \vec {v_1},\vec {v_2} leżą w płaszczyźnie wyznaczonej przez r_1, r_2, a jeśli nie leżą, to jak duże są ich składowe prostopadłe do tej płaszczyzny.

Masz dane te dwa punkty r w znanym odstępie czasu: dt;
Masz też prędkość: v w obu punktach...

w sumie masz tu całe setki takich punktów, bo liczymy to stale podczas symulacji numerycznej.. i całego układu wielu ciał, np. Słonecznego: Merkury, Wenus, ... itd.

Niemniej dwa kolejne punkty mi wystarczą, bo do wychwycenia minimum |r| wystarczy sprawdzać iloczyny skalarne typu:
v(t)*r(t) < 0, a w kolejnym kroku masz już: v(t+dt)*r(t+dt) > 0,

czyli minimum |r|, które nas interesuje siedzie pomiędzy tymi punktami, więc sobie interpolujemy i wyliczamy.

A ponieważ dwa punkty, r i v każdy, dają krzywą kubiczną, więc należy wyliczyć ekstremum kubicznej.

Dodatkowo wiemy że prędkość radialna czyli |r|' to takie coś:
v_r = v * r / |r|, więc można sobie to użyć w wyznaczaniu minimum.

0
bogdans napisał(a):

chciałbym zobaczyć czy wektory \vec {v_1},\vec {v_2} leżą w płaszczyźnie wyznaczonej przez r_1, r_2, a jeśli nie leżą, to jak duże są ich składowe prostopadłe do tej płaszczyzny.

W praktyce one nie leżą w płaszczyźnie, a jedynie prawie, co znaczy, że i ta krzywa kubiczna: r = r(x,y,z,t), nie będzie płaska, lecz prawie płaska.

A dzieje się tak z powodu tzw. precesji płaszczyzny orbity... co w astronomii nazywają precesją węzłów.

user image

0

W sumie to kicha z tego wychodzi.

Biorąc symetryczne dane, jak jest na elipsie, czy dowolnej krzywej symetrycznej czyli:
y0 = y1, oraz: v0 = -v1;

wówczas interpolacja kubiczna redukuje się do kwadratowej.. bowiem wtedy: a = 0.

y(t) = (2t<sup>3-3t</sup>2+1)r0 + (t<sup>3-2t</sup>2+t)v0 + (-2t<sup>3+3t</sup>2)r1 +(t<sup>3-t</sup>2)v1&lt;/tex]</p> <p>co w wersji typu:<br> <tex>y = at<sup>3 + bt</sup>2 +ct + d

daje:
a = -2(r1-r0) + (v0+v1)
b = 3(r1-r0) -2v0-v1
c = v0; d = r0

zatem mając te symetryczne dane, otrzymamy:
a = 0, b = -(2v0+v1) = -v0 - (v0+v1) = -v0;

czyli zamiast kubicznej będzie:
y = -v0t<sup>2 +v0t + r0 = v0(-t</sup>2+t) + r0
kwadratowe tylko...

pochodna z tego będzie tylko liniowa, czyli kubiczne nic tu nie dają...

0

Niemniej wsadzając t = 1 otrzymamy:
y = r0, zamiast r1, czy to jest kompletnie zjeb***...

Prawdopodobnie idioci produkują te wzory.

0

Oj! Sorry... OK jest, przecież r0 = r1. :)

1 użytkowników online, w tym zalogowanych: 0, gości: 1