Cálculo das Órbitas Planetárias (versão português)

Andrés Esteban de la Plaza

1

Cálculo das Órbitas Planetarias por Ing. Andrés Esteban de la Plaza © De livre publicação desde que o autor seja mencionado, direitos reservados. Última revisão 29 junho 2016 – Rio de Janeiro, Brasil. [email protected]

Homines dum docent discunt. - Séneca o Moço (Córdoba 4 a.C.- Roma 65 d.C.)

ÍNDICE: 1. REPRESENTAÇÃO DE VETORES 2. ACELERAÇÃO EM COORDENADAS POLARES 3. FORÇAS CENTRAIS E GRAVITAÇÃO 4. RESOLUÇÃO DE EQUAÇÃO DIFERENCIAL DO PROBLEMA 5. CONSIDERAÇÕES PARA AS ÓRBITAS ELÍPTICAS 6. CÁLCULO DA POSIÇÃO NA ÓRBITA ELÍPTICA 7. CONSIDERAÇÕES SOBRE A CONSTANTE DE GAUSS, O MOVIMENTO DIURNO MÉDIO E A EQUAÇÃO DE KEPLER 8. CÁLCULO DE ÓRBITAS NÃO ELÍPTICAS: ÓRBITAS PARABÓLICAS E HIPERBÓLICAS

1

Andrés Esteban de la Plaza

2

1. REPRESENTAÇÃO DE VETORES 

A posição de um ponto P no espaço pode ser definida por um vetor posição P com origem num centro O de coordenadas ortogonais retangulares. Definindo nesta terna xyz, os versores  (vetores unitários) dos eixos, obtemos respectivamente: o versor i como versor da direção x, o   versor j como o versor da direção y, e o versor k como o versor da direção z. Assim, qualquer vetor definido num destes eixos estará representado pelo produto de um escalar igual ao     módulo e do versor da direção  X  X  i  x  i .



Portanto, podemos decompor o vetor P em    seus componentes ortogonais x P , y P e zP no



sistema de centro O: ( projeções de P nos eixos    x,  y, e z ) de forma que P  xP  yP  zP  xP i  yP j  zP k . Também podemos representar este vetor de forma matricial: x    P P  yP   zP 

 r P 

 x    y    z  2

P

2

P

2

P

Porém esta não é única  representação possível para o ponto P. Podemos imaginar um plano que contém o vetor P e o centro O, de maneira que, embora ainda continuemos no espaço tridimensional, a representação de P é feita num plano, entretanto, este plano será um plano complexo ou seja, o eixo x será real, e o eixo y será imaginário, ou seja: A vantagem de procedermos assim está dada pela possibilidade de representarmos um número complexo na sua forma exponencial, como um vetor. Ë fácil observar  que se definimos o versor i para a direção x, então a  direção de r estará definida por um versor r tal que

  r  e i  i

onde i =  1 , logo:     i r  e i  cos   i sen    i  i  cos   i i sen     e como i i  j pois multiplicar por i equivale a girar i 90 :    r  i  cos   j  sen  de forma que:

    P = r  r  r  r  e i  i que é a expressão do raio vetor em coordenadas polares.

2

Andrés Esteban de la Plaza

3

2. ACELERAÇÃO EM COORDENADAS POLARES







i



Podemos agora proceder a derivar P = r  r  r  r  e  i respeito de tempo para  encontrar a aceleração. Porém antes vamos calcular a derivada respeito do tempo de r , pois é   uma função de , ou seja que r = r ():

 i  d e  i i d e  i d      d r      r      ie  i i    n      n onde n  versor normal a r dt dt d dt    logo: r    n  i      d n d n d  d  ie i  d   n      i 2 e  i i    e  i i    r   dt d dt d dt    logo: n  r     vamos calcular então a primeira derivada de r , ou seja r :

   dr  d r d rr  d r    r   r  r   r  r  r   n dt dt dt dt    r  r  r  r   n   vamos calcular agora a segunda derivada de r , ou seja r :

     2  d r d r  d   d n  d r d r d  r  r  r  n d r  r    r  r    n  r   n  r   2  dt dt dt dt dt dt dt dt      n  r    r    r  r + r   n  r   n  r         rr  2r  n  r n  r 2 r  r  r  r 2 r  2r   r  n

    n r   r  r 2 r  2r  r onde:

a r  r  r 2  aceler. normal: a n  2r   r    aceler. total en coordenadas polares: a  a r r  a n n

aceler. radial :

3

Andrés Esteban de la Plaza

4

3. FORÇAS CENTRAIS E GRAVITAÇÃO Inicialmente vamos deduzir as fórmulas que determinam as trajetórias dos corpos que se encontram sujeitos à atração de forças centrais, no caso particular da a lei da gravitação universal, ou seja

 M m F  G S2 r . De acordo com a figura 1 temos: r Onde o ponto m está submetido à força :

 Mm  F  G 2 r r  r é o versor (vetor unitário) da direção de F e:

  r  r( t )   ( t )

t = tempo MS = massa corpo central m = massa objeto A resolução deste problema deve proporcionar as seguintes funções:

r  r(  )   ( t ) As expressões para a aceleração em função das componentes normal e radial da aceleração são: 2     d r   ar 2  a rr  a nn dt

onde

d2 r d 2 2  a r   r  r   r ( ) d t2 dt dr d d2    a n = 2r   r  2  r dt dt d t2

Mas nos sistemas de forças centrais sabemos que a aceleração normal é nula  an=0, portanto a aceleração será: 2  d r    a  r = 2  arr dt 2 d r d 2 2   a r  r - r  r ( ) dt d t2

e lembrando que a equação diferencial geral do movimento central neste caso é:

4

Andrés Esteban de la Plaza

5

 MSm  d2 r m 2  G 2 r dt r 4. RESOLUÇÃO DA EQUAÇÃO DIFERENCIAL DO PROBLEMA A equação diferencial geral do problema pode ser decomposta em duas particulares:

M m   mar r  an n    G s2 r r

e

an  0

temos então: mar =  G man=

M sm r2

0

 ar  G

 an 

Ms r2

0

(pois m  0)

[2]

então a equação diferencial que temos que resolver será:

d2 r d 2 2  ar  r  r   r ( ) d t2 dt

[1]

porém não conhecemos as relações r  r( t ),  = ( t) , mas lembrando da equação [2] obtemos que:

an = 0  2r  r  0 ou seja

dr d d2  2  r 2 0 dt dt dt De

  0 2r   r

entretanto

dr dr d   dt d dt



r  r 

 , ou seja: podemos obter 

     2r   = d   e como  r dt    = d    2r   mas r  r    dt r  2r  2r  2 dr    e como r   r r d  d dr 1    2    [2a] dt d r fazendo passagem de términos na [2 a] obtemos: 5

Andrés Esteban de la Plaza

6

d d r 2  d   2   dt r d r 2   d   2   r  

dr  d   2   r d  dr   2  r

[3]

d  dr    2 r

integrando agora a [3]:

ln   2 ln r  constante 

ln   2 ln r  C1  constante ln   ln r 2  C  constante 1

ln(   r 2 )  C1  const. introduzindo a constante 1/ 2

portanto 

1 ln(   r 2 )  C 2  const. 2

1 2   r    e c2  constante 2

que não é outra coisa senão a velocidade areolar do ponto m na sua trajetória, como podemos observar na figura 2:

6

Andrés Esteban de la Plaza

7

 Se a velocidade areolar é constante, então: 

“O raio vetor que une o planeta ao sol varre áreas iguais em tempos iguais”.

A área do triângulo será:

1 rrd  2 1 d A   r2  d  2

dA=

Vamos fazer agora algumas considerações energéticas. Para isto definimos:

 r  raio vetor m  massa   p  quantidade de movimento LINEAR = m  v    L  quantidade de movimento angular ou momento da quantidade de movimento linear = r  p

 L  chamando a l m

   1   1     temos que l  r  v  ( r  p )  ( r  mv )  r  v m m  1   1 2 d A  mas a metade do modulo de l vale rv  r   A 2 2 dt 1   1 r 2 logo l  A 2 2 L ou  r 2 que como sabemos = constante, portanto m

L L   constante   = m mr 2

L   r m 2

[3a]

podemos determinar agora a partir da [3a] o valor de r 2 pois

r     2

r

2

r3 L2 2  r  2 3 mr

2

 L    m  r3

2

[4]

Logo, para resolver a equação diferencial [1], já temos o valor r 2 porém falta conhecermos o valor de r . Para isto, e voltando agora à aceleração radial:

7

Andrés Esteban de la Plaza

8

M r  r 2  G 2s r

[5]

observamos que:

r 

dr dr d dr       dt d dt d dr L   d  mr 2 dr L  2  r md   1 L  d    r md  

L d  1    m d   r

logo r 

d r d r d    dt d dt d r     d d  L d  1   L          d   m d   r    mr 2  r = 

L2 d 2  1    m 2r 2 d  2  r 

[6]

Agora, substituímos [6] e [4] na eq. diferencial [5] e obtemos então:

r



r 2 

L2 d 2  1 L2  2 2  2   -r 2 4  m r d   r mr L2 m2

ou seja Fazendo z 

 d  1 1   2      d   r r 

Ms r2 M  G 2s r G

GM s 2

d 2  1 1 M sm2    G d 2  r  r L2 1 e substituindo, chegamos a uma equação diferencial [7] de fácil resolução: r M sm2 d 2  1 1   G   d 2  r r L2 M sm2 d2 z  z  G d 2 L2

8

Andrés Esteban de la Plaza

9

M sm 2 z  z  G L2

[7]

A solução geral da equação [7] compõe-se de duas soluções; a solução homogênea zh que GM s m 2 provém de resolver z  z  0 , e a solução particular zp  . L2  Solução homogênea:

z h  e  (C 1 cos   C 2 sen  ) ou seja

onde  = 0

e

 =1

z h  C 1 cos   C 2 sen 

introduzindo agora as constantes  0 e Q de forma que:

C1  Q cos  0

e

C 2  Q sen  0

chegamos a

tan  0 

C2 C1

z h  Q cos(    0 )

 Solução particular:

Msm2 zp  G L2  Solução completa:

Msm2 z  zh  zp  Q cos(   0 ) + G L2 de forma que

Msm2 1  Q cos(   0 ) + G r L2 operando chegamos a

1  r



Msm2 L2 Q cos(   0 ) + 1G GM s m 2 L2

ou seja

9

Andrés Esteban de la Plaza

10

L2

GM s m 2 r L2 Q 1 2 cos( -  0 ) GM s m

[8]

que obviamente corresponde à equação de uma cônica (Figura 3), cuja fórmula geral é:

r

p

1  e cos   0 

p  1  e  f

[8aa]

ou seja que: 2 L 1  e f 

GM s m 2

[8 a]

10

Andrés Esteban de la Plaza

11

5. CONSIDERAÇÕES PARA AS ÓRBITAS ELÍPTICAS Figura 4

Os parâmetros da elipse na Figura 4 são: a: semi-eixo maior b: semi-eixo menor c=a2-b2 e: excentricidade = c/a f: distância focal = a( 1 e ) A: área da elipse = ab= a 2 1 e 2 T: período para uma revolução

A equação da elipse em coordenadas cartesianas é:

A equação da elipse em coordenadas polares é:

( x  c )2 y 2  2 1 a2 b

p r 1  e cos 

b2 p  a( 1  e 2 ) a

A velocidade com que o raio vetor r varre a área da órbita vale:

ou seja

  d A  1 r 2  A dt 2   L A [9] 2m

mas

1 2 L r   constante 2 2m

se agora integramos a [9] no tempo T obteremos a superfície total A: T

T

T

 dt L dt L dt L T AA 0 2m 2m 0 2m 0

11

Andrés Esteban de la Plaza

12

A  a

LT 1 e  2m

2

2

agora vamos isolar T:

2ma 2 1  e 2  T L L2 agora vamos calcular T2 lembrando que f ( 1  e )  GM s m 2 , ou seja:

4 2 a 3 T  GM s 2

[10]

porém para o sistema solar, Ms=massa solar, a=distância média do planeta ao sol, T=seu

4 2 período de revolução. Olhando a eq. [10] é óbvio que a quantidade GM s é uma constante T12 T22 4 2   para todos os planetas do sistema solar. Assim temos que onde os a 13 a 32 GM s subíndices 1,2, representam os diferentes planetas. Tal é a 3 a. Lei de Kepler: “Os quadrados dos períodos de revolução dos planetas são inversamente proporcionais às distâncias medias do sol”.

6. CÁLCULO DA POSIÇÃO NA ÓRBITA ELÍPTICA Até agora descobrimos a função r=r(), de forma que sabemos qual será o tipo de curva descrita pelo ponto m porém, falta conhecermos qual a função que determina as variações de  em função do tempo t, isto é, a função =(t). Vamos desvendar o mistério... Sabemos que para a elipse, em coordenadas polares, a posição do ponto de massa m está definida pela relação r(), porém não conhecemos a dependência dos parâmetros r e  em função do tempo t. Vamos deduzi-las então:

  d A  1 r 2   L A dt 2 2m

a velocidade areolar é:

portanto, após ter descrito um ângulo  na trajetória, em um tempo t, o ponto m terá varrido a área A que será a resultante de integrarmos: A

t

t

t

L L L A    d A   A d t   dt d t  t  2 m 2 m 2 m 0 0 0 0 porém para termos os limites de integração em função de  :

12

Andrés Esteban de la Plaza

13

  d A  1 r 2  1 r 2 d  A dt 2 2 dt dA 1 2 d = r dt 2 dt

ou seja e eliminando d t :

 dA

1 2 r d 2

que agora podemos integrar



A

1 A   d A   r2 d  2 0 0

portanto

Substituindo r 

a( 1  e 2 ) 1  e cos 

 L   t  m

=

Lt 2m

[11]

na integral acima e multiplicando por 2, temos: 

r





2

d

a 2 (1  e 2 ) 2 0 (1  e cos ) 2 d  

=

o

então 

 a (1  e ) 2

2

2

1

 1  e cos  

2

d

0

 e sin  1   a (1  e )  2  2    e  11  e cos    e  1 2

2

2

 1  d  0 1  e cos      

  a 2 ( 1  e 2 )2  e sin  1   d    e 2  1 1  e cos  0 1  e cos   



porém para a elipse  e < 1  e2 < 1 logo a solução da

1

 1  ecos  d  será: 0



1 0 1  e cos  d  

 1 e   arctan tan    2  2  1  e2  1 e 2

