Otimização de Empuxo de um Foguete Movendo-se na Atmosfera

sidinpebrMTCm1380200508241850TDI OTIMIZAÇÃO DE EMPUXO DE UM FOGUETE MOVENDOSE NA ATMOSFERA Maxime Koffi Dissertação de Mestrado do Curso de PósGraduação em Engenharia e Tecnologia EspaciaisCombustão e Propulsão orientada pelo Dr Jerzy Tadeusz Sielawa aprovada em 02 de junho de 2005 URL do documento original httpurlibnet6qtX3pFwXQZGivnJSYHfKwM INPE São José dos Campos 2005 PUBLICADO POR Instituto Nacional de Pesquisas Espaciais INPE Gabinete do Diretor GB Serviço de Informação e Documentação SID Caixa Postal 515 CEP 12245970 São José dos Campos SP Brasil Tel012 320869236921 Fax 012 32086919 Email pubtcsidinpebr COMISSÃO DO CONSELHO DE EDITORAÇÃO E PRESERVAÇÃO DA PRODUÇÃO INTELECTUAL DO INPE DEDIR544 Presidente Marciana Leite Ribeiro Serviço de Informação e Documentação SID Membros Dr Gerald Jean Francis Banon Coordenação Observação da Terra OBT Dr Amauri Silva Montes Coordenação Engenharia e Tecnologia Espaciais ETE Dr André de Castro Milone Coordenação Ciências Espaciais e Atmosféricas CEA Dr Joaquim José Barroso de Castro Centro de Tecnologias Espaciais CTE Dr Manoel Alonso Gan Centro de Previsão de Tempo e Estudos Climáticos CPT Dra Maria do Carmo de Andrade Nono Conselho de PósGraduação Dr Plínio Carlos Alvalá Centro de Ciência do Sistema Terrestre CST BIBLIOTECA DIGITAL Dr Gerald Jean Francis Banon Coordenação de Observação da Terra OBT Clayton Martins Pereira Serviço de Informação e Documentação SID REVISÃO E NORMALIZAÇÃO DOCUMENTÁRIA Simone Angélica Del Ducca Barbedo Serviço de Informação e Documentação SID Yolanda Ribeiro da Silva Souza Serviço de Informação e Documentação SID EDITORAÇÃO ELETRÔNICA Marcelo de Castro Pazos Serviço de Informação e Documentação SID André Luis Dias Fernandes Serviço de Informação e Documentação SID sidinpebrMTCm1380200508241850TDI OTIMIZAÇÃO DE EMPUXO DE UM FOGUETE MOVENDOSE NA ATMOSFERA Maxime Koffi Dissertação de Mestrado do Curso de PósGraduação em Engenharia e Tecnologia EspaciaisCombustão e Propulsão orientada pelo Dr Jerzy Tadeusz Sielawa aprovada em 02 de junho de 2005 URL do documento original httpurlibnet6qtX3pFwXQZGivnJSYHfKwM INPE São José dos Campos 2005 Dados Internacionais de Catalogação na Publicação CIP Koffi Maxime K821o Otimização de empuxo de um foguete movendose na atmosfera Maxime Koffi São José dos Campos INPE 2005 xxii 67 p sidinpebrMTCm1380200508241850TDI Dissertação Mestrado em Engenharia e Tecnologia EspaciaisCombustão e Propulsão Instituto Nacional de Pesquisas Espaciais São José dos Campos 2005 Orientador Dr Jerzy Tadeusz Sielawa 1 Trajetória do foguete 2 Otimização 3 Míssil 4 Empuxo 5 Cálculo variacional 6 Combustível ITítulo CDU 6297 Esta obra foi licenciada sob uma Licença Creative Commons AtribuiçãoNãoComercial 30 Não Adaptada This work is licensed under a Creative Commons AttributionNonCommercial 30 Unported License ii Aprovado a pela Banca Examinadora em cumprimento ao requisito exigido para obtenção do Título de Mestre em Engenharia e Tecnologia EspaciaisCombustão e Propulsão iv v Com muito amor a meus pais Koffi Kan Goerges e Oka Amenan Odette vi vii AGRADECIMENTOS Ao meu orientador Prof Dr Jerzy Tadeusz Sielawa a quem tenho uma grande admiração e respeito agradeço pela predisposição compreensão e amizade que teve comigo Aos meus professores do Instituto Nacional de Pesquisas Espaciais INPE em especial ao Dr Demétrio Bastos Netto e ao Dr Fernando de Souza Costa pelo conhecimento compartilhado durante toda a realização do curso Aos membros da banca examinadora pelas sugestões recebidas e ao Governo Brasileiro através da Coordenação de Aperfeiçoamento de Pessoal de Nível Superior CAPES pela bolsa de mestrado concedida Aos meus colegas de trabalho e aos funcionários agradeço pela amizade e confiança Não poderia ainda esquecer da minha família Maxime Koffi Junior e da minha esposa Lúcia Helena Agostinho que sempre estiveram presentes nas horas difíceis viii ix RESUMO O empuxo é tratado como uma variável de controle na otimização da altitude de um foguete movendose na atmosfera São determinadas as características de empuxo para obter a altitude máxima de um foguete lançado verticalmente na atmosfera com uma dada quantidade de combustível Neste trabalho é mostrada passo a passo a seqüência de controle ótima a ser aplicada desde o instante inicial do lançamento do foguete até o fim da queima considerandose a superfície da Terra como um sistema de referência inercial Em seguida é feita também uma análise do alcance máximo do foguete quando a superfície da Terra é considerada curvada ou seja admitindose os efeitos da rotação da Terra Ao final apresentase um exemplo de resolução numérica para um foguete hipotético O presente estudo utiliza o cálculo variacional para a análise da trajetória de um foguete usando o empuxo como variável de controle no entanto o método pode ser aplicado a outros tipos de problemas considerando outras variáveis de controle x xi OPTIMIZATION OF THE PUSH OF A ROCKET MOVING IN THE ATMOSPHERE ABSTRACT In this work the thrust is considered as a control variable in order to maximize the altitude attained by a rocket Initially the Earth Surface is taken as an inertial reference frame The rocket is launched vertically into the atmosphere with a given quantity of fuel It is shown how optimal control can be applied from the initial launch point up to the burnout An example with a numerical solution is presented for a hypothetical rocket The rocket trajectory launched vertically is also analyzed considering the effects of the Earth rotation The present study makes possible the analysis of a rocket trajectory from a mathematical point of view using Variational Calculus with thrust as a control variable obviously the method can be applied using other control variables xii xiii LISTA DE SÍMBOLOS Latinos a h velocidade do som P20 A constante de integração eq4244 e A área de saída da tubeira eq2113 B constante de integração eq4243 C constante de integração eq4234 D C coeficiente de arrasto eq2122 L C coeficiente de sustentação P20 D força de arrasto eq212 Det determinante eq3152 ETot energia total do foguete eq5124 EmTot energia mecânica total eq4134 e excentricidade eq4149 neˆ versor normal eq4312 reˆ versor radial das bases do mesmo sistema eq4111 teˆ versor tangencial eq4313 F força de tração eq211 F r força do campo eq4111 gi t função genérica eq3114 0 g aceleração gravitacional na superfície da Terra eq2133 g h aceleração gravitacional no ponto onde se encontra o corpo eq2133 G constante gravitacional eq2132 h t altitude do foguete no instante t eq5119 hq tq altitude do foguete no final da queima eq5122 xiv hmax altitude máxima do foguete eq5123 H função Hamiltoniana eq31211 h distancia entre o centro de massa do foguete e a superfície da Terra eq2133 h quantidade do movimento angular eq4132 J funcional eq311 k variável de controle eq2214c kmax variável de controle ótima eq5120 S k variável de controle singular eq3149 L força de sustentação P20 m massa do foguete eq2131 0 m massa inicial do foguete massa da estrutura massa do propelente eq313c m t massa instantânea do foguete eq5111 c m massa do combustível eq5120 f m massa final eq5111 m taxa mássica do foguete eq2214c T M massa da Terraeq2131 M numero de Mach P20 M parâmetro eq4256 N parâmetro eq4257 P h força provocada no sentido do movimento devido à diferença entre a pressão dos gases na saída da tubeira do motor e a pressão atmosférica no ponto onde se encontra o foguete eq2113 e P pressão dos gases na saída da tubeira do motor eq2113 P1 h pressão atmosférica no ponto onde se encontra o foguete eq2113 p semi latus rectumeq4148 xv Q parâmetro eq4258 T R raio da Terra eq2133 r coordenada radial do sistema cilíndrico eq4147 or posição inicial do foguete em relação ao centro da Terra no burn out eq4221 0r componente radial da velocidade eq4222 S t função Comutadora eq3134 S fator de proporcionalidade constante eq4112 T empuxoeq2112 qt tempo de queimaeq5120 ft tempo final eq3129 0t tempo inicial eq3129 cT energia cinética eq4131 t velocidade instantânea do foguete eq518 q tq velocidade no final da queima eq5121 E velocidade de ejeção dos gaseseq2214 0 velocidade resultante inercial no lançamento eq4313 V0 velocidade do foguete em relação à superfície da Terra eq4311 tV velocidade linear na direção tangencial eq4323 V energia potencial eq4114 u variável adimensional eq4144 W força gravitacional terrestre eq2132 0 W força gravitacional na superfície da Terra que age num corpo a uma altitude h eq2131 xi t componente vetorial das variáveis de estado eq3111 xvi Gregos parâmetro eq4145 ângulo de lançamentoeq4221 parâmetro eq4145 0 ângulo inicial da trajetória em relação á linha horizontal eq4221 max alcance da trajetória eq4251 módulo da velocidade angular eq4324 densidade atmosférica eq2122 0 variável adimensional eq4326 i componentes vetoriais linha dos multiplicadores de Lagrange 3 2 1 eq3124 i T matriz transposta aos multiplicadores de lagrange eq3124 componente polar da velocidade eq4223 0 latitude do sítio do lançamento eq4324 velocidade angular eq4211 variaçãoeq31210 t intervalo de tempo eq4328 denominador da variável de controle singular eq31410 xvii LISTA DE SIGLAS E ABREVIATURAS CTA Centro Técnico Aeroespacial DRAG Força de Arrasto EUA Estados Unidos da América INPE Instituto Nacional de Pesquisas Espaciais xviii xix LISTA DE FIGURAS Pág FIGURA 21 Foguete lançado verticalmente na atmosfera 5 FIGURA 41 Configuração da curva cônica nas vizinhanças do foco ativo 37 FIGURA 42 Esquema da trajetória de um míssil 39 FIGURA 43 Trajetória elíptica inercial de um míssil lançado da superfície da 45 FIGURA 44 Ângulo 0 de saída do foguete 48 FIGURA 51Taxa de variação de massa versus tempo de queima 54 FIGURA 52 Velocidade do foguete versus tempo de queima 55 FIGURA 53 Altitude do foguete versus tempo de queima 56 FIGURA 54 Altitude máxima do foguete versus tempo de queima 57 FIGURA 55 Energia total gasta pelo foguete versus tempo de queima 58 xx xxi SUMÁRIO Pág 1 INTRODUÇÃO 1 11 Objetivos 2 12 Hipóteses Admitidas 3 2 ASPECTOS FÍSICO E MATEMÁTICO 5 21 Aspecto Físico 5 211 Força de Tração 6 212 Força Aerodinâmica 7 213 Força Gravitacional 8 22 Aspecto Matemático 9 221 Equação Geral do Movimento 9 222 Observação 10 3 OTIMIZAÇÃO DA ALTITUDE DE UM FOGUETE 11 31 Otimização da Altitude 11 311 Formulação Matricial do Problema 12 312 Condição Necessária para Mínimo 14 3121 Função Hamiltoniana 16 3122 Equação de Euler para Variável de Controle 18 3123 Equação de Euler para as Variáveis de Estados 20 313 Condição de Weierstrass 20 314 Condição de Controle Singular 21 315 Determinação da Seqüência de Variável de Controle Ótimo 25 4 EFEITOS DA TERRA CURVADA EM ROTAÇÃO 31 41 Equacionamento de Trajetórias no Campo Gravitacional Newtoniano 31 411 Algumas Características do Campo Newtoniano 31 412 Leis de Kepler 33 413 Equação de Energia 34 414 Equação da Trajetória 35 415 Interpretação dos Parâmetros p e e 37 416 Forma Alternativa da Equação de Trajetória 38 xxii 42 Equação da Trajetória de um Míssil em relação à Superfície da Terra 38 421 Considerações Gerais 38 422 Condições Iniciais do Lançamento 39 423 Equação Diferencial e sua Solução Geral 40 424 Solução Particular 41 425 Alcance da Trajetória 43 43 Lançamento Vertical Considerando a Rotação da Terra 44 431 Considerações Gerais 44 432 O Ponto de Impacto considerando a Rotação da Terra 46 5 RESOLUÇÃO NUMÉRICA 49 51 Exemplo 49 52 Conclusão 60 53 Sugestão para Trabalho Futuro 61 REFERÊNCIAS BIBLIOGRÁFICAS 63 1 1 INTRODUÇÃO A explicação de como os corpos se movem em suas órbitas em torno da Terra foi dada pela lei da gravitação universal de Newton no final do século XVII Segundo esta lei se um corpo fosse lançado para cima com uma velocidade suficientemente elevada ele poderia entrar em órbita da Terra Tal princípio permeou os primeiros textos de ficção científica principalmente os de Júlio Verne 18281905 que inspiraram todos os pioneiros da Astronáutica O cientista russo Konstantin Tsiolkwoski 18571935 demonstrou matematicamente a possibilidade de vôos interplanetários por meio de foguetes no artigo Exploração do Espaço Cósmico por meio de um Engenho Reativo publicado em 1896 Neste trabalho Tsiolkwoski estabeleceu a equação clássica da velocidade de um foguete em termos da razão das massas inicial e final do foguete Os primeiros foguetes utilizavam propelentes sólidos como a pólvora negra ou pólvoras de base simples nitrocelulose posteriormente passaram a empregar pólvoras de base dupla nitroglicerina e nitrocelulose e composites explosivos contendo um ligante plástico Em 1926 o cientista americano Robert H Goddard 18821945 lançou o primeiro foguete empregando propelentes líquidos Os propelentes líquidos apresentam maiores velocidades de ejeção que os sólidos A segunda guerra mundial pôs em evidência o alemão Wernher Von Braun 19121977 chefe técnico da construção dos mísseis alemães das séries V1 e V2 Von Braun continuou seu trabalho nos Estados Unidos da América nos anos iniciais da guerra fria comandando o desenvolvimento dos foguetes Saturno Nos anos 1970 Von Braun esteve no Brasil e visitou o Centro Técnico Aeroespacial e o Instituto Nacional de Pesquisas Espaciais 2 11 Objetivos O objetivo deste trabalho é determinar a altitude máxima alcançada por um foguete lançado verticalmente na atmosfera com características fixas utilizandose a taxa mássica dm dt k m na saída da tubeira do motor como variável de controle Admitindose que a velocidade de ejeção seja constante então o empuxo pode ser considerado também como a variável de controle O problema é chamado de Problema de Goddard e foi proposto no ano de 1919 Vários pesquisadores como G Hamel em 1927 Tsien e Evans 1951 Miele 1958 e Garfinked 1963 trataram deste assunto bem como Fonseca em 1977 no ITABrasil O problema aqui tratado envolve dois aspectos diferentes Aspecto Físico Consiste no estudo e análise do movimento dos foguetes na atmosfera o que corresponde fundamentalmente à Mecânica de Vôo Atmosférico Aspecto Matemático Consiste na aplicação da otimização do empuxo recorrendose aos resultados clássicos do Cálculo Variacional eou do Cálculo da Otimização dos Sistemas Dinâmicos incluindo o problema de controle 3 12 Hipóteses Admitidas A Terra é considerada como um sistema de referência inercial plano No caso dos foguetes que se movem nas vizinhanças imediatas da Terra por exemplo dentro da atmosfera terrestre a influência de outros corpos celestes é negligenciável em relação à força gravitacional terrestre São desprezadas as forças menores como a força centrífuga e a de Coriolis comparadas com as forças de Tração Aerodinâmica e Gravitacional O foguete sendo lançado verticalmente na atmosfera tem o ângulo de ataque 0 devido à sua simetria axial Logo a força de sustentação é 0 L h No conjunto das forças não foram incluídas forças elásticas devidas às oscilações de foguete e também as forças devidas a mudanças da posição do centro de massa de foguete A taxa mássica que sai da tubeira do motor do foguete será representada pela variável de controle dm dt k m A força P h provocada pela diferença de pressão entre a pressão eP na saída da tubeira do motor de foguete e a pressão atmosférica 1 h P na altitude h é uma parcela de caráter corretivo 4 5 2 ASPECTOS FÍSICO E MATEMÁTICO 21 Aspecto Físico Queremos determinar a altitude máxima possível que um foguete lançado verticalmente pode alcançar controlando seu consumo de combustível através da variável dm dt k m que representa a taxa mássica que sai da tubeira do motor do foguete com parâmetros como altitude h t pressão P h densidade h temperatura T h velocidade t v arrasto CD h e a massa m t As variáveis serão estabelecidas de acordo com o modelo conveniente da atmosfera A Figura 21 mostra um esquema de um foguete lançado verticalmente na atmosfera FIGURA 21 Foguete lançado verticalmente na atmosfera 6 As forças que agem sobre o foguete em movimento lançado verticalmente na atmosfera são F força de tração 211 D força de arrasto 212 W força gravitacional 213 As outras forças como as centrifugas e as de Coriolis são muito menores portanto podem ser desprezadas se comparadas ás três primeiras 211 Força de Tração É a força produzida pelo motor do foguete que atua na direção do movimento P h T F 2111 com m E T 2112 e e P h A P P h 1 2113 Substituindo 2111 e 2112 em 2113 vem e e E P h A P m F 1 2114 onde sin 0 1 de ejeçãodos gases velocidade é decorrente do fato da massa decaircom o tempo negativo al o taxa de va riação da massa do foguete m m saidadatubeiradomotor da área transversal a de A altitudeh acompanhad na pressão atmosférica motor de foguete e P h do pressãona saida da tubeira diferençaentrea P à força provocadano se ntidodomovimento devido h P força doempuxodevidoao fluxodos gases T E e e 7 212 Força Aerodinâmica É a força exercida pelo ar sobre o foguete em movimento Em nosso problema a única força aerodinâmica que age sobre o foguete em movimento é a força de arrasto Esta força é a componente da força aerodinâmica na direção oposta à da trajetória do vôo representada pela letra D Drag 2 1 2 M SC h h D D 2121 e 2 2 1 S h D CD M 2122 Onde velocidadedo somnoar h a númerodeMach do foguete a h M coeficientedearasto é funçãodonúmerodeMach M C velocidadedo fogueteem relação aoar t cte seção dereferência do foguete S densidade atmosférica h força dearasto h D D e Notese que para M fixo e valores do ângulo de ataque pequenos C L é uma função aproximadamente linear de sendo que no caso de foguetes devido à simetria axial 0 conseqüentemente C L 0 e a força de sustentação 0 L h 8 213 Força Gravitacional É a força exercida sobre o foguete por todos os corpos celestes sol planetas lua etc Ela age na direção do centro de massa do corpo em questão Esta força é chamada de Peso e é designada por W A força gravitacional na superfície da Terra é dada por T T R GM m mg W 2 0 0 2131 A força gravitacional na superfície da terra que age num corpo a uma altitude h é da forma h 2 R GM m mg h W T T 2132 dividindo 2132 por 213 1 temos 2 0 h R R g h g T T 2133 onde 0 0 distância entre ocentro de massa do foguete e a Superficie daTerra h aceleração da gravidadena Superficie da Terra g aceleração da gravidade no ponto onde se encontra ocorpo h g força gravitacional terrestreno ponto onde seencontrao corpo W força gravitacional na Superficie daTerra W massa da Terra M constaante gravitacional G raio daTerra R massa do foguete m T T 9 22 Aspecto Matemático O estudo deste trabalho possibilitará a apresentação da análise da trajetória do foguete sob o ponto de vista matemático recorrendose ao chamado Cálculo Variacional O problema de cálculo variacional é determinar os valores extremos máximos e mínimos de uma função 221 Equação Geral do Movimento Com base na 2º Lei de Newton a Equação 21 pode ser escrita na forma F m 2211 onde D W P h T F 2212 Igualando as Equações 2212 e 2211 temos W D P h m m E 1 2213 Admitindo que W m h W h D D P h P de m depende dt dm k m E 2214 Finalmente 10 4 122 4 122 4 122 1 c dt dm k m b h a W m h h D P h k m E As Equações acima formam o sistema de equações diferenciais do movimento do foguete 222 Observação A tubeira é um bocal que tem por objetivo acelerar os gases produzidos na câmara de combustão até velocidades supersônicas A tubeira é dividida em três partes o convergente a garganta e o divergente O escoamento dos gases de combustão ocorre da seguinte forma Câmara de Combustão Convergente Garganta Divergente Na garganta da tubeira a velocidade dos gases não pode ultrapassar a velocidade sônica portanto a variável de controle k da Equação 2214c deve ser limitada fisicamente logo max 0 k k 2221 11 3 OTIMIZAÇÃO DA ALTITUDE DE UM FOGUETE 31 Otimização da Altitude O problema de maximizar ft h nos leva a minimizar o funcional J expresso através das Equações 221 4a b e c h ft h m k J 311 sujeito aos seguintes vínculos 0 h h m m k i 312 onde 2 13 0 2 13 0 2 13 0 1 3 2 1 c k m h m m k h b h h m m k h a W m h h D P h m k h h m m k E As condições iniciais 3 13 0 3 13 0 0 3 13 0 0 0 c m m b h a e final 12 4 13 4 13 b valora ser máximizado t h a valor conhecido m t m f f f Sendo que as variáveis de estado hm e a variável de controle k são variáveis dependentes funções do tempo e o tempo estático t a variável independente 311 Formulação Matricial do Problema Sejam 3 2 1 e x t x t x t as coordenadas das variáveis de estados hm representadas da seguinte forma 3 2 1 t m t h t v t x t x t x xi t 3111 onde i 1 23 Logo temse que t k t x x t W x t x t D x t P x t k t t x t m t h t t x t x t x t x E i 1 1 3 2 2 1 2 3 3 2 1 3112 13 Admitindo t k t x x t W x t x t D x t P x t k t t x t g E i 1 1 3 2 2 1 2 3 3113 temos g t t x i i 3114 ou ainda x t x t x t k t g x t x t x t k t g x t x t x t k t g t x t x t x t xi 3 2 1 3 3 2 1 2 3 2 1 1 3 2 1 3115 Tendo em vista a Equação 312 a representação matricial tem a forma 0 g t x t h h m m k i i i 3116 ou seja 0 0 0 3 2 1 3 3 2 1 2 3 2 1 1 3 2 1 3 3 2 2 1 1 x t x t x t k t g x t x t x t k t g x t x t x t k t g t x t x t x t x t x t x t x t x t k t x 3117 14 312 Condição Necessária para Mínimo Queremos determinar a função k kt que minimiza fh ou maximiza fh a altitude do foguete Usando a Equação 311 e os Multiplicadores de Lagrange temos ft t i i T f dt h t x t x t x t k t t J 0 3 2 1 3121 onde f f x t h t 2 3122 assim f f t t i T i t t i T i f g t dt x t dt x t x t x t x t k t t J 0 0 2 3 2 1 3123 Usando a matriz transposta aos Multiplicadores de Lagrange temos 3 2 1 T 3 2 1 iT 3124 sendo que integrando por parte vem f f f t t i T i t t i T i i t t iT x t dt x t t dt x 0 0 0 3125 e substituindo a Equação 3125 em 3123 temos f f t t i T i i T i t i t i T f x t dt g t x t x x t x t x t k t t J 0 0 2 3 2 1 3126 A função hamiltoniana do sistema acima é dada pela relação g t H i T i i 3127 assim substituindo a Equação 3127 em 3126 obtémse 15 f f t t i T i i t i t i T f x t dt H x t x x t x t x t k t t J 0 0 2 3 2 1 3128 Para que o funcional J seja um mínimo é necessário que a primeira variação seja nula isto é J 0 Admitindo que f f f f f m t x livre t x t x livre t fixo t 3 2 1 0 0 3129 temos 0 x t dt t x J dt d t x J J f 0 t t i i i 31210 0 t dt k k t H t x x t H 1 x t t t t H t x t J f 0 t t i i T i i f 2 f 2 f f f 2 f 31211 uma vez que xi t e k t são independentes Anulando cada termo da Equação 31211 em relação a zero temos 15 213 0 k t H 14 213 t x H 13 213 1 t 12 213 t x t t H T i i f 2 f f 2 f 16 3121 Função Hamiltoniana Para resolver nosso problema recorrermos ao elemento mais importante da otimização a função Hamiltoniana descrita por Pontryagin 7 cuja sua expressão é da forma k W m h h D P h k m g k t h m H E m h v i iT 1 31211 logo h m E k mg h D P h k m k t h m H 31212 Para saber o comportamento da função hamiltoniana ao longo da trajetória ótima basta derivar a Equação 31212 em relação ao tempo Assim temos t k t k g x t x t g g g g dt dH T i i i i T i i T i i T i i iT 31213 Sabemos que t x g x t H k t g k t H i 31214 17 Levando as expressões 31214 na Equação 31213 vem t k t k H x t x t H x t g H dt dH i i i 31215 Tendo em vista que 0 dk t dH g t x i i 31216 e substituindo 31216 em 31215 temos dt 0 dH 31217 ou ainda 0 2 f f f f t t x t H t 31218 como 0 ft H 31219 então H 0 312110 Concluiremos que ao longo da trajetória ótima a função Hamiltoniana é nula 18 3122 Equação de Euler para Variável de Controle Desejase transformar o intervalo fechado da variável de controle k da Equação 2215 em um intervalo aberto para que a nova variável de controle assuma qualquer valor real Admitimos 2 k kmax cos R 31221 e substituindo 31221 em 2215 vem max 2 max cos 0 k k 31222 onde 1 cos 0 2 31223 Derivando parcialmente em relação à nova variável de controle vem 0 k H k H 31224 Assim 0 0 2 sen 1 max H k k 31225 19 Substituindo 31221 em 31212 temos 2 max 2 max cos cos k g m h D m h P m k k t h m H m h E 31226 e derivando parcialmente 31226 em relação a vem 0 2 sen 0 0 max m E m k k H 31227 Analisando a Equação 31227 são definidos os possíveis valores ótimos da variável de controle k 3 2 1 0 2 0 sen 2 0 m E m 31228 Substituindo a Equação 31228 em 31221 obtemse uma solução da forma 10 2213 0 531 9 2213 420 2 cos max 2 max k para k k para k k ímpares pares 0 m E m 312211 O caso 312211 não pode ser determinado pois não fornece nenhuma informação de como obter a variável de controle k Este caso será objeto de estudo no próximo Capítulo 20 3123 Equação de Euler para as Variáveis de Estados Derivando a função Hamiltoniana em relação a cada variável de estado temos 3 3213 2 3213 1 3213 2 3 2 1 h D P h k m mg D P m m D x H x H x H E m h h h h h 313 Condição de Weierstrass A condição necessária do mínimo para satisfazer o controle ótimo é dada pela seguinte equação 0 t H x k t H x k H 3131 0 0 1 1 k k W D P m k W D P k m H E E m h v 3132 Para h 0 temos 0 k k m H m E 3133 O termo m E m é chamado de função comutadora Switching Function e é denotado como S t m E m S t 3134 21 A função comutadora S t é o processo que permite determinar a altitude máxima do foguete em três casos possíveis 7 313 0 0 3 6 313 0 2 0 5 313 0 1 max max S t se k S t se k k S t se k k 8 Notese que os casos 3135 e 3137 definem os subarcos não singulares embora seja importante ressaltar que quando a variável de controle ótima passa bruscamente do valor máximo para o valor mínimo sem passar pela variável de controle ótima singular tratase de Problema do tipo Bang Bang O caso 3136 chamase de controle singular ou subarcos singulares 314 Condição de Controle Singular Queremos determinar a solução da Equação 3136 calculando a primeira derivada da função comutadora S t até que se encontre explicitamente a variável de controle ótima procurada Derivando a função comutadora t S em relação ao tempo temos m E E m E m m m m dt d S t 2 3141 Introduzindo as expressões e m na equação vem 22 h D D P h m m S t E h E 2 3142 se S t 0 e m 0 temos h D D h P m E h E 3143 Podese ver que na primeira derivada de t S nenhum termo da variável de controle k foi encontrada Logo calcularemos em seguida a segunda derivada de t S e assim por diante até que seja explicitamente encontrada a variável de controle ótima de k procurada Assim h D D P h m dt d S t dt d S t E E h 2 2 3144 ou ainda h D D P h h D D P h m m S t E E E h E h 3145 mas D h D D h 3146a D h D h D h 3146b P h P h h 3146c Substituindo as Equações 3146a b c e 3143 em 3145 temos 23 h D D P h k h D D h P D D P m h D D h P W m h h D P h k D D h D D h P D h D D P h mg D P m t S E E h E E h h h E h E E E E h E E E h h h E h 2 2 3147 Para h 0 e S t 0 temos P h h D D D h D h P P h h D D D W m h h D P h D D P h h D D D W m h h D P h D D D P m mg D P m k E E E E E E E E E E h h h h h h E 2 2 2 2 2 2 2 2 2 2 2 3148 Finalmente manipulando a Equação acima a variável de controle procurada na Equação 3136 apareceu explicitamente em 3148 e será chamada de variável de controle singular S k P h h D D D h D h P P h h D D D W m h h D P h D D P h h D D D W m h h D P h D D D P m mg D P m k E E E E E E E E E E h h h h h h E S 2 2 2 2 2 2 2 2 2 2 2 3149 Admitindo que 24 P h h D D D E E 2 2 31410 e levando a Equação 31410 em 3149 vem W m h h D P h D D W m h h D P h D h D h P h D D m mD m P g m D m P h m k E E E h E h h E h E E S 2 2 1 1 2 2 2 31411 Resumido 13 413 0 0 3 6 313 0 2 12 413 0 1 max S t se k vide eq controle singular S t se k k S t se k k S 25 315 Determinação da Seqüência de Variável de Controle Ótimo Para que possamos calcular a altitude máxima de um foguete lançado verticalmente na atmosfera devese primeiramente determinar a seqüência de variável de controle ótima a ser aplicada no instante inicial do movimento do vôo Desejase determinar os multiplicadores de Lagrange das Equações anteriores 31212 313 4 e 3142 que formam um sistema de três equações com três incógnitas conforme mostrado abaixo Resolvendo o sistema temos 1 513 0 0 1 513 0 0 1 513 0 c m b m h D D h P a mk m mg h D P h k m E E h E m h E e onde m m h D D h P mk m mg h D P h k Det E E E E 0 0 3152 0 2 2 2 2 E E E E m k h D D P h m mg h D P h k m Det 26 Para m 0 e Det 0 temos h D D h P mg h D h P E E 3153 Substituindo a Equação 3153 no sistema de equação acima vem 0 E h E m D h D P h 3154 Introduzindo a Equação 3154 no sistema de equações acima chegase às seguintes 1 h E E m E E D h D h P D h D h P m 3155 levando as expressões m e h da Equação 3155 nas Equações 31212 313 4 e 3142 Obtemse 0 0 0 t H t S t S 3156 Convem salientar mais uma vez que quando tivermos a condição 3153 teremos as condições da Equação 3156 27 Então ao longo da trajetória ótima 0 0 t S t H 3157 Neste caso como já foi dito antes a variável de controle a ser aplicada será obviamente k S k vide eq 3136 Para m 0 e Det 0 temos h D D h P mg h D h P E E 3158 aplicaremos as variáveis de controle seguinte max k k e k 0 vide eq 31412 e 31414 É interessante notar que para determinar a seqüência de variável de controle a ser aplicada ao movimento do foguete no instante inicial t 0 analisaremos a Equação 3153 recorrendo às Equações 313 a 313b e 313c Aplicando as condições iniciais na Equação 3153 vem 0 0 0 0 0 0 0 0 0 0 0 h D h D P h W m h D P E E 3159 Como 00 0 0 h D 00 0 0 h D 31510 28 então 0 0 0 0 0 P W m E P 31511 Recorrendo à Equação 2111 vem e e E P h A P m P h T F 1 1 31512 onde 0 força gravitacional W m h a e muito menor que m que a parcela estático T do muito menor t se ndo no in st ante inicial corretivo à diferença de pressão é uma parcela à caráter devido força causada na saida da tubeira do motor do foguete h P E e finalmente 0 0 0 0 E W m P 31513 Portanto para que haja movimento no instante inicial t 0 do vôo a taxa mássica que sai da tubeira do motor tem que ser diferente de zero 0 k m o que confirma a condição da Equação 3158 ou seja 0 0 0 0 0 P W m E P 31515 Logo 29 max max 0 0 0 0 0 0 k k k k k k S 31516 Finalmente da Equação 31516 concluiremos que a variável de controle a ser aplicada no inicio do movimento é max k k k 30 31 4 EFEITOS DA TERRA CURVADA EM ROTAÇÃO Neste Capítulo são estudados os efeitos da rotação da Terra sobre a trajetória de um foguete lançado verticalmente 41 Equacionamento de Trajetórias no Campo Gravitacional Newtoniano 411 Algumas Características do Campo Newtoniano Como se sabe o equacionamento do movimento em campo central de forças do tipo newtoniano é dado pela equação re r S F r ˆ 2 4111 onde tan ˆ te fator de proporcionalidade cons S versor radial dasbases do mesmo sistema e coordenada radial do sistema cilíndrico r força do campo r F r A Equação 4111 pode ser escrita na forma r2 S F r 4112 com F r r F 4113 32 Nota 1 A origem do sistema de coordenadas é colocada no centro de massa Terra míssil Como a massa do míssil é m MT massa da Terra para todos os fins práticos podemos admitir que a origem esteja colocada no centro da Terra Veremos mais adiante que para as velocidades iniciais que não ultrapassam a velocidade de escape do míssil este vaise mover ao longo de uma elipse cujo foco ativo vai coincidir com o centro da Terra É fácil ver que o campo representado pela Equação 4111 é potencial De fato o potencial é r S V 4114 pois F r r S r V 2 4115 Recorrendo à variável u r r u 1 1 4116 e substituindo a Equação 4114 em 4112 temos 2 1 Su F u 4117 No campo Newtoniano atrativo tal como o campo gravitacional temse S 0 Comparandose as Equações 2131 e 5112 temos GM m S r S r GM m F r T T 2 2 4118 33 412 Leis de Kepler Como se sabe da Mecânica Celeste a trajetória no campo do tipo da Equação 4117 satisfaz as chamadas leis de Kepler 1 Todos os planetas movemse em torno do Sol em órbita elíptica estando o Sol em um dos dois focos chamado de ativo 2 O vetor posição do planeta em relação ao Sol varre a mesma área em intervalos iguais de tempo 3 O quadrado do período de revolução do planeta é diretamente proporcional ao cubo da distancia média do Sol Como uma das conseqüências das leis de Kepler temos As órbitas dos planetas são planares assim o equacionamento do movimento pode ser expresso num sistema de coordenadas polares com a origem no foco ativo do campo gravitacional A chamada velocidade areal isto é referida à área varrida pelo vetorposição é constante Ela pode ser expressa pela chamada quantidade de movimento angular h Nota 2 As leis de Kepler valem para qualquer campo gravitacional Por exemplo no caso de um míssil movendose no campo terrestre a palavra Sol deve ser substituída pela Terra e a palavra Planeta por Míssil por exemplo 34 413 Equação de Energia Sendo o campo gravitacional potencial a energia mecânica global é preservada Temos 1 A equação da energia cinética é dada por 2 2 2 2 1 r m r Tc 4131 onde pela segunda lei de Kepler temos h r2 constante 4132 Assim temos 2 2 2 2 1 r h m r Tc 5133 2 A energia potencial V é dada pela Equação 4114 Assim a energia mecânica total é r S r h m r V T E c mTot 2 2 2 2 1 4134 35 414 Equação da Trajetória como r r h r dr d 2 4141 associando as Equações 4134 e 4141 temos 2 2 2 2 r h r S E m hmr dr d 4142 ou então h m mSr mEr dv hmr 2 2 1 0 2 2 4143 introduzindo a variável r u 1 4144 resulta que 2 2 2 0 2 2 u mh Su mh E du 4145 Chamando 2 2 mh E 2 2 mh S 1 e recorrendo às tabelas de integrais temos 4 2 arccos 1 2 0 u 36 ou seja 2 2 2 2 2 0 2 arccos mh E h m S mh S u 4146 onde 0 é uma constante de integração que pode ser determinada pelas condições iniciais Voltando à variável r vide eq 4144 obtemse 0 cos 1 e p r 4147 onde S mh p 2 4148 2 2 2 1 S Eh e 4149 O parâmetro p é chamado de semi latus rectum e e de excentricidade em latim latus rectum significa lado reto Como se sabe da Geometria Analítica a Equação 4147 representa as chamadas curvas cônicas elípticas incluindo circunferência parábola ou hipérbole com a origem colocada num dos focos chamado ativo da curva 37 415 Interpretação dos Parâmetros p e e Observemos que a Equação 4147 é expressa em coordenadas polares r Notemos o valor mínimo de r corresponde ao valor 0 neste caso e p r r 1 min 4151 Assim e p 1 é a distancia mais curta entre o foco e a curva Notemos ainda que se 2 0 e e 1 temse p r r Lr semi latus rectum 4152 Desta maneira podemos admitir 0 0 e o eixo x a direção da distância mais curta A Figura 41 mostra um esquema da configuração da curva cônica nas vizinhanças do foco ativo FIGURA 41 Configuração da curva cônica nas vizinhanças do foco ativo A expressão e excentricidade denota o tipo de curva cônica Assim temos 1 0 e Elipse quando e 0 circunferência e 1 Parábola e 1 Hipérbole 38 Em nosso caso só interessam trajetórias fechadas isto é com 1 0 e 416 Forma Alternativa da Equação de Trajetória Diferenciando duas vezes a Equação 4145 em relação a u obtemse uma equação diferencial linear de ordem 2 extremamente simples C u d d u 2 2 constante 4161 A solução da Equação 4161 é imediata C B Asin r u cos 1 4162 da qual recorrendose às condições iniciais e calculando ABC obtemse a solução 4161 Embora a solução geral seja simples o cálculo das constantes de integração apresenta certa dificuldade 42 Equação da Trajetória de um Míssil em relação à Superfície da Terra Curvada sem Movimento 421 Considerações Gerais Qualquer míssil lançado da superfície da Terra movese num campo gravitacional Terrestre e portanto após o instante de burnout não levando em conta a resistência do ar movese ao longo da trajetória cônica com um dos focos o ativo colocado no centro da Terra 8 O movimento parabólico só poderia ser considerado se a superfície da Terra fosse sem curvatura A trajetória do tipo em questão sofrerá a ação dos fatores diferentes 1 A própria curvatura da superfície como já foi mencionado 39 2 A rotação Terrestre com a velocidade angular h rad 24 2 4211 que agirá impondo uma componente tangencial ao movimento do míssil em relação à superfície fixa 422 Condições Iniciais do Lançamento A Figura 42 mostra um esquema da trajetória de um míssil FIGURA 42 Esquema da trajetória de um míssil Como na Seção anterior introduziremos um sistema polar de coordenadas r com a origem no centro da Terra Sejam 40 sup 0 0 inicial datrajetória em relação àlinha horizontal local ângulo erficie altitude do ponto inicial da trajetória em cima da H raio da Terra R out posição inicial do foguete em relação ao centro da Terra no burn r As expressões t r t e são as mesmas que 0r e 0 só que num instante genérico t Escolhendo um certo instante como inicial t0 é fácil verificar que no instante inicial temos t 0 0 0 0 r r 4 221 3 224 cos 2 224 sin 0 0 0 0 0 0 polar da velocidade componente r componente radial da velocidade r O valor de h pode ser também determinado das condições iniciais embora permaneça invariante durante toda a trajetória 0 0 0 2 0 0 0 0 2 cos cos r r r r h constante 4224 423 Equação Diferencial e sua Solução Geral Levando em conta a Equação 4161 temos C u d d u 2 2 4231 41 onde 0 2 2 0 2 0 2 cos r GM mh GMm C constante 4232 Introduzindo uma variável adimensional GM r 2 4233 temos da Equação 4232 0 2 0 0 0 2 0 2 0 0 cos 1 1 cos 1 r r r GMm C 4234 A Equação 4231 é obviamente linear não homogênea de ordem 2 com coeficientes constantes Como é fácil verificar a sua solução geral é C B Asin r u cos 1 4235 onde A e B são constantes de integrações e C é dada por 4232 424 Solução Particular Diferenciando 4235 em relação ao tempo temos B sin A r r u cos 1 2 4241 Substituindo as Equações 4221 em 4235 temos 0 0 1 cos 0 0 1 r C B C B Asin r u 4242 ou então 42 0 2 0 0 0 2 0 0 0 0 cos 1 1 cos 1 1 1 1 Br r r C r B 4243 Substituindo as Equações 4222 em 4241 temos 0 cos cos 0 cos 0 0 0 0 0 0 2 0 0 0 sin B r r A r sin u 0 0 0 0 0 0 0 0 0 0 cos sin cos sin tg Ar r A r A 4244 Assim a Equação 4235 fica 0 2 0 0 0 2 0 0 0 0 cos 1 1 cos cos 1 1 1 1 r r sin r tg r u 4245 ou então 0 2 0 0 0 0 2 0 0 0 2 0 0 0 0 cos cos 1 cos cos cos 1 1 cos cos 1 1 1 sin sin r r sin tg r r 0 2 0 0 0 0 0 cos cos 1 cos cos cos sin sin r r 4246 assim obtemse 0 2 0 0 0 0 cos cos 1 cos cos r r 4247 Nota 3 Como se sabe da Mecânica Celeste esta é a equação de qualquer curva cônica com foco ativo no ponto r 0 centro da Terra Curvas cônicas não 43 degeneradas são elipse incluindo o círculo parábola e hipérbole As duas ultimas constituem trajetórias que se estendem até o infinito Neste trabalho estudaremos somente trajetórias fechadas elipse ou no caso particular a circunferência 425 Alcance da Trajetória O alcance da trajetória max é obtido quando o míssil tocar na superfície terrestre isto é quando r R Assim a relação entre max e as condições iniciais é dada por 0 2 0 max 0 max 0 0 cos cos 1 cos cos r r 4251 onde as grandezas 0 0 0 r e denotam as condições iniciais Se a trajetória começar na superfície terrestre a relação entre max e os valores iniciais é 1 cos cos 1 cos cos 0 2 0 max 0 max 0 4252 A relação acima pode fornecer uma dependência explícita entre o alcance max e as condições iniciais Realmente 4252 pode ser escrito da seguinte forma 0 2 0 max max 0 max 0 0 0 cos cos 1 cos cos cos sin sin ou então 0 2 0 max 0 0 0 max 0 2 0 cos 1 1 cos 1 cos cos sin sin 4253 44 Observemos que a equação acima tem a forma 0 0 max 0 0 max 0 0 cos Q sin N M 0 e assim max 2 2 2 2 max max 2 max cos cos 0 cos 1 cos N N Q M Q N M 0 cos cos 2 cos max 2 2 2 max max 2 2 N N MQ M 4254 e finalmente 0 cos 2 cos 2 max max 2 2 2 N MQ N M 4255 A Equação 4255 é uma equação quadrática da qual se pode calcular cosmax ou termos das expressões M N Q que dependem apenas das condições iniciais Nota 4 Os valores de M N e Q são 8 524 cos 1 7 524 sin cos 6 524 cos 0 2 0 0 0 0 0 2 0 Q N M 43 Lançamento Vertical Considerando a Rotação da Terra 431 Considerações Gerais 45 Notemos que o movimento rotacional da Terra não influencia a trajetória do míssil descrita nas Seções 431 e 432 no plano inercial da força central pois o centro da Terra não se move em relação à superfície Assim no plano inercial o míssil movese ao longo de uma cônica com foco ativo no centro da Terra Ademais se a trajetória for uma curva fechada ela tem que ser elíptica A Figura 43 mostra um esquema da trajetória elíptica inercial de um míssil lançado da superfície da Terra Curvada esférica FIGURA 43 Trajetória elíptica inercial de um míssil lançado da superfície da Terra Curvada esférica Na Seção anterior a velocidade inicial do míssil foi determinada em relação à superfície da Terra tendo duas componentes 1 A vertical 0 V 2 A tangencial tV Notemos que as duas componentes determinam também a O ângulo 0 pois 46 V V tg 0 0 4311 b O versor neˆ normal ao plano inercial isto é plano aonde é aplicado a 2a Lei de Newton t t r n V e V e e ˆ ˆ ˆ 0 4312 Conhecendose estes dois elementos e a velocidade inercial temos t t r V e V e ˆ 0 ˆ 0 4313 bem como as coordenadas geográficos do sitio do lançamento a equação da trajetória pode ser imediatamente determinada incluindo o alcance max 432 O Ponto de Impacto considerando a Rotação da Terra A velocidade linear na direção tangencial é dada pela equação t R 4321 onde módulo da velocidade angular da Terra da Equação 4211 Notemos que para um observador no plano inercial valem os resultados das seções anteriores deste Capítulo Porém um observador fixo na superfície vê o míssil sujeito às forças de Coriolis e às centrífugas A velocidade resultante 0 inercial no lançamento é composta da velocidade de V0 do foguete em relação à superfície da Terra mais a velocidade tV imposta pelo movimento da Terra em relação ao plano inercial onde obviamente temos 47 0r Vt 4323 com 0 0 r cos r 4324 onde 0 é a latitude do sítio do lançamento Assim 0 0 0 r V 4325 A variável adimensional inicial 0 fica GM r V r GM r 2 0 0 0 2 0 0 0 4326 Para lançamento vertical em relação ao observador na superfície o ângulo 0 pode ser determinado da seguinte equação 0 tg 0 4327 A Figura 44 mostra um esquema do ângulo 0 de saída do foguete 48 FIGURA 44 Ângulo 0 de saída do foguete O versor que normaliza o plano inercial pode ser determinado pela Equação 4312 O alcance max pode ser determinado da Equação 4255 Quanto ao ponto de impacto isto é latitude e longitude temos que levar em conta que durante o percurso da trajetória o intervalo t a posição da Terra em relação ao plano inercial fica deslocada de ângulo t 4328 Para determinar o tempo t decorrido durante a trajetória podese recorrer à Equação 4132 49 5 RESOLUÇÃO NUMÉRICA 51 Exemplo Queremos através de um exemplo resolver o problema de otimização da altitude de um foguete hipotético lançado verticalmente na atmosfera com os seguintes dados Massa inicial do foguete kg m 0 1000 Massa de combustível kg mc 750 Massa final kg m m f 250 4 1 0 Tempo de queima s tq 60 Considerar constante m s de escape dos gases velocidade s m g gravitacional aceleração k k demassa fluxo E 2000 9 81 2 0 max Determinar nos diferentes casos considerados 1 O fluxo de massa 2 A velocidade q q t do foguete na fase final 1 t qt 3 A altitude q hq t no instante t qt 4 A altitude máxima do foguete 5 A energia total do foguete q q q q f f gh t t m E 2 2 1 6 Apresentar a tabela dos resultados para todos os casos considerados 50 Resolução Pela 2a Lei de Newton F m 511 Como sabemos um foguete movendose no espaço livre isto é no vácuo sem a pressão Ph e sem ângulo de ataque os resultantes das forças durante o vôo são dados pela seguinte equação mg0 m W T F E 512 Igualando as Equações 511 e 512 temos g0 dt dm m dt d E 513 A Equação 513 pode integrada imediatamente entre t 0 e t genérico t t E t g dt m dm d 0 0 0 0 514 As condições iniciais são 0 0 0 0 m m 515 assim g t m t t E 0 0 ln 0 516 51 logo g t m t m t E 0 0 ln 517 A velocidade t do foguete no instante t arbitrário da fase1 0 t qt g t m t m t E 0 0 ln 518 onde massa combustivel m burn out te final da fase propulsiva ins t massa final do foguete istoé sem propelente m m m c q f c tan 0 A altitude h t do Foguete no Instante t arbitrário da fase 1 0 t qt O caminho t h percorrido pelo foguete é obviamente t t dt t h 0 519 ou então t t E g tdt m t dt m t h 0 0 0 0 ln 5110 É obvio que c f q t m m m t m m dt m t m 0 0 0 5111 52 Substituindo 5111 em 5110 vem t t E g tdt mtdt m m t h 0 0 0 0 0 ln 5112 Como sabemos a velocidade q q t não depende da distribuição temporal da função m t e sim da quantidade do propelente gasto f c m m m 0 assim o caminho hq tq é dada pela seguinte equação q q t t E q q g tdt m t dt m t h 0 0 0 0 ln 5113 Fazendose uma mudança de variáveis temos t q t E q g tdt t dt m mt t t h 0 0 0 0 ln 1 5114 q q q q dXt dt Xt t t dt dX t t X 5115 Substituindo 5115 em 5114 resulta t t q E q g tdt dX m mXt t t h 0 0 0 0 ln 1 5116 53 Admitindo que 0 0 0 1 1 1 0 1 m mt Y X Y X dX m mt dY m mXt Y q q q 5117 Temos t m t m E g tdt YdY m m t h q 0 0 1 1 0 0 ln 5118 usando a tabela de integração dos logaritmos temos t q q q E g tdt m mt t m m t t h 0 0 0 0 ln 1 5119 54 A equação da taxa de variação de massa m é q c t m m k max 5120 A Figura 51 mostra um gráfico da taxa de variação de massa versus tempo de queima 0 10 20 30 40 50 60 0 100 200 300 400 500 600 700 800 taxa de variação de massa kgs tempo de queima s FIGURA 51Taxa de variação de massa versus tempo de queima 55 A equação da velocidade do foguete q q t no final da queima é q f E q q g t m m t 0 0 ln 5121 A Figura 52 mostra um gráfico da velocidade do foguete versus tempo de queima 0 10 20 30 40 50 60 2100 2200 2300 2400 2500 2600 2700 2800 velocidade do foguete ms tempo de queima s FIGURA 52 Velocidade do foguete versus tempo de queima 56 A equação da altitude do foguete no final da queima q hq t é 2 ln 1 2 0 0 q f c f E q q q g t m m m m t h t 5122 A Figura 53 mostra um gráfico da altitude do foguete versus tempo de queima 0 10 20 30 40 50 60 0 10k 20k 30k 40k 50k altitude do foguete km tempo de queima s FIGURA 53 Altitude do foguete versus tempo de queima 57 A equação da altitude máxima do foguete hmax é q q q q h t g t h 2 2 max 5123 A Figura 54 mostra um gráfico da altitude máxima do foguete versus tempo de queima 0 10 20 30 40 50 60 2800k 3000k 3200k 3400k 3600k 3800k 4000k altitude máxima km tempo de queima s FIGURA 54 Altitude máxima do foguete versus tempo de queima 58 A equação da energia total gasto do foguete é q q q f T gh t m E 2 2 1 5124 A Figura 55 mostra um gráfico da energia total gasto do foguete versus tempo de queima 0 10 20 30 40 50 60 750E008 800E008 850E008 900E008 950E008 energia total gasta J tempo de queima s FIGURA 55 Energia total gasta pelo foguete versus tempo de queima 59 A Tabela 51 a seguir mostra os resultados de todos os casos considerados TABELA 51 Resultados dos casos considerados Tempos s kmax s kg q tq kms hq tq km hmax km ETot 106 x J 0 0 0 0 0 1 750 2763 1071 3901 9568 5 150 2724 5256 3833 9404 10 75 2674 1027 3748 9204 15 50 2625 1503 3664 9009 20 375 2576 1955 3579 8820 25 30 2527 2383 3494 8636 30 25 2478 2786 3409 8458 35 214 2429 3164 3324 8285 40 188 2380 3518 3239 8117 45 167 2331 3848 3155 7955 50 15 2282 4153 3007 7798 55 136 2233 4433 2985 7647 60 125 2184 4689 290 7501 5125 60 52 Conclusão No início do movimento do foguete o valor da variável de controle a ser aplicado deve ser k kmax Este valor deve ser mantido até que a função comutadora t S e a velocidade atinjam a igualdade ou seja S t Quando acontecer a condição da Equação 31413 aplicaremos então a lei de variação singular de controle ao sistema k S k dt S t d dt t dS 2 2 0 0 Quando a função comutadora é nula S t 0 a variável de controle singular deverá permanecer atuante até que toda a quantidade total de combustível se esgote ou até que seu valor se anule Caso a variável de controle singular se anular 0 k S antes de combustível se esgotar não teria sentido físico a aplicação da variável de controle singular Após o instante em que 0 k S a aplicar k 0 o foguete perderá energia cinética sob a ação das forças gravitacional e de arrasto Portanto não é viável a aplicação de k 0 após que 0 k S uma vez que isso fugira da solução ótima Mas se o combustível se esgotar durante a aplicação da variável de controle singular obviamente o próximo e último valor ótimo que a variável poderia assumir é k 0 neste caso a seqüência de controle ótima a ser aplicado será 0 max k k k k k S No final da queima t qt toda a quantidade de combustível está esgotada portanto o foguete vai voar com sua inércia ou seja a energia armazenada ou adquirida até o final da queima 61 No limite quando o tempo de queima qt 0 a aceleração isto deixa claro que o tempo de queima qt não pode ser arbitrariamente curto Para se obter um vôo rigorosamente vertical de um míssil será necessário introduzir continuamente correções na variação de massa dos gases de saída o que não seria prático Assim um lançamento vertical pode ser considerado apenas como uma aproximação muito boa para altitude de ordem de até 10 km Para ter noção do desvio do lançamento vertical podese recorrer às equações válidas para sistemas de coordenadas fixas na superfície da Terra introduzindose as forças de Coriolis e centrífugas 53 Sugestão para Trabalho Futuro Minimizar o tempo gasto para se obter o alcance máximo de um foguete a partir de uma dada quantidade de combustível 62 63 REFERÊNCIAS BIBLIOGRÁFICAS 1Azimov MD Intermediatethrust arcs in mayers variation problem Journal of Applied Mathematics and Mechanics v 64 n 1 p 8795 2000 2Azimov DM Analytic solutions for intermediatetrust arcs of rocket trajectories in a newtonians field PMM Journal of Applied Mathematics and Mechanics v 60 n 3 p 421427 1996 3Azizov G A Korshunova NA An analytical solution of the optimum trajectory problem in gravitational field Celestial Mechanics v 38 n 4 p 297306 Apr 1986 4Bate RR Muelles D D While J E Fundamentals of astrodynamics New York Dover 1971 5Bliss G A Lectures on the calculus of variations Chicago University of Chicago 1946 6Boas L M Mathematical methods in the physical sciences New York NY John Wiley 1966 7Boltyanskii VG Gamkrelidze RV Pontryagin LS The mathematical theory of optimal processes New York Intersciences 1962 8Bryson Jr A E Ho YC Applied optimal control New York NY John Wiley 1975 9Chern JS Hong Z C Optimal trajectory for vertical ascent to geosynchronous Earth orbit Acta Astronautica v 29 n 7 p 497502 July 1993 10Cornelisse WJ Rocket propulsion and space flight dynamics London England Pitman 1979 64 11Elsgolts L Differential equations and the calculus of variations Moscow MIR 1973 12Fonseca M A F Otimização da altitude de foguete lançado verticalmente numa atmosfera com densidade variável 1977 Dissertação Mestrado em Física Instituto Tecnológico de Aeronáutica São José dos Campos 1977 13Garfinkel B Solution of goddard problem Journal of the Society for Industrial and Applied Mathematics Series A Control v 1 n 3 p 319348 1962 14Gelfand I M Fomin SV Calculus of variations New Jersey Prentice Hall 1963 15Goldstein H Classical mechanics Reading MA AddisonWesley 1980 16Hestenes R M Calculus of variations and optimal control theory New York NY Wiley 1966 17Hildebrand F Methods of applied mathematics New Delhi India Prentice Hall 1968 18Jacoby S L S Kowalik J S Pizzo T J Iterative methods for nonlinear optimization problems New Jersey PrendiceHall 1972 19Kelley J H Uzzell B R McKay S S Rocket trajectory optimization by a secondorder numerical technique AIAA Journal v 7 n 5 p 879884 1969 20Kelley J H Tsiotras P Draglaw effect in the goddard problem Acta Astronáutica v 27 n 3 p 481490 May 1991 65 21Kelley J H A second variation test for singular extremals AIAA Journal v 2 n 8 p 13801382 1964 22Kulakowsi J L Stancil Rocket boots trajectories for maximum burn out velocity Texas ConvairFort Worth Sd 23Lawden F D Rocket trajectory optimization Journal of Guidance Control and Dynamics v 14 n 4 p 705711 JulyAug 1991 24Lawden F D Optimal trajectories for space navigation London Butterworth 1963 25Lawden F D Optimal rocket trajectories Jet Propulsion v 27 n 12 p 12631264 1957 26Lawden F D Necessary conditions for optimal trajectories New Zealand University of Canterbury 1958 27Lawden H G Trajectory optimization for a rocket with a generalized thrust characteristic Astronautica Acta v 10 n 56 p 279281 1964 28Leitmann G Optimization techniques New York NY Academic Press 1962 29Marec P J Optimal space trajectories Amsterdam Netherlands Elsevier 1979 30McDanell J P Powers W F Necessary conditions for joining optimal singular and nonsingular subarcs Siam Journal of Control v 9 n 2 p 161173 May 1971 31McDanell J P Powers W F New jacobitype necessary and sufficient conditions for singular optimization problems AIAA Journal v 8 n 8 p 14161420 1970 66 32Meirovitch L Methods of analytical dynamics New York NY McGraw Hill 1970 33Miele A Flight mechanics and variation problems of a linear type Journal of the Aerospace Sciences v 25 n 9 p 581590 1958 34Miele A Mathematical programming and calculus of variations Houston Texas Rice University 19691970 Manuscripts of courses 35Nash S G Sofer A Linear and nonlinear programming New York NY McGrawHill 1996 36National Aeronautics and Space Administration NASA USStandard Atmosphere Supplements Washington DC NASA 19661976 37Pierson L B Ong S Y Minimumfuel aircraft transition trajectories Advances in Space Research v 11 n 12 p 1927 1991 38Ross I M An analysis of firstorder singular thrustarcs in rocket trajectory optimization Acta Astronautica v 39 n 6 p 417422 Sept 1996 39Shen Y C Formulation of the optimum rocket trajectory ARS Journal v 31 n 1 p 890 1961 40Sielawa J T Apontamentos de aulas sobre cálculo variacional São José dos Campos Instituto Tecnológico de Aeronáutica 1987 Anotações de aulas 41Sielawa J T Baêta Segundo J A Mecânica analítica São José dos Campos ITA 1987 Apostila 42Sielawa J T Luz M S Mecânica Sl Os autores 2000 67 43Singe J L Griffith B A Principles of mechanics New York NY McGrawHill 1959 44Struble R A Stewart C E Granton J The trajectory of a rocket with thrust Jet Propulsion v 28 n 7 p 472478 1958 45Tsien S H Evans R E Optimum thrust programming for a sounding rocket Journal of the American Rocket Society v 21 n 5 p 99107 1951 46Well K H Tandon S R Rocket ascent trajectory optimization via recursive quadratic programming Journal of the Astronautically Sciences v 30 n 2 p 101116 1982 47William TT Introduction to space dynamicsNew York NY Dover 1986