Soluções analítica e numérica para o movimento vertical de um minifoguete

Analytical and numerical solutions for the vertical motion of a model rocket

Resumos

Neste trabalho é descrito um modelo para a solução analítica da equação diferencial de movimento de um minifoguete que descreve uma trajetória vertical num movimento de subida. As soluções para a velocidade, v(t), e altitude, y(t), são comparadas com uma solução numérica utilizando o método de Euler. Para testar estes modelos e comparar os métodos para a obtenção de v(t) e y(t), foi utilizado o voo experimental do minifoguete Epsilon-8 com motor a propelente sólido, do grupo de foguetemodelismo Carl Sagan da Universidade Federal do Paraná (UFPR). Na solução analítica é considerado um modelo aproximado no qual a massa total do minifoguete é constante durante o tempo de queima do propelente. As previsões numérica e analítica para v(t) e y(t) podem ser obtidas utilizando parâmetros físicos relacionados ao minifoguete: as curvas experimentais para o empuxo, T(t), e a força de arrasto. Foram utilizadas curvas T(t) de quatro motores fabricados pela empresa Boa Vista Modelismo e testados pelo grupo de Carl Sagan, que são do mesmo lote do motor usado no voo experimental do Epsilon-8. Embora o método analítico seja satisfatório, o método numérico retorna valores de velocidade e altitudes mais próximos ao voo experimental do Epsilon-8.

Palavras-chave:
Método analítico; método de Euler; minifoguete; equações de trajetória

In this work, a model for the analytical solution of the differential equation of motion of a model rocket that describes a vertical trajectory in an upward movement is described. The solutions for velocity, v(t), and altitude, y(t), are compared with a numerical solution using Euler’s method. To test these models and compare the methods for obtaining v(t) and y(t), the experimental flight of the Epsilon-8 model rocket, with a solid propellant engine, from the Carl Sagan rocketry group of the Federal University of Paraná (UFPR). In the analytical solution, an approximate model is considered in which the total mass of the mini rocket is constant during the propellant burn time. Numerical and analytical predictions for v(t) and y(t) can be obtained, using physical parameters related to the mini rocket, the drag force and the experimental curves for the thrust, T(t). T(t) curves of four engines manufactured by the Carl Sagan group were used, similar to the one used in the experimental flight of the Epsilon-8. Although the analytical method is satisfactory, the numerical method returns values of speed and altitudes closer to the experimental flight of the Epsilon-8.

Keywords:
Analytical method; Euler’s method; model rocket; trajectory equations

Introdução

A corrida espacial no século XX trouxe muitos avanços científicos relacionados à construção de potentes telescópios e foguetes. Telescópios, como o Hubble, revolucionaram a astronomia tornando possível a observação detalhada de nebulosas e de milhares de galáxias, a descoberta de planetas fora do sistema solar e a detecção de buracos negros [¹, ²]. Estas observações necessitam ser realizadas na escuridão do espaço e a uma distância apropriada do nosso planeta. Para isto, torna-se necessário o transporte de telescópios e sondas espaciais via foguetes. As longas distâncias que os foguetes podem alcançar são devidas à eficiência dos seus motores, que utilizam propelentes do tipo sólido, líquido ou a gás, permitindo a exaustão de gases a velocidades hipersônicas, conferindo empuxos elevados por longos períodos de tempo [³].

No ensino da educação básica e superior, temas relacionados a foguetes e às descobertas espaciais levam fascínio para os estudantes, pois envolvem missões audaciosas, “aventureiras”, permitindo-se desbravar o universo e compreender sua origem. Com relação aos foguetes, atividades como o lançamento, determinação de parâmetros físicos relacionados ao desempenho e estimativa de altitudes, são temas que podem ser abordados em sala de aula. Neste contexto, [⁴] obtiveram experimentalmente o coeficiente de arrasto do ar, utilizando um lançamento oblíquo e queda vertical de um “foguete” de garrafa PET. Queiroz & Cuzinatto [⁵] apresentam um trabalho teórico no qual obtiveram as equações da trajetória para “foguetes” balísticos, sem propulsão, que se movem num movimento oblíquo a baixas velocidades, com coeficiente de arrasto proporcional a velocidade, v, e para foguetes com forças de arrasto proporcionais a v². Nesta mesma linha, [⁶] utilizaram um sensor eletrônico, denominado como altímetro [⁷], que permite a obtenção de dados experimentais para v(t) e y(t). Equações teóricas foram utilizadas para se ajustarem às curvas v(t) e y(t), possibilitando obter parâmetros importantes como a taxa de variação de massa, dm/dt, a velocidade de exaustão, v_e, e o coeficiente de arrasto do ar, α. Em outro trabalho, diversos modelos teóricos foram empregados para predizer as curvas v(t) e y(t) num movimento vertical de subida [⁸].

O minifoguete é um pequeno foguete aerodinâmico lançado no ar por alguma fonte. A mais comum é um motor que carrega propelente do tipo sólido como exemplo o KNSu, abreviação do inglês de Potassium (K) Nitrate (N) and Sucrose (Su) altímetro [⁷]. Os modelos podem variar em tamanho de aproximadamente 20 cm até bem mais de 1 m de comprimento e podem ser montados a partir de kits ou feitos do zero. Além da aplicação didática, eles podem ser utilizados para fotografias aéreas de pequeno formato [⁹]. A vantagem de se utilizar o minifoguete no ensino é que ele se assemelha muito aos foguetes a propelente sólido de grande porte, levando os estudantes para algo que se aproxima da engenharia espacial.