Logo

Lt a 2 ( 1  e 2 ) 2  m  e 2  1

  1 e 2  e sin        arctan tan       2 2   1  e cos  2  1 e  1 e   

13

Andrés Esteban de la Plaza

14

 L   t = - a 2 (1  e 2  m

  1 e  2  e sin    )  arctan tan   2   1 e2  1 e  1  e cos  

 

[12]

Esta fórmula é muito bonita porém, o resultado obtido é o inverso que desejávamos, na verdade nós queríamos =(t) e acabamos obtendo t=f(). É fácil observar que despejar  nos leva a uma equação muito complexa, ou seja, devemos encontrar outro método. Introduzimos agora uma quantidade auxiliar E denominada anomalia excêntrica, que cumpre a seguinte condição:

 

 

1 e tan E 2  tan  2 1 e

[12a]

que equivale a dizer que

 1+ e   = 2arctan  tan E 2   1 e 

 

então introduzindo a [12a] na [12], temos:

 e s en   L   t  - a 2 1  e 2     m 1  e cos 

  1  e2  E

[13]

entretanto se pode demonstrar que para a quantidade auxiliar E introduzida:

1  e2 sen E   sen  1  e cos  de forma que

sen E 1  e2



sen  1  e cos 

substituindo esta última na [13] chegamos a:

 e sen  E   L    t  a 2 1  e2     m 1  e2  1  e cos   e sen E a 2 1  e2  E   L 2 2  e sen E  E t   a 1  e        2 2  2  m 1  e 1  e 1  e   Logo

 L   t  a 2 1  e 2  E  sen E  m

[14]

14

Andrés Esteban de la Plaza

15

 Vamos definir agora as constantes para o cálculo prático, como:

L 2ab 2a 2 1  e 2   constante = m T T para o planeta Terra  T=Tt

[15]

e at= 1 ua (unidade astronômica),

 logo para um planeta qualquer com T, a, temos (aplicando [10]) que:

a 3t



a3

3

Tt2

a 3t a  3 a 2 1   T TT TT2



T2

e substituindo a [16] na [15] obtemos:

L 2 1  e 2 a 2 a  m TT

3

[16]

2

ou seja que a [14] fica:

 L   t  a 2 1  e 2 E  e sen E  m 2 1  e 2 a 2 a TT 3

2a TT

2

3

2

t  a 2 1  e 2 E  e sen E

t  E  e sen E [17]

Esta última fórmula é básica, pois é a empregada no cálculo das efemérides.

As seguintes denominações são as usuais:

 1+ e  tan E 2   1 e 

 

 = anomalia verdadeira [radianos] = 2arctan 

E = anomalia excêntrica =

2 arctan

a = distancia média ao sol [ua]

[18]

 

1 e tan  2 [radianos] 1 e

(unidades astronômicas)

15

Andrés Esteban de la Plaza

16 3

n = movimento diurno médio = n =k a

3

2

onde

2a TT

2

[radianos/dia]

k=constante de Gauss =

M = anomalia média = n t

-2  2  radianos  ua 3  [18a]  365.256  dia 

[19]

M = E - e sen (E) [radianos]

[20]

(fórmula de Kepler)

L ogo as fórmulas reduzem-se a: M=nt

e

M = E - e sen (E)

sendo t = tempo após passo pelo periélio, em dias. TT = ano trópico = 365.256374 dias

Cálculo da Anomalia Media M O processo de cálculo se reduz a determinar a anomalia média M [19], e de posse desta procedermos a calcular a anomalia excêntrica E que satisfaz a equação de Kepler [20], logo calculamos  segundo a [18]. Com  calculamos r , o que determina a posição do ponto m na órbita. A resolução da Equação de Kepler é realizada por iterações sucessivas, de forma que o algoritmo de cálculo é o seguinte: t = T-T0 onde T0 e a data da última passagem pelo periélio a = distância média ao sol Nota: Normalmente o tempo da passagem T0 está dado em dias Julianos JD, de forma que a data de calcula também está em dias julianos e o cálculo de t se reduz apenas a encontrar a diferença (em dias julianos) entre T e T0.

16

Andrés Esteban de la Plaza

17

O diagrama de fluxo do algoritmo é:

Vamos fazer um exemplo de aplicação: caso Júpiter a = 5.208174 ua ; e = 0.049284 T = 2450896.510556 dj =24 Março 1998 T0 = 2446966.84378 dj =24 Junho 1987 (periélio) t = 3929.666776 dj n = 0.00144728573606 M = n t = 5.687350672374 radianos ( =325.861190138) E=M= 5.687350672374 17

Andrés Esteban de la Plaza

18

