Message Boards Message Boards

Geodesics equation - solving and plotting via Mathematica

Recently I've found Mathematica software so I tried to use it for solving some equations. I went to library and found entire book wrote as manual for this software. And there was chapter dedicated to geodesics. I tried to run:

  Module[{x,P,z},
    b = 2 Pi;
    x = torus[3,1,1];
    P={Pi,Pi};
    sg= solvegeo[x,P,b,30];
    z[2]=ParametricPlot[Evaluate[{u[t],v[t]}/.sg],{t,0,b},
        PlotRange\[Rule]{-5,10},PlotStyle\[Rule]Hue[0]];
    z[1] = ParametricPlot3D[Evaluate[x[u[t],v[t]]/.sg], {t,0,b}];
    S @ GraphicsArray[Array[z,2],GraphicsSpacing\[Rule] 0]
    ]
  Show[Graphics[Table[complete[{u[t],v[t]}/.sg,Hue[0]],
        {t,0,b,1/2}]],AspectRatio->Automatic]

But unfortunately instead of solution soft. gave me:

\!(* RowBox[{(NDSolve::"ndnum"), ((:)(\ )), "\<\"Encountered non- numerical value for a derivative at \!\(s\) == \!\(0.`\). \!\(\*ButtonBox[\\"MoreĀ…\\", \ ButtonStyle->\\"RefGuideLinkText\\", ButtonFrame->None, \ ButtonData:>\\"NDSolve::ndnum\\"]\)\"\>"}])

And a lot of other's.

Any ideas what is wrong?

Thanks Michal Buday

POSTED BY: michal buday
5 Replies

For surfaces z=f[x,y] the formula is already in Wolfram Language in the New Kind of Science book in a note on page 1049

That page also has differential geometry formulas, including Christoffel symbols.

POSTED BY: Todd Rowland
solvegeo[x_, {u0_, v0_}, b_, m_] :=
  Module[{su, e, so},
    christoff2[x][u, v];
    su = {u -> u[s], v -> v[s], p -> u'[s], q -> v'[s]};
    e[j_] := e[j] = G[j, 1, 1]p^2 + 2G[j, 1, 2]p q + G[j, 2, 2]q^2 /. su;
    eqic := {u''[s] + e[1] == 0, e[2] + v''[s] == 0,
        u[0] == u0, v[0] == v0, u'[0] == Cos[Theta], v'[0] == Sin[Theta]}; 
    so := NDSolve[eqic, {u, v}, {s, 0, b}];
    Flatten[Table[so, {Theta, 2 Pi/m, 2 pi, 2 pi/m}], 1]
    ]

Ok so this is function that computes the geodesics i hope. I'm new user of Mathematica. previously i used Scilab for numerical computations. So now i have function that computes the geodesics. Last step is defining the objects like sphere or thorus.

to := torus[3, 1, 1]

But now the question is, if i want to run the whole code together should i order it in the same notebook? First should be the function solvegeo, then defining the torus and last my first code posted here?

Frank: I know that the easiest way how to define this curve is to use the arc length parametrization, but many times is not so easy to find that parametrization e.g. rotation ellipsoid.

POSTED BY: michal buday

If i can get back to the point i can show you a little math at this moment. Using Christoffel symbols of second kind i can write system of partial differencial equations which describes geodesics.

So (is there any chance that i can input latex code?):

$$u^{\prime \prime} + \Gamma^{1}_{11} (u^\prime)^2 + 2\Gamma^1_{12} u^\prime v^\prime + \Gamma^{1}_{22} (v^\prime)^2 = 0$$ $$v^{\prime \prime} + \Gamma^{1}_{21} (u^\prime)^2 + 2\Gamma^2_{12} u^\prime v^\prime + \Gamma^{2}_{22} (v^\prime)^2 = 0$$

There are a few ways how to compute Christoffel symbols easiest one is from metric of a surface. When i want to solve this type of eqaution systems a need to define some conditions (sorry i am not familiar with english therm). In this case they are:

u[0] == u0, v[0] == v0, u'[0] == Cos[Theta], v'[0] == Sin[Theta]}; 

So i think that problem could be undefined the angle [Theta] that says the direction, because it's not in the input. I would like to ask if anyone could help me and try to run this on his computer.

First i define function:

solvegeo[x_, {u0_, v0_}, b_, m_] :=
  Module[{su, e, so},
    christoff2[x][u, v];
    su = {u -> u[s], v -> v[s], p -> u'[s], q -> v'[s]};
    e[j_] := e[j] = G[j, 1, 1]p^2 + 2G[j, 1, 2]p q + G[j, 2, 2]q^2 /. su;
    eqic := {u''[s] + e[1] == 0, e[2] + v''[s] == 0,
        u[0] == u0, v[0] == v0, u'[0] == Cos[Theta], v'[0] == Sin[Theta]}; 
    so := NDSolve[eqic, {u, v}, {s, 0, b}];
    Flatten[Table[so, {Theta, 2 Pi/m, 2 pi, 2 pi/m}], 1]
    ]

and apply it on:

hp:=paraboloid[-1,1], b:=2,m:=50 ,sg:=solvegeo[hp,{1,-1}/2,b,m]

or, you can try any of imported notebooks

Thank you for any good ideas

Attachments:
POSTED BY: michal buday

Calculating a geodesic on a surface is a minimization problem. You need to express the distance between two points in terms of the parameters that describe the surface and minimize that distance.

POSTED BY: Frank Kampas

As far as I can tell from the code you posted the example contains a number of undefined functions and parameters: e.g., torus, solvegeo, and so on.... Without these the code will not run.

Also, use the code sample icon in the formatting bar to post code rather than the blockquote icon.

POSTED BY: David Reiss
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract