## Líneas geodésicas

Supongamos la superficie

$\left\{\begin{array}{l}x=x\left(u,v\right)\\ y=y\left(u,v\right)\\ z=z\left(u,v\right)\end{array}$

La longitud de la línea que une dos puntos A y B y que está contenida en dicha superficie es

$\begin{array}{l}L=\underset{A}{\overset{B}{\int }}ds=\underset{A}{\overset{B}{\int }}\sqrt{d{x}^{2}+d{y}^{2}+d{z}^{2}}\\ dx=\frac{\partial x}{\partial u}du+\frac{\partial x}{\partial v}dv,\text{ }dy=\frac{\partial y}{\partial u}du+\frac{\partial y}{\partial v}dv,\text{ }dz=\frac{\partial z}{\partial u}du+\frac{\partial z}{\partial v}dv\\ d{x}^{2}={\left(\frac{\partial x}{\partial u}\right)}^{2}d{u}^{2}+{\left(\frac{\partial x}{\partial v}\right)}^{2}d{v}^{2}+2\frac{\partial x}{\partial u}\frac{\partial x}{\partial v}dudv\\ d{y}^{2}={\left(\frac{\partial y}{\partial u}\right)}^{2}d{u}^{2}+{\left(\frac{\partial y}{\partial v}\right)}^{2}d{v}^{2}+2\frac{\partial y}{\partial u}\frac{\partial y}{\partial v}dudv\\ d{z}^{2}={\left(\frac{\partial z}{\partial u}\right)}^{2}d{u}^{2}+{\left(\frac{\partial z}{\partial v}\right)}^{2}d{v}^{2}+2\frac{\partial z}{\partial u}\frac{\partial z}{\partial v}dudv\\ L=\underset{A}{\overset{B}{\int }}\sqrt{\left[{\left(\frac{\partial x}{\partial u}\right)}^{2}+{\left(\frac{\partial y}{\partial u}\right)}^{2}+{\left(\frac{\partial z}{\partial u}\right)}^{2}\right]d{u}^{2}+\left[{\left(\frac{\partial x}{\partial v}\right)}^{2}+{\left(\frac{\partial y}{\partial v}\right)}^{2}+{\left(\frac{\partial z}{\partial v}\right)}^{2}\right]d{v}^{2}+2\left[\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}\right]dudv}\\ L=\underset{A}{\overset{B}{\int }}\sqrt{P+R{\left(\frac{dv}{du}\right)}^{2}+2Q\frac{dv}{du}}·du=\underset{A}{\overset{B}{\int }}\sqrt{P+R{\stackrel{˙}{v}}^{2}+2Q\stackrel{˙}{v}}·du\end{array}$

Emplemos el cálculo de variaciones para determinar la curva que pasa por los puntos A y B contenida en dicha superficie y cuya longitud L sea mínima. Vamos a considerar dos casos particulares:

1. P y R son funciones de u y Q=0

2. La función integrando, $f(u, v ˙ )$, no depende de v. La ecuación de Euler-Lagrange se escribe

$\begin{array}{l}\frac{\partial f}{\partial \stackrel{˙}{v}}={c}_{1}\\ \frac{R\stackrel{˙}{v}}{\sqrt{P+R{\stackrel{˙}{v}}^{2}}}={c}_{1}\\ {R}^{2}{\left(\frac{dv}{du}\right)}^{2}={c}_{1}^{2}\left(P+R{\left(\frac{dv}{du}\right)}^{2}\right)\\ v={c}_{1}\int \frac{\sqrt{P}}{\sqrt{{R}^{2}-{c}_{1}^{2}R}}du+{c}_{2}\end{array}$

3. P y R son funciones de v y Q=0

4. La función integrando, $f(v, v ˙ )$, no depende de u. La ecuación de Euler-Lagrange se escribe