E1=M+e senE= 5.659692503366 E1-E= -0.0276581690088 E=E1= 5.659692503366 E1=M+e senE=5.658575010 E1-E=-0.001117493 E=E1=5.658575010 E1=M+e senE=5.658530316 E1-E=-0.000044694 E=E1=5.688530316 E1=M+e senE=5.658528529 E1-E=-0.000001787 E=E1=5.658528529 2 tan E 2 e como senE = =-0.584818898  0 0 sen  < 0 sen  > 0 sen  > 0 sen  < 0 sen  < 0

 = 0  = 0  = 0  = 180 + 0  = 180 + 0  = 360 + 0

     

0 = 90 0 = 270 0 > 0 0 < 0 0 > 0 0 < 0

0 <  < 90 90 <  < 180 180 <  < 270 270 <  < 360

Também, na hora dos cálculos, devemos observar o seguinte:  O movimento diurno médio n = k a-1.5 pode ser expressado tanto em [ radianos / dia ] como em [  / dia]. Ou seja, a constante de Gauss k, pode ser expressada em: k = 2 / 365.2564

[ radianos / dia ]

k = 360 / 365.2564 [  / dia ]  A anomalia média M = n t poderá então estar expressada em [ radianos ] ou [  ]. Normalmente as tabelas ou efemérides proporcionam M em [  ].  Para o cálculo da Equação de Kepler, é obrigatório o emprego de M em [ radianos ], e o valor de E obtido também está em radianos. Se é desejado se pode converter E para [] e calcular os valores de seno, coseno e tangente; entretanto se E não estiver dentro de uma função trigonométrica, então E deverá ser tomada em radianos!  Se a última passagem pelo periélio ocorreu há mais de um período do planeta, então o cômputo da anomalia média certamente dará um valor superior a 360 ou 2 radianos (6.2831853...). Neste caso é melhor reduzirmos a anomalia média para o valor fracionário de um giro completo. Expressando M em radianos ou graus sexagesimais, esto se reduz à seguinte fórmula:  M   M  M calculo  radianos     IP    2  2    2 

 M   M   M calculo         IP     360     360   360

onde IP representa a função que nos dá a parte inteira do argumento contido. Vejamos um exemplo: Se M = 9.28 radianos = 531.7048339, então Mcálculo{rad} = [ (9.28/2) - IP(9.28/2)]2 = [1.4769579 - 1]2 = = [0.4769579]2 = 2.9968147 rad Mcálculo {} = [ (531.7048339/360) - IP(531.7048339/360)]360 = [ 1.4769579 - 1]360= = [0.4769579] 360 = 171.7048339  logo é evidente que 2.9968147 radianos = 171.7048339

19

Andrés Esteban de la Plaza

20

7. CONSIDERAÇÕES SOBRE A CONSTANTE DE GAUSS, O MOVIMENTO DIURNO MÉDIO E A FÓRMULA DE KEPLER  CONSTANTE DE GAUSS E MOVIMENTO DIURNO MÉDIO L Dada a equação [14] :   t  a 2 1  e2  E  sen E , e como:  m

M = E - sen(E) =

2 t T

Aelipse = ab =  a 2 1  e 2  L (sendo At a área varrida até o instante t )   t = 2At  m resulta que a equação [14] eqüivale a:

2A t 

A elipse 

M

[20a]

que encontramos foi a expressão do dobro da área orbital varrida em um tempo t desde o periélio. Por outro lado, se na eq. [20a] isolamos M obtemos a eq. [20aa], logo podemos deduzir a eq. [20aaa]:

 At   M  2    A elipse 

[20aa]



 t A t  A elipse     T

[20aaa]

 At   não é outra coisa Prestando atenção na equação [20aa] vemos que o quociente   A elipse  senão a expressão em forma de % da área varrida respeito da área total; porém, M é um ângulo dependente dp tempo t, de forma que o produto 2.% dá um ângulo proporcional à área varrida, ou seja ao tempo t. Isto faz pensar num movimento circular. Na eq. [20aaa] verificamos a proporcionalidade entre o tempo t do movimento circular e a área varrida At (pois a velocidade angular é constante na circunferência), ou seja, se falamos de um setor circular varrido em um tempo t igual a 25% do período T obtemos 25% da área do círculo, então o ângulo do setor será 25%.2 = /4 radianos ou 25%.360 = 90. Estas considerações levam a pensar que talvez M esteja relacionado com algum tipo de circunferência auxiliar; ou seja, a órbita e circular e a velocidade angular é constante. Quando tratarmos da Equação de Kepler veremos que esta idéia é absolutamente válida.

Analisemos agora as dimensões da constante L/m. O momento da quantidade de movimento m [L] tem a dimensão { m.kg }, a massa [m] ={ kg }, de forma que a dimensão do quociente s 2 {L/m} é {m /s}, que multiplicado pelo tempo t nos dá uma superfície. Isto resulta bastante claro se observamos o membro direito da equação [14], onde o único fator não adimensional é o semi-eixo maior da órbita a elevado ao quadrado. Ë claro que pouco interessa se a está em metros ou unidades astronômicas, a superfície será sempre {m2} ou {ua2}.

20

Andrés Esteban de la Plaza

