Indice | Previo: Fundamentos.Conicas.Areas | Siguiente: Fundamentos.ProblemasSeleccionados
4.2.16. Cónicas en el espacio¶
En las secciones anteriores desarrollamos todas las posibles representaciones algebraicas de las curvas cónicas sobre el plano en el que están definidas. En esta sesión daremos el salto a tres dimensiones y resolveremos la pregunta de ¿cuál es la representación algebraica o geométrica más general de las cónicas en el espacio?
En la Sección Rotación de las cónicas en el plano habíamos visto que es posible, partiendo de la descripción algebraica de la cónica en un sistema de coordenadas referido a su eje de simetría y que llamaremos en lo sucesivo el sistema de referencia natural de la cónica, aplicar una rotación sobre el plano para obtener la ecuación más general de la cónica sobre ese mismo plano. El salto al espacio de tres dimensiones es simplemente una generalización de este procedimiento.
4.2.17. Ángulos de Euler¶
Partiendo del sistema de referencia natural de la cónica, es posible orientar de forma arbitraria la curva en el espacio realizando en total tres rotaciones independientes (ver Figura (angulos_euler).) Los ángulos en los que se realizan esas rotaciones, y que llamaremos en este texto $(\Omega,i,\omega)$, se conocen universalmente como los ángulos de Euler.
La secuencia de rotaciones mostradas en la Figura (aungulos_euler) se puede describir cualitativamente como:
Rotación del sistema original $x-y-z$ (sistema natural) en un ángulo $\Omega$ alrededor del eje z, para obtener un nuevo sistema de ejes $x'-y'-z'$
Rotación del sistema $x'-y'-z'$ en un ángulo $i$ alrededor del eje $z'$, para obtener un nuevo sistema de ejes $x''-y''-z''$
Rotación del sistema $x''-y''-z''$ en un ángulo $\omega$ alrededor del eje $z''$, para obtener un nuevo sistema de ejes $x'''-y'''-z'''$.
Al sistema de ejes final lo llamamos el sistema de coordenadas del observador.
Usando la representación matricial de las rotaciones en el plano de la Ec. (rotacion2d_matricial), la relación entre las coordenadas del observador y las coordenadas naturales de las cónicas será:
\begin{equation} \label{eq:rotacion3d_observador_R} \left( \begin{array}{c} x''' \\ y''' \\ z''' \end{array} \right) = R_z(\omega) R_x(i) R_z(\Omega) \left( \begin{array}{c} x \\ y \\ z \end{array} \right) \end{equation}Si usamos la definición de las matrices de rotación de la Ec. (matriz_rotacion2d):
$$ \label{matrix-inverted} \left( \begin{array}{c} x''' \\ y''' \\ z''' \end{array} \right)= \left( \begin{array}{ccc} \mathrm{c} \omega & \mathrm{s} \omega & 0\\ -\mathrm{s} \omega & \mathrm{c} \omega & 0\\ 0 & 0 & 1 \end{array} \right) \left( \begin{array}{ccc} 1 & 0 & 0\\ 0 & \mathrm{c} i & \mathrm{s} i\\ 0 & -\mathrm{s} i & \mathrm{c} i \end{array} \right) \left( \begin{array}{ccc} \mathrm{c} \Omega & \mathrm{s} \Omega & 0\\ -\mathrm{s} \Omega & \mathrm{c} \Omega & 0\\ 0 & 0 & 1 \end{array} \right) \left( \begin{array}{c} x \\ y \\ z \end{array} \right) $$donde se ha abreviado $\mathrm{c}\theta\equiv\cos\theta$ y $\mathrm{s}\theta\equiv\sin\theta$.
Al realizar las multiplicaciones matriciales explícitas queda:
\begin{equation} \label{eq:rotacion3d_observador_M} \left( \begin{array}{c} x''' \\ y''' \\ z''' \end{array} \right)= {\cal M}(\Omega,i,\omega) \left( \begin{array}{c} x \\ y \\ z \end{array} \right) \end{equation}donde la matriz de rotación en tres dimensiones,
\begin{equation} \label{eq:matrizM_definicion} {\cal M}(\omega,i,\Omega)\equiv R_z(\omega) R_x(i) R_z(\Omega) \end{equation}se escribe explícitamente como:
\begin{equation} \label{eq:matrizM_explicita} {\cal M}(\omega,i,\Omega)= \left( \begin{array}{ccc} \mathrm{c} \omega\;\mathrm{c} \Omega - \mathrm{c} i\;\mathrm{s} \omega\; \mathrm{s} \Omega & \mathrm{c} \omega\; \mathrm{s} \Omega + \mathrm{c} \Omega\; \mathrm{c} i\; \mathrm{s} \omega & \mathrm{s} i\; \mathrm{s} \omega \\ -\mathrm{c} \Omega\; \mathrm{s} \omega - \mathrm{s} \Omega\; \mathrm{c} i\; \mathrm{c} \omega & -\mathrm{s} \Omega\; \mathrm{s} \omega + \mathrm{c} \Omega\; \mathrm{c} i\; \mathrm{c} \omega & \mathrm{s} i\; \mathrm{c} \Omega \\ \mathrm{s} \Omega\; \mathrm{s} i & -\mathrm{c} \Omega\; \mathrm{s} i & \mathrm{c} i \\ \end{array} \right) \end{equation}Por ser ${\cal M}$ el producto de matrices de unitarias (Ec. matrizM_definicion), ella es en sí misma una matriz unitaria, es decir $\mathrm det{\cal M}=1$, pero más importante:
$$ {\cal M}^{-1}={\cal M}^\mathrm{T} $$explícitamente,
\begin{equation} \label{eq:matrizM_explicita_transpuesta} {\cal M}(\omega,i,\Omega)^\mathrm{T}= \left( \begin{array}{ccc} \mathrm{c} \omega \;\mathrm{c} \Omega - \mathrm{c} i \;\mathrm{s} \Omega \;\mathrm{s} \omega & -\mathrm{c} \Omega \;\mathrm{s} \omega - \;\mathrm{c} \omega \;\mathrm{c} i \;\mathrm{s} \Omega & \mathrm{s} i \;\mathrm{s} \Omega \\ \mathrm{c} \omega \;\mathrm{s} \Omega + \mathrm{s} \omega \;\mathrm{c} i \;\mathrm{c} \Omega & -\mathrm{s} \Omega \;\mathrm{s} \omega + \mathrm{c} \omega \;\mathrm{c} i \;\mathrm{c} \Omega & \mathrm{s} i \;\mathrm{c} \omega \\ \mathrm{s} \omega \;\mathrm{s} i & \mathrm{c} \omega \;\mathrm{s} i & \mathrm{c} i \\ \end{array} \right) \end{equation}A la inversa, las coordenadas naturales de la cónica $(x,y,z)$ se pueden obtener de las coordenadas del observador, invirtiendo la Ec. (rotacion3d_observador):
\begin{equation} \label{eq:rotacion3d_natural_R} \left( \begin{array}{c} x \\ y \\ z \end{array} \right) = R_z(-\Omega) R_x(-i) R_y(-\omega) \left( \begin{array}{c} x''' \\ y''' \\ z''' \end{array} \right) \end{equation}o equivalentemente, usando la forma explícita en la Ec. (rotacion3d_observador_M) y aprovechando las propiedades de la matriz de rotación en tres dimensiones $\cal M$:
\begin{equation} \label{eq:rotacion3d_natural_M} \left( \begin{array}{c} x \\ y \\ z \end{array} \right)= {\cal M}^\mathrm{T} \left( \begin{array}{c} x''' \\ y''' \\ z''' \end{array} \right) \end{equation}4.2.18. Matrices de rotación generales¶
Usando SPICE podemos construir la matriz de rotación en tres dimensiones ${\cal M}$ por dos medios distintos. El primero es usar la rutina rotate que habíamos introducido en el Alg. (spice_rotate):
#Angulos
from numpy import pi
Omega=pi/6
omega=pi/3
i=pi/4
#Matrices indiduales
from spiceypy import rotate
RzOmega=rotate(Omega,3)
Rxi=rotate(i,1)
Rzomega=rotate(omega,3)
La matriz de rotación en tres dimensiones ${\cal M}$, de acuerdo con su definición en la Ec. (matrizM_definicion), puede obtenerse aplicando sucesivamente la rutina de multiplicación de matrices mxm:
from spiceypy import mxm
M=mxm(Rzomega,mxm(Rxi,RzOmega))
Existe sin embargo una rutina compacta y mucho más general en SPICE que permite calcular la matriz rotación de cualquier sucesión de rotaciones en el espacio, no solamente la que usamos aquí.
Para ello podemos definir una matriz de rotación general ${\cal R}$ resultante de aplicar rotaciones sucesivas $\theta_i$ alrededor del eje $\hat{e}_k$, $\theta_j$ alrededor del eje $\hat{e}_j$ y $\theta_i$ alrededor del eje $\hat{e}_i$ como:
$$ {\cal R}(\theta_i,\theta_j,\theta_k,i,j,k)\equiv R_{i}(\theta_i)R_{j}(\theta_j)R_{k}(\theta_k) $$En términos de esta matriz general, la matriz de rotación canónica ${\cal M}$ presentada en la Sección Ángulos de Euler se escribe:
$$ {\cal M}={\cal R}(\omega,i,\Omega,3,1,3) $$En el paquete SPICE la rotación general ${\cal R}$ se implementa con la rutina eul2m que podemos usar para obtener $\cal M$ como:
from spiceypy import eul2m
M=eul2m(omega,i,Omega,3,1,3)
Naturalmente el resultado coincide con el obtenido en el Alg. (calculo_M_explicito). Podemos ahora verificar la unitariedad de $\cal M$ calculando su determinante, su inversa (calculada numérica y usando la propiedad en la Ec. rotacion3d_natural_R) y su transpuesta:
#Determinante
from numpy.linalg import det
detM=det(M)
#Inversa
from numpy.linalg import inv
Minv=inv(M)
#Inversa por definicion
Minv_def=eul2m(-Omega,-i,-omega,3,1,3)
#Transpuesta
MT=M.transpose()
print(f"det(M) = {detM:g}")
print(f"inversa M (numérica) = \n{Minv}")
print(f"inversa M (definición) = \n{Minv_def}")
print(f"transpuesta M = \n{MT}")
4.2.19. Gráfico de una cónica rotada en el espacio¶
Con estos elementos a la mano podemos escribir algoritmos para, usando las ecuaciones encontradas en la Sección Conicas en coordenadas cilíndricas y los algoritmos de la Sección Gráfico de una cónica en el plano, representar gráficamente una cónica en el espacio.
Comenzamos por obtener los puntos de la cónica en su sistema natural de coordenadas (plano de la cónica, origen en el foco, semieje $x+$ apuntando hacia el periapsis) usando para ello la rutina puntos_conica que habíamos introducido en el Alg. (puntos_conica):
from pymcel.export import puntos_conica
p=10.0
e=0.8
xppps,yppps,zppps=puntos_conica(p,e)
Ahora podemos rotar los puntos de la cónica usando la matriz M calculada en el Alg. (calculo_M_general). Para hacerlo nos apoyaremos nuevamente en la rutina general rota_puntos que habíamos introducido antes en el Alg. (\ref{cod:rota_puntos}):
from pymcel.export import rota_puntos
xs,ys,zs=rota_puntos(M,xppps,yppps,zppps)
Una gráfica en tres dimensiones de la cónica en su plano natural y en los ejes rotados se puede realizar usando el siguiente código:
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
fig=plt.figure()
ax=fig.gca(projection='3d')
#Gráfica de los puntos originales
ax.plot(xs,ys,zs,'k--')
#Gráfica de la cónica rotada
ax.plot(xppps,yppps,zppps,'b-')
#Decoración
ax.set_xlabel("$x$")
ax.set_ylabel("$y$")
ax.set_zlabel("$z$")
from pymcel.plot import fija_ejes3d_proporcionales
fija_ejes3d_proporcionales(ax);
fig.tight_layout();
4.2.2. Elementos orbitales¶
Partiendo de las ecuaciones paramétricas de la cónica en su plano natural (Ecs. conica_parametricas_f y conica_ecuacion_cilindricas):
\begin{eqnarray} \nonumber x''' & = & \frac{p\cos f}{1+e\cos f}\\ \nonumber y''' & = & \frac{p\sin f}{1+e\cos f}\\ \nonumber z''' & = & 0 \end{eqnarray}y usando las expresiones explícitas para la rotación en tres dimensiones dadas por las Ecs. (rotacion3d_natural_M) y (matrizM_explicita_transpuesta), podemos escribir las ecuaciones paramétricas generales de una cónica en el espacio:
\begin{equation} \label{eq:elementos_estado_f} \begin{array}{rcl} x & = & r[\cos \Omega \cos(\omega+f) - \cos i \sin \Omega \sin(\omega+f)]\\ y & = & r[\sin \Omega \cos(\omega+f) + \cos i \cos \Omega \sin(\omega+f)]\\ z & = & r[\cos f\sin \omega \sin i + \sin f \cos \omega \sin i]\\ \end{array} \end{equation}donde $r=p/(1+e\cos f)$
Vemos aquí entonces que para específicar la posición de cualquier punto sobre una cónica, independiente de su orientación espacial, hace falta indicar el valor de 6 parámetros: $p$, $e$, $i$, $\Omega$, $\omega$ y $f$. A estas cantidades las llamamos en mecánica celeste los elementos orbitales clásicos y volveremos sobre ellos en el Capítulo El problema de los dos cuerpos.
En realidad, de los 6 elementos orbitales clásicos, 5 de ellos $(p,e,i,\Omega,\omega)$ permiten especificar el tamaño, forma y orientación de la cónica y son compartidos por todos los puntos que definen la curva. El último elemento, $f$ permite especificar la posición de un punto específico.
Usando los elementos orbitales clásicos podemos dibujar una cónica en el espacio usando un algoritmo más directo que el que usamos en la Sección Gráfico de una cónica rotada en el espacio. Para ello hemos diseñado la rutina conica_de_elementos:
def conica_de_elementos(p=10.0,e=0.8,i=0.0,Omega=0.0,omega=0.0,
df=0.1,
elev=30,azim=60,
figreturn=False):
#Convierte elementos angulares en radianes
from numpy import pi
p=float(p)
e=float(e)
i=float(i)*pi/180
Omega=float(Omega)*pi/180
omega=float(omega)*pi/180
#Compute fmin,fmax
if e<1:
fmin=-pi
fmax=pi
elif e>1:
from numpy import arccos
psi=arccos(1/e)
fmin=-pi+psi+df
fmax=pi-psi-df
else:
fmin=-pi+df
fmax=pi-df
#Valores del ángulo
from numpy import linspace,pi
fs=linspace(fmin,fmax,500)
#Distancia al periapsis
q=p/(1+e)
#Distancia al foco
from numpy import sin,cos
rs=p/(1+e*cos(fs))
#Coordenadas
xs=rs*(cos(Omega)*cos(omega+fs)-cos(i)*sin(Omega)*sin(omega+fs))
ys=rs*(sin(Omega)*cos(omega+fs)+cos(i)*cos(Omega)*sin(omega+fs))
zs=rs*(cos(fs)*sin(omega)*sin(i)+sin(fs)*cos(omega)*sin(i))
#Posición del periapsis (f=0)
xp=q*(cos(Omega)*cos(omega)-cos(i)*sin(Omega)*sin(omega))
yp=q*(sin(Omega)*cos(omega)+cos(i)*cos(Omega)*sin(omega))
zp=q*sin(omega)*sin(i)
#Posición del nodo ascendente
rn=p/(1+e*cos(omega))
xn=rn*cos(Omega)
yn=rn*sin(Omega)
zn=0
#Gráfico
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
plt.close("all")
fig=plt.figure()
ax=fig.gca(projection='3d')
#Gráfica de los puntos originales
ax.plot(xs,ys,zs,'b-')
#Posición del periapsis
ax.plot([0,xp],[0,yp],[0,zp],'r-')
#Posición del nodo ascendente
ax.plot([0,xn],[0,yn],[0,zn],'g-')
#Fija punto de vista
ax.view_init(elev=elev,azim=azim)
#Decoración
from pymcel.plot import fija_ejes3d_proporcionales
xrange,yrange,zrange=fija_ejes3d_proporcionales(ax);
ax.set_title(f"Cónica con:"+\
f"$p={p:.2f}$, $e={e:.2f}$, "+\
f"$i={i*180/pi:.2f}$, "+\
f"$\Omega={Omega*180/pi:.1f}$, "+\
f"$\omega={Omega*180/pi:.1f}$"
)
#Dibuja Ejes
ax.plot([0,xrange[1]],[0,0],[0,0],'k-')
ax.plot([0,0],[0,yrange[1]],[0,0],'k-')
ax.plot([0,0],[0,0],[0,zrange[1]],'k-')
ax.text(xrange[1],0,0,"$x$",ha='left',va='top')
ax.text(0,yrange[1],0,"$y$",ha='left',va='top')
ax.text(0,0,zrange[1],"$z$",ha='left',va='bottom')
fig.tight_layout();
if figreturn:return fig
El gráfico de la cónica sería:
from numpy import pi
conica_de_elementos(p=10.0,e=0.5,Omega=60.0,i=30.0,omega=20.0);
Para ver una versión interactiva de esta gráfica vaya a las libertas disponibles en la versión electrónica el libro.
Indice | Previo: Fundamentos.Conicas.Areas | Siguiente: Fundamentos.ProblemasSeleccionados