INPE-12970-TDI/1018
ESTUDO COMPARATIVO DE TÉCNICAS DE CONTROLE DE
ATITUDE EM TRÊS EIXOS PARA SATÉLITES ARTIFICIAIS
Gilberto Arantes Júnior
Dissertação de Mestrado do Curso de Pós-Graduação em Engenharia e Tecnologia
Espaciais/Mecânica Espacial e Controle, orientada pelo Dr. Ijar Milagre da Fonseca,
aprovada em 23 de fevereiro de 2005.
INPE
São José dos Campos
2005
629.7.062.2
ARANTES JR, G.
Estudo comparativo de técnicas de controle de atitude
em três eixos para satélites artificiais / G. Arantes Jr. – São
José dos Campos: INPE, 2005.
201p. – (INPE-12970-TDI/1018).
1.Estabilização de veículos espaciais. 2.Controle de
atitude. 3.Estabilização em três eixos. 4.rodas de reação.
5.Controle magnético. I.Título.
Aprovado (a) pela Banca Examinadora em
cumprimento ao requisito exigido para
obtenção do Título de Mestrado em
Engenharia e Tecnologia Espaciais/Mecânica
espacial e Controle
Aluno (a): Gilberto Arantes Junior
São José dos Campos, 23 de fevereiro de 2005
”No meio de qualquer dificuldade encontra-se a oportunidade”
ALBERT EINSTEIN
“A maior recompensa do trabalho do Ser humano não é o que ele(a) ganha com
isso, mas sim o que ele(a) se torna com isso.”
DESCONHECIDO
A meus pais, GILBERTO MOURA ARANTES e
CÉLIA DIAS ARANTES,
pelo apoio e paciência.
E ao tio e amigo
MIGUEL BERNARDES DE CASTRO (in memorian).
AGRADECIMENTOS
A meus pais, pelo amor, compreensão, paciência e por incentivar e acreditar na
importância de ir em busca dos sonhos.
À minha famı́lia pela certeza de sempre poder contar com seu eterno apoio e incentivo, em especial aos meus Avós, José, Terezinha, Jerônimo e Cacilda.
À querida tia Maria José e ao inesquecı́vel tio Miguel (in memorian) por todos os
conselhos, incentivos e valiosa torcida.
Aos tios Aluisio pela amizade, generosidade e todos os “ensinamentos filosóficos”, e
Branca pelo amor de mãe que me foi dedicado, serei eternamente grato.
Ao orientador e amigo Prof. Dr. Ijar Milagre da Fonseca pela orientação, diretrizes,
conselhos e à valiosa confiança e por acreditar que eu pudesse realizar este trabalho.
Ao apoio financeiro dos meus pais e amigos.
Ao grande amigo e colega de curso Rolf Vargas por poder dividir minhas dificuldades
e pelo seu bom humor, nos momentos difı́ceis.
A todos os colegas do curso, pela amizade e companheirismo demonstrados, em
especial a Leandro e Cecı́lia, por toda a colaboração.
Ao Laboratório LABSIM e todos os técnicos, pela oportunidade de estudos e utilização de suas instalações.
Ao Instituto Nacional de Pesquisas Espaciais (INPE), pela oportunidade.
Aos professores do INPE em especial a André Fenili, Marcelo Lopes, Waldemar de
Castro, Mario Ricci e Evandro Rocco, pelo conhecimento compartilhado.
Aos professores Hélio Koiti Kuga e Roberto Vieira da Fonseca Lopes pelas, sugestões
na elaboração desse trabalho.
A meus amigos da República, mestre Clério, à simpática mexicana Nora Trelles e
ao disciplinado Marcos Timbo por toda a paciência e stress compartilhado.
À amiga Rose e ao amigo Fred, pelo incentivo.
À Coordenação de Aperfeiçoamento de Pessoal de Nı́vel Superior (CAPES), por um
ano de bolsa concedida.
A todas aquelas pessoas que de uma maneira ou de outra contribuı́ram, e por omissão
não constam nessa lista, peço desculpas e agradeço.
RESUMO
Neste trabalho propõe-se um estudo sobre técnicas de controle de atitude para
satélites estabilizados em três eixos e, através da modelagem e simulação
computacional, analisar, comparar e fazer um estudo de alternativas/viabilidade
de sistemas de controle de atitude (ACS) em três eixos, com base nos requisitos de
missões espaciais. Para o estudo de alternativas/viabilidade dos sistemas de
controle de atitude foi realizado um estudo comparativo de diferentes técnicas de
controle, utilizadas para estabilização de satélites em três eixos, utilizando-se
diferentes atuadores, tais como: 1) rodas (de reação e volantes de inércia); 2)
bobinas magnéticas. Os procedimentos de estabilização estudados foram: 1)
controle de atitude em três eixos utilizando rodas (de reação e volantes de inércia)
e bobinas magnéticas; 2) controle de atitude em três eixos utilizando apenas
bobinas magnéticas. Foram utilizados a teoria do Regulador Linear Quadrático
(LQR) e Regulador Quadrático Gaussiano (LQG) e controladores não lineares
baseados em energia (energy based control ) para o desemvolvimento do projeto de
controle no modo de estabilização. A teoria do LQR rastreio (tracking) e o
controlador Proporcional Derivativo (PD) foram usados no modo de aquisição de
atitude. Para a fase de redução da velocidade angular (detumbling) foi utilizado o
controlador de Wisniewski ou Bdot. As diferentes configurações dos ACS são
discutidas, comparadas e analisadas, visando avaliar o desempenho dos
procedimentos de controle aqui desenvolvidos para os modos de detumble,
estabilização e aquisição da atitude. A formulação obtida nesse trabalho foi
aplicada no controle de atitude do satélite brasileiro EQUARS (Equatorial
Atmosphere Research Satellite), que motivou este trabalho. Os resultados obtidos
atendem as especificações e os requisitos do satélite.
COMPARATIVE STUDY OF THREE AXES ATTITUDE CONTROL
TECHNIQUES FOR ARTIFICIAL SATELLITES
ABSTRACT
This work deals with the study of attitude control techniques for three axes
stabilized satellite, and through modeling and computational simulation, analyze,
compare and develop a feasibility study for attitude control system based on space
missions requirements. The attitude control system feasibility study is carried out
by using different control techniques for satellite three-axis stabilization, and
different actuators such as 1) reaction wheels and momentum wheels; 2) torque
coils. The procedures used for stabilization were: 1) Three axis attitude control by
using reaction wheels and momentum wheels combined with torque coils; 2) Three
axis attitude control by using torque coils only. Each of these technique was
implemented in computer for the attitude control simulations. The control law was
based on the Linear Regulator Quadratic (LQR) and Linear Quadratic Gaussian
(LQG) for linear systems and on energy based control for no linear systems. These
techniques were used for the stabilization mode. The LQR tracking and the
proportional derivative (PD) techniques were used for the acquisition mode. The
Wisniewski or Bdot approach has been used for detumblig phase. The different
configuration results for the control modes are analyzed and discussed in terms of
the performance of the control procedures associated with the attitude control
modes (detumble, stabilization and acquisition). The control formulation has been
applied for the brazilian satellite EQUARS (Equatorial Atmosphere Research
Satellite). The results comply with the satellite specification and requirements.
SUMÁRIO
LISTA DE FIGURAS
LISTA DE TABELAS
LISTA DE SÍMBOLOS
LISTA DE SIGLAS E ABREVIATURAS
CAPÍTULO 1– INTRODUÇÃO
31
CAPÍTULO 2– OBJETIVO
37
2.1
Meios e Recursos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
2.2
Metodologia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
CAPÍTULO 3– REVISÃO BIBLIOGRÁFICA
41
3.1
Atitude . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.2
Controle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
3.3
Literatura do INPE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
CAPÍTULO 4– DEFINIÇÕES DE NOTAÇÕES
4.1
51
Sistemas de Referência . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.1.1
Referencial Inercial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.1.2
Referencial Orbital . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.1.3
Referencial do Satélite . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.2
Representação da Atitude . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4.2.1
Matriz de Rotação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
4.2.2
Parâmetros de Euler . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
4.2.3
Ângulos de Euler . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
4.2.4
Ângulo Eixo Equivalente . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.2.5
Rotações Infinitesimais . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
4.3
Sumário de Notações . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
CAPÍTULO 5– FORMULAÇÃO DO PROBLEMA
5.1
5.1.1
63
Modelagem Matemática: Cinemática e Dinâmica . . . . . . . . . . . . . . 63
Equações da Cinemática . . . . . . . . . . . . . . . . . . . . . . . . . . 64
5.1.2
Equações da Dinâmica . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
5.1.3
Torque devido ao Gradiente de Gravidade . . . . . . . . . . . . . . . . 71
5.2
Campo Magnético Terrestre . . . . . . . . . . . . . . . . . . . . . . . . . 71
5.3
Tratamento das Equações da Dinâmica para os Casos Estudados nesse
Trabalho . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
5.4
Linearização . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
5.4.1
Linearização do Modelo do Satélite Equipado com Rodas . . . . . . . . 76
5.4.2
Linearização do Modelo do Satélite Equipado com Bobinas . . . . . . . 78
5.5
5.5.1
5.6
Considerações Sobre os Atuadores, os Modelos Matemáticos e Sistemas
de Referênica Adotados . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
Modelo dos Atuadores . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
Torques Ambientais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
5.6.1
Torque Devido ao Gradiente de Gravidade . . . . . . . . . . . . . . . . 85
5.6.2
Torque Aerodinâmico . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
5.6.3
Torque de Pressão de Radiação Solar . . . . . . . . . . . . . . . . . . . 88
5.6.4
Torque Devido ao Dipolo Residual . . . . . . . . . . . . . . . . . . . . 88
5.7
Perturbações internas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89
CAPÍTULO 6– PROJETO DE CONTROLE
6.1
6.1.1
6.2
91
Modo de Detumble . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91
Controlador Bdot . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91
Modo de Estabilização . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
6.2.1
Método LQR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
6.2.2
Método LQG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97
6.2.3
Controladores Baseados em Energia . . . . . . . . . . . . . . . . . . . . 101
6.3
Modo de Aquisição . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102
6.3.1
LQR Tracking
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
6.3.2
Proporcional Derivativo . . . . . . . . . . . . . . . . . . . . . . . . . . 104
CAPÍTULO 7– SIMULAÇÕES
109
7.1
Modo de Detumble . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109
7.2
Modo de Estabilização . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113
7.2.1
Método LQR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114
7.2.2
Controladores Baseados em Energia . . . . . . . . . . . . . . . . . . . . 124
7.2.3
Método LQG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130
7.3
Modo de Aquisição . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133
7.3.1
LQR Tracking
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
7.3.2
Controle PD . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
CAPÍTULO 8– CONCLUSÃO
8.1
141
Sugestões para Tabalhos Futuros . . . . . . . . . . . . . . . . . . . . . . 144
REFERÊNCIAS BIBLIOGRÁFICAS
147
APÊNDICE A–CÁLCULO DO GRADIENTE DE GRAVIDADE
157
APÊNDICE B–PROPAGAÇÃO DA ÓRBITA E TRANSFORMAÇÕES
163
B.1 Cálculo da Órbita . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163
B.1.1 Posicionamento de Satélites - Problema Direto . . . . . . . . . . . . . . 163
B.1.2 Equação de Kepler . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 164
B.1.3 Matriz de Rotação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166
B.2 Matriz de Rotação do Sistema Inercial - Sistema do Sátelite . . . . . . . 169
OF
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 170
B.2.1 Matriz RPECI
B.2.2 Matriz ROF
P OF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172
B.2.3 Matriz RBF
OF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173
B.3 Cálculo da Latitude, Longitude e Altura . . . . . . . . . . . . . . . . . . 175
APÊNDICE C–IMPLEMENTAÇÃO EM SIMULINK
183
APÊNDICE D–PROGRAMAS EM MATLAB
193
D.1 Projeto LQR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193
D.2 Projeto LQG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193
APÊNDICE E– TOOLBOX ATITUDE
199
E.1 Exemplos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 200
E.1.1 Equações do Movimento . . . . . . . . . . . . . . . . . . . . . . . . . . 200
E.1.2 Ambiente . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 200
LISTA DE FIGURAS
4.1
Sistemas de referência, inercial (ECI), orbital (OF) e do satélite (BF)
4.2
Construção dos ângulos de Euler . . . . . . . . . . . . . . . . . . . . 57
5.1
Seqüência de rotações 3(ψ) − 2(θ) − 1(φ)
5.2
Ilustração do satélite com rodas de reação, bobinas e os referenciais
orbital OF (xo , yo , zo ) e do corpo BF (x, y, z) . . . . . . . . . . . . . . 66
5.3
Configuração das bobinas magnéticas no satélite . . . . . . . . . . . . 70
5.4
Campo magnético local Bo usando o modelo IGRF 2000 . . . . . . . 72
6.1
Configuração do controle LQR . . . . . . . . . . . . . . . . . . . . . . 94
6.2
Sistema planta mais controlador . . . . . . . . . . . . . . . . . . . . . 97
6.3
Estrutura básica do sistema de controle LQG
6.4
Estrutura do filtro de Kalman-Bucy . . . . . . . . . . . . . . . . . . . 100
6.5
Configuração do controle LQR tracking . . . . . . . . . . . . . . . . . 104
6.6
Configuração do controlador PD . . . . . . . . . . . . . . . . . . . . . 108
7.1
Velocidade angular do satélite ω bib . . . . . . . . . . . . . . . . . . . . 111
7.2
velocidade angular ω bib . . . . . . . . . . . . . . . . . . . . . . . . . . 112
7.3
Dipolo magnético, torque magnético e potência . . . . . . . . . . . . 112
7.4
Dipolo magnético, torque magnético e potência . . . . . . . . . . . . 113
7.5
Ângulos de atitude roll, pitch e yaw em função do tempo para o caso 1116
7.6
Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas para o caso 1 . . . . . . . . . . . . . . . . . . . . . . . . 116
7.7
Velocidades angulares ω bib e ω bob e rotações por minuto das rodas de
reação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
7.8
Ângulos de atitude roll, pitch e yaw em função do tempo para o caso 2118
7.9
Torque τ bw , torque de acoplamento e quantidade de movimento angular das rodas para o caso 2 . . . . . . . . . . . . . . . . . . . . . . . . 119
52
. . . . . . . . . . . . . . . 65
. . . . . . . . . . . . . 99
7.10 Velocidades angulares ω bib e ω bob e rotações por minuto das rodas para
o caso 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119
7.11 Ângulos de atitude roll, pitch e yaw e tempo para o caso 3 . . . . . . 120
7.12 Velocidades angulares ω bib e ω bob e rotações por minuto das rodas para
o caso 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121
7.13 Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas para o caso 3 . . . . . . . . . . . . . . . . . . . . . . . . 121
7.14 Ganho do controlador em função das órbitas Kc (1, 1) . . . . . . . . . 122
7.15 Ângulos de atitude roll, pitch e yaw em função do número de órbitas,
utilizando bobinas . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
7.16 Dipolo magnético, torque de controle (τ bm ) e potência das bobinas em
função do número de órbitas, utilizando bobinas . . . . . . . . . . . . 124
7.17 Ângulos de atitude em função do número de órbitas . . . . . . . . . . 125
7.18 Dipolo magnético, torque de controle e potência das bobinas em
função do número de órbitas . . . . . . . . . . . . . . . . . . . . . . . 126
7.19 Velocidades angulares ω bib e ω bob e campo magnético terrestre local Bb 127
7.20 Ângulos de atitude . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128
7.21 Dipolo magnético, torque de controle e potência das bobinas em
função do número de órbitas . . . . . . . . . . . . . . . . . . . . . . . 129
7.22 Velocidades angulares ω bib e ω bob e campo magnético local Bb em
função do número de órbitas . . . . . . . . . . . . . . . . . . . . . . . 129
7.23 Ângulos de atitude estimados em função do tempo . . . . . . . . . . . 131
7.24 Ângulos de atitude simulados . . . . . . . . . . . . . . . . . . . . . . 131
7.25 Ângulos de atitude simulados e estimados . . . . . . . . . . . . . . . 132
7.26 Erro dos ângulos de atitude e variação da atitude . . . . . . . . . . . 132
7.27 Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133
7.28 Ângulos de atitude roll, pitch e yaw em função do tempo para a
aquisição de atitude . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135
7.29 Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas em função do tempo . . . . . . . . . . . . . . . . . . . . 136
7.30 Rotações por minuto das rodas de reação e velocidades angulares ω bib
e ω bob . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 137
7.31 Ângulos de Euler em função do tempo . . . . . . . . . . . . . . . . . 138
7.32 Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas em função do tempo . . . . . . . . . . . . . . . . . . . . 139
7.33 Rotações por minuto das rodas e velocidades angulares ω bib e ω bob em
função do tempo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140
B.1 Referenciais (inercial, da órbita) e os elementos Keplerianos (i, ω, Ω) . 164
B.2 Referencial inercial (ECI) e pseudo orbital (POF) . . . . . . . . . . . 170
B.3 Orientação do referencial orbital (OF) e do referencial pseudo orbital
(POF) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 172
B.4 Referencial inercial e cartesiano geocêntrico . . . . . . . . . . . . . . . 179
B.5 Longitude (λ), latitude (φ) e altura (h) do satélite . . . . . . . . . . . 179
C.1 Implemetação em SIMULINK do Modo de detumble . . . . . . . . . 183
C.2 Lei de controle Bdot . . . . . . . . . . . . . . . . . . . . . . . . . . . 184
C.3 Implementação em SIMULINK do modo de estabilização: satélite
equipado com bobinas . . . . . . . . . . . . . . . . . . . . . . . . . . 184
C.4 Controle LQR e controladores baseados em energia . . . . . . . . . . 185
C.5 Modo de estabilização: satélite equipado com rodas, usando o método
LQR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 185
C.6 Modo de estabilização: satélite equipado com rodas, usando o método
LQG . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186
C.7 Modo de aquisição usando o controlador PD . . . . . . . . . . . . . . 186
C.8 Lei de controle LQR tracking . . . . . . . . . . . . . . . . . . . . . . 187
C.9 Lei de controle PD . . . . . . . . . . . . . . . . . . . . . . . . . . . . 187
C.10 Modelo dinâmico do satélite equipado com bobinas . . . . . . . . . . 187
C.11 Modelo dinâmico do satélite equipado com rodas . . . . . . . . . . . . 188
C.12 Modelo das rodas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 188
C.13 Modelo das bobinas . . . . . . . . . . . . . . . . . . . . . . . . . . . . 189
C.14 Modelo do torque de gradiente de gravidade . . . . . . . . . . . . . . 189
C.15 Modelo do campo magnético . . . . . . . . . . . . . . . . . . . . . . . 190
C.16 Propagação da órbita . . . . . . . . . . . . . . . . . . . . . . . . . . . 190
C.17 Cálculo da latitude, longitude e altura . . . . . . . . . . . . . . . . . 191
C.18 Matriz de transformação entre os sistemas de referência ECI e BF . . 191
E.1 Bloco Gyrostat . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 200
E.2 Bloco gradiente de gravidade . . . . . . . . . . . . . . . . . . . . . . . 201
LISTA DE TABELAS
5.1
Elementos orbitais do satélite EQUARS . . . . . . . . . . . . . . . . 73
7.1
Parâmetros de simulação . . . . . . . . . . . . . . . . . . . . . . . . . 110
7.2
Condições iniciais da simulação para o modo de estabilização utilizando a metodologia LQR . . . . . . . . . . . . . . . . . . . . . . . . 114
7.3
Parâmetros de simulação . . . . . . . . . . . . . . . . . . . . . . . . . 130
7.4
Condições iniciais da simulação para o modo de aquisição utilizando
rodas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
8.1
Alternativas para o ACS . . . . . . . . . . . . . . . . . . . . . . . . . 141
B.1 Programa MATLAB para o cálculo da órbita circular . . . . . . . . . 169
BF
B.2 Programa MATLAB para o cálculo de c1i , (i = 1, 2, 3) da matriz RECI
176
BF
B.3 Programa MATLAB para o cálculo de c2i , (i = 1, 2, 3) da matriz RECI
177
BF
B.4 Programa MATLAB para o cálculo de c3i , (i = 1, 2, 3) da matriz RECI
178
B.5 Programa MATLAB para o cálculo da longitude, latitude do ponto
sub-satélite e altura do satélite . . . . . . . . . . . . . . . . . . . . . . 181
D.1 Programa MATLAB para o projeto LQR do modelo satélite equipado
com rodas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 194
D.2 Programa MATLAB para o projeto LQR do modelo satélite equipado
com rodas (cont.) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 195
D.3 Programa MATLAB para o projeto LQR do modelo satélite equipado
com bobinas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 196
D.4 Programa MATLAB para o projeto LQG do modelo satélite equipado
com rodas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 197
D.5 Programa MATLAB para o projeto LQG do modelo satélite equipado
com rodas (cont.) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 198
LISTA DE SÍMBOLOS
Bb
Bo
In
J
Ji
Jo
Js
Jw
Kd
Kp
L
Lb
Qc
Qf
Rab
Rc
Rf
Sa
T
A
B
C
G(s)
Ho
Jx , Jy , Jz
K(s)
Kc
Kf
Pc
Pf
To
V
X, Y, Z
x̂
cb1
cb2
cb3
fw
hs
Vetor campo geomagnético expresso no sistema BF
Vetor campo geomagnético expresso no sistema OF
Matriz identidade de ordem n
Tensor de inercial em relação ao sistema BF
Índice de desempenho
Tensor de inercial em relação ao sistema o
Tensor de inercial do satélite sem as rodas
Tensor de inercial das rodas
Matriz de ganhos derivativos do controlador PD
Matriz de ganhos proporcionais do controlador PD
Quantidade de movimento angular do satélite referido e expresso no sistema inercial
Quantidade de movimento angular do satélite referido ao sistema inercial e expresso no sistema BF
Matriz de ponderação de estados
Matriz de covariânça do ruı́do nas medidas
Matriz de rotação do sistema b para o sistema a
Matriz de ponderação da lei de controle
Matriz de covariânça do ruı́do da dinâmica
Operador anti-simétrico
Torque interno
Matriz do sistema
Matriz dos atuadores
Matriz de sensores
Matriz função de transferência de um sistema
Quantidade de movimento angular nominal no eixo de pitch
Momentos principais de inércia
Função de transferência do controlador
Matriz de ganhos do regulador LQR
Matriz de ganhos do filtro de Kalman
Matriz solução da equação de Riccati no regime estacionário
para o caso do regulador
Matriz de covariânças dos estados estimados
Perı́odo orbital
Função de Lyapunov do sistema
Coordenadas do sistema inercial
Estimativa do estado x
Cosenos diretores de xo em relação aos eixos do corpo x, y, z
Cosenos diretores de yo em relação aos eixos do corpo x, y, z
Cosenos diretores de zo em relação aos eixos do corpo x, y, z
Sinal de comando
Quantida de movimento angular do satélite
hw
lw
mb
mc
q
r
A
D
e
e
n
gm
, hnm
h
h
i
ik
k
N
u
v
w
x, y, z
xo , yo , zo
y
Ω
Φ
β
ω bbw
ω bib
ω biw
ω bob
τ bbc
τ bcw
τ bc
τ bd
τ bg
τ bm
τa
τ ext
τs
Quantida de movimento angular das rodas relativo ao satélite
Quantida de movimento angular das rodas
Momento dipolo magnéico das bobinas
Dipolo magnético residual para o detumble
Quaternion
Sinal de referência
Área da secção da bobina magnética
Matriz de influênça do controle na saı́da
Excentricidade da órbita
Sinal de erro
Coeficiente de Gauss
Altura do satélite
Ganho do controlador EBC com realimentação de velocidade
angular
Inclinação da órbita
Corrente que passa pela bobina magnética disposta na direção
k
Ganho do controlador Bdot
Número de espiras
Sinal de contole
Ruı́do branco na medida
Ruı́do branco na dinâmica
Coordenadas do sistema do satélite
Coordenadas do sistema orbital
Sinal de saı́da do sistema
Ascensão do nodo ascendente
Ângulo de rotação em torno do vetor unitário λ
Ganho do controlador EBC com realimentação de atitude
Vetor do quaternion
Vetor velocidade angular das rodas referido ao sistema BF, expresso no sistema BF
Vetor velocidade angular do satélite referido ao sistema ECI,
expresso no sistema BF
Vetor velocidade angular das rodas referido ao sistema ECI,
expresso no sistema BF
Vetor velocidade angular do satélite referido ao sistema OF,
expresso no sistema BF
Torques de controle das bobinas magnéticas ou magneto torques
Torques de controle das rodas
Torques de controle expresso no sistema BF
Torques extenos que agem sobre o satélite expresso no sistema
BF
Toque de gradiente de gravidade
Torque magnético devido ao dipolo residual
Toque de arraste aerodinâmico
Torque(s) exteno(s) que agem sobre o satélite
Toque de pressão de radiação solar
η
ωo
φr , θr , ψr
φ
ψ
θ
Parte real do quaternion
Velocidade orbital média
Atitude de referência
Ângulo de rotação em torno do eixo de roll
Ângulo de rotação em torno do eixo de yaw
Ângulo de rotação em torno do eixo de pitch
LISTA DE ABREVIATURAS E SIGLAS
ACS
BF
EBC
ECI
EQUARS
FG
IAGA
IGRF
INPE
LQG
LQR
MIMO
OF
OOF
PD
PMM
POF
SISO
VLHL
Sistema de Controle de Atitude
Referencial do Satélite
Controladores Baseados em Energia
Referencial Inercial
Equatorial Atmosphere Research Satellite
Rereferencial Cartesiano Terrestre
Association of Geomagnetism and Aeronomy
International Geomagnetic Reference Field
Istituto Nacional de Pesquisas Espaciais
Regulador Quadrático Gaussiano
Regulador Linear Quadrático
Multiplas Entradas Multiplas Saı́das
Referencial Orbital
Referencial da Órbita
Proporcional Derivativo
Plataforma Multi-Missão
Referencial Pseudo Orbital
Única Entrada Única Saı́da
Vertical Local e Horizontal Local
CAPÍTULO 1
INTRODUÇÃO
A orientação de um satélite, em relação a um sistema de referência conhecido, é
denominada atitude e o movimento de rotação em torno do seu centro de massa
é denominado movimento de atitude. Assim, a atitude e o movimento de atitude
especificam a orientação espacial e o movimento rotacional em torno do centro de
massa do satélite (Wertz, 1978 e Hughes, 1986). Para determinar a atitude de um
satélite em relação a um sistema de referência, esse deve estar equipado com sensores
que possam fornecer a sua orientação em relação ao Sol, a Terra, e/ou alguma estrela
fixa bem como em relação ao vetor campo magnético terrestre (direção e magnitude).
Em certos casos é necessário considerar elementos orbitais do satélite para que seja
possı́vel determinar completamente a atitude e o movimento de atitude do veı́culo
espacial. A análise de atitude pode ser dividida em determinação, previsão e controle
de atitude (Wertz, 1978).
Este trabalho trata da previsão e o controle de atitude. A previsão é o processo de
prever a orientação do satélite pelo uso de modelos que permitam extrapolar sua atitude, sendo necessário conhecer as forças perturbadoras que agem sobre o satélite,
portanto, necessário modelá-las. O controle é o processo de orientar o satélite de
maneira que esse adquira ou mantenha a atitude prefixada pela missão. A implementação do modelo cinemático e dinâmico do satélite é feita neste trabalho utilizando o software MATLAB e SIMULINK, bem como a implementação das leis e
estratégias de controle de atitude em três eixos.
Satélites artificiais, em uma grande variedade de missões espaciais, seja para fins
meteorológicos, de telecomunicações, de sensoriamento remoto e cientı́fico, devem
estar com uma ou mais faces/equipamentos orientada(o)s para direções especı́ficas,
tais como para a Terra (BrasilSAT, Syncom, etc), para o Sol (Soho) ou para outras
estrelas (Hubble). É evidente a necessidade de se conhecer a atitude do satélite, de
modo que se possa estabilizar seu movimento de acordo com a atitude nominal especificada. Isto é feito através de um sistema de controle, projetado de acordo com
os requisitos da missão. O controle de atitude é, portanto, o processo de orientar o
satélite de maneira que este adquira ou mantenha a atitude nominal especificada.
Este processo pode ser visto sob dois modos: modo de aquisição de atitude, que
31
consiste em levar o satélite para a atitude nominal a partir de uma atitude qualquer e modo normal de operação, que consiste em estabilizar/manter a atitude, de
acordo com a atitude nominal, fazendo correções, quando necessárias. Durante a
fase de aquisição de atitude são, em geral, necessárias manobras de grandes ângulos,
uma vez que o sistema de controle de atitude (ACS) deve ser capaz de, a partir de
uma atitude qualquer, levar o satélite para a atitude nominal. Já a fase de estabilização/manutenção, em geral, requer manobras de pequenos ângulos objetivando
correções na atitude. Eventualmente os requisitos do ACS podem impor manobras
de grandes ângulos para reorientar o satélite em determinada direção. Um requisito
deste tipo requer manobra de reaquisição da atitude.
As duas configurações básicas para estabilizar satélites são aquelas em um eixo
(exemplos do SCD1 e SCD2) e em dois eixos, essa conhecida como estabilização
em três eixos (exemplos do CBERS1 e CBERS2), uma vez que estabilizando-se
dois eixos garante-se a estabilização do terceiro eixo. O estudo da estabilização
em um eixo é descrito em Kuga L. D. e Guedes (1987), Kuga e Guedes (1987),
Quirelli (2002) e Zanardi et al. (2003a). Esse trabalho é dedicado ao estudo da
estabilização em três eixos de satélites artificiais. A estabilização em três eixos pode
ser descrita em relação ao sistema orbital. Nesse referencial o movimento em torno
da direção da velocidade orbital é denominado roll (rolamento). O movimento em
torno da direção normal à órbita é denominado pitch (arfagem), e finalmente o
movimento em torno da direção Nadir/Zênite é denominado yaw (guinada). Satélites
estabilizados em três eixos podem ser vistos como um sistema com quantidade de
movimento angular nulo (zero momentum system) ou então como um sistema em
que a quantidade de movimento angular é não nula (bias momentum system). Os
sistemas com quantidade de movimento angular não nulo contém, geralmente, um
volante de inércia com quantidade de movimento angular nominal diferente de zero
(Wertz, 1978 e Wie, 1998). Este equipamento difere das rodas de reação pelo fato de
ter velocidade angular nominal diferente de zero em certa direção (Larson e Wertz,
1992).
Atualmente o programa espacial brasileiro desenvolve, entre outros, projetos de dois
satélites; o Satélite Equatorial Upper Atmosphere Research Satellite de aplicações
cientı́ficas (EQUARS) e a plataforma multi-missão (PMM), visando uma plataforma única para satélites com diferentes missões. Os requisitos dessas missões
impõem um sistema de controle em três eixos. Essas duas missões espaciais consti32
tuem a principal motivação para o desenvolvimento do trabalho.
Diversas técnicas de estabilização podem ser utilizadas quando se deseja estabilizar
satélites em três eixos. Dentre elas destacam-se nesse trabalho, aquelas que empregam como atuadores: rodas de momentum e bobinas magnéticas. Nesse trabalho são
estudados/analisados três modos de operação:
• Detumble: reduzir a velocidade angular do satélite.
• Estabilização: adquirir e manter o referencial do satélite alinhado com o
referencial orbital (vertical local e horizontal local).
• Aquisição: aquisição de uma orientação arbitrária.
Os modelos matemáticos do veı́culo espacial usados nos procedimentos de controle
de atitude avaliados para o modo de detumble e o modo de estabilização
são: rodas de reação combinadas com bobinas magnéticas e um sistema empregando
somente bobinas magnéticas. As bobinas são utilizadas para a redução da velocidade
angular, ou seja, no modo de detumble e desaturação. Entretanto, nesse estudo a
estratégia de dessaturação não é estudada. Para o modo de aquisição é avaliado
o emprego de rodas com objetivo de analizar a exequibilidade do procedimento a
partir das especificações das rodas para o satélite EQUARS.
As leis e estratégias de controle empregadas nesse trabalho são associadas à técnica
do Regulador Linear Quadrático (LQR), LQR tracking, Regulador Quadrático Gaussiano (LQG), controle proporcional derivativo (PD), controladores baseados em
energia (attitude feedback e angular velocity feedback ) usando a teoria de Lyapunov
e o controlador Bdot que também é baseado na teoria de Lyapunov. As estratégias
de controle adotadas para os três modos de operação são:
• O controlador Bdot, usando apenas medidas dos magnetômetros, para o descapotamento/detumbling do satélite. Os atuadores empregados são as bobinas magnéticas.
• Os controladores LQR/LQG e os baseados em energia para o modo de estabilização. A metodologia LQR é usada para rodas e bobinas e a metodologia
33
LQG apenas para rodas. Os controladores baseados em energia são usados com
o emprego de bobinas magnéticas.
• A metodologia LQR tracking e o controlador PD para o modo de aquisição,
com o objetivo de avaliar a exequibilidade/factibilidade do emprego das rodas
especificadas para o satélite brasileiro EQUARS.
A dissertação está organizada de acordo com a seqüência descrita a seguir.
No Capı́tulo 1 é feita uma breve introdução sobre o trabalho, descrevendo brevemente
os procedimentos e estratégias de controle adotadas. Nesse Capı́tulo é incluı́da a
motivação do estudo, que se refere à aplicação desse estudo ao satélite brasileiro
EQUARS.
O Capı́tulo 2 apresenta os objetivos do trabalho.
O Capı́tulo 3 apresenta a revisão bibliográfica sobre o tema.
No Capı́tulo 4 é realizada uma revisão das definições e fundamentos matemáticos
envolvendo atitude e movimento de atitude de veı́culos espaciais.
No Capı́tulo 5 são apresentados os modelos matemáticos da cinemática e dinâmica do
veı́culo espacial. É mostrado o modelo do campo magnético terrestre e o modelo dos
principais torques ambientais (gradiente de gravidade, pressão de radiação, arrasto
atmosférico e dipolo residual). No Capı́tulo são discutidos brevemente os torques de
origem interna e uma aproximação do seu modelo. É apresentado ainda a linearização
dos modelos dinâmicos e cinemáticos da atitude para os diferentes procedimentos
adotados.
No Capı́tulo 6 são apresentados as estratégias de controle, os fundamentos e metodologias usadas na formulação das leis de controle empregadas nos sistemas de
controle de atitude (ACS).
No Capı́tulo 7 são apresentadas as simulações da dinâmica e controle do veı́culo para
as diferentes configurações do ACS, envolvendo os diferentes modos de operação. O
estudo é aplicado ao satélite brasileiro EQUARS. É feita a análise dos resultados
34
e comparações das estratégias de controle empregadas.
Capı́tulo 8 encerra-se com a análise e discussões dos resultados, sendo feitas sugestões
para trabalhos futuros.
35
36
CAPÍTULO 2
OBJETIVO
Este trabalho visa modelar, simular, analisar, comparar resultados e apresentar um
estudo de alternativas de sistemas de controle de atitude (ACS) para satélites estabilizados em três eixos em diferentes modos de operação, utilizando como atuadores:
• Rodas de reação/volantes de inércia (reaction wheel /momentum wheel );
• Bobinas magnéticas.
Os procedimentos adotados para os modos de detumble e estabilização
usando os atuadores citados acima são: 1) satélite equipado com rodas e bobinas
magnéticas; 2) satélite equipado apenas com bobinas magnéticas. O projeto de controle usado para o modo de estabilição utiliza diferentes reguladores: 1) os baseados nas técnicas modernas de controle multivariável: Regulador Linear Quadrático
(LQR) e Regulador Quadrático Gaussiano (LQG); 2) os controladores baseados em
energia (energy based control ) com realimentação de velocidade angular (angular
velocity feedback ) e realimentação de atitude (attitude feedback ), ambos utilizando
a teoria de Lyapunov. O projeto de controle usado para o modo de detumble usa
apenas medidas dos magnetômetros e é baseado na teoria de Lyapunov. O projeto
de controle para o modo de aquisição utiliza: 1) controlador Proporcional Derivativo (PD); 2) LQR rastreio/tracking. Esse projeto tem o objetivo de avaliar e
discutir a factibilidade do emprego das rodas especificadas para o satélite brasileiro
EQUARS, em adquirir uma atitude arbitrária, e ainda comparar os controladores
usados.
O objetivo principal desse trabalho é realizar um estudo de alternativas e um
estudo comparativo em função do desempenho e análise de custo & benefı́cio
dos diferentes sistemas de controle de atitude. As configurações dos ACS(s) estudados serão aplicados ao satélite brasileiro EQUARS (Equatorial Upper Atmosphere Research Satellite). O estudo desenvolvido nesse trabalho (estudo de alternativas/viabilidade) pode se extendido por exemplo, a plataforma Multi-Missão
(PMM).
37
2.1
Meios e Recursos
Os principais meios e recursos utilizados para se alcançar os objetivos propostos no
trabalho são:
• pesquisa bibliográfica utilizando os recursos da biblioteca do INPE, artigos
publicados na literatura da área, e pesquisa no portal de periódicos da CAPES;
• uso dos recursos da divisão de Mecânica Espacial e Controle, no que se refere
a meios computacionais para simulação, laboratório LABSIM, incluindo os
ambientes MATLAB e SIMULINK para o desenvolvimento das simulações.
2.2
Metodologia
O desenvolvimento do trabalho proposto está fundamentado nas seguintes metodologias:
• Modelagem matemática da dinâmica de atitude para um veı́culo espacial contendo como atuadores: rodas (de reação e/ou volantes de inércia) e bobinas
magnéticas, incluindo os efeitos do gradiente de gravidade.
• Análise dinâmica do modelo matemático, incluindo a linearização do modelo,
usados no desenvolvimento do projeto de controle (reguladores), para os controladores lineares (PD, LQR e LQG).
• Modelagem e implementação, no ambiente SIMULINK, da órbita e do campo
magnético terrestre, modelo IGRF (International Geomagnetic Reference Field ).
• Estudo do problema do regulador linear quadrático (LQR) e regulador
quadrático gaussiano (LQG).
• Estudo dos controladores baseados em energia, utilizando a teoria de Lyapunov.
• Estudo do controlador proposto por Wisniewski (1996) ou Bdot.
38
• Implementação dos sistemas dinâmicos no ambiente MATLAB e SIMULINK.
• Simulação dos sistemas dinâmicos (não lineares).
• Análise e comparação das estratégias de controle empregadas nas diferentes
configurações dos ACS(s), desenvolvidos para satélites estabilizados em três eixos, nos diversos modos de operação. Aplicação do estudo ao satélite brasileiro
EQUARS.
39
40
CAPÍTULO 3
REVISÃO BIBLIOGRÁFICA
Este Capı́tulo apresenta uma revisão da literatura relacionada ao principais assuntos
envolvidos nesse trabalho: controle de atitude de satélites, técnicas de estabilização
em três eixos, modelamento matemático dos sistemas (satélite, atuadores, meio) e
estratégias de controle. As metodologias modernas de controle multivariável LQR,
LQR/tracking e LQG são revistas, controladores PD e reguladores não lineares baseados em energia também são apresentados.
3.1
Atitude
As referências Wertz (1978), Kaplan (1976), Wie (1998) e Hughes (1986) apresentam
um estudo do movimento de atitude de satélites artı́ficiais, mostrando o histórico,
os fundamentos e os conceitos fı́sicos fundamentais utilizados na previsão, controle
e determinação de atitude. A modelagem, análise e controle de sistemas dinâmicos
são revistos, as técnicas de controle de atitude para satélites artificiais estabilizados
em três eixos também são discutidas. Os tipos e descrição dos sensores e atuadores
empregados no movimento de atitude são apresentados em Wertz (1978) e em Pilchowski (2001). Uma descrição do uso de sensores em satélites é mostrada em Wright
e Wong (1989). As várias formas de representação/parametrizações da atitude (quaternions, ângulo eixo equivalente, ângulos de Euler, cossenos diretores, variáveis de
Andoyer, parâmetros de Gibbs) são apresentadas em Wertz (1978), Rodrigues e Zanardi (2004), Wie (1998), Fauske (2003) e Junkins e Turner (1986). O problema de
estabilização em três eixos é também discutido em Martins Neto (2001).
A análise de missão é uma fase na qual define-se o que deve ser feito sem necessariamente definir como fazer. A análise de missão envolvendo especificações do sistema
de controle de atitude para estabilização em três eixos é descrita e discutida em
Larson e Wertz (1992). Em Marteau e Rogers (1996) é feita uma análise de sistemas de controle de atitude (ACS), abordando aspectos gerais sobre as especificações
(atuadores, sensores, computadores de bordo) de um sistema de controle de atitude
(ACS). O Autor discute um ACS com custo “razoável”, para pequenos satélites
(< 500kg) que requeiram apontamento de “razoável”precisão (20 arcsec), precisão
tı́pica requerida para satélites já operacionais como IRAS (Infrared Astronomical
41
Satellite) (Beichman et al., 2004), ASCA (Advanced Satellite for Cosmology and
Astrophysics) (Tamura, 1998) e SOHO (Solar & Heliospheric Observatory) (Gurman, 2004).
A modelagem da dinâmica de atitude utilizando-se rodas (de reação e volantes de
inércia) para estabilização em três eixos é apresentada nas referências Wie (1998)
(onde se utilizam duas rodas de reação e um volante de inércia), Yairi (1994) e Fichter
e Zentcraf (1996). Yairi (1994) e Fichter e Zentcraf (1996) apresentam a mesma
modelagem de Wie (1998), mas utilizam quatro rodas de reação, sendo uma delas
disposta na diagonal (skew symmetric), redundante. A roda redundante é disposta
de tal forma que forneça uma quantidade de movimento angular nas três direções
principais de inércia, o que pode, eventualmente, em caso de falhas substituir uma ou
mais rodas ao longo dessas direções. Essas referências mostram o desenvolvimento
das equações do movimento rotacional de um corpo rı́gido equipado com rodas.
O problema ótimo no procedimento de estabilização de satélites artifı́ciais com a
utilização de rodas é apresentado em El-Gohary (2003). Este autor apresenta um
estudo de estabilidade segundo Lyapunov. Em Varatharajoo e Fasoulas (1975) é feita
a análise do problema de atitude empregando rotores para satélites de observação
da Terra. Spindler (2000) aborda o problema de controle para estabilização em
três eixos empregando N volantes de inércia, exemplificando o procedimento para o
satélite MONS-ballerina.
O problema de estabilização de atitude em três eixos usando apenas atuadores eletromagnéticos é discutida em Kaplan (1976), Wertz (1978), Bushenkov e Smirnov
(2002), Psiaki (2001), Wang e Shtessel (1998), Wisniewski e Blanke (1999) e Musser
e Ebert (1989). Em Wertz (1978) e Carrara (1982) é encontrada a modelagem das
perturbações ambientais (forças e torques) que atuam sobre o satélite no espaço,
devido ao campo magnético, o campo gravitacional e a radiação proveniente do
Sol e da Terra. A modelagem do campo magnético, modelo IGRF (International
Geomagnetic Reference Field ) usado nesse trabalho é encontrado em Macmillan e
Quinn (2000). Outros modelos do campo magnético (dipolo e quadripolo) podem
ser encontrados em Zanardi et al. (2003b) e Zanardi et al. (2004).
O Desenvolvimento das equações da dinâmica de atitude usando atuadores eletromagnéticos/bobinas magnéticas é apresentado em Psiaki (2001), Musser e Ebert
(1989), Wisniewski (1997) e Marteau e Psiaki (1988). Outras referências também
42
trazem as equações da dinâmica e cinemática de atitude, como Bushenkov e Smirnov (2002), Fauske (2002), Grassi e Moccia (1995) e Wang e Shtessel (1998). Essas referências discutem a utilização de bobinas magnéticas para o controle para
pequenos satélites visando menor custo para precisão requerida de 1o a 2o , sem requerer alta fonte de energia. Wisniewski e Blanke (1999) analisam a técnica para o
satélite dinamarquês Orsted (Clausen, 2004) (satélite aplicado no estudo do campo
magnético da Terra). O autor mostra que é factı́vel obter estabilização em três eixos
utilizando apenas torques magnéticos para satélites de órbita baixa (LEO), sujeitos
ao gradiente de gravidade. O ASRI (Australian Space Research Institute) (Ardebil,
2004) confirma a viabilidade de se usar apenas atuadores magnéticos para o controle de atitude de pequenos satélites, através da missão cientı́fica TechSAT. Silani
e Lovera (2003) apresentam uma revisão do problema de estabilização de atitude em
três eixos para pequenos satélites, usando atuadores magnéticos (bobinas), baseado
na teoria de controle linear e não linear. Junkins e Carrington (1980) realizam um
estudo de otimização para manobras de atitude usando atuadores magnéticos.
A referência Kim e Choi (1999) trata da utilização de um sistema de controle de
atitude (ACS), utilizando três rodas de reação combinadas com bobinas magnéticas
para o microsatélite Koreano KITSAT-3 (satélite de telecomunicações), visando estabilização em três eixos, em órbitas baixas (LEO). Nos requisitos de payload /carga
útil do KITSAT-3 são exigidas alta precisão de apontamento (0.05o ) e estabilidade
(0.014o rad/s) da plataforma. Os atuadores magnéticos também podem ser utilizados
para a desaturação das rodas. Bang e Choi (2003) analisam a desaturação de rodas
utilizadas para manobras de grandes ângulos. Outra técnica de controle para estabilização em três eixos, desenvolvida primeiramente para os satélites TIROS (Television Infrared Observation Satellite), por Harold Perkel, utiliza uma configuração
particular de um volante de inércia (momentum wheel ), ao longo do eixo de arfagem (pitch), combinando com dois atuadores eletromagnéticos, dispostos ao longo
do eixo de rolamento (roll ) e guinada (yaw ). Esta técnica é denominada stabilite
e é encontrada na referencia Perkel (1966). Em Hamzah e Hashida (1999) é descrito
o desenvolvimento do sistema de controle de atitude para o satélite TiungSAT-1,
utilizando-se a configuração de Perkel obtendo uma precisão de ampontamento de
±1o .
Whitford e Forrest (1998) descrevem o desenvolvimento de um sistema de controle
de atitude (ACS) para estabilização em três eixos também baseado na combinação
43
de rodas (de reação e volantes de inércia) e bobinas magnéticas. Os vários modos
de operação são avaliados para sistema de controle de atitude (ACS). O ACS é aplicado ao satélite CATSAT (Co-operative Astrophysical and Technology Satellite) do
programa STEDI (Student Explorer Demonstration Initiative). Um dos principais
problemas inerentes à utilização de rodas é o da saturação no processo de controle,
devido a manobras de grandes ângulos (Bang e Choi, 2003) e/ou torques seculares.
Em Buckingham e Smirnov (1972) é discutida a desaturação de rodas de reação
empregando atuadores magnéticos. Em Gökçev e Meerkov (2001) é apresentada a
metodologia para sistemas com saturação de atuadores utilizados na estabilização
de satélites sujeito a perturbações seculares.
Da vasta literatura que apresenta o desenvolvimento das equações da dinâmica não
linear, para satélites equipados com bobinas, destacam-se os trabalhos de Fauske
(2002), Cohen (1973), Spencer (1977), Shigehara (1972) e Alfriend (). As equações
do sistema tornam-se variantes no tempo devido ao campo geomagnético. A análise
da controlabilidade de satélites equipados apenas com bobinas magnéticas é encontrada em Bhat e Dham (2003). Wisniewski e Markley (1999) desenvolvem uma
metologia de controle ótimo aplicada ao controle de atitude em três eixos. Propriedades elétricas de materiais utilizados em bobinas magneticas são apresentados e
discutidos em Legg (2003). A descrição, especificações e precisões de rodas (de reação
e volantes de inércia) da TELDIX utilizadas para estabilização em três eixos são
encontrados em Auer (1983) e Heidelberg (2004). Em Auer (1983) são também apresentadas repostas para os comandos de controle.
3.2
Controle
Nesta seção é feita um revisão bibliográfica das metodologias LQR, LQR rastreio
(tracking) e LQG. Uma revisão de controladores baseados em energia, também reguladores, obtidos a partir da teoria de Lyapunov é apresentada. O controlador usado
para o modo de detumble, conhecido como Bdot é revisto.
Metodologia LQR e LQG
A teoria do regulador linear quadrático (LQR) e do regulador linear gaussiano (LQG)
aplicadas em controle multivariável é encontrada em uma vasta literatura. As prin44
cipais referências estudadas pelo autor para elaboração e implementação dessas metodologias no sistema de controle de atitude foram Dorato e Cerone (1995), Maciejowski (1989), Kwakernaak e Sivan (1972), Moore e Anderson (1990), Brown (1997)
e Kirk (1970). Essas referências fornecem os métodos para formular a lei de controle
para sistemas dinâmicos. O projeto de controle LQR e LQG é baseado na linearização dos sistemas dinâmicos, definindo uma função objetivo a ser minimizada, e
na obtenção de uma matriz de ganhos (variantes no tempo ou não) usada na realimentação (Overby, 2004), o método se aplica a sistemas variantes no tempo, como
o caso de satélites equipados com bobinas magnéticas, tornando-se uma alternativa
atraente, devido a robustez e fácil implementação, essa metodologia é usada por
Wisniewski e Markley (1999), Wisniewski (1996) e Marteau e Psiaki (1988) para
o sistema de controle de atitude de satélites equipados com atuadores magnéticos.
Yairi (1994) usa a formulação LQR para o sistema de controle de atitude de satélites
estabilizados em três eixos equipados com rodas.
Nesse trabalho a lei de controle, formulada a partir do sistema linear, é aplicada
ao sistema dinâmico completo, não linear, assim como feito por Musser e Ebert
(1989) e Fauske (2002). No controle de atitude para satélites estabilizados em três
eixos é necessário manter a atitude em relação a um sistema de referência, inercial ou
não inercial, como é o caso de satélites Terra apontados em que o sistema do satélite
deve estar alinhado com a vertical local e horizontal local (VLHL) ou referencial
orbital, um sistema de referência móvel (Wertz, 1978). Portanto, no problema de
estabilização a referência para o controlador é zero ou constante, tratando-se de um
regulador (Distefano e Willians, 1972). Em Kirk (1970) e Dorato e Cerone (1995) é
apresentada a formulação das lei de controle para o caso em que a referência difere
de zero, essa formulação é conhecida como LQR rastreio/tracking (Dorato e Cerone,
1995).
Em Overby (2004) é discutido várias estratégias de controle e reguladores usados
para a estabilização do satélite Norueguês Ncube. Entre essas estratégias a metodologia LQR é apresentada, mostrando boas propridades de estabilidade. Overby
(2004) e Kristiansen (2000) apresentam um algoritmo para a seleção dos ganhos das
matrizes de ajuste usadas no LQR, baseado na saturação dos atuadores e intervalos de operação para os estados. Em Dorato e Cerone (1995), Kwakernaak e Sivan
(1972) e Souza (1987) é discutida a escolha/seleção das matrizes dos ganhos com
base nos requisistos e compromissos/trade-off do processo. Na metologia LQR não
45
é levada em conta as incertezas do problema, ou seja, a metodologia LQR é determinı́stica. Através da utilização da metodologia LQG, um controlador baseado
em observador, essas incertezas (estocásticas) podem ser atribuı́das à planta ou processo. A solução do problema do LQG é obtida pelo uso do princı́pio da separação
que possibilita a separação do problema original (estocástico) em dois problemas:
1) estimação ótima do estado, de modo que sua covariância seja minimizada.
Esse sub-problema é resolvido pelo uso de um filtro de Kalman, ignorando-se completamente o problema de controle; 2) formulação da lei de controle usando a
metodologia LQR, fazendo uso da estimativa como se fosse a medida exata do estado, ignorando-se completamente os aspectos ou natureza estocásticas do processo
(Flora, , Moore e Anderson, 1990 e Dorato e Cerone, 1995).
A solução do problema 1 consiste na determinação da matriz ganho ou ganho de Kalman, a estrutura do filtro é encontrada em Maciejowski (1989), Kwakernaak e Sivan
(1972) e Athans (1971). A solução do problema 2 consiste na derminação da matriz
ganho do controlador usado na realimentação. Para a solução do problema LQG
ser assintoticamente estável é necessário que o sistema seja completamente observável, para a solução do problema 1 e completamente controlável, para a solução
do problema 2 (Brown, 1997 e Dorato e Cerone, 1995).
Controladores baseados em energia
Uma importante ferramenta em teoria de controle é o uso de controladores não lineares baseados em energia (energy based control ), desenvolvidos através da teoria de
Lyapunov (Fauske, 2002 e Khalil, 2000). Fauske (2002) sugere diferentes controladores baseados em energia para satélites equipados apenas com bobinas magnéticas;
as leis de controle avaliadas para o modo de estabilização são: 1) realimentação de
velocidade angular (angular velocity feedback ); 2) realimentação de atitude (attitude
feedback ). Overby (2004) realiza um estudo comparativo dos controladores 1 e 2 com
o controlador obtido pela teoria linear ou seja usando a metodologia LQR. Essas estratégias de controle são aplicadas ao satélite Norueguês NCUBE. Hegrenes (2004)
realiza um estudo comparativo entre um controlador usando o método LQR e um
controlador PD, para o satélite NCUBE, no seu modo de estabilização.
A teoria de Lyapunov é desenvolvida em detalhes nas referências Khalil (2000) e
46
Overby (2004). Para o modo de operação de descapotamento/detumble, a lei de
controle é sugerida por Wisniewski (1996) e é baseada na teoria de Lypunov. Essa
lei de controle também é conhecida como Bdot (Silani e Lovera, 2003). Wells e
Jeans (2002) utilizam essa estratégia de controle para o satélite Canadense CanX1. Nesse trabalho a estratégia de controle de Wisniewski (1996) é usada para o
modo de operação de detumble. Portanto nos dois sistemas de controle de atitude
(ACS(s)), o detumble é feito com o uso de bobinas magnéticas.
3.3
Literatura do INPE
No acervo do INPE pesquisado, destacam-se o trabalho de Souza (1981), que apresenta o estudo e desenvolvimento de um sistema de controle ativo de atitude em
três eixos para satélites artificiais, usando atuadores pneumáticos a gás frio e rodas
de reação. O estudo é feito através de simulação digital; Souza (1987) apresenta
o estudo de um sistema de controle de atitude em três eixos utilizando rodas de
reação. O autor ainda apresenta o modelo dinâmico do satélite equipado com rodas,
utilizando as equações de corpo rı́gido ou equações de Euler, assim como desenvolvido nesse trabalho; Carrara (1982) faz o estudo das principais forças e torques
atuantes em satélites artificiais, como gradiente de gravidade, pressão de radiação
(solar e terrestre) e arrasto atmosférico, desenvolvendo um programa/software para
o computo numérico dessas influências na atitude e órbita do satélite, a partir da
sua geometria.
Souza (1987) apresenta o estudo e o projeto experimental de um atuador do tipo
roda de reação. Sua aplicação está no controle fino da atitude de satélites artificiais
atráves da troca de quantidade de movimento angular entre a roda e o satélite,
realiza, ainda testes de desempenho e confrontação com os resultados obtidos em
simulação digital. Trivelato (1988) apresenta o desenvolvimento de controladores
digitais de torques de rodas de reação, usadas no controle de atitude de satélites
artificiais, o autor realiza simulações digitais e simulação com hardware na malha de
controle, qualificando o controlador da roda para seu uso em sistemas de controle de
atitude. O estudo detalhado de um sistema de controle por referência para operação
de rodas é apresentado em Trivelato (1988).
O projeto de controle de atitude de satélites estabilizados em três eixos utilizando a
47
metodologia LQR e LQG é apresentada em Moscati (1992) mostrando ser uma opção
viável para o sistema de controle de atitude, atravês dos resultados obtidos. Moscati
(1992) utiliza o modelo linear em suas simulações. Outros trabalhos se destacam no
estudo de atitude de satélite que incluem apêndices flexı́veis usando metodologias
modernas e distintas de controle (LQG/LTR (Loop Transfer Recovery) e H-infinito),
Roma (1991) desenvolve um software de um manipulador simbólico para obtenção
das equações do movimento de atitude. Flora () realiza o projeto de um sistema de
controle de atitude em três eixos usando os métodos LQG/LTR e H-infinito. Souza
(1987) apresenta uma extensão na teoria do regulador linear quadrático (LQR) para
o caso de uma dinâmica não linear, sugerindo uma nova lei de controle de atitude
para satélites artificiais.
A partir da revisão da literatura do INPE nota-se uma grande diversidade e qualidade nos trabalhos voltados para o estudo de atitude em três eixos, e implementações
de metodologias de controle. Entretanto, os trabalhos revistos tratam de técnicas
que empregam apenas rodas como atuadores, incluindo até mesmo o projeto da
mesma. As análises e simulações se restringem a um único modo de operação, o de
estabilização, envolvendo apenas manobras de pequenos ângulos. O presente trabalho que inclue o estudo de atitude empregando rodas, especificadas para o satélite
EQUARS, é usado para outros modos de operação, estabilização e aquisição. O modelo do satélite equipado com atuadores magnéticos também é desenvolvido. Esse
trabalho aborda o problema do uso de rodas usando as metodologias usadas em
Souza (1987), Flora () e Moscati (1992) para a obtenção da lei de controle.
Na simulação digital foi adotada o modelo completo da dinâmica (não linear e acoplada), assim como feito por Souza (1987). Na literatura disponı́vel do INPE não é
encotrada nenhuma referência que trata da utilização de bobinas magnéticas como
atuadores para a concepção do sistema de controle de atitude para satélites estabilizados em três eixos. Constitue, portanto, uma importante contribuição e desafio
desse trabalho, o estudo da estabilização em três eixos de satélites artificiais através
de atuadores magnéticos. Cabe salientar que esse estudo foi motivado pelos projetos dos satélites Brasileiros EQUARS e a Plataforma Multi-Missão (PMM)
que estão sendo desenvolvidos pelo INPE. Esses satélites empregam como atuadores
rodas e bobinas magnéticas para o controle de atitude. Como contribuição desse estudo desenvolveu-se um toolbox em SIMULINK, voltados para dinâmica e controle
de atitude.
48
Neste trabalho as simulações digitais são feitas no ambiente MATLAB e SIMULINK. Nesses recursos computacionais são encontrados nas referências Hanselman e
Benjamin (1994), Levine e Leonard (1995), Hanselman e Littlefield (2003) e Bishop
(1997) e na documentação disponibilizada pela Mathworks. Levine e Leonard (1995).
Ogata (1998) e Bishop (1997) apresentam uma introdução do uso do MATLAB em
engenharia de controle e automação. Em Hanselman e Benjamin (1994) e Wie (1998)
é apresentado a implementação das metodologias LQR e LQG através das functions do control toolbox do MATLAB. Hanselman e Littlefield (2003) mostra como
usar os recursos do MATLAB e apresentam todas as caracterı́sticas e benefı́cios da
programação em MATLAB. Hanselman e Littlefield (2003) mostram como realizar
o desevolvimento de um toolbox, usando SIMULINK. Para a geração do software de
atitude em SIMULINK, foram usados além desses recursos a documentação disponibilizada pela Mathworks (http://www.mathworks.com/support/).
49
50
CAPÍTULO 4
DEFINIÇÕES DE NOTAÇÕES
A posição e a atitude de um satélite podem ser representadas de várias formas
(Wertz, 1978, Larson e Wertz, 1992). Essa Seção introduz definições e notações
necessárias para o desenvolvimento do trabalho proposto.
4.1
Sistemas de Referência
Os sistemas de referência utilizados nesse trabalho são definidos a seguir, a Figura
(4.1) ilustra os referenciais definidos.
4.1.1
Referencial Inercial
O referencial ECI, no qual as leis de Newton do movimento são formuladas, é considerado um referencial inercial (Fauske, 2002). O referencial inercial (X, Y, Z) definido
tem origem no centro da Terra. O eixo Z aponta na direção do polo norte geográfico,
o eixo X aponta na direção do ponto vernal γ. O eixo Y é encontrado usando a regra
da mão direita, completando o sistema dextrogiro.
Nota: De fato, o referencial acima é aproximadamente inercial, devido ao movimento
da Terra; para as escalas de tempo adotadas, entretanto, o referencial ECI constitui
uma boa aproximação de um referencial inercial (Chobotov, 1996).
4.1.2
Referencial Orbital
O referencial orbital (OF) definido por (xo , yo , zo ) é um sistema de coordenadas com
origem no centro de massa do satélite. O eixo zo aponta na direção do centro da
Terra (direção Nadir) e o eixo yo aponta na direção normal a órbita. O eixo xo é
obtido pela regra da mão direita, e coincide com a direção do vetor velocidade orbital
linear, para uma órbita circular (Wertz, 1978).
Nota: Esse referencial também é referido como Vertical Local e Horizontal Local
(VLHL).
51
4.1.3
Referencial do Satélite
O referencial do corpo (BF), ou do satélite definido por (x, y, z) é um sistema de
coordenadas com origem no centro de massa do satélite. Os eixos são escolhidos como
sendo coincidentes com os eixos dos momentos principais de inércia. Para estudos
de satélite estabilizados em três eixos, Terra-apontado, é prático definir os eixos de
roll, pitch e yaw como sendo (Wie, 1998, Moscati, 1992)
• eixo de roll x, nominalmente alinhado com xo
• eixo de pitch y, nominalmente alinhado com yo
• eixo de yaw z, nominalmente alinhado com zo
Figura 4.1- Sistemas de referência, inercial (ECI), orbital (OF) e do satélite (BF).
52
4.2
Representação da Atitude
A atitude em três eixos é convenientemente representada através de uma transformação de coordenadas, na qual transforma-se um conjunto de coordenadas no
espaço inercial em um conjunto de coordenadas fixadas ao satélite (Wertz, 1978).
Existem parametrizações/representações alternativas para estas transformações. A
principais representações são:
• Matriz de rotação ou matriz de atitude;
• Ângulos de Euler;
• Quaternions ou parâmetros de Euler;
• Ângulo eixo equivalente;
Cada um desses conjuntos de parametrizações apresentam vantagens e desvantagens. Em Wertz (1978), Hughes (1986), Fauske (2003) e Wie (1998) são discutidas
pormenorizadamente cada uma dessas parametrizações.
4.2.1
Matriz de Rotação
A matriz de rotação, também chamada de matriz de cossenos diretores, tem as
interpretações (Fauske, 2003):
• Descreve a orientação mútua entre dois referênciais, onde cada vetor coluna
são os cossenos dos ângulos entre os dois referenciais;
• Transforma um vetor representado em um referencial para outro;
A matriz de rotação R de um referencial a para um referencial b denotada por Rba
é um elemento do conjunto SO (3), definido como
SO (3) = R | R ∈ <3×3 , RT R = I e detR = 1 ,
53
(4.1)
onde I é a matriz identidade 3 × 3.
A notação seguinte é usada para transformar um vetor r de um referencial para
outro:
f rom
rto = Rto
f rom r
(4.2)
onde o ı́ndice superior descreve em qual referencial o vetor está expresso.
Devido a propriedade de ortogonalidade, RT R = I, pode-se mostrar que a derivada
no tempo da matriz de rotação é dada por
Ṙab = S(ω aab )Rba
(4.3)
onde a notação usada ω aab representa a velocidade angular do referencial b em relação
ao referencial a, expressa no referencial a. Essa notação é usada em todo o trabalho.
A velocidade angular tem a propriedade ω aab = −ω aba (Fauske, 2003), e S(ω) é o
operador anti-simétrico (skew-symmetric):


0 −ωz ωy


S (ω) =  ωz
0 −ωx  ,
−ωy ωx
0
 
ωx
 
ω = ωy 
ωz
(4.4)
reescrevemos (4.3), temos
Ṙab = S(ω aab )Rba = −S(ω aba )Rba
Denotaremos a matriz de rotação por
54
(4.5)
Rba
=
h
cb2
cb1
cb3
i
(4.6)
T
onde cbi = cbix cbiy cbiz (i = 1, 2, 3) são matrizes coluna, tendo como componentes
os cossenos diretores dos eixos de a em relação aos eixos do referencial b. Daqui para
frente o ı́ndice b será usado para representar o referencial do corpo e o ı́ndice o
será usado para representar o referencial orbital. A ausência de ı́ndice representa o
referencial inercial.
4.2.2
Parâmetros de Euler
Os parâmetros de Euler, são também chamados quaternions (Fauske, 2003), e são
uma representação atrativa devido a parametrização sem singularidades e equações
diferenciais da cinemática serem lineares (Fauske, 2003, Wie, 1998). A representação
em quaternion requer menos tempo computacional que a representação por ângulos
de Euler (Arantes e Fonseca, 2004a), e portanto é usado em aplicações onde os
recursos computacionais são limitados (Wertz, 1978).
Um quaternion q é definido como uma grandeza hiper-imaginária composta por uma
parte real η e um vetor dada por
η = cos
Φ
,
2
= λ sin
Φ
2
(4.7)
que, representa uma rotação de Φ em torno do vetor unitário λ. Um quaternion
satisfaz o vı́nculo qT q = 1, ou seja
η 2 + 21 + 22 + 23 = 1
(4.8)
As equações diferenciais da cinemática são dadas por (Wertz, 1978, Wie e Arapostathis, 1989, Wie, 1985):
55
1
q̇ = Ωq
2
(4.9)

0
ωz −ωy ωx


−ωz

0
ω
ω
x
y

Ω=
 ω

−ω
0
ω
x
z
 y
−ωx −ωy −ωz 0
(4.10)
onde

A matriz de atitude (4.1) obtida a partir dos quaternions é dada por (Wie, 1998)
Rba = η 2 − qT q I + 2qqT − 2ηQ ()
(4.11)
onde Q () é o operador anti-simétrico dado por


0 −3 2


Q () =  3
0 −1 
−2 1
0
(4.12)
Dada a matriz de atitude (4.1), podemos determinar η e pelas relações
1
η = (1 + c11 + c22 + c33 ) 2
para 0 ≤


c23 − c32
1 

=
c31 − c13 
4η
c12 − c21
Se η = 0 escolhe-se outra sequência de rotações.
56
se η 6= 0
Φ
≤π
2
(4.13)
(4.14)
4.2.3
Ângulos de Euler
A orientação de um corpo também pode ser descrita por três ângulos (três
parâmetros independentes) denominados ângulos de Euler. A Figura (4.2) ilustra
a seqüência de rotações necessárias para levar um referencial X̄ ≡ (X, Y, Z) a outro
x̄ ≡ (x, y, z), que são listadas como.
0
0
0
• Rotação em torno do eixo Z de um ângulo φ que leva ao sistema ξ , θ , ζ .
0
• Rotação em torno do eixo ξ de um ângulo θ que leva ao sistema ξ, θ, ζ.
• Rotação em torno do eixo ζ de um ângulo ψ que leva ao sistema x, y, z.
Para essa seqüência de rotações 3 − 1 − 3 a matriz de atitude é dada por (Kaplan,
1976)


cos ψ cos φ − sin ψ cos θ sin φ
cos ψ sin φ + sin ψ cos θ cos φ sin ψ sin θ


RX
x = − cos ψ cos φ − cos ψ cos θ sin φ − sin ψ sin φ + cos ψ cos θ cos φ cos ψ sin θ 
sin θ sin φ
− sin θ cos φ
cos θ
Figura 4.2- Construção dos ângulos de Euler.
Em problemas de estabilização de atitude em três eixos é comum definir:
• ângulo de roll (φ) é o ângulo de rotação em torno do eixo de roll
57
• ângulo de pitch (θ) é o ângulo de rotação em torno do eixo de pitch
• ângulo de yaw (ψ) é o ângulo de rotação em torno do eixo de yaw
Note que os ângulos são definidos para rotações em torno de eixos distintos. Usados
comumente para relacionar o referêncial orbital (OF) ou LVHL com o referêncial
do satélite (BF), adequado para tais aplicações (Moscati, 1992). Existem doze (12)
conjuntos possı́veis de ângulos de Euler para descrever um referencial em relação a
outro (Wertz, 1978). Esses conjuntos/seqüências são dividas em dois tipos:
Tipo 1 (anti-simétrica): Nesse caso as rotações são feitas, sucessivamente, em
cada um dos três eixos. Esse tipo apresenta uma singularidade em θ = ± π2 . As
seqüências são 1 − 2 − 3, 2 − 1 − 3, 3 − 2 − 1, 2 − 3 − 1, 3 − 1 − 2, 1 − 3 − 2.
Tipo 2 (simétrica): Nesse caso a primeira e terceira rotação são feitas sobre o
mesmo eixo. Esse tipo apresenta uma singularidade em θ = π e θ = 0. As seqüências
são 3 − 1 − 3, 2 − 1 − 2, 1 − 2 − 1, 2 − 3 − 2, 1 − 3 − 1, 3 − 2 − 3.
As matrizes de rotação e as equações cinemáticas para cada uma dessas seqüências
de rotação (simétricas e anti-simétricas) são encontradas nas referência Wertz (1978)
e Hughes (1986). Apesar das equações da cinemática apresentarem singularidades do
tipo 1 e tipo 2, os ângulos de Euler tem um clara interpretação fı́sica e são utilizadas
na entrada e saı́da das simulações.
Dada a matriz de rotação (4.11) escrita em termos de quaternions, os ângulos de
roll (φ), pitch (θ) e yaw (ψ) para a seqüência de rotação 3 − 2 − 1, ou seja, primeira
rotação em yaw, segunda rotação em pitch e terceira rotação em roll, são dados por
(Wie, 1998, Hughes, 1986):
φ = arctan
c23
c33
0 ≤ φ ≤ 2π
π
π
θ = arcsin (−c13 ) − ≤ θ ≤
2
2
c12
ψ = arctan
0 ≤ ψ ≤ 2π
c11
58
(4.15)
(4.16)
(4.17)
Nota: A seqüencia de rotações escolhida para ser utilizada nesse trabalho é 3−2−1.
4.2.4
Ângulo Eixo Equivalente
O ângulo eixo equivalente é definido a partir do teorema de Euler que diz: um
corpo rı́gido pode atingir qualquer orientação/atitude a partir de uma orientação
arbitrária (de referência), através de uma rotação Φ sobre um eixo unitário fixo λ
em relação ao corpo e a referência. Esse eixo é chamado de eixo equivalente ou eixo
de Euler.
Definindo
 
λ1
 
λ = λ2 
λ3

e

0 −λ3 λ2


Λ =  λ3
0 −λ1 
−λ2 λ2
0
(4.18)
Pode-se mostrar que a matriz de atitude é dada por (Wie, 1998, Fauske, 2002):
Rba = cos ΦI + (1 − cos Φ) λλT − sin ΦΛ
(4.19)
onde I é a matriz identidade. A expressão (4.19) é também conhecida como fórmula
de Rodrigues.
Dada a matriz de atitude (4.1) Φ pode ser obtido como
cos Φ =
1
(c11 + c22 + c33 )
2
e pode se mostrar (Wie, 1998) que λ é dado por
59
(4.20)


 
c23 − c32
λ1
1 

 
λ = λ2  =
c31 − c13 
2 sin Φ
c12 − c21
λ3
para Φ 6= 0, ±π, ±2π, . . .
(4.21)
Usando a parametrização (4.19) e a identidade trigonométrica (4.7), a matriz de
rotação pode ser escrita em função dos quaternions como
Rba = 1 + 2ηS () + 2S2 ()
4.2.5
(4.22)
Rotações Infinitesimais
No modo de estabilização durante a manutenção da atitude nominal, em geral, são
realizadas manobras de pequenos ângulos para correções da atitude, justificando
aproximar a matriz de atitude por (Moscati, 1992)


1
ψ −θ


Rbo = −ψ 1
φ
θ −φ 1
(4.23)
que é a matriz de transformação do sistema orbital (VLHL) para o sistema do
sátelite. Os ângulos φ, θ, ψ são tomados no sentido anti-horário. Para rotações infinitesimais a seqüência de rotações não é importante (Wie, 1998, Moscati, 1992).
A representação em quaternions, para pequenos ângulos pode ser aproximado por
η∼
=1
φ
1 ∼
=
2
θ
2 ∼
=
2
ψ
3 ∼
=
2
60
(4.24)
(4.25)
(4.26)
(4.27)
sendo necessário a normalização dos quaternions para satisfazer a relação
η 2 + 21 + 22 + 23 = 1
4.3
(4.28)
Sumário de Notações
As principais notações usadas nesse trabalho são:
Velocidade angular do referencial b em relação ao referencial e exω neb
presso/escrito no referencial n
ra
Vetor r expresso no referencial a
a
Matriz de rotação do referencial b para o referencial a. ra = Rab rb
Rb
S (k) Operador anti-simétrico
r̂ Versor ou vetor unitário
61
62
CAPÍTULO 5
FORMULAÇÃO DO PROBLEMA
Neste Capı́tulo o problema objeto desta dissertação é formulado. A formulação aqui
apresentada envolve:
• os modelos matemáticos da cinemática e da dinâmica de um satélite contendo
3 rodas de reação e 3 bobinas magnéticas;
• a modelagem do gradiente de gravidade e dos torques magnéticos;
• a modelagem do campo geomagnético, modelo IGRF;
• o tratamento das equações diferenciais, incluindo a linearização dos modelos
matemáticos, para serem utilizadas no projeto de controle linear nos modos
de controle de estabilização e aquisição;
• uma discussão do emprego de atuadores usados neste trabalho com enfoque
nas suas aplicações, factibilidade, vantagens e desvantagens;
• uma discussão dos torques ambientais incluindo o modelo de pior caso.
5.1
Modelagem Matemática: Cinemática e Dinâmica
Esta sessão apresenta e discute os modelos matemáticos da cinemática e da dinâmica
de um satélite equipado com três rodas de reação e bobinas magnéticas. As equações
da cinemática são obtidas a partir da representação da atitude em ângulos de Euler e
quaternions. Essas parametrizações são amplamente utilizadas na área de dinâmica
de atitude uma vez que permitem representar a dinâmica e a atitude dos veı́culos
espaciais em diferentes sistemas de coordenadas. As equações da cinemática não
envolvem as forças ou torques associados ao movimento. Elas descrevem a velocidade
do veı́culo espacial em termos de sua orientação em relação a um ou mais sistemas
de coordenadas. A orientação do veı́culo em relação a um referencial conhecido
é chamada atitude. A atitude pode ser representada por diferentes conjuntos de
parâmetros, como é mostrado no Capı́tulo (4).
63
5.1.1
Equações da Cinemática
A cinemática descreve a orientação de um veı́culo em relação a um sistema de eixos
conhecido. O modelo matemático da cinemática é escrito na forma de um sistema de
equações diferenciais de primeira ordem que, em conjunção com as equações de Euler,
permite escrever as equações do movimento na forma de estado e, por integração
numérica, obter os ângulos de atitude e sua derivada temporal em função do tempo.
Em geral as equações da cinemática podem ser representadas por quaternions, matriz
de rotação, ângulos de Euler, etc. Cada uma das formas de representar a atitude tem
vantagens e desvantagens. A representação por quaternions tem a vantagem de não
conter funções trigonométricas e nem conter singularidades, que são tı́picas quando
se representa a cinemática em ângulos de Euler. Wertz (1978), Kaplan (1976) e
Hughes (1986) apresentam as equações da cinemática em quaternions na forma
1 T b
ω ob
2
1 b
1
˙ =
ηω ob − ω bob × 2
2
η̇ =
(5.1)
(5.2)
Usando o operador anti-simétrico S (ω), a Equação (5.2) pode ser reescrita como
˙ =
1
[ηI + S ()] ω bob
2
(5.3)
Os parâmetros de Euler ou quaternions dados nas equações (5.2) são apresentados
na Seção (4.2). ω bob é o vetor velocidade angular do satélite escrito no referencial
do corpo, referido ao sistema orbital. As equações da cinemática podem também
ser escritas em função dos ângulos de Euler. A forma do sistema de equações da
cinemática em termos dos ângulos de Euler depende da seqüência de rotações escolhidas para se passar de um sistema de referencia para outro. A forma mostrada a
seguir pode ser vista em Wertz (1978) e Wie e Arapostathis (1989) e é obtida para
a seqüência de rotações 3 − 2 − 1.
64
 




φ̇
cos θ sin φ sin θ cos φ sin θ
sin ψ
1 
ωo 
 


cos φ cos θ − sin φ cos θ ω bib +
 θ̇  =
 0
cos θ cos ψ 
cos θ
cos θ
ψ̇
0
sin φ
cos φ
sin θ sin ψ
(5.4)
onde ωo é a velocidade orbital, que para uma órbita circular é constante.
A Figura (5.1) mostra a seqüencia de rotações associadas à Equação (5.4) para
representar o sistema do corpo em relação ao sistema orbital, através de rotações
subseqüentes dos ângulos ψ, θ e φ.
Figura 5.1- Seqüência de rotações 3(ψ) − 2(θ) − 1(φ) .
Os eixos x, y e z e os demais sistemas de eixos formados pelas rotações são ortogonais.
Entretanto, as rotações em torno dos eixos que são escritas as velocidades angulares
não são ortogonais . Por esta razão tem-se seno ou co-seno no denominador (dependendo da seqüência de rotações escolhida) quando se representa o vetor [φ̇ θ̇ ψ̇]T
em função dos respectivos ângulos e das componentes do vetor velocidade angular ω bob . Tal caracterı́stica explica a singularidade sempre presente na representação
da atitude via ângulos de Euler, uma vez que para cos ±n π2 , n = 1, 3, 5, · · · e
sin (±nπ), n = 0, 1, 2, · · · tem-se denominador nulo. A Equação (5.4) ilustra o caso
do co-seno.
65
5.1.2
Equações da Dinâmica
O modelo real do satélite é aproximado por um modelo fı́sico constituı́do de um
corpo rı́gido principal (main bus) contendo internamente 3 rodas de reação, como
mostrado na Figura (5.2), também representadas por corpos rı́gidos. O sistema,
entretanto, não constitui de fato um corpo rı́gido uma vez que contém partes moveis
internas, representadas pelas rodas de reação. O modelo matemático da dinâmica é
descrito pelas equações de Euler, expressas no sistema BF (x, y, z) do corpo, referida
ao sistema inercial. A origem do sistema de eixos do veı́culo é definida no centro de
massa, CM , do satélite. Os eixos x, y e z são coincidentes com os eixos principais de
inércia da espaçonave. A Figura ilustra o conceito do satélite objeto da modelagem
matemática bem como os sistemas de referências orbital OF(xo , yo , zo ) e do corpo
BF (x, y, z). O referencial orbital foi apresentado e comentado na Seção (4.1.2).
As equações de Euler utilizadas aqui para representar a dinâmica do satélite, são
obtidas a partir da derivada temporal da quantidade de movimento angular, como
é mostrado a seguir.
Figura 5.2- Ilustração do satélite com rodas de reação, bobinas e os referenciais
orbital OF (xo , yo , zo ) e do corpo BF (x, y, z).
Seja Lb a quantidade de movimento angular do satélite, expressa no sistema do
66
corpo. O torque externo τext é dado pela taxa de variação no tempo da quantidade
de movimento angular. Matematicamente:
b
dL
dLb
=
+ ω bib × Lb = τ ext
dt
dt b
(5.5)
onde
dLb
=
dt
L̇b
= Js ω˙bib + Jw ω˙bib + ω˙bbw
(5.6)
b
onde Js é a matriz de inércia do satélite (sem considerar as rodas), Jw é a matriz de
inércia das rodas, ω bbw é a velocidade angular das rodas relativa ao satélite e ω bib é
a velocidade angular do satélite em relação ao sistema inercial expresso no sistema
do satélite.
Lb pode ser escrito com a soma dos vetores quantidades de movimento angular do
corpo principal do veı́culo (main bus) mais o vetor quantidade de movimento angular
das rodas de reação. Matematicamente
Lb = hs + lw
(5.7)
temos que hs = Js ω bib é a quantidade de movimento angular do satélite e lw =
Jw ω bib + ω bbw é a quantidade de movimento angular das rodas. Substituindo na
Equação (5.5)
Escrevendo Lb dessa forma pode-se reescrever a Equação (5.5) na forma:
dLb
b
= Js ω̇ib
+ Jw ω̇ bib + ω̇ bbw + ω bib × Js ω bib + Jw ω bib + ω bbw = τ ext
dt
67
(5.8)
A matriz de inércia Js do corpo principal do satélite (main bus) é dada por


Jx −Jxy −Jxz


Js =  −Jxy Jy −Jyz 
−Jxz −Jyz
Jz
Z
Jx =
2
y +z
B
Jxy
Z
Z
x2 + y 2 dm
dm, Jy =
x + z dm, Jz =
B
B
Z
Z
Z
=
(xy) dm, Jxz =
(xz) dm, Jyz =
(yz) dm
2
2
B
2
B
(5.9)
(5.10)
B
Como considerou-se o sistema BF (x, y, z) coincidente com os eixos principais de
inércia, essas matrizes são diagonais e os elementos Jx , Jy e Jz são os momentos
principais de inércia. Matematicamente:


Jx 0 0


Js =  0 J y 0 
0 0 Jz
(5.11)
Jw é matriz de inércia associada às três rodas de reação. Como cada rotor tem seu
eixo principal de inércia alinhado com cada um dos eixos principais de inércia do
veı́culo, a matriz de inércia Jw , associada aos três rotores, é também diagonal e
os elementos da diagonal são os momentos de inércia axiais associados ao eixo de
rotação de cada roda de reação. Matematicamente:


Jwx 0
0


Jw =  0 Jwy 0 
0
0 Jwz
(5.12)
Reescrevendo a Equação (5.8) em termos das velocidade angulares e sua taxa de
variação no tempo na forma:
68
dLb
b
= (Js + Jw ) ω̇ib
+ Jw ω̇ bbw + ω bib × (Js + Jw ) ω bib + Jw ω bbw ) = τ ext
dt
(5.13)
A matriz de inércia total do veı́culo é J = Js + Jw . Matematicamente


Jx + Jwx
0
0


J = Js + Jw = 
0
Jy + Jwy
0

0
0
Jz + Jwz
(5.14)
Substituindo a Equação (5.14) na Equação (5.13) e usando o operador anti-simétrico,
reescrevemos a Equação (5.13)
b
Jω̇ib
+ Jw ω̇ bbw + S ω bib
Jω bib + Jw ω bbw = τ ext
(5.15)
Pode-se identificar a expressão correspondente ao torque de controle gerado pelas
rodas de reação:
Jw ω̇ bbw + S(ω bib )Jw ω bbw = τ bw
(5.16)
ou seja, o torque das rodas de reação aplicados no satélite para fins de controle e
estabilização. Usando esta definição pode-se reescrever a expressão do torque (5.15)
como:
b
Jω̇ib
+ S ω bib Jω bib = τ ext − τ bw
(5.17)
O torque externo τ ext representa o somatório de todos os torques externos que
atuam no veı́culo espacial. Neste trabalho considerou-se apenas o torque devido ao
gradiente de gravidade e o torque de controle gerado pela interação das bobinas
magnéticas com o campo magnético da Terra. Considerou-se 3 bobinas, dispostas
69
segundo os eixos x, y e z, como mostrado na Figura (5.3). O torque externo pode,
então, ser escrito como na forma:
Figura 5.3- Configuração das bobinas magnéticas no satélite.
τ ext = τ bd + τ bm
(5.18)
onde τ bd representa todos os torques de perturbação externa (gradiente de gravidade,
pressão de radiação, arrasto atmosféricos, etc) e τ bm representa o torque gerado pelas
bobinas, dado por (Fauske, 2002)
τ bm = mb × Bb
(5.19)
onde mb é o momento dipolo magnético gerado pelas bobinas, cujas caracterı́sticas
depende de projeto (magnitude e direção do torque) e pode ser representado na
forma

  
N x i x Ax
mx

  
b
b
b
b
m = mx + my + mz =  Ny iy Ay  = my  ,
N z iz Az
mz
(5.20)
Nessa Equação x, y e z são os eixos do sistema BF (x, y, z) e N , i e A são o número
de espiras, a corrente elétrica e a área da espira (parâmetros de projeto), respectivamente, para as três bobinas, dispostas segundo os eixos x, y e z.
b
B =
h
Bxb
Byb
Bzb
iT
é o vetor campo magnético local, escrito no referencial do
70
satélite.
Portanto o modelo matemático que representa o satélite equipado com bobinas é
dado por:
b
Jω̇ib
+ S ω bib Jω bib = τ bd + τ bm
5.1.3
(5.21)
Torque devido ao Gradiente de Gravidade
O torque de gradiente de gravidade é causado devido a variação da força gravitacional ao longo do corpo do satélite. O torque depende da altitude, da geometria e da
atitude do satélite (Carrara, 1982). Esse efeito é usado para estabilização passiva de
satélites O gradiente de gravidade alinha o eixo de menor momento de inércia com
a vertical local (Fauske, 2002 e Arantes e Fonseca, 2004a). Esse fato requer que um
dos momentos principais de inércia do satélite seja muito menor dos que os outros.
Em geral essa configuração é obtida através de mastros extensı́veis colocados no
satélite. Uma vez em órbita o mastro é estendido segundo o eixo conveniente a ser
alinhado com a vertical local. Exemplos de satélites que usam a estabilização por
gradiente de gravidade são DODGE, GEOS-3 e RAE-2.
Como deduzido no apêndice A, o torque devido ao gradiente de gravidade, expresso
no referencial do satélite, é dado por
τ bg = 3ωo2 cb3 × Jcb3
(5.22)
onde cb3 é o vetor coluna com os cosenos diretores da vertical local ou eixo zo nas
direções dos eixos x, y e z e ωo é a velocidade orbital.
5.2
Campo Magnético Terrestre
O campo magnético terrestre é calculado como o gradiente de uma função potencial
escalar, que é modelada por um série de harmônicos esféricos, dada por (Fauske,
71
2002)
V (r, θ, φ) = a
n
∞ X
a n+1 X
n=1
r
m
[gnm cos mφ + hm
n sin mφ] Pn (θ)
(5.23)
m=0
A função potencial V é uma função da posição do satélite em coordenadas esfêricas
n
no sistema ECI definido em 4.1.1. Os coeficientes gm
e hnm são chamados coeficientes
de Gauss, que são constantes e obtidos empiricamente. A cada cinco anos a IAGA
(Association of Geomagnetism and Aeronomy) publica um novo conjunto de coeficientes que constituem o modelo IGRF (International Geomagnetic Reference Field ).
Discussões e detalhes desse modelo são encontrados em Macmillan e Quinn (2000).
O modelo (código fonte nas linguagens C e FORTRAN) é fornecido em Mclean, .
Nesse trabalho a implementação do modelo IGRF é feita no ambiente MATLAB e
SIMULINK. A Figura (5.4) mostra as componentes do vetor campo magnético local
T
Bo = Bxo Byo Bzo , expresso no referencial da órbita (OF) ou VLHL, as componentes do campo magnético são dadas para três revoluções do satélite em torno da
Terra.
Figura 5.4- Campo magnético local Bo usando o modelo IGRF 2000.
A propagação da órbita é descrita no apêndice B. Foi usada uma órbita circular kepleriana. A órbita propagada para o cálculo do campo magnético é a órbita nominal
72
do satélite EQUARS. Os elementos orbitais são mostrados na Tabela (5.1).
Tabela 5.1- Elementos orbitais do satélite EQUARS.
altura (h)
5.3
excentricidade (e)
750 km
∼
=0
ascenção do nodo ascendente (Ω)
30 graus
inclinação (i)
20 graus
Tratamento das Equações da Dinâmica para os Casos Estudados nesse
Trabalho
O controle de atitude de satélites artificiais tem inı́cio na separação foguete-satélite.
Até então, qualquer que seja o controle, sua atuação é sobre o veı́culo lançador.
Por exemplo, o lançador pode ter um sistema de controle que reduza a rotação do
último estágio, para nı́veis de projeto, imediatamente antes da separação. Conhecida a atitude final do último estágio do foguete, tem-se uma estimativa razoável
da atitude do satélite imediatamente após a separação. Por outro lado, o último
estágio do foguete pode garantir apenas a inserção do veı́culo na órbita correta, não
fornecendo dados sobre a atitude do satélite imediatamente após a separação. Nesse
caso, o sistema de controle do satélite deve ser capaz de controlar o movimento de
atitude do mesmo a partir de qualquer configuração de atitude imediatamente após
a separação. Portanto, se o lançador não fornecer, a priori, uma atitude final que
possa ser estendida ao satélite imediatamente após a separação, o sitema de controle
de atitude deverá, de acordo com requisitos de missões:
• Reduzir para valores especificados ou eliminar, de acordo também com especificações de projeto, a rotação do satélite imediatamente após a separação.
Esta fase constitui um dos modos de estabilização, comumente conhecido como
detumble, traduzido aqui como descapotamento. É a operação de eliminar o
movimento de rotação arbitrária do satélite imediatamente após a separação.
Este aspecto é tratado nesta Seção porque o controle para este modo pode requerer manobras de grandes ângulos e, neste caso, a linearização do sistema de
equações diferenciais e as técnicas de controle linear não se aplicam. Portanto
neste modo de operação utiliza-se técnicas de controle não linear;
73
• Uma vez que o sistema de controle de atitude faça o detumbling e tenha controle sobre a velocidade angular e atitude do veı́culo, o satélite entra no modo
de estabilização. Uma vez reduzida a velocidade angular do satélite imediatamente após a separação, o satélite deve ser manobrado para adquirir a
atitude nominal. Nesse caso, a rotação residual pode sobrecarregar o sistema
de controle e, portanto, é desejável que, no modo anterior, a velocidade angular
residual não seja crı́tica, ou seja, não comprometa a eficiência do sistema de
controle. A fase de estabilização pode requerer, também, manobras de grandes
ângulos, uma vez que, no modo detumble, o objetivo é a redução da velocidade
angular e não manobras angulares para levar o veı́culo para uma configuração
nominal em atitude. Novamente, neste caso, deve se ter um certo conhecimento da atitude do satélite na fase final do detumbling para se decidir sobre
a aplicação de leis de controle linear ou não linear.
• Uma vez adquirida a condição normal de operação o sistema de controle deve
manter o veı́culo estável em torno da atitude nominal, de acordo com especificações de projeto para atender as missões. Nesse trabalho, esse processo
também faz parte do modo de estabilização. Nesse modo o controle deve fazer
sempre manobras de pequenos ângulos (correções) visando manter o veı́culo
na atitude nominal especificada. Neste modo de operação as equações do movimento podem ser linearizadas em torno da atitude nominal e as técnicas de
controle linear se aplicam. Neste estudo as técnicas de controle linear aplicadas
são as do PD, LQR e LQG.
Eventualmente, dependendo da missão, pode-se requerer a mudança de uma atitude
nominal para outra, arbitrária, por um certo intervalo de tempo. Nesse trabalho esse
modo de contigência é chamado modo de aquisição. Quando isto ocorre novamente
requer-se manobras de grandes ângulos e posterior reaquisição da nova atitude nominal.
A visão do sistema de controle em modos de operação é mais prática sob o ponto
de vista de que o projeto de sistemas de controle reais é implementado por modos
de operação. São necessárias diferentes estratégias e técnicas de controle para cada
modo.
74
5.4
Linearização
Para o projeto do controlador linear usando as metodologias LQR e LQG, que são
baseadas em uma planta linear, é feita a linearização das equações do movimento
em torno de um ponto de operação. As equações linearizadas também podem ser
usadas para estudos de estabilização/manutenção da atitude nominal, para os casos
onde são realizadas manobras de pequenos ângulos (< 15o ). As equações usadas são
escritas em termos dos ângulos de Euler, já que a sigularidade do tipo 1 (±π/2),
não acontece para a faixa de manobras realizadas (pequenos ângulos). Outra justificativa é a clara interpretação fı́sica que os ângulos de Euler fornecem. São feitas as
linearizações para as versões:
• Satélite equipado com rodas, Equacão (5.17);
• Satélite equipado com bobinas, Equação (5.21).
A Equação da cinemática (5.4) pode ser escrita como (Bryson, 1994)



 
cos θ sin ψ
1
0
− sin θ
φ̇



 
ω bib = 0 cos φ sin φ cos θ   θ̇  − ωo sin φ sin θ sin ψ + cos φ cos ψ 
cos φ sin θ sin ψ − sin φ cos ψ
ψ̇
0 − sin φ cos φ cos θ
(5.24)
Para pequenos ângulos a expressão (5.24) pode ser aproximada por (Bryson, 1994)
ω1 ≡ φ̇ − ωo ψ − ψ̇θ
ω2 ≡ θ̇ − ωo − ψ̇φ
ω3 ≡ ψ̇ + ωo φ − θ̇φ
(5.25)
O valor médio, por órbita, dos termos não lineares ψ̇θ, ψ̇φ, θ̇φ são pequenos comparados com as magnitude média dos termos lineares em (5.25). Se φ e θ oscilam em
torno de zero ou se a magnitude de ψ̇ e θ̇ são pequenos comparados com a velocidade
orbital ωo . Nesses casos a Equação (5.25) pode ser aproximada por
75
φ̇ ≡ ω1 + ωo ψ
θ̇ ≡ ω2 + ωo
ψ̇ ≡ ω3 − ωo φ
Usando a notação ω bib = [ω1
ω2
(5.26)
ω3 ]T temos
 
 
φ̇
−ψ
 
 
b
ω ib = I  θ̇  + ωo  −1 
ψ̇
φ
(5.27)
onde I é a matriz identidade de ordem 3.
5.4.1
Linearização do Modelo do Satélite Equipado com Rodas
Definimos o vetor de variáveis de estado como
h
iT
x = φ θ ψ φ̇ θ̇ ψ̇
(5.28)
Usando a notação hw = [hw1 hw2 hw3 ]T onde hwi , (i = 1, 2, 3) são as componentes, ao longo das direções (x, y, z), da quantidade de movimento angular total das
rodas. Para a configuração adotada para as rodas, cada componente da quantidade
de movimento angular corresponde a uma roda. Substituindo a Equação (5.27) na
Equação (5.17) e incluindo o torque de gradiente de gravidade dado pela Equação
(5.22), obtemos
Jx φ̈ = [4ωo2 (Jz − Jy ) + ωo hw2 ] φ + ωo h3 + [ωo (Jx − Jy + Jz ) + hw2 ] ψ̇
+ [ωo (Jy − Jz ) φ + hw3 + (Jz − Jy )] θ̇ − ḣw1
76
(5.29)
Jy θ̈ = 3ωo2 (Jz − Jx ) θ h− ωo φh1 − [ωo2 (Jz − Jx ) φ + ωo hw3 ]iψ − [ωo (Jz − Jx ) ψ + hw1 ] ψ̇
− ωo (Jx − Jz ) φ − hw3 + (Jx − Jz ) ψ̇ φ̇ − ḣw2
(5.30)
Jz ψ̈ = [ωo2 (Jxh− Jy ) + ωo hw2 ] ψ − ωo hw1 − [ωo (Jxi− Jy + Iz ) + hw2 ] φ̇
+ ωo (Jy − Jx ) ψ + hw1 + (Jx − Jy ) φ̇ θ̇ − ḣw3
(5.31)
onde Jx , Jy e Jz são os momentos principais de inércia.
Transformando as equações (5.29), (5.30) e (5.31), de segunda ordem em equações
de primeira ordem e fazendo a expansão em série de Taylor em (x = 0), ou seja, a
vertical e horizontal local (VLHL), obtemos
x∗ =
h∗w
=
0
h
0 Ho
iT
0
(5.32)
Onde o termo Ho é a quantidade de movimento angular na direção do eixo de pitch
gerado pelo volante de inércia. Esse procedimento é utilizado para produzir rı́gidez
giroscópica durante a fase de manutenção da atitude.
Escrevendo na forma de estados as equações linearizadas, temos
ẋ = Ax + Bu
com
77
(5.33)

