9

This question has been asked twice but not been answered on this site, and I really need to understand this concept.

All I want is, based on 2 arbitrary points A and B on a given ellipsoid, to find the geodesic containing those 2 points.

Makogan
  • 3,329
  • Since some ellipsoids are surfaces of revoution, https://en.wikipedia.org/wiki/Clairaut%27s_relation helps really much in such cases. – Jack D'Aurizio Oct 22 '17 at 02:08
  • 1
    @Jack D'Aurizio: two of the axes should be equal for that; but it's a good start, past the sphere case. – orangeskid Oct 22 '17 at 02:15
  • I wrote a full-on honors research paper on this freshman year. I'll go through it now and write as in-depth of an answer as I can, please don't delete this question. – Andrew Tawfeek Oct 22 '17 at 02:16
  • Andrew Tawfeek, thank you. I won't. – Makogan Oct 22 '17 at 02:20
  • 1
    It looks like the first solution to the triaxial case came from Jacobi: https://en.wikipedia.org/wiki/Geodesics_on_an_ellipsoid – Jack D'Aurizio Oct 22 '17 at 02:32
  • If you want to look at open source software that does this well, please see GeographicLib. The author, Charles Karney, is a major contributor to the related Wikipedia pages, eg Geodesics on an ellipsoid. FWIW, in this answer I show how to do simple geodesic triangulation on an ellipsoid of revolution, in particular the WGS84 ellipsoid; there's also a link to some relevant Python code that uses GeographicLib. – PM 2Ring Oct 22 '17 at 07:51
  • Unfortunately I cannot rely on simply having a computer do it for me, I need to understand at least some of the theory. I am doing research on computer graphics related to this topic – Makogan Oct 23 '17 at 03:04

2 Answers2

13

So a geodesic is essentially the minimization of a functional, that is, a function of functions. They are the main study of Calculus of Variations. Now, if we have a functional of the form $$J[y]=\int_a^b F(x,y,y') \ dx$$ defined on the set of functions $y(x)$ which have continuous first derivatives in $[a,b]$, and satisfy the boundary conditions $y(a)=A$, $y(b)=B$, then a neccesary condition for $J[y]$ to have an extremum for a function $y(x)$ is that $y(x)$ satisfy the Euler-Lagrange equation $$F_y-\frac{d}{dx}F_{y'}=0.$$

In essence, we need to solve that above differential equation.

Without getting anymore nitty gritty on things, the functional we care about here is the arc length integral.


Firstly, suppose we have a surface $\mathcal{S}$ defined by the vector equation $$\vec{r}(u,v)=x(u,v) \hat{i} + y(u,v) \hat{j} + z(u,v) \hat{k}.$$ The geodesic curve lying on surface $\mathcal{S}$ can be specified by the equations $$u=u(t), \qquad v=v(t),$$ and can be found by minimizing the arc length integeral $$L=\int_{t_0}^{t_1} \sqrt{\left( \frac{dx}{dt} \right)^2+\left( \frac{dy}{dt} \right)^2+\left( \frac{dz}{dt} \right)^2} \ dt$$ But since $x,y,z$ are each functions of more than one variable, the chain rule is going to turn this integral into a monstrosity. Thankfully, it's a monstrosity with a pretty neat pattern so we can make things look a bit prettier.

Letting $u'=\frac{du}{dt}$ and $v'=\frac{dv}{dt}$, the integeral above can just be rewritten as

$$J[u,v]=\int_{t_0}^{t_1} \sqrt{Eu'^2+2Fu'v'+Gv'^2} \ dt,$$

where $E$, $F$, and $G$ are the coefficients of the first fundamental (quadradic) form of the surface, i.e., $$\begin{cases} \displaystyle{E=\vec{r}_u\cdot \vec{r}_u = \left( \frac{\partial x}{\partial u} \right) ^2 +\left( \frac{\partial y}{\partial u} \right) ^2 + \left( \frac{\partial z}{\partial u} \right) ^2} \\ \\ \displaystyle{F=\vec{r}_u \cdot \vec{r}_v = \frac{\partial x}{\partial u} \frac{\partial x}{\partial v} + \frac{\partial y}{\partial u} \frac{\partial y}{\partial v} + \frac{\partial z}{\partial u} \frac{\partial z}{\partial v}}\\ \\ \displaystyle{G=\vec{r}_v\cdot \vec{r}_v=\left( \frac{\partial x}{\partial v} \right) ^2 +\left( \frac{\partial y}{\partial v} \right) ^2 + \left( \frac{\partial z}{\partial v} \right) ^2}. \end{cases}$$

The Euler-Lagrange equation in this case corresponds to the two different equations $$\displaystyle{ F_u- \frac{d}{dt} F_{u'} = 0,} \ \ \ \ \ \ \ \displaystyle{F_v- \frac{d}{dt} F_{v'} = 0,}$$ hence, we obtain \begin{equation}\label{geodesiceq1} \frac{E_u u'^2+2F_u u' v' + G_u v'^2}{\sqrt{Eu'^2+2Fu'v'+Gv'^2}} - \frac{d}{dt} \frac{2(Eu'+Fv')}{\sqrt{Eu'^2+2Fu'v'+Gv'^2}} = 0, \end{equation} \begin{equation}\label{geodesiceq2} \frac{E_v u'^2+2F_v u' v' + G_v v'^2}{\sqrt{Eu'^2+2Fu'v'+Gv'^2}} - \frac{d}{dt} \frac{2(Fu'+Gv')}{\sqrt{Eu'^2+2Fu'v'+Gv'^2}} = 0. \end{equation}

