Full text
Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Trabajo de Fin de Grado Ingeniería Aeroespacial Estudio de los puntos de cruce para órbitas de pequeña excentricidad y con traza repetida Autor: Alejandro Montero Miñán Tutor: Rafael Vázquez Valenzuela Dpto. Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2018
Trabajo de Fin de Grado Ingeniería Aeroespacial Estudio de los puntos de cruce para órbitas de pequeña excentricidad y con traza repetida Autor: Alejandro Montero Miñán Tutor: Rafael Vázquez Valenzuela Profesor Titular Dpto. Ingeniería Aeroespacial Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2018
Trabajo de Fin de Grado: Estudio de los puntos de cruce para órbitas de pequeña excentricidad y con traza repetida Autor: Alejandro Montero Miñán Tutor: Rafael Vázquez Valenzuela El tribunal nombrado para juzgar el trabajo arriba indicado, compuesto por los siguientes profesores: Presidente: Vocal/es: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha:
Resumen Uno de los mayores problemas a la hora de controlar y recuperar los datos obtenidos de un satélite se debe a que el número de oportunidades viene limitado por la localización del centro desde donde se recuperan y de los parámetros de la propia órbita. Por esta razón, las órbitas con traza repetida, las cuales tienen la peculiaridad de que, tras un número de revoluciones, vuelven a sobrevolar los mismos puntos de la Tierra, son muy ventajosas, ya que permiten controlar de forma precisa los puntos de observación y recuperación de datos. Además, este tipo de trazas pueden presentar puntos de cruce, es decir, puntos que el satélite sobrevuela dos veces antes de que la traza se repita. Estos puntos son de especial interés para la colocación de estos centros de recuperación de datos y observación de satélites, ya que permiten duplicar el número de ventanas de observación. Este documento se encarga de caracterizar con el mayor grado de detalle posible cómo los distintos elementos orbitales afectan a la disposición y cantidad de estos puntos de cruce en órbitas geocéntricas de excentricidad pequeña. Para ello, primero se desarrollan brevemente algunos conceptos necesarios de Mecánica orbital, además de realizar un pequeño resumen del caso circular ya estudiado en otros documentos para que sirva de introducción. Tras esto, se obtendrán las ecuaciones que permiten obtener la posición de los puntos de corte, se explicará su resolución y se analizará cómo afecta la variación de algunos elementos orbitales, por ejemplo, la excentricidad, el argumento de perigeo o la inclinación. Después, se ilustrará con ejemplos la precisión y las limitaciones del modelo utilizado para resolver el problema. Por último, se incluirán algunas conclusiones generales y se discutirán algunas líneas de trabajo que podrían tomarse en un futuro a partir de este trabajo. I
Abstract One of the greatest problems when monitoring and retrieving data obtained from a satellite is that the number of actual chances is limited by the location of the retrieving site and the parameters of the orbit. For this reason, repeated groundtrack orbits, whose main characteristic is the repetition of the points that this satellite flies by after a given number of complete periods, are really important, as they allow a precise control of the observing and retrieving points. Also, this type of groundtracks can have crossover points, meaning, points by which the satellite flies by twice before the groundtrack repeats. These points are very interesting for the location of an observation site, because, as mentioned before, they provide twice the number of observation windows. This document tries to characterize as deeply as possible how the different orbital elements influence the position and number of these crossover points in geocentric low excentricity orbits. For this, first various orbital mechanics’concepts are faintly described, followed by a short sum up of the previous studied circular case, which serves as a introduction. Then, equations to calculate the crossover points’location are obtained, their resolution process explained, as well as how the variations of some orbital elements such as inclination, excentricity or perigee argument affects the position or ammount of these crossover points. Also, the precision and limitations of the model used to solve this problem is discussed. Finally, some conlusions and some new research paths are included. III
XLista de imágenes 4.12 Zoom 1 de la gráfica 4.11 42 4.13 Zoom 2 de la gráfica 4.11 42 4.14 Zoom 3 de la gráfica 4.11 43 4.15 Zoom 4 de la gráfica 4.11 43 4.16 Zoom 5 de la gráfica 4.11 44 4.17 Evolución de la posición de los puntos de cruce para órbitas con k=7, m=3, i=85º y ω=−15 º 44 4.18 Zoom 1 de la gráfica 4.17 45 4.19 Zoom 2 de la gráfica 4.17 45 4.20 Zoom 3 de la gráfica 4.17 46 4.21 Evolución de la posición de los puntos de cruce para órbitas con k=3, m=1, i=80º y e=0.05 48 4.22 Zoom 1 de la gráfica 4.21 48 4.23 Zoom 2 de la gráfica 4.21 49 4.24 Zoom 3 de la gráfica 4.21 49 4.25 Zoom 2 de la gráfica 4.21 50 4.26 Evolución de la posición de los puntos de cruce para órbitas con k=3, m=2, i=85º y e=0.05 50 4.27 Zoom 1 de la gráfica 4.26 51 4.28 Zoom 2 de la gráfica 4.26 51 4.29 Evolución de la posición de los puntos de cruce para órbitas con k=3, m=2, i=85º y e=0.05 52 4.30 Zoom 1 de la gráfica 4.29 52 4.31 Zoom 2 de la gráfica 4.29 53 4.32 Zoom 3 de la gráfica 4.29 53 4.33 Zoom 4 de la gráfica 4.29 54 4.34 Zoom 5 de la gráfica 4.29 54 4.35 Zoom 6 de la gráfica 4.29 55 4.36 Evolución de la posición de los puntos de cruce para órbitas con k=7, m=3, i=85º y e=0.05 55 4.37 Zoom 1 de la gráfica 4.36 56 4.38 Zoom 2 de la gráfica 4.36 56 4.39 Zoom 3 de la gráfica 4.36 57 4.40 Zoom 4 de la gráfica 4.36 57 4.41 Zoom 5 de la gráfica 4.36 58 4.42 Variación de la posición del punto de corte situado en el hemisferio Norte con la excentricidad y el argumento del perigeo para órbitas con i=83º k=3 y m=1 59 4.43 Zoom de la gráfica 4.42 59 4.44 Variación de la posición de uno de los puntos de corte situados en el hemisferio Sur con la excentricidad y el argumento del perigeo para órbitas con i=87º k=3 y m=2 60 4.45 Variación de la posición de uno de los puntos de corte situados en el hemisferio Norte con la excentricidad y el argumento del perigeo para órbitas con i=85º k=3 y m=1 60 4.46 Variación de la posición del punto de corte situado sobre el Ecuador con la excentricidad y el argumento del perigeo para órbitas con i=85º k=7 y m=3 61 4.47 Zoom de la gráfica 4.46 61 4.48 Ejemplo de órbita con traza repetida cada 5 revoluciones (k=5) y cada tres días (m=3), de excentricidad e=0.05, argumento del perigeo ω=30 º e inclinación crítica. 62 4.49 Ejemplo de órbita con traza repetida cada 9 revoluciones (k=9) y cada 7 días (m=7), de excentricidad e=0.025, argumento del perigeo ω=70 º e inclinación crítica. 63 4.50 Primer par de inclinaciones críticas 64 4.51 Segundo par de inclinaciones críticas 65 4.52 Quinta inclinación crítica 65 4.53 Primer par de inclinaciones críticas 65 4.54 Segundo par de inclinaciones críticas 66 4.55 Primera inclinación crítica 66 4.56 Último par de inclinaciones críticas 66 4.57 Variación con la excentricidad de las inclinaciones críticas para órbitas del caso impar 68 4.58 Variación con la excentricidad de las inclinaciones críticas para órbitas del caso par con mimpar 68 4.59 Variación con la excentricidad de las inclinaciones críticas para órbitas del caso par con mpar 69 4.60 Variación con el argumento del perigeo de las inclinaciones críticas para órbitas del caso impar 69 4.61 Variación con el argumento del perigeo de las inclinaciones críticas para órbitas del caso impar (zoom) 69 4.62 Variación con el argumento del perigeo de las inclinaciones críticas para órbitas del caso par con mimpar 70
Lista de imágenes XI 4.63 Variación con el argumento del perigeo de las inclinaciones críticas para órbitas del caso par con mpar 70 4.64 Variación con el argumento del perigeo de las inclinaciones críticas para órbitas del caso par con mpar (zoom) 70 4.65 Ejemplo de tangencia en el hemisferio opuesto para órbita con k=7, m=5, e=0.03, segunda inclinación crítica y argumento del perigeo ω=−70 º. 72 4.66 Generación de nuevos puntos de corte para inclinaciones cercanas a cero con excentricidades pequeñas y kymgrandes 75 4.67 Función cos(i)−m kcos2(i∗)con m kgenérico 80 4.68 Función cos(i)−m kcos2(i∗)con m k≈180 4.69 Valor de excentricidad que hace que se anule el término Ccon u∗ 2=−3π 4yω=π 281 4.70 Comparación entre las variaciones de la latitud y longitud del punto de corte situado en el hemisferio Norte con k=3, m=1 e i=83º obtenidas por la ecuación aproximada y el sistema de ecuaciones completo. 84 4.71 Zoom de la gráfica 4.70 85 4.72 Comparación entre las variaciones de la latitud y longitud de uno de los puntos de corte situado en el hemisferio Sur con k=3, m=2 e i=87º obtenidas por la ecuación aproximada y el sistema de ecuaciones completo. 85 4.73 Zoom de la gráfica 4.72 86 4.74 Comparación entre las variaciones de la latitud y longitud de uno de los puntos de corte situado en el hemisferio Norte con k=6, m=5 e i=83º obtenidas por la ecuación aproximada y el sistema de ecuaciones completo. 86 4.75 Zoom de la gráfica 4.74 87 4.76 Comparación entre las variaciones de la latitud y longitud del punto de corte situado sobre el Ecuador con k=7, m=3 e i=85º obtenidas por la ecuación aproximada y el sistema de ecuaciones completo. 87 4.77 Zoom de la gráfica 4.76 88 4.78 Creación de nuevos puntos de corte para excentricidades altas en una órbita con e=0.9, k=3, m=1, i=0.1 rad y ω=0º 88 4.79 Zoom de la gráfica 4.78 89 4.80 Creación de nuevos puntos de corte para excentricidades altas en una órbita con e=0.77, k=7, m=3, i=0.1 rad y ω=0º 89 4.81 Zoom de la gráfica 4.80 90
Índice de Tablas 3.1 Soluciones de la ecuación (3.25) para órbita circular del caso impar con i=83º, k=5ym=322 3.2 Soluciones de la ecuación (3.37) para órbita circular del caso par con i=85º, k=4ym=324 4.1 Comparativa entre las soluciones obtenidas por la ecuación aproximada y el sistema de ecuaciones no lineales completo para una órbita cuyas características son k=3, m=2, ω=25 º, i=85º y e=0.15 33 4.2 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=3, m=1, ω=20 º, i=80º y e=0.05 35 4.3 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=5, m=1, ω=20 º, i=85º y e=0.05 35 4.4 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=7, m=6, ω=20 º, i=85º y e=0.05 35 4.5 Variación de la posición de algunos puntos de corte de la órbita con traza repetida cuyos elementos orbitales son k=3, m=1, i=80º y ω=50 º 37 4.6 Inclinaciones críticas para órbita con para k=7, m=5, e=0.03, ω=70 º 64 4.7 Inclinaciones críticas para órbita con para k=7, m=4, e=0.03, ω=70 º 65 4.8 Inclinaciones críticas para órbita con para k=8, m=3, e=0.05, ω=−20 º 66 4.9 Inclinaciones críticas para órbita con para k=7, m=5 y e=0 67 4.10 Inclinaciones críticas para órbita con para k=7, m=4 y e=0 67 4.11 Inclinaciones críticas para órbita con para k=8, m=3 y e=0 67 4.12 Inclinaciones críticas para órbita con para k=7, m=5, e=0.03, ω=70 º 82 4.13 Inclinaciones críticas para órbita con para k=7, m=4, e=0.03, ω=70 º 82 4.14 Inclinaciones críticas para órbita con para k=8, m=3, e=0.05, ω=−20 º 82 4.15 Inclinaciones críticas para órbita con para k=7, m=5 y e=0 82 4.16 Inclinaciones críticas para órbita con para k=7, m=4 y e=0 82 4.17 Inclinaciones críticas para órbita con para k=8, m=3 y e=0 82 4.18 Tabla comparativa para órbitas con k=3, m=2, i=83º y ω=35 º con e=0.25 83 4.19 Tabla comparativa para órbitas con k=3, m=2, i=83º y ω=35 º con e=0.35 83 4.20 Tabla comparativa para órbitas con k=3, m=2, i=83º y ω=35 º con e=0.45 83 XIII
Notación ~ RiVector de posición de la masa i miMasa i GConstante de gravitación universal ¨xSegunda derivada respecto al tiempo de la variable x µParámetro gravitacional de masa combinadas ∈Energía específica ~ hMomento cinético específico pParámetro de la cónica ~rVector de posición del satélite desde el foco no vacante θAnomalía verdadera rpRadio de periápsis raRadio de apoápsis a Semieje mayor de la elipse b Semieje menor de la elipse c Distancia del origen de coordenadas a uno de los focos de la elipse eExcentricidad de la órbita nVelocidad angular TPeriodo de la órbita ∆tTiempo trancurrido entre dos posiciones de la órbita MAnomalía media EAnomalía excéntrica Primer punto de Aries Nodo ascendente ωArgumento de perigeo/periápsis iInclinación de la órbita RAAN oΩAscensión recta del nodo ascendente ¯ ωLongitud del perigeo uArgumento de latitud φLatitud λLongitud λTLongitud verdadera GST0Ascensión recta del meridiano de Greenwich en t=0 ⊕Referencia a la Tierra SAT Referencia al satélite ∆λRetraso nodal k Número de revoluciones completas del satélite que tarda en repetirse la traza mNúmero de días sidéreos que tarda en repetirse la traza λuÁngulo utilizado para el cálculo de trazas XV
XVI Notación r Número de revoluciones necesarias para encontrar un punto de corte δi,jCoeficiente (i,j) de la delta de Kronecker d c Parte entera sin Función seno cos Función coseno tan Función tangente arcsin Función arco seno arc cos Función arco coseno arctan Función arco tangente ⩽Menor o igual ⩾Mayor o igual >Mayor <Menor ∂ ∂Derivada parcial rad.Radianes revs.Revoluciones Sol.Solución aprox.Aproximadamente N−SNorte-Sur
1 Introducción El uso de satélites artificiales para la investigación, la comunicación o incluso para misiones bélicas es ya un hecho cotidiano en nuestras vidas. Para diseñar la órbita un satélite y hacer que cumpla una cierta misión, una de las herramientas de mayor utilidad es la traza del satélite, que no es más que el conjunto de puntos que este sobrevuela a lo largo de su órbita. En concreto, cuando esta traza se cierra, se dice que el satélite se encuentra en una órbita de traza repetida. Este tipo concreto de trazas, puede presentar puntos de corte o ”crossover points” que son aquellos puntos de la superficie sobre los que pasa el satélite dos veces sin que su traza se haya repetido aún. Estos puntos son de especial interés para ciertas misiones, ya que suponen una doble oportunidad para la observación, para la toma de datos o para la transferencia de los mismos. Se incluyen a continuación un par de ejemplos de trazas repetidas en las que además pueden observarse ”crossover points” o puntos de cruce: Imagen 1.1 Ejemplo de órbita circular con traza repetida y puntos de cruce. En 1.1 puede verse cómo la traza del satélite (línea verde) se cierra sobre sí misma, de forma que sea cual sea el tiempo que lleve el satélite en órbita, se encontrará siempre sobre la misma curva (idealmente). Además, marcados en rojo pueden verse los denominados anteriormente puntos de cruce o ”crossover points”, que 1
2Capítulo 1. Introducción son el objetivo de estudio de este documento. En concreto, 1.1 corresponde a la traza de una órbita circular. El caso circular ha sido estudiado por Eva Tesón Muñoz en ”Estudio detallado de los puntos de cruce de las trazas de satélites geocéntricos” (ref. [1]). En este documento, se estudia cómo afectan los distintos parámetros de la órbita (elementos orbitales) a la localización y cantidad de estos puntos de cruce para órbitas circulares, tanto de forma analítica como mediante cálculo numérico. Además se realiza una breve introducción a este mismo estudio generalizado para órbitas excéntricas. El estudio completo del caso circular se encuentra resumido en el capítulo 3, que servirá como introducción al tema propio de este documento. La extensión de este estudio a órbitas excéntricas es posible ya que la introducción de esta excentricidad no evita que la traza pueda cerrarse sobre sí misma, y producir además puntos de cruce, tal y como puede verse en 1.2: Imagen 1.2 Ejemplo de órbita con excentricidad, traza repetida y puntos de cruce. Este documento constituye una extensión del estudio realizado en [1] para órbitas circulares, particularizado para órbitas geocéntricas de pequeña excentricidad (a lo largo del documento se detallará qué órbitas engloba este adjetivo ”pequeña”), en concreto, se estudia la localización de estos puntos de corte en función de los elementos orbitales del satélite y la influencia de la variación de estos parámetros sobre los puntos de corte. Al poderse calcular el número de puntos de corte y su situación sobre un mapa terrestre, es posible emplazar los centros de observación y recopilación de datos en estos puntos de corte, o diseñar una órbita de forma que ciertos puntos específicos de la superficie sean observados más de una vez por el satélite. La principal diferencia con respecto a [1] es que, debido a la complejidad de las ecuaciones al introducir la excentricidad, en este documento se va a hacer uso de una aproximación de la ecuación para excentricidad pequeña, la cual se explicará con detalle en el capítulo 4, de forma que los resultados deberán compararse con los obtenidos al resolver el problema real, para poder comprobar su exactitud. Además, en el caso circular no puede estudiarse la variación de la posición de estos puntos de corte respecto a la excentricidad o al argumento del perigeo, ya que el primer elemento orbital es nulo para el caso circular, y el segundo solo puede definirse en órbitas excéntricas. Todos estos conceptos, como el de traza o el de elemento orbital, se definirán con más detalle en el capítulo 2.
1.1 Antecedentes históricos 3 1.1 Antecedentes históricos Durante varios siglos, las personas han estado interesadas en el cielo y en los cuerpos que en él se encuentran. Esto los condujo a intentar estudiar el comportamiento de los mismos, y a elaborar teorías de cómo y por qué estos cuerpos celestes se movían de la forma que lo hacían. La primera teoría que se considera correcta a la hora de describir el movimiento de los cuerpos en el espacio fue la formulada por Kepler en el siglo XVII, la cual se compone principalmente de tres leyes, las cuales son la base de la mecánica orbital. A finales del mismo siglo, Newton desarrolla la teoría de gravitación universal, así como el cálculo diferencial, los cuales son de suma importancia para todos los desarrollos matemáticos realizados en este documento. Estos dos grandes acontecimientos forman la base para todos los desarrollos y estudios relacionados con la mecánica celeste posteriores. Además, durante la Guerra Fría, Estados Unidos y la Unión Soviética intentan hacer uso de estos conocimientos para poner en órbita el primer satélite artificial de la historia. El primero en conseguir este gran hito fue la URSS, cuando el 4 de octubre de 1957 puso en órbita el Sputnik-1. Desde ese momento, el número de satélites artificiales que orbitan alrededor de la Tierra ha crecido exponencialmente, así como sus campos de aplicación. Para poder controlar, monitorizar y corregir las trayectorias de estos satélites artificiales, se desarrolló la Astrodinámica, la ciencia que aplica la mecánica celeste a los satélites artificiales, la cual es la base para los estudios realizados en este documento. 1.2 Descripción y clasificación de satélites Un satélite es todo aquel cuerpo que orbita alrededor de un planeta. En concreto, este documento se centra en el estudio de los satélites artificiales, y más concretamente, de los geocéntricos, es decir, de aquellos que orbitan alrededor de la Tierra (aunque puede realizarse fácilmente una extrapolación a cualquier otro planeta tomando los parámetros correspondientes al mismo). Actualmente, los satélites tienen formas muy variadas para poder desempeñar sus funciones concretas. Van desde telescopios de más de 13 toneladas hasta micro-satélites enfocados a la toma de muestras o a la transmisión de información. Además, a excepción de estos micro-satélites que están entrando ahora en auge, la producción en serie de satélites completamente similares está muy limitada, debido al elevado coste de lanzamiento. Por esta misma razón, usualmente se han construido satélites grandes que sean capaces de desarrollar varias tareas, y su diseño se hacía de forma casi individual. La clasificación principal de los satélites se basa en las características de su órbita. Debido a la amplia variedad de órbitas existentes, se describen aquí únicamente las más utilizadas. •Órbitas LEO (Low Earth Orbit): son aquellas que no se extienden más allá de los 2.000 km. Además, se suele evitar las órbitas a alturas inferiores a 300 km debido a los problemas que genera la resistencia atmosférica. Estas órbitas son las más próximas a la Tierra, por lo que son aquellas para las cuales es necesario una menor energía para poner en órbita al satélite. Por otra parte, debido a la atmósfera residual y a la forma no esférica de la Tierra, este tipo de órbita sufre perturbaciones, haciendo necesario el uso de combustible en el satélite para mantener la órbita nominal. Un ejemplo de cuerpo artificial que orbita en este tipo de órbita es la Estación Espacial Internacional. •Órbitas geoestacionarias GEO: son aquellas cuya principal característica es que el satélite siempre sobrevuela el mismo punto terrestre. Esto se consigue al hacer que el satélite tenga una órbita síncrona, la cual es un tipo de órbita cuyo periodo coincide con el de rotación terrestre, circular y ecuatorial. Este tipo de órbitas son muy utilizadas por satélites de comunicaciones, ya que permanecen ”fijos” sobre el mismo punto de la Tierra. La desventaja de este tipo de órbitas es que, para que el satélite tenga el mismo periodo que el de rotación de la Tierra, su distancia a la misma debe ser 35.786 km, lo que supone que es necesario una gran cantidad de energía para llevar el satélite hasta su órbita destino. Junto con el anterior, estos dos tipos de órbitas son las más utilizadas actualmente. •Órbitas MEO (Medium Earth Orbit): son las que se encuentran entre las órbitas LEO y GEO. Proporcionan una mayor cobertura que las órbitas LEO y necesitan un menor coste energético para su uso.
10 Capítulo 2. Conceptos básicos de mecánica orbital Haciendo uso de (2.13) , (2.12) y (2.9) puede relacionarse directamente la anomalía verdadera θ y el tiempo transcurrido desde el último paso por periápsis mediante el siguiente método: A partir de la anomalía verdadera de un cierto punto de la órbita, se calcula el tiempo que debe transucrrir desde el paso por periápsis: •Primero, a partir de θy con (2.12), se calcula la anomalía excéntrica E asociada. • Después, con la anomalía excéntrica correspondiente y la ecuación de Kepler (2.13) , se calcula la anomalía media M (esta ecuación debe resolverse de forma numérica). • Por último, a partir de la definición de la anomalía media (2.9) , se obtiene el tiempo transcurrido desde el último paso por periápsis. En caso de que se quisiera resolver el problema contrario, es decir, la determinación de la anomalía verdadera que ocuparía el cuerpo transcurrido un tiempo desde el paso por periápsis, se seguirá el mismo proceso pero de forma inversa. Aclaración con respecto a la nomenclatura: debido a que este documento se centra en el estudio de órbitas geocéntricas, los puntos de periápsis y apoápsis también pueden denominarse perigeo y apogeo, respectivamente. Por esta razón se utilizará en adelante esta nomenclatura en lugar de la genérica utilizada hasta ahora. Esta ley horaria se simplificará haciendo uso de un desarrollo en serie de Taylor de segundo orden con la excentricidad como variable, para poder simplificar los cálculos asociados a la determinación de la posición de los puntos de corte en una órbita geocéntrica de excentricidad pequeña. El desarrollo de la aproximación de la ley horaria puede verse en la referencia [3], y su uso en este documento en el capítulo 4. 2.3 Elementos orbitales Los elementos orbitales son parámetros que definen geométricamente por completo la órbita. Estos elementos (para órbitas geocéntricas) se definen en un sistema de referencia geocéntrico ecuatorial, cuyo plano XY contiene el Ecuador y cuyo eje X se identifica con el primer punto de Aries , el cual, aunque se llama punto, realmente define una dirección, la que se obtiene al unir la Tierra con el Sol durante el equinoccio de primavera. El ángulo que forma la proyección sobre el Ecuador de la recta que une el satélite con el origen de coordenadas y el eje X se conoce como ascensión recta , mientras que el ángulo formado por esta misma recta que une el satélite con el origen y el Ecuador se conoce como declinación. Imagen 2.5 (2018). Sistema de referencia ecuatorial inercial [Imagen]. Obtenido de https://lco.global/spacebook/equatorialcoordinate-system/. Antes de definir los elementos como tal, es necesario introducir una serie de conceptos que facilitarán esta definición: • Existen dos planos clave en la definición de una órbita: el plano de referencia y el plano orbital. El primero es, como su nombre indica, el plano a partir del cual se toman medidas y direcciones. Puede ser, por ejemplo, la eclíptica para el estudio de planetas, el plano que contiene el Ecuador, etc... Por otra parte, el plano orbital es aquel en el que está contenida la propia órbita. • La recta en la que se cortan estos dos planos determina la línea de nodos. El punto en el cual el satélite cruza el plano de referencia de Sur a Norte se conoce como nodo ascendente , mientras que el punto en el que lo hace de Norte a Sur es el nodo descendente . • El vector que une el centro de la Tierra con el nodo ascendente se denomina~n, y será útil a la hora de medir ángulos. Con esto, pueden describirse los elementos orbitales que definirán geométricamente la órbita. Para órbitas no degeneradas (las cuales se aclararán más tarde) son los siguientes:
2.3 Elementos orbitales 11 Descripción de la órbita en su propio plano: •Semieje mayor a: ya se definió cuando se explicaron los distintos elementos de una elipse. •Excentricidad e: ya definida. •Argumento de periápsis ω : es el ángulo entre el nodo ascendente y el punto de periápsis en el plano orbital, medido en la dirección del movimiento. Imagen 2.6 Descripción de la órbita en su propio plano. Descripción de la ubicación del plano orbital: •Inclinación i : es el ángulo entre los planos de referencia y orbital. Se debe medir ”mirando” el plano orbital desde el nodo ascendente, es decir, ”enfrentándose” al vector ~n . Esto hace que i∈[0,π] . Si i∈[0,π 2) , la órbita será directa, mientras que en caso contrario, esta se denominará retrógrada. Para el caso i=π 2 , las órbitas son polares. Otra forma más fácil de medir la inclinación es obteniendo el ángulo que forman el vector momento cinético de la órbita ( ~rx~v ) y el vector normal al plano del Ecuador en dirección al polo Norte. •Ascensión recta del nodo ascendente RAAN o Ω : como indica su nombre, es el ángulo que forma el eje X del sistema de referencia (que viene definido por el primer punto de Aries) y el vector~n. Imagen 2.7 Descripción del plano orbital. Para caracterizar la posición del propio cuerpo sobre su órbita, será necesario definir su anomalía verdadera θ , quedando completamente definida entonces la posición del cuerpo en cuestión mediante la siguiente terna de datos: {a,e,ω,RAAN,i,θ} . Cabe aclarar que la anomalía verdadera puede sustituirse por la excéntrica E, la media M o incluso el tiempo transcurrido desde el último paso por periápsis. Por último, se explican los casos degenerados mencionados anteriormente, en los que uno o varios elementos orbitales no están bien definidos, con lo que es necesario utilizar otros: •Órbitas elípticas ecuatoriales: en este caso ni Ω ni ω están bien definidos, al no existir línea de nodos. Se sustituye por el ángulo longitud del perigeo ω=Ω+ω , el cual indica el ángulo entre el primer punto de Aries y el perigeo de la órbita medido en sentido antihorario.
12 Capítulo 2. Conceptos básicos de mecánica orbital •Órbitas circulares: en este caso, al no existir una variación del radio, la línea de ápsides no puede definirse, y por tanto tampoco ω . Por ello, se utiliza el ángulo argumento de latitud u=θ+ω para medir la posición del satélite desde en el sentido del movimiento. Cabe aclarar que esta variable u no tiene limitado su uso a órbitas circulares, ya que aunque ω esté definida, puede calcularse el argumento de latitud para órbitas excéntricas mediante su definición. •Órbitas circulares ecuatoriales: en este caso no existe ni línea de nodos ni línea de ápsides, de forma que ni Ω , ni ω ni θ están bien definidos. Por ello se usa la longitud verdadera λT=Ω+ω+θ , que es el ángulo formado entre el satélite y el primer punto de Aries . 2.4 Definición y determinación de trazas Se define traza como la curva proyectada sobre la superficie terrestre (o del planeta central) formada por todos los puntos que sobrevuela el satélite (puntos subsatélite). Imagen 2.8 Sistema de referencia ecuatorial no inercial. Debido a que es necesario tener en cuenta la rotación terrestre, el sistema de referencia utilizado no será el mismo, sino que en este caso rotará con la Tierra. En este caso, el plano XY seguirá coincidiendo con el Ecuador, pero el eje X será en todo momento la recta que une el origen del sistema y el punto de corte del meridiano de Greenwich con el Ecuador. De esta forma, el eje X del sistema rotará a la misma velocidad que la Tierra. Las coordenadas de este sistema de referencia son la longitud λ de un punto, siendo esta el ángulo que forman el meridiano de Greenwich y el meridiano que contiene el punto en cuestión y la latitud φ , que es el ángulo que forma la recta que une el punto de estudio con el Ecuador. Dicho esto, las trazas se representarán como gráficas del tipo (φ(t),λ(t)) . Para representar la superficie terrestre sobre la que se incluirá esta curva, pueden utilizarse dos tipos distintos de proyecciones. Las dos más comunes son las siguientes: •Proyección cilíndrica equidistante: es un tipo de proyección en la que las distancias iguales abarcan ángulos iguales. Por tanto, no es una transformación conforme, pero permite representar la Tierra en una superficie finita. Debido al uso de esta escala lineal para latitud y longitud, los puntos cercanos a los polos quedarán deformados, de forma que las distancias estarán más deformadas cuanto más alejadas estén del Ecuador. (Imagen 2.9a) •Proyección esterográfica: proyección conforme sobre un plano tangente a uno de los polos. Esta proyección respeta los ángulos pero no las áreas, pero es muy útil para estudiar trazas con inclinaciones cercanas a 90º. (Imagen 2.9b) (a) Proyección cilíndrica equidistante. (b) Proyección estereográfica. Imagen 2.9 Ejemplos de proyecciones de la superficie terrestre. Para obtener ecuaciones con las que determinar las trazas, es necesario hacer uso de la trigonometría esférica, por lo que se va a hacer una pequeña introducción de la misma:
2.4 Definición y determinación de trazas 13 2.4.1 Trigonometría esférica Se define como parte de la geometría esférica que estudia los polígonos que se forman sobre la superficie de una esfera, en especial, los triángulos. Un triángulo en geometría esférica se forma mediante la intersección de 3 rectas esféricas, las cuales se definen como la intersección de un plano que contiene el origen de la esfera y la superficie de la misma. En un triángulo esférico existen seis ángulos: los formados entre las rectas en los vértices, que se denominan α , β y γ (los cuales no suman 180 º ) y tres ángulos interiores, a, b y c que se oponen a los anteriores. Si el radio de la esfera es unitario, entonces a, b y c (en radianes) corresponden a las longitudes de los arcos que forman el triángulo; por ello se denominan ”lados”, mientras que α,βyγson los ”ángulos” del triángulo. Imagen 2.10 Ejemplo de triángulo esférico. Imagen 2.11 Triángulo esférico (zoom). A la derecha puede verse el triángulo esférico por sí solo. En esta geometría y utilizando la nomenclatura de 2.11 se cumplen las siguientes relaciones: sin a sin α=sin b sin β=sin c sin γ(2.14) cos α=−cos βcos γ+sin βsin γcos a(2.15) cos β=−cos αcos γ+sin αsin γcos b(2.16) cos γ=−cos βcos α+sin βsin αcos c(2.17) cos a=cos bcos c+sin bsin ccos α(2.18) cos b=cos acos c+sin asin ccos β(2.19) cos c=cos bcos a+sin bsin acos γ(2.20) Haciendo uso de las relaciones trigonométricas (2.14) - (2.20) pueden calcularse expresiones con las que se puede representar la traza de un satélite.
14 Capítulo 2. Conceptos básicos de mecánica orbital 2.4.2 Cálculo de la traza Puede obtenerse una expresión analítica de la traza siempre que se conozca la ascensión recta del meridiano de Greenwich ( GST0 ) y los elementos orbitales del satélite en el instante de referencia t=0. En la imagen 2.12 se representa una órbita genérica en un instante cualquiera. Si al triángulo esférico que se forma mediante el ecuador, el meridiano sobre el que se encuentra el satélite y la proyección de la propia órbita sobre la superficie terrestre se le aplican las relaciones trigonométricas descritas anteriormente, se llega a: Imagen 2.12 Órbita genérica. sin(φ(t)) = sin(u(t))sin(i)(2.21) tan(GST0+ω⊕t+λ(t)−Ω) = tan(u)cos(i)(2.22) Donde ω⊕es la velocidad angular de rotación de la Tierra De esta forma, las ecuaciones que determinarán la traza serán: φ(t) = arcsin(sin(u(t))sin(i)) (2.23) λ(t) = arctan(tan(u(t))cos(i))−GST0−ω⊕t+Ω(2.24) Para representarla, sería necesario contar con una función de la forma u=u(t) , la cual puede obtenerse de forma muy sencilla para el caso circular, pero no para el caso elíptico. Al representarse normalmente estas trazas de forma numérica, se barrerán distintos valores de t, con los que se calcularán los valores correspondientes de θ mediante el proceso explicado en el apartado ”Leyes horarias”. Con estos valores de θ , puede calcularse el valor correspondiente de u conocido el argumento del perigeo, y por tanto, los valores de φyλcorrespondientes a ese instante de tiempo. Una peculiaridad de las trazas es que, a pesar de que la órbita (sin considerar perturbaciones) permanece inalterable en el tiempo, su proyección sobre la superficie terrestre no lo hace, debido a la propia rotación de la Tierra. Esta rotación hace que tras un periodo del satélite, la traza se haya visto desplazada hacia el Oeste una cierta cantidad de grados, la cual será proporcional al propio periodo del satélite y a la velocidad angular de la Tierra como ∆λ=−ω⊕TSAT y que se denomina retraso nodal . Esto puede verse en la siguiente traza de ejemplo, en la que únicamente se representa un periodo orbital:
2.4 Definición y determinación de trazas 15 Imagen 2.13 Ejemplo de representación de traza (con su retraso nodal). 2.4.3 Órbitas con traza repetida El estudio en el que se centra este documento es la determinación de los puntos de corte en órbitas cuya traza es repetida pero, ¿qué significa que la traza es repetida? Esto quiere decir que tras un cierto número de periodos, la traza se cierra sobre sí misma, generando una única curva sobre la superficie terrestre independientemente del número de periodos que se representen. Para que esto sea posible, el retraso nodal ∆λ debe ser múltiplo de 2π . Esto se traduce en que el periodo orbital del satélite debe ser múltiplo o submúltiplo del periodo de rotación de la Tierra, de forma que se cumpla: TSAT =m kT⊕(2.25) En (2.25) k y m son dos números enteros, coprimos, que representan el número de periodos orbitales que deben transcurrir para que la traza se repita y el número de días sidéreos (no son exactamente 24h sino un poco menos) necesarios que para que esto ocurra, respectivamente. Esta nomenclatura será la que se utilice en el resto del documento. Se incluye a continuación una traza de ejemplo de una órbita cuya traza se repite cada 5 días y 6 revoluciones (k=6ym=5): Imagen 2.14 Ejemplo de representación de traza repetida.
16 Capítulo 2. Conceptos básicos de mecánica orbital 2.5 Definición de λu Debido a que para los posteriores desarrollos matemáticos se hará uso de la variable λu , incluiremos aquí su definición: esta representa el ángulo que recorre la traza proyectada sobre el Ecuador, a partir del nodo ascendente, hacia el Este y sin tener en cuenta la rotación terrestre. Se puede calcular del triángulo de la imagen 2.12 de varias formas, por ejemplo tan(λu) = tan(u)cos(i) . Es importante tener en cuenta que, en función de si la órbita es directa o retrógrada, será necesario ajustar la solución para que λu se encuentre en el mismo cuadrante que u para órbitas directas, y en cuadrantes opuestos para órbitas retrógradas (el primer cuadrante es el opuesto del cuarto y el segundo lo es del tercero). Si la órbita es polar, λu solo puede tomar dos valores, λu=0cuando u∈(−π 2,π 2)yλu=180 cuando u∈(π 2,3π 2). Para ilustrar esta necesidad de ajuste, se incluyen dos esquemas para cada uno de los casos mencionados no degenerados: una órbita directa y otra retrógrada: Imagen 2.15 Relación uyλupara órbita directa. Imagen 2.16 Relación uyλupara órbita retrógrada.
3 Resumen del problema circular Para analizar el problema de los puntos de corte en órbitas con trazas repetidas, se comienza estudiando el caso más sencillo, que son aquellas órbitas cuya excentricidad es nula, es decir, las circulares. En este capítulo se obtendrá una ecuación maestra general para determinar los argumentos de latitud y el número de periodos orbitales que han de transcurrir para que se produzcan estos cruces, además de analizar cómo solucionar esta ecuación para los dos casos que pueden darse en función de los parámetros de la órbita: kym, todo ello de manera sucinta, ya que este estudio ya ha sido realizado en [1]. 3.1 Obtención de la ecuación maestra Se denomina ecuación maestra a aquella que permite la obtención de los argumentos de latitud para los que se obtienen puntos de cruce, así como el número de periodos orbitales que han de transcurrir desde la vuelta de referencia para que se produzcan estos cruces. Las ecuaciones de partida son las siguientes: sin(u)sin(i) = sin(φ)(3.1) tan(λu) = tan(u)cos(i)(3.2) Ω+λu=GST0+ω⊕t+λ(3.3) Los puntos de cruce serán aquellos cuyas latitud φ y longitud λ sean iguales sin tener el mismo argumento de latitud u, ya que el satélite no puede pasar dos veces por un mismo punto en el mismo sentido sin haberse repetido la traza, es decir, se cumplirá que (diferenciando los puntos mediante los subíndices 1y2): φ1=φ2(3.4) λ1=λ2(3.5) u16=u2(3.6) λu,16=λu,2(3.7) Como la órbita es la misma durante todo el análisis y no se están considerando en este documento las perturbaciones existentes, los parámetros Ω e i se mantienen constantes durante todo el estudio. Además, los parámetros propios del planeta alrededor del que orbitan los satélites en cuestión (en este caso la Tierra), GST0yω⊕, también se mantienen constantes. Con estas consideraciones, y con las ecuaciones (3.1) - (3.7) , puede obtenerse la ecuación maestra para la resolución del problema circular: Primero se particulariza (3.1) para los puntos 1 y 2 , para después dividir estas ecuaciones entre sí, obteniendo: sin(u1) sin(u2)=1(3.8) 17
18 Capítulo 3. Resumen del problema circular Teniendo en cuenta la expresión (3.6) , la única relación entre los argumentos de latitud que cumple la ecuación (3.8) es: u2=π−u1(3.9) Si ahora se siguen los mismos pasos con la ecuación (3.2), se llega a: tan(λu,1) tan(λu,2)=tan(u1) tan(u2)=−1(3.10) Donde se ha tenido en cuenta que tan(u2) = −tan(u1) , debido a la relación (3.9) que se cumple entre los argumentos de latitud. Dada la ecuación (3.10) , las posibles relaciones entre las variables λu,1 y λu,2 son λu,1=−λu,2 o λu,2=π−λu,1 . Para comprobar cuál de las dos es la correcta, se hace uso de la relación entre las variables u y λu , (3.2) , la cual además depende de si la inclinación es mayor o menor de π 2 . Al hacer uso de esta relación, puede verse que, ya sea la órbita directa ( i<π 2 ) o retrógrada ( i>π 2 ), la relación que se obtiene entre las variables λues: λu,2=π−λu,1(3.11) Una vez obtenidas las relaciones entre las variables u1 , u2 , λu,1 y λu,2 , se procede a obtener la ecuación maestra como tal. Para ello, se particulariza la ecuación (3.3) para ambos puntos y se restan entre sí, de forma que se obtiene otra relación entre λu,1 y λu,2 , que viene dada a través del tiempo que transcurre desde que pasa por primera vez por la latitud y longitud del punto de cruce, hasta que vuelve a pasar. λu,2−λu,1=ω⊕(t2−t1)(3.12) Para no trabajar con las unidades temporales, estas se relacionan con los argumentos de longitud a través de la ley horaria para órbitas circulares: ui=u0+nsat (ti−t0)−2πri(3.13) Donde nsat =qµ⊕ a3 es la velocidad angular del satélite, a el semieje mayor de la órbita, µ⊕≈Gm⊕ es el parámetro gravitacional de la Tierra, el subíndice i hace referencia a cualquiera de los dos puntos que se analizan, el subíndice 0 al instante de referencia y r representa el número de órbitas completas que ha dado el satélite desde el momento t0 . Para simplificar el problema, se asumirá que las variables 1 corresponden a parámetros de la primera vuelta, es decir, r1=0 , de forma que los puntos de corte que se extraerán serán los correspondientes a esta órbita ”inicial” con el resto. Después bastará aplicar las propiedades de repetición de la traza para obtener el resto de puntos de corte, o simplemente tener en cuenta que esta órbita ”inicial” puede ser cualquiera. Con esto, si se particulariza la ecuación (3.13) para ambos puntos y se restan entre sí, se llega a la siguiente relación: u2−u1=nsat (t2−t1)−2πr(3.14) Donde ahora r representa el número de vueltas completas realizadas desde la vuelta de referencia hasta que se obtiene el punto de cruce. Despejando de (3.14) el tiempo transcurrido entre ambos sobrevuelos del mismo punto sobre la Tierra y, teniendo en cuenta que, por ser una órbita de traza repetida, se cumple que ω⊕ nsat =m k y sustituyendo la relación obtenida en (3.12) , se consigue otra expresión que relaciona las variables uyλu, además de la ecuación (3.2): λu,2−λu,1=m k(u2−u1+2πr)(3.15) Aplicando las relaciones (3.9) y (3.11) para dejar la ecuación (3.15) únicamente en función de variables 2 , se obtiene la relación que han de cumplir u2yλu,2para que pertenezcan a un punto de corte de la traza: λu,2=π 2+m k(u2−π 2+πr)(3.16) Sustituyendo (3.16) en (3.2) se obtiene una ecuación en u2 que aporta los argumentos de latitud para los cuales se da un punto de corte en la órbita número r(desde la referencia): tan(π 2−m k[π 2−u2−πr]) = tan(u2)cos(i)(3.17)
3.2 Resolución de la ecuación 19 En esta ecuación, las posibles soluciones de u2 se encontrarán en el intervalo [−π,π] , mientras que r podrá encontrarse en el intervalo [0,k−1] , ya que si el satélite completa k periodos, su traza se estará repitiendo, y lógicamente sobrevolará de nuevo los mismos puntos. 3.2 Resolución de la ecuación En esta sección se explica cómo ha de modificarse la ecuación maestra previamente obtenida para su resolución analítica, aunque no es competencia de este documento incluir esta resolución analítica completa, porque ya ha sido llevada a cabo en [1]. Lo que sí se va a explicar con detalle es cómo se obtienen las soluciones de esta ecuación maestra de manera numérica, ya que estas soluciones pueden ser la base para la obtención de las correspondientes a esta misma ecuación para órbitas excéntricas, como se explicará más adelante. Para resolver la ecuación maestra, primero se reescribe haciendo uso de la propiedad de la tangente tan(α) = tan(α+pπ)con pun número entero real, de forma que queda: tan(π 2+(u2−π 2)m k−πr k) = tan(u2)cos(i)(3.18) Donde ahora r es otra variable distinta, que sigue representando el número de periodos completos que ha realizado el satélite desde el primer sobrevuelo hasta el segundo, pero que sigue un orden distinto al de la primera variable, todo ello para simplificar la ecuación. Por comodidad, y debido a que más adelante se realizarán otros cambios de variable en r , todas estas variables se llamarán igual, pero es necesario remarcar que cada vez que se modifique la ecuación maestra, la variable r será distinta, debiendo, una vez resuelta la ecuación pertinente, deshacer todos los cambios de variable para obtener el número de vueltas reales que el satélite ha completado entre ambos sobrevuelos. 3.2.1 Región de soluciones La ecuación (3.18) presenta simetrías respecto a la variable u2 si se modifica el valor de r que se analiza. Por ejemplo, se comprueba que si u2=u∗ 2 es solución de (3.18) para un r=r∗ , entonces u2=−u∗ 2 es solución para otro rdistinto: tan(π 2+(−u∗ 2−π 2)m k−πr0 k) = tan(−u∗ 2)cos(i)(3.19) Teniendo en cuenta que tan(−α) = −tan(α) , y que u∗ 2 es solución de la ecuación (3.18) , debe cumplirse que: −π 2+m ku∗ 2+m k π 2+πr0 k=π 2+(u∗ 2−π 2)m k−πr∗ k(3.20) De donde se obtiene la relación entre la r0 para la cual −u∗ 2 es solución de la ecuación (3.18) : r0=k−m−r∗ . Gracias a esta simetría, puede reducirse el intervalo de las soluciones u2de [−π,π]a[0,π]. Repitiendo este proceso con una solución del tipo u2=π−u∗ 2 , se llega a otra relación de este tipo, en la que existe otro valor de la variable r0 para el cual π−u∗ 2 es solución de (3.18) : r0=k−r∗ . Gracias a la simetría de las soluciones con respecto a π 2 , puede reducirse el intervalo de r que se analiza a la mitad, es decir, de [0,k−1]a[0,k−1 2]si kes impar o a [0,k 2]si kes par. Con estas simplificaciones, el problema a resolver será:
26 Capítulo 3. Resumen del problema circular Por ilustrar también el caso impar, se incluyen también las siguientes gráficas: (a) i=53.34º.(b) i=82º. Imagen 3.5 Trazas de una órbita circular con k=3ym=1y varias inclinaciones. De nuevo, puede verse cómo entre i=53º e i=82º existirá una inclinación crítica, ya que el número de puntos de cruce se ha visto aumentado de una órbita a otra. En [1] se realiza un estudio detallado de estas inclinaciones para las órbitas circulares, así como la influencia de estas en el número de puntos de corte. En el presente documento, se estudiará también qué influencia tiene la inclinación de la órbita en el número de puntos de cruce, y cómo la inclusión de la excentricidad afecta a la cantidad y al valor de estas inclinaciones críticas.
4 Problema excéntrico El desarrollo del problema circular servirá de base para el objetivo principal de este documento: la determinación de los puntos de corte de órbitas geocéntricas con excentricidad y traza repetida. En este capítulo se determinarán tres posibles métodos para la obtención de estos puntos de corte en los casos en los que la órbita sea excéntrica, además de realizarse varios estudios de sensibilidad, es decir, se analizará cómo la variación de algunos elementos orbitales afecta a la posición y cantidad de puntos de corte. 4.1 Obtención de la ecuación maestra En el caso excéntrico, si se quiere determinar con completa precisión la localización de los puntos de corte, será necesario hacer uso de la ley horaria para órbitas excéntricas, de forma que la ”ecuación maestra” a obtener será un sistema de 3 ecuaciones no lineales con 3 incógnitas. Debido a que la resolución de este sistema es muy costosa en recursos y tiempo, después se realizará una aproximación de esta ley horaria mediante Taylor, de forma que podrá reducirse este sistema a una única ecuación muy similar a la obtenida en el caso circular. Como ya se ha comentado, primero se va a determinar el sistema de ecuaciones que proporciona de manera exacta la posición de los puntos de corte para órbitas excéntricas de traza repetida. Para ello, se parte de las mismas ecuaciones que en el caso circular: sin(u)sin(i) = sin(φ)(4.1) tan(λu) = tan(u)cos(i)(4.2) Ω+λu=GST0+ω⊕t+λ(4.3) Las condiciones para que exista un punto de corte es que se cumplan las siguiente relaciones (de nuevo las mismas que en el caso circular): φ1=φ2(4.4) λ1=λ2(4.5) u16=u2(4.6) λu,16=λu,2(4.7) Al igual que en el capítulo dedicado a las órbitas circulares, no se están considerando las perturbaciones que afectan a la órbita, de forma que, los elementos orbitales de la misma permanecen invariables durante su estudio: Ω , i , e y ω . De forma similar al caso circular, se particularizan las ecuaciones (4.1) y (4.2) para ambos puntos y se dividen entre sí, llegando de nuevo a las mismas relaciones: sin(u1) = sin(u2)(4.8) tan(λu,1) = −tan(λu,2)(4.9) 27
28 Capítulo 4. Problema excéntrico Atendiendo a las relaciones (4.6) y (4.7) , así como por similitud con el caso circular, se llega a las siguientes relaciones entre los argumentos de latitud y entre ambas λu: u2=π−u1(4.10) λu,2=π−λu,1(4.11) Una vez obtenidas estas relaciones, se particulariza la ecuación (4.3) para ambos puntos y se restan entre sí, para obtener una relación entre el tiempo transcurrido entre ambos sobrevuelos del mismo punto y la variación del parámetro λuen ambos periodos: λu,2−λu,1=ω⊕(t2−t1)(4.12) Tras esto, en el caso circular, se aplicaba la ley horaria para obtener otra relación entre las variables u y λu , además de (4.2) . En este caso también se va a aplicar la ley horaria, pero al tener excentricidad la órbita, esta ley horaria no relaciona directamente los argumentos de latitud y los incrementos de tiempo, sino que se relacionan las anomalías medias Micon los incrementos de tiempo: Mi=M0+nsat (ti−t0)−2πri(4.13) De nuevo el subíndice i hace referencia a cualquiera de los dos puntos, el subíndice 0 se corresponde con el instante de referencia, el parámetro ri indica el número de periodos completos que el satélite ha realizado desde la órbita de referencia y la velocidad angular del satélite sigue teniendo la misma definición: nsat =qµ⊕ a3, con µ⊕el parámetro gravitacional de la Tierra y ael semieje mayor de la órbita. De nuevo, a lo largo de todo el desarrollo, se tomará como órbita de referencia la inicial, de forma que todos los puntos de corte que se obtengan del sistema de ecuaciones (y posteriormente de la ecuación maestra aproximada) serán los correspondientes a esta primera órbita. Esto implica que r1=0 y r2=r . Para obtener los puntos de corte del resto de órbitas entre sí, basta con propagar los de la primera órbita mediante el retraso nodal o tener en cuenta que esta órbita de referencia podría ser cualquiera. Particularizando (4.13) para ambos puntos y restando ambas ecuaciones se llega a una relación entre los incrementos de tiempo y de anomalía media muy similar al del caso sin excentricidad: M2−M1=nsat (t2−t1)−2πr(4.14) La anomalía media se puede relacionar directamente con la anomalía excéntrica Ei mediante la ecuación de Kepler: Mi=Ei−esin(Ei)(4.15) Particularizando (4.15) para ambos puntos y sustituyendo las relaciones en (4.14) se llega a una ecuación a partir de la cual puede obtenerse el tiempo transcurrido entre ambos sobrevuelos si se conocen las anomalías excéntricas del punto de corte para ambos periodos: E2−esin(E2)−E1+esin(E1)+ 2πr=nsat (t2−t1)(4.16) Si se despeja el incremento de tiempo t2−t1 de (4.16) y se introduce en (4.12) teniendo en cuenta que por tratarse de un satélite con traza repetida la relación que existe entre las velocidades angulares del satélite y de la Tierra ω⊕ y nsat es ω⊕ nsat =m k , se llega a la siguiente ecuación que relaciona las anomalías excéntricas E con las variables λu: 2λu,2−π=m k(E2−esin(E2)−E1+esin(E1)+2πr)(4.17) Es necesario tener en cuenta que para obtener (4.17) se ha utilizado la relación entre ambas λu que proporciona la ecuación (4.11) . Si se despeja λu,2 de (4.17) y se sustituye en (4.2) particularizada para los puntos 2se puede obtener la primera ecuación del sistema: tan(π 2+m 2k(E2−esin(E2)−E1+esin(E1)+2πr)) = tan(u2)cos(i)(4.18)
4.1 Obtención de la ecuación maestra 29 Para completar el sistema de ecuaciones es necesario establecer las relaciones entre las anomalías excéntricas E1 y E2 con el argumento de latitud u2 . Para ello, se hace uso de la definición de la anomalía excéntrica (4.19) , que la relaciona directamente con la anomalía verdadera y la definición del argumento de latitud (4.20), que lo relaciona con la anomalía verdadera: tan(Ei 2) = r1−e 1+etan(θi 2)(4.19) ui=θi+ω(4.20) Particularizando estas ecuaciones para ambos puntos y teniendo en cuenta la relación obtenida entre u2 y u1, pueden obtenerse las dos ecuaciones restantes del sistema, de forma que este queda como: tan(π 2+m 2k(E2−esin(E2)−E1+esin(E1)+2πr)) = tan(u2)cos(i) tan(E1 2) = q1−e 1+etan(π−u2−ω 2) tan(E2 2) = q1−e 1+etan(u2−ω 2) (4.21) Como ya se ha comentado, este sistema de ecuaciones no lineales es complejo de resolver, pero en caso de hacerlo, aporta con exactitud el argumento de latitud u2 , las anomalías excéntricas E1 y E2 , así como el número de periodos completos que ha realizado el satélite rpara los cuales se da un punto de corte. Para poder trabajar de forma similar al caso circular, se va a aproximar la ley horaria de la órbita mediante varios desarrollos de Taylor de orden 2, tal y como viene explicado en [3]. Los pasos para obtener la ecuación maestra aproximada son los mismos que para obtener el sistema completo hasta (4.14) . A partir de ahí, se hará uso de la relación aproximada (4.22) entre la anomalía media y la anomalía verdadera para llegar a una relación entre λuyu. Mi≈θi−2esin(θi)+ 3e2 4sin(2θi)(4.22) Si se particulariza (4.22) para ambos puntos y se sustituyen ambas relaciones en (4.14) se llega a: nsat (t2−t1) = 2πr+θ2−θ1−2esin(θ2)+2esin(θ1) + 3e2 4sin(2θ2)−3e2 4sin(2θ1)(4.23) Haciendo uso de sin(α)−sin(β) = 2 cos(α+β 2)sin(α−β 2) y teniendo en cuenta la relación existente entre θ2 y θ1θ2+2ω=π−θ1 (obtenida de restar a ambos miembros de la relación (4.10) el argumento del perigeo ω) puede llegarse a la siguiente ecuación: nsat (t2−t1) = θ2−θ1−4ecos(π 2−ω)sin(θ2−θ1 2)+ 3e2 4cos(π−ω)sin(θ2−θ1)+ 2πr(4.24) Teniendo en cuenta que θ2−θ1=u2−u1 , puede usarse (4.24) para obtener una relación entre u2 , u1 , λu,2 yλu,1, si se despeja de (4.24) el tiempo transcurrido entre ambos sobrevuelos y se introduce en (4.12): λu,2−λu,1=ω⊕ nsat u2−u1−4ecos(π 2−ω)sin(u2−u1 2)+ 3e2 4cos(π−2ω)sin(u2−u1)+ 2πr(4.25) Haciendo uso de las relaciones previamente establecidas entre las variables 1 y 2 para los puntos de corte y teniendo en cuenta la relación existente entre las velocidades angulares de la Tierra y del satélite por tratarse de una órbita con traza repetida, puede obtenerse una relación directa entre λu,2yu2además de (4.2): λu,2=π 2+m ku2−π 2−2ecos(π 2−ω)sin(u2−π 2)+ 3e2 4cos(π−2ω)sin(2u2−π)+ πr(4.26)
30 Capítulo 4. Problema excéntrico Con esta relación puede obtenerse la ecuación maestra aproximada para la determinación de los puntos de corte en órbitas excéntricas, sustituyendo (4.26) en (4.2) una vez particularizada para los puntos 2 . La ecuación maestra aproximada es la siguiente: tan π 2−m kπ 2+2ecos(π 2−ω)sin(u2−π 2)−u2−3e2 4cos(π−2ω)sin(2u2−π)−πr= =tan(u2)cos(i)(4.27) Puede comprobarse que si se hace e=0 , la ecuación maestra aproximada (4.27) que se ha obtenido coincide con la que se obtuvo para el caso circular (3.17). Por último, aplicando relaciones trigonométricas, puede simplificarse (4.27) de la siguiente forma: tan π 2−m kπ 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2)−πr= =tan(u2)cos(i)(4.28) 4.2 Resolución de la ecuación maestra aproximada Debido a la no linealidad de la ecuación (4.28) , la forma más rápida y eficiente de resolverla es de forma numérica. Para ello, se proponen dos métodos distintos: -La utilización de la función fzero de MAT LAB , la cual calcula el valor de una variable de referencia que hace cero la función que se implemente en el comando. Para que el código funcione correctamente, la función fzero necesita un iterante inicial, es decir, un valor de la variable con la que resuelve la función lo más cercano posible a la solución real. Este valor inicial de u2 se obtiene al resolver el problema circular, razón por la cual era necesario identificar para qué valores de r cada solución u∗ 2 resuelve la ecuación (3.17) , ya que r será un parámetro que también habrá que introducir en la nueva ecuación (4.28) para resolverla. Este método obtiene de forma precisa los puntos de corte asociados a las soluciones correspondientes circulares, pero presenta el inconveniente de que, en caso de generarse nuevas soluciones no asociadas a ningún punto de corte circular, el programa no las calculará. -La aproximación de la ecuación maestra mediante polinomios, y el cálculo de todas sus soluciones para cada uno de los valores posibles de r mediante el comando roots . Esta aproximación es posible gracias al paquete de MAT LAB chebfun , el cual permite una aproximación muy precisa y rápida de funciones con mucha variabilidad. Para obtener más información o su descarga puede consultarse [4]. Debido a la robustez y sencillez de este método, será el utilizado principalmente. Dicho esto, se explicarán ambos. La ecuación (4.28) presenta discontinuidades, por lo que será necesario modificarla para poder trabajar con ella cómodamente en MAT LAB . Para eliminar la discontinuidades, se multiplican ambos miembros por los cosenos de los argumentos de ambas tangentes, con el fin de eliminar los denominadores de la ecuación: cos(u2)sin π 2−m kπ 2−2esin(ω)cos(u2)−u2 3e2 4cos(2ω)sin(2u2)−πr= =cos π 2−m kπ 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2)−πrsin(u2)cos(i)(4.29) Al tratar con esta nueva ecuación (4.29) , es necesario comprobar si se introducen nuevas ”soluciones” que antes no estaban, las cuales será necesario descartar. Sabiendo que la función tan(α) se anula para los mismos valores que lo hace sin(α) y diverge cuando cos(α) = 0 (y por tanto aquellos valores que anulen el coseno no serán soluciones de la ecuación tan(α) = 0 ), será necesario comprobar si existen valores de u2 que anulan a la vez ambos cosenos de la ecuación (4.29) y no tenerlos en cuenta en caso de que los comandos fzero o roots los tome como soluciones del problema. Para ello, y sabiendo que cos(u2) se anula para u2=(2g+1)π 2 siendo g un número entero, será necesario comprobar para qué valores de g y r se anula el coseno del miembro derecho.
4.2 Resolución de la ecuación maestra aproximada 31 Si se sustituye u2=(2g+1)π 2 en el argumento del coseno del miembro derecho, puede comprobarse que los términos multiplicados por la excentricidad se anulan, ya que en ambos casos van multiplicados por cos((2g+1)π 2) y sin(2gπ) respectivamente, siendo g un entero, de forma que ambos términos serán iguales a 0. Teniendo esto en cuenta, será necesario encontrar los valores de g y r que hacen que se cumpla la siguiente igualdad: cosπ 2−m kπ 2−(2g+1)π 2−πr=0(4.30) Los valores para los que se anula la función cos son aquellos múltiplos impares de π 2 , es decir, se deberá cumplir (4.31) siendo jcualquier número entero: π 2−m kπ 2−(2g+1)π 2−πr=(2j+1)π 2(4.31) Operando con (4.31) se obtiene la siguiente relación entre los parámetros g,ryj: m k(g+r) = j(4.32) Con todo esto, el método para discernir entre las soluciones reales y las que no corresponden a la ecuación maestra original es el siguiente: 1.- Comprobar si la solución de los comandos fzero o roots puede expresarse de la forma u∗ 2=(2g+1)π 2 , siendo gun número entero. 2.- En caso de que 1.- sea posible, comprobar si el número g+r k es un entero, siendo r el número de periodos completos para el cual u∗ 2 es solución. En caso de serlo, será necesario descartar la solución e introducir un nuevo iterante inicial o simplemente no tener en cuenta esa solución, ya que, como se verá más adelante, algunos puntos de cruce pueden desaparecer al introducir la excentricidad. Una vez calculadas todas las soluciones u∗ 2 y descartadas todas las ”soluciones falsas” tal y como se acaba de explicar, pueden utilizarse las ecuaciones (4.1) - (4.3) así como la ley horaria aproximada para calcular las latitudes y longitudes de los puntos de corte. Esta ley horaria aproximada puede encontrarse en [3] desarrollada, pero por comodidad se incluye a continuación: ∆t≈ui−ω+2πr−2esin(ui−ω+2πr)+ 3e2 4sin(2(ui−ω+2πr))−M0 nsat (4.33) Para resolver la ecuación (4.33) es necesario tener en cuenta que ∆t es el tiempo transcurrido desde el instante de referencia t0 hasta el momento en el que se tiene el argumento de latitud ui en la vuelta número r desde la referencia. En caso de que el argumento de latitud del instante de referencia t0 sea distinto de 0, será necesario calcular la anomalía media en ese instante M0 para tenerla en cuenta en la ecuación a través de las relaciones (4.15),(4.19) y (4.20).
32 Capítulo 4. Problema excéntrico Se incluyen a continuación ejemplos resueltos para distintos valores de k,m,i,ωye: Imagen 4.1 Ejemplo de órbita con traza repetida cada 3 revoluciones (k=3) y cada día (m=1), de excentricidad e=0.1, inclinación 80ºy argumento de perigeo ω=50 º.. Imagen 4.2 Ejemplo de órbita con traza repetida cada 5 revoluciones (k=5) y cada 4 días (m=4), de excentricidad e=0.05, inclinación 120ºy argumento de perigeo ω=50 º.. Ambas imágenes representan un caso resuelto del conocido como ”caso impar”, con una inclinación de órbita directa, y otro caso resuelto del tipo ”impar”, en este caso con una órbita retrógrada. Una vez obtenidas las soluciones de la ecuación aproximada, pueden usarse para calcular los iterantes iniciales necesarios para la resolución del sistema de ecuaciones completo, para así obtener la solución exacta. Para ello, se hace uso del comando fsolve de MAT LAB y de las ecuaciones (4.19) , (4.20) y (4.10) para calcular estos iterantes iniciales con los que empezará a resolver el sistema para las variables E1 , E2 y u2 . Se incluye a continuación una tabla a modo de comparación con las soluciones que se obtienen resolviendo la ecuación aproximada y el sistema de ecuaciones completo:
4.2 Resolución de la ecuación maestra aproximada 33 Tabla 4.1 Comparativa entre las soluciones obtenidas por la ecuación aproximada y el sistema de ecuaciones no lineales completo para una órbita cuyas características son k=3, m=2, ω=25 º , i=85 º y e=0.15. k=3, m=2, e=0.15, ω=25 º, i=85º,Ω=0u∗ 2[rad] r[revs] φ[º]λ[º] Ecuación maestra aproximada 1.37370 1 77.6628 116.3952 0.8325 1 47.4597 114.3015 1.7679 2 77.6628 -3.6048 -1.33783 0 -75.7624 36.1206 -1.80376 2 -75.7624 156.1206 2.30909 2 47.4597 -5.6985 -1.07742 0 -61.3287 37.2475 -2.06417 2 -61.3287 157.2475 Sistema de ecuaciones completo 1.37340 1 77.6476 116.3857 0.83385 1 47.5367 114.3133 1.76819 2 77.6476 -3.6143 4.94467 2 -75.7991 36.1274 -1.80307 2 -75.7991 156.1274 2.30774 2 47.5367 -5.6867 5.20789 2 -61.2088 37.2583 -2.06630 2 -61.2088 157.2583 En la tabla puede comprobarse como ambas soluciones son muy similares. En concreto, para el argumento de latitud u2 , la mayor diferencia se da en la milésima de radián, para la latitud φ , el mayor error es una diferencia de décimas de grado y para la longitud λ un error de centésimas de grado. Un comentario que ha de hacerse es que, por simplicidad a la hora de resolver el sistema de ecuaciones y calcular el resto de parámetros, las soluciones contenidas en el intervalo [−π,0) que resolvían la ecuación maestra aproximada para r=0 se han pasado al intervalo [π,2π) , de forma que ahora serán soluciones para r=k−1=2 (esto puede hacerse porque un argumento de latitud negativo implica que el satélite se encuentra aún en el periodo anterior). Puede comprobarse en las columnas correspondientes a la latitud y la longitud que este cambio no les afecta, y se siguen obteniendo resultados muy similares. La conclusión que puede obtenerse de esta tabla es que la aproximación hecha de la ley horaria no es burda, y funciona con precisión de décimas de grado para excentricidades del orden de 0.15. Se incluye a continuación un gráfico con la traza de la órbita resuelta en la tabla, para completar así el ejemplo comparativo: Imagen 4.3 Traza de la órbita resuelta en la tabla 4.1.
34 Capítulo 4. Problema excéntrico 4.3 Obtención de función de simetría de las soluciones u2 En el caso circular, siempre que u∗ 2 fuese solución de la ecuación maestra, −u∗ 2 lo era para otro valor de r . Por ello, cabe preguntarse si existe alguna función del tipo u∗ 2,2=−u∗ 2,1+f(e,ω,i...) que permita una relación similar entre las soluciones del caso excéntrico. Primero se comprueba que no se mantiene la simetría de las soluciones con respecto a u∗ 2=0 como ocurría en el caso circular. Para ello, suponiendo que para unos ciertos elementos orbitales e, ω , i, k y m, u2=u∗ 2 es solución de la ecuación maestra aproximada para r=r∗ , se comprueba la existencia de un valor entero, positivo y menor que k−1de rque hace que u2=−u∗ 2también sea solución: tan π 2−m kπ 2+2ecos(π 2−ω)sin(u∗ 2−π 2)−u∗ 2−3e2 2cos(π−2ω)sin(2u∗ 2−π)−πr∗= =tan(u∗ 2)cos(i)(4.34) Para que sea más fácil comprobar el signo de los distintos elementos de los argumentos, se reescribe la ecuación (4.34): tan π 2−m kπ 2−2ecos(π 2−ω)cos(u∗ 2)−u∗ 2+3e2 2cos(π−2ω)sin(2u∗ 2)−πr∗= =tan(u∗ 2)cos(i)(4.35) Al sustituir u∗ 2 por −u∗ 2 , el segundo miembro de la ecuación (4.35) cambiará de signo, con lo que, llamando f(u2) al argumento de la primera tangente, deberá comprobarse que tan(f(u∗ 2)) = −tan(f(−u∗ 2)) . Para que esta igualdad sea cierta, debe cumplirse que f(u∗ 2) = pπ−f(−u∗ 2) (con p un número entero y positivo), es decir, que debe comprobarse que a partir de la igualdad (4.36) puede obtenerse un valor de r=r0 que sea entero, positivo y menor que k−1: π 2−m kπ 2−2ecos(π 2−ω)cos(u∗ 2)−u∗ 2+3e2 2cos(π−2ω)sin(2u∗ 2)−πr∗= =pπ−π 2+m 2π 2−2ecos(π 2−ω)cos(u∗ 2)+ u∗ 2−3e2 2cos(π−2ω)sin(2u∗ 2)−πr0(4.36) Puede comprobarse que en (4.36) los términos asociados a e2 y los propios términos u∗ 2 se cancelan entre sí, pero no ocurre lo mismo con los asociados a e . Esto implica que, a la hora de despejar r0 , estos términos influirán, siendo normalmente números no enteros debido a los valores de excentricidad con los que se está trabajando y a la propia naturaleza de las funciones trigonométricas. Para poder eliminar este término asociado a e , se buscan valores posibles de r’que vayan asociados a una solución de la forma −u∗ 2+4ecos(u∗ 2)cos(π 2−ω), para lo cual se hacen las siguientes aproximaciones: cos(−u∗ 2+4ecos(u∗ 2)cos(π 2−ω)) ≈cos(−u∗ 2)(4.37) sin(−2u∗ 2+4ecos(u∗ 2)cos(π 2−ω)) ≈sin(−2u∗ 2)(4.38) Estas aproximaciones están justificadas, ya que se está trabajando con excentricidades pequeñas, y los términos trigonométricos como mucho valdrán 1. Sustituyendo ahora −u∗ 2+4ecos(u∗ 2)cos(π 2−ω) en (4.35) , en lugar de −u∗ 2 y realizando las aproximaciones descritas anteriormente, se llega a la siguiente relación: π 2−m k(π 2−πr∗) = pπ−π 2+m k(π 2−πr0)(4.39) Despejando de esta igualdad, se llega a la siguiente relación entre r0yr∗: r0=1−r∗+k m(p−1)(4.40)
4.3 Obtención de función de simetría de las soluciones u235 En (4.40) , tomando correctamente un valor de p en función de los valores de k y m , puede obtenerse un valor de r0que cumpla con todas las condiciones de ser entero, positivo y menor que k−1. Con esto, puede concluirse que una primera aproximación de la fórmula que relaciona una de las soluciones de la ecuación maestra aproximada con la simétrica respecto a u2=0de otra de ellas es la siguiente: u∗ 2,2≈−u∗ 2,1+4ecos(u∗ 2,1)cos(π 2−ω)(4.41) Para comprobar la exactitud de esta aproximación, se comparan distintos casos de órbitas con distintos elementos orbitales: Tabla 4.2 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=3, m=1, ω=20 º, i=80ºy e=0.05. k=3, m=1, ω=20 º, i=80ºy e=0.05 u∗ 2,1[rad] u∗ 2,2[rad] r∗[revs] r0[revs] (sol. numérica) u∗ 2,2aprox. [rad] pr0(por (4.40)) -1.9765 1.9629 0 1 1.9495 1 1 -1.1651 1.1787 2 2 1.1921 2 2 -3.0707 -0.0709 0 2 3.0024 x x En este primer ejemplo puede comprobarse cómo la aproximación es exacta hasta la centésima de radián para las dos primeras parejas de soluciones. Sin embargo, existe una tercera pareja de soluciones que se diferencian aproximadamente en π, para la cual la relación (4.41) no se cumple. Tabla 4.3 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=5, m=1, ω=20 º, i=85ºy e=0.05. k=5, m=1, ω=20 º, i=85ºy e=0.05 u∗ 2,1[rad] u∗ 2,2[rad] r∗[revs] r0[revs] (sol. numérica) u∗ 2,2aprox. [rad] pr0(por (4.40)) -1.6358 1.6357 0 1 1.6314 1 1 -1.5058 1.5059 2 4 1.5102 2 4 1.8984 -1.9062 2 4 -1.9205 2 4 -1.2354 1.2431 3 3 1.2579 2 3 -0.0603 -3.0813 3 4 0.12855 x x De nuevo se repite el patrón del ejemplo anterior: la relación (4.41) aporta valores exactos hasta la centésima de radián en todas las parejas de soluciones excepto en una, en la que, de nuevo, la diferencia es un factor de π. Tabla 4.4 Comparativa entre las soluciones obtenidas por la ecuación (4.41) y la ecuación maestra aproximada para una órbita cuyas características son k=7, m=6, ω=20 º, i=85ºy e=0.05. k=7, m=6, ω=20 º, i=85ºy e=0.05 u∗ 2,1[rad] u∗ 2,2[rad] r∗[revs] r0[revs] (sol. numérica) u∗ 2,2aprox. [rad] pr0(por (4.40)) -1.5244 1.5247 0 1 1.5275 1 1 -0.3279 0.2544 0 1 0.3926 1 1 -2.8137 2.8872 2 6 2.7490 7 6 -1.6172 1.6169 2 6 1.6141 7 6 0.8866 -0.9438 2 6 -0.8433 7 6 1.4323 -1.4289 2 6 -1.4229 7 6 -2.1978 2.2550 3 5 2.1577 7 5 -1.7127 1.7093 3 5 1.7030 7 5
42 Capítulo 4. Problema excéntrico Imagen 4.12 Zoom 1 de la gráfica 4.11. Imagen 4.13 Zoom 2 de la gráfica 4.11.
4.4 Análisis de sensibilidad 43 Imagen 4.14 Zoom 3 de la gráfica 4.11. Imagen 4.15 Zoom 4 de la gráfica 4.11.
44 Capítulo 4. Problema excéntrico Imagen 4.16 Zoom 5 de la gráfica 4.11. Imagen 4.17 Evolución de la posición de los puntos de cruce para órbitas con k=7, m=3, i=85 º y ω=−15 º .
4.4 Análisis de sensibilidad 45 Imagen 4.18 Zoom 1 de la gráfica 4.17. Imagen 4.19 Zoom 2 de la gráfica 4.17.
46 Capítulo 4. Problema excéntrico Imagen 4.20 Zoom 3 de la gráfica 4.17. Observando todas las gráficas, pueden sacarse dos conclusiones respecto a la variación de la latitud y la longitud de los puntos de cruce al aumentar la excentricidad: • Cuanto mayor es la latitud inicial del punto de cruce, menor será su variación de latitud con respecto a la excentricidad. • En función del signo del argumento del perigeo, los puntos de un hemisferio se acercarán, mientras que los del otro hemisferio se alejarán. Esto puede comprobarse si se comparan las gráficas 4.4 y 4.8.
4.4 Análisis de sensibilidad 47 4.4.2 Variación con el argumento del perigeo ω: Para ver cómo afecta la posición del perigeo dentro de la propia órbita se realiza un proceso completamente análogo al realizado con la excentricidad, en el que se utiliza la derivación implícita para calcular las derivadas ∂u2 ∂ω ,∂φ ∂ω y∂λ ∂ω . Debido a que el proceso es muy similar al previamente desarrollado, en este caso los desarrollos serán menos concretos. Para obtener ∂u2 ∂ω , se aplica la derivación implícita a la ecuación maestra aproximada (4.42) y se despeja de la misma, obteniéndose: C=∂u2 ∂ω D,Donde: C=−m k−2ecos(ω)cos(u2)+ 3e2 2sin(2ω)sin(2u2)1+tan2(f(ω)) D=cos(i)1+tan2(f(ω))−m k1+3e2 2cos(2ω)cos(2u2)−2esin(ω)sin(u2)1+tan2(f(ω)) (4.49) En estas expresiones, f(ω) hace referencia a la misma función que en el apartado anterior, es decir, a (4.43) , tomando en este caso como variable el argumento del perigeo ω en lugar de la excentricidad e . Una vez obtenida ∂u2 ∂ω , pueden calcularse las otras dos derivadas que se necesitan, haciendo uso de las mismas relaciones que en el apartado anterior: (4.1)-(4.3) y (4.33). Volviendo a realizar derivaciones implícitas y despejando las derivadas en las que estamos interesados, se llega a: ∂φ ∂ω =∂φ ∂u2 ∂u2 ∂ω =cos(u2)sin(i) cos(φ) ∂u2 ∂ω =cos(u2)sin(i) p1−(sin(u2)sin(i))2 ∂u2 ∂ω (4.50) ∂λ ∂ω =∂λu ∂u2 ∂u2 ∂ω −ω⊕ ∂∆t ∂ω =1+tan2(u2) 1+tan2(λu) ∂u2 ∂ω −ω⊕ ∂∆t ∂ω = =1+tan2(u2) 1+(tan(u2)cos(i))2 ∂u2 ∂ω −ω⊕ nsat −1+2ecos(u2−ω)(1−∂u2 ∂ω ) +3e2 2cos(2(u2−ω))(∂u2 ∂ω −1)(4.51) Al igual que en el caso de la excentricidad, estas expresiones permiten predecir el comportamiento de los puntos de corte, en cuanto a su localización se refiere, conforme se varía el argumento del perigeo ω , siempre que el argumento de latitud u2de ese punto sea conocido de antemano. Se incluyen en la siguiente página varias gráficas en las que puede verse la evolución de la posición de los puntos de corte al variar el argumento de perigeo de -180ºa 180º.
48 Capítulo 4. Problema excéntrico Imagen 4.21 Evolución de la posición de los puntos de cruce para órbitas con k=3, m=1, i=80ºy e=0.05. Imagen 4.22 Zoom 1 de la gráfica 4.21.
4.4 Análisis de sensibilidad 49 Imagen 4.23 Zoom 2 de la gráfica 4.21. Imagen 4.24 Zoom 3 de la gráfica 4.21.
50 Capítulo 4. Problema excéntrico En estas gráficas se incluyen tres trazas de órbitas con los mismos elementos orbitales (los descritos al pie de 4.21) a excepción del argumento del perigeo, que es el parámetro que se varía para obtener las distintas posiciones de los puntos de corte. En las gráficas, las curvas cián corresponden a un argumento del perigeo ω=−π 2 , las curvas verdes a ω=π 2 y las curvas rojas a ω=0 . Por último, las curvas negras representan las distintas posiciones que los diferentes puntos de corte toman para todos los argumentos del perigeo analizados, mientras que los números son de nuevo para identificar los distintos zooms realizados. Para que el zoom 4.23 incluyese toda la variación de posición del punto, no puede aumentarse mucho la imagen, de forma que esta no se aprecia bien. Por ello, se incluye a continuación esa misma variación de posición en latitud y longitud no representada sobre un mapa de la Tierra: Imagen 4.25 Zoom 2 de la gráfica 4.21. Al igual que para el caso de variación con la excentricidad, se incluyen a continuación tres ejemplos más, para poder apreciar distintas variaciones con los elementos orbitales: Imagen 4.26 Evolución de la posición de los puntos de cruce para órbitas con k=3, m=2, i=85ºy e=0.05.
4.4 Análisis de sensibilidad 51 Imagen 4.27 Zoom 1 de la gráfica 4.26. Imagen 4.28 Zoom 2 de la gráfica 4.26.
58 Capítulo 4. Problema excéntrico Imagen 4.41 Zoom 5 de la gráfica 4.36. Por último, se incluyen a continuación una serie de conclusiones que pueden extraerse de las gráficas presentadas anteriormente: • Debido a que los extremos del intervalo ω=−π y ω=π son en realidad el mismo punto, la variación del argumento del perigeo entre estos dos valores da lugar a curvas de variación de posición cerradas, como puede observarse en todas las gráficas. • Al igual que ocurría con la variación respecto a la excentricidad, cuanto mayor es la latitud inicial del punto de corte a estudiar, menor será la variación de la misma. Esto, llevado al extremo, puede a dar lugar a curvas que parezcan completamente planas, en lugar de las curvas cerradas que se ven en puntos más cercanos al Ecuador. Un ejemplo de este supuesto puede apreciarse en la imagen 4.41.
4.4 Análisis de sensibilidad 59 4.4.3 Variación con la excentricidad y el argumento del perigeo conjuntamente Para terminar de ilustrar la variación de la posición de los puntos de corte con la excentricidad y el argumento del perigeo de la órbita, se incluyen a continuación varios casos en los que se analiza la variación de uno de los puntos de corte al modificar conjuntamente ambos elementos orbitales: Imagen 4.42 Variación de la posición del punto de corte situado en el hemisferio Norte con la excentricidad y el argumento del perigeo para órbitas con i=83ºk=3 y m=1. Imagen 4.43 Zoom de la gráfica 4.42.
60 Capítulo 4. Problema excéntrico Imagen 4.44 Variación de la posición de uno de los puntos de corte situados en el hemisferio Sur con la excentricidad y el argumento del perigeo para órbitas con i=87ºk=3 y m=2. Imagen 4.45 Variación de la posición de uno de los puntos de corte situados en el hemisferio Norte con la excentricidad y el argumento del perigeo para órbitas con i=85ºk=3 y m=1.
4.4 Análisis de sensibilidad 61 Imagen 4.46 Variación de la posición del punto de corte situado sobre el Ecuador con la excentricidad y el argumento del perigeo para órbitas con i=85ºk=7 y m=3. Imagen 4.47 Zoom de la gráfica 4.46.
62 Capítulo 4. Problema excéntrico Observando las gráficas obtenidas mediante la variación conjunta de ambos elementos orbitales se llega a las siguientes conclusiones: • Las curvas de variación con el argumento del perigeo ω son todas cerradas y no se cortan entre sí al variar la excentricidad, excepto en el caso de que la latitud del punto llegue a ser 0. En este caso, la curva de variación con el argumento del perigeo se corta consigo misma y con el resto de curvas en φ=0º. • La variación de posición de los puntos con el argumento del perigeo ω aumenta al aumentar la excentricidad de la órbita. • La propiedad observada para la variación por separado en la que, cuanto mayor fuese la latitud inicial del punto de corte, menor sería la variación de esta se conserva al variar ambos elementos orbitales a la vez. 4.4.4 Inclinaciones críticas Para continuar con el análisis de sensibilidad de los puntos de corte con respecto a las variaciones de los elementos orbitales, es necesario estudiar qué influencia tiene la inclinación en los mismos. Para ello se va a recordar la definición realizada en el capítulo 3 de un concepto clave para el estudio: las inclinaciones críticas . Se denominan críticas a aquellas inclinaciones que generan trazas que son tangentes en ciertos puntos. Esto implica que para inclinaciones ligeramente mayores o menores, se crearán o desaparecerán ciertos puntos de corte, y por tanto, son esenciales para la determinación del número de puntos de corte. Se incluyen a continuación un par de trazas ejemplo en las que se produce la tangencia, de forma que las inclinaciones indicadas serán críticas para esa combinación de elementos orbitales: Imagen 4.48 Ejemplo de órbita con traza repetida cada 5 revoluciones (k=5) y cada tres días (m=3), de excentricidad e=0.05, argumento del perigeo ω=30 ºe inclinación crítica..
4.4 Análisis de sensibilidad 63 Imagen 4.49 Ejemplo de órbita con traza repetida cada 9 revoluciones (k=9) y cada 7 días (m=7), de excentricidad e=0.025, argumento del perigeo ω=70 ºe inclinación crítica.. Los puntos de tangencia vienen marcados mediante un recuadro verde en cada una de las trazas. Cabe destacar que el número de inclinaciones críticas para unos elementos orbitales dados no siempre es uno, es decir, que puede existir más de una. En la primera parte de este capítulo se van a obtener las ecuaciones que permiten el cálculo de estas inclinaciones mediante la ecuación maestra aproximada. La condición de tangencia entre las trazas puede traducirse en una condición de tangencia de la siguiente función representada frente a la variable u2con el eje horizontal: F=tan(π 2−m kπ 2+2ecos(π 2−ω)sin(u2−π 2)−u2−3e2 4cos(π−2ω)sin(2u2−π))−tan(u2)cos(i) (4.52) Llamando de nuevo f al argumento de la primera tangente, la condición de tangencia entre las trazas se puede expresar como sigue: (tan(f)−tan(u2)cos(i) = 0 ∂f ∂u2(1+tan2(f))−cos(i)(1+tan2(u2)) = 0(4.53) El sistema de ecuaciones (4.53) debe resolverse en las variables u2 e i , de forma que sus soluciones serán los argumentos de latitud para los cuales se dan la tangencia en la vuelta r y la inclinación crítica i∗ para la cual se da la tangencia. Para dejar completamente determinado el sistema de ecuaciones (4.53) es necesario calcular ∂f ∂u2, lo cual ya se ha hecho anteriormente, y cuyo resultado es: ∂f ∂u2 =−m k2ecos(π 2−ω)cos(u2−π 2)−1−3e2 2cos(π−2ω)cos(2u2−π)(4.54) Para poder trabajar de una forma más sencilla, se despeja cos(i) de la segunda ecuación del sistema, para sustituirla en la primera, de forma que el sistema quedaría:
64 Capítulo 4. Problema excéntrico tan(f)−sin(u2)cos(u2)∂f ∂u2 1 cos2(f)=0 cos(i) = ∂f ∂u2 cos2(u2) cos2(f) (4.55) A la hora de despejar se ha tenido en cuenta que 1 cos2(α)=1+tan2(α) . Para resolver el sistema, bastaría con obtener las soluciones u2 de la primera ecuación, y sustituirlas en la segunda para obtener las inclinaciones críticas. Es necesario realizar una serie de aclaraciones: -Debido a las posibles discontinuidades que podrían aparecer dentro del dominio para la ecuación en u2 , es necesario operar con la misma para eliminar el denominador. Al hacer esto, es posible que se estén introduciendo nuevas soluciones, con lo que cuando se obtengan las soluciones u2 para la ecuación (4.56) , será necesario comprobar que no coinciden con ninguna de las soluciones de cos(f) = 0 . Eliminado el denominador, la ecuación a partir de la cual se obtienen los argumentos de latitud para los cuales se producen las tangencias es la siguiente: cos(f)sin(f)−sin(u2)cos(u2)∂f ∂u2 (4.56) -Cuando se obtengan las soluciones de la ecuación (4.56) habrá que calcular el coseno de la inclinación crítica asociada, y tras esto la propia inclinación. La cuestión a aclarar es que, a través de la segunda ecuación del sistema (4.55) , es posible obtener valores mayores que 1, de forma que será necesario descartarlos, al no tener asociada una inclinación real que pueda tener como coseno ese valor. Para resolver la ecuación (4.56) se vuelve a hacer uso de la aproximación polinómica que ofrece el paquete chebfun y del comando roots para obtener todas las soluciones para cada una de las r . Tras esto se calculan todos los valores de cos(i) y se descartan aquellos cuyo valor absoluto sea mayor que 1, para finalmente obtener todas las inclinaciones críticas posibles. Debido a que para una misma inclinación crítica se producen varias tangencias (lo cual puede verse en ambos ejemplos mostrados anteriormente), varias soluciones de (4.56) llevarán asociada la misma inclinación crítica, con lo que, para que el programa sea eficiente, será necesario evitar que se repitan. Para ilustrar el método se incluyen a continuación dos ejemplos, uno de cada tipo de caso y sus correspondientes trazas: Caso impar: k=7, m=5, ω=70 º, e=0.03 Tabla 4.6 Inclinaciones críticas para órbita con para k=7, m=5, e=0.03, ω=70 º. i∗ 1=52.3688 ºi∗ 2=79.0059 ºi∗ 3=80.3424 ºi∗ 4=88.9167 ºi∗ 5=89.0336 º (a) Primera inclinación crítica. (b) Segunda inclinación crítica. Imagen 4.50 Primer par de inclinaciones críticas.
4.4 Análisis de sensibilidad 65 (a) Tercera inclinación crítica. (b) Cuarta inclinación crítica. Imagen 4.51 Segundo par de inclinaciones críticas. Imagen 4.52 Quinta inclinación crítica. Si se observan todas las trazas, puede verse cómo van creándose puntos de corte al pasar de una inclinación crítica a la siguiente, además de las propias tangencias asociadas a la inclinación crítica correspondiente. Caso par: k=7, m=4, ω=70 º, e=0.03 Tabla 4.7 Inclinaciones críticas para órbita con para k=7, m=4, e=0.03, ω=70 º. i∗ 1=75.2376 ºi∗ 2=77.2457 ºi∗ 3=88.6388 ºi∗ 4=88.7867 º (a) Primera inclinación crítica. (b) Segunda inclinación crítica. Imagen 4.53 Primer par de inclinaciones críticas.
66 Capítulo 4. Problema excéntrico (a) Tercera inclinación crítica. (b) Cuarta inclinación crítica. Imagen 4.54 Segundo par de inclinaciones críticas. Caso par: k=8, m=3, ω=−20 º, e=0.05 Tabla 4.8 Inclinaciones críticas para órbita con para k=8, m=3, e=0.05, ω=−20 º. i∗ 1=83.0339 ºi∗ 2=83.5845 ºi∗ 3=90 º Imagen 4.55 Primera inclinación crítica. (a) Segunda inclinación crítica. (b) Tercera inclinación crítica. Imagen 4.56 Último par de inclinaciones críticas.
4.4 Análisis de sensibilidad 67 Para comparar estos resultados con los que se obtienen en el caso circular, se incluyen a continuación las inclinaciones críticas para los mismos valores de kym, pero sin excentricidad: Caso impar: k=7, m=5 Tabla 4.9 Inclinaciones críticas para órbita con para k=7, m=5 y e=0. i∗ 1=44.4153 ºi∗ 2=79.7077 ºi∗ 3=88.9774 º Caso par: k=7, m=4 Tabla 4.10 Inclinaciones críticas para órbita con para k=7, m=4 y e=0. i∗ 1=76.3061 ºi∗ 2=88.7157 º Caso par: k=8, m=3 Tabla 4.11 Inclinaciones críticas para órbita con para k=8, m=3 y e=0. i∗ 1=83.3402 ºi∗ 2=90 º La primera característica que salta a la vista es que el número de inclinaciones críticas es menor. De hecho, se cumple para órbitas del caso impar ( m−k=2l ) o para órbitas del caso par con m impar que, si j es el número de inclinaciones críticas de la órbita circular, para valores de excentricidad pequeños se obtendrán 2j−1 inclinaciones críticas. Por otra parte, si el problema es del tipo par ( k o m par) y j el número de inclinaciones críticas para la órbita circular, al incluir la excentricidad se obtendrán 2j inclinaciones críticas para la órbita excéntrica. Tanto el aumento como la discrepancia entre las inclinaciones que se duplican en ambos casos se debe a la ruptura de la simetría N-S de la traza al incluir la excentricidad y un argumento de perigeo distinto de 0 º y 180 º (es decir, que para aquellas órbitas con excentricidad pero con ω=0º o ω=180 º , el número de inclinaciones críticas seguirá siendo j ). Pero, ¿por qué el número de inclinaciones que se duplican es distinto? La respuesta está en la disposición de los puntos de corte: para las órbitas de los casos par, cada punto de corte del hemisferio norte lleva asociado un punto de corte en el hemisferio sur, de forma que si es posible conseguir una tangencia en el hemisferio norte, también lo será en el sur. En cambio, para las órbitas del caso impar, existen puntos cerca del ecuador, los cuales no tienen ningún punto asociado en el otro hemisferio. Debido a esto, aunque se rompa la simetría, ese punto tiene asociada una única inclinación crítica, la cual no se duplica.
74 Capítulo 4. Problema excéntrico 4.5 Determinación del número de inclinaciones críticas de forma analítica Para poder determinar analíticamente el número de inclinaciones críticas, primero se determinan el número de puntos de cruce que se obtienen para i=0º , después se realizará la misma operación para i=180 º . Tras esto, se determina cuántos puntos de cruce se introducen por cada inclinación crítica, y en combinación con lo descrito anteriormente, puede obtenerse de manera analítica el número de inclinaciones críticas que tiene la órbita: 4.5.1 Número de soluciones para i=0º La ecuación maestra aproximada particularizada para i=0ºes la siguiente: F(u2) = tan(π 2−m kπ 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2)−πr)−tan(u2)(4.76) En una ecuación del tipo tan(x)−tan(y) = 0 , las soluciones son x=y+pπ siendo p un número entero. Teniendo esto en cuenta, para encontrar las soluciones de (4.76) se buscan las soluciones de la siguiente ecuación: u2+2p−1 2π=m k−π 2+2esin(ω)cos(u2)+ u2+3e2 4cos(2ω)sin(2u2)+ πr(4.77) O lo que es lo mismo: u2 k−m k=(1−2p)π 2−m kπ 2−2esin(ω)cos(u2)−3e2 4cos(2ω)sin(2u2)−πr(4.78) Primero se comprueba el número de soluciones para e=0. Para ello, se deben encontrar las soluciones de la ecuación: u2=π 2(k−m)m(2r−1)−k(2p−1)(4.79) El número de soluciones de la ecuación (4.79) será igual a la cantidad de números enteros que pueda tomar m(2r−1)−k(2p−1) k−m en el intervalo (−2,2) , o lo que es lo mismo, la cantidad de números enteros que puede tomar m(2r−1)−k(2p−1) en el intervalo (−2(k−m),2(k−m)) de forma que m(2r−1)−k(2p−1) 2(k−m)6=±π 2 . Este número de soluciones no se corresponde estrictamente con el número de puntos de cruce por órbita, ya que la primera órbita corta a varias, dando lugar a puntos de cruce que pueden obtenerse mediante un retraso nodal de otro que también se ha obtenido como solución. Esto podía no ocurrir en el problema circular debido a la creación de la variable z , la cual evitaba las repeticiones. Dicho esto, del caso circular sabemos que (4.79) lleva asociadas m−k−1 soluciones con distinta latitud (para que no sean reproducibles mediante un retraso nodal) sean cuales sean los valores de k y m , gracias al estudio realizado en [1]. Teniendo esto en cuenta, lo que se va a tratar en este documento es de demostrar que la introducción de la excentricidad no genera nuevas soluciones, sino que solo desplaza las existentes (al menos para valores de excentricidad pequeños). Para ello, se debe comprobar que la ecuación (4.77) sigue siendo continua y derivable al introducir los términos asociados a la excentricidad y que no se introducen nuevas soluciones al incluirse los mismos. Existen tres formas de que se introduzcan nuevas soluciones: •A partir de un punto crítico , es decir, que se de una combinación de valores de e y u2 tal que se anule (4.78) y su derivada con respecto a u2. • Que se genere una nueva solución desde las asíntotas de las funciones tangente que componen la ecuación (4.76). •Que se genere una nueva solución a partir de los extremos del intervalo de análisis [−π,π]. Existencia de puntos críticos Para comprobar la existencia de este tipo de puntos, se calcula la derivada de la función (4.78) respecto a u2 , y se comprueba si existe alguna combinación de elementos orbitales para la cual se anula. Esta derivada queda de la forma:
4.5 Determinación del número de inclinaciones críticas de forma analítica 75 k−m=−m2esin(ω)sin(u2)−3e2 2cos(2ω)cos(2u2)(4.80) k−m=−m(−2e−3e2 2)(4.81) Observando la ecuación puede verse que, para valores de k y m pequeños, no existirá una excentricidad tal que permita que se cumpla la igualdad, pero tomando las funciones trigonométricas de forma que se llegue a la situación más desfavorable ( (4.81) con ω=−90 º y u2=90 º ) y analizando la expresión de (4.80) , puede verse que cuanto mayores sean k y m , y más cercanos estén, podrían dar lugar a que existiese una excentricidad pequeña que resolviese esta ecuación. Para ilustrar la posible existencia de excentricidades pequeñas asociadas a k y m grandes que producen nuevos puntos de corte, se analiza ese caso más desfavorable con m=k−1,ω=−90 ºyu2=90 º: e=−2±q4+6 m 3(4.82) Tomando +q4+6 m , puede obtenerse un límite superior a partir del cual es posible que se tengan nuevos puntos de corte aumentando únicamente la excentricidad, en función de la variable m . Por ejemplo, para m=1 se obtiene e=0.387 , mientras que para m=20 se obtiene e=0.02455 . Puede comprobarse en la siguiente gráfica cómo, efectivamente, para inclinaciones muy cercanas a cero, existe un punto de corte para órbitas con k=21 , m=20 , e=0.05 y ω=90 º , cuando, según el análisis del caso circular, deberían existir 21 −20 −1=0: Imagen 4.66 Generación de nuevos puntos de corte para inclinaciones cercanas a cero con excentricidades pequeñas y kymgrandes. Por debajo de la excentricidad que se obtiene de (4.82) para cada m (o k y m si se resuelve (4.81) completa) puede asegurarse que no se introducen nuevos puntos de corte a partir de tangencias en las trazas al introducir la excentricidad. Introducción de soluciones a través de las asíntotas Para que puedan introducirse nuevas soluciones a través de las asíntotas debe cumplirse que (4.76) sea igual a 0 para valores de u2 de la forma (2p−1)π 2 . Para valores genéricos de r , ω y e , esta función toma el valor ” ±∞ ”, al tratarse de unos valores asociados a las asíntotas de la segunda tangente. Para que esta función pudiese valer 0, debe cumplirse que la primera tangente también tienda a ∞ , pero con signo opuesto, de forma que deben encontrarse valores de pyrdentro de los rangos aceptables que cumplan: π 2−m kπ 2−(2p−1)π 2−πr= (2g−1)π 2(4.83)
76 Capítulo 4. Problema excéntrico Debe tenerse en cuenta que, para que realmente se llegue a una indeterminación del tipo ∞−∞ debe cumplirse que 2p−1 2 y 2g−1 2 puedan expresarse de la forma 1 2+2q o 2w−1 2 , pero ambos de la misma forma, todo ello siendo tanto pcomo q,wygnúmeros enteros y con pcontenido en ∈[0,1]. Resolviendo (4.83) se llega a: r=1−p+k(g−1) m(4.84) Para que el valor de r sea entero, los dos valores que puede tomar g son g=m+1 o g=−m+1 , y tomando p de forma que pueda llegarse a una indeterminación del tipo ∞−∞ en función del valor de m se llegan a las siguientes expresiones de r: r=k(4.85) r=1−k(4.86) r=k+1(4.87) r=−k(4.88) El único caso que permitiría obtener un valor de r dentro del rango válido sería k=1 y m=1 a partir de la ecuación (4.86) . Tomando p=0 se cumple que ambos números 2p−1 2 y 2g−1 2 pueden expresarse de la misma forma, de modo que se procede a estudiar este caso. Si se analiza la indeterminación del tipo ∞−∞ para este caso mediante L’Hopital, se llega a: (∂f ∂u2−1)(sin(u2)sin(f(u2))+cos(u2)cos(f(u2)) −sin(f(u2))cos(u2)∂f ∂u2−cos(f(u2))sin(u2)(4.89) Sustituyendo los valores para los cuales se está realizando el estudio, es decir, aquellos que consiguen anular los denominadores de ambas funciones trigonométricas, se llega a: ∂f ∂u2−1 0(4.90) Para que (4.90) siga siendo una indeterminación debe cumplirse que ∂f ∂u2=1 . Teniendo en cuenta que se están analizando casos del tipo u2=(2p−1)π 2 , la expresión de la derivada (4.91) y que para que r se encuentre en un rango válido, tanto kcomo mdeben valer 1, se llega a la ecuación (4.92): ∂f ∂u2 =−m k2esin(ω)sin(u2)−1−3e2 4cos(2ω)cos(2u2)(4.91) 2esin(ω)+ 3e2 2cos(2ω) = 0(4.92) La ecuación (4.92) lleva asociada a la solución e=0 y puede llevar asociada otra solución válida de e en función del valor de ω . La primera solución e=0 se debe a que para una órbita geosíncrona ( k=1 y m=1 ) y circular, la primera inclinación crítica es i=0º (puede verse en [1]), de forma que ya se han introducido nuevos puntos de corte para esta inclinación. Esto mismo ocurre para aquellos valores de ω que proporcionan un valor de excentricidad menor que 1, es decir, que para la combinación de elementos orbitales e y ω tales que se resuelve (4.92) , se da una tangencia en la traza para i=0º , es decir, que i=0º es una inclinación crítica. Además de este caso de tangencia, para otras combinaciones de e y ω pueden darse las siguientes situaciones: • Que a partir del método explicado en la sección correspondiente, se obtenga una inclinación crítica i∗>0º , de forma que para i=0º , se tengan 0 puntos de cruce. Un ejemplo de esta combinación puede verse en una órbita con k=m=1,e=0.05 yω=20. • El caso contrario, que a partir del método previamente descrito, no se obtenga ninguna inclinación crítica viable (todos los valores de cos i que se obtienen están fuera del rango [0,1] ). Esto quiere decir que la traza puede tener un punto de corte desde i=0º , sin que este punto de corte sea una tangencia. Un ejemplo de esto se tiene para k=m=1,e=0.2yω=−3rads.
4.5 Determinación del número de inclinaciones críticas de forma analítica 77 Por tratarse de un caso muy concreto, el estudio más detallado de las órbitas geosíncronas se propone como trabajo futuro. Con esto, puede asegurarse que la ecuación (4.76) no toma el valor 0 para valores de u2 de la forma (2p−1)π 2 para valores de k y m distintos a 1, y para las órbitas geosíncronas debe tenerse en cuenta que su primera inclinación crítica puede ser 0. Introducción de soluciones a través de los extremos −π,π Ahora debe comprobarse si (4.76) se iguala a 0 para valores de u2=±π . Para ello se buscan valores enteros de gyrdentro de sus respectivos rangos de forma que se cumpla: π 2−m kπ 2−2esin(ω)±π−πr=gπ(4.93) Para que (4.93) lleve asociado un valor de r entero, debe cumplirse que e=0 o ω=0/ω=±180 º , es decir, elementos orbitales que conserven la simetría N-S. Despejando r, y tomando g=m+1 2(para lo cual m debe ser impar), se llega a: r=±1+k+1 2(4.94) Para que (4.94) lleve asociado un valor de r entero, deberá cumplirse que k sea impar, de forma que tanto m como k lo sean. Esto nos lleva a que los únicos casos para los cuales se tiene una solución en los extremos del intervalo es para los casos impares con simetría N-S, soluciones que se encontrarán además sobre el Ecuador. Este tipo de órbitas llevarán asociadas otra solución con u2=0 y un r distinto, de forma que no será necesario contabilizarlas, ya que llevarán asociada la misma latitud ( φ=0º ), y podrán obtenerse mediante un retraso nodal. Teniendo todos los casos en cuenta, se ha llegado a que la introducción de la excentricidad (por debajo de los valores obtenidos en la sección de puntos críticos) no aumenta el número de soluciones de la ecuación (4.76), de forma que se siguen teniendo k-m-1 puntos de corte para el caso impar e ídem para el caso par. 4.5.2 Particularización para i=180º En este caso, la ecuación que se analiza es: F(u2) = tan(π 2−m kπ 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2)−πr)+ tan(u2)(4.95) En una ecuación del tipo tan(x)+tan(y) = 0 la relación que cumplen las variables x e y es x=pπ−y , de forma que siguiendo los pasos descritos anteriormente, se llega a una ecuación del tipo: u2 k+m k=(2p−1)π 2+m kπ 2−2esin(ω)cos(u2)−3e2 4cos(2ω)sin(2u2)−πr(4.96) Al igual que para el caso i=0º , primero se comprueba que el número de soluciones de (4.96) coincide con las del caso circular cuando se hace e=0: u2=π 2(k+m)k(2p−1)−m(2r−1)(4.97) Se sabe por el estudio realizado en [1], que el número de soluciones distintas que la ecuación (4.97) lleva asociadas son k+m−1 , por lo que de nuevo deben comprobarse los 3 casos expuestos anteriormente, para comprobar que la introducción de excentricidad no genera nuevas soluciones en (4.95): Existencia de puntos críticos De forma análoga al caso i=0º , se estudia cuando se anula la derivada de (4.96) respecto a u2 , cuya expresión se corresponde con (4.98), pero se simplifican de nuevo las funciones trigonométricas para llegar a (4.99): k+m=m2esin(ω)sin(u2)−3e2 2cos(2ω)cos(2u2)(4.98) 3e2 2−2e+k m+1=0(4.99)
78 Capítulo 4. Problema excéntrico Para pasar de (4.98) a (4.99) se tiene en cuenta el caso más desfavorable al ser e<1, es decir, que el término asociado a esea positivo y el asociado a e2negativo (ya que ambos no pueden ser positivos). A partir de (4.99) se llega a la expresión (4.100) , la cual lleva asociados valores de excentricidad imaginarios, por ser siempre k m+1>1, y por tanto, valores no válidos: e= 2±q4−6(k m+1) 3(4.100) Con esto se demuestra que no podrán darse puntos críticos tras la introducción de la excentricidad. Introducción de soluciones a través de las asíntotas Siguiendo un procedimiento completamente análogo al comentado en el caso de i=0º , pero teniendo en cuenta que ahora los términos 2p−1 2 y 2g−1 2 deben poder expresarse como 1 2+2q y 2w−1 2 o viceversa, pero no de la misma forma, para que realmente pueda darse la indeterminación del tipo ∞−∞ y que los valores que puede tomar gson, de nuevo g=m+1og=−m+1, se llegan a las mismas expresiones de r: r=−k(4.101) r=k+1(4.102) r=k(4.103) r=1−k(4.104) Puede comprobarse que en ningún caso, las ecuaciones (4.101) - (4.104) aportan valores de r dentro del intervalo válido, ya que el único valor posible sigue siendo r=1−k con k=m=1 , pero en este caso, para que los términos 2p−1 2 y 2g−1 2 puedan expresarse de forma distinta, al tomar g=−m+1 se deberá tomar p=1 , con lo que la expresión de r a la que se llegaría sería r=−k , en lugar de r=1−k . Con esto, puede asegurarse que no se generan nuevas soluciones a través de las asíntotas al introducir la excentricidad en la ecuación (4.95). Introducción de soluciones a través de los extremos −π,π De nuevo, el argumento de la primera tangente es completamente análogo, con lo que se llegan a los mismos resultados que para el caso i=0º : los únicos casos en los que se tiene una solución en los extremos del intervalo son los casos impares que conservan su simetría N-S, soluciones que no están siendo tenidas en cuenta en el cómputo. De nuevo, puede asegurarse que no se introducen nuevas soluciones al tenerse una órbita ligeramente excéntrica, con lo que el número de soluciones para i=180 ºse mantiene en k+m−1. 4.5.3 Determinación del número de soluciones introducidas por cada inclinación crítica Una vez demostrado que el número de puntos de corte al introducir la excentricidad no varía, será necesario obtener el número de soluciones que introduce cada inclinación crítica. Para hacer esto se realiza un desarrollo en serie de Taylor alrededor de la solución u∗ 2 asociada a esta inclinación crítica i∗de la ecuación maestra aproximada: F(u2,i) = tanπ 2−m k(π 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2) −πr)−tan(u2)cos(i)(4.105) Se calculan los tres primeros coeficientes para el desarrollo de la función F con respecto a la variable u2 , teniendo en cuenta que se realiza la siguiente denominación: f(u2) = π 2−m k(π 2−2esin(ω)cos(u2)−u2−3e2 4cos(2ω)sin(2u2)−πr)(4.106)
4.5 Determinación del número de inclinaciones críticas de forma analítica 79 F0=tan(f)−tan(u2)cos(i)(4.107) ∂F ∂u2 =∂f ∂u2 (1+tan2(f))−(1+tan2(u2))cos(i)(4.108) ∂2F ∂u2 2 =∂2f ∂u2 2 (1+tan2(f))+ 2(∂f ∂u2 )2tan(f)(1+tan2(f))−2 tan(u2)(1+tan2(u2))cos(i)(4.109) Estas expresiones pueden simplificarse teniendo en cuenta que para que i∗ y u∗ 2 estén relacionadas, deben cumplir el siguiente sistema de ecuaciones: F=tan π 2−m k(π 2−2esin(ω)cos(u∗ 2)−u2−3e2 4cos(2ω)sin(2u∗ 2)−πr∗)−tan(u∗ 2)cos(i∗) = 0 ∂F ∂u2=∂f ∂u2|u∗ 2(1+tan2(f(u∗ 2)))−(1+tan2(u∗ 2))cos(i∗) = 0(4.110) Con esto, se reescriben las expresiones de las derivadas anteriores: F0=tan(u∗ 2)(cos(i∗)−cos(i)) (4.111) ∂F ∂u2|u∗ 2= (1+tan2(u∗ 2))(cos(i∗)−cos(i)) (4.112) ∂2F ∂u2 2|u∗ 2=∂2f ∂u2 2 (1+tan2(u∗ 2)cos2(i∗))+ ∂f ∂u2 2 tan(u∗ 2)cos2(i∗)(1+tan2(u∗ 2)) −2 tan(u2)(1+tan2(u2))cos(i)(4.113) Una vez calculadas estas derivadas y particularizadas para el punto central del desarrollo de Taylor u2=u∗ 2 , puede describirse la función F en la cercanía de este punto como F≈A+B(u2−u∗ 2)+C(u2−u∗ 2)2 , donde: A=tan(u∗ 2)(cos(i∗)−cos(i))(4.114) B= (1+tan2(u∗ 2))(cos(i∗)−cos(i))(4.115) C=−m k2esin(ω)cos(u∗ 2)+ 3e2cos(2ω)sin(2u∗ 2)(1+tan2(u∗ 2)cos2(i∗)) − m k2esin(ω)sin(u∗ 2)−1−3e2 2cos(2ω)cos(2u∗ 2)cos2(i∗)+ cos(i)!2 tan(u∗ 2)(1+tan2(u∗ 2))(4.116) Para ver cuántas soluciones se introducen al superarse una de las inclinaciones críticas, se estudian las soluciones de la ecuación A+B(u2−u∗ 2)+C(u2−u∗ 2)2=0 para i>i∗ y i<i∗ . Observando el caso circular, lo que se espera es que para las órbitas con i>i∗ las soluciones de esta ecuación sean dos reales, mientras que para i<i∗no exista ninguna solución real. Resolviendo la ecuación cuadrática, se tiene: u2−u∗ 2=−B±√B2−4AC 2C(4.117) Teniendo en cuenta que se va a resolver la ecuación (4.117) para valores de inclinación de la forma i=i∗±di siendo di un valor infinitesimal, el término asociado a B2 puede despreciarse, al ir acompañando a un término de orden menor di2. Para comprobar la exactitud de los cálculos, primero se verifica que la ecuación (4.117) tiene 2 soluciones reales para i>i∗y ninguna para i<i∗con e=0. Para ello, se particulariza 4AC para e=0: 4AC =8 tan2(u∗ 2)(cos(i)−m kcos2(i∗))(1+tan2(u∗ 2))(cos(i)−cos(i∗)) (4.118)
80 Capítulo 4. Problema excéntrico El término cos(i)−m kcos2(i∗) es siempre positivo en las inmediaciones de i=i∗ , como puede comprobarse en las siguientes gráficas: (a) Centrado en i<i∗.(b) Centrado en i>i∗. Imagen 4.67 Función cos(i)−m kcos2(i∗)con m kgenérico. (a) Centrado en i<i∗.(b) Centrado en i>i∗. Imagen 4.68 Función cos(i)−m kcos2(i∗)con m k≈1. Puede comprobarse que, incluso en el caso más desfavorable con m k≈1 , la función cos(i)−m kcos2(i∗) es positiva en las inmediaciones de i=i∗(a excepción de para i=π 2, donde vale 0). Dicho esto, el término cos(i)−cos(i∗) deberá ser negativo si i>i∗ y positivo si i<i∗ , lo cual se cumple, ya que la función cos es decreciente con el ángulo en el intervalo [0,π 2] . Con esto, queda demostrado que las expresiones (4.114) - (4.116) son válidas cuando se particularizan para el caso circular. Ahora sería necesario comprobar que el término −4AC cumple la condición de ser positivo para i>i∗ y negativo en caso contrario, al introducir un valor de excentricidad cualquiera. Escribiendo de forma contraída el término −4AC: −4 tan(u∗ 2)(cos(i∗)−cos(i))C(4.119) El término cos(i∗)−cos(i) ya se ha comprobado que es menor que 0 para i<i∗ y mayor que 0 en caso contrario. Para que existieran soluciones reales, deberá cumplirse que C sea menor que cero para i>i∗ siempre que tan(u∗ 2)>0, y menor que cero en caso contrario. Debido a la complejidad del término C , se analiza un caso concreto, que se ha considerado el más desfavorable: el caso en el que u∗ 2=−3π 4yω=π 2. En este caso, la expresión del término C queda: C=m k√2(1+5 cos(i∗))e+4(m kcos2i∗−cos(i)) (4.120)
4.5 Determinación del número de inclinaciones críticas de forma analítica 81 Si se calcula la excentricidad que hace que la expresión (4.120) se anule con i=i∗se llega a la siguiente distribución, la cual ha sido truncada en e=1: Imagen 4.69 Valor de excentricidad que hace que se anule el término Ccon u∗ 2=−3π 4yω=π 2. Observando 4.69 puede llegarse a una serie de conclusiones: • Para valores de m k pequeños, la excentricidad necesaria para anular C es mayor que 1, con lo que se puede asegurar que al introducir la excentricidad, cada inclinación crítica llevará asociadas dos soluciones. • Conforme el valor de m k aumenta, la excentricidad necesaria para que pueda anularse C disminuye, pero solo para valores de m k muy cercanos a 1 se llegan a valores de excentricidad que podrían considerarse dentro del rango de ”excentricidad pequeña”, región para la cual es válida la ecuación a partir de la cual se ha extraído la expresión de C. • Cuanto menor sea la inclinación crítica, menor será la excentricidad necesaria. Como puede verse a lo largo de todo el apartado dedicado a las inclinaciones críticas, estas se encuentran en el intervalo [45 º ,90 º ] para los valores de m y k estudiados, de forma que observando la gráfica incluso para valores de m k cercanos a 1, puede asegurarse que las excentricidades necesarias para anular el término C son lo suficientemente grandes como para no considerarlas aceptables en nuestro modelo. Para valores de m y k tales que m k≈1 y tales que se obtenga una inclinación crítica cercana a 0, la excentricidad necesaria para anular el término C estaría dentro del rango considerado como ”de excentricidad pequeña”, de forma que, al igual que con el número de puntos de cruce para i=0º , existirían valores de m y kgrandes para los cuales el modelo tendría problemas. Visto que para uno de los casos más desfavorables la excentricidad necesaria para que los resultados del caso circular puedan no cumplirse es muy elevada, se asume, que al utilizar valores de excentricidad pequeños, cada inclinación crítica introducirá dos soluciones, de forma completamente análoga al caso circular. 4.5.4 Cálculo del numero de inclinaciones críticas Una vez realizado todo lo anteriormente expuesto, pueden obtenerse el número de inclinaciones críticas para cada órbita: teniendo en cuenta que para i=0ºse tienen k−m−1soluciones, que para i=180 ºse tienen k+m−1 , y que por cada inclinación crítica se introducen 2 nuevas soluciones, es fácil obtener la conclusión de que el número de inclinaciones críticas es m . La veracidad de este resultado puede verse a través de los siguientes ejemplos, calculados numéricamente (y ya expuestos anteriormente):
82 Capítulo 4. Problema excéntrico Caso excéntrico: Caso impar: k=7, m=5, ω=70 º, e=0.03 Tabla 4.12 Inclinaciones críticas para órbita con para k=7, m=5, e=0.03, ω=70 º. i∗ 1=52.3688 ºi∗ 2=79.0059 ºi∗ 3=80.3424 ºi∗ 4=88.9167 ºi∗ 5=89.0336 º Caso par: k=7, m=4, ω=70 º, e=0.03 Tabla 4.13 Inclinaciones críticas para órbita con para k=7, m=4, e=0.03, ω=70 º. i∗ 1=75.2376 ºi∗ 2=77.2457 ºi∗ 3=88.6388 ºi∗ 4=88.7867 º Caso par: k=8, m=3, ω=−20 º, e=0.05 Tabla 4.14 Inclinaciones críticas para órbita con para k=8, m=3, e=0.05, ω=−20 º. i∗ 1=83.0339 ºi∗ 2=83.5845 ºi∗ 3=90 º Caso circular: Caso impar: k=7, m=5 Tabla 4.15 Inclinaciones críticas para órbita con para k=7, m=5 y e=0. i∗ 1=44.4153 ºi∗ 2=79.7077 ºi∗ 3=88.9774 º Caso par: k=7, m=4 Tabla 4.16 Inclinaciones críticas para órbita con para k=7, m=4 y e=0. i∗ 1=76.3061 ºi∗ 2=88.7157 º Caso par: k=8, m=3 Tabla 4.17 Inclinaciones críticas para órbita con para k=8, m=3 y e=0. i∗ 1=83.3402 ºi∗ 2=90 º Visto que los resultados analíticos se corresponden con los obtenidos numéricamente en apartados anteriores, puede asumirse que todo el desarrollo es válido para órbitas de excentricidad pequeña (sujeto a las condiciones de mykgrandes). 4.6 Precisión y limitaciones del modelo Un apunte importante que debe hacerse es la delimitación que supone las llamadas ”órbitas de excentricidad pequeña”. A lo largo de todo el documento, el mayor valor de excentricidad utilizado está en torno a e=0.3, y únicamente para obtener las gráficas de la evolución de las inclinaciones críticas de las órbitas. Además, anteriormente, se había comprobado la exactitud de la solución de la ecuación maestra en comparación con la obtenida del sistema de ecuaciones completo, viendo que para valores menores de e=0.15, los resultados obtenidos eran muy precisos. Dicho esto, cabe preguntarse, ¿hasta dónde es válida la aproximación realizada para excentricidad pequeña? Es decir, ¿cuándo las soluciones obtenidas por la ecuación maestra aproximada comienzan a distanciarse de las que se obtienen resolviendo el sistema de ecuaciones completo? Para responder
4.6 Precisión y limitaciones del modelo 83 a estas cuestiones se incluyen unas tablas comparativas para órbitas con k=3, m=2, i=83 º y ω=35 º . La elección de los valores de estos parámetros simplemente se basa en tomar una órbita que presente pocos puntos de corte, para que sea más sencillo comparar los valores que se obtienen como solución. Tabla 4.18 Tabla comparativa para órbitas con k=3, m=2, i=83ºyω=35 ºcon e=0.25. k=3, m=2, e=0.25, ω=35 º, i=83º,Ω=0u∗ 2[rad] r[revs] φ[º]λ[º] Ecuación maestra aproximada 0.7073 1 40.1584 123.7603 1.2845 1 72.2042 126.9271 1.8571 2 72.2042 6.9271 2.4343 2 40.1584 3.7603 Sistema de ecuaciones completo 0.7127 1 40.4660 123.7512 1.2818 1 72.0609 126.9450 1.8598 2 72.0609 6.9450 2.4289 2 40.4660 3.7512 En esta primera tabla puede verse cómo, a pesar de que la excentricidad se ha tomado como e=0.25, los resultados de ambos métodos solo comienzan a diferenciarse en el segundo diferencial, es decir, se tienen errores en la centésima de grado o de radián, lo que implica resultados aún precisos. Tabla 4.19 Tabla comparativa para órbitas con k=3, m=2, i=83ºyω=35 ºcon e=0.35. k=3, m=2, e=0.35, ω=35 º, i=83º,Ω=0u∗ 2[rad] r[revs] φ[º]λ[º] Ecuación maestra aproximada 0.5105 1 29.0114 125.7585 1.3107 1 73.5698 131.6377 1.8309 2 73.5698 11.6377 2.6311 2 29.0114 5.7585 Sistema de ecuaciones completo 0.5116 1 29.0736 125.6259 1.3058 1 73.3152 131.6951 1.8358 2 73.3152 11.6951 2.6300 2 29.0736 5.6259 Al aumentar la excentricidad a e=0.35, puede verse cómo los resultados pasan a discrepar en la décima de grado en latitud y longitud, lo que los hace menos precisos en comparación con el caso anterior, como era de esperar. Tabla 4.20 Tabla comparativa para órbitas con k=3, m=2, i=83ºyω=35 ºcon e=0.45. k=3, m=2, e=0.45, ω=35 º, i=83º,Ω=0u∗ 2[rad] r[revs] φ[º]λ[º] Ecuación maestra aproximada 0.3269 1 18.5830 126.7022 1.3303 1 74.5776 135.8505 1.8112 2 74.5776 15.8505 2.8147 2 18.5830 6.7022 Sistema de ecuaciones completo 0.3062 1 17.4101 126.1864 1.3227 1 74.1891 135.9799 1.8189 2 74.1891 15.9799 2.8354 2 17.4101 6.1864 En este caso, si se sigue aumentando la excentricidad, siguen creciendo las discrepancias entre ambas soluciones: los valores de los argumentos de latitud u2 siguen discrepando en los segundos decimales, pero en este caso las diferencias son mayores. Esto hace que, en el caso de la latitud y la longitud, se estén cometiendo errores cercanos al grado, errores que pueden no ser aceptables en ciertas aplicaciones de observación geológica o geográfica, por ejemplo.
90 Capítulo 4. Problema excéntrico Imagen 4.81 Zoom de la gráfica 4.80. Con estas gráficas puede verse las límitaciones propias de la aproximación de excentricidad pequeña, la cual no considera esta creación de nuevos puntos únicamente por el aumento de excentricidad (concepto similar al de ”excentricidad crítica”).
5 Conclusiones y trabajo futuro 5.1 Conclusiones De forma general, puede afirmarse que, al agregar el factor de la excentricidad a la órbita no se producen cambios extremos en cuanto a la cantidad y disposición de los puntos de corte, siempre que esta sea ”pequeña”. Dicho esto, la introducción de la excentricidad complica la resolución del problema al no poder definir la variable z , de forma que el desarrollo analítico de las ecuaciones es más complejo, como ha podido comprobarse por ejemplo, al estudiar la variación de la posición de los puntos de cruce con la excentricidad o el argumento del perigeo. En cuanto a esta variación con los elementos orbitales, cabe destacar que, en general, la variación de la latitud y la longitud viene influenciada por el argumento del perigeo de la órbita y, además, también por la posición inicial del punto de corte (se constató que cuanto mayor fuese la latitud inicial, menor sería la variación de la misma). Por otra parte, respecto a las inclinaciones críticas, se ha visto cómo el número de estas aumenta siempre que se rompa la simetría N-S de la traza (es decir, ω6=0º o ω6=180º ), debido a que para las órbitas circulares, cuando se producía una tangencia en uno de los hemisferios, debido a la simetría de la traza, también se producía en el otro. Esto se ha demostrado con varios ejemplos, así como en un estudio de los valores de las inclinaciones críticas conforme se varía el argumento del perigeo y la excentricidad de las órbitas. Por último, cabe destacar las limitaciones del propio modelo, derivadas del uso de una ecuación aproximada mediante un desarrollo de Taylor alrededor de e=0 de la ley horaria de la órbita. Como ha podido comprobarse en la sección correspondiente, existen combinaciones de elementos orbitales para las cuales, la solución obtenida discrepa notablemente de la que se debería obtener (la cual puede verse al representar la traza) o incluso, fallar a la hora de obtener el conjunto completo de soluciones si la excentricidad de la órbita es lo suficientemente grande. Dicho esto, también es necesario indicar que, tal y como también se ha expuesto, esta aproximación ofrece muy buenos resultados para órbitas de pequeña excentricidad, quedando reflejada esta precisión en las tablas y gráficas de la sección 4.5. Además, y para el desarrollo analítico de las propias ecuaciones, el modelo se ve limitado a valores de m y k no excesivamente grandes, ya que estos modifican el número de puntos de corte para i=0 y pueden propiciar que la excentricidad necesaria para que el término C del desarrollo de Taylor de la ecuación maestra aproximada alrededor del punto de tangencia u∗ 2 tome un valor que no garantice la existencia de soluciones reales para i>i∗ , como ha podido verse a lo largo del documento. 5.2 Trabajo futuro A pesar de que este documento ilustra en gran medida la dependencia de los puntos de corte de la traza de las órbitas de pequeña excentricidad con respecto a los elementos orbitales de la misma, es posible seguir profundizando en este tema, para lo que se proponen posibles líneas de investigación futura: • Como ya se ha comentado, este estudio se centra en las órbitas de pequeña excentricidad, con lo que una de las líneas más interesantes podría ser la generalización de este estudio a cualquier tipo de órbita, incluidas las de alta excentricidad. Dentro de este estudio cabría destacar el análisis de lo que en este documento se ha denominado ”excentricidades críticas”, es decir, valores de excentricidad que hacen 91
92 Capítulo 5. Conclusiones y trabajo futuro que aparezcan nuevos puntos de corte sin variar el resto de elementos orbitales, ya que estas no pueden estudiarse con la ecuación maestra aproximada. • Por otra parte, durante los desarrollos realizados en el documento, no se ha tenido en cuenta ninguna de las perturbaciones a las que se ven sometidas las órbitas, las cuales hacen que el satélite se aleje de su órbita ”ideal” o inicial. Por esto, otra línea de investigación podría ser la inclusión de estas perturbaciones en los desarrollos, de forma que cuantas más perturbaciones se incluyan, más preciso será el cálculo de los puntos de corte. Algunos ejemplos relevantes de estas perturbaciones a incluir podrían ser la resistencia atmosférica para órbitas bajas, o incluso, para órbitas más alejadas de la Tierra, la influencia de un tercer cuerpo en la órbita, como la Luna o el Sol. • El estudio en profundidad del número de inclinaciones críticas de forma analítica, teniendo en cuenta otros casos de u∗ 2y resolviendo la posible discrepancia obtenida para kymgrandes. • El estudio de las órbitas geosíncronas para inclinaciones cercanas a i=0º , para determinar la dependencia del número de puntos de corte y de las inclinaciones críticas con la excentricidad e y el argumento del perigeo ω. • Por último, una línea de investigación muy interesante podría ser la resolución del problema inverso, es decir, la determinación de los elementos orbitales necesarios para que la traza tenga un número concreto de puntos de corte y unas ubicaciones específicas sobre la superficie terrestre.
6 Referencias [1] Eva Tesón Muñoz. ”Estudio detallado de los puntos de cruce de las trazas de satélites geocéntricos”. Trabajo de fin de grado. Departamento de Ingeniería aeroespacial, Universidad de Sevilla, 2015. [2] John P. Boyd. ”Computing the zeros of a Forier series or a Chebyshev series or general orthogonal polynomial series with parity symmetries”. ScienceDirect. Computers and Mathematic with Applicaions 54, 2007. Department of Atmospheric, Oceanic and Space Science. University of Michigan. [3] Rafael Vázquez Valenzuela. Problemas y cuestiones de la asignatura Mecánica Orbital y Vehículos Espaciales, problemas de los tema 3 y 4, problema 13. Departamento de ingeniería aeroespacial, Universidad de Sevilla, 2016. http://www.aero.us.es/move/index.html [4] Chebfun website. MathWorks, OCCAM, EPSRC, erc. Package downloaded on 06/04/2018 from http://www.chebfun.org/download/ [5] M.C. Kim. ”Theory of satellite ground-track crossovers”, Journal of Geodesy, Volume 71, 1997. [6] Davis A. Vallado. Fundaments of astrodynamics and applications. Springer Science & Business Media, 2007. 93