0
0
0





A =  4ωo2 (Jz −Jy )−ωo Ho

Jx


0

0
0
0
0
0
3ωo2 (Jz −Jx )
Jy
0
0
0
0
0
0
1
0
0
0
0
ωo2 (Jx −Jy )−ωo Ho
Jz
ωo (−Jx +Jy −Jz )+Ho
Jz
0
1
0
0
0
0
0
1





ωo (Jx −Jy +Jz )−Ho 

Jx


0

0
0
(5.34)
e

0
0

 0
0

 0
0

B= 1
− J
0
 x
 0 −1

Jy
0
0
0
0
0
0
0
− J1z











(5.35)
O vetor de controle é dado por
iT
h
u = u1 u 2 u3
(5.36)
u1 = ḣw1 − ωo hw3
u2 =
ḣw2
u1 = ḣw3 + ωo hw1
(5.37)
onde
5.4.2
Linearização do Modelo do Satélite Equipado com Bobinas
Definimos o vetor de variáveis de estado e o vetor de controle como
78

h
x = φ θ ψ φ̇ θ̇ ψ̇
iT
(5.38)
e
h
u = mx my mz
iT
(5.39)
Usando a notação Bb = Bxb Byb Bzb , onde Bib ,(i = x, y, z) são as componentes do
vetor campo magnético local expresso no sistema do satélite. Substituindo a Equação
(5.27) na Equação (5.21) e incluindo o torque de gradiente de gravidade dado pela
Equação (5.22), obtemos
Jx φ̈ =
h
4ωo2
i
(Jz − Jy ) − ωo (Jz − Jy ) θ̇ φ + ωo (Jx − Jy + Jz ) ψ̇ + (Jy − Jz ) θ̇ψ̇
+my Bzb (t) − mz Byb (t)
(5.40)
i
h
Jy θ̈ = 3ωo2 (Jx − Jz ) θ − ωo2 (Jz − Jx ) ψ + ωo (Jx − Jz ) φ̇ φ + ωo (Jx − Jz ) ψ ψ̇
+ (Jz − Jx ) φ̇ψ̇ + mz Bxb (t) − mx Bzb (t)
(5.41)
h
i
h
i
Jy ψ̈ = ωo2 (Jx − Jy ) + ωo (Jy − Jx ) θ̇ ψ − ωo (Jx − Jy + Jz ) + (Jy − Jx ) θ̇ φ̇
+mx Byb (t) − my Bxb (t)
(5.42)
Transformando as equações (5.40), (5.41) e (5.42), de segunda ordem, em equações
de primeira ordem e fazendo a expansão em série de Taylor em torno de (x = 0), ou
seja, a vertical e horizontal local (VLHL), obtemos
ẋ = Ax + B (t) u
79
(5.43)
com