21

Examinando agora a equação de Kepler, seja na sua forma dada pela equação [17], segundo 2 a  3 2 t  E  e sen E , ou na sua forma dada pela equação [20], a qual: TT segundo a qual: M = E - sen(E) , vemos que não existe fator dimensional, exceto o valor 2, que está expressado em radianos, ou como já vimos antes, em graus sexagesimais (não esquecer que para operar no lado direito da expressão devemos trabalhar em radianos!). É interessante entender o que significa a constante de Gauss:

 L  t  a 2 1  e2  E  sen E , com a [20]: n t = E - sen(E),  m

Vamos igualar a eq. [14]: 

para isto tomamos a eq. [14] e isolamos E - sen(E), de modo que agora temos:  L  E - sen(E) =  2 m 2   t  a 1 e 

comparando esta última com a eq. [20] chegamos a:

  L     m   a 2 1  e2  

    t   

=

  L     m   a 2 1  e2  

     

n =

n  t

k a

=

3

2

 t

ou seja:

=

2

porém da eq. [8 a]: 1  e f  L

L2

m

2

 GMS (1  e)f

então:

   

a

GM sm2

k a

2

[21]

deduzimos que :

L  GMS  (1+ e)f m



3

e como

(1 + e)f  a (1  e2 )

L  GMS  a (1 + e 2 ) , de forma que introduzindo este valor na equação [21] temos: m

GM S   2  1 e 

1  e2  a

2

a

3

2



=

n =

GM S =

k a

n =

GMS = k

3

2

k a

3

2

[22]

[23]

21

Andrés Esteban de la Plaza

22

3 4  2a 3 2 a por outro lado, da equação [10]: T  deduzimos que GM S  4 2 de forma que temos: GM S T 2

3

a 2 , portanto a equação [22] GM S  2  T forma:

ou seja:

3

n

=

a

n



a

2

3

2

do movimento diurno médio fica da seguinte



GM S = a



GM S

3

3

2



a 2  2 T 2 T

= =

k a

3

2

2  3 2 a TT

[24]

3

do qual:

a 2 k  GM S  2 T

[25]

portanto, uma maneira de encontrar k e utilizar os valores de a e T terrestres, de forma que se a = 1 ua, e T = TT, (atenção que para qualquer outro planeta, k se calcula pela [25]),então: k  GMS 

2 TT

[26]

tal como foi feito na [18a] quando introduzimos a constante de Gauss. Calculemos então o valor de k segundo a eq. [23], para isto devemos passar as unidades de [G] m3 m3 m3 ua 3 -11 ={ } para { }, ou seja que se G [ ] = 6.67 10 então: kg dia 2 kg s2 kg s2 kg s2 ua 3 m3 -11 G[ ] = 6.67 10 (1ua/1.495979 1011 m)3 (86400 s/dia)2 = kg dia 2 kg s2 ua 3 ua 3 -34 G[ ] = 1.4872 10 kg dia 2 kg dia 2

e como MS = 1.991 1030 kg então:



k=

GMS = 0.0172076  0.0172

o que concorda com o valor introduzido pela equação [26] pois: 2 = 2 / 365.2564 dia = 0.017202 k TT Resumindo:  O movimento diurno médio (n) não é senão uma expressão da velocidade angular média do planeta (  ), que resulta de imaginar a sua órbita circular, de forma que   2 T . A anomalia média (M), dá a idéia do ângulo do setor circular varrido até o instante t, se a órbita fosse circular.  A constante de Gauss (k) é um valor constante para o sistema solar, portanto de igual valor para todos os planetas. Digamos que k quantifica a “grandeza gravitacional” do sol, 22

Andrés Esteban de la Plaza

23

pois k = GMS , ou seja sua capacidade de influir gravitacionalmente sobre outras massas. Como dato interessante, esta “grandeza gravitacional”, multiplicada pela distância média do planeta ao sol (elevada a 1.5), nos dá o movimento diurno médio, ou seja como o planeta se comporta a essa distância do sol.  A EQUAÇÃO DE KEPLER Quando resolvemos a equação [12], introduzimos o parâmetro E tal que se cumpria a relação 1 e dada pela equação [12a] : tan E 2  tan  2 . Resta perguntar agora o que a eq. [12a] 1 e significa.

 

 

Vamos dar uma olhada na Figura 6:

Figura 6

Nesta figura observamos o sol, dado pelo ponto S, a posição do planeta, dado pelo ponto P, na sua órbita elíptica de periélio SQ¸excentricidade e, e dimensões orbitais a e b. Agora traçamos uma órbita circular de raio a e centro O, que obviamente não coincide com S, nesta órbita circular temos um planeta imaginário, de idênticas caraterísticas a P, denominado P’. Este planeta P’ tem o mesmo período orbital T que P porém, movimenta-se com velocidade angular uniforme (já que a órbita é circular). Resulta evidente que a velocidade angular do planeta P’ vale 2/T, mas isto não é outra coisa senão o movimento diurno médio do planeta P.

A posição do planeta P em coordenadas polares de centro S está definida pela anomalia verdadeira (), e o raio vetor (r). A posição do planeta P em coordenadas polares de centro O (excêntrico respeito de S) está definida pela anomalia excêntrica (E), e o raio vetor (a). Vamos demonstrar agora a correspondência entre as posições de P e P’ tal como mostradas na Figura 6. Analisando a figura observamos que a componente horizontal do raio vetor a de P’, vale acosE. Este valor é igual à distância c + x, mas c = ae, e x = rcos. Lembrando que r a1  e2  vale (elipse com o centro de coordenadas polares centrado no foco direito), então 1  e cos a1  e2  cos  podemos dizer que: . Após operarmos nesta igualdade a cos E = ae + 1  e cos  podemos encontrar o valor de cosE, ou seja:

cos E =

e  cos  1  e cos

[27]

23

Andrés Esteban de la Plaza

24

 

1  cos E Lembrando das funções trigonométricas que tan E 2  , podemos introduzir a 1  cos E equação [27] nesta última, de forma que:

 

1  cos E 1  e cos   e  cos  tan E 2    1  cos E 1  e cos   e  cos  mas

 

1  cos   tan  2 assim: 1  cos 

1  e1  cos   1  e 1  cos   1  e1  cos   1  e 1  cos 

 

tan E 2



 

1 e tan  2 1 e

=

que é precisamente a equação [12a], tal como queríamos demonstrar. Agora sabemos o real significado do conceito anomalia excêntrica.

8. CÁLCULO DE ÓRBITAS NÃO ELÍPTICAS: ÓRBITAS PARABÓLICAS E HIPERBÓLICAS  CÁLCULO DAS ÓRBITAS PARABÓLICAS Quando chegamos à equação [8aa] que nos dava a fórmula geral da cônica, continuamos nosso raciocínio assumindo que 0 < e < 1, logo, para chegar até a equação de Kepler, introduzimos o valor de r no cálculo da integral dada pela equação [11]. Para o caso das órbitas parabólicas, sabemos que e = 1, porém vamos refazer a integração da [11] considerando a equação da parábola (que resulta de introduzir e=1 na eq. [8aa]). Vejamos a Figura 7: a distância focal SQ do sol à passagem pelo periélio, foi chamada de f; a excentricidade da elipse vale e = 1, assim colocando estes dados na equação [8aa]: 2f 1  ef r , chegamos a equação da órbita parabólica: r  . Agora, apenas por 1  cos  1  e cos  uma questão de nomenclatura, denominamos à distância focal f de q, ou seja q = f. Portanto a equação da parábola é:

r

Agora vamos fazer a integração:

2q 1  cos 

 L 2A =    t   m

[28]



r

2

d  , introduzindo o valor de r dado pela

o

eq. [28].

 L   t =  m



r o



2

d =

4q 2

 1  cos  0

 

2

 

 

1 1  d  = 4q 2   tan  2  tan 3  2  6 2 

 

1   = 2q 2 tan  2  tan 3  2  3   porém como

24

Andrés Esteban de la Plaza

25

 L     m

GMS 

 L     m

GMS  2  q

1  ef e na parabola e  1 e f  q resulta que:

então:

1

2

   

   

1  L      t = 2q 2 tan  2  tan 3  2   m 3   1 1   GM S  2  q 2  t  2q 2 tan  2  tan 3  2  3  

isolando agora o termo das tangentes:

 

 

1   3   tan 2  3 tan 2  

1

2 q 2 GM S   t  2 q2

GM S 

q

3

2

2

t

ou seja que:

 

 

1 tan  2  tan 3  2 3



GM S 

q

3

2

2

t

[29]

equação análoga a [17] no sentido que na sua resolução obtemos o valor de . Neste caso temos que resolver uma equação cúbica em tan(/2). Pode-se demonstrar que esta equação tem uma raiz real (a que nos interessa) e duas imaginárias. Existem bastantes métodos de resolução. Propomos o Método de Newton, que diz: Se x0 é um valor aproximado da raiz da função f(x)=0 então como aproximação mais exata se toma x1  x 0 

f (x0 ) f  (x0 )

Substituindo agora x1 por x0 obtemos uma melhor aproximação x2.

Logo, para nosso problema (resolução da equação [29]):  Para facilitar os cálculos: multiplicamos a eq. [29] por três e chamamos a tan(/2) = E, ou seja a equação [29] passa agora a ser: E3 + 3E - M = 0  definimos M = n .t,

onde n = 3  GM S 

q

3

2

2

 3 k 

q

3

2

2

k  constante de Gauss

 aplicamos o método de Newton, de forma que denominando E0 ao valor aproximado da raiz, o próximo valor será: E1  E 0 

logo:

E1  E0 

f (E 0 ) f (E 0 )

e como f ( E )  0  f ( E )  E 3  3E  M  0

