Witam,
Wiem, że wątek jest stary, ale jest pierwszy pod względem ilości wyświetleń dla tagu "Algorytmy i struktury danych" i nie jest do końca wyjaśniony. A więc zacznę od tego, że przedmiotem pytania jest tzw. "odwrotny problem geodezyjny" ( można to znaleźć w angielskiej wikipedii: https://en.wikipedia.org/wiki/Geodesy#Geodetic_problems ). Zwykły problem geodezyjny nie jest aż tak trudny do rozwiązania i jak się zna analizę tensorową to rozwiązanie jest w miarę proste. A więc, żeby rozwiązać problem "w przód" ( czyli znamy współrzędne startu i początkowy kierunek poruszania się ) to wystarczy zapisać równania ruchu po liniii geodezyjnej na elipsoidzie obrotowej i wycałkować je numerycznie. Równanie linii geodezyjnej ( opisane tu: https://en.wikipedia.org/wiki/Geodesic ) można traktować w pewnym sensie jako uogólnienie drugiej zasady dynamiki z zerowym przyspieszeniem dla zakrzywionych przestrzeni i dowolnych układów współrzędnych krzywoliniowych. W wzorze na linię geodezyjną pojawiają się symbole Christoffela, które zależą od zakrzywienia przestrzeni i tego jakiego układu współrzędnych się używa. Tutaj wydaje mi się, że najłatwiej jest użyć współrzędnych geocentrycznych ( nie są to te same współrzędne, których się na codzień używa! jednak przeliczenie jednych w drugie jest łatwe, różnica między nimi jest opisana tutaj: https://en.wikipedia.org/wiki/Geodetic_coordinates#Geodetic_vs._geocentric_coordinates ). Żeby obliczyć symbole Christoffela ( tutaj jest opis: https://en.wikipedia.org/wiki/Christoffel_symbols#General_definition ) trzeba najpierw obliczyć tensor metryczny, którego składowe to odpowiednie iloczyny skalarne kowariantnych wektorów bazowych ( opis co to jest baza kowariantna: https://en.wikipedia.org/wiki/Covariance_and_contravariance_of_vectors ). Jak się wszystko złoży do kupy to otrzymuje się dwa równania różniczkowe drugiego stopnia jedno na szerokość a drugie na długość geocentryczną. Do tego się dopisuje dwa równania na pozycję początkową i dwa równania na prędkość początkową i voila - całkuje się to numerycznie. Co do 2 równań na prędkość to w rzeczywistości wystarczy tylko jeden parametr - kierunek początkowy, wartość prędkości można przyjąć równą np. 1 i nie zmieni ona kształtu samej linii geodezyjnej. Opis całej procedury wyprowadzenia równań opisałem tutaj: http://www.zenker.pl/fizyka/analiza-tensorowa/linie-geodezyjne-na-elipsoidzie-obrotowej.pdf ( tekst nie jest skończony, ale wystarczy doczytać do zdania "Rozwiązanie drugiego równania" ). Poniżej wstawiam przykładowe linie geodezyjne wygenerowane w wolframie:
geodesic1.png
geodesic2.png
geodesic3.png
Odwrotny problem geodezyjny wydaje się trudniejszy i niedokończony opis jak się to wylicza zostawiłem w dalszej części pdfa.
Link do notebooka wolframa ( jeżeli ktoś chce sobie samemu zmieniać parametry ): http://www.zenker.pl/fizyka/analiza-tensorowa/geodesic.nb