I want to solve (numerically) the geodesics of the upper half plane and plot. The results are quite known (i) straight lines parallel to $y$-axis and (ii) semicircles centered on the $x$-axis. Now the problem is how to solve the coupled ODEs. I tried with NDSolve but the results are not coming. Can you please write a code and help me.

Here are the coupled ODEs:

$$\begin{align*} \ddot{x} &= \frac{2\dot{x}\dot{y}}{y}\\ \ddot{y} &= \frac{-\dot{x}^2+\dot{y}^2}{y} \end{align*}$$

Thanking you in advance.

  • $\begingroup$ dear chris can you please explain why the t in the NDSolve runs for 0 to 1 and then in parametric plot from 1 to 10. Also If I can put range of t 1 to 10 in NDSolve then two curves are coming ??? can you please explain $\endgroup$
    – ricci1729
    Commented May 25, 2013 at 10:01

2 Answers 2


If you use NDSolve as follows

  NDSolve[{ x''[t] == 2 x'[t] y'[t]/y[t], 
          y''[t] == (y'[t]^2 - x'[t]^2)/y[t], 
          x[0] == 1, y[0] == 1, 
          x'[0] == 1, y'[0] == 1}, {x[t], y[t]}, 
           {t, 0, 8}]; 


 ParametricPlot[{x[t], y[t]} /. sol, {t, 0, 8}]


Note that if/when the solution is driven towards 'y=0' you will run into problems.

You might also want to have a look at EquationTrekker

<< EquationTrekker`
  • $\begingroup$ dear chris; that's very true . y=0 is the singularity that one can immediately see from the ODE's. But my question is, if you change the limit of t from 0 to 9 you see two curves. but why is that is not clear to me??? $\endgroup$
    – ricci1729
    Commented May 25, 2013 at 10:18
  • $\begingroup$ its a bit strange... $\endgroup$
    – chris
    Commented May 25, 2013 at 10:28

I think you can solve the equations analytically. The change of variables y[t] -> Exp[z[t]] changes your system into

x''[t] == 2 x'[t] z'[t]
Exp[2 z[t]] z''[t] == - x'[t]^2

The first equation is easy and you get

z[t] == 1/2 Log[x'[t]] + a

with a a constant; the second equation, assuming x'[t]!=0 gives

1/2 Exp[2 a] x''[t]/x'[t] == -x[t] + b

with b a constant. This can be handled by Mathematica :

solx[t_] = x[t] /. First@
  DSolve[1/2 Exp[2 a] x''[t]/x'[t] == -x[t] + b, x[t], t] /. {C[1] -> c, C[2] -> d}

from which we can retrieve the solution for z[t] :

solz[t_] = 1/2 Log[D[solx[t], t]] + a

In order to check against @chris let's solve for the constants {a,b,c,d} by using the initial conditions :

const = FindRoot[
          {solx[0] == 1, Derivative[1][solx][0] == 1, 
           solz[0] == 0, Derivative[1][solz][0] == 1}, 
           {{a, 1}, {b, 2}, {c, 1}, {d, 1}}] //Chop ;

ParametricPlot[{solx[t] /. const, Exp[solz[t] /. const]}, {t, 0, 5}]


  • $\begingroup$ thank you very much but i know it is exactly solvable model. the question is how to generate this graph itself from the geodesic equation by NDSolve. $\endgroup$
    – ricci1729
    Commented May 28, 2013 at 8:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.