$\begin{array}{l}f-\frac{\partial f}{\partial \stackrel{˙}{v}}\stackrel{˙}{v}={c}_{1}\\ \sqrt{P+R{\stackrel{˙}{v}}^{2}}-\frac{R{\stackrel{˙}{v}}^{2}}{\sqrt{P+R{\stackrel{˙}{v}}^{2}}}={c}_{1}\\ {P}^{2}-P{c}_{1}^{2}=R{\left(\frac{dv}{du}\right)}^{2}\\ u={c}_{1}\int \frac{\sqrt{R}}{\sqrt{{P}^{2}-{c}_{1}^{2}P}}dv+{c}_{2}\end{array}$

## Superficie esférica

La ecuación de una superficie esférica de radio r es

$\left\{\begin{array}{l}x=r\mathrm{cos}\phi \mathrm{sin}\theta \\ y=r\mathrm{sin}\phi \mathrm{sin}\theta \\ z=r\mathrm{cos}\theta \end{array}$

Para esta superficie, u es φ y v es θ

$\begin{array}{l}P={\left(\frac{\partial x}{\partial \phi }\right)}^{2}+{\left(\frac{\partial y}{\partial \phi }\right)}^{2}+{\left(\frac{\partial z}{\partial \phi }\right)}^{2}={r}^{2}{\mathrm{sin}}^{2}\phi {\mathrm{sin}}^{2}\theta +{r}^{2}{\mathrm{cos}}^{2}\phi {\mathrm{sin}}^{2}\theta ={r}^{2}{\mathrm{sin}}^{2}\theta \\ R={\left(\frac{\partial x}{\partial \theta }\right)}^{2}+{\left(\frac{\partial y}{\partial \theta }\right)}^{2}+{\left(\frac{\partial z}{\partial \theta }\right)}^{2}={r}^{2}{\mathrm{cos}}^{2}\phi {\mathrm{cos}}^{2}\theta +{r}^{2}{\mathrm{sin}}^{2}\phi {\mathrm{cos}}^{2}\theta +{r}^{2}{\mathrm{sin}}^{2}\theta ={r}^{2}\\ Q=\frac{\partial x}{\partial \phi }\frac{\partial x}{\partial \theta }+\frac{\partial y}{\partial \phi }\frac{\partial y}{\partial \theta }+\frac{\partial z}{\partial \phi }\frac{\partial z}{\partial \theta }=-{r}^{2}\mathrm{sin}\phi \mathrm{sin}\theta ·\mathrm{cos}\phi \mathrm{cos}\theta +{r}^{2}\mathrm{cos}\phi \mathrm{sin}\theta ·\mathrm{sin}\phi \mathrm{cos}\theta =0\end{array}$

La longitud de la línea que une dos puntos A y B de la superficie esférica y que está contenida en dicha superficie es

$L=r\underset{A}{\overset{B}{\int }}\sqrt{{\mathrm{sin}}^{2}\theta +{\left(\frac{d\theta }{d\phi }\right)}^{2}}·d\phi$

Q=0, P y R son funciones de v=θ. Estamos en el segundo caso con u=φ

$\phi ={c}_{1}\int \frac{\sqrt{R}}{\sqrt{{P}^{2}-{c}_{1}^{2}P}}d\theta +{c}_{2}$

Resolvemos la integral

$\begin{array}{l}{c}_{1}\int \frac{1}{\sqrt{{r}^{2}{\mathrm{sin}}^{4}\theta -{c}_{1}^{2}{\mathrm{sin}}^{2}\theta }}d\theta =\int \frac{\frac{1}{{\mathrm{sin}}^{2}\theta }}{\sqrt{{\left(\frac{r}{{c}_{1}}\right)}^{2}-\frac{1}{{\mathrm{sin}}^{2}\theta }}}d\theta =\int \frac{\frac{1}{{\mathrm{sin}}^{2}\theta }}{\sqrt{{\left(\frac{r}{{c}_{1}}\right)}^{2}-\frac{{\mathrm{sin}}^{2}\theta +{\mathrm{cos}}^{2}\theta }{{\mathrm{sin}}^{2}\theta }}}d\theta =\\ \int \frac{\frac{1}{{\mathrm{sin}}^{2}\theta }}{\sqrt{\left({\left(\frac{r}{{c}_{1}}\right)}^{2}-1\right)-\frac{{\mathrm{cos}}^{2}\theta }{{\mathrm{sin}}^{2}\theta }}}d\theta =\int \frac{-dw}{\sqrt{\left({\left(\frac{r}{{c}_{1}}\right)}^{2}-1\right)-{w}^{2}}}=\\ -\mathrm{arcsin}\left(\frac{w}{\sqrt{\left({\left(\frac{r}{{c}_{1}}\right)}^{2}-1\right)}}\right)=-\mathrm{arcsin}\left(\frac{\mathrm{cot}\theta }{\sqrt{\left({\left(\frac{r}{{c}_{1}}\right)}^{2}-1\right)}}\right)\end{array}$

Obtenemos la ecuación del plano que pasa por dos puntos y el origen. Su intersección con la superficie esférica de radio r es una circunferencia máxima

$\begin{array}{l}\phi =-\mathrm{arcsin}\frac{\mathrm{cot}\theta }{{c}_{3}}+{c}_{2}\\ \mathrm{sin}\left({c}_{2}-\phi \right)-\frac{\frac{\mathrm{cos}\theta }{\mathrm{sin}\theta }}{{c}_{3}}=0\\ \mathrm{sin}{c}_{2}·r·\mathrm{sin}\theta \mathrm{cos}\phi -\mathrm{cos}{c}_{2}·r\mathrm{sin}\theta \mathrm{sin}\phi -\frac{r\mathrm{cos}\theta }{{c}_{3}}=0\\ x\mathrm{sin}{c}_{2}-y\mathrm{cos}{c}_{2}-\frac{z}{{c}_{3}}=0\\ \mathrm{sin}A·x-\mathrm{cos}A·y-Bz=0\end{array}$

Elegimos dos puntos de la superficie esférica de radio unidad, r=1, especificando sus ángulos (φ, θ).

${P}_{1}\left\{\begin{array}{l}{x}_{1}=\mathrm{cos}{\phi }_{1}\mathrm{sin}{\theta }_{1}\\ {y}_{1}=\mathrm{sin}{\phi }_{1}\mathrm{sin}{\theta }_{1}\\ {z}_{1}=\mathrm{cos}{\theta }_{1}\end{array}\text{ }{P}_{2}\left\{\begin{array}{l}{x}_{2}=\mathrm{cos}{\phi }_{2}\mathrm{sin}{\theta }_{2}\\ {y}_{2}=\mathrm{sin}{\phi }_{2}\mathrm{sin}{\theta }_{2}\\ {z}_{2}=\mathrm{cos}{\theta }_{2}\end{array}$

Determinamos la ecuación del plano que pasa por el origen y por los puntos P1 y P2.

$\begin{array}{l}\left\{\begin{array}{l}\mathrm{sin}A·{x}_{1}-\mathrm{cos}A·{y}_{1}-B·{z}_{1}=0\\ \mathrm{sin}A·{x}_{2}-\mathrm{cos}A·{y}_{2}-B·{z}_{2}=0\end{array}\\ A=\mathrm{arctan}\left(\frac{{y}_{1}{z}_{2}-{y}_{2}{z}_{1}}{{x}_{1}{z}_{2}-{x}_{2}{z}_{1}}\right)\\ B=\frac{\mathrm{sin}A·{x}_{1}-\mathrm{cos}A·{y}_{1}}{{z}_{1}}\end{array}$

Para trazar la circunferencia máxima, en la ecuación del plano, sustituimos x, y, z en coordenadas polares de radio r=1

$\begin{array}{l}\mathrm{sin}A·x-\mathrm{cos}A·y-Bz=0\\ \mathrm{sin}A·\mathrm{cos}\phi \mathrm{sin}\theta -\mathrm{cos}A·\mathrm{sin}\phi \mathrm{sin}\theta -B\mathrm{cos}\theta =0\\ \mathrm{sin}\theta \left(\mathrm{sin}A·\mathrm{cos}\phi -\mathrm{cos}A·\mathrm{sin}\phi \right)=B\mathrm{cos}\theta \\ \mathrm{tan}\theta =\frac{B}{\mathrm{sin}\left(A-\phi \right)}\end{array}$

Damos valores a φ entre 0 y 2π, obtenemos el ángulo θ, dibujamos el punto de coordenadas (φ,θ), unimos los puntos para trazar la circunferencia máxima que pasa por P1 y P2

Ejemplo

• Punto P1: φ=0, θ=π/4 (45°)
• Punto P2: φ=π/2, θ=π/2 (90°) en el ecuador
phi=linspace(0,pi,30);
theta=linspace(0,2*pi,40);
[phi,theta]=meshgrid(phi,theta);
x=sin(phi).*cos(theta);
y=sin(phi).*sin(theta);
z=cos(phi);
hold on
h1=mesh(x,y,z);
set(h1,'EdgeColor',[0.6,0.6,0.6]);

%punto P1
th=pi/4;
phi=0;
x1=cos(phi)*sin(th);
y1=sin(phi)*sin(th);
z1=cos(th);
plot3(x1,y1,z1,'o','markersize',4,'markeredgecolor','b','markerfacecolor','b')

%punto P2
th=pi/2;
phi=pi/2;
x2=cos(phi)*sin(th);
y2=sin(phi)*sin(th);
z2=cos(th);
plot3(x2,y2,z2,'o','markersize',4,'markeredgecolor','b','markerfacecolor','b')

%plano
A=atan2((y1*z2-y2*z1),(x1*z2-x2*z1));
B=(x1*sin(A)-y1*cos(A))/z1;
[x,y] = meshgrid(-1:0.1:1);
z = (x*sin(A)-y*cos(A))/B;
mesh(x,y,z)

%geodésica
phi=(0:360)*pi/180;
th=atan2(B,sin(A-phi));
plot3(cos(phi).*sin(th),sin(phi).*sin(th),cos(th), 'r','lineWidth',1.5)

axis equal
axis off
view(100,5)
hold off
xlabel('x'); ylabel('y'); zlabel('z')
title('Superficie esférica')

En la figura, vemos los dos puntos en color azul, el plano que pasa por los dos puntos y el origen y su intersección con la superficie esférica de radio unidad, línea gruesa de color rojo

Eliminaremos el plano y dibujaremos los puntos P1 y P2 y la circunferencia que pasa por dichos puntos

Los caso particulares más simples se producen cuando los puntos P1 y P2 están en el mismo meridiano φ=0, la circunferencia máxima es el meridiano. O cuando ambos están en el ecuador, θ=π/2

Otro caso particular, se produce cuando los P1 y P2 están alineados con el origen. P1 y P2 son antípodas. Hay infinitas circunferencias máximas que pasan por los dos puntos

phi=linspace(0,pi,30);
theta=linspace(0,2*pi,40);
[phi,theta]=meshgrid(phi,theta);
x=sin(phi).*cos(theta);
y=sin(phi).*sin(theta);
z=cos(phi);
hold on
h1=mesh(x,y,z);
set(h1,'EdgeColor',[0.6,0.6,0.6]);

%punto P1
th=pi/4;
phi=0;
x1=cos(phi)*sin(th);
y1=sin(phi)*sin(th);
z1=cos(th);
plot3(x1,y1,z1,'o','markersize',4,'markeredgecolor','b','markerfacecolor','b')

%punto P2
th=3*pi/4;
phi=pi;
x2=cos(phi)*sin(th);
y2=sin(phi)*sin(th);
z2=cos(th);
plot3(x2,y2,z2,'o','markersize',4,'markeredgecolor','b','markerfacecolor','b')

%plano
A=atan2((y1*z2-y2*z1),(x1*z2-x2*z1));
B=(x1*sin(A)-y1*cos(A))/z1;
%circunferencia máxima
phi=(0:360)*pi/180;
th=atan2(B,sin(A-phi));
plot3(cos(phi).*sin(th),sin(phi).*sin(th),cos(th), 'r','lineWidth',1.5)

axis equal
axis off
view(100,5)
hold off
xlabel('x'); ylabel('y'); zlabel('z')
title('Superficie esférica')

### Longitud del arco que une dos puntos

Los vectores que unen el origen con los puntos P1 y P2 de la superficie esférica, son

$\begin{array}{l}\stackrel{\to }{O{P}_{1}}=\mathrm{cos}{\phi }_{1}\mathrm{sin}{\theta }_{1}·\stackrel{^}{i}+\mathrm{sin}{\phi }_{1}\mathrm{sin}{\theta }_{1}·\stackrel{^}{j}+\mathrm{cos}{\theta }_{1}·\stackrel{^}{k}\\ \stackrel{\to }{O{P}_{2}}=\mathrm{cos}{\phi }_{2}\mathrm{sin}{\theta }_{2}·\stackrel{^}{i}+\mathrm{sin}{\phi }_{2}\mathrm{sin}{\theta }_{2}·\stackrel{^}{j}+\mathrm{cos}{\theta }_{2}·\stackrel{^}{k}\end{array}$

El producto escalar de estos dos vectores de módulo unidad es

$\stackrel{\to }{O{P}_{1}}·\stackrel{\to }{O{P}_{2}}=\mathrm{cos}\alpha =\mathrm{sin}{\theta }_{1}\mathrm{sin}{\theta }_{2}\mathrm{cos}\left({\phi }_{1}-{\phi }_{2}\right)+\mathrm{cos}{\theta }_{1}\mathrm{cos}{\theta }_{2}$

Conocido el ángulo α, la longitud del arco es s=Rα. Siendo R el radio de la superficie terrestre

## Superficie cilíndrica

La ecuación de una superficie lateral cilíndrica de radio r es

$\left\{\begin{array}{l}x=r\mathrm{cos}\phi \\ y=r\mathrm{sin}\phi \\ z=z\end{array}$

Para esta superficie, u es z y v es φ

$\begin{array}{l}P={\left(\frac{\partial x}{\partial z}\right)}^{2}+{\left(\frac{\partial y}{\partial z}\right)}^{2}+{\left(\frac{\partial z}{\partial z}\right)}^{2}=1\\ R={\left(\frac{\partial x}{\partial \phi }\right)}^{2}+{\left(\frac{\partial y}{\partial \phi }\right)}^{2}+{\left(\frac{\partial z}{\partial \phi }\right)}^{2}={r}^{2}{\mathrm{sin}}^{2}\phi +{r}^{2}{\mathrm{cos}}^{2}\phi ={r}^{2}\\ Q=\frac{\partial x}{\partial z}\frac{\partial x}{\partial \phi }+\frac{\partial y}{\partial z}\frac{\partial y}{\partial \phi }+\frac{\partial z}{\partial z}\frac{\partial z}{\partial \phi }=0\end{array}$

La longitud de la línea que une dos puntos A y B de la superficie lateral cilíndrica y que está contenida en dicha superficie es

$L=\underset{A}{\overset{B}{\int }}\sqrt{1+{r}^{2}{\left(\frac{d\theta }{dz}\right)}^{2}}·dz$

Q=0, P y R no son funciones de v=φ. Estamos en el primer caso con u=z

$\begin{array}{l}\phi ={c}_{1}\int \frac{1}{\sqrt{{r}^{2}-{c}_{1}^{2}r}}dz+{c}_{2}\\ \phi ={c}_{3}z+{c}_{2}\end{array}$

Dados dos puntos P1 y P2 de la superficie lateral cilíndrica, calculamos las constantes c2 y c3. Dibujamos una superficie cilíndrica de radio r=1 y la espiral φ=0.05z

r=1*ones(30,1);
phi=linspace(0,2*pi,30);
[r,phi]=meshgrid(r,phi);
x=r.*cos(phi);
y=r.*sin(phi);
z=repmat(linspace(0,2,30),30,1);
hold on
surfl(x,y,z);
colormap(gray);
%espiral
th=(0:10:6*360)*pi/180;
plot3(cos(th),sin(th),0.05*th, 'r','lineWidth',1.5)

hold off
axis equal
axis off
xlabel('x'); ylabel('y'); zlabel('z')
title('Superficie cilíndrica')
view(60,15)

## Superficie cónica

La ecuación de una superficie cónica de ángulo α es

$\left\{\begin{array}{l}x=r\mathrm{sin}\alpha \mathrm{cos}\phi \\ y=r\mathrm{sin}\alpha \mathrm{sin}\phi \\ z=r\mathrm{cos}\alpha \end{array}$

Para esta superficie, u es r y v es φ

$\begin{array}{l}P={\left(\frac{\partial x}{\partial r}\right)}^{2}+{\left(\frac{\partial y}{\partial r}\right)}^{2}+{\left(\frac{\partial z}{\partial r}\right)}^{2}={\mathrm{sin}}^{2}\alpha {\mathrm{cos}}^{2}\phi +{\mathrm{sin}}^{2}\alpha {\mathrm{sin}}^{2}\phi +{\mathrm{cos}}^{2}\alpha =1\\ R={\left(\frac{\partial x}{\partial \phi }\right)}^{2}+{\left(\frac{\partial y}{\partial \phi }\right)}^{2}+{\left(\frac{\partial z}{\partial \phi }\right)}^{2}={r}^{2}{\mathrm{sin}}^{2}\phi {\mathrm{sin}}^{2}\alpha +{r}^{2}{\mathrm{cos}}^{2}\phi {\mathrm{sin}}^{2}\alpha ={r}^{2}{\mathrm{sin}}^{2}\alpha \\ Q=\frac{\partial x}{\partial r}\frac{\partial x}{\partial \phi }+\frac{\partial y}{\partial r}\frac{\partial y}{\partial \phi }+\frac{\partial z}{\partial r}\frac{\partial z}{\partial \phi }=-\mathrm{sin}\alpha \mathrm{cos}\phi ·r\mathrm{sin}\alpha \mathrm{sin}\phi +\mathrm{sin}\alpha \mathrm{sin}\phi ·r\mathrm{sin}\alpha \mathrm{cos}\phi =0\end{array}$

La longitud de la línea que une dos puntos A y B de la superficie cónica y que está contenida en dicha superficie es

$L=\underset{A}{\overset{B}{\int }}\sqrt{1+{r}^{2}{\mathrm{sin}}^{2}\alpha {\left(\frac{d\phi }{dr}\right)}^{2}}·dr$

Q=0, P y R no son funciones de v=φ. Estamos en el primer caso con u=r

$\phi ={c}_{1}\int \frac{1}{\sqrt{{r}^{4}{\mathrm{sin}}^{4}\alpha -{c}_{1}^{2}{r}^{2}{\mathrm{sin}}^{2}\alpha }}dr+{c}_{2}$

Resolvemos la integral haciendo el cambio de variable, s=rsinα/c1

$\int \frac{dr}{r\mathrm{sin}\alpha \sqrt{\frac{{r}^{2}{\mathrm{sin}}^{2}\alpha }{{c}_{1}^{2}}-1}}=\frac{1}{\mathrm{sin}\alpha }\int \frac{ds}{s\sqrt{{s}^{2}-1}}$

Hacemos el cambio, s=sect, s=1/cost, ds=sint/cos2t

$\frac{1}{\mathrm{sin}\alpha }\int dt=\frac{1}{\mathrm{sin}\alpha }t$

Deshacemos los cambios y despejamos r

$\begin{array}{l}\frac{1}{\mathrm{sin}\alpha }\text{arcsec(s)}\\ \phi =\frac{1}{\mathrm{sin}\alpha }\text{arcsec}\left(\frac{r\mathrm{sin}\alpha }{{c}_{1}}\right)-{c}_{2}\\ \frac{r\mathrm{sin}\alpha }{{c}_{1}}=\frac{1}{\mathrm{cos}\left(\phi \mathrm{sin}\alpha +{c}_{2}\right)}\\ r=\frac{{c}_{1}}{\mathrm{sin}\alpha ·\mathrm{cos}\left(\phi \mathrm{sin}\alpha +{c}_{2}\right)}\end{array}$

Para una geodésica que pase por los puntos P1 (φ1=0, r1) y P2 (φ2, r2), las constantes c1 y c2, valen

$\begin{array}{l}{r}_{1}=\frac{{c}_{1}}{\mathrm{sin}\alpha ·\mathrm{cos}{c}_{2}}\\ {r}_{2}=\frac{{c}_{1}}{\mathrm{sin}\alpha ·\mathrm{cos}\left({\phi }_{2}\mathrm{sin}\alpha +{c}_{2}\right)}\\ \frac{{r}_{1}}{{r}_{2}}=\frac{\mathrm{cos}\left({\phi }_{2}\mathrm{sin}\alpha +{c}_{2}\right)}{\mathrm{cos}{c}_{2}}\\ \mathrm{tan}{c}_{2}=\frac{\mathrm{cos}\left({\phi }_{2}\mathrm{sin}\alpha \right)-\frac{{r}_{1}}{{r}_{2}}}{\mathrm{sin}\left({\phi }_{2}\mathrm{sin}\alpha \right)}\\ {c}_{1}={r}_{1}\mathrm{sin}\alpha ·\mathrm{cos}{c}_{2}\end{array}$

Ejemplo

• Angulo del cono, α=π/3 (60°)
• Punto P1: r1=1, φ1=0
• Punto P2: r2=1, φ2=5π/6 (150°)
alfa=pi/3;
r=linspace(0,1,30);
phi=linspace(0,2*pi,30);
[r,phi]=meshgrid(r,phi);
x=r.*cos(phi)*sin(alfa);
y=r.*sin(phi)*sin(alfa);
z=r*cos(alfa);
hold on
surfl(x,y,z);
colormap(gray);

%puntos
r1=1;
r2=1;
th_2=5*pi/6;
c2=atan2((cos(th_2*sin(alfa))-r1/r2),sin(th_2*sin(alfa)));
c1=sin(alfa)*cos(c2)*r1;
%geodésica
th=linspace(0,th_2,100);
r=c1./(sin(alfa)*cos(th*sin(alfa)+c2));
x=r.*cos(th)*sin(alfa);
y=r.*sin(th)*sin(alfa);
z=r*cos(alfa);
plot3(x,y,z,'r','lineWidth',1.5)
hold off
axis off
xlabel('x'); ylabel('y'); zlabel('z')
title('Superficie cónica')
view(47,50)

Comprobamos que para obtener la geodésica, el ángulo φ2<π/sinα. Si disminuimos el ángulo α del cono podemos aumentar el ángulo φ2 por encima de 2π, como en este ejemplo

• Angulo del cono, α=π/7 (25.7°)
• Punto P1: r1=1, φ1=0
• Punto P2: r2=1, φ2=6.46 (370°)

## Referencias

J Villanueva. Geodesics on surfaces by Variational Calculus. Florida Memorial University, Miami