Recentemente, [¹⁰] demonstrou que a construção e o lançamento de minifoguetes a propelente sólido no ensino levam a uma quantificação dos ganhos de aprendizado do aluno e melhoram significativamente suas habilidades de raciocínio quantitativo e científico. Os resultados dos questionários dos alunos mostraram um aumento significativo no interesse e confiança dos alunos em abordar questões científicas, bem como uma compreensão de como as missões planetárias são projetadas e conduzidas.

No trabalho aqui desenvolvido foi utilizado o voo experimental do minifoguete Epsilon-8 realizado pelo Grupo de Foguetemodelismo Carl Sagan da UFPR [¹¹]. Ele é o atual recordista nacional de apogeu para minifoguetes com motores da classe E [¹²], atingindo em 2015 o apogeu de 723 m [¹³]. Dois métodos foram utilizados para confrontar teoria e experimento: (i) um método puramente teórico, obtido da solução analítica da equação diferencial de movimento e (ii) o método numérico de Euler que utiliza um programa computacional para a obtenção de gráficos v(t) e y(t) [¹⁴, ¹⁵]. Conforme será descrito, com os dados do motor do foguete e parâmetros físicos associados ao voo pode-se obter: (i) as curvas v(t) e y(t), (ii) a velocidade máxima e (iii) a máxima altitude atingida pelo minifoguete.

Adicionalmente, quatro motores, semelhantes ao inserido no Epsilon-8 para seu lançamento, foram utilizados. Com eles foi possível obter curvas de empuxo em função do tempo, T(t). Os dados das curvas T(t) puderam ser inseridos diretamente na equação diferencial de movimento para a obtenção numérica das curvas v(t) e y(t). Nas simulações teóricas de voo, utilizou-se um empuxo médio constante, para a obtenção de v(t) e y(t).

Este trabalho serve como um guia para que estudantes e professores possam simular o voo de minifoguetes como uma atividade teórico-experimental.

Embasamento teórico

Um modelo analítico para o movimento vertical de subida do minifoguete

A Equação diferencial para o voo de um minifoguete, que se move com aceleração gravitacional constante, g, e força de arrasto, F_r, pode ser escrita como [¹⁶]:

(1)$$m\frac{dv}{dt}=-v_e\frac{dm}{dt}-mg+F_r.$$

Nesta equação, m = m(t) é a massa remanescente do minifoguete. Ela varia a uma taxa dm/dt, no qual os gases são expelidos do motor a uma velocidade v_e. O membro esquerdo na Equação (1) é interpretado como o somatório das forças que atuam no minifoguete. Já o primeiro termo no membro direito desta equação é definido como o empuxo, T. Tendo em vista que dm/dt é negativo, o empuxo é positivo na Equação (1). A força de arrasto quadrática do ar, depende da velocidade, v, de um objeto, sendo escrita como:

(2)$$F_r=-\alpha v^2,$$

no qual α é uma constante positiva que especifica o poder da força de arrasto [¹⁷].

Na prática, T é função do tempo e um gráfico T(t) pode ser obtido experimentalmente a partir de testes estáticos do motor do foguete [³]. Entretanto, para algumas estimativas teóricas para a solução da Equação (1), pode-se utilizar T(t) = constante, o que é válido para um caso no regime estacionário do motor quando os regimes transientes (início e término do funcionamento do motor) são breves [¹⁸]. Em uma simplificação, considera-se também na Equação (1) que dm/dt é constante.