0
0
0





A =  4ωo2 (Jz −Jy )

Jx


0

0
0
0
0
0
3ωo2 (Jz −Jx )
Jy
0
0
0
0
0
0
1
0
0
0
0
0
1
0
0
0
ωo2 (Jx −Jy )
Jz
ωo (−Jx +Jy −Jz )
Jz
0
0
1






ωo (Jx −Jy +Jz ) 

Jx


0

0
0
(5.44)
e

0
0
0





B=
 0
 b
− Bz (t)
 Jy
Byb (t)
Jz
0
0
0
Bzb (t)
Jx
0
B b (t)
− Jyz
0
0
0






Byb (t) 
− Jx 

Bxb (t) 

Jy
0
(5.45)
O sistema linear dado pela Equação (5.43) é variante no tempo (LTV), devido
ao campo magnético terrestre não ser constante. O campo é aproximadamente
periódico. No sistema linearizado, o vetor campo magnético expresso no referencial do satélite (Bb ) coincide com o vetor campo magnético expresso no referencial
orbital (Bo ). Portanto a matriz B pode ser escrita como

0
0
0





B=
 0
 o
− Bz (t)
 Jy
Byo (t)
Jz
0
0
0
Bzo (t)
Jx
0
B o (t)
− Jyz
80
0
0
0






o
By (t) 
− Jx 

Bxo (t) 

Jy
0
(5.46)
As equações (5.43) e (5.33) representam o modelo matemático linearizado do satélite
equipado com rodas (de reação e volante de inércia) e bobinas magnéticas, que
constituem uma boa aproximação do modelo para manobras de pequenos ângulos
(< 15o ). Os projetos do controladores multivariável LQR e LQG foram feitos sobre o
modelo linear dado pelas equações (5.43) e (5.33). Entretanto, a simulação digital é
feita sobre o modelo completo, não linear, equações (5.4), (5.21) e (5.17). O sistema
de controle no modo normal de operação (que é a fase de manutenção da atitude
nominal), constitue um problema de regulador, ou seja, a referência é constante e
igual a zero (x = 0). O que corresponde ao satélite estar alinhado com a VLHL.
O modelo linear do satélite equipado com bobinas, a atitude é descrita em quaternions e é dado em Arantes e Fonseca (2004b) e Overby (2004). No modelo do satélite
parametrizado em quaternions o estado é dado por
h
x = 1 2 3 ˙1 ˙2 ˙3
iT
(5.47)
Escrevendo o modelo na forma da Equação (5.43), as matrizes A e B são dadas por

0
0
0





A =  4ωo2 (Jz −Jy )

Jx


0

0
0
0
0
0
3ωo2 (Jz −Jx )
Jy
0
0
0
0
0
0
1
0
0
0
0
ωo2 (Jx −Jy )
Jz
ωo (−Jx +Jy −Jz )
Jz
e
81
0
1
0
0
0
0
0
1






ωo (Jx −Jy +Jz ) 

Jx


0

0
0
(5.48)

0
0
0





B=
 0
 o
− Bz (t)
 2Jy
Byo (t)
2Jz
0
0
0
Bzo (t)
2Jx
0
Byo (t)
− 2J
z
0
0
0






Byo (t) 
− 2Jx 

Bxo (t) 
2Jy 
0
(5.49)
Nota: As linearizações das equações da dinâmica foram feitas com o auxı́lio do
manipulador simbólico do programa MATHEMATICA e as equações obtidas foram
aferidas com as encontradas na literatura.
A próxima Seção discute algumas considerações sobre os modelos matemáticos e
sistemas de referencias associados, e os atuadores escolhidos para compor os sistemas
de controle nesse trabalho.
5.5
Considerações Sobre os Atuadores, os Modelos Matemáticos e Sistemas
de Referênica Adotados
As considerações sobre os atuadores são mostradas a seguir.
5.5.1
Modelo dos Atuadores
Existem vários tipos de atuadores que podem ser usados para o controle de atitude
de satélites artificiais. Os atuadores podem ser dividos em três categorias: 1) propulsores; e 2) dispositivos de troca de quantidade de movimento angular, como rodas
e 3) atuadores magnéticos. É comum o emprego de mais de um tipo de atuador
em satélites, dependendo dos requisitos da missão, e das caracterı́sticas do projeto
de controle, Wertz (1978) e Larson e Wertz (1992) apresentam uma descrição do
emprego de atuadores, empregados em diferentes missões.
Roda de Reação e Volante de Inércia
As rodas de reação e volantes de inércia podem ser modelados como motores DC,
essa modelagem pode ser encontrada em Bryson (1994), Wilson (2000) e Trivelato
82
(1988). Uma aproximação para a função de transferência desses atuadores é dada
em Souza (1981). Define-se (Larson e Wertz, 1992)
• Rodas de reação: rodas com quantidade de movimento angular nominal nulo
• Volante de inércia: são rodas de reação com quantidade de movimento
angular nominal diferente de zero
Satélites equipados com volantes de inércia são referidos como satélites com quantidade de movimento angular embarcado (momentum bias system), e satélites equipados com rodas de reação são referidos como satélites com quantidade de movimento
angular nulo (zero momentum system). Basicamente, o funcionamento consiste na
geração de torques devido a aceleração de uma roda, ligada ao rotor de um motor
elétrico, em relação a seu estator que é fixado à estrutura do satélite (Souza, 1987).
A dinâmica da roda é função do seu atrito, sua velocidade, e voltagem aplicada pelo
motor. Considera-se nesse trabalho que exista um sistema de controle fino para a
operação das rodas, que seja eficiente, fornecendo o torque requisitado pelo sistema
de controle de atitude. Um estudo detalhado de um sistema de controle por referência
para operação de rodas pode ser encontrado em Trivelato (1988). Em Souza (1987)
é apresentado o estudo e o projeto de um modelo experimental de uma roda de
reação.
O uso de volantes de inércia, que apresentam quantidade de movimento angular
nominal diferente de zero, fornece rı́gidez giroscópica ao satélite. Para as rodas de
reação espera-se que operem em torno de velocidade (relativa) nula ao longo de um
perı́odo orbital, adequada para absorver torques externos cı́clicos. A presença de
torques seculares, pode levar à velocidade angular máxima da roda, exigindo a sua
dessaturação. Isso pode ser feito por bobinas magnéticas e/ou propulsores.
Nesse trabalho considerou-se somente rodas de reação combinadas com atuadores
magnéticos para compor o sistema de controle. Propulsores não foram considerados.
83
Bobinas Magnéticas
As bobinas magnéticas podem ser modeladas como um circuito RL (resistência e
indutância) usando as leis de Kirchoff. As bobinas tem uma dinâmica muito rápida
comparada com a dinâmica do satélite. Incluindo a dinâmica das bobinas na simulação, o sistema torna-se muito lento, pois envolve duas escalas de tempo distintas, portanto, a dinâmica das bobinas é ignorada na simulação. O erro indroduzido
devido a isso é negligenciável (Fauske, 2002). Logo, a dinâmica total do sistema é
relativamente lenta.
Bobinas magnéticas são comumente usadas em aplicações espaciais desde que os
requisitos de tempo não sejam muito exigentes. As manobras de grandes ângulos via
bobinas magnéticas podem levar horas, dependendo das especificações da bobina.
As bobinas são também muito usadas para a dessaturação de rodas de reação. A
saturação das rodas ocorre quando elas atingem sua quantidade máxima de quantidade de movimento especificada. Quando isto ocorre não se consegue mais acelerar
as rodas e, portanto, não se consegue mais usá-las para finalidade de torque objetivando o controle da atitude. As bobinas podem também ser usadas para manobras de
pequenos ângulos. Uma aplicação muito comum para manobras de grandes ângulos
se refere as operações do controle no modo detumble.
O estudo de dessaturação das rodas pode ser encontrado nas referências Bryson
(1994) e Bang e Choi (2003). Para manobras de pequenos e grandes ângulos, ou seja,
estabilização e aquisição, respectivamente, as bobinas são, em geral, usadas para
pequenos satélites (< 500kg) (Wisniewski e Blanke, 1999,Wisniewski e Markley,
1999), como é o caso do satélite brasileiro EQUARS (∼
= 150kg). O controle de
pequenos satélites usando somente bobinas magnéticas só é possı́vel para órbitas
em que a variação do campo magnético seja suficiente para garantir o controle do
satélite (Wisniewski, 1997, Musser e Ebert, 1989).
O uso das bobinas magnéticas para operações em órbita faz uso da interação entre
o momento dipolo magnético gerado pelas bobinas e do campo magnético terrestre. Por isso, para efeitos de simulações para avaliação do desempenho dos torques
magnéticos é muito importante um bom modelo do campo magnético terrestre.
Em termos de projeto é sempre possı́vel, dentro das limitações de espaço, massa e
energia, alterar o momento das bobinas magnéticas, via parâmetros de projeto tais
84
como número de espiras, corrente disponı́vel e área das bobinas. Entretanto não se
tem controle do campo geomagnético, que é parte do ambiente espacial. O modelo
do campo geomagnético e sua eficiência para os sistemas de controle com bobinas
magnéticas depende dos parâmetros orbitais. Portanto, o estudo do controle via bobinas exige conhecimento do modelo do campo bem como do modelo da órbita do
satélite em questão.
Esse Capı́tulo é concluı́do com a apresentação dos modelos de pior caso para os
principais torques ambientais e um modelo de torques de perturbação internos.
5.6
Torques Ambientais
Um satélite está sujeito a pequenos mas a persistentes torques externos devido a
perturbações ambientais. Sem resistência a essas perturbações o satélite pode perder
a atitude nominal. Os principais torques ambientais são:
• torque devido ao gradiente de gravidade;
• torque aerodinâmico;
• torque devido a pressão de radiação (solar e terrestre);
• torque devido a interação do dipolo residual com o campo magnético.
Essa Seção discute brevemente esses torques ambientais, e apresenta o modelo de pior
caso, conveniente para dimensionalização dos atuadores (torque de controle), para
satélites em geral. Uma discussão pormenorizada é encontrada em Carrara (1982),
Wertz (1978), Hughes (1986) e Fauske (2002). Nesse trabalho somente o torque
devido ao gradiente de gravidade é considerado. Esse torque é considerado, nesse
estudo, como o mais importante em termos de perturbação, desde que se considere
um bom balanceamento magnético. Os outros torques não são considerados.
5.6.1
Torque Devido ao Gradiente de Gravidade
Este torque depende da altitude do veı́culo, das suas propriedades de inércia e de sua
atitude em relação à vertical local. Uma expressão que leva em conta estes aspectos e
85
que é apropriada para se calcular em primeira aproximação a magnitude do torque
do gradiente de gravidade para o pior caso é dada em Larson e Wertz (1992) e
reproduzida aqui
τgb =
3ωo2
|Jx − Jz | sin 2θ
2
(5.50)
onde θ é o ângulo que descreve o deslocamento angular na direção da vertical local
ou nadir. O gradiente de gravidade tem a propriedade de alinhar o eixo de menor
momento de inércia com a vertical local, em uma configuração de estabilidade chamada estabilização por gradiente de gravidade. Por esta razão satélites estabilizados
por gradiente de gravidade requerem o projeto estrutural de tal forma a permitir
que o eixo de apontamento para a Terra seja o eixo de menor momento de inércia.
Um exemplo interessante refere-se à primeira versão do satélite brasileiro de coleta de dados, SCD-1, que deveria ser estabilizado por gradiente de gravidade. O
veı́culo deveria conter um mastro a ser destendido em órbita. O mastro teria um
comprimento de 10 m, com uma massa de 3 kg na ponta. Isto aumentaria de aproximadamente dez vezes os momentos de inércia transversais do satélite enquanto
manteria o momento de inércia axial (eixo de apontamento para a Terra) aproximadamente igual. Posteriormente a configuração de estabilização foi reconsiderada
e o projeto modificado para a estabilização por spin (direção de ω coincidente com
eixo axial do veı́culo). Esta configuração de estabilização requer apenas que o eixo
de rotação do veı́culo seja o de maior momento de inércia. A mudança diminuiu
os riscos associados as operações de abertura e captura do satélite (SCD-1) em
torno da vertical local. Normalmente os satélites estabilizados por gradiente de gravidade requerem um mastro a ser aberto em órbita para adequar as propriedades
de inércia do satélite a estabilização por gradiente de gravidade. Portanto a complicada dinâmica na fase transitória da abertura dos mastros representam um aspecto
delicado para tal tipo de estabilização. Como existem restrição de espaço e volume
nos veı́culos lançadores não há como lançar um satélite estabilizado por gradiente
de gravidade já com a configuração adequada em termos de propriedades de inércia.
A geometria em forma de lápis (pencil-shaped ) que caracteriza os satélites estabilizados por gradiente de gravidade é, de fato, inadequada para os veı́culos lançadores.
Daı́ a necessidade de mastros com massa adicional na ponta para dar ao satélite a
forma e as propriedades de inércia apropriadas para a estabilização por gradiente de
86
gravidade.
Neste trabalho considerou-se apenas o torque do gradiente de gravidade como perturbação ambiental. A justificativa é que assumi-se um satélite de órbita baixa da
Terra, LEO (Low Earth Orbit). Assumiu-se também que o satélite objeto de simulações, neste caso o EQUARS, tenha um balanço magnético adequado de tal
forma que a perturbação predominante seja associada ao gradiente de gravidade.
Entretanto outros torques podem perturbar o movimento do satélite em órbita da
Terra. A seguir ao feitas algumas considerações sobre outros torques importantes na
análise da dinâmica e controle de atitude de satélites artificiais da Terra.
5.6.2
Torque Aerodinâmico
Para satélites de órbitas baixas (< 900km) a densidade do ar é suficiente para
influenciar a dinâmica de atitude do satélite. Além da altitude o arrasto atmosférico
também depende das dimensões do satélite, da sua geometria e velocidade relativa
(Fauske, 2002). O torque aerodinâmico pode ser escrito como
1
τ a = ρV 2 Cd A (uv × scp )
2
(5.51)
onde ρ é a densidade da atmosfera em (kg/m3 ), A em (m2 ) é a área perpendicular a
uv , uv é o versor na direção da velocidade, Cd é o coeficiente de arrasto aerodinâmico,
V é a velocidade em (m/s) e scp é o vetor distância do centro de massa ao centro
de pressão.
Uma expressão para o pior caso é dada por Fauske (2002)
τa = F (cpa − cg )
F = 0.5 ρCd AV
(5.52)
2
onde cg é o centro de gravidade e cpa é o centro de pressão.
87
(5.53)
5.6.3
Torque de Pressão de Radiação Solar
A radiação e partı́culas do Sol afetam a dinâmica do satélite. Para baixas órbitas o
efeito pode ser negligenciado se comparado com outros torques ambientais (Fauske,
2002). Uma expressão para o pior caso é dada por Larson e Wertz (1992)
τs = F (cps − cg )
Fs
F =
As (1 + q) cos i
c
(5.54)
(5.55)
onde Fs é uma constante solar (1367W/m2 ), As é a área superficial irradiada em
(m2 ), cg é o centro de gravidade (m), i é o ângulo de incidência do Sol, c é a
velocidade da luz na vácuo em (m/s), cps é o centro de pressão solar (m) e q é o
fator de reflectância, que varia de 0 a 1.
5.6.4
Torque Devido ao Dipolo Residual
Os sistemas eletrônicos no satélite podem gerar um dipolo magnético residual. Esse
dipolo residual irá interagir com o campo magnético da Terra. O dipolo residual
também mascara as medidas de Bb , feitas pelos magnetômetros. O torque resultante
pode ser expresso como (Fauske, 2002)
τm = D · B
(5.56)
onde D é o dipolo magnético residual em (Am2 ) varia de 0.1 à 20Am2 , dependendo
do tamanho e desbalanceamento magnético do satélite, e B é o campo magnético
da Terra medido em Tesla. Para uma órbita polar B pode ser aproximado por
B=
2M
Rs2
88
(5.57)
onde M = 7.95 · 1015 T · m3 é o momento magnético da Terra e Rs é distância do
centro de massa da Terra ao centro de massa do satélite.
Nesta Seção faz-se também algumas considerações sobre torques internos ao sistema.
Esses torques não alteram a quantidade de movimento angular do sistema, implicando torque externo nulo. Entretanto, os efeitos dos torques internos podem causar
problemas no apontamento de antenas, painéis, telescópios (tipo Hubble). Por outro
lados efeitos internos dissipadores podem ser úteis, principalmente quando se tratam
dos amortecedores de nutação, quase sempre requeridos quando se usa estabilização
por spin (SCD-1 e SCD-2).
5.7
Perturbações internas
Perturbações internas são torques internos exercidos sobre o corpo do satélite devido
a partes móveis. O efeito de torques internos é a dissipação de energia cinética e a
redistribuição das componentes da quantidade de movimento angular no sistema do
satélite (Wertz, 1978), como acontece na presença de rodas. Entretanto, a presença
de perturbações internas não altera a quantidade de movimento angular referido ao
espaço inercial. As principais fontes de perturbações internas são:
• Movimento de painéis, rodas, mastros, antenas, etc;
• slosh devido ao movimento de lı́quidos (combustı́vel);
• Deformação na estrutura (snap) devido a dilatação térmica;
• Choque de gases no corpo do satélite propelidos pelo sistema de propulsão;
• etc.
Um exemplo bem conhecido onde a influência de perturbações internas, provocado
pelo movimento de antenas, é o do satélite explorer 1, lançado em 1958. O efeito
de dissipação de energia, resultou na mudança do eixo de rotação, eixo de menor
momento de inércia, para o eixo de maior momento de inércia (Kaplan, 1976).
Perturbações internas são difı́ceis de serem modeladas e podem gerar jitter. Uma
aproximação para esse efeito é fornecida por Wilson (2000)
89
T = a0 +
X
ak cos (kωo t + β) + bk sin (kωo t + β)
(5.58)
k
no qual ωo é o velocidade orbital, ao , ak e bk são coeficientes constantes e β a fase. A
influência desse efeito é conveniente para estudos de longa duração, particularmente
para sistema de controle com quantidade de movimento angular embarcado (Wilson,
2000).
90
CAPÍTULO 6
PROJETO DE CONTROLE
Este Capı́tulo é dedicado a apresentar as estratégias de controle usadas para os
modos de operação do satélite: 1) detumble, fase em que o satélite é injetado em
órbita; 2) estabilização, fase em que o satélite adquire a vertical local e realiza a
manutenção da atitude nominal e 3) aquisição, fase em que o satélite faz uma reorientação da atitude. É apresentado o controlador Bdot, os controladores baseados
em energia e as metodologias do Regulador Linear Quadrático (LQR) e do Regulador Quadrático Gaussiano (LQG) usados no modo de estabilização. O método
LQR é analisado primeiramente, formando a base para um melhor entendimento do
método LQG. Essas estratégias são discutidas pormenorizadamente nas referências
Dorato e Cerone (1995), Maciejowski (1989), Kwakernaak e Sivan (1972), Moore e
Anderson (1990), Brown (1997), Kirk (1970), Fauske (2002) e Overby (2004).
6.1
Modo de Detumble
O objetivo do controlador usado no modo de detumble é dissipar a energia cinética
do satélite, reduzindo sua velocidade de rotação.
6.1.1
Controlador Bdot
O controlador usado nesse trabalho para o detumble é proposto por Wisniewski
(1996) e é conhecido como controlador Bdot (Silani e Lovera, 2003). O controlador
Bdot ou de Wisniewski (1996) usa a taxa de variação das medidas feitas pelos
magnetômetros embarcados no satélite. Sua viabilidade é confirmada pela aplicação
no satélite Canadense CanX-1 (Wells e Jeans, 2002).
A lei de controle é dada por
mb = −k Ḃb − mc
onde
91
(6.1)
mc = [0 0 mc ]
(6.2)
no qual k é uma constante escalar maior que zero e mc o momento dipolo nominal
da bobina disposta na direção z.
A demostração analı́tica de que a lei de controle Bdot dissipa a energia cinética de
rotação do satélite, é dada usando a teoria de Lyapunov. Essa demostração é encontrada em Fauske (2002) e Wisniewski (1996). Nos dois procedimentos ou sistemas de
controle de atitude propostos, o detumbling é feito com o uso do controlador acima,
através dos atuadores magnéticos.
6.2
Modo de Estabilização
Para modo de estabilização definido no Capı́tulo 1 as leis de controle são baseadas:
• na metodologia LQR e LQG;
• nos controladores baseados em energia.
Essas estratégias são apresentadas a seguir.
6.2.1
Método LQR
O método LQR é baseado na linearização de sistemas dinâmicos, pois a metodologia
é formulada para sistemas lineares. O problema de controle ótimo consiste em minimizar uma função custo quadrática e gerar uma matriz de ganhos para realimentação
(Dorato e Cerone, 1995 e Maciejowski, 1989). A metodologia LQR é formulada para
os dois sistemas dinâmicos envolvidos: 1) satélite equipado com rodas e 2) satélite
equipado com bobinas. As equações linearizadas dos modelos, Equação (5.33), invariante no tempo e Equação (5.43) variante no tempo são utilizadas para o projeto
de controle. A formulação do problema é feita no domı́nio do tempo e é descrita a
seguir.
92
Seja um sistema linear descrito por
ẋ = A(t)x + B(t)u
(6.3)
O problema de otimização consiste em encontrar uma lei de controle linear do tipo
u = −Kc (t)x
(6.4)
que minimize o ı́ndice de desempenho quadrático
Z
Jp =
T
T
x Qc x + uT Rc u
(6.5)
0
onde Qc é uma matriz positiva semi-definida (Qc ≥ 0), Rc positiva definida (Rc >
0), x é o vetor de estado de dimensão n × 1 e u é o vetor de controle de dimensão
m × 1. As matrizes Qc e Rc são as ponderações no vetor de estado e no vetor de
controle, respectivamente.
Nota: para a existência e estabilidade da solução do problema LQR, a condição
necessária e suficiente é que o sistema seja completamente controlável (Dorato e
Cerone, 1995 e Maciejowski, 1989). As análises de controlabilidade, mostradas em
Kwakernaak e Sivan (1972), feitas para os sistema envolvidos. As equações (5.33) e
(5.43), garantem a solução do problema LQR.
A solução do problema LQR, ou seja, o cálculo do ganho de controle é obtido resolvendo a Equação de Ricatti, dada por (Moore e Anderson, 1990)
T
Ṗc = −Pc A(t) − A(t)T Pc + Pc B(t)R−1
c B Pc − Qc
(6.6)
Resolvendo a Equação diferencial matricial de Ricatti obtemos o controlador variante
no tempo:
93
u = −R−1
c B(t)Pc (t)x
(6.7)
Kc = R−1
c B(t)Pc (t)
(6.8)
Da Equação (6.4), temos
A Figura (6.1) mostra a configuração do LQR para o sistema em malha fechada.
Figura 6.1- Configuração do controle LQR.
O controlador (6.7) exige um grande esforço computacional, devido ao fato de que
a cada passo de integração Pc (t) deve ser computado, nos casos em que o sistema
for variante no tempo, Equação (5.43).
Método LQR Aplicado para o Satélite Equipado com Rodas
Devido ao fato do sistema (5.33) ser invariante no tempo e considerando o intervalo
de otimização infinito, ou seja, T da Equação (6.5) tendendo ao infinito. A solução
do problema LQR é dada para o estado estacionário. Esse problema é referido como
problema do LQR estacionário (steady-state LQR) (Dorato e Cerone, 1995). Nessas
condições a Equação diferencial matricial de Ricatti torna-se uma Equação matricial
algébrica, dada por:
T
0 = −Pc A − AT Pc + Pc BR−1
c B Pc − Qc
94
(6.9)
e a lei de controle para o sistema dado pela Equação (5.33) resulta em
u = −R−1
c BPc x
(6.10)
onde as matrizes A e B são dadas em (5.34) e (5.35), respectivamente.
Método LQR Aplicado para o Satélite Equipado com Bobinas
Para a formulação da lei de controle usada no sistema (5.43) utilizando a metodologia LQR, pode ser feita a seguinte aproximação: consideramos o campo magnético
apresentado na seção 5.2, periódico, e tomamos um valor médio de Bo . Dessa forma
assumimos que o modelo (5.43) seja invariante no tempo. Com essas considerações
e tomando (T → ∞) o problema se reduz ao problema LQR estacionário. Temos
que a Equação de Ricatti resulta
T
0 = −Pc A − AT Pc + Pc BR−1
c B Pc − Qc
(6.11)
e a lei de controle para o sistema dado pela Equação (5.33) resulta em
u = −R−1
c BPc x
(6.12)
Overby (2004) e Silani e Lovera (2003) apresentam resultados satisfatórios usando
o procedimento acima. Entretanto, nesse trabalho, a Equação matricial algébrica
de Ricatti é calculada para cada passo de integração, considerando-se a matriz de
entrada B(t) variante no tempo, obtida a partir das componentes do campo geomagnético, dada pela Equação (5.49). A solução da Equação algébrica matricial de
Ricatti resulta em um ganho variante no tempo para o controlador, dado por
u = −R−1
c B(t)Pc (t)x
95
(6.13)
ou
u = −Kc (t)x
(6.14)
onde Kc (t) = R−1
c B(t)P(t).
As matrices de ponderação Qc e Rc são definidas como:
Rc = diag ([r1 , r2 , · · · , rna ])
(6.15)
Qc = diag ([q1 , q2 , · · · , qns ])
(6.16)
onde na é o número de atuadores no sistema de controle e ns é o número de estados
de interesse. Para o EQUARS nos dois procedimentos/modelos temos na = 3 e
ns = 6. O desempenho desejado do sistema é obtido pelo ajuste das ponderações.
Como sugerido por Overby (2004) e Kristiansen (2000) as ponderações são escolhidas
como:
qi =
1
(∆xi )2
e ri =
1
(∆ui )2
(6.17)
Os valores de ∆ui são baseados no máximo esforço de controle ou valor máximo
de operação dos atuadores, ou seja, máximo torque no caso de rodas e máximo
momento dipolo no caso de bobinas. Os valores de ∆xi são baseados na faixa/intervalo de operação dos estados. Para os sistemas dinâmicos avaliados as escolhas dos
parâmetros usados são baseadas nas especificações do satélite EQUARS, temos
∆xi = 10 graus,
∆ui = 75 mN m (i = 1, 2, 3) e ∆ẋi = 1 graus/s (i = 4, 5, 6)
(6.18)
96
para o modelo do satélite equipado com rodas
e
∆xi = 10 graus,
∆ui = 1 Am2
(i = 1, 2, 3) e ∆ẋi = 1 graus/s (i = 4, 5, 6)
(6.19)
para o modelo do satélite equipado com bobinas.
6.2.2
Método LQG
Nessa seção o problema do LQG é apresentado. A metodologia será implementada
ao sistema (5.43). A solução é obtida no estado estacionário (T → ∞). As medidas
são corrompidas por ruı́dos, modelados como distribuições aleatórias gaussianas. É
apresentado a estrutura do filtro de Kalman-Bucy e o princı́pio da separação.
O problema LQG é apresentado como segue. Dado o sistema da Figura (6.2), onde
G é a função de transferência da planta. O problema LQG pode ser colocado como
sendo o de calcular uma lei de controle que matenha o sistema estável e minimize um
critério de erros quadráticos (Maciejowski, 1989). O sinal de referência r na Figura
é tomado como nulo (r = 0).
Figura 6.2- Sistema planta mais controlador.
97
Representado o sistema da Figura (6.2) na forma de variáveis de estado, temos
ẋ = Ax + Bu + Γw
(6.20)
y = Cx + v
(6.21)
onde x é o vetor de estado, u é o vetor de controle e y é o vetor de saı́das corrompidas por v; w e v são modelados como ruı́dos brancos, caracterizando processos
estocásticos gaussianos com média zero. Considera-se que w e v não são correlacionadas no tempo e entre si, tendo as covariâncias:
E wwT = Rf ≥ 0
E νν T = Qf > 0
E wν T = 0
(6.22)
No problema LQG deseja-se minimizar a função custo:
Z
Jp =
T
T
x Qc x + uT Rc u
(6.23)
0
Onde as matrizes Qc e Rc são as mesmas definidas em 6.2.1. Considera-se que o
vetor de saı́da seja os estados do sistema, nesse caso a matriz de saı́da é uma matriz
identidade, ou seja:
y = Cx
(6.24)
onde C = I, I é a matriz indentidade de ordem 6.
O sistema 6.21 está sujeito a perturbações da planta e ruı́dos na observação da saı́da,
a filosofia do projeto do controlador Kc pode ser estruturada (Athans, 1971) em três
passos:
98
• projeto/análise determinı́stica do problema do controle;
• projeto/análise do problema de estimação estocástica do estado;
• projeto de um sistema de controle estocástico.
A Figura (6.3) mostra a configuração do sistema de controle utilizado, note que o
controle LQG é baseado em observador.
Figura 6.3- Estrutura básica do sistema de controle LQG.
A solução do problema LQG é conseguida pelo princı́pio da separação que possibilita a separação do problema original em dois problemas:
Problema
1: Obter
n
o uma estimativa ótima x̂ do estado x de modo que
T
E (x̂ − x) (x̂ − x) seja minimizado. Esse sub-problema é resolvido pelo uso de
um filtro de Kalman-Bucy, ignorando-se completamente o problema de controle.
Problema 2: Obter o controlador para o problema linear quadrático determinı́stico
(LQR considerando a referência zero (r = 0)), fazendo uso da estimativa x̂ como se
fosse a medida exata do estado, ignorando-se completamente os aspectos estocásticos
do problema. O problema do LQR é resolvido na seção 6.2.1.
Nota: para a existência e estabilidade da solução do problema de estimação e/ou
observação, a condição necessária e suficiente é que o sistema seja completamente
99
observável (Dorato e Cerone, 1995 e Maciejowski, 1989). As análises de observabilidade, mostradas em Kwakernaak e Sivan (1972), feitas para os sistema envolvido,
Equação (5.33), garantem a solução do problema 2.
A solução do problema 1 consiste na determinação da matriz ganho Kf . A estrutura
do filtro de Kalman-Bucy é mostrado na Figura (6.4).
Figura 6.4- Estrutura do filtro de Kalman-Bucy.
A matriz ganho de Kalman Kf é dada por (Athans, 1971 e Maciejowski, 1989)
Kf = Pf CT R−1
f
(6.25)
Onde Pf é dada pela Equação matricial algébrica de Ricatti
0 = −Pf A − AT Pf + Pf CT R−1
f Pf − Qf
(6.26)
É possı́vel mostrar que os autovalores do sistema completo, filtro mais controlador,
são compostos pelo soma dos autovalores do filtro e do LQR (Kwakernaak e Sivan,
1972). As referências Flora (), Dorato e Cerone (1995) e Maciejowski (1989) mos100
tram que os autovalores do filtro e do LQR estão no semi-plano esquerdo (SPE)
do plano complexo, ou seja os dois projetos são assintoticamente estáveis, então o
sistema completo é, também assintoticamente estável. Na análise de robustez os dois
projetos quando considerados isoladamente apresentam boas caracterı́sticas de robustez. Entretanto, como mostrado por Doyle (1978), o sistema completo, formado
pelo filtro de Kalman-Bucy e pelo controlador determinı́stico não mantém as boas
caracterı́sticas de robustez apresentadas, isoladamente, pelo filtro e pelo controlador.
A recuperação das propriedades de robustez é conseguida pelo método LQG/LTR
(loop transfer recovery), esse método é discutido em Maciejowski (1989) e Kwakernaak e Sivan (1972).
6.2.3
Controladores Baseados em Energia
Nesta seção os controladores com realimentação de velocidade angular e atitude
são investigados para o modo de estabilização. Ambos são baseados em energia e
usam os atuadores magnéticos para o controle. Diz-se que esses controladores são
baseados em energia pois a partir do cálculo da energia do sistema é escolhido uma
lei de controle que satisfaz os critérios de estabilidade, usando para isso a teoria de
Lyapunov. A função de Lyapunov V é escolhida como sendo a energia do sistema. As
discussões são baseadas nas referencias Wisniewski (1996) e Fauske (2002). A prova
da estabilidade dos controladores propostos e a obtenção dos pontos de equilı́brio
estáveis do sistema são obtidos a partir da teoria de Lyapunov, e não fazem parte da
proposta desse trabalho. Essas provas são fornecidas em Wisniewski (1996), Fauske
(2002) e Overby (2004).
Realimentação de Velocidade Angular
A lei de controle com realimentação de velocidade angular é dada por (Fauske, 2002)
mb = hω bob × Bb
(6.27)
onde h é uma constante escalar maior que zero (h > 0), que faz com que o sistema
satélite equipado com bobinas, seja assintoticamente estável sobre quatro pontos de
equilı́brio.
101
ω bob , cb3 , cb2 : (0, ±co3 , ±co2 )
(6.28)
onde coi i = 1, 2, 3 são matrizes coluna tendo como componentes os cosenos diretores
dos eixos do referencial do corpo em relação aos eixos do referencial orbital. Os
quatro pontos de equilı́brio correspondem a orientação dos eixos x, y e z na mesma
direção dos eixos xo , yo e zo , respectivamente.
Realimentação de Atitude
A lei de controle com realimentação de atitude é dada por (Fauske, 2002)
mb = hω bob × Bb − α × Bb
(6.29)
onde é a componente vetorial do quaternion, Equação (4.7), e α é uma constante
escalar, essa lei de controle faz com que o sistema, satélite equipado com bobinas,
seja assintoticamente estável sobre a referência ou a VLHL.
6.3
ω bob , cb3 , cb2 : (0, co3 , co2 )
(6.30)
Modo de Aquisição
Para o estudo de aquisição de uma atitude arbitrária, definido nesse trabalho como
modo de aquisição é feita a implementação da metodologia LQR rastreio (tracking)
e um controlador proporcional derivativo (PD). O objetivo é avaliar, através da formulação desenvolvida nesse trabalho, a factibilidade do uso das rodas da TELDIX
especificadas para o satélite EQUARS, para manobras de atitude e também comparar e discutir as leis de controle formuladas. Nas simulações os controladores são
aplicados à dinâmica completa, não linear.
102
6.3.1
LQR Tracking
Na aquisição de uma atitude arbitrária, onde a referência difere de zero (x 6= 0),
o problema de formular uma lei de controle utilizando a metodologia LQR é conhecido como LQR rastreio/tracking. No LQR traking a lei de controle envolve um
termo antecipativo (feedfoward ) adicionado ao termo de realimentação do estado
(state-feedback ) (Dorato e Cerone, 1995 e Kirk, 1970), apresentado na seção 6.2.1. O
problema de LQR tracking consiste em minimizar um ı́ndice de desempenho definido
como
Z
Jp =
T
T
x̃ Qc x̃ + uT Rc u
(6.31)
0
onde
x̃(t) = x(t) − xd (t)
(6.32)
onde xd é a trajetória de referência ou o estado desejado, x é obtido através do
sistema dinâmico
ẋ = A(t)x + B(t)u
(6.33)
No problema assumimos que o estado desejado (xd ) é conhecido e que o todos os
estados x são disponı́veis para realimentação. As matrizes Qc e Rc são as mesmas
definidas na seção 6.2.1. A lei de controle ótima para o problema de tracking é dada
por (Kirk, 1970)
u = −Kc (t)x + fw (t)
(6.34)
onde Kc (t) é a matriz de ganhos para realimentação dos estados, obtido pela Equação
(6.8), obtido através da Equação de Ricatti, Equação (6.6), definida na seção 6.2.1:
103
Ṗ = −PA(t) − A(t)T P + PB(t)R−1
c P − Qc
(6.35)
e o termo antecipativo (fw (t)) é denominado sinal de comando e é dado por
fw (t) = −Rc BT s(t)
(6.36)
onde s(t) é cálculado a partir da Equação diferencial linear (Kirk, 1970)
−1
s(t) + Qc xd
ṡ = − A(t) − PB(t)R−1
c B(t)
(6.37)
Note que o sinal de comando (fw (t)) depende dos parâmetros do sistema e do sinal
de referência. A Figura (6.5) mostra a configuração do LQR tracking.
Figura 6.5- Configuração do controle LQR tracking.
6.3.2
Proporcional Derivativo
O controlador proporcional derivativo (PD), também é projetado a partir da
planta linearizada, Equação (5.33). Nesse trabalho o projeto do controlador PD
é baseado em Kaplan (1976). Uma discussão pormenorizada dos controladores industriais PD e PID é encontrada em Ogata (1998). Primeiro é obtida uma lei de
controle para o movimento em torno do eixo de pitch, como se fosse um sistema
104
SISO, pois o movimento em pitch (θ) é desacoplado dos outros movimentos, roll e
yaw. Em seguida é projetada uma lei de controle para o sistema de equações diferenciais acopladas, que descrevem o movimento de roll e yaw, e constitue um sistema
MIMO com duas entradas e duas saı́das.
Projeto em Pitch
A partir da Equação (5.33), obtemos a Equação linearizada do movimento em pitch
dada por
J2 θ̈ + 3ωo2 (J1 − J3 ) θ + u2 = 0
(6.38)
onde u2 representa o torque de controle da roda de reação disposta na direção
do eixo de pitch. Note que a Equação (6.38) indica que o sistema não apresenta
amortecimento, isso pode ser obtido através do sinal de controle u2 . Uma forma
satisfatória para a lei de controle, discutida em Kaplan (1976) é dada por
u2 = Kp2 Td2 θ̇ + θ − θr
(6.39)
onde o termo θ̇ introduz amortecimento ao sistema, Kp2 é o ganho proporcional
(P) e Td2 a constante de tempo. O ganho derivativo é dado por Kd2 = Kp2 Td2 , θr
é o ângulo de pitch de referência, considerado constante. Com a lei de controle o
movimento em pitch torna-se o clássico sistema amortecido de segunda ordem.
J2 θ̈ + Kp2 Td2 θ̇ + Kp2 + 3ωo2 (J1 − J3 ) θ − Kp2 θr = 0
(6.40)
Os ganhos do controlador PD (Kp2 e Kd2 ), são obtidos nesse trabalho a partir da
metodologia LQR, ou seja, a partir da matriz Kc tomando-se os ganhos da linha 2,
que correspondem ao sinal de controle u2 . Uma outra alternativa para escolha dos
ganhos é ajustar/sintonizar os ganhos a partir do diagrama de lugar das raı́zes, esse
procedimento é usado em Kaplan (1976).
105
Projeto em Roll e Yaw
A partir da Equação (5.33), obtemos as equações linearizadas do movimento em roll
e yaw dadas por
J1 φ̈ + a1 φ + a2 ψ̇ + u1 = 0
(6.41)
J3 ψ̈ + c1 ψ + c2 φ̇ + u3 = 0
(6.42)
onde os coeficientes a1 , a2 , c1 e c2 , são
a1 = 4ωo2 (J2 − J3 ) + ωHo
(6.43)
a2 = −ωo (J1 − J2 + J3 ) + Ho
(6.44)
c1 = ωo2 (J2 − J1 ) + ωHo
(6.45)
c2 = ωo (J1 − J2 + J3 ) − Ho
(6.46)
os sinais de controle u1 e u3 são os torques de controle gerados pelas rodas de reação
dispostas na direção do eixo de roll e yaw, respectivamente. Devido ao acoplamento
dos movimentos as leis de controle são escolhidas como sendo
u1 = Kp1 Td1 φ̇ + Kp1 (φ − φr ) + Kc1 ψ̇
(6.47)
é lei de controle para o movimento em roll. O termo φ̇ introduz amortecimento
ao sistema, Kp1 é o ganho proporcional (P), Td1 a constante de tempo. O ganho
derivativo é dado por Kd1 = Kp1 Td1 , φr é o ângulo de roll de referência, considerado constante. O termo Kc1 é inserido devido ao acoplamento entre os movimentos
em roll e yaw. Com a lei de controle o movimento em roll torna-se um sistema
amortecido de segunda ordem.
106
J1 φ̈ + a1 φ + a2 ψ̇ + Kp1 Td1 φ̇ + Kp1 (φ − φr ) + Kc1 ψ̇ = 0
(6.48)
Anologamente, a lei de controle para o movimento em yaw é escolhida como
u3 = Kp3 Td3 ψ̇ + Kp3 (ψ − ψr ) + Kc3 φ̇
(6.49)
O termo ψ̇ introduz amortecimento ao sistema, Kp3 é o ganho proporcional (P),
Td3 a constante de tempo, o ganho derivativo é dado por Kd3 = Kp3 Td3 , ψr é o
ângulo de yaw de referência, considerado constante. O termo Kc3 é inserido devido
ao acoplamento entre os movimentos. Com a lei de controle o movimento em yaw
torna-se em um sistema amortecido de segunda ordem.
J3 ψ̈ + c1 ψ + c2 φ̇ + Kp3 Td3 ψ̇ + Kp3 (ψ − ψr ) + Kc3 φ̇ = 0
(6.50)
Reescrevendo a equações (6.41) e (6.42) temos
J1 φ̈ + Kp1 Td1 φ̇ + (a1 + Kp1 ) φ − Kp1 φr + (Kc1 + a2 ) ψ̇ = 0
(6.51)
J3 ψ̈ + Kp3 Td3 ψ̇ + (c1 + Kp3 ) ψ − Kp3 ψr + (Kc3 + c2 ) φ̇ = 0
(6.52)
A Figura (6.6) mostra a configuração do controlador PD.
Os ganhos são obtidos a partir da metodogia LQR (matriz Kc ), assim como feito
para o projeto de controle em pitch. A partir da matriz Kc toma-se os ganhos das
linhas 1 e 3, que correspondem ao sinal de controle para u1 e u3 , respectivamente.
Os valores de referência (φr , θr , ψr ) correspondem a uma atitude nominal arbitrária.
A vantagem do uso do controlador PD é que a realimentação é apenas dos ângulos
de atitude, não sendo necessário medir as taxas de variação dos ângulos, entretanto,
o procedimento de derivação dos ângulos gera erros.
107
Figura 6.6- Configuração do controlador PD.
108
CAPÍTULO 7
SIMULAÇÕES
Neste Capı́tulo, são feitas as simulações utilizando as estratégias de controle para os
modos de operação:
• Detumble;
• Estabilização;
• Aquisição.
As configurações propostas para o ACS utilizam como atuadores: rodas de reação
combinadas com bobinas magnéticas e apenas bobinas magnéticas para o controle
autônomo do satélite, nos dois primeiros modos de operação (detumble e estabilização). O controlador Bdot é utilizado para o descapotamento/detumbling do
satélite. O controlador proposto para estabilização e manutenção da atitude nominal
é obtido usando a metodologia LQR e LQG. Outra alternativa apresentada para fase
de estabilização e manutenção da atitude nominal, usando bobinas magnéticas, são
os controladores baseados em energia: realimentação da velocidade angular (angular
velocity feedback ) e realimentação da atitude ()(attitude feedback ). Nas simulações
são usados os parâmetros do satélite brasileiro EQUARS. Para a fase de aquisição
é analisada a exeqüibilidade da utilização das rodas de reação onde é proposto o
controle PD e LQR rastreio/tracking.
Os parâmetros usados na simulação são mostrados na Tabela (7.1). Os parâmetros
e especificações obtidos para o satélite EQUARS são baseados em Carvalho (2003)
e Heidelberg (2004). As três rodas de reação são iguais e as três bobinas também
são iguais.
7.1
Modo de Detumble
O controlador usado para o modo de detumble foi o controlador Bdot. As condições
iniciais e os parâmetros do controlador são listados a seguir.
109
Tabela 7.1- Parâmetros de simulação.
Parâmetros
Valores
Jx = 13.31, Jy = 14.22, Jz = 11.20 kgm2
1o sobre cada eixo (φ, θ, ψ)
Tensor de inércia do satélite (J)
Precisão de apontamento
Bobinas magnéticas:
Máximo momento dipolo
Seção da área
Número de espiras
Resistência
2 Am2
0.075 m2
100
20 Ω
Rodas de reação:
Momento de inércia
Máxima rotação por minuto
Máxima quantidade de movimento angular
Máximo torque
Parâmetros orbitais:
Altura (h)
Excentricidade (e)
Ascenção do nodo ascendente (Ω)
Inclinação (i)
Velocidade orbital média
Perı́odo orbital
Data da simulação
Jw = 0.015 kgm2
7500 rpm
12 N ms
75 mN m
750 km
∼
=0
30 graus
20 graus
ωo = 0.001 rad/s
To ∼
= 1h40mim
31 − 12 − 2004
110
Condições iniciais
Parâmetros de controlador
b
ωib
= [10 10 10]T graus/s
k = 2 · 105 ,
mc = −0.1 Am2
onde a lei de controle Bdot é dada pela Equação (6.1)
mb = −k Ḃb − mc
(7.1)
A escolha da velocidade angular inicial foi baseada nas simulações feitas para o
modo de estabilização, Seção (6.2). Para essa condição inicial de velocidade angular
b
b
b
(ωix
, ωiy
, ωiz
) (≈ 2 rpm) os atuadores especificados saturam rapidamente no modo
de estabilização e não são capazes de estabilizar ou controlar a atitude do veı́culo,
definindo uma situação ou velocidade angular crı́tica. Os parâmetros do controlador
foram obtidos por tentativa. A Figura (7.1) mostra a redução da velocidade angular
do satélite. Depois de quatro revoluções em volta da Terra, o veı́culo já está apto a
executar manobras de atitude, ou seja, entrar no modo de estabilização.
Figura 7.1- Velocidade angular do satélite ω bib .
A Figura (7.2) mostra o amortecimento da velocidade angular do satélitea partir
da quarta órbita. Reduzindo a velocidade de rotação do satélite para abaixo de (2
graus/s) .
111
Figura 7.2- velocidade angular ω bib .
A Figura (7.3) mostra o dipolo magnético das bobinas, o torque de controle gerado e
a potência de cada bobina. Os valores do dipolo magnético atendem as especificações
dos atuadores magnéticos (mx , my e mz < 2 Am2 ).
Figura 7.3- Dipolo magnético, torque magnético e potência.
112
A Figura (7.4) mostra os mesmos resultados mas tomados a partir da quarta órbita.
Figura 7.4- Dipolo magnético, torque magnético e potência.
A lei de controle apresentada (Bdot) é de simples implentação e não requer informações da atitude do veı́culo. Apesar das simulações terem avaliado uma condição
inicial crı́tica de velocidade angular, essa pode ser maior, dependendo de como o foguete lança o satélite na sua órbita, ou seja, depende do veı́culo lancador. Ainda
assim o controlador Bdot pode ser utilizado, entretanto, exigindo um tempo maior
para a redução da velocidade angular.
7.2
Modo de Estabilização
As estratégias de controle usada no modo de estabilização são baseadas nas metodologias LQR e LQG, apresentadas no Capı́tulo 6. A metodologia LQR é implementada
nos dois modelos: satélite equipado com bobinas e satélite equipado com rodas. A
metodologia LQG é usada no modelo do satélite equipado com rodas. Os controladores baseados em energia são usados no modelo do satélite equipado com bobinas.
113
7.2.1
Método LQR
Nessa Seção o método LQR é usado na formulação da lei de controle do modelo do
satélite equipado com rodas. São investigados três cenários:
• condições iniciais crı́ticas (ωx , ωy , ωz > 1 rpm) ;
• condições iniciais quase crı́ticas (ωx , ωy , ωz = 1 rpm);
• velocidade angular (ωx , ωy , ωz ∼
= 0 rpm).
O objetivo é avaliar situações nas quais o projeto LQR ainda é válido no controle de
atitude em três eixos. A Tabela (7.2) mostra as condições iniciais para os três casos.
Tabela 7.2- Condições iniciais da simulação para o modo de estabilização utilizando a metodologia LQR.
casos
1
2
3
velocidade angular inicial
atitude inicial
ω bib = [0 0 0]T graus/s
ω bib = [6 6 6]T graus/s
ω bib = [10 10 10]T graus/s
(φ, θ, ψ) = (30, −20, 25) graus
(φ, θ, ψ) = (30, −20, 25) graus
(φ, θ, ψ) = (30, −20, 25) graus
Em seguida será usada a formulação da lei de controle do modelo do satélite equipado
com bobinas. O objetivo e estudar/analisar os dois procedimentos para o modo de
estabilização
Satélite Equipado com Rodas
Caso 1: Velocidade Angular
Para o primeiro caso considera-se que a velocidade angular do satélite seja nula,
conseguida através do controlador Bdot. A partir dessa condição inicial temos que
a lei de controle para o modo de estabilização é dada pela Equação (6.7), reescrita
aqui
114
u = −Kc x
(7.2)
Kc = R−1
c BPc
(7.3)
onde
onde Pc é dado pela solução da Equação de Ricatti, Equação (6.6). Os valores
númericos das matrizes A e B, usadas no projeto de controle, são obtidos a partir
das equações (5.34) e (5.35). As matrizes de ponderação Rc e Qc selecionadas são
mostradas na Seção (6.2.1). Resolvendo a Equação de Ricatti, obtemos o ganho
usado na lei de controle