Which are the two differential equations whose solutions provide the geodesic on surface $\mathcal{S}$.

Now let me note one important case. The case where $E$ and $G$ are explicit functions of $v$ only and $F=0$, we have $$\frac{E_v+v'^2G_v}{2\sqrt{E+Gv'^2}}-\frac{d}{du} \left( \frac{Gv'}{\sqrt{E+Gv'^2}} \right) =0,$$ so $$ E_v + v'^2G_v - 2 G \sqrt{E+Gv'^2} \bigg[ \frac{v''}{\sqrt{E+Gv'^2}} + \left( \frac{1}{2} \right) \frac{v'(2Gv'v'')}{(E+Gv'^2)^{\frac{3}{2}}} \bigg] =0 $$ $$E_v + v'^2 G_v - 2Gv'' + \frac{2G^2v'^2v''}{E+Gv'^2} = 0 $$ \begin{equation} \frac{Gv'^2}{\sqrt{E+Gv'^2}} - \sqrt{E+Gv'^2} = c_1 \end{equation} which can be made even more helpful by noting $\displaystyle{v'= \frac{dv}{du}}$, giving us $$Gv'^2-(E+Gv'^2)=c_1 \sqrt{E+Gv'^2}$$ $$\left( - \frac{E}{c_1} \right) ^2 = E + Gv'^2$$ $$\frac{E^2-{c_1}^2E}{G{c_1}^2} = v'^2,$$ finally providing \begin{equation} u={c_1} \int_{v_1=v(u_1)}^{v_2=v(u_2)} \sqrt{ \frac{G}{E^2-{c_1}^2E }} \ dv. \end{equation}


Now lets talk about the ellipsoid. An ellipsoid with semi-axis lengths $a,b,c$ is given by the vector equation $$\vec{r}(u,v)= \langle a \cos(u) \sin(v), b \sin(u) \sin(v), c \cos(v) \rangle$$ with $u \in [0,2\pi)$ and $v\in [0, \pi]$.

The first fundamental form of the ellipsoid which I sadly found by hand before I realized I could've went to Google for is

$$\begin{cases} E= b^2 \cos^2(u) \sin^2(v)+a^2 \sin^2(u) \sin^2(v) \\ F= \left( b^2-a^2 \right) \sin(v) \cos(v) \sin(u) \cos(u)\\ G= \cos^2(v) \left[ (a^2-b^2) \cos^2(u)+b^2 \right] +c^2 \sin^2(v) \end{cases}$$

Now listen, I'm actually procrasinating the hell out of a transfer essay by being here and looking at the fundamental form above is pretty worrying, because it's far from elegant and being easy to deal with and would entail solving that frightening differential equation I mentioned earlier.

So I'm going to let $a=b$. Playing around with $c$ is ellipsoid-y enough for me, and if you want to go for the full $a,b,c$ then be my guest, haha! So doing what I said, things are a lot prettier:

$$\begin{cases} E= a^2 \sin^2(v) \\ F= 0 \\ G=b^2 \cos^2(v) + c^2 \sin^2(v) \end{cases}$$

and this is where we can apply that last remark I made earlier regarding the first fundamental form looking like this. We thus have the geodesic on this surface given by

$$u={c_1} \int \sqrt{ \frac{G}{E^2-{c_1}^2E }} \ dv = c_1 \int \sqrt{\frac{b^2 \cos^2(v) + c^2 \sin^2(v)}{a^4 \sin^4(v) - c_1^2 a^2 \sin^2(v)}} \ dv.$$

Now in the case of a sphere, solving the integeral and rearranging gives a plane -- and that planes interesection with the sphere represents geodesics on a sphere. I myself couldn't get anywhere with integeral, but I'm sure someone else can so I'll just leave things off here.

... Now back to my essay.


On a very final note, here is my research paper from my first semester. I prove every step up to proving the very Euler-Lagrange equations themselves and then show some examples of it in use, like finding the geodesic on a sphere and geodesic of a function revolved around an axis.

Like Jack mentioned in the comments to your post, if you consider just the ellipsoids that can be obtained through revolving a function $f(x)$, the integeral is much easier. I prove it in the above paper, but here it is again:

$$v=c_1 \int_{u_1}^{u_2} \frac{\sqrt{1+(f'(u))^2}}{f(u) \sqrt{(f(u))^2-{c_1}^2}} \ du$$

Note, $x=u$ in the above integeral.

  • From your calculation, your $G$ should be (after setting $a=b$) $G=b^2\cos^2(v)+c^2\sin^2(v)$? –  Oct 22 '17 at 04:38
  • @JohnMa Ahh yes, you're correct! Sorry about that, let me edit it now. – Andrew Tawfeek Oct 22 '17 at 04:48
  • Thank you so much, I am an undergrad in computer science, and this is hard to swallow right now so i am just going to stare at it for a while doing my best to understand it. Thank you very much.

    As a side question, please tell me you are at least a Grad student so that my ego doesn't plummet to the ground :p

    – Makogan Oct 22 '17 at 04:55
  • @Makogan there will always be students stronger than you, and for many undergraduate math majors, differential geometry is a standard topic. Don't be disheartened, even if many know this material well. – Andres Mejia Oct 22 '17 at 04:58
  • Andres Mejia, no the disheartening thing here would be, that if this is being taught in FRESHMAN years at his uni at the undergraduate level, I feel I went to a bad university :p – Makogan Oct 22 '17 at 05:00
  • @Makogan I dont want to turn this into an extended conversation, but this probably is a result from the author's freshman year undergrad. It was a "research topic," evidently, so just something that was interesting to the student. I'm sure that with some effort students at your own university could produce work as "complicated" or "difficult." – Andres Mejia Oct 22 '17 at 05:03
  • @AndresMejia Was actually right on point. If it makes you feel better, I actually just go to a community college (in the US)! I took multivariable/multivariate Calculus my first semester when I enrolled and as part of an honors program at this college I had the wonderful opportunity to do a semester-long research project on a topic of my choosing alongside a professor and I went with geodesics! So it took me an entire semester of work to get that 13-page paper, don't feel too bad ahaha! Also, there's no class on Differential Geometry here, it only goes up to Differential Equations :) – Andrew Tawfeek Oct 22 '17 at 05:15
  • If you do not understand the answer, you don't need to accept it. @Makogan –  Oct 22 '17 at 05:15
  • Understanding has multiple levels, I have stared at wolframs definition on geodesics enough to see that his answer is correct. And he actually gave enough information such that If I work out the math I feel I'll understand what is happening. I won't be able to read this in depth for now because I have 5 midterms 2 assignments 1 group project... Comming up all of next week, so i accepted it for now, since on a high level analysis it seems exactly what I need to solve my problem. – Makogan Oct 22 '17 at 05:19
  • Let me just be clear: I mean you don't have to accept it right away. You can wait a bit until you have time, or wait till there are more answers @Makogan . –  Oct 22 '17 at 05:35
  • I'll give a hint of the punch line: you won't be able to avoid dealing with (incomplete) elliptic integrals. – J. M. ain't a mathematician Oct 22 '17 at 06:05
  • @J.M ELI5 please? – Makogan Oct 22 '17 at 06:17
  • 2
    @Mak, I will say that explaining elliptic integrals to a five-year-old looks like an unreasonable request. ;) Nevertheless: these are nonelementary functions that will inevitably show up when you're integrating something that came out of an ellipse or an ellipsoid. Most modern computing environments should have a stable implementation of them. – J. M. ain't a mathematician Oct 22 '17 at 06:26
  • @Andrew Tawfeek: Which suggests you must have had a ton of prior maths knowledge before coming to that college, well in advance of your typical freshman to college. Where/how did you get that? Since that's at least a full course of lin alg at the beginning of that paper plus multivariable calc (which I just see you mentioned) and some intro Analysis in there as well. – The_Sympathizer Oct 22 '17 at 15:27
  • @The_Sympathizer Ahaha I appreciate the comment! I actually learned all that stuff as I did the research on geodesics/studied Calculus of Variations, so that was my first exposure to spaces. I actually still haven't learned Linear Algebra or Analysis yet, I've been spending all my time since then with Aluffi's amazing Chapter 0 Algebra book :) – Andrew Tawfeek Oct 22 '17 at 15:35
  • @Makogan If you feel I've sufficiently answered your question, I would appreciate marking it as the answer :) – Andrew Tawfeek May 02 '20 at 01:53
  • 1
    @AndrewTawfeek woops, this one slipped under the radar – Makogan May 02 '20 at 23:57
2

The following addresses only initial values only.

If all geodesics are tangent to circle

$$ x^2+y^2= r_{min}^2 @ \,r=r_{min}$$

as the Clairaut constant then for an ellipsoid of revolution

$$ (z/c)^2+(r/a)^2 = (z/c)^2+(x^2+y^2)/a ^2 = 1 $$

we differentiate to find meridian slope

$$ \tan \phi=-(a/c) \sqrt{(a/r)^2-1}, $$

and $r=f(\theta) $ is found by integration (to get elliptic integral result)

$$ \frac{dr}{r \, d\theta}= \sin\phi \, \sqrt{1-(r/r_{min})^2 }.$$

similarly $z$ can be set up. We have $ (r,z)$ as functions of $\theta$

With given $ (r_B,r_A), \theta_{AB}, (z_B,z_A) $ we cannot directly use this, however..

we can start with initial conditions $(r_B, \theta_A =0, z_A ) $ but to arrive at last point we can vary /adjust Clairaut constant $r_{min}$ several small increments in a numerical procedure.

Narasimham
  • 40,495