Em um modelo “simples” para o movimento de subida na vertical de um minifoguete, considera-se um valor médio para a massa do foguete (m_o, que é tratada como constante, mesmo durante a queima do combustível [¹]. Esta aproximação é válida se o tempo de queima é pequeno quando comparado ao tempo, t_a, para o foguete atingir o apogeu. Neste caso, m_o é dado pela média:

(3)$$m_o=\frac{m_{inic}+m_f}{2}$$

na qual m_inic e m_f são as respectivas massas inicial e final do minifoguete. Reagrupando os termos da Equação (1) pode-se escrever:

(4)$$\frac{m_0}{\beta -\alpha v^2}dv=dt,$$

no qual: $\beta =-v_e\frac{dm}{dt}-m_{0}g$. A partir da Equação (4) pode-se obter as Equações para v(t) e y(t) [¹⁹].

Durante o tempo de queima, t_q, do propelente as equações de movimento são dadas por:

(5)$$v_I\left(t\right)=\sqrt{\frac{\beta }{\alpha }}tgh\left(\frac{\sqrt{\alpha \beta }}{m_0}t\right),$$

(6)$$y_I\left(t\right)=\frac{m_0}{\alpha }\ln \left[cosh\left(\frac{\sqrt{\alpha \beta }}{m_0}t\right)\right],0s\le t\le t_q.$$

. O subíndice “I” se refere a todo instante no qual t ≤ t_q. Após t_q, o minifoguete tem um movimento balístico de subida na vertical, nas quais apenas as forças de arrasto e gravitacional existem. O tempo necessário para o minifoguete atingir o apogeu, t_a (Equação (7)), é utilizado nas deduções das equações para v(t) e y(t) após t_q. Esse tempo é calculado como:

(7)$$t_a=\sqrt{\frac{m_o}{\alpha g}}\cdot t{g}^{-1}\left(v_q\sqrt{\frac{\alpha }{m_{o}g}}\right)+t_q.$$

As Equações (8) e (9) que se seguem representam v(t) e y(t), respectivamente, após t_q:

(8)$$v_{II}\left(t\right)=\sqrt{\frac{m_{o}g}{\alpha }}\cdot tg\left[\sqrt{\frac{\alpha g}{m_o}}\left(t_a-t\right)\right]$$

(9)$$y_{II}\left(t\right)=y\left(t_q\right)+\frac{m_o}{\alpha }\ln \left\{\frac{cos\left[\sqrt{\frac{\alpha g}{m_o}}\left(t_a-t\right)\right]}{cos\left[\sqrt{\frac{\alpha g}{m_o}}\left(t_a-t_q\right)\right]}\right\},t_q\le t\le t_a.$$

Utiliza-se o subíndice “II” para se referir a todo instante em que t_q ≤ t ≤ t_a. O estudante pode verificar a continuidade da equação para v(t = t_q). Para isso, basta obter “t_a − t_q” na Equação (7), substituir na Equação (8) e realizar alguns cálculos algébricos. Por outro lado, a observação da continuidade de y(t = t_q) é imediata.

Para t = t_a, é possível estimar o apogeu, y_máx. A Equação(9) toma a forma:

(10)$$y_{máx}=y_I\left(t_q\right)-\frac{m_o}{\alpha }\ln \left[\text{cos}\left(\sqrt{\frac{\alpha g}{m_o}}\left(t_a-t_q\right)\right)\right].$$

Substituindo as Equações (5),(6) e (7) para v(t_q), y(t_q), t_a, respectivamente, obtém-se, após algumas manipulações trigonométricas:

(11)$$y_{máx}=\frac{m_o}{\alpha }\ln \left[cosh\left(\frac{\sqrt{\alpha \beta }}{m_o}{t}_q\right)\cdot \sqrt{1+\frac{\beta }{m_{o}g}tan{h}^2\left(\frac{\sqrt{\alpha \beta }}{m_o}{t}_q\right)}\right],$$

que é função apenas de t_q. Esse tempo pode ser obtido da curva experimental de empuxo, de cada um dos quatro motores utilizados.

Na estimativa do apogeu pela Equação (10), pode-se utilizar um valor de t_a teórico, Equação (7), ou um valor experimental, obtido dos dados de um altímetro. Por outro lado, t_q pode ser obtido da curva T(t) ou com os dados de um altímetro. A vantagem de se utilizar a Equação (11) é que se realiza apenas um cálculo para estimar o apogeu.

Modelos analíticos mais sofisticados consideram que na Equação (1): m = m_i − |dm/dt|.t, com dm/dt constante ou variando exponencialmente no tempo. As soluções para v(t) envolvem funções de Bessel de primeira e segunda ordens, necessitando de uma matemática mais avançada [¹⁶, ²⁰]. Mesmo as previsões destes modelos para v(t) e y(t) sendo melhores, os resultados não são muito divergentes das soluções analíticas deduzidas nessa seção, considerando voos em que t_q é relativamente pequeno quando comparado com t_a [⁵].

Método numérico de Euler

O método numérico de Euler é um método do tipo explícito, utilizado em soluções de equações diferenciais ordinárias de primeira ordem na forma:

(12)$$\frac{dy}{dx}=f\left(x,y\right).$$

Este método utiliza o chamado passo simples que fornece uma solução aproximada no ponto (x_i+1, y_i+1) partindo da solução conhecida no ponto (x_i, y_i) [²¹]. Neste caso, utilizam-se as Equações:

(13)$$x_{i+1}=x_i+h,$$

(14)$$y_{i+1}=y_i+\left({\left.\frac{dy}{dx}\right|}_{x=x_i}\right)h,$$

no qual h é o tamanho do passo para a variável independente e ${\left.\frac{dy}{dx}\right|}_{x=x_i}$ é a inclinação da reta no ponto (x_i, y_i).

Para aplicar o método de Euler e obter uma função para v(t), no movimento do minifoguete, primeiro escreve-se a Equação (1) como:

(15)$$\frac{dv}{dt}=-\frac{1}{m}\left(v_e\frac{dm}{dt}\right)-g-\frac{\alpha v^2}{m},$$

Lembrando que o termo entre parênteses no membro direito dessa equação é o empuxo, T(t). A Equação (15) escrita pelo método de passo simples fica [¹⁴]:

(16)$${\left.\frac{dv}{dt}\right|}_{t=t_i}=\frac{T_i}{m_i}-g-\frac{\alpha v_i{}^2}{m_i}.$$

Utilizando a Equação (14), tem-se para a velocidade:

(17)$$v_{i+1}=v_i+\left(\frac{T_i}{m_i}-g-\frac{\alpha v_i{}^2}{m_i}\right)\Delta t,$$

em que:

(18)$$m_i=m_{inic}-\left|\frac{dm}{dt}\right|t.$$

Nas Equações (17) e (18) as interações em Δt devem ser realizadas até o tempo t_a. Para a altitude, função y(t), tem-se que:

(19)$$y_{i+1}=y_i+v_i.\Delta t$$

Na solução numérica, os valores de α, m_inic, g, dm/dt e T(t) são conhecidos, sendo possível plotar gráficos de v(t) e y(t). Neste trabalho, utilizou-se a programação em Python [²²] para a obtenção de tabelas de dados. Posteriormente utilizou-se o programa Origin [²³] para gerar os gráficos com os dados das tabelas.

Metodologia

As soluções analíticas e numéricas para v(t) e y(t) do movimento de subida foram confrontadas com os respectivos gráficos experimentais de v(t) e y(t). Os valores previamente calculados do coeficiente de arrasto (α), da taxa de variação de massa expelida (dm/dt) e da velocidade de exaustão (v_e), dentre outros, que entram nas equações da trajetória, podem ser utilizados como parâmetros durante as simulações computacionais. Vale lembrar que no modelo analítico aqui apresentado, o empuxo, T(t), do minifoguete foi considerado constante e por isso os parâmetros dm/dt e v_e, serão mantidos constantes na construção dos gráficos analíticos. Por outro lado, na simulação numérica, podem ser utilizados dados experimentais do empuxo e o método numérico de Euler para simular v(t) e y(t). Pode-se também utilizar um empuxo constante como uma aproximação para a simulação dos gráficos, quando não se tem a curva T(t).

Resultados e discussões

Parâmetros teóricos/experimentais relacionados aos motores

Como descrito da seção 1 foi utilizado o voo do minifoguete Epsilon-8. O motor acoplado a este foi do tipo propelente sólido da classe E7 [¹⁷]. Parâmetros de importância sobre as características desse motor podem ser obtidos experimentalmente através de testes estáticos. Nestes testes, o motor é acoplado a um anteparo que ligado a um sistema eletrônico fornece os dados de T(t) quando o motor é acionado [³]. A Figura 1 representa os testes estáticos realizados com 4 motores [²⁴]. Estes motores são semelhantes em suas características físicas (tamanho, peso e quantidade de propelente) quando comparado com o motor empregado no voo do Epsilon-8; são todos do mesmo lote de fabricação [²¹]. As escalas vertical e horizontal são as mesmas para todos os gráficos. A linha vermelha horizontal no motor 1 simboliza o empuxo médio, T_m. Representações semelhantes são válidas para os demais motores. A grandeza t_q está indicada em cada figura por uma seta.

Figura 1
Curvas experimentais de empuxo como função do tempo, T(t), para 4 motores de classificação E7. A linha horizontal em vermelho simboliza o empuxo médio. O tempo de queima dos motores são indicados por t_q em cada gráfico. Fonte: O autor.

Nos Gráficos da Figura 1 a magnitude do impulso total, I_tot, é definida de acordo com a Equação (20) como a integral de T(t), durante o tempo de funcionamento do motor:

(20)$$I_{tot}={\displaystyle {\int }_0^{t_q}T\left(t\right)dt.}$$

A magnitude do empuxo médio, T_m, que o motor exerce sobre o minifoguete é calculado pela Equação:

(21)$$T_m=\frac{I_{tot}}{t_q}.$$

As equações analíticas de v(t) e y(t), descritas na seção ^2.1, consideram um empuxo constante, dado por T(t) = T_m. Assim, a magnitude do impulso total pode ser calculada como:

(22)$$I_{tot}=v_e\frac{dm}{dt}{t}_q.$$

Admitindo que v_e e dm/dt são ambos constantes, pode-se estimar um valor médio para essas grandezas. Medindo a massa do propelente, m_P, obtém-se para dm/dt um valor médio dado por:

(23)$$dm/dt=\frac{m_p}{t_q},$$

tornando possível estimar v_e como [²¹]:

(23)$$v_e=\frac{I_{tot}}{m_p}.$$

O impulso específico, I_sp, obtido da relação:

(24)$$I_{sp}=\frac{I_{tot}}{m_{p}g},$$

é outro parâmetro relevante. Com ele mede-se a eficiência do propelente, já que, quanto maior for o impulso específico, maior será a variação da quantidade de movimento, por unidade de peso do propelente [³].

Na Tabela 1 estão apresentados os parâmetros relevantes obtidos dos testes estáticos dos quatro motores utilizados, cujas curvas T(t) estão apresentadas na Figura 1. Dessas curvas extraem-se os valores de t_q e os valores de I_tot, obtidos através da Equação (20). Para a obtenção de dm/dt de cada motor utilizou-se os valores de m_P, medido em uma balança, e o tempo t_q. Os valores de v_e foram obtidos com o uso das Equação (24).

Tabela 1
Parâmetros obtidos dos quatro motores similares ao utilizado no voo do Epsilon-8.

Motores	mp (g)	tq (s)	Itot (N.s)	dm/dt (g/s)	ve (m/s)
Motor 1	47,12 ± 0,01	6,51 ± 0,01	36,60 ± 0,01	7,23 ± 0,01	776,7 ± 0,4
Motor 2	42,61 ± 0,01	4,97 ± 0,01	33,53 ± 0,01	8,58 ± 0,02	786,9 ± 0,4
Motor 3	42,47 ± 0,01	5,02 ± 0,01	34,62 ± 0,01	8,47 ± 0,02	815,2 ± 0,4
Motor 4	42,37 ± 0,01	4,93 ± 0,01	34,82 ± 0,01	8,59 ± 0,02	821,8 ± 0,4

Observa-se nessa tabela que os valores de dm/dt e v_e apresentam valores próximos. O motor 1 por queimar propelente mais lentamente, apresenta menor valor de v_e. Adicionalmente, esse motor por apresentar maior valor de m_P possui menor valor de v_e. Observa-se que os valores de v_e compreendem a região de velocidades supersônicas.

Estimativas teóricas/experimentais extraídas do motor e do voo do Epsilon-8

Na Tabela 2 são apresentados os parâmetros medidos/ estimados do motor do Epsilon-8. Nessa tabela, m_inic é a massa inicial do minifoguete, “motor + minifoguete”; m_f é a massa final; e m_P é a massa do propelente. Com os valores de m_inic e m_f, utiliza-se a Equação (3), para se obter m_o. Esse valor será utilizado na solução analítica para estimativas teóricas do voo do Epsilon-8. Para a estimativa de v_e, com a Equação (24), foi utilizado uma média obtida a partir dos quatro valores de I_tot, apresentados na Tabela 1. Isso foi necessário pois não foi possível obter a curva experimental T(t) para o motor acoplado ao Epsilon-8. O valor médio calculado para o impulso total foi: $I_{tot}^{méd}=\left(35\pm 1\right)N.s$. Adicionalmente, com os valores de m_P e t_q foi possível estimar, utilizando a Equação (23), o parâmetro dm/dt. A grandeza t_q foi obtida com o altímetro acoplado ao Epsilon-8, durante seu voo, conforme se discute na próxima seção.

Tabela 2
Alguns parâmetros experimentais de voo e do motor acoplado ao Epsilon-8.

minic (g)	mf (g)	mp (g)	mo (g)	tq (s)	ve (m/s)	dm/dt (g/s)
131,20 ± 0,01	86,50 ± 0,01	44,70 ± 0,01	108,85 ± 0,02	5,35 ± 0,01	780,8 ± 0,1	8,34 ± 0,01

Na Figura 2 estão apresentadas as respectivas curvas experimentais para v(t) (círculos sólidos) e y(t) (círculos abertos), obtidas do altímetro acoplado ao Epsilon-8, em seu lançamento em 2015. As grandezas experimentais, relevantes do voo, estão indicadas na própria figura. A linha tracejada em vermelho indica o instante t_q. Nesse instante a velocidade passa por um máximo, v_máx. A curva v(t) apresenta uma cúspide pois antes de t_q a resultante das forças aponta no mesmo sentido do movimento; após t_q a resultante das forças apontam no sentido oposto ao movimento, pois a força propulsora, o empuxo, vai a zero rapidamente. Outras grandezas experimentais relevantes são: a altitude no instante em que o propelente é queimado, y(t_q), o apogeu, y_máx, o tempo para que o mesmo seja atingido, t_a.

Figura 2
Curvas experimentais de v(t) e y(t) para o voo real do Epsilon-8.

Comparação dos gráficos obtidos das soluções analítica e numérica

Nesta seção comparam-se os gráficos obtidos via solução analítica e numérica discutidos nas seções ^2.1 e ^2.2, respectivamente. Os gráficos serão comparados entre si e com o voo experimental do Epsilon-8. Vale lembrar que alguns resultados da solução analítica já foram realizados por [⁸] em um trabalho recente. É de importância apresentá-los aqui para que se possa realizar uma comparação da solução analítica com a numérica. O autor mostrou que, para o voo do Epsilon-8, o valor do coeficiente de arrasto α da Equação (2), obtido através de ajustes das curvas teóricas para v(t) e y(t) às respectivas curvas experimentais, foi de α = 0, 000306 kg/m. Com o impulso médio $I_{tot}^{méd}=\left(35\pm 1\right)N.s$, estimou-se por meio da Equação (21) um empuxo médio T_m = (6, 52 ± 0, 02) N. Os valores de α e T_m, juntamente com aqueles apresentados na Tabela 2, foram utilizados para a construção dos gráficos obtidos via soluções analíticas e numéricas. Na Figura 3(a) e (b) estão apresentados os respectivos gráficos para v(t) e y(t), obtidos da solução analítica (marrom) e numérica (azul); os pontos experimentais estão simbolizados por “o”.

Figura 3
Comparação de gráficos de (a) v(t) e (b) y(t) obtidos analiticamente e numericamente com os respectivos gráficos experimentais obtidos via altímetro.

Como exemplo para o estudante, no Apêndice I foi deixado um script da programação em Phyton para v(t), no caso em que o empuxo é variável e dado pela curva T(t) do motor 1, Figura 1. Foi utilizado o método numérico de Euler para a obtenção de v(t). Para y(t) utilizou-se integração a partir dos dados de v(t).

Para a construção dos gráficos obtidos via solução analítica, são utilizadas as Equações (5) e (6) para v_I(t) e y_I(t), respectivamente, e as Equações (8) e (9) para v_II(t) e y_II(t), respectivamente. Adicionalmente, para a construção dos gráficos via método numérico, são usadas as Equações (17) e (19) descritas na seção 2.2 e o tempo de apogeu t_a = 11, 70 s, obtido via altímetro durante o voo do Epsilon-8. O passo utilizado no tempo foi de Δt = 0, 005 s. Observa-se no gráfico v(t) que a curva obtida da solução numérica se encontra abaixo da curva obtida analiticamente, se ajustando melhor aos gráficos experimentais, principalmente nos instantes em que t ≥ t_q. Este melhor ajuste também é observado na curva y(t) no qual a previsão para a altitude máxima é mais satisfatória do que o método analítico. Na Tabela 3 estão apresentadas grandezas físicas de importância obtidas das soluções analítica e numérica.

Tabela 3
Grandezas teóricas de importância estimadas da solução analítica e numérica.

Método	vmáx (m/s)	yI (tq)(m)	ymáx (m)	ymáx * (m)	ta (s)
Analítico	129	474	785	785	12,20
Numérico	130	452	730	–	11,70

• Nota: _máx* foi o valor obtido com a Equação (11)

Para os dois métodos utilizados o valor obtido para v_máx está em excelente concordância com seu valor experimental, não apresentando diferença percentual significante. Os demais valores obtidos analiticamente também estão em boa concordância e com pequenas diferenças percentuais relativas aos valores experimentais: para as grandezas y_I(t_q) e y_máx elas são, respectivamente 5,8% e 8,6%; para o tempo t_a essa diferença é de 4,2%. Esses resultados demonstram que as equações teóricas apresentadas na seção ^2.1 podem ser utilizadas no ensino do movimento vertical de um minifoguete, obtendo resultados satisfatórios. Por outro lado, as previsões obtidas com o método numérico são ainda mais satisfatórias; as diferenças percentuais relativas para y_I(t_q) e y_máx são ambas da ordem de 1%, relativamente menores que na previsão analítica. Conclui-se que embora a solução analítica seja satisfatória para a previsão de velocidades e altitudes, a solução numérica retorna valores mais próximos do voo real do Epsilon-8.

A aceleração do minifoguete

Para uma análise mais detalhada da curva v(t), está apresentado na Figura 4 o gráfico da aceleração em função do tempo, a(t). Ela pode ser obtida com auxílio de ferramentas computacionais. No presente trabalho, utilizou-se a ferramenta “differentiate” do programa Origin. A curva representada pelo símbolo “o” na Figura 4 é a derivada da curva v(t) experimental, enquanto a curva contínua na cor verde representa a derivada da curva v(t) teórica (veja a Figura (1)). A reta em vermelho possui apenas uma representação simbólica.

Figura 4
Curva a(t) obtida da derivada da curva v(t).

Observa-se que a curva teórica se ajusta bem aos pontos experimentais, passando aproximadamente pela região média destes. A aceleração experimental é positiva para t ≤ t_q. Nos instantes iniciais a velocidade do Epsilon-8 vai de 0 a 128 m/s num intervalo de tempo de 5,35 s, com uma aceleração média, a ≈ 24 m/s², equivalente a ≈ 2, 5 g. Isto ocorre devido aos baixos valores da força de arrasto para velocidades relativamente baixas. Essa aceleração leva o foguete a altitudes elevadas, em intervalos de tempo relativamente curtos.

Com o decorrer do tempo, a aceleração diminui principalmente devido ao aumento da força de arrasto com a velocidade, que é do tipo quadrática. O estudante também pode verificar esse decréscimo na aceleração, observando que a inclinação da reta tangente à curva da v(t) (Figura 3(a)) diminui nos instantes em que t ≤ t_q. Em t ≈ t_q, o minifoguete passa a experimentar uma alta variação na aceleração: ela passa por valores médios da ordem de 15 m/s² para valores da ordem de 50 m/s² num curto intervalo de tempo ≈ 1s (veja a inclinação da reta em vermelho). Essa abrupta mudança na aceleração ocorre devido ao cancelamento da força de empuxo, que coincide com o máximo da força de arrasto (a velocidade do minifoguete é máxima em t ≈ t_q).

É importante que o leitor perceba que em t ≈ t_q a aceleração passa brevemente por valores da ordem de −9, 8 m/s² em instantes t ≈ t_q. Durante uma aula, o professor pode levantar este ponto importante em um comentário apenas ilustrativo: “se houvessem astronautas tripulando o Epsilon-8 estes sentiriam um peso aparente nulo, semelhante a uma flutuação. Essa flutuação é muito breve e a aceleração atinge valores negativos maiores rapidamente. Neste instante os astronautas sentiriam uma elevada força retardadora proveniente dos cintos de segurança”.

Experimentalmente, a aceleração possui valores bem definidos em t ≈ t_q. Entretanto, a curva teórica, obtida analiticamente, não possui valores bem definidos nestes instantes. Isso ocorre devido ao fato de que, teoricamente, a força propulsora se anula instantaneamente. É importante que o estudante perceba que isto reflete na abrupta variação da aceleração experimentada pelo minifoguete, resultado da cúspide na curva v(t) em t = t_q.

Previsões teóricas de voo do Epsilon-8 acoplado a cada motor

Nesta seção realizam-se previsões para o voo vertical do minifoguete, utilizando os métodos analítico e numérico. Ou seja, se cada um dos motores apresentados na seção ^4.1 fossem utilizados no voo do Epsilon-8, como seriam as curvas v(t) e y(t)? Valores máximos de velocidade e o apogeu do minifoguete podem ser analisados com relação ao tempo de queima e do impulso médio ou especifico de cada motor.

Na Figura 5, motor 1–4, estão apresentadas as curvas obtidas na solução analítica (em marrom) e da solução numérica (em azul). As curvas simbolizadas por “o” correspondem ao voo experimental do Epsilon-8. Estas curvas foram repetidas em todos os gráficos de v(t) e y(t), para se ter uma comparação com as curvas obtidas das soluções analítica e numérica. A linha vertical pontilhada é a penas um guia para os olhos, para se observar comparativamente os “tempos de queima”.

Figura 5
Estimativas para um voo teórico do Epsilon-8 acoplado a cada um dos motores. A esquerda gráficos v(t) e a direita os respectivos gráficos y(t).

É importante ressaltar que nos gráficos da solução analítica, a força de empuxo e a massa foram mantidas constantes. Por exemplo, na estimativa analítica de voo para o motor 2 (veja a Tabela 1) acoplado ao Epsilon-8, tem-se que v(t) e y(t). Com isso, o empuxo médio a ser utilizado na simulação vale T_m = 6, 75 N. O valor da massa do Epsilon-8 sem o motor é de 86, 50 g (Tabela 2). Com a massa de propelente do motor 2, m_p = 42, 61 g (Tabela 1), tem-se que o valor a ser considerado para a massa é m_o = 107, 81 g. Os demais valores a serem utilizados como v_e e dm/dt estão apresentados na Tabela 1. O tempo de apogeu para cada voo simulado analiticamente foi calculado de acordo com a Equação (7). Para os demais motores, cálculos semelhantes podem ser utilizados.

Por outro lado, na solução numérica o empuxo tem valores bem definidos em todos os instantes de tempo (veja a Figura (1)). Além disso, a massa varia linearmente com o tempo, de acordo com a Equação (18). Por exemplo, para o Epsilon-8 acoplado ao motor 2, tem-se que a massa inicial é $m_{epsilon}+m_p^{motor2}=129$, 11 g, e termina em m_epsilon = 86, 50 g quando t = t_q. Na simulação numérica a condição de que T(t) = 0 para t ＞ t_q ocorre naturalmente devido às características das curvas T(t) (Figura 1). A simulação foi finalizada quando v(t) = 0, em t = t_a.

São observadas semelhanças em suas formas, em todas as curvas v(t) e y(t). Nas simulações numéricas nota-se, relativamente à solução analítica, que: (i) a velocidade máxima ocorre em instantes antes e com valores de velocidades ligeiramente superiores, (ii) após t_q a curva v(t) da simulação numérica se encontra abaixo da curva analítica, (iii) o tempo de apogeu obtido na simulação numérica é inferior ao observado para a solução analítica, (iv) a altitude máxima prevista na solução numérica é menor do que a observada para a solução analítica.

A observação “(i)” é justificada pelo fato de que o termo entre parênteses na Equação (17) começa a apresentar valores negativos num instante antes do propelente se esgotar. A partir deste instante, a velocidade terá acréscimos negativos, começando a decrescer com o tempo. Devido a isso, o minifoguete também apresentará um apogeu instantes antes do simulado, o que justifica a observação “(iii)”. Devido a observação “(ii)” a curva y(t) numérica se encontra abaixo da analítica e apresenta um valor menor para o apogeu. Para uma melhor comparação das curvas obtidas da simulação numérica, estão apresentadas na Tabela 4 o impulso específico, I_sp, e algumas grandezas numéricas importantes no voo de um minifoguete: velocidade máxima v_máx, o instante em que v_máx ocorre, t_vmáx (para diferenciar de t_q), a altitude neste instante, Y_Imáx, e o apogeu, y_IImáx. Para o cálculo do I_s,p, usou-se g = 9, 81 m/s².

Tabela 4
Parâmetros importantes do motor e do voo teórico do Epsilon-8.

Motores	Isp(s)	tvmáx(s)	vmáx(m/s)	yImáx(m)	yIImáx(m)
Motor 1	79,1	6,05	129	465	765
Motor 2	80,3	4,65	135	400	692
Motor 3	83,1	4,72	137	410	711
Motor 4	83,8	4,67	139	413	710

Observa-se na Tabela 4 que v_máx está diretamente relacionado com I^sp: quanto maior I^sp maior foi o valor observado para v_máx. O I^sp está diretamente relacionado com a força média, por peso de propelente, que será exercida num minifoguete acoplado. Isto confere maior aceleração e portanto maior valor de v_máx. Por outro lado, t_vmáx tem um papel importante para a obtenção de y_Imáx e y_IImáx. Mesmo com I^sp relativamente menor, o motor 1 com maior t_vmáx apresenta maiores altitudes y_Imáx e y_IImáx. Observa-se na linha pontilhada da Figura 5 que o motor 1, por apresentar maior t_vmáx apresentará valores de velocidade superiores aos demais motores em uma faixa de tempo entre t_vmáx de cada motor e t_vmáx do motor 1. Teoricamente, isso faz com que o minifoguete acoplado atinja altitudes maiores.

As grandezas y_Imáx e y_IImáx estão relacionada com t_vmáx e I^sp, com um “peso” maior em t_vmáx, desde que os valores de I^sp não apresentem variações consideráveis. Dessa forma, o motor 1 com maior valor de t_vmáx foi o que apresentou maiores valores de y_Imáx e y_IImáx. O motor 2, com menor valor de t_vmáx, foi o que apresentou menores valores para y_Imáx e y_IImáx. Os motores 3 e 4, por terem valores mais próximos de I^sp e t_vmáx são os que apresentaram valores semelhantes de v_máx, y_Imáx e y_IImáx.

Quando se compara os gráficos de altitude para os motores 3 e 4, simulados numericamente com o voo experimental do Epsilon-8, observa-se uma boa proximidade entre as curvas, indicando que estes motores apresentam características semelhantes ao motor utilizado no voo experimental do Epsilon-8.

Considerações finais

Neste trabalho foi possível aplicar um modelo analítico e outro numérico para predizer a Equação de movimento de um minifoguete a propelente sólido. Em um primeiro momento, utilizou-se uma força de empuxo constante em ambos os modelos para a predição das equações de v(t) e y(t) comparadas às respectivas equações experimentais. Em ambos os modelos, as formas das curvas para v(t) e y(t) se assemelham entre si e ao voo experimental do minifoguete, mostrando ao estudante como as teorias são importantes na previsão de movimentos. Ambos os modelos são satisfatórios, entretanto, o numérico se destaca por melhores proximidades com as curvas experimentais e na predição da altitude máxima do minifoguete, com um erro percentual de apenas 1%. O modelo analítico aqui empregado, por utilizar uma aproximação de massa constante durante o funcionamento do motor, prevê velocidades e altitudes maiores que as obtidas experimentalmente no voo do minifoguete.

Já em um segundo momento as curvas de empuxo de 4 motores foguetes, obtidas experimentalmente, foram utilizadas para se comparar as previsões analíticas e numéricas com o voo do Epsilon-8. Para todas as simulações, considerou-se que no modelo analítico a força de empuxo fosse constante. A vantagem da abordagem numérica é que a massa do minifoguete, bem como a força de empuxo tem valores bem definidos para cada instante de tempo, possibilitando obter valores de velocidade e altitudes mais condizentes com o voo experimental do minifoguete. Observou-se que o impulso específico está relacionado com a máxima velocidade atingida pelo minifoguete enquanto que o tempo de queima está relacionado com o apogeu do minifoguete.

Enfim, foi possível confrontar medidas teóricas com as experimentais e relacionar grandezas físicas com o voo de um minifoguete.

Agradecimentos

Agradecemos à Prograd/UFES e a FAPES por conceder apoio financeiro a este projeto. Agradecemos ao professor Carlos Henrique Marchi, bolsista PQ-2/CNPq, pelos dados concedidos do minifoguete Epsilon-8.

Material suplementar

O seguinte material suplementar está disponível online:

Apêndice I: Código da Simulação em Python.

Referências

• ^[1]
  M. Abate, E. Anandapadmanaban, L. Bao, S. Challani, J. Gaughan, A. Jiang, A. Lingineni, A. Vora, C. Yang, D. Zhao, Correlation Between Simulated, Calculated, and Measured Model Rocket Flight, disponível em:http://ftp.demec.ufpr.br/foguete/bibliografia/Abate_et_al_2014.pdf, acessado em: 19/09/2023.
  » http://ftp.demec.ufpr.br/foguete/bibliografia/Abate_et_al_2014.pdf
• ^[2]
  C.E. Aguiar e G.A. Rubini, Revista Brasileira de Ensino de Física 26, 297 (2004).
• ^[3]
  Altus Metrum, Introducing Micro Peak by Altus Metrum, disponível em:https://altusmetrum.org/MicroPeak/, acessado em 29/03/2023.
  » https://altusmetrum.org/MicroPeak/
• ^[4]
  A.L. Alves, A.L. Paneto, K.A. Littike, S.S. Bento e C.H. Marchi, Revista Brasileira de Ensino de Física 42, e20200390 (2020).
• ^[5]
  A.L. Alves, S.S. Bento e C.H. Marchi, Revista Brasileira de Ensino de Física 43, e20200479 (2021).
• ^[6]
  B.N. Biswas, S.P. Chetterjee, P. Mukherjee e S. Paul, Eletronic Journal of Mathematical Analysis and Applications 1, 294 (2013).
• ^[7]
  Boa Vista Modelismo (BVM), disponível em:http://www.boavistamodelismo.com.br/, acessado em 30/04/2023.
  » http://www.boavistamodelismo.com.br/
• ^[8]
  Grupo de Foguetemodelismo Carl Sagan (GFCS), disponível em:http://ftp.demec.ufpr.br/foguete/TRABALHO_TESTES_ESTATICOS/, acessado em 30/04/2023.
  » http://ftp.demec.ufpr.br/foguete/TRABALHO_TESTES_ESTATICOS/
• ^[9]
  J.O. Chang, F.L.Q. Costa, I.C. Dutra, G.A.F. Nery, L.H. Nery, E.R.R. Rodrigues, A.M.G. Selin, F.T. Zambello, D.J. Silva e R.R. Cuzinatto, Caderno Brasileiro de Ensino de Física 36, 529 (2019).
• ^[10]
  L. Hennemann, Journal of Aeronautical Sciences 7,1 (2016).
• ^[11]
  R.R.S Júnior, Da Guerra Fria à Corrida Espacial: Um Recorte Interdisciplinar para a Educação Básica. Dissertação de Mestrado, Universidade Estadual de Feira de Santana, Feira de Santana, Bahia (2021).
• ^[12]
  N.P. Lang, Open Journal of Statistics 10, 139 (2020).
• ^[13]
  O.B. Loureda: In: Anais do 12º Encontro de Iniciação Científica e Pós-Graduação do ITA – XII ENCITA (São José dos Campos, 2006).
• ^[14]
  C.H. Marchi, Foguete UFPR, disponível em:https://fogueteufpr.blogspot.com/2015/12/723-metros-novo-recorde-de-apogeu-para.html, acessado em 30/04/2023.
  » https://fogueteufpr.blogspot.com/2015/12/723-metros-novo-recorde-de-apogeu-para.html
• ^[15]
  C.H. Marchi, Recordes Brasileiros de Minifoguetes, 19ª edição, disponível em:http://ftp.demec.ufpr.br/foguete/Recordes/2022-09-16_Recordes-BAR-19_resumo.pdf, acessado em 30/04/2023.
  » http://ftp.demec.ufpr.br/foguete/Recordes/2022-09-16_Recordes-BAR-19_resumo.pdf
• ^[16]
  R. Nakka, Richard Nakka’s Experimental Rocketry Web Site, disponível em:https://www.nakka-rocketry.net/articles/nakka_theory_pages.pdf, acessado em 28/03/2023.
  » https://www.nakka-rocketry.net/articles/nakka_theory_pages.pdf
• ^[17]
  NATIONAL ASSOCIATION OF ROCKETRY (NAR), disponível em:https://www.nar.org/standards-and-testing-committee/, acessado em 30/04/2023.
  » https://www.nar.org/standards-and-testing-committee/
• ^[18]
  OriginLAB CORPORATION, Origin 8.1 Getting Started Booklet, disponível em:https://www.originlab.com/pdfs/Origin_8.1_Getting_Started_Booklet.pdf, acessado em 29/03/2020.
  » https://www.originlab.com/pdfs/Origin_8.1_Getting_Started_Booklet.pdf
• ^[19]
  PHYTON, disponível em:https://www.python.org/doc/essays/everybody/, acessado em 29/03/2020.
  » https://www.python.org/doc/essays/everybody/
• ^[20]
  A.M. Queiroz e R.R. Cuzinatto, Revista Brasileira de Iniciação Científica 4, 82 (2017).
• ^[21]
  M.V.F Ribeiro, Metodologia de Projeto e Validação de Motores Foguete a Propelente sólido. Dissertação de Mestrado, Universidade de São Paulo, São Carlos, (2013).
• ^[22]
  H. Rodrigues, M.O. Pinho, D. Portes e A.J. Santiago, European Journal of Physics 29, 527 (2008).
• ^[23]
  H. Rodrigues, M.O. Pinho, D. Portes e A.J. Santiago, European Journal of Physics, 30, 185 (2009).
• ^[24]
  A.M. Santos, Lançamento de Foguetes: uma Análise Introdutória da Física, Matemática e Solução Numérica. Trabalho de Conclusão de Curso, Universidade Estadual da Paraíba, Monteiro (2018).
• ^[25]
  J. Waters, Undergraduate Journal of Mathematical: Modeling one + Two 6, 1 (2014).

Datas de Publicação

• Publicação nesta coleção
  27 Out 2023
• Data do Fascículo
  2023

Histórico

• Recebido
  14 Ago 2023
• Revisado
  27 Ago 2023
• Aceito
  04 Set 2023