−0.0286 0.0000
0.0004 −0.8476 0.0000 0.0000


Kc =  0.0000 −0.0287 0.0000
0.000 −0.9180 0.0000
−0.0004 0.0000 −0.0286 0.000
0.0000 0.9510
(7.4)
É importante resaltar que em todas as simulações aplicou-se a lei de controle linear
para o modelo completo, não linear. A Figura (7.5), mostra os ângulos de atitude
em função do tempo, durante o modo de estabilização. O resultado mostra um bom
desempenho do sistema de controle. Note que a partir de 100 segundos o satélite
já atinge a precisão requerida (< 1o ) para a atitude nominal especificada, (VLHL),
realizando em seguida a manutenção da atitude nominal.
A Figura (7.6) mostra o torque gerado pelas três rodas de reação. O primeiro gráfico
mostra o torque (τ bw ) gerado pelas rodas dispostas nas três direções principais de
inércia. O segundo gráfico mostra o torque de acoplamento devido ao movimento de
rotação do satélite, sendo da ordem de 100 vezes menor que o torque τ bw gerado pelas
rodas. O torque efetivo de controle é a soma dos dois torques. Para efeito de projeto
da roda o efeito de acoplamento é considerado uma perturbação. O último gráfico
mostra a quantidade de movimento angular das rodas, note que a roda ao longo do
eixo de pitch (normal a órbita) tem uma quantidade de movimento angular nominal
diferente de zero, fornecendo rı́gidez giroscópica ao sistema durante a manutenção
115
Figura 7.5- Ângulos de atitude roll, pitch e yaw em função do tempo para o caso
1.
da atitude nominal.
Figura 7.6- Torque τ bw , torque de acoplamente e quantidade de movimento angular
das rodas para o caso 1.
A Figura (7.7) mostra as rotações por minuto de cada uma das rodas e as velocidades
116
angulares ω bib e ω bob , respectivamente. Como era de se esperar a velocidade relativa
entre os referenciais do satélite (BF) e orbital (OF) é nula.
Figura 7.7- Velocidades angulares ω bib e ω bob e rotações por minuto das rodas de
reação.
A simulações mostram que, para a condição inicial de velocidade angular assumida,
Tabela (7.2), todos os requisitos de precisão e operação são atendidos, não havendo
saturação dos atuadores no modo de estabilização.
Caso 2: Velocidade Angular Quase Crı́tica
Para o segundo caso considera-se que a velocidade angular residual do satélite seja
de 1 rpm, que equivale a 6 graus/s. A partir dessa condição inicial usa-se a lei de
controle dada pela Equação (6.7), como no caso 1. O objetivo é avaliar se mesmo
sob uma condição de velocidade angular não nula. Esta configuração poderia, por
exemplo, ser causada por imprecisão de medidas ou falha no sistema de controle
para o modo de detumble. Os resultados mostram que ainda assim, o sistema de
controle ainda é capaz de estabilizar e controlar o movimento de atitude, a partir
dos atuadores (rodas) especı́ficadas para o satélite EQUARS. Os valores de Kc
obtidos e as matrizes de ponderação Rc e Qc selecionadas são as mesmas utilizadas
no caso 1.
117
A Figura (7.8), mostra os ângulos de atitude durante o modo de estabilização. O
resultado mostra um desempenho pobre do sistema de controle. Note que, apesar do
desempenho, a partir de 200 segundos o satélite atinge a precisão requerida (< 1o )
para a atitude nominal especificada (VLHL), realizando em seguida manutenção da
atitude nominal sem maiores problemas.
Figura 7.8- Ângulos de atitude roll, pitch e yaw em função do tempo para o caso
2.
A Figura (7.9) mostra o torque gerado pelas três rodas de reação, o primeiro gráfico
mostra o torque (τ bw ) gerado pelas rodas, o gráfico mostra, claramente, a saturação
dos atuadores, todas as rodas em algum momento sofrem saturação o segundo gráfico
mostra o torque de acoplamento devido ao movimento de rotação do satélite, sendo
bem maior a influência que no caso 1. O último gráfico mostra a quantidade de movimento angular das rodas, note que apesar da saturação em torques (> 75mN m),
nenhuma das rodas apresentam saturação na quantidade de movimento angular
(< 12N ms).
A Figura (7.10) mostra as rotações por minuto de cada uma das rodas e as velocidades angulares ω bib e ω bob .
A simulações mostram que para a condição inicial de velocidade angular assumida
(Tabela (7.2)) mesmo sob saturação das rodas, todos os requisitos de precisão são
118
Figura 7.9- Torque τ bw , torque de acoplamento e quantidade de movimento angular
das rodas para o caso 2.
Figura 7.10- Velocidades angulares ω bib e ω bob e rotações por minuto das rodas para
o caso 2 .
atendidos pelo sistema de controle.
119
Caso 3: Velocidade Angular Crı́tica
Para o terceiro caso analisado usamos o mesmo projeto de controle dos casos 1 e 2,
e as condições iniciais mostradas na Tabela (7.2). A Figura (7.11) mostra os ângulos
de atitude em relação ao tempo. Nota-se que o sistema de controle não é capaz de
atender os requisitos no modo de estabilização, devido a uma velocidade angular
crı́tica de aproximamente 1.5 rotações por minuto (rpm). Portanto, para o projeto
do sistema de controle avaliado é necessário garantir que no modo de estabilização
a velocidade angular inicial seja igual ou menor que 1.5 rpm. Essa redução é obtida
no modo de detumble como mostrado nas simulações, Seção (7.1).
Figura 7.11- Ângulos de atitude roll, pitch e yaw e tempo para o caso 3.
A Figura (7.12) mostra os gráficos das velocidades angulares ω bib e ω bob e as rotações
por minuto das rodas, nota-se que as rodas dispostas ao longo dos eixos de roll e yaw
saturam (> 7500rpm). A roda em roll atinge rapidamente sua velocidade máxima
de operação.
A Figura (7.13) mostra a saturação em torque e quantidade de movimento angular
das rodas.
Os resultados da simulação mostram que o sistema de controle não é capaz de
120
Figura 7.12- Velocidades angulares ω bib e ω bob e rotações por minuto das rodas para
o caso 3.
Figura 7.13- Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas para o caso 3.
estabilizar e controlar a atitude. Isso se deve, fundamentalmente, a limitação dos
atuadores e a dinâmica se tornar altamente não linear.
121
Satélite Equipado com Bobinas
A formulação da lei de controle utilizando a metodologia LQR para o modelo do
satélite equipado com bobinas foi discutida na Seção (6.2.1). O ganho de controle
Kc é variante no tempo, dada pela Equação (6.13)
u = −Kc (t)x
(7.5)
A matriz solução da Equação de Ricatti Pc , Equação (6.6), é calculada para
cada passo de integração, sendo aproximadamente periódica devido ao campo geomagnético (Bo ) ser quase periódico. Os valores númericos das matrizes A, B, usadas
no projeto de controle, são obtidos a partir das equações (5.48) e (5.49). A Figura
(7.14) mostra a variação do elemento Kc (1, 1) da matriz de ganho do controlador,
mostrando claramente que os ganhos são aproximadamente periódicos. Isso permite
reduzir o esforço computacional introduzindo uma matriz de ganhos previamente
calculada para uma órbita nominal.
Figura 7.14- Ganho do controlador em função das órbitas Kc (1, 1).
As matrizes de ponderação Rc e Qc selecionadas são mostradas na Seção (6.2.1). As
condições iniciais são listadas a seguir
122
Condições iniciais:
Velocidade angular
Atitude inicial
ω bib = [0 0 0]T graus/s
(φ, θ, ψ) = (15, −10, 10)
Novamente utilizou-se nas simulações o modelo completo, não linear e no caso do
satélite equipado com bobinas o modelo também é variante no tempo. A Figura
(7.15) mostra os ângulos de atitude durante o modo de estabilização, o resultado
mostra um desempenho razoável do sistema de controle, a partir dos parâmetros
selecionados para o controlador. Note que é necessário duas revoluções do satélite
em torno da Terra (≈ 3h20min) para que esse atinga a precisão requerida (< 1o )
realizando em seguida a manutenção da atitude nominal, ou seja, alinhado com a
vertical local (VLHL).
Figura 7.15- Ângulos de atitude roll, pitch e yaw em função do número de órbitas,
utilizando bobinas.
A Figura (7.16) mostra o dipolo magnético das bobinas, o torque de controle gerado,
e a potência de cada bobina, os valores do dipolo magnético atendem as especificações dos atuadores magnéticos (mx , my e mz < 2 Am2 ), ficando bem abaixo do
dimensionado para o sistema de controle de atitude.
Apesar das bobinas apresentarem um tempo de resposta bem maior, a potência
exigida pelo sistema de controle é muito menor do que aquela usada no sistema de
123
Figura 7.16- Dipolo magnético, torque de controle (τ bm ) e potência das bobinas
em função do número de órbitas, utilizando bobinas.
controle de atitude do satélite equipado com rodas; tornando-se uma opção atrativa,
devido a seu custo e baixo consumo de energia.
7.2.2
Controladores Baseados em Energia
Uma alternativa para o sistema de controle de atitude empregando bobinas
magnéticas são os controladores não lineares propostos na Seção (6.2.3). Sua principal vantagem é a de que as leis de controle são formuladas a partir do sistema
não linear, baseado nos crı́terios de estabilidade de Lyapunov. Isso garante a estabilização, não sendo válida somente em torno de um ponto de operação, como no caso
de sistema lineares. São analisados duas leis de controle: realimentação de velocidade
angular e realimentação de atitude.
Realimentação de Velocidade Angular
A lei de controle é dada pela Equação (6.27)
124
mb = hω bob × Bb
(7.6)
As condições iniciais e os parâmetros do controlador obtidos por tentativa são listados a seguir
Condições iniciais:
Velocidade angular
ω bib = [0 0 0]T graus/s
Atitude inicial
(φ, θ, ψ) = (60, −30, 40)graus
Parâmetros do controlador
h = 3.25 · 106
A Figura (7.17) mostra os ângulos de atitude no modo de estabilização usando a
realimentação de velocidade angular. O resultado mostra um desempenho razoável
para os ângulos de roll e pitch. Para o eixo de yaw é obtido um pobre desempenho do
sistema de controle. Esse resultado ocorre devido a falta de informação da atitude.
É necessário quatro revoluções do satélite em torno da Terra (≈ 6h40min) para que
esse atinga a precisão requerida (< 1o ) nos eixos de roll e pitch, sendo necessário
quase o dobro de revoluções para o eixo de yaw. Em seguida realiza-se a manutenção
da atitude nominal, ou seja, alinhado com a vertical local (VLHL).
Figura 7.17- Ângulos de atitude em função do número de órbitas.
A Figura (7.18) mostra o dipolo magnético das bobinas, o torque de controle gerado,
125
e a potência de cada bobina, os valores do dipolo magnético atendem as especificações dos atuadores magnéticos, ficando abaixo do dimensionado para o sistema
de controle de atitude. O terceiro gráfico mostra que a potência exigida pelas bobinas
é bem baixa.
Figura 7.18- Dipolo magnético, torque de controle e potência das bobinas em
função do número de órbitas.
A Figura (7.19) mostra as velocidades angulares e o campo geomagnético expresso
no sistema do satélite.
É possı́vel mostrar (Fauske, 2002 e Overby, 2004) que a lei de controle apresentada
não garante, sobre outras condições iniciais, a estabilização e o controle em três eixos
de interesse, ou seja, alinhado com a VLHL. Isso se deve, ao sistema dinâmico ter
mais de um ponto de equilı́brio, sendo necessário informações da atitude para a estabilização de interesse. A lei de controle apresentada é uma alternativa que deve ser
considerada, para os casos em que nenhuma informação da atitude é disponı́vel para
o sistema de controle, devido a possı́veis falhas nos sensores de atitude. Entretanto,
não se pode garantir, em alguns casos, o controle em para a orientação de interesse
(VLHL). Novamente o uso de bobinas apresenta um baixo consumo de energia.
126
Figura 7.19- Velocidades angulares ω bib e ω bob e campo magnético terrestre local
Bb .
Realimentação de Atitude
A lei de controle é dada pela Equação (6.29)
mb = hω bob × Bb − α × Bb
(7.7)
As condições iniciais e os parâmetros do controlador são listados a seguir
Condições iniciais:
Velocidade angular
ω bib = [0 0 0]T graus/s
Atitude inicial
(φ, θ, ψ) = (60, −30, 40)graus
Parâmetros do controlador
h = 3.25 · 106
α = −3000
A Figura (7.20) mostra os ângulos de atitude no modo de estabilização usando a
realimentação de velocidade angular e atitude. O resultado mostra um desempenho
razoável para os ângulos de roll, pitch e yaw. É necessário duas revoluções do satélite
em torno da Terra (≈ 3h20min) para que esse atinja a precisão requerida (< 1o )
nos eixos de roll, pitch e yaw. O desempenho do sistema de controle é próximo
127
ao obtido com a metodologia LQR, entretanto, a atitude inicial assumida é bem
maior. O resultado mostra que a informação da atitude é importante para um bom
desempenho do sistema de controle, na estabilização da atitude em três eixos. A
partir da segunda revolução do satélite em torno da Terra ocorre a manutenção
da atitude nominal (VLHL), ou seja, o sistema se estabiliza em torno da atitude
nominal.
Figura 7.20- Ângulos de atitude.
A Figura (7.21) mostra o dipolo magnético das bobinas, o torque de controle gerado, e a potência de cada bobina. Os valores do dipolo magnético atendem as
especificações dos atuadores magnéticos, Tabela (7.1), ficando abaixo do dimensionado para o sistema de controle de atitude. O terceiro gráfico mostra que a potência
exigida pelas bobinas é bem baixa.
A Figura (7.22) mostra as velocidades angulares e o campo geomagnético expresso
no sistema do satélite.
A lei de controle apresentada é uma boa alternativa para o sistema de controle de
atitude em três eixos, pois apresenta um desempenho razoável. A principal vantagem
em relação a metodologia LQR é que mesmo sob as não linearidades, a eficiência do
sistema de controle de atitude na estabilização do sistema é garantida. Novamente o
128
Figura 7.21- Dipolo magnético, torque de controle e potência das bobinas em
função do número de órbitas.
Figura 7.22- Velocidades angulares ω bib e ω bob e campo magnético local Bb em
função do número de órbitas.
uso de bobinas apresenta uma opção atrativa devido a baixa potência, o que equivale
a um baixo consumo de energia.
129
7.2.3
Método LQG
Como descrito na Seção (6.2.2), metodogia LQG é usada no modo de estabilização
para o satélite equipado com rodas. As condições iniciais, e os parâmetros do filtro
de Kalman são listados na Tabela (7.3)
Tabela 7.3- Parâmetros de simulação.
Condições iniciais:
Velocidade angular ω bib
Atitude inicial (φ, θ, ψ)
Parâmetros do filtro:
Matriz Rf
Matriz Qf
[0 0 0]T graus/s
(30, −20, 25)graus
diag[(1, 1, 1, 0.1, 0.1, 0.1)]0.52 × (π/180o )
diag[(1, 1, 1, 1, 1, 1)]152 × (π/180o )
A matriz Rf é obtida a partir do desvio padrão de sensores solares (0.5o ) e de
girômetros (0.05graus/s) tı́picos. A matriz Q é obtida através de ajustes por tentativa. A partir dos parâmetros do filtro obtemos o matriz de ganho de Kalman