E0 3  3E0  M 3E0 2  1



3E0 E0 2  1  E0 3  3E0  M 3E0 2  1



f ( E )  3E 2  3  3E 2  1

3E30  3E0  E20  3E0  M 3E0 2  1



2E30  M

3E02  1

25

Andrés Esteban de la Plaza

26

E i 1 =

ou seja

2E 3i  M

[29a]

3E 2i  1

Este valor será iterado até que f( Ei )  0 , por exemplo até que  f( Ei )  ≤ 10-8.  A primeira aproximação, isto é E0 , vale E0 = 0  iteramos então de acordo ao seguinte algoritmo em QBASIC da Microsoft: 10 20 30 40 50 60 70

M = nt E=0 E1 = 2E3 + M / [3(E2 + 1)] IF E1 3 + 3E1 - M > 10-8 THEN E = E1 : GOTO 30 END IF  = 2 . Atan ( E )

 CÁLCULO DAS ÓRBITAS HIPERBÓLICAS

ae2  1 Neste caso se procede como no anterior, lembrando que se e > 1 então r  pois 1  e cos segundo a definição na hipérbole, c > a  f = c - a  (1+e)f = (1+e)(e-1)a = a(e2-1), logo:  L   t =  m



r d 2



=

o

 L    t = a 2 e 2  m

 0

a 2 ( e 2  1 )2 d  ( 1  e cos  )2

e sen  1 2   1  2  2  e  1    e  11  e cos  



 0

 d 1  e cos     

então:   e sen    L  2 2     t = a  e  1   m 1  e cos       

   

  e  1 tan  2   ln  e2  1    e  1 tan 2  1

 e2  1     [30] e 2  1   

Vamos agora fazer algumas considerações a respeito das funções trigonométricas hiperbólicas, para isto tomamos o fator dentro do logaritmo natural na equação [30]:

  e  1 tan 2  

e  1 tan  2  e 2  1 e2  1



 2   1  2   1

e 1 tan e 1 e 1 tan e 1

e se

 

e 1 tan  2 > 1 e+1

então temos que:

   ln   

 2  1   2 Arctanh   2  1

e 1 tan e 1 e 1 tan e 1

 

 e 1 tan  2  e 1 

  

= 2 Arctanh tanh E 2  E 26

Andrés Esteban de la Plaza

27

logo chamando a

  ln    

resulta que:

e e e e

 

 

e -1 tan  2  tanh E 2 e +1

   

 1 tan  2  1 1  E 1  tan  2  1 1 

[31]

podemos agora calcular o valor de senh(E), ou seja:

   

   

 

2   2 tanh E 2 e  1 sen 2  e  1 cos 2 e2  1 sen  senh E  2    1  e cos 1  e cos  e  1 cos  1  tanh 2 E 2 2

logo se deduz que



sin   1  e cos 

1 e2  1

senh E

[32]

assim, introduzindo as equações [31] e [32] na equação [30] obtemos:    sen    L   2 2    t = a  e  1 e      m  1  e cos       

  senh E    L      t = a 2  e 2  1 e    2  m     e  1

   

  e  1 tan  2   ln  e2  1    e  1 tan 2  1

 2 E   a 2 e  1   1

 L    t  a 2 e2  1  e senh E  E  m

e

2

 1

e2  1

 e2  1     e 2  1   

e senh E  E

[33]

lembrando agora que na equação geral das cônicas [8], o numerador representa o parâmetro p que vale (1+e)f, então, de acordo com a eq. [8a] obtemos, como principio geral para as órbitas cônicas: L  GM s  1  e f m L  GM s  1  e q m

e tinhamos feito f  q entao

porém na definição inicial de hipérbole tínhamos visto que 1  e f  1  e q  ae 2  1 , ou seja que:

27

Andrés Esteban de la Plaza

28 L  GM s  1  e q  GM s  a e2  1 m L  GM s   e2  1a [34] m

de forma que introduzindo a eq. [34] na equação [33]:  L    t  a2  m GM s

3

2

[33]

e 2  1e senh E  E

e2  1  t  a 2

a

GM s  a

e 2  1  e senh E  E

 e senh E  E

t

[35]

e como a = q / (e-1) para a hipérbole, então a eq. [35]: GM s  a

GM s 

3

q

3

2

 e  1

 e senhE  E

t

2

3

2

t

= e senh E  E

[35]

[36]

Equação parecida a obtida no cálculo das órbitas elípticas, neste caso a resolução da eq. [36] é idêntico ao da equação [17] ou [20], ou seja: E = e . senh(E) - M M = n.t 3

n= k

q

2

 e  1

3

2

onde k  GMS

De livre publicação desde que o autor seja mencionado, direitos reservados. Última revisão 29 junho 2016 – Rio de Janeiro, Brasil. [email protected]

Homines dum docent discunt. - Séneca o Moço (Córdoba 4 a.C.- Roma 65 d.C.)

28