30.0002 0.0000 0.0000
0.9091
0.0000
0.0000



 0.0000 30.0002 0.0000
0.000
0.9109
0.0000


 0.0000 0.0000 30.0002

0.000
0.0000
0.9509


Kf = 

 0.0000 0.0000 0.0000 300.0000 0.0000

0.0001


 0.0000 0.0091 0.0000

0.0000
300.0000
0.0000


0.0000 0.0000 0.0091
0.001
0.0000 300.0000
(7.8)
Os parâmetros do controlador usado são os mesmos utilizados no método LQR.
A Figura (7.23) mostra os ângulos de atitude ou de Euler estimados com o filtro de
Kalman. O desempenho do controlador é muito próximo ao obtido com o projeto
LQR (caso 1).
A Figura (7.24), mostra os ângulos de Euler simulados, esses são corrompidas por
ruı́dos de distribuição aleatória e uniforme.
130
Figura 7.23- Ângulos de atitude estimados em função do tempo.
Figura 7.24- Ângulos de atitude simulados.
A Figura (7.25) confronta os ângulos de atitude estimados e simulados, mostrando
claramente a suavização do sinal estimado, que é usado na realimentação do sistema
de controle. A Figura (7.26) mostra os gráficos do erro (valor real (simulado) −
valor estimado). O primeiro gráfico mostra o erro nos ângulos de atitude e o segundo
131
gráfico mostra o erro nas taxas de variação dos ângulos de atitude. Em ambos os
casos o erro é menor que a covariância em ∼
= 70% dos casos, representada pela linha
vermelha.
Figura 7.25- Ângulos de atitude simulados e estimados.
Figura 7.26- Erro dos ângulos de atitude e variação da atitude.
A Figura (7.27) mostra o torque gerado pelas três rodas de reação, o torque de
132
acoplamento e a quantidade de movimento angular. O primeiro gráfico mostra o
torque (τ bw ) gerado por cada roda, o segundo gráfico mostra o torque de acoplamento
devido ao movimento de rotação do satélite. As magnitudes ou esforço de controle
são próximas às obtida no projeto LQR (caso 1), como era de se esperar, pois o
controlador é o mesmo desenvolvido no caso 1. O último gráfico mostra a quantidade
de movimento angular das rodas. A roda ao longo do eixo de pitch (normal a órbita)
tem uma quantidade de movimento angular nominal diferente de zero, fornecendo
rı́gidez giroscópica ao sistema durante a manutenção da atitude.
Figura 7.27- Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas.
7.3
Modo de Aquisição
Para o sistema de controle no modo de aquisição são avaliados duas alternativas
para o projeto de controle: 1) controlador Proporcional Derivativo (PD); 2) LQR
rastreio/tracking. Esse projeto tem o objetivo de avaliar a factibilidade do emprego
das rodas especificadas para o satélite brasileiro EQUARS (ver Heidelberg, 2004),
na execusão de manobras de grandes ângulos. No modo de aquisição, assumese que a condição inicial o satélite está alinhado com a vertical e horizontal local
(VLHL), ou seja, o satélite parte do modo de estabilização. As condições iniciais e
a atitude de refêrencia são mostradas na Tabela (7.4).
133
Tabela 7.4- Condições iniciais da simulação para o modo de aquisição utilizando
rodas.
Lei de controle
velocidade angular inicial
atitude inicial
(φ, θ, ψ)
atitude de referência
(φr , θr , ψr )
LQR tracking
PD
ω bib = [0 0 0]T graus/s
ω bib = [0 0 0]T graus/s
(0, 0, 0) graus
(0, 0, 0) graus
(−60, 70, 130) graus
(−60, 70, 130) graus
7.3.1
LQR Tracking
A lei de controle ao método LQR tracking no modo de aquisição, apresentada na
Seção (6.3.1), é dada pela Equação (6.34),
u = −Kc (t)x + fw (t)
(7.9)
onde a matriz Kc é obtida através da solução da Equação de Ricatti. As matrizes de
ponderação Qc e Rc são as mesma utilizadas no projeto LQR. O sinal de comando
é dado pela Equação (6.36).
fw (t) = −Rc BT s(t)
(7.10)
onde o vetor s(t) é função da referência, atitude (φr , θr , ψr ) e velocidades angulares
(φ̇r , θ̇r , ψ̇r ), calculada pela Equação (6.37). As matrizes obtidas no projeto LQR
tracking são:


−0.0286 0.0000
0.0004 −0.8476 0.0000
0.0000


Kc =  0.0000 −0.0286 0.0000
0.0000 −0.9180 0.0000 
0.0004
0.0000 −0.0286 0.0000
0.0000 −0.9510
e
134
(7.11)


0.1015


−0.1285


−0.2472


s=

 1.3734 


−1.8581


−3.7070
(7.12)
A Figura (7.28) mostra os ângulos de atitude durante o modo de aquisição. O resultado mostra um bom desempenho do sistema do controle, usando o método LQR
tracking. A partir de 200 segundos o satélite atingi a precisão requerida (< 1o ) para
a atitude especificada ou de referência (φr , θr , ψr ) = (−60, 70, 130). O modelo do
satélite usado nas simulações é o modelo completo, não linear.
Figura 7.28- Ângulos de atitude roll, pitch e yaw em função do tempo para a
aquisição de atitude.
A Figura (7.29) mostra o torque gerado pelas três rodas de reação, o primeiro gráfico
mostra o torque (τ bw ) gerado pelas rodas dispostas nas três direções principais de
inércia o segundo gráfico mostra o torque de acoplamento devido ao movimento
de rotação do satélite. O torque efetivo de controle é a soma dos dois torques. O
resultado mostra que mesmo para manobras de grandes ângulos, não há saturação
da quantidade de movimento das rodas (hw < 12N ms), ficando abaixo do valor
135
de saturação. O primeiro grafico mostra que não ocorre saturação do esforço ou
torque de controle (τ bw < 75mN m). Esse resultado mostra a viabilidade da utilização
das rodas especificadas para o satélite EQUARS para a aquisição de uma atitude
arbitrária.
Figura 7.29- Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas em função do tempo.
A Figura (7.30) mostra as rotações em rpm de cada uma das rodas e as velocidades
angulares ω bib e ω bob . Na manutenção da atitude de referência as velocidades das
rodas são nominalmente diferentes de zero.
Nas simulações, apesar do projeto do controlador LQR tracking se basear na planta
linear, o controle empregado sobre a planta não linear apresenta boas caracterı́sticas
de desempenho. O resultado se deve em parte as boas propriedades de robustez do
controle LQR. O controle LQR é uma alternativa atraente devido a otimalidade. A
escolha de ponderações para o estado e controle, permite ao projetista um projeto
de controle com base em informações e especificações dos atuadores e respectivas
faixas de operação. Um menor esforço para o ajuste dos ganhos Qc e Rc também é
conseguido, através do algoritmo apresentado, Equação (6.17).
136
Figura 7.30- Rotações por minuto das rodas de reação e velocidades angulares ω bib
e ω bob .
7.3.2
Controle PD
Os ganhos para o controlador PD foram obtidos usando a metodologia LQ, Seção
(6.2.1). A lei de controle PD, obtida na Seção (6.3.2), equações (6.47), (6.39) e (6.49)
podem ser reescritas como
 

φ − φr
φ̇
 


u = Kp  θ − θr  + Kd  θ̇ 
ψ̇
ψ − ψr

(7.13)
Onde as matrizes Kp e Kd são as matrizes com os ganhos proporcionais e derivativos,
respectivamente. Os ganhos do controlador PD obtidos são


−0.0286 0.0000
0.0000


Kp =  0.0000 −0.0286 0.0000 
0.0000
0.0000 −0.0286
137
(7.14)
e


−0.8476
0.0000 6.2714 · 10−5


Kd =  0.0000
−0.9180
0.0000 
8.1064 · 10−5 0.0000
−0.9510
(7.15)
A Figura (7.31) mostra os ângulos de atitude durante o modo de aquisição. O
resultado mostra um bom desempenho do sistema do controle, a partir dos ganhos Kp e Kd usados para o controlador PD. A partir de 200 segundos o satélite
já atinge a precisão requerida (< 1o ) para a atitude especificada ou de referência
(φr , θr , ψr ) = (−60, 70, 130). Novamente foi usado o modelo completo (não linear),
nas simulações.
Figura 7.31- Ângulos de Euler em função do tempo.
A Figura (7.32) mostra o torque gerado pelas três rodas de reação. O primeiro gráfico
mostra o torque (τ bw ) gerado pelas rodas dispostas nas três direções principais de
inércia. O segundo gráfico mostra o torque de acoplamento devido ao movimento
de rotação do satélite. O torque efetivo de controle é a soma dos dois torques. O
último gráfico mostra a quantidade de movimento angular das rodas. Note que as
três rodas ao longo dos eixos de roll, pitch e yaw apresentam uma quantidade de
138
movimento angular nominal diferente de zero, fornecendo rı́gidez giroscópica ao sistema durante a manutenção da atitude de referência. O resultado mostra que mesmo
para manobras de grandes ângulos não há saturação de quantidade de movimento
angular das rodas (hw < 12N ms), ficando abaixo desse valor. O primeiro grafico
mostra que não ocorre saturação do esforço de controle (τ bw < 75mN m), chegando a
um valor máximo no eixo de yaw de ∼
= 60mN m. Esse resultado mostra que as rodas
especificadas para o satélite EQUARS podem ser usadas para realizar manobras
de grandes ângulos, ou seja, ser usadas para aquisição de uma atitude arbitrária.
Figura 7.32- Torque τ bw , torque de acoplamente e quantidade de movimento angular das rodas em função do tempo.
A Figura (7.33) mostra as rotações por minuto de cada uma das rodas e as velocidades angulares ω bib e ω bob . Na manutenção da atitude de referência as velocidades
das rodas são nominalmente diferentes de zero.
O resultados obtidos mostram que é factı́vel o uso das rodas especificadas (Tabela
7.1), para o modo de aquisição do satélite EQUARS. Uma das vantagens no uso
do controlador PD em relação ao LQR é usar como realimentação apenas a atitude.
Entretanto, tem como incoveniente a derivação dos sinais. A metodologia LQ usada
para seleção/ajuste dos ganhos mostra-se muito prática para o controle PD.
139
Figura 7.33- Rotações por minuto das rodas e velocidades angulares ω bib e ω bob em
função do tempo.
140
CAPÍTULO 8
CONCLUSÃO
Nesse Capı́tulo apresentam-se as principais conclusões relacionadas aos resultados
obtidos, encerrando com sujestões para trabalhos futuros.
Nesse trabalhos apresentou-se um estudo de alternativas de sistemas de controle
(ACSs) para satélites artificiais estabilizados em três eixos. Os procedimentos adotados foram: 1) satélite utilizando rodas e bobinas como atuadores e 2) satélite
utilizando apenas bobinas magnéticas. Para cada um dos procedimentos analisou-se
a exeqüibilidade dos ACS(s) para os dois principais modos de operação estudados:
1) detumble, 2) estabilização. As leis de controle empregadas para esses modos foram
LQR, LQG e controladores baseados em energia (Bdot, realimentação de atitude e
velocidade angular). Os métodos LQR e LQG resultam em projetos lineares e os reguladores baseados em energia resultam em projetos não lineares. Apesar do projeto
linear os métodos LQR e LQG foram aplicados à dinâmica não linear, apresentando
bons resultados, no que se refere ao desempenho e robustez, também fornecendo
boas margens de ganho e de fase.
A Tabela (8.1) resume as diferentes configurações dos ACS(s) estudados e as diferentes estratégias de controle.
Tabela 8.1- Alternativas para o ACS.
Modo
ACS
Lei de controle
Detumble
bobinas
Rodas e bobinas
Bdot
Bdot
Estabilização
Rodas e bobinas
Rodas e bobinas
Bobinas
Bobinas
Bobinas
LQR
LQG
LQR
Realimentação de velocidade angular
Realimentação de atitude
Aquisição
Rodas e Bobinas PD
Rodas e Bobinas LQR tracking
A lei de controle Bdot utilizada no modo de detumble, apresentou um bom resul141
tado em termos de tempo, reduzindo a velocidade de rotação do satélite para um
nı́vel tı́pico admissı́vel (< 0.01rpm), colocando o satélite em condições de inicializar
o modo de estabilização. A partir de uma velocidade crı́tica ∼
= 2rpm(definida nesse
trabalho como a velocidade em que o ACS no modo de estabilização não é apropriado, devido a saturação dos atuadores) a estratégia de controle leva aproximadamente quatro órbitas para atingir uma faixa de velocidades angulares admissı́veis.
A estratégia de controle é factı́vel de ser empregada, tomando como base a consumo
de energia dos atuadores e valores de saturação. A análise de factibilidade do sistema de controle foi feita com base na potência especificada para satélites tı́picos do
EQUARS. A lei de controle é fácil de ser implementada e não requer informações
ou conhecimento da atitude do veı́culo, mas apenas dados dos magnetômetros. Apesar das simulações terem avaliado uma condição inicial de velocidade angular crı́tica
de ∼
= 2rpm, ela pode ser maior dependendo do veı́culo lançador. Entretanto, o controlador Bdot ainda pode ser usado, se o requisito de tempo não for muito estreito,
pois exigiria um tempo maior para a redução da velocidade angular.
Para o modo de estabilização foram comparados o uso de bobinas e rodas. Na
utilização de rodas foram analisadas três cenários: 1) com velocidade residual de
rotação nula (ω bib = [0, 0, 0]T graus/s), 2) com velocidade quase crı́tica (ω bib =
[6, 6, 6]T graus/s) e 3) com velocidade crı́tica (ω bib = [10, 10, 10]T graus/s). Os resultados mostram que para o primeiro caso o sistema de controle consegue atender
os requisitos de operação muito bem, não ocorrendo saturação na velocidade da rodas nem em esforço de controle ou torque. Para o segundo caso mesmo ocorrendo
saturação no torque o sistema de controle consegue estabilizar a atitude e realizar a manutenção ou o alinhamento com a VLHL. Para o terceiro caso o sistema
de controle desestabiliza a atitude devido a saturação nos atuadores. Uma solução
para esse problema seria ajustar os ganhos do controlador através das matrizes de
ponderação do estado e controle (Qc e Rc ), Diminuindo o esforço de controle e aumentado o tempo para estabilização que nos casos 1 e 2 estudados é relativamente
rápido (< 200s).
O método LQR desenvolvido para os dois modelos: satélite equipado com rodas e
satélite equipado com bobinas, atende as especificações do sistema de controle em
ambos os casos. As principais caracterı́sticas apresentadas pelo método são: controle
ótimo, robustez, fácil implementação, ponderações no estado e controle permitindo
ao projetista um maior sentimento fı́sico do problema. A principal desvantagem do
142
método é o esforço computacional para o caso de sistema variantes no tempo, como
é o caso do satélite equipado com bobinas. Entretanto isso pode ser contornado
estocando-se um ganho previamente calculado que, devido a quase periodicidade
do campo geomagnético, o ganho é quase periódico, com mostrado nas simulações.
Outro importante resultado é que, devido a robustez do método ele se mostra eficiente quando aplicado para o sistema não linear. No método LQR é necessário
conhecer todos os estados para a realimentação, caso não seja possı́vel medir todos
os estados podemos usar um controlador baseado em observador como é o caso do
método LQG. Para o método LQG implentado para o satélite equipado com rodas, considerou-se um modelo com incertezas estocásticas na saı́da e na dinâmica,
tornando o estudo mais realista. As medidas dos sensores foram corrompidas com
ruı́dos de distribuição gaussiana, a partir dos ajustes do filtro de Kalman obtemos
uma boa estimação do estado, como pode se ver nos resultados das simulações. Apesar da perda de robustez que a inserção do filtro pode causar os resultados mostram
um desempenho próximo do obtido pelo método LQR, o que evidencia que a seleção
das ponderações do controle sugeridas por Kristiansen (2000) é conveniente.
Os controladores baseados em energia usando realimentação de atitude e velocidade
angular, para o modo de estabilização, mostram-se eficientes, apresentando como
vantagem a simplicidade na implementação. A lei de controle utilizando apenas
realimentação de velocidade angular não garante, para outras condições iniciais, a
estabilização e o controle em três eixos de interesse, sendo necessário informações da
atitude. Entretanto ainda é uma alternativa que deve ser considerada, para os casos
em que nenhuma informação da atitude é disponı́vel para o sistema de controle. Os
resultados mostram que a lei de controle apresentada é uma boa alternativa para
o sistema de controle com realimentação de velocidade angular e atitude, obtendo
um desempenho razoável para o controle em três eixos. Em ambos os casos o uso
de bobinas apresenta uma opção atrativa devido ao baixo consumo de energia. Uma
vantagem em relação a metodologia LQR é que mesmo sob manobras de grandes
ângulos, a eficiência do sistema de controle de atitude na estabilização do sistema é
garantida.
Realizou-se o estudo de viabilidade do uso de rodas para manobras de grandes
ângulos durante o modo de aquisição. Para esse modo de operação as leis de controle empregadas foram: 1) LQR tracking e 2) controle PD. Para o controle PD
foi utilizado a metodologia LQ para o cálculo dos ganhos do controlador. As si143
mulações mostram desempenhos muito próximos dos dois controladores devido ao
uso da metodologia LQ para o controlador PD. O método LQR utiliza o estado como
realimentação e o PD apenas informações da atitude, entretanto tem a desvantagem
de necessitar derivar esses sinais. Os resultados mostram que é viável e factı́vel o uso
das rodas especificadas para manobras de aquisição. E mesmo para as manobras de
grandes ângulos os controladores que são baseados na planta linear apresentam um
bom desempenho.
Finalizando esse trabalho contribuiu para a área de dinâmica e controle de atitude
no sentido de que:
• apresenta um estudo de referenciais e parâmetros mais usados para representar
a atitude;
• apresenta o modelo matemático da dinâmica de atitude para veı́culos espaciais
contendo rodas de reação e bobinas magnéticas;
• discute e apresenta as técnicas de controle ótimo (LQR, LQR tracking é LQG)
bem como as leis de controle desenvolvindas pelo método de energia;
• desenvolve as leis de controle utilizando as técnicas referidas no ı́tem anterior
para aplicação de dois conceitos de sistemas de controle para satélites estabilizados em três eixos. Um procedimento que utiliza uma combinação de rodas
de reação e bobinas magnéticas como atuadores e outro que utiliza somente
bobinas magnéticas com atuadores;
• simula em ambiente MATLAB/SIMULINK o controle de satélites estabilizados em três eixos, utilizando diferentes leis de controle para as concepções de
sistemas de controle;
• disponibiliza um pacote de sofware desenvolvido na plataforma MATLAB/SIMULINK que poderá ser utilizado no futuro para o estudo de
dinâmica e controle de atitude.
8.1
Sugestões para Tabalhos Futuros
Referente ao modelo seria propõe-se além do modelo do torque de gradiente de
gravidade os modelos dos torques de pressão de radiação, arrasto atmosférico, dipolo
144
residual e o efeito jitter, associado a fontes de perturbação interna. Outra sugestão
seria a implementação do método LQG para o satélite equipado com bobinas e
a apresentação rigorosa das provas de estabilidade para os controladores baseados
em energia usando a teoria de Lyapunov. Sugere-se também inserir um filtro para
estimação dos dados dos magnetômetros para o projeto do ACS que utiliza bobinas.
Para uma simulação mais realista sugere-se:
• incluir o modelo da roda de reação;
• incluir o modelo do magnetômetro;
• testes de robustez através de Monte Carlo.
145
146
REFERÊNCIAS BIBLIOGRÁFICAS
Alfriend, T. K. Magnetic attitude control system for dual-spin satellite. AIAA
Journal, v. 13, n. 6, p. 817–822.
Arantes, J. G.; Fonseca, I. M. da. A comparasion between quaternions and euler
angles for satellite atittude dynamics. In: Colóquio Brasileiro de Dinâmica
Orbital, 12, 2004. Ubatuba, SP. Anais... Ubatuba, SP: INPE, 2004.
Arantes, J. G.; Fonseca, I. M. da. Three-axis attitude dynamics by using torque
coils only. In: Colóquio Brasileiro de Dinâmica Orbital, 12, 2004. Ubatuba, SP.
Anais... Ubatuba, SP: INPE, 2004.
Ardebil.
ASRI
program.
2004.
Disponı́vel
<http://www.asri.org.au/ASRI/index.xml>. Acesso em: out 2004.
em:
Athans, M. The role and use of stochastic linear quadratic gaussian problem in
control system design. IEEE Transactions on Automatic Control, v. 6, p.
529–552, 1971.
Auer, W. A. A double gimballed momentum wheel for precision three-axis
attitude control. In: AGARD Symposium Guidance and Control Panel, 37,
1983. Florence, Italy. Proceedings... Florence, Italy: AGARD, 1983.
Bang, M. J. T. H.; Choi, H. D. Large angle attitude control of spacecraft with
actuator saturation. Control engineering practice, n. 11, p. 989–997, 2003.
Beichman, C.; Neugebauer, G.; Chester, T. IRAS explanatory supplement.
2004. Disponı́vel em: <http://irsa.ipac.caltech.edu/IRASdocs/iras.html>.
Acesso em: 10 Aug 2004.
Bhat, S. P.; Dham, A. S. Controllability of spacecraft attitude under magnetic
actuation. In: Conference on Decision and Control IEEE, 42, 2003. Maui,
Hawai, USA. Proceedings... Maui, Hawai, USA: IEEE, 2003.
Bishop, R. Modern control system analysis & design using Matlab &
Simulink. Mento Park, California: Eddison Wesley, 1997.
Brown, C. D. Spacecraft mission design. Reston, Virginia: Education Series
AIAA, 1998. 182 p.
147
Brown, R. G. Introduction to randon signals and applied Kalman filtering.
New York: John Wiley & Sons, 1997.
Bryson, A. E. Control of spacecraft and aircraft. New jersey: Princeton
University Press, 1994. 378 p.
Buckingham, O. V. A.; Smirnov, G. V. Magnetic torques for momentum
desaturation of space station control moment gyros. Journal Spacecraft,
v. 9, n. 6, p. 324–330, 1972.
Bushenkov, M. Y. O. A. V.; Smirnov, V. G. Attitude stabilization of a satellite by
magnetic coils. Acta Astronautica, v. 50, p. 721–728, 2002.
Carrara, V. Modelagem das forças e torques atuantes em satélites. 1982.
153 p. (INPE–2454–TDL/094). Dissertação (Mestrado em Ciência Espacial),
INPE, São José dos Campos. 1982.
Carvalho, H. de C. EQUARS mission analysis. São José dos Campos, INPE,
2003. 40 p. (E2000-TRP-001v00).
Chobotov, V. A. Orbital mechanics. Reston, Virginia: Education Series AIAA,
1996. 447 p.
Clausen, L. T. The Orsted satellite project. 2004. Disponı́vel em:
<http://web.dmi.dk/fsweb/projects/oersted/homepage.html>. Acesso em: jan
2004.
Cohen, V. D. Attitude dynamics of an orbiting electromagnet. Journal
Spacecraft, v. 11, p. 252–256, 1973.
Distefano, A. R. S. J. J.; Willians, I. J. Sistemas de retroação de controle. São
Paulo: McGRAW-HILL, 1972. 478 p.
Dorato, C. A. P.; Cerone, V. Linear quadratic control an introdution.
Englewood Cliffs: New Jersey: Prentice Hall, 1995. 205 p.
Doyle, J. C. Guaranteed margins for lqg regulators. IEEE Transactions on
Automatic Control, v. 4, n. 23, p. 756–757, 1978.
El-Gohary, A. Optimal stabilization of a rigid body motion using rotors system.
Applied Mathematics and Computation, v. 136, p. 229–239, 2003.
148
Escobal, P. R. Methods of orbit determination. New York: John Wiley &
Sons, 1965.
Fauske, K. M. NCUBE attitude control. Norwegian, Trondhein: Departament
of Enginneering Cybernetics, NTNU, dec 2002.
Fauske, K. M. Attitude stabilization of an underactuated rigid spacecraft.
2003. 55 p. SIV.ING Thesis, (Departament of Engineering Cybernetics).
Norwegian University of Technology and Science, Trondheim, Norwegian. Jan
2003.
Fichter, M. S. W.; Zentcraf, P. Control design for generalized normal mode
operation of bias momentum satellites. Control Engeneering Practice,
v. 4, p. 1355–1360, 1996.
Flora, A. L. Projeto de um sistema de controle de atitude de um satélite
com apêndices flexı́veis pelos métodos LQG/LTR e H-Infinito. 1995.
223 p. Dissertação (Mestrado em Ciência Espacial), INPE, São José dos
Campos.
Gökçev, P. T. M. C.; Meerkov, S. M. An lqr/lqg theory for system with saturation
actuators. Transactions on Automatic Control, v. 46, n. 10, p. 1529–1542,
2001.
Grassi, S. V. M.; Moccia, A. Preliminary design of the attitude control system of a
microsatellite for earth observation. Space Technological, Biarritz, France,
v. 15, p. 223–230, 1995.
Gurman, J. B. The SOHO solar cycle mission. 2004. Disponı́vel em:
<http://sohowww.nascom.nasa.gov/publications>. Acesso em: sep 2004.
Hamzah, N.; Hashida, Y. Tiungsat-1 momentum wheel commissioning. In: World
Engineering Congress, 1, 1999. Kuala Lumpur. Proceedings... Kuala
Lumpur: World Eng. Cong., 1999.
Hanselman, D.; Littlefield, B. Matlab 6 curso completo. São Paulo: Prentice
Hall, 2003. 676 p.
Hanselman, D. C.; Benjamin, C. K. MATLAB tools for control system
analysis and design. Englewood Cliffs: New Jersey: Princeton-Hall, 1994.
149
Hegrenes, O. Attitude control by means of explicit model predictive
control, via multi-parametric quadratic programming. 2004. 120 p.
Thesis (Engineering Cybernetics), Trondheim Norwegian. 2004.
Heidelberg. Computers, displays and space products. 2004. Disponı́vel em:
<http://www.teldix.de/>. Acesso em: ago 2004.
Hughes, P. C. Spacecraft attitude dyanmics. New York: John Wiley & Sons,
1986. 564 p.
Junkins, J. L.; Carrington, C. K. Time optimal magnetic attitude maneuvers. In:
IAA/AAS Astrodynamics conference, 11, 1980. USA. Proceedings... USA:
IAA/ASS, 1980.
Junkins, J. L.; Turner, J. D. Optimal spacecraft rotational maneuvers. New
York: Elsevier Science Publishers B. V., 1986. 515 p.
Kaplan, M. H. Modern spacecraft dynamic & control. New York: John Wiley
& Sons, 1976. 415 p.
Khalil, H. K. Nonlinear system. [S.l.: s.n.], 2000.
Kim, H. L. B. J.; Choi, S. D. Three axis reaction wheel attitude control
system for KITSAT-3 microsatellite. Taejon, Korea: Sattelite Technology
Reserch Center, KAIST, 1999.
Kirk, D. E. Optimal control theory: an introduction. Englewood Cliffs: New
Jersey: Princeton Hall, 1970.
Kristiansen. Norwegiann micro salellite.
Cybernetics), Trondheim Norwegian. 2000.
2000.
Thesis
(Engineering
Kuga, H. K.; Guedes, U. T. V. Dinâmica de atitude para satélites
estabilizados por rotação. São José dos Campos, SP, 1987.
(INPE-4403-TVTE/275).
Kuga, H. K.; Kondapalli, R. R. Introdução à mecânica orbital. São José dos
Campos, 1995. 73 p. (INPE-5615-PUD/064).
Kuga L. D., F. H. K.; Guedes, U. T. V. Simulação de atitude de manobras
para o satélite brasileiro estabilizado por rotação. São José dos
Campos, SP, 1987. (INPE-4271-PRE/1143).
150
Kwakernaak, H.; Sivan, R. Linear optimal control system. New York: John
Wiley & Sons, 1972. 564 p.
Larson, W. J.; Wertz, J. R. Space mission analysis and design. Torance,
California: Space Technology Series, 1992. 865 p.
Legg, V. E. Survey of magnetic material and applications in the
telephone
system.
2003.
Disponı́vel
em:
<http://www.telepsystem.com/reporter/index.xml>. Acesso em: dez 2003.
Levine, S. W.; Leonard, N. E. Using Matlab to analyse and control system.
New York: Addison Wesley Publishing, 1995.
Maciejowski, J. M. Multivariable feedback design. New York: Addison Wesley
Publishing, 1989.
Macmillan, S.; Quinn, J. M. The derivation of world magnetic model 2000.
London: British Geological Survey, 2000. 278 p. (Geomagnetism Series –
WM/00/17R).
Marteau, P. P. F.; Psiaki, M. Active magnetic control system for gravity gradient
stabilized spacecraft. In: Annual AIAA/USU conference for small satellites, 2,
1988. Logan (Utah), USA. Proceedings... Logan (Utah), USA: AIAA, 1988.
Marteau, S. B. G. F.; Rogers, E. Attitude determination and control for small
spacecraft. In: UKACC International Conference on Control, 13, 1996. U.K.
Proceedings... U.K.: IEE, 1996.
Martins Neto, A. F. Atitude e seu controle. In: Kuga, A. F. B. de A. P. . H. K.
(Ed.). Fundamentos de tecnologia espacial. São José dos Campos, INPE:
INPE, 2001. v. 1, p. 65–79.
Mclean, S. The world magnetic model and associated software.
Washington, D.C., National Oceanic & Atmospheric Administration (NOAA).
Disponı́vel em: <http://www.ngdc.noaa.gov>. Acesso em: out 2004.
Moore, J. B.; Anderson, D. O. B. Optimal control linear quadratic methods.
Englewood Cliffs: New Jersey: Princeton Hall, 1990. 380 p.
Moscati, N. R. Projeto de um sistema de controle de atitude (três eixos)
de satélites utilizando a metodologia LQG/LTR. 1992. 218 p.
151
(INPE–5473–TDI/504). Dissertação (Mestrado em Ciência Espacial), INPE,
São José dos Campos. 1992.
Musser, L. K.; Ebert, L. W. Autonomous spacecraft attitude control using
magnetic torquing only. In: AIAA Guidance, Navigation, and Control
Conference, 12, 1989. NASA Goddard Space Flight Center, Greenbelt.
Proceedings... NASA Goddard Space Flight Center, Greenbelt: NASA,
1989. p. 23–38.
Ogata, K. Engenharia de controle moderno. Rio de Janeiro: LTC - Livros
Técnicos e Cientı́ficos Editora, 1998. 812 p.
Overby, E. J. Attitude control for the norwegian student satellite nCube.
2004. 74 p. Thesis (Engineering Cybernetics), Trondheim Norwegian. 2004.
Perkel, H. Stabilite - three axis attitude control system utilizing a single reaction
wheel. In: AIAA Comunication Satellite System Conference, 1, 1966.
Washington, D.C. Proceedings... Washington, D.C.: AIAA, 1966. p.
375–400.
Pilchowski, C. C. W. S. H. U.; Ferreira, D. D. L. Introdução à mecânica
celeste. São José dos Campos, INPE, 1981. (INPE-COM.4/RPE
C.D.U.:521.3).
Pilchowski, H. U. Sensores de atuadores. In: Kuga, A. F. B. de A. P. . H. K. (Ed.).
Fundamentos de tecnologia espacial. São José dos Campos: INPE, 2001.
v. 1, p. 50–64.
Psiaki, L. M. Magnetic torquer attitude control via asymptotic periodic linear
quadratic regulation. Journal of Guidance, Control, and Dynamics,
v. 24, p. 386–394, 2001.
Quirelli, I. M. P. Spin stabilized satellite attitude propagation. 2002. Master
Thesis, UNESP, Guaratinguetá, SP. 2002.
Rodrigues, D. S. S.; Zanardi, M. C. Spacecraft attitude propagation with different
representation. In: Kuga, A. F. B. de A. P. . H. K. (Ed.). Advances in space
dynamics 4: celestial mechanics and astronautics. São José dos Campos:
INPE, 2004. v. 1, p. 143–150.
152
Roma, A. M. Análise dinâmica e controle de um satélite artificial com
paı́neis solares flexı́veis. 1991. 177 p. (INPE–5220–TDL/436). Dissertação
(Mestrado em Ciência Espacial), São José dos Campos. 1991.
Shigehara, M. Geomagnetic attitude control of an axisymemtric spinning satellite.
Journal Spacecraft, v. 9, n. 6, p. 623–635, 1972.
Silani, E.; Lovera, M. Magnetic spacecraft attitude control survey and some new
results. Control Engineering Practice, v. 1, p. 1, 2003.
Souza, L. C. G. Controle de atitude de um satélite artificial através da
extensão da teoria do regulador linear quadrático. 1987. 60 p.
(INPE–4407–TDL/304). Dissertação (Mestrado em Ciência Espacial), INPE,
São José dos Campos. 1987.
Souza, M. L. O. Estudo e desenvolvimento de um sistema de controle de
atitude ativo em três eixos para satélites artificiais usando atuadores
pneumáticos
a
gás
frio
e
volantes
a
reação.
1981.
(INPE–2000–TDL/042). Dissertação (Mestrado em Ciência Espacial), INPE,
São José dos Campos. 1981.
Souza, P. N. Análise, projeto, construção e testes de um modelo de roda
de reação para aplicações espaciais. 1987. 185 p. (INPE–4358–TDL/299).
Dissertação (Mestrado em Ciência Espacial), INPE, São José dos Campos.
1987.
Spencer, T. M. Automatic magnetic control of a momentum-biased observatory in
equatorial orbit. Journal Spacecraft, v. 14, n. 4, p. 211–218, 1977.
Spindler, K. Single attitude control laws for momentum-wheel actuated spacecraft.
In: International Symposium Space Dynamics, 14, 2000. Biarritz, France.
Proceedings... Biarritz, France: AIAA, 2000.
Tamura, T. ASCA Measurements of the Gravitational Potential Profile in
the Central Region of Galaxy Clusters. Tese (Ph.D Thesis) — Graduate
School of Science, University of Tokyo, Tokyo. 1998.
Trivelato, G. da C. Controle de rodas de reação através de técnicas digitais
usando modelos de referência. 1988. 209 p. (INPE–4618–TDL/335).
Dissertação (Mestrado em Ciência Espacial), INPE, São José dos Campos.
1988.
153
Ulrich, H. Satélites Artificiais movimento de atitude. Notas de aula. 2004.
Varatharajoo, R.; Fasoulas, S. The combined energy and attitude control
system for small satellites earth observation missions. Boston, USA:
Sattelite Technology Reserch Center, 1975.
Wang, P.; Shtessel, B. Y. B. Satellite attitude control using only magnetictorques.
In: AIAA Guidance, Navigation, and Control Conference, 12, 1998. Boston.
Proceedings... Boston: AIAA, 1998.
Wells, L. S. J. G.; Jeans, T. Canada’s smallest satellite: The canadian advanced
nanospace experiment (canx-1). In: Annual AIAA/USU Conference on Small
Satellites, 16, 2002. Logan, Utah. Proceedings... Logan, Utah, 2002.
Wertz, J. R. Spacecraft attitude determination and control. London,
England: D. Reideil Publishing Company, 1978. 861 p.
Whitford, C.; Forrest, D. The catsat attitude control system. In: AIAA/USU
Conference on Small Satellites, 12, 1998. USA. Proceedings... USA: AIAA,
1998.
Wie, B. Space vehicle dynamics and control. Reston, Virginia: AIAA
Education Series, 1998. 661 p.
Wie, H. W. B.; Arapostathis, A. Quaternion feedback regulator for spacecraft
eigenaxis rotaions. Journal of Guidance, Control and Dynamics, v. 12,
n. 3, p. 375–380, 1989.
Wie, P. M. B. B. Quaternion feedback for spacecraft large angle maneuvers.
Journal of Guidance, Control and Dynamics, v. 3, n. 3, p. 360–365,
1985.
Wilson, D. Attitude and Orbit Control Using the Spacecraft Control
Toolbox v4.6. Princeton - New Jersey, Dec 2000.
Wisniewski, R. Satellite attitude control using only electromagnetic
actuation. Tese (Ph.D Thesis) — Aalborg University, Departament of
Control Engineering, Danish. 1996.
Wisniewski, R. Linear time varying approach to satellite control using only
electromagnetic actuation. AIAA Guidance, navigation control, v. 11, p.
243–251, 1997.
154
Wisniewski, R.; Blanke, M. Fully magnetic attitude control for spacecraft subject
to gravity gradient. Automatica, v. 35, p. 1201–1214, 1999.
Wisniewski, R.; Markley, F. L. Optimal magnetic attitude control. In: IAFC
World Congress, 14, 1999. USA. Proceedings... USA: IAFC, 1999. p.
313–318.
Wright, P. S.; Wong, H. S. An overview of sensors in spacecraft engineering.
New York, USA: Sattelite Technology Reserch Center, 1989.
Yairi, T. On-board reconfigurable attitude control system with optimazation. In:
International Symposium on Space Technology and Science,19,1994.
Yokohama, Japan. Proceedings... Yokohama, Japan: ISTS, 1994.
Zanardi, M. C.; Assis, S. C. de; Kuga, H. K. Torque residual médio com modelo
de quadripolo. In: Congresso Temático de Dinâmica, Controle e Aplicaçòes, 3,
2004. Ilha Solteira - SP. Anais... Ilha Solteira - SP: Série Arquimedes, 2004.
Zanardi, M. C.; Quirelli, I. M. P.; Kuga, H. K. Analytical attitude propagation of
the spin stabilized earth artificial satellite. In: International Symposium of
Space Flight Dynamics, 17, 2003. Moscou - Rússia. Proceedings... Moscou Rússia: CD-ROM, 2003.
Zanardi, M. C.; Quirelli, I. M. P.; Kuga, H. K. Torques magnéticos: Aplicações à
satélites estabilizados por rotação. In: Congresso Temático de Dinâmica,
Controle e Aplicaçòes, 2, 2003. São José dos Campos - SP. Anais... São José
dos Campos - SP: Série Arquimedes, 2003. p. 3167–3176.
155
156
APÊNDICE A
CÁLCULO DO GRADIENTE DE GRAVIDADE
No desenvolvimento do modelo do torque devido ao gradiente de gravidade
consideram-se as seguintes aproximações, válidas para satélites em órbita da Terra
em sua maioria.
• somente o campo gravitacional da Terra é considerado;
• a Terra é considerada esférica com distribuição de massa uniforme;
• as dimensões do satélite são pequenas em comparação com a distância Terra
- satélite;
• o satélite é um corpo simples .
A partir das aproximações 1, 2 e 3 o torque devido ao gradiente de gravidade pode
ser expresso por
Z
τ g = −µ
sat
r×R
dm
R3
(A.1)
onde µ é a constante geo-gravitacional, R é o vetor distância do centro da Terra
ao elemento de massa dm do satélite e r é o vetor distância do centro de massa do
satélite ao elemento de massa dm. Tem-se que
R = |Rs + r|
(A.2)
onde Rs é o vetor distância do centro de massa do satélite ao centro de massa da
Terra.
Substituindo a Equação (A.2) em (A.1), expandindo o termo |Rs + r|−3 em série
2
|r|
de Taylor e desprezando os termos de ordem igual e superior a dois, O |R
,
s|
|r|
admissı́vel com base na aproximação 4 ( |R
<< 1) obtemos
s|
157
µ
τg = 5
Rs
Z
(Rs · r) (r × Rs ) dm
(A.3)
sat
Manipulando a Equação (A.3) como segue (Wie, 1998)
τg
Z
µ
= −3 5 Rs ×
r (r · Rs ) dm
Rs
Zsat
µ
= −3 5 Rs ×
rrdm · Rs
Rs
sat
(A.4)
(A.5)
Usando a notação dyadic, pormenorizadas em Hughes (1986) e Wie (1998), o tensor
de inércia pode ser escrito como
Z
J=
r2 Î − rr
(A.6)
sat
onde Î = x̂x̂ + ŷ ŷ + ẑ ẑ é o dyadic unitário, e x̂, ŷ e ẑ são os vetores unitários/versores
que formam a base do sistema do satélite.
Substituindo (A.6) na Equação (A.5), temos
τg
Z
µ
2
= −3 5 Rs ×
r Î − J · Rs
Rs
sat
Z
µ
µ
= −3 5 Rs ×
r2 Îdm · Rs + 3 5 Rs × J · Rs
Rs
Rs
sat
(A.7)
(A.8)
de (A.6) e da relação Î · Rs = Rs obtemos
τg = 3
µ
Rs × J · Rs
Rs5
158
(A.9)
O torque devido ao gradiente de gravidade pode ser expresso por
Rs
µ Rs
×
J
·
Rs3 Rs
Rs
µ
= 3 3 R̂s × J · R̂s
Rs
τg = 3
(A.10)
(A.11)
onde R̂s é o vetor unitário na direção Rs . Escrevendo o torque de gradiente de
gravidade no referencial orbital ou VLHL, resulta
τ og = 3
µ
ẑo × J · ẑo
Rs3
(A.12)
onde ẑo é o vetor intário na direção nadir/vertical local.
A velocidade orbital média é dada por
r
ωo =
µ
Rs3
(A.13)
Substituindo na Equação (A.12) obtemos
τ og = 3ωo2 ẑo × J · ẑo
(A.14)
Usando a notação do operador anti-simétrico para o produto vetorial, reescrevemos
a Equação (A.14)
τ og = 3ωo2 S (ẑo ) J · ẑo
Expressando o vetor unitário ẑo no referencial do satélite (x̂, ŷ, ẑ) temos
159
(A.15)
ẑo = cb3x x̂ + cb3y ŷ + cb3z ẑ
(A.16)
onde cb3i , (i = x, y, z) são os cosenos diretores do vetor initário ẑo nas direções x̂, ŷ, ẑ.
Desenvolvendo a Equação (A.15)
  b 
c3x
0
−cb3z cb3y
  b 
b
2 b
b
τ g = 3ωo  c3z
0
−c3x  J c3y 
−cb3y cb3x
0
cb3z

(A.17)
Como os eixos do sistema de referência do satélite coincide com os eixo principais
de inércia, temos que os produtos de inércia são zeros, temos que a Equação (A.17)
resulta
 b 
c3x Jx
0
−cb3z cb3y
 b 
b
2 b
b
τ g = 3ωo  c3z
0
−c3x  c3y Jy 
b
b
0
cb3z Jz
−c3y c3x

(A.18)
Desenvolvendo a Equação (A.18) obtemos o modelo do torque de gradiente de gravidade expresso no referencial do satélite


(Jz − Jy ) cb3y cb3z


τ bg = 3ωo2 (Jx − Jz ) cb3x cb3z 
(Jy − Jx ) cb3x cb3y
(A.19)
Podemos escrever o gradiente de gravidade em termos dos ângulos de Euler/de
atitude; roll (ψ), pitch (θ) e yaw (φ) como segue
Para a seqüência de rotação 3 − 2 − 1 os cosenos diretores cb3 são dados por (Wertz,
1978)
160
cb3x = − sin θ
(A.20)
cb3y = sin φ cos ψ
(A.21)
cb3z = cos φ cos ψ
(A.22)
Substituindo na Equação (A.19) obtemos o modelo do torque de gradiente de gravidade em função dos ângulos de Euler


(Jz − Jy ) sin φ cos φ cos θ2


τ bg = 3ωo2  (Jx − Jz ) cos φ cos θ sin θ 
(Jy − Jx ) sin φ sin θ cos θ
(A.23)
Usando a Equação (4.11) do Capı́tulo (4), Seção (4.2.2), podemos expressar o gradiente de gravidade em função dos parâmetos de Euler/quaternions como

(Jz − Jy ) (1 η + 2 3 ) (1 − 221 − 222 )


τ bg = 6ωo2 (Jx − Jz ) (1 3 − 2 η) (1 − 221 − 222 )
2 (Jy − Jx ) (1 3 − 2 η) (1 η + 2 3 )

(A.24)
Linearização da Expressão do Torque de Gradiente de Gravidade
No modo de operação onde são realizadas manobras de pequenos ângulos, é razoável
aproximar as equações por equações lineares, válidas em torno de um ponto de
operação. Considerando θ, φ e ψ menores que até 15 graus, podemos fazer as seguintes aproximações sin φ ≡ φ, sin θ ≡ θ, cos θ ≡ 1 e θφ ≡ 0. Fazendo essas
aproximações, o modelo do torque de gradiente de gravidade, dado pela Equação
(A.23) é aproximado por (Wie, 1998, Kaplan, 1976)


(Jz − Jy ) φ


τ bg = 3ωo2  (Jx − Jz ) θ 
0
161
(A.25)
Em termos dos parâmetros de Euler o modelo do torque de gradiente de gravidade
pode ser aproximado por


(Jz − Jy ) 21


τ bg = 3ωo2 (Jx − Jz ) 22 
0
162
(A.26)
APÊNDICE B
PROPAGAÇÃO DA ÓRBITA E TRANSFORMAÇÕES
Esse apêndice descreve de maneira detalhada a teoria necessária para o cálculo da
propagação da órbita Kepleriana (modelo de dois corpos). São descritas ainda as
transformações entre os diferentes referênciais envolvidos e o cálculo da longitude,
latitude e altura do ponto sub-sátelite. Todas as transformações e o cálculo da órbita
são implementados em MATLAB.
Para a propagação do campo magnético terrestre (modelo IGRF) é necessário incluir:
• cálculo/propagação da órbita (Kepleriana);
• matriz de transformação do sistema inercial para o sistema do satélite;
• obtenção da latitude, longitude e a altura do ponto sub-satélite;
A partir desses cálculos podemos obter as componentes do campo magnético terrestre, espressas no referencial do satélite.
Nota: O modelo IGRF, disponı́vel nas linguagens C e FORTRAN foi convertido para ambiente MATLAB e SIMULINK, usando S-functions. Detalhes
desse procedimento é encotrado na documentação da mathworks, disponı́vel em
www.mathworks.com.
B.1
B.1.1
Cálculo da Órbita
Posicionamento de Satélites - Problema Direto
O problemo direto consiste em dados os elementos Keplerianos (a, e, i, ω, Ω, M ) determinar à posição do satélite em coordenadas cartesianas X = [Xi Zi Zi ]T (Kuga
e Kondapalli, 1995). Onde X é o vetor posição no referencial inercial (ECI).
163
Figura B.1- Referenciais (inercial, da órbita) e os elementos Keplerianos (i, ω, Ω).
Pela Figura (B.1) pode-se notar que os elementos Keplerianos que representam os
ângulos de Euler são:
• Ω é a ascensão do nodo ascendente 0o ≤ Ω ≤ 360o ;
• i é a inclinação da órbita em relação ao equador −90o ≤ i ≤ 90o ;
• ω é o argumento do perigeu 0o ≤ ω ≤ 360o .
Esses elementos Keplerianos definem a orientação da órbita. Os elementos Keplerianos (a, e) definem o tamanho e tipo (elı́ptica, circular, hiperbólica) da órbita.
Sistema da órbita (OOF): O sistema da órbita definido (Xo , Yo , Zo ) é um sistema
de coordenadas com origem no centro de massa da Terra. O eixo Xo aponta na
direção do perigeu Π, o eixo Zo aponta na direção normal ao plano da órbita. O eixo
Yo é obtido pela regra da mão direita.
B.1.2
Equação de Kepler
As coordenadas do satélite em relação ao sistema da órbita são dadas por
164
Xo = a (cos u − e)
√
Yo = a 1 − e2 (sin u)
(B.1)
Zo = 0
(B.3)
(B.2)
onde a é o semi-eixo maior, u é a anomalia excêntrica e e é a excentricidade da
órbita.
A Equação de Kepler é dada por (Brown, 1998)
M = u − e sin u
(B.4)
onde M é a anomalia média dada por
M = n (t − to )
(B.5)
n é o movimento médio, que é dado por
r
µ
a3
(B.6)
µ = GM
(B.7)
n=
e
onde G é a constante gravitacional e M a massa do corpo central. Para a Terra
3
temos µ = 398600.4 km
s2
165
B.1.3
Matriz de Rotação
Dados os elementos orbitais (i, Ω, ω) temos que a matriz de transformação/rotação
que fornece as coordenadas do sistema inercial em função das coordenadas do sistema
da órbita é dada por (Kuga e Kondapalli, 1995).
R (i, Ω, ω) = RZo (−Ω) RXo (−i) RZo (−ω)
(B.8)
Sendo a relação/transformação dos sistemas expressa por
 
 
Xi
Xo
 
 
 Yi  = R (i, Ω, ω)  Yo 
Zi
Zo
(B.9)
Temos, para cada seqüência de rotações, as matrizes


cos ω − sin ω 0


RZo (−ω) =  sin ω cos ω 0
0
0
1
(B.10)


cos Ω − sin Ω 0


RZo (−Ω) =  sin Ω cos Ω 0
0
0
1
(B.11)


1 0
0


RXo (−i) = 0 cos i − sin i
0 sin i cos i
(B.12)
Substituindo na Equação (B.8) obtemos.
166


cos Ω cos ω − sin Ω cos i sin ω − cos Ω sin ω − sin Ω cos i cos ω sin Ω sin i


R (i, Ω, ω) = sin Ω cos ω + cos Ω cos i sin ω − sin Ω sin ω + cos Ω cos i cos ω − cos Ω sin i
sin i sin ω
sin i cos ω
cos i
(B.13)
Portanto, de (B.1), (B.2), (B.3) e (B.13) em (B.9), obtemos as coordendas do satélite
no sistema inercial.
  

Xi
(cos Ω cos ω − sin Ω cos i sin ω) Xo − (cos Ω sin ω − sin Ω cos i cos ω) Yo
  

 Yi  = (sin Ω cos ω + cos Ω cos i sin ω) Xo − (sin Ω sin ω + cos Ω cos i cos ω) Yo 
Zi
sin i sin ωXo + sin i cos ωYo
(B.14)
Órbita Circular
Para uma órbita circular temos que:
a = R⊗ + h
e ∼
= 0
(B.15)
ω = 0
(B.17)
(B.16)
Note que a referência para o inı́cio do cálculo da órbita circular é tomado sobre o
nodo. R⊗ é o raio médio da Terra e h a altitude do satélite. Da Equação de Kepler
(B.4) resulta.
M =u=f
logo de (B.5) temos
167
(B.18)
u = n (t − to )
(B.19)
onde n = ωo = velocidade orbital constante.
As equações (B.1), (B.2) e (B.3) se reduzem a
Xo = a cos u
(B.20)
Yo = a sin u
(B.21)
Zo = 0
(B.22)
e a matriz de transformação (B.13) se reduz a


cos Ω − sin Ω cos i sin Ω sin i


R (i, Ω, ω) =  sin Ω cos Ω cos i − cos Ω sin i
0
sin i
cos i
(B.23)
Substituindo na Equação (B.9), obtemos as coordenadas dos satélite no sistema
inercial para o caso de uma órbita circular.
  

Xi
a cos Ω cos u − a sin Ω cos i sin u
  

 Yi  = a sin Ω cos u + a cos Ω cos i sin u
Zi
a sin i sin u
(B.24)
A Tabela (B.1) apresenta o programa em MATLAB para a propagação da órbita
circular.
168
Tabela B.1- Programa MATLAB para o cálculo da órbita circular.
function [X] = SimOrbitaECI(par)
% by Gilberto Arantes 2004 INPE
%
% [X] = SimOrbitaECI(par) Propaga uma orbita circular
%
% Calculo da posicao do satelite no referencial ECI
%
% Input
% par = [ h i asc nodo u]
% h = altura do satelite (km)
% i = inclinacao da orbita (graus)
% asc nodo = ascensao do nodo ascendente (graus)
% u = anomalia excentrica (rad)
%
% Output
% X = vetor posicao no referencial inercial (km)
%
h = par(1);
i = par(2);
asc nodo = par(3);
u = par(4);
r = 6378 + h; % a = distancia do centro da Terra ao satelite em km
% Calculo da matriz de tranformacao do sistema orbital para o inercial (ref.: Kuga, 1995)
fator = pi/180;
a = 0; % Orbita circular: argumento do perigeu e zero
b = fator*i;
c = fator*asc nodo;
R = [(cos(c)*cos(a)-sin(c)*cos(b)*sin(a)) (-cos(c)*sin(a)-sin(c)*cos(b)*cos(a)) (sin(c)*sin(b))
(sin(c)*cos(a)+cos(c)*cos(b)*sin(a)) (-sin(c)*sin(a)+cos(c)*cos(b)*cos(a)) (-cos(c)*sin(b))
sin(b)*sin(a)
sin(b)*cos(a)
cos(b) ];
% Posicao no referencial orbital
xo = r*cos(u);
yo = r*sin(u);
zo = 0;
Xo=[xo;yo;zo]; % Coordenadas no Sistema da Orbita
% Posicao no referencial inercial ECI
X = R*Xo; % X = X [X Y Z]’
% end
B.2
Matriz de Rotação do Sistema Inercial - Sistema do Sátelite
O cálculo da matriz de transformação do sistema inercial ECI (definido em 4.1.1)
para o sistema do satélite é mostrado nessa Seção. As etapas do procedimento para
o cálculo da matriz de transformação são:
• obtemos a matriz de transformação do referencial inercial (ECI) para o pseudo
169
OF
orbital (POF) (RPECI
);
• obtemos a matriz de transformação do referencial pseudo orbital (POF) para
o orbital (OF)(ROF
P OF );
• obtemos a matriz de transformação do referencial orbital (OF) para o sistema
do satélite (BF) (RBF
OF ).
Obtida cada uma das matrizes de rotação, temos que a matriz de transformação do
sistema inercial para o sistema do satélite será dada por
BF OF
P OF
RBF
ECI = ROF RP OF RECI
B.2.1
(B.25)
OF
Matriz RPECI
Referencial pseudo orbital: O referencial pseudo orbital (POF) definido por
0
0
0
Xo , Yo , Zo é um sistema de coordenadas com origem no centro da Terra. O eixo
0
0
Xo aponta na direção do centro de massa do satélite, o eixo Zo aponta na direção
0
normal a plano da órbita. O eixo Yo é obtido pela regra da mão direita (Ulrich,
2004). Note que o referencial pseudo orbital é um referencial girante, assim como o
referencial orbital, definido em 4.1.2. A Figura (B.2) ilustra o referido referencial.
Figura B.2- Referencial inercial (ECI) e pseudo orbital (POF).
170
A partir da Figura (B.2) podemos ver que a seqüência de rotações que leva o sistema
inercial para o sistema pseudo orbital é dada por RZi (u) ← RXi (i) ← RZi (Ω), ou
seja, primeira rotação de um ângulo Ω em torno do eixo Zi , segunda rotação de um
ângulo i em torno do eixo Xi e terceira rotação de um ângulo u em torno do eixo
Zi , essa seqüencia é válida para uma órbita circular.
Temos para cada seqüência de rotação


cos u sin u 0


RZi (u) = − sin u cos u 0
0
0
1
(B.26)


cos Ω sin Ω 0


RZi (Ω) = − sin Ω cos Ω 0
0
0
1
(B.27)


1
0
0


RXi (i) = 0 cos i sin i 
0 − sin i cos i
(B.28)
Temos que a matriz de transformação é dada por
R (i, Ω, u) = RZi (u) RXi (i) RZi (Ω)
(B.29)
 0
 
Xo
Xi
 0
 
 Yo  = R (i, Ω, u)  Yi 
0
Zo
Zi
(B.30)
onde
Substituindo a Equação (B.26), (B.27) e (B.28) em (B.29), obtemos a matriz de
171
rotação que leva o sistema inercial (ECI) para o sistema pseudo orbital (POF)

OF
RPECI
B.2.2

cos u cos Ω − cos i sin u sin Ω cos Ω cos i sin u + cos u cos Ω sin i sin u


= − cos Ω sin u − cos i cos u cos Ω cos Ω cos i cos u − sin u sin Ω sin i cos u
sin i sin Ω
− cos Ω sin i
cos i
(B.31)
Matriz ROF
P OF
O referencial orbital foi definido em (4.1.2), a Figura (B.3) ilustra a orientação do
referencial orbital em relação ao referencial pseudo orbital, note que o referencial
pseudo orbital é representado no centro de massa do satélite.
Figura B.3- Orientação do referencial orbital (OF) e do referencial pseudo orbital
(POF).
A seqüência de rotações que descreve o referencial pseudo orbital em relação ao
referencial orbital é RXo0 (−90o ) ← RZo0 (90o ), onde


0 1 0


RZo0 (90o ) = −1 0 0
0 0 1
172
(B.32)


1 0 0


RXo0 (−90o ) = 0 0 −1
0 1 0
(B.33)
Temos que a matriz de transformação é dada por
o
o
ROF
P OF = RXo0 (−90 ) RZo0 (90 )
(B.34)
onde
 
xo
 
OF
 yo  = RP OF
zo
 0
Xo
 0
 Yo 
0
Zo
(B.35)
Substituindo as equações (B.32) e (B.33) em (B.34), obtemos a matriz de transformação do referencial pseudo orbital para o referencial orbital

ROF
P OF
B.2.3

0 1 0


=  0 0 −1
−1 0 0
(B.36)
Matriz RBF
OF
A matriz de rotação do sistema orbital para o sistema do corpo é obtida a partir
dos ângulos de atitude ou ângulos de Euler. A seqüência de rotações escolhida foi
3 − 2 − 1, cuja matriz de atitude é dada por (Wertz, 1978)
173

RBF
OF

cos θ cos ψ
cos θ sin ψ
− sin θ


= − cos φ sin ψ + sin φ sin θ cos ψ cos ψ cos φ + sin ψ sin θ sin φ sin φ cos θ 
sin ψ sin φ + cos φ sin θ cos ψ − sin φ cos ψ + cos φ sin θ sin ψ cos φ cos θ
(B.37)
Matriz RBF
ECI
Substituindo as equações (B.31), (B.36) e (B.37) na Equação (B.25), e com o auxı́lio
do manipulador simbólico do programa MATHEMATICA, obtemos finalmente a
matriz de transformação do sitema inercial para o sistema do satélite.
RBF
ECI


c11 c12 c13


= c21 c22 c23 
c31 c32 c33
(B.38)
com
c11 = −cθsisψsΩ − cθcψ (cΩsu − cicusΩ) − sθ (cisusΩ − cucΩ)
c12 = cθsisψcΩ + cθcψ (−sΩsu − cicucΩ) − sθ (−cisucΩ − cusΩ)
c13 = cucθsψsi + sisusθ − cicθsψ
(B.39)
que são os cosenos diretores de Xi em relação aos eixos do satélite (x, y, z)
c21 = −sisΩ (cφ + sθsφsψ) + (cψsθsφ − cφsψ) (−cΩsu − cicusΩ) + cθsφ (−cucΩ + cisusΩ)
c22 = sicΩ (cφ + sθsφsψ) + (cψsθsφ − cφsψ) (−sΩsu + cicucΩ) + cθsφ (−cusΩ − cisucΩ)
c23 = −cθsisusφ + cusi (cψsθsφ − cφsψ) − ci (cφcψ + sθsφsψ)
(B.40)
que são os cosenos diretores de Yi em relação aos eixos do satélite (x, y, z)
174
c31 = −sisΩ (−cψsφ + cφsθsψ) + (cφcψsθ + sφsψ) (−cΩsu − cicusΩ) + cθcφ (cisusΩ − cucΩ)
c22 = sicΩ (−cψsφ + cφsθsψ) + (cφcψsθ + sφsψ) (−sΩsu + cicucΩ) + cθcφ (−cisucΩ − cusΩ)
c23 = −cθcφsisu − ci (cφsθsψ − cψsφ) + cusi (cφcψsθ + sφsψ)
(B.41)
que são os cosenos diretores de Zi em relação aos eixos do satélite (x, y, z)
A matriz RBF
ECI é usada para expressar o vetor campo magnético, fornecido em
coordenadas inerciais, no sistema do satélite.
As Tabelas (B.2),(B.3) e (B.4) apresentam o programa em MATLAB para o cálculo
BF
da matriz RECI
.
B.3
Cálculo da Latitude, Longitude e Altura
Para o cálculo da latitude, longitude e altura, primeiro expressamos a posição do
satélite no referencial cartesiano geocêntrico.
Referencial cartesiano terrestre ou geocêntrico: (FG) (Xg , Yg , Zg ) é um sistema de coordenadas com origem no centro de massa da Terra. O eixo Xo aponta
na direção do meridiano de Greenwich, o eixo Zo aponta na direção normal ao plano
do equador. O eixo Yo é obtido pela regra da mão direita. A Figura (B.4) ilustra o
referencial inercial e o referencial cartesiano geocêntrico (Kuga e Kondapalli, 1995).
Na Figura (B.4) θg é o tempo sideral de Greeenwich. O cálculo do tempo sideral
de Greenwich é mostrado em Escobal (1965). A matriz de transformação do sistema
inercial para o sistema cartesiano geocêntrico é dada por (Pilchowski e Ferreira,
1981)


cos θ sin θ 0


G
RFECI
(θg ) = − sin θ cos θ 0
0
0
1
175
(B.42)
BF .
Tabela B.2- Programa MATLAB para o cálculo de c1i , (i = 1, 2, 3) da matriz RECI
function [R1] = ECI2BF X(par)
% by Gilberto Arantes 2004 INPE
%
% [R1] = ECI2BF X(par) Cosenos diretores X.x X.y X.z
%
% Calculo dos cosenos diretores
%
% Input
% par = [orb atitude]
%
% atitude [phi teta psi] = angulos de atitude (rad)
% orb = [u i asc nodo] = elementos orbitais/Keplerianos (rad, graus, graus]
%
% Output
% R1 = [r11,r21,r31]
% r11 = cos(X.x)
% r21 = cos(X.y)
% r31 = cos(X.z)
orb(1) = par(1); %u
orb(2) = par(2); %i
orb(3) = par(3); %asc nodo
%
f = pi/180; %graus->radianos
at(1) = par(4); %rad phi
at(2) = par(5); %rad teta
at(3) = par(6); %rad psi
%
% Calculo dos cosenos diretores
%
% cos(X.x)
r11 = -cos(at(2))*sin(orb(2)*f)*sin(at(3))*sin(orb(3)*f)+...
cos(at(2))*cos(at(3))*(-cos(orb(3)*f)*sin(orb(1))-cos(orb(2)*f)*cos(orb(1))*sin(orb(3)*f))-...
sin(at(2))*(-cos(orb(1))*cos(orb(3)*f)+cos(orb(2)*f)*sin(orb(1))*sin(orb(3)*f));
% cos(X.y)
r21 = cos(at(2))*sin(orb(2)*f)*sin(at(3))*cos(orb(3)*f)+...
cos(at(2))*cos(at(3))*(-sin(orb(3)*f)*sin(orb(1))+cos(orb(2)*f)*cos(orb(1))*cos(orb(3)*f)-...
sin(at(2))*(-cos(orb(1))*sin(orb(3)*f)-cos(orb(2)*f)*sin(orb(1))*cos(orb(3)*f));
% cos(X.z)
r31 = cos(at(2))*sin(orb(2)*f)*cos(at(3))*cos(orb(1))+...
sin(orb(2)*f)*sin(orb(1))*sin(at(2))-cos(orb(2)*f)*cos(at(2))*sin(at(3));
%
R1 = [r11 r21 r31];
%
% end
Obtido as coordenadas do satélite no referencial cartesiano geocêntrico podemos
obter a latitude φ a longitude λ e a altura h do satélite, como segue.
176
BF .
Tabela B.3- Programa MATLAB para o cálculo de c2i , (i = 1, 2, 3) da matriz RECI
function [R2] = ECI2BF Y(par)
% by Gilberto Arantes 2004 INPE
%
% [R2] = ECI2BF Y(par) Cosenos diretores Y.x Y.y Y.z
%
% Calculo dos cosenos diretores
%
% Input
% par = [orb atitude]
% atitude [phi teta psi] = angulos de atitude (rad)
% obr = [u i asc nodo] = elementos orbitais/Keplerianos (rad, graus, graus]
%
% Output
% R1 = [r12,r22,r32]
% r12 = cos(Y.x)
% r22 = cos(Y.y)
% r32 = cos(Y.z)
%
f = pi/180; %graus->radianos
u = par(1); % rad
i = par(2); % graus
nodo = par(3); % graus
%
phi = par(4); %rad phi
teta = par(5); %rad teta
psi = par(6); %rad psi
%
% Calculo dos cosenos diretores
%
% cos(Y.x)
r12 = -sin(i*f)*sin(nodo*f)*(cos(phi)*cos(psi)+sin(teta)*sin(phi)*sin(psi))...
+(cos(psi)*sin(teta)*sin(phi)-cos(phi)*sin(psi))*(-cos(nodo*f)*sin(u)-cos(i*f)*cos(u)*sin(nodo*f))+...
cos(teta)*sin(phi)*(-cos(u)*cos(nodo*f)+cos(i*f)*sin(u)*sin(nodo*f));
% cos(Y.y)
r22 = cos(nodo*f)*sin(i*f)*(cos(phi)*cos(psi)+sin(teta)*sin(phi)*sin(psi))+...
cos(teta)*sin(phi)*(-cos(u)*sin(nodo*f)-cos(i*f)*sin(u)*cos(nodo*f))+...
(cos(psi)*sin(teta)*sin(phi)-cos(phi)*sin(psi))*(-sin(nodo*f)*sin(i*f)+...
cos(i*f)*cos(u)*cos(nodo*f));
% cos(Y.z)
r32 = -cos(teta)*sin(i*f)*sin(u)*sin(phi)+...
cos(u)*sin(i*f)*(cos(psi)*sin(teta)*sin(phi)-cos(phi)*sin(psi))-...
cos(i*f)*(cos(phi)*cos(psi)+sin(teta)*sin(phi)*sin(psi));
%
R2 = [r12 r22 r32];
%
% end
Observando a Figura (B.5), temos que
177
BF .
Tabela B.4- Programa MATLAB para o cálculo de c3i , (i = 1, 2, 3) da matriz RECI
function [R3] = ECI2BF Z(par)
% by Gilberto Arantes 2004 INPE
%
% [R3] = ECI2BF Z(par) Cosenos diretores Z.x Z.y Z.z
%
% Calculo dos cosenos diretores
%
% Input
% par = [orb atitude]
% atitude [phi teta psi] = angulos de atitude (rad)
% obr = [u i asc nodo] = elementos orbitais/Keplerianos (rad, graus, graus]
%
% Output
% R3 = [r13,r23,r33]
% r13 = cos(Z.x)
% r23 = cos(Z.y)
% r33 = cos(Z.z)
%
f = pi/180; %graus->radianos
u = par(1); % rad
i = par(2); %rad
nodo = par(3); %rad
%
phi = par(4); %rad phi
teta = par(5); %rad teta
psi = par(6); %rad psi
%
% Calculo dos cosenos diretores
%
% cos(Z.x)
r13 = -sin(i*f)*sin(nodo*f)*(-sin(phi)*cos(psi)+sin(teta)*cos(phi)*sin(psi))+...
(cos(psi)*sin(teta)*cos(phi)+sin(phi)*sin(psi))*(-cos(nodo*f)*sin(u)-cos(i*f)*cos(u)*sin(nodo*f))+...
cos(teta)*cos(phi)*(-cos(u)*cos(nodo*f)+cos(i*f)*sin(u)*sin(nodo*f));
% cos(Z.y)
r23 = cos(nodo*f)*sin(i*f)*(-sin(phi)*cos(psi)+sin(teta)*cos(phi)*sin(psi))+...
cos(teta)*cos(phi)*(-cos(u)*sin(nodo*f)-cos(i*f)*sin(u)*cos(nodo*f))+...
(cos(psi)*sin(teta)*cos(phi)+sin(phi)*sin(psi))*(-sin(nodo*f)*sin(u)+cos(i*f)*cos(u)*cos(nodo*f));
% cos(Z.z)
r33 = -cos(teta)*sin(i*f)*sin(u)*cos(phi)+...
cos(u)*sin(i*f)*(cos(psi)*sin(teta)*cos(phi)+sin(phi)*sin(psi))-...
cos(i*f)*(-sin(phi)*cos(psi)+sin(teta)*cos(phi)*sin(psi));
%
R3 = [r13 r23 r33];
%
% end
 


Xg
cos φ cos λ
 


 Yg  = (R⊗ + h)  cos φ sin λ 
Zg
sin φ
178
(B.43)
Figura B.4- Referencial inercial e cartesiano geocêntrico.
Figura B.5- Longitude (λ), latitude (φ) e altura (h) do satélite.
Logo, a longitude será dada por:
λ = arctan
Yg
Xg
onde 0 ≤ λ ≤ 2π
(B.44)
A latitude será dada por:
φ = arcsin
Zg
R⊗ + h
onde
179
−
π
π
≤φ≤
2
2
(B.45)
e a altura dada por:
h=
q
Xg2 + Yg2 + Zg2 − R⊗
(B.46)
A Tabela (B.5) apresenta o programa em MATLAB para o cálculo da longitude,
latitude do ponto sub-satélite e altura do satélite.
180
Tabela B.5- Programa MATLAB para o cálculo da longitude, latitude do ponto sub-satélite e altura do satélite.
function [LLA] = ECI2LLA(par)
% by Gilberto Arantes 2004 INPE
%
% [LLA] = ECI2LAA(par) para uma orbita circular
%
% Calculo da latitude, longitude do ponto sub-satelite e altitude do satelite
%
% Input
% par = [ teta g X ]
% teta g = Tempo sideral em graus
% X = [X Y Z] coordenadas no referencial inercial em km
%
% Output
% LLA = [latitude longitude altitude]
% latitude em rad
% Longitude em rad
% altitude em km
%
teta g = par(1);
X(1) = par(2)*1000; % transformando para metros
X(2) = par(3)*1000;
X(3) = par(4)*1000;
%
Raio t = 6378000; % R = Raio da Terra em m
%
% Calculo da matriz de tranformacao do sistema inercial para o terresrte (ref.: Silva, 1981)
%
fator = pi/180;
a = teta g*fator;
%
R = [ cos(a) sin(a) 0
-sin(a) cos(a) 0
0 0 1];
% Posicao no referencial terrestre
x g = R*transpose(X); % x g = [x g y g z g]’
% Calculo da altura
h = ((sqrt( x g(1)^
2 + x g(2)^
2 + x g(3)^
2 ) - Raio t))/1000;
% Calculo da latitude
lat = asin(x g(3)/(Raio t+h)); % Sul (-) Norte (+)
% Calculo da longitude
long = atan2(x g(2),x g(1)); % Leste (+) Oeste (-)
% conversao da saida (0 <= long <= 360)
if sign(long) < 0;
long = 2*pi+long;
else sign(long) >= 0;
long = long;
end
%
LLA = [lat long h];
%
% end
181
182
APÊNDICE C
IMPLEMENTAÇÃO EM SIMULINK
Os modelos apresentados no Capı́tulo (5) e as leis de controle apresentadas no
Capı́tulo (6) foram implementedas em SIMULINK. O sistema foi dividido nos blocos
de: 1) dinâmica que contêm o modo completo do satélite equipado com bobinas
e rodas; 2) ambiente que contem o modelo do gradiente de gravidade e do campo
magnético; 3) órbita com o cálculo da órbita (modelo kepleriano) do satélite e 3)
com os controladores avaliados para os três modos de operação. A implentação é
feita para os três modos de operação:
• Detumble;
• Estabilização;
• aquisição.
A Figura (C.1) mostra a implemtação em SIMLINK do modo de detumble. A
Figura (C.2) mostra o controlador usado para o detumbling.
Figura C.1- Implemetação em SIMULINK do Modo de detumble.
183
Figura C.2- Lei de controle Bdot.
A Figura (C.3) mostra o modo de estabilização, para o procedimento de controle
caracterizado por bobinas, os diferentes controladores, LQR e baseados em energia,
que são implentados usando-se multiport switches, Figura (C.4), o que torna fácil
selecionar a lei de controle para o modo de estabilização.
Figura C.3- Implementação em SIMULINK do modo de estabilização: satélite
equipado com bobinas.
As Figuras (C.5) e (C.6) mostram o modo de estabilização, com comtrole via rodas
de reação. A primeira utilizando o método LQR e a segunda utilizando o método
LQG. São mostrados o modelo linear e o completo (não linear) da dinâmica, para
efeito de comparação dos resultados relativos à aplicação das leis de controle para
ambos os modelos.
184
Figura C.4- Controle LQR e controladores baseados em energia.
Figura C.5- Modo de estabilização: satélite equipado com rodas, usando o método
LQR.
185
Figura C.6- Modo de estabilização: satélite equipado com rodas, usando o método
LQG.
A Figura (C.7) mostra o sistema de controle para o modo de aquisição utilizandose um controlador PD. A Figura (C.8) mostra a lei de controle do método LQR
tracking implementada para a fase de aquisição.
Figura C.7- Modo de aquisição usando o controlador PD.
186
Figura C.8- Lei de controle LQR tracking.
Figura C.9- Lei de controle PD.
As Figuras (C.10) e (C.11) mostram a implementação em SIMULINK dos modelos:
satélite equipado com bobinas e satélite equipado com rodas, respectivamente.
Figura C.10- Modelo dinâmico do satélite equipado com bobinas.
A Figura (C.12) mostra o modelo das rodas de reação. Note que as saturações são
incluı́das. A Figura (C.13) mostra o modelo das bobinas magnéticas.
187
Figura C.11- Modelo dinâmico do satélite equipado com rodas.
Figura C.12- Modelo das rodas.
A Figura (C.14) mostra o modelo do gradiente de gravidade e a Figura (C.15) mostra
o modelo do campo magnético.
As Figuras (C.16), (C.17) e C.18) mostram a implemetação da propagação da órbita,
do cálculo da latitude e longitude do ponto sub-satélite e das transformações dos
referenciais envolvidos, respectivamente. Note que a implemtentação para órbita e as
transformações utilizam os arquivos mostrados em MATLAB mostrados na Seção
(B).
188
Figura C.13- Modelo das bobinas.
Figura C.14- Modelo do torque de gradiente de gravidade.
189
Figura C.15- Modelo do campo magnético.
Figura C.16- Propagação da órbita.
190
Figura C.17- Cálculo da latitude, longitude e altura.
Figura C.18- Matriz de transformação entre os sistemas de referência ECI e BF.
191
192
APÊNDICE D
PROGRAMAS EM MATLAB
A seguir são apresentados os arquivos de inicialização com o projeto LQR para os
modelos do satélite equipado com rodas e bobinas e o projeto LQG.
D.1
Projeto LQR
As Tabelas (D.1) e (D.3) mostram os programas em MATLAB para o projeto LQR
dos dois modelos analisados.
D.2
Projeto LQG
A Tabela (D.4) mostra o programa em MATLAB para o projeto LQG, usado no
modelo do satélite equipado com rodas.
193
Tabela D.1- Programa MATLAB para o projeto LQR do modelo satélite equipado com rodas.
% Por Gilberto Arantes @INPE 2003
%
% Gera as matrizes para o projeto do sistema de controle LQR e LQR rastreio
% Calculo da equacao matricial de Ricatti e implementacao do LQR
% Matrizes A, B, C e D do sistema: satélite equipado com rodas
% d/dt(x) = Ax + Bu
% y = Cx +du
% Inicializacao clc
clear all
close all
%
deg2rad = pi/180;
rad2deg = 180/pi;
mi=3.986e14; % moviemnto medio em (m)
r=(6378+750)1000; % raio da orbita circular em (m)
n=sqrt(mir^
3)); % movimento orbital medio em (rad/s)
H o = n*0.01; % Momento angular nominal kgm2/s
% Momentos e produtos de inercia (14/10/2003 de Sebastiao)
% Ixx = 13.31 kg.m2 Ixy = 0.37 kg.m2
% Iyy = 14.22 kg.m2 Ixz = 0.39 kg.m2
% Izz = 11.2 kg.m2 Iyz = -0.22 kg.m2
I = [13.31 0.37 0.39
0.37 14.22 -0.22
0.39 -0.22 11.2];
% Momentos principais de inercia
[v, e] = eig(I);
I1 = e(1,1); I2 = e(2,2); I3 = e(3,3);
% inercia da roda
Iw = 12/(45000*(pi/180)) ;
hs = 12; %saturacao
ws = 7500*6 ;%7500rpm*6graus/s
a2rpm = 1/6;
%--------------------------------------------------------------------% Dinamica linearizada: satelite + rodas
%---------------------------------------------------------------------% Matriz A
A =[0 0 0 1 0 0;0 0 0 0 1 0;0 0 0 0 0 1;...
(-4*n^
2*(I2-I3)-n*H o)/I1 0 0 0 0 (n*(-I2+I3+I1)-H o)/I1;...
0 3*n^
2*(-I1+I3)/I2 0 0 0 0;...
0 0 (n^
2*(-I2+I1)-n*H o)/I3 (-n*(I1-I2+I3)+H o)/I3 0 0];
% Matriz B
B=[0 0 0
0 0 0
0 0 0
-1/I1 0 0
0 -1/I2 0
0 0 -1/I3];
%Matriz C e D
C=eye(6);
D=zeros(6,3);
194
Tabela D.2- Programa MATLAB para o projeto LQR do modelo satélite equipado com rodas (cont.).
% Calculo da matriz de controlabilidade
co = ctrb(A,B);
posto co = rank(co);
posto A = rank(A);
%------------------------------------------------% Lei LQ Regulador
%-----------------------------------------------% Matrizes R e Q (ganho) ref. Overby, 2004
rc = 1/((5e-3)^
(2)); % u = 5mNm
qc1 = 1/((10*(pi/180))^
(2)); % x = 10 graus
qc2 = 1/((1*(pi/180))^
(2)); % xdot = 1 grau/s
Rc = rc*diag([1,1,1]); %peso do controle
Qc = diag([qc1,qc1,qc1,qc2,qc2,qc2]); %peso do estado
% Matriz Kc: ganho do controlador
[Kc, S, E] = LQR(A, B, Qc, Rc);
%----------------------------------------------% Tracking r(t)
%----------------------------------------------% u = Kc*x + v
% v = -inv(R)B’s
% ds/dt = -[A’-K*B*inv(R)*B’]s + Q*r(t)
% r(t) vector reference
At = - [A’-S*B*inv(Rc)*B’];
Bt = Qc;
Ct = C;
Dt = zeros(size(A));
% simulink
run(’SimLQEquarsRodasSat’)
% end
195
Tabela D.3- Programa MATLAB para o projeto LQR do modelo satélite equipado com bobinas.
% Por Gilberto Arantes INPE 2004
%
% sintaxe [K] = LQcontrolNcube(B)
% Output
% K = control matrix gain
% Input
% B = magnetic field
% -----------------------------% LQ Initialization
%------------------------------function [K] = LQcontrolNcube(B)
% Initial Values
Ix = 18.941;
Iy = 21.726;
Iz = 16.983;
n = 1.083e-3;
kx = (Iy - Iz)/Ix;
ky = (Ix - Iz)/Iy;
kz = (Iy - Ix)/Iz;
% The Goeomagnetic field
Bx o = B(1);
By o = B(2);
Bz o = B(3);
A = [ 0 0 0 1 0 0;
0 0 0 0 1 0;
0 0 0 0 0 1;
-4*kx*n^
2 0 0 0 0 (1-kx)*n;
0 -3*ky*n^
2 0 0 0 0;
0 0 -kz*n^
2 -(1-kz)*n^
2 0 0];
% The input matrix for the linearied system
B = [zeros(3,3);
0 Bz o/(2*Ix) -By o/(2*Ix);
-Bz o/(2*Iy) 0 Bx o/(2*Iy);
By o/(2*Iz) -Bx o/(2*Iz) 0];
% LQ weighting matrices
Q = diag([1 1 1 0 0 0])*inv(10*pi/180)^
2;
P = diag([1 1 1])*inv(.1)^
2;
K = -lqr(A,B,Q,P);
% Coil Parameters
Nx = 100; Ny = 100; Nz = 100; % nunber of coil windings
Ax = 0.075; Ay = 0.075; Az = 0.075; % coil area [m^
2]
Rx = 20; Ry = 20; Rz = 20; % Coil resistance [Ohm]
M = Nx*Ax;
%end
196
Tabela D.4- Programa MATLAB para o projeto LQG do modelo satélite equipado com rodas.
% Por Gilberto Arantes @INPE 2003
%
% Gera as matrizes para o projeto do sistema de controle LQG
% Calculo da equacao matricial de Ricatti e implementacao do LQG
% Matrizes A, B, C e D do sistema d/dt(x) = Ax + Bu
% y = Cx +du
% Inicializacao
clc
clear all
close all
deg2rad = pi/180;
rad2deg = 180/pi;
mi=3.986e14; % moviemnto medio em (m)
r=(6378+750)*1000; % raio da orbita circular em (m)
n=sqrt(mi/(r^
3)); % movimento orbital medio em (rad/s)
H o = n*0.01; % Momento angular nominal kgm2/s
% Momentos e produtos de inercia (14/10/2003 de Sebastiao)
% Ixx = 13.31 kg.m2 Ixy = 0.37 kg.m2
% Iyy = 14.22 kg.m2 Ixz = 0.39 kg.m2
% Izz = 11.2 kg.m2 Iyz = -0.22 kg.m2
I = [13.31 0.37 0.39
0.37 14.22 -0.22
0.39 -0.22 11.2];
% Momentos principais de inercia
[v,e] = eig(I);
I1 = e(1,1); I2 = e(2,2); I3 = e(3,3);
% inercia da roda
Iw = 12/(45000*(pi/180)) ;
hs = 12; %saturacao
ws = 7500*6 ;%7500rpm*6graus/s
a2rpm = 1/6;
%-------------------------------------------------------% Dinamica linearizada: satelite + rodas
%-------------------------------------------------------% matriz A
A = [0 0 0 1 0 0;0 0 0 0 1 0;0 0 0 0 0 1;...
(-4*n^
2*(I2-I3)-n*H o)/I1 0 0 0 0 (n*(-I2+I3+I1)-H o)/I1;...
0 3*n^
2*(-I1+I3)/I2 0 0 0 0;...
0 0 (n^
2*(-I2+I1)-n*H o)/I3 (-n*(I1-I2+I3)+H o)/I3 0 0];
% Matriz B
B=[0 0 0
0 0 0
0 0 0
-1/I1 0 0
0 -1/I2 0
0 0 -1/I3];
%Matriz C e D
C=eye(6);
D=zeros(6,3);
197
Tabela D.5- Programa MATLAB para o projeto LQG do modelo satélite equipado com rodas (cont.).
% Calculo da matriz de controlabilidade
co = ctrb(A,B);
posto co = rank(co);
posto A = rank(A);
%--------------------------------------------------% Lei LQ Regulador
%--------------------------------------------------Matria R e Q (ganho)
rc = 1/((5e-3)^
(2)); % u = 5mNm
qc1 = 1/((10*(pi/180))^
(2)); % x = 10 graus
qc2 = 1/((1*(pi/180))^
(2)); % xdot = 1 grau/s
Rc = rc*diag([1,1,1]); %peso do controle
Qc = diag([qc1,qc1,qc1,qc2,qc2,qc2]); %peso do estado
% Matriz Kc: ganho do controlador
[Kc,S,E] = LQR(A,B,Qc,Rc);
%--------------------------------------------------% Tracking r(t)
%--------------------------------------------------% u = Kc*x + v
% v = -inv(R)B’s
% ds/dt = -[A’-K*B*inv(R)*B’]s + Q*r(t)
% r(t) vector reference
At = - [A’-S*B*inv(Rc)*B’];
Bt = Qc;
Ct = C;
Dt =zeros(size(A));
%-------------------------------------------------%--------------------------------------------------% Implementation Kalman filter
% dynamic d/dt(x) = Ax + Bu + Gw
% y = Cx + du + v
% state estimation
% d/dt(xe) = Axe + Bu + Kf(delta - Cxe - Du) delta = yn-ym
% d/dt(xe) = Axe + B(-kcxe) + Kf(delta - Cxe - Du)
% covariance noise sensor: Evv’ = Rk
% covariance noise dynamic Eww’ = Qk
G = eye(6,6);
sigma angle = (0.5*deg2rad)^
2; % desvio padrao dos sensor de orientacao
sigma rate = (0.05*deg2rad)^
2; % desvio padrao dos giros
Qk = eye(6)*(15*deg2rad)^
2; % angle +/-(0.5 graus )
Rk = [eye(3)*sigma angle zeros(3); zeros(3) eye(3)*sigma rate];
[Kf,P,E2] = lqe(A,G,C,Qk,Rk);
% Structure of compensator
% Ac = A - B*Kc - Kf*C;
% Bc = Kf;
% Simulink
run(’SimLQGEquarsRodas’)
%end
198
APÊNDICE E
TOOLBOX ATITUDE
Durante o desenvolvimento do trabalho foi desenvolvido um toolbox para atitude,
no ambiente MATLAB e SIMULINK, constituindo uma importante contribuição
desse trabalho. Esse apêndice tem por objetivo divulgar o uso do toolbox, apresentando alguns exemplos. A versão final do toolbox será disponibilizado para uso do
INPE e instituições interessadas em pesquisas envolvendo dinâmica e controle de
atitude. A documetação do toolbox, estará disponı́vel na bibiloteca do INPE.
O toolbox é divido nos blocos:
• Equações do movimento: Cinemática e dinâmica
• Órbita
• Ambiente
• Controle
• Visualização
O bloco Equações do movimento contêm as equações da cinemática em diferentes
representações: ângulos de Euler e quaternions. O modelo dinâmico da atitude (corpo
rı́gido) e o modelo dinâmico do gyrostat (extensão das equações do corpo rı́gido
envolvendo movimento interno de rodas).
O bloco Órbita contém modelo de órbita Kepleriana, cálculo da latitude e longitude
sub-satélite e altura do veı́culo e transformações entre os referenciais envolvidos no
estudo de atitude e órbita.
O bloco Ambiente contém o modelo do torque de gradiente de gravidade, modelo
do campo magnético IGRF, e modelo do torque devido ao momento dipolo residual
do satélite.
O bloco Controle contém os controladores desenvolvidos nesse trabalho: LQR, LQR
tracking, LQG, PID e controladores baseados em energia (Energy based control ).
199
O bloco Visualização contém a visualização/animação da atitude do veı́culo.
E.1
E.1.1
Exemplos
Equações do Movimento
A Figura (E.1) mostra um dos blocos que fazem parte do bloco equações do movimento, o bloco mostrado integra as equações do movimento de um corpo rı́gido
contendo rodas (gyrostat), modelo completo, não linear. No ambiente SIMULINK
existem vários integradores disponı́veis de passo fixo e/ou variável. A janela mostra
os parâmetros de entrada.
(a) Modelo do Gyrostat
(b) Parâmetros de entrada
Figura E.1- Bloco Gyrostat.
E.1.2
Ambiente
O bloco mostrado na Figura (E.2) calcula o modelo completo, não linear, do torque
de gradiente de gravidade. A janela mostra os parâmetros de entrada.
200
(a) Modelo do gradiente
de gravidade
(b) Parâmetros de entrada
Figura E.2- Bloco gradiente de gravidade.
201
Livros Grátis
( http://www.livrosgratis.com.br )
Milhares de Livros para Download:
Baixar livros de Administração
Baixar livros de Agronomia
Baixar livros de Arquitetura
Baixar livros de Artes
Baixar livros de Astronomia
Baixar livros de Biologia Geral
Baixar livros de Ciência da Computação
Baixar livros de Ciência da Informação
Baixar livros de Ciência Política
Baixar livros de Ciências da Saúde
Baixar livros de Comunicação
Baixar livros do Conselho Nacional de Educação - CNE
Baixar livros de Defesa civil
Baixar livros de Direito
Baixar livros de Direitos humanos
Baixar livros de Economia
Baixar livros de Economia Doméstica
Baixar livros de Educação
Baixar livros de Educação - Trânsito
Baixar livros de Educação Física
Baixar livros de Engenharia Aeroespacial
Baixar livros de Farmácia
Baixar livros de Filosofia
Baixar livros de Física
Baixar livros de Geociências
Baixar livros de Geografia
Baixar livros de História
Baixar livros de Línguas
Baixar livros de Literatura
Baixar livros de Literatura de Cordel
Baixar livros de Literatura Infantil
Baixar livros de Matemática
Baixar livros de Medicina
Baixar livros de Medicina Veterinária
Baixar livros de Meio Ambiente
Baixar livros de Meteorologia
Baixar Monografias e TCC
Baixar livros Multidisciplinar
Baixar livros de Música
Baixar livros de Psicologia
Baixar livros de Química
Baixar livros de Saúde Coletiva
Baixar livros de Serviço Social
Baixar livros de Sociologia
Baixar livros de Teologia
Baixar livros de Trabalho
Baixar livros de Turismo
Download

estudo comparativo de técnicas de controle de atitude em três eixos