REGISTRO DOI: 10.70773/revistatopicos/788061573
RESUMO
Este trabalho aborda a modelagem dinâmica e a simulação computacional de um veículo aéreo não tripulado do tipo quadrotor. O objetivo principal é obter uma representação matemática capaz de descrever o comportamento dinâmico do veículo a partir das forças e torques gerados pelos quatro rotores, sem a aplicação de uma estratégia de controle. Para isso, a modelagem é desenvolvida com base no formalismo de Newton–Euler. A partir dessa formulação, são obtidas as equações não lineares de movimento, relacionando o empuxo total, os torques de atitude, a matriz de rotação entre os referenciais inercial e fixo ao corpo e os parâmetros físicos do sistema. Em seguida, o modelo é implementado computacionalmente em espaço de estados e integrado numericamente pelo método de Runge–Kutta de quarta ordem, permitindo avaliar a evolução temporal das variáveis de posição, velocidade, orientação e velocidade angular. As simulações são realizadas por meio da aplicação de entradas prescritas, como variações no empuxo total e nos torques de rolagem, arfagem e guinada, com o objetivo de analisar a resposta natural do modelo dinâmico. Os resultados obtidos permitem verificar a coerência física das equações implementadas, evidenciando o acoplamento entre inclinação e deslocamento translacional, a influência do empuxo na altitude e o comportamento instável do quadrotor na ausência de controle. Dessa forma, o trabalho fornece uma base consistente para estudos posteriores envolvendo controle, identificação de parâmetros, rejeição a perturbações e validação experimental.
Palavras-chave: Quadrotor; Modelagem Dinâmica; Newton–Euler; Runge–Kutta de Quarta Ordem; Simulação Computacional.
ABSTRACT
This work addresses the dynamic modeling and computational simulation of an unmanned aerial vehicle of the quadrotor type. The main objective is to obtain a mathematical representation capable of describing the dynamic behavior of the vehicle based on the forces and torques gen-erated by its four rotors, without the application of a control strategy. To this end, the model is developed using the Newton–Euler formalism. From this formulation, the nonlinear equations of motion are derived, relating the total thrust, attitude torques, the rotation matrix between the inertial and body-fixed reference frames, and the physical parameters of the system. The model is then computationally implemented in state-space form and numerically integrated using the fourth-order Runge–Kutta method, allowing the time evolution of position, velocity, orienta-tion, and angular velocity variables to be evaluated. The simulations are performed by applying prescribed inputs, such as variations in total thrust and roll, pitch, and yaw torques, in order to analyze the natural response of the dynamic model. The obtained results make it possible to ver-ify the physical consistency of the implemented equations, highlighting the coupling between attitude and translational motion, the influence of thrust on altitude, and the unstable behavior of the quadrotor in the absence of control. Thus, this work provides a consistent basis for fur-ther studies involving control, parameter identification, disturbance rejection, and experimental validation.
Keywords: Quadrotor; Dynamic Modeling; Newton–Euler; Fourth-Order Runge–Kutta; Computational Simulation.
1. INTRODUÇÃO
Os veículos aéreos não tripulados do tipo quadrotor têm se destacado em estudos de ro-bótica aérea por sua capacidade de realizar decolagem e pouso vertical, pairar em determinada altitude e executar manobras em espaços reduzidos. Essas características tornam essa configuração atrativa para aplicações como inspeção, monitoramento, mapeamento, vigilância, busca e resgate, além de atividades de ensino e pesquisa em dinâmica de sistemas (Chovancová et al., 2014; Okasha; Kralev; Islam, 2022). Em comparação com outras aeronaves, o quadrotor apresenta simplicidade mecânica, pois seus movimentos são obtidos pela variação das velocidades angulares dos quatro rotores, sem a necessidade de mecanismos complexos para alteração do passo das hélices (Beard, 2008; Bouabdallah, 2007).
Apesar dessa simplicidade construtiva, a descrição matemática do quadrotor envolve aspectos dinâmicos relevantes. O sistema possui seis graus de liberdade, sendo três associados ao movimento translacional e três à orientação da aeronave, mas conta com apenas quatro entradas principais de atuação. Além disso, apresenta comportamento não linear, acoplamento entre os movimentos de translação e rotação, subatuação e instabilidade natural, principalmente quando analisado sem a presença de uma estrutura de controle (Lima et al., 2014; Lima, 2015). Dessa forma, a obtenção de um modelo dinâmico coerente é uma etapa fundamental para compreender o comportamento físico do veículo e servir como base para estudos posteriores.
Neste trabalho, o foco está na modelagem dinâmica e na simulação computacional de um quadrotor, sem o projeto de controladores. A proposta consiste em representar o veículo como um corpo rígido no espaço, sujeito ao empuxo total produzido pelos rotores, aos torques de rolagem, arfagem e guinada e à ação da gravidade. Para isso, adota-se o formalismo de Newton–Euler, abordagem também utilizada por Sá (Sá, 2012), por permitir relacionar diretamente as forças e torques aplicados ao veículo com suas acelerações translacionais e rotacionais.
A modelagem considera dois sistemas de referência: um referencial inercial, fixo em relação ao solo, e um referencial solidário ao corpo do quadrotor, localizado em seu centro de massa. A orientação entre esses referenciais é descrita por meio dos ângulos de Euler, permitindo obter a matriz de rotação responsável por transformar o empuxo gerado no eixo do corpo para o referencial inercial. Com isso, são obtidas as equações não lineares de movimento que descrevem a evolução das posições, velocidades, ângulos de atitude e velocidades angulares do sistema.
Após a formulação matemática, o modelo é implementado computacionalmente em espaço de estados, de modo que a dinâmica do quadrotor possa ser escrita na forma x˙ = f (x, u). A integração temporal das equações diferenciais é realizada pelo método de Runge–Kutta de quarta ordem, escolhido por apresentar boa precisão numérica e implementação direta para sistemas não lineares. Dessa forma, torna-se possível avaliar a resposta temporal do modelo a partir de diferentes condições iniciais e entradas prescritas.
As simulações realizadas têm como objetivo analisar o comportamento natural do quadrotor sem controle. Para isso, são aplicadas variações no empuxo total e nos torques de atitude, observando-se seus efeitos sobre a altitude, os ângulos de rolagem, arfagem e guinada, bem como sobre os deslocamentos nos eixos x, y e z. Essa abordagem permite verificar a coerência física do modelo implementado, evidenciando como pequenas inclinações podem gerar deslocamentos translacionais e como a ausência de realimentação impede a estabilização automática do sistema.
Assim, este trabalho contribui com uma base de modelagem e simulação para o estudo de quadrotores. A estrutura desenvolvida permite compreender a dinâmica não linear do veículo, avaliar a influência das entradas de força e torque e fornecer suporte para trabalhos futuros envolvendo controle, identificação de parâmetros, inclusão da dinâmica dos motores, perturbações externas e validação experimental.
2. CONCEITOS PRELIMINARES
Os veículos aéreos não tripulados, conhecidos como VANTs, são aeronaves capazes de operar sem piloto embarcado, podendo ser controladas remotamente ou por sistemas computacionais. Entre as diferentes configurações existentes, o quadrotor se destaca por sua simplicidade mecânica, capacidade de decolagem e pouso vertical, possibilidade de voo pairado e boa ma-nobrabilidade em espaços reduzidos (Chovancová et al., 2014; Okasha; Kralev; Islam, 2022). Essas características tornam esse tipo de veículo uma plataforma recorrente em estudos de ro-bótica aérea, modelagem dinâmica, simulação computacional e desenvolvimento de estratégias de navegação.
O quadrotor é formado por quatro conjuntos motor-hélice distribuídos de forma simétrica em relação ao centro de massa da estrutura. Dois rotores giram em um sentido e os outros dois no sentido oposto, de modo a reduzir o torque resultante sobre o corpo da aeronave. Diferentemente de helicópteros convencionais, o quadrotor não utiliza mecanismos complexos para alteração do passo das hélices; seus movimentos são produzidos pela variação das velocidades angulares dos rotores (Beard, 2008; Bouabdallah, 2007).
Do ponto de vista dinâmico, o quadrotor é um sistema de seis graus de liberdade, sendo três associados à posição no espaço e três relacionados à orientação da aeronave. As coordenadas translacionais são representadas por x, y e z, enquanto a orientação é descrita pelos ângulos de Euler: rolagem ϕ, arfagem θ e guinada ψ, também conhecidos como roll, pitch e yaw. Apesar de possuir seis graus de liberdade, o veículo conta com quatro entradas principais de atuação, associadas ao empuxo total e aos torques em torno dos três eixos do corpo. Por isso, trata-se de um sistema subatuado, não linear e com acoplamento entre movimentos translacionais e rotacionais (Lima et al., 2014; Lima, 2015).
A descrição do movimento exige a definição de dois sistemas de referência. O primeiro é o referencial inercial, considerado fixo em relação ao solo, no qual são descritas as posições e velocidades translacionais do veículo. O segundo é o referencial fixo ao corpo, com origem no centro de massa do quadrotor, no qual são usualmente representadas as forças, torques e velocidades angulares. A relação entre esses referenciais é feita por meio de matrizes de rotação, construídas a partir dos ângulos de Euler (Beard, 2008; Sá, 2012).
Os movimentos principais do quadrotor são obtidos pela combinação das forças geradas pelos rotores. O deslocamento vertical ocorre quando o empuxo total é alterado simultaneamente nos quatro rotores. A rolagem é produzida pela diferença de empuxo entre rotores laterais, enquanto a arfagem ocorre pela diferença entre rotores dianteiro e traseiro. Já a guinada é gerada pela diferença entre os torques de reação dos pares de hélices que giram em sentidos opostos. Dessa forma, as ações dos rotores podem ser agrupadas em quatro entradas equivalentes: o empuxo total U1 e os torques U2, U3 e U4, associados respectivamente aos movimentos de rolagem, arfagem e guinada (Beard, 2008; Bouabdallah, 2007).
A modelagem dinâmica de um quadrotor pode ser desenvolvida por diferentes abordagens, sendo comuns as formulações de Newton–Euler e Euler–Lagrange. A formulação de Newton–Euler parte do equilíbrio de forças e momentos aplicados ao corpo rígido, enquanto a formulação de Euler–Lagrange utiliza uma abordagem energética, baseada nas energias cinética e potencial do sistema. Embora partam de princípios diferentes, ambas podem descrever a dinâmica do quadrotor quando aplicadas de forma consistente (Martini et al., 2024). Neste trabalho, adota-se o formalismo de Newton–Euler, por permitir relacionar diretamente as forças e os torques produzidos pelos rotores com as acelerações translacionais e rotacionais do veículo (Sá, 2012; Beard, 2008).
Na modelagem adotada, o quadrotor é considerado um corpo rígido e simétrico, com centro de massa coincidente com a origem do referencial fixo ao corpo. As forças principais consideradas são o empuxo total, aplicado ao longo do eixo vertical do corpo, e a força peso, aplicada no referencial inercial. Já a dinâmica rotacional é determinada pelos torques de controle, pelos momentos de inércia da estrutura e, quando considerado, pelo efeito giroscópico associado à rotação dos rotores. Essas hipóteses simplificam a formulação sem retirar os principais efeitos necessários para a análise dinâmica do sistema.
A partir das equações de Newton–Euler, obtém-se um conjunto de equações diferenciais não lineares que descreve a evolução temporal das posições, velocidades, ângulos de orientação e velocidades angulares do quadrotor. Essas equações podem ser organizadas em espaço de estados, na forma x˙ = f (x, u), em que x representa o vetor de estados do sistema e u representa o vetor de entradas. Essa representação é adequada para implementação computacional, pois permite calcular a derivada dos estados a cada instante de tempo a partir da condição atual do veículo e dos comandos aplicados.
Como as equações obtidas não possuem, em geral, solução analítica simples para entradas arbitrárias, a simulação computacional é realizada por integração numérica. Neste trabalho, utiliza-se o método de Runge–Kutta de quarta ordem, que estima a evolução do sistema a partir de avaliações intermediárias da função dinâmica dentro de cada passo de integração. Esse método apresenta boa precisão para sistemas não lineares e permite simular o comportamento temporal do quadrotor sem a necessidade de linearizar previamente o modelo.
Dessa forma, os conceitos apresentados neste capítulo servem como base para a formulação matemática e para a implementação computacional do modelo dinâmico do quadrotor. A partir deles, torna-se possível desenvolver as equações de movimento, organizar o sistema em espaço de estados e realizar simulações em malha aberta para analisar a resposta natural da aeronave diante de diferentes entradas de empuxo e torque.
3. METODOLOGIA
A metodologia deste trabalho foi estruturada de forma sequencial, partindo da definição das hipóteses de modelagem até a implementação computacional e análise das simulações do quadrotor. Inicialmente, são estabelecidos os parâmetros físicos utilizados, os sistemas de referência adotados e as simplificações necessárias para representar o veículo como um corpo rígido no espaço. Em seguida, a dinâmica do sistema é desenvolvida com base no formalismo de Newton–Euler, relacionando as forças e torques produzidos pelos rotores com os movimentos translacionais e rotacionais da aeronave.
Após a obtenção das equações dinâmicas, o modelo não linear é organizado em espaço de estados, de modo que a evolução temporal do sistema possa ser representada pela forma x˙ = f (x, u). Nessa etapa, são definidos o vetor de estados, composto pelas posições, velocidades, ângulos de Euler e velocidades angulares, e o vetor de entradas, formado pelo empuxo total e pelos torques de rolagem, arfagem e guinada.
Com o modelo matemático definido, a simulação computacional é realizada por integração numérica utilizando o método de Runge–Kutta de quarta ordem. Essa escolha permite simular diretamente as equações diferenciais não lineares do quadrotor, sem a necessidade de linearização ou aplicação de controladores. Assim, o comportamento do sistema é analisado em malha aberta, a partir de diferentes condições iniciais e entradas prescritas.
Por fim, as respostas obtidas nas simulações são avaliadas observando a evolução das variáveis de posição, velocidade, orientação e velocidade angular. A análise busca verificar a coerência física do modelo implementado, destacando efeitos como a variação de altitude causada pelo empuxo total, os movimentos angulares provocados pelos torques aplicados e o acoplamento entre a inclinação da aeronave e seus deslocamentos translacionais.
3.1. Configuração Física do Quadrotor
Um drone quadrotor pode ser definido como uma aeronave composta por quatro rotores com asas rotativas, distribuídos de maneira simétrica em relação ao seu centro de massa. Na sua forma mais simples e comum, os motores são fixados nas extremidades de sua estrutura, formando uma cruz, de modo que cada rotor contribui para a geração de empuxo e também para a geração dos torques responsáveis pelos movimentos da aeronave. Essa configuração permite que a aeronave realize decolagem e pouso vertical, mantenha-se pairando no ar e execute manobras por meio do controle das velocidades angulares dos motores (Beard, 2008; Bouabdallah, 2007).
A Figura 1 apresenta um esquemático dessa configuração física do quadrotor, destacando a simetria dos quatro rotores em torno do centro de massa. Cada hélice é acionada por um motor independente, por consequência, cada uma é responsável por gerar uma força de sustentação perpendicular ao plano da estrutura. A intensidade dessa força é modelada como proporcional ao quadrado da velocidade angular do respectivo rotor, conforme usualmente adotado na modelagem dinâmica de quadrirotores (Chovancová et al., 2014; Belkheiri et al., 2012; Delgado-Reyes et al., 2024). Desse modo, pode-se afirmar que o controle de um drone com essa configuração é realizado por meio da alteração das velocidades angulares dos motores.
Para reduzir o efeito do torque de reação causado pelas rotações das hélices, os rotores são configurados em pares que giram em sentidos opostos. Dessa forma, enquanto um par gira no sentido horário, o outro par gira no sentido anti-horário. Em condições ideais de voo pairado, esses torques de reação se anulam devido à configuração em pares, permitindo que o drone mantenha sua orientação constante no ar. Esse mesmo efeito é utilizado para realizar o controle da guinada.
A numeração e a definição do sentido de rotação dos rotores são importantes para a construção do modelo matemático, pois determinam como as velocidades angulares dos rotores se relacionam com o empuxo total e com os torques de rolagem, arfagem e guinada. Assim, a configuração física adotada neste trabalho será utilizada como base para a definição das entradas de controle do sistema e para a formulação das equações dinâmicas que descrevem os movimentos translacionais e rotacionais do quadrotor (Sá, 2012).
3.2. Sistemas de Referência e Eixos de Movimento
Para ser possível descrever os movimentos da aeronave, é necessário definir os sistemas de referência utilizados na modelagem. Neste trabalho, são considerados um sistema de coordenadas inercial, ou seja, fixo em relação à Terra, e um sistema de referência móvel, fixo ao corpo da aeronave. O referencial inercial é utilizado para representar a posição da aeronave no espaço, enquanto o referencial fixo ao corpo é utilizado para descrever as grandezas relacionadas à aeronave, como forças, torques e velocidades angulares dos rotores (Beard, 2008; Sá, 2012).
O referencial inercial é representado pelos eixos x, y e z, sendo utilizado para descrever a posição translacional do centro de massa do quadrotor. Dessa forma, as coordenadas x e y estão associadas ao deslocamento horizontal da aeronave, enquanto a coordenada z está relacionada ao deslocamento vertical, ou seja, à altitude. O referencial fixo ao corpo, por sua vez, possui origem no centro de massa da aeronave e acompanha sua orientação durante os movimentos.
A Figura 2 demonstra os referenciais utilizados para a descrição dos movimentos da aeronave. A relação entre os referenciais é obtida por meio de uma matriz de rotação, construída a partir dos ângulos de Euler. Esses ângulos permitem representar a orientação da aeronave no espaço tridimensional e são definidos como rolagem, arfagem e guinada, estando associados às variáveis ϕ, θ e ψ, respectivamente (Lima, 2015; Beard, 2008).
O ângulo de rolagem, ϕ, representa a rotação do quadrotor em torno do eixo x do referencial fixo ao corpo, estando associado à inclinação lateral da aeronave. O ângulo de arfagem, θ, representa a rotação em torno do eixo y do corpo, estando relacionado à inclinação para frente ou para trás. Por fim, o ângulo de guinada, ψ, representa a rotação em torno do eixo z, alterando a orientação horizontal do veículo. Esses três ângulos são fundamentais para descrever a atitude do quadrotor e serão utilizados na formulação das equações dinâmicas e no projeto das malhas de controle de atitude.
Além da posição e da orientação, também é necessário distinguir as velocidades lineares e angulares do sistema. As velocidades lineares descrevem a variação temporal das coordenadas de posição do centro de massa, enquanto as velocidades angulares descrevem a rotação do corpo da aeronave em torno dos eixos do referencial fixo ao corpo. Sendo assim, as variáveis p, q e r representam as componentes da velocidade angular do quadrotor em torno dos eixos e x, y e z do corpo.
A relação entre as velocidades angulares do corpo e as derivadas dos ângulos de Euler é apresentada por (Diebel, 2006) e pode ser escrita conforme a Equação 1.
Observa-se que essa relação depende dos ângulos de rolagem e arfagem, apresentando termos trigonométricos associados a ϕ e θ. Entretanto, para a condição de voo pairado, considera-se que os ângulos de inclinação permanecem próximos de zero. Assim, adotando as aproximações sin ϕ ≈ 0, cos ϕ ≈ 1, tan θ ≈ 0 e sec θ ≈ 1, a matriz de transformação se aproxima da matriz identidade. Dessa forma, nas proximidades do ponto de operação em hover, as derivadas dos ângulos de Euler podem ser aproximadas pelas velocidades angulares do corpo, conforme indicado na Equação 2.
3.3. Matrizes de Rotação
Como o quadrotor possui um referencial fixo ao corpo e um referencial inercial, é necessário estabelecer uma relação matemática entre esses dois sistemas de coordenadas. Essa relação é feita por meio de matrizes de rotação, responsáveis por transformar vetores descritos no referencial do corpo para o referencial inercial, ou no sentido inverso, dependendo da convenção adotada. No caso do quadrotor, essa transformação é fundamental, pois o empuxo gerado pelos rotores atua inicialmente na direção definida pelo corpo da aeronave, mas o movimento translacional deve ser descrito em relação ao referencial inercial (Beard, 2008; Sá, 2012).
A rotação em torno do eixo x, associada ao ângulo de rolagem ϕ, é dada por:
A rotação em torno do eixo y, associada ao ângulo de arfagem θ, é dada por:
A rotação em torno do eixo z, associada ao ângulo de guinada ψ, é dada por:
Considerando a sequência de rotações em guinada, arfagem e rolagem, a matriz de rotação completa é obtida pelo produto:
Essa matriz será utilizada posteriormente na formulação da dinâmica translacional do quadrotor, permitindo expressar, no referencial inercial, a força de empuxo gerada no referencial fixo ao corpo. Além disso, sua definição é essencial para relacionar a atitude da aeronave com os deslocamentos no espaço, especialmente na análise do voo pairado, na linearização do modelo e no projeto das malhas de controle.
3.4. Forças, torques e equações de Newton-Euler
Após a definição dos referenciais e da matriz de rotação, a dinâmica do quadrotor é descrita a partir das forças e torques atuantes sobre a aeronave. Neste trabalho, o quadrotor é considerado como um corpo rígido de massa m, com centro de massa localizado na origem do referencial fixo ao corpo. Além disso, é assumido que a estrutura apresenta simetria em relação aos seus eixos principais, de modo que a matriz de inércia pode ser representada por uma matriz diagonal, conforme indicado na Equação 8.
Cada rotor i gera uma força de empuxo Fi proporcional ao quadrado de sua velocidade angular Ωi. Essa relação pode ser escrita como:
em que b representa a constante de empuxo da hélice. Além da força de empuxo, cada rotor também produz um torque de reação associado ao arrasto aerodinâmico da hélice, dado por:
em que d representa a constante de arrasto ou constante de torque da hélice.
A partir das forças geradas pelos rotores, é definido as entradas de controle do sistema. A entrada U1 representa o empuxo total produzido pelos quatro rotores, sendo responsável pelo movimento vertical da aeronave. As entradas U2, U3 e U4 representam os torques de rolagem, arfagem e guinada. Considerando a numeração dos rotores adotada neste trabalho, essas entradas são definidas pelas Equações 11 a 12.
Nas equações anteriores, l representa a distância entre o centro de massa do quadrotor e cada rotor. A entrada U1 possui unidade de força, enquanto U2, U3 e U4 possuem unidade de torque. Portanto, o vetor de entradas de controle pode ser escrito como:
As Equações de 11 a 12 mostra que as entradas de controle não correspondem diretamente às velocidades individuais dos motores, mas sim a combinações dessas velocidades. Essa relação é importante porque os controladores são projetados para atuar sobre grandezas físicas de interesse, como empuxo total e torques de atitude, enquanto os motores recebem comandos individuais de velocidade angular.
Para se obter a dinâmica translacional e rotacional do quadrotor, aplica-se o formalismo de Newton-Euler para corpos rígidos. Considerando que a origem do referencial fixo ao corpo da aeronave coincide com o centro de massa do veículo, o que simplifica a formulação dinâmica, pois elimina termos associados ao deslocamento entre a origem do referencial do corpo e o centro de massa. Dessa forma, as equações de movimento podem ser escritas em forma matricial conforme a Equação 14 (Beard, 2008; Sá, 2012).
Na Equação 14, F representa a força total atuando sobre o centro de massa, m é a massa do quadrotor, I ∗ 3 é a matriz identidade de ordem três, a ∗ cm é a aceleração do centro de massa, τ representa o torque total atuando em torno do centro de massa, Icm é a matriz de inércia em torno do centro de massa, ω é a velocidade angular do corpo e α é a aceleração angular do corpo.
A força de empuxo gerada pelos rotores atua no referencial fixo ao corpo da aeronave, na direção perpendicular ao plano do quadrotor. Como os rotores produzem uma força resultante na direção do eixo z do corpo, essa força pode ser representada pela matriz FB, dado pela Equação 15.
Como esse vetor está descrito no referencial do corpo, é necessário transformá-lo para o referencial inercial por meio da matriz de rotação R. Considerando vetores coluna, essa transformação é dada por obtém
Substituindo a matriz de rotação completa e o vetor de força no referencial do corpo,
Além do empuxo gerado pelos rotores, atua sobre o quadrotor a força peso. Considerando o eixo z positivo para cima, pode representar a contribuição gravitacional a ser subtraída da força de empuxo por meio do vetor Fg, dado por
Pela segunda lei de Newton, a força resultante sobre o centro de massa é igual ao produto entre a massa e a aceleração. Assim, a aceleração translacional do quadrotor no referencial inercial pode ser escrita como
Substituindo as forças na Equação 19, obtêm as acelerações translacionais nos eixos x,
Para a dinâmica rotacional, é considerado a matriz de velocidades angulares do quadrotor no referencial fixo ao corpo. Essa matriz é composta pelas velocidades angulares em torno dos eixos x, y e z do corpo, sendo representado por
em que p, q e r representam, respectivamente, as velocidades angulares em torno dos eixos x, y e z do referencial fixo ao corpo. Os torques aplicados ao quadrotor são associados às entradas de controle U2, U3 e U4, a matriz de torques de controle pode ser escrito como
Além dos torques de controle produzidos pela diferença de empuxo e pelos torques de reação dos rotores, é considerado o efeito giroscópico associado à rotação das hélices. Esse efeito ocorre porque os rotores possuem momento de inércia próprio e estão girando com velocidades angulares elevadas. Assim, mudanças na orientação do corpo do quadrotor produzem torques adicionais devido ao acoplamento entre a rotação do corpo e a rotação dos rotores (Sá, 2012; Lima, 2015).
Esse efeito pode ser representado por uma matriz de torque giroscópico dos rotores, dado por
em que Jr é o momento de inércia do rotor em torno do seu próprio eixo de rotação e Ωr representa a velocidade angular resultante dos rotores, considerando seus sentidos de rotação. Para a convenção adotada neste trabalho, essa velocidade resultante pode ser escrita como
Dessa forma, o vetor de torques total aplicado ao corpo do quadrotor é obtido pela soma entre os torques de controle e o torque giroscópico dos rotores
A dinâmica rotacional é então obtida a partir da parte rotacional do formalismo de Newton-Euler
em que α é o vetor coluna de acelerações angulares do corpo, dado por
Considerando a matriz de inércia diagonal, definida em 8.
Antes de aplicar a aproximação associada ao voo pairado, a dinâmica rotacional é mantida em termos das velocidades angulares do corpo, dadas por p, q e r. Assim, substituindo os torques na equação rotacional de Newton–Euler, obtém-se
As Equações 30, 31 e 32 representam a dinâmica rotacional sem assumir diretamente que as velocidades angulares do corpo são iguais às derivadas dos ângulos de Euler. Nesse caso, p, q e r são as componentes da velocidade angular no referencial fixo ao corpo do quadrotor. Para obter a evolução dos ângulos ϕ, θ e ψ, deve-se utilizar a relação cinemática entre as velocidades angulares do corpo e as derivadas dos ângulos de Euler.
Nas proximidades do voo pairado, considerando pequenas inclinações e velocidades angulares reduzidas, pode-se adotar a aproximação
Com isso, as acelerações angulares podem ser escritas diretamente em função dos ângulos de Euler:
Portanto, fora da aproximação de voo pairado, a forma mais adequada é manter a dinâmica rotacional em termos de p, q e r. Já nas proximidades do hover, a aproximação p ≈ ϕ˙, q ≈ θ˙ e r ≈ ψ˙ simplifica a representação e permite escrever diretamente as acelerações angulares em função dos ângulos de Euler. Em ambos os casos, os torques U2, U3 e U4 atuam sobre os movimentos de rolagem, arfagem e guinada, enquanto os termos associados a JrΩr representam o efeito giroscópico dos rotores, introduzindo acoplamentos adicionais principalmente nas dinâmicas de rolagem e arfagem.
Os efeitos giroscópicos dos rotores ocorrem devido ao momento angular gerado pela rotação das hélices e dos motores. Como cada rotor possui uma inércia própria e gira com determinada velocidade angular, mudanças na orientação do corpo do quadrotor produzem torques adicionais associados à conservação do momento angular. No caso ideal de voo pairado, parte desses efeitos tende a se cancelar devido aos pares de rotores girando em sentidos opostos. Entretanto, quando há diferenças entre as velocidades dos rotores ou movimentos angulares do corpo, surge uma velocidade residual Ωr, que introduz acoplamentos na dinâmica rotacional, principalmente entre rolagem e arfagem. Por isso, os termos giroscópicos tornam o modelo mais fiel ao comportamento físico do quadrotor, embora possam ser pequenos em condições próximas ao hover.
3.5. Dinâmica dos Motores
Nas equações anteriores, as velocidades angulares dos rotores Ωi foram utilizadas diretamente para calcular o empuxo total e os torques aplicados ao quadrotor. Entretanto, em uma representação mais completa, essas velocidades não devem ser tratadas como entradas instan-tâneas do sistema, pois cada conjunto motor–hélice possui dinâmica própria. Assim, a tensão aplicada ao motor não produz imediatamente uma velocidade angular, mas sim uma aceleração no rotor, determinada pela relação entre torque eletromagnético, torque resistente, inércia e perdas mecânicas.
Seguindo a modelagem apresentada por Sá (Sá, 2012), o motor sem escovas pode ser descrito por uma equação de movimento que relaciona o torque desenvolvido pelo motor, o torque da carga, a inércia equivalente e o amortecimento viscoso. Essa relação é dada por
em que τe representa o torque eletromagnético gerado pelo motor, τl é o torque de carga, Jr é a inércia associada ao rotor, A é o coeficiente de atrito viscoso e Ω é a velocidade angular do motor. O torque eletromagnético é proporcional à corrente elétrica que circula no motor, sendo escrito como
em que kt é a constante de torque do motor e i(t) é a corrente elétrica.
Para incluir o conjunto motor–hélice na simulação, considera-se também a influência da força contraeletromotriz e da resistência interna do motor. Desconsiderando os efeitos elétricos transitórios, como a indutância da armadura, a corrente pode ser aproximada por
em que V (t) é a tensão aplicada ao motor, Km é a constante de força contraeletromotriz e Rm é a resistência interna do motor.
Além disso, o torque resistente associado à hélice é considerado proporcional ao quadrado da velocidade angular do rotor. Dessa forma,
em que d representa a constante de arrasto da hélice. Substituindo as Equações 38, 39 e 40 na Equação 37, e considerando a inércia equivalente do conjunto motor–rotor como Jm + Jr, obtém-se a equação diferencial da velocidade angular do motor:
Na Equação 41, o primeiro termo representa a redução da velocidade angular devido à força contraeletromotriz do motor, o segundo termo representa o efeito resistente associado ao arrasto aerodinâmico da hélice e o terceiro termo representa a contribuição da tensão aplicada ao motor. Assim, a velocidade angular de cada rotor passa a ser obtida por integração numérica, e não imposta diretamente como uma entrada instantânea.
Para simplificar a implementação computacional, a Equação 41 pode ser escrita na forma
com
Dessa forma, o modelo completo do quadrotor passa a possuir também quatro estados associados às velocidades angulares reais dos rotores. A cada passo da simulação, as tensões Vi aplicadas aos motores são usadas para atualizar Ωi pela dinâmica dos motores. Em seguida, essas velocidades são utilizadas no cálculo das forças e torques equivalentes U1, U2, U3 e U4, que alimentam as equações translacionais e rotacionais do corpo rígido. Com isso, a simulação representa de forma mais fiel o comportamento físico do sistema, pois inclui o atraso e as limitações naturais entre o comando elétrico aplicado aos motores e a força efetivamente gerada pelas hélices.
3.6. Simulação Computacional
A simulação computacional foi desenvolvida a partir do modelo dinâmico não linear apresentado nas seções anteriores. O objetivo é avaliar a evolução temporal do quadrotor considerando simultaneamente a dinâmica do corpo rígido e a dinâmica dos motores. Dessa forma, as velocidades angulares dos rotores não são impostas diretamente, mas calculadas a partir das tensões aplicadas aos motores, conforme a equação diferencial apresentada na Equação 42.
Com a inclusão da dinâmica dos motores, as entradas externas da simulação passam a ser as tensões aplicadas aos quatro motores. Assim, o vetor de entrada é definido por
O vetor de estados reúne as variáveis necessárias para representar o movimento translacional, a atitude, as velocidades angulares do corpo e as velocidades reais dos rotores. Neste trabalho, adota-se
Inicialmente, a dinâmica completa pode ser escrita de forma compacta como
em que a função f representa o conjunto de equações diferenciais do sistema. Essa função inclui: a dinâmica dos motores, o cálculo das velocidades angulares dos rotores, a obtenção das entradas equivalentes U1, U2, U3 e U4, o cálculo da velocidade residual Ωr, as equações translacionais, as equações rotacionais e a relação cinemática entre p, q, r e os ângulos de Euler.
Em cada instante de simulação, o procedimento computacional segue uma sequência definida. Primeiro, as tensões Vi são aplicadas ao modelo dos motores, permitindo calcular a variação das velocidades Ωi. Em seguida, essas velocidades são utilizadas nas equações de empuxo e torque para determinar as entradas equivalentes do corpo rígido. Com isso, são calculadas as acelerações translacionais no referencial inercial e as acelerações angulares no referencial fixo ao corpo. Por fim, as derivadas dos ângulos de Euler são obtidas a partir da relação cinemática apresentada na Equação 1.
Como o modelo resultante é não linear e acoplado, a solução temporal é obtida por integração numérica. Neste trabalho, utiliza-se o método de Runge–Kutta de quarta ordem, conhecido como RK4, por apresentar boa precisão para sistemas dinâmicos não lineares e por permitir implementação direta em ambiente computacional.
Para uma equação diferencial genérica dada por
considerando um passo de integração h, o método RK4 calcula quatro estimativas intermediárias da derivada do sistema:
A partir dessas estimativas, o estado no próximo instante é calculado por
No contexto deste trabalho, cada avaliação da função f dentro do RK4 representa uma avaliação completa do modelo do quadrotor. Portanto, em cada chamada de f , o algoritmo recalcula a resposta dos motores, os empuxos, os torques, os termos giroscópicos e as derivadas dos estados. Isso é importante porque o método RK4 não atualiza apenas as posições e ângulos, mas também as velocidades dos motores, que fazem parte do vetor de estados.
As simulações são realizadas em malha aberta. Portanto, não há realimentação nem controlador atuando sobre o sistema. Os sinais V1(t), V2(t), V3(t) e V4(t) são definidos previamente para cada cenário simulado. Essa escolha permite observar a resposta natural do modelo diante de diferentes condições de atuação e verificar como o quadrotor se comporta quando não existe correção automática de erro.
Além dos ensaios individuais de empuxo e atitude, foram definidas trajetórias geométricas simples no espaço, como uma reta tridimensional e uma trajetória circular em altitude aproximadamente constante. Essas trajetórias são utilizadas como referências planejadas para comparação com a resposta obtida pelo modelo. Entretanto, é importante destacar que tais referências não são utilizadas por um controlador em tempo real. Elas servem apenas para gerar comandos prescritos de tensão e para permitir a comparação entre o caminho desejado e o movimento efetivamente realizado pelo quadrotor.
Assim, para cada trajetória simulada, define-se a referência de posição como
enquanto a posição obtida na simulação é dada por
O erro de trajetória é calculado apenas para análise posterior, sendo definido por
Como não há realimentação, esse erro não modifica as tensões aplicadas aos motores. Dessa forma, caso o veículo se afaste da referência, o sistema não realiza nenhuma correção automática. Essa característica é importante para evidenciar a limitação da operação em malha aberta, principalmente em trajetórias no plano horizontal, nas quais pequenos erros de atitude podem gerar desvios acumulados de posição.
Para avaliar o efeito de perturbações externas, também é incluído um cenário com vento. Nesse caso, considera-se uma força externa aplicada ao centro de massa do quadrotor, expressa no referencial inercial por
Com a presença dessa perturbação, a dinâmica do sistema pode ser escrita de forma compacta como
No modelo implementado, o vento atua diretamente na dinâmica translacional, alterando a força resultante sobre o veículo. Assim, a aceleração translacional passa a depender não apenas do empuxo rotacionado para o referencial inercial e da gravidade, mas também da força de perturbação externa. Essa abordagem permite observar como uma rajada ou variação aleatória de vento pode deslocar o quadrotor de sua trajetória planejada quando não existe controle em malha fechada.
A análise dos resultados é feita a partir da evolução temporal das posições x, y e z, das velocidades lineares, dos ângulos ϕ, θ e ψ, das velocidades angulares p, q e r, das velocidades dos motores Ωi e das tensões aplicadas aos motores. Além disso, para os cenários com trajetória, compara-se graficamente a referência planejada com a trajetória realizada pelo modelo, tanto no espaço tridimensional quanto em cada eixo separadamente.
Dessa forma, é possível verificar se o modelo apresenta comportamento coerente com a física esperada: aumento de altitude quando o empuxo supera o peso, inclinação da aeronave quando há diferença de empuxo entre rotores opostos, deslocamento horizontal associado aos ângulos de atitude e resposta gradual dos motores devido à dinâmica própria dos atuadores. Nos cenários com vento, também é possível observar que perturbações externas geram desvios acumulados na trajetória, reforçando a necessidade de uma estratégia de controle para rejeição de distúrbios e correção de erro.
Assim, a simulação computacional proposta permite estudar o quadrotor como um sistema dinâmico completo, no qual a resposta do corpo rígido depende diretamente da resposta dos motores. A utilização de trajetórias de referência simples, como reta e círculo, torna mais clara a análise visual do comportamento do veículo em malha aberta. Já a inclusão de perturbações de vento evidencia que, sem realimentação, o quadrotor não possui capacidade de corrigir automaticamente desvios de posição, atitude ou altitude.
3.7. Parâmetros Utilizados na Simulação
Para a implementação computacional do modelo, foram definidos os parâmetros físicos do quadrotor, os coeficientes associados aos rotores, os parâmetros da dinâmica simplificada dos motores e as configurações numéricas das simulações. A definição desses parâmetros foi feita com base em valores típicos utilizados em modelos de quadrirotores de pequeno porte e em trabalhos de modelagem dinâmica, simulação e controle de VANTs do tipo quadrotor (Beard, 2008; Bouabdallah, 2007; Sá, 2012; Lima, 2015; Silva, 2015; Chovancová et al., 2014).
As relações entre empuxo, torque e velocidade angular dos rotores seguem a formulação clássica adotada na literatura, em que a força produzida por cada rotor é proporcional ao quadrado de sua velocidade angular e os torques de atitude são obtidos por combinações das velocidades dos quatro rotores (Beard, 2008; Bouabdallah, 2007; Sá, 2012; Lima et al., 2014). A inclusão da dinâmica dos motores foi baseada principalmente na formulação apresentada por Sá (Sá, 2012), com apoio de trabalhos que também tratam da identificação ou simulação de motores e atuadores em quadrirotores (Silva, 2015; Iyer; Bansal, 2021). Já as trajetórias geométricas simples utilizadas nos ensaios, como reta e círculo, foram adotadas com o objetivo de avaliar o comportamento do modelo em malha aberta, de forma semelhante a estudos que utilizam trajetórias planejadas para análise dinâmica e controle de quadrirotores (Muñoz, 2012; Lima, 2015; Okasha; Kralev; Islam, 2022).
É importante destacar que os valores utilizados nesta etapa possuem caráter didático e foram escolhidos de forma coerente com plataformas de pequeno porte. Assim, eles permitem avaliar o comportamento dinâmico do sistema e verificar a consistência da implementação computacional. Entretanto, em uma aplicação experimental, esses valores devem ser substituídos por parâmetros medidos diretamente na plataforma física ou identificados a partir de ensaios específicos, como feito em trabalhos de identificação de parâmetros de quadrirotores (Silva, 2015; Chovancová et al., 2014).
A Tabela 1 apresenta os parâmetros associados ao corpo rígido do quadrotor. Esses parâmetros são utilizados diretamente nas equações translacionais e rotacionais obtidas pelo formalismo de Newton–Euler.
Tabela 1. Parâmetros físicos do quadrotor utilizados na simulação. Fonte: Adaptado a partir de valores típicos empregados em modelos de quadrirotores (Beard, 2008; Bouabdallah, 2007; Sá, 2012; Lima, 2015).
Parâmetro | Símbolo | Valor | Unidade |
Massa do quadrotor | m | 0, 468 | kg |
Aceleração da gravidade | g | 9, 81 | m/s2 |
Momento de inércia em x | Ixx | 4, 856 × 10−3 | kg m2 |
Momento de inércia em y | Iyy | 4, 856 × 10−3 | kg m2 |
Momento de inércia em z | Izz | 8, 801 × 10−3 | kg m2 |
Distância do rotor ao centro de massa | l | 0, 225 | m |
Os coeficientes associados aos rotores e hélices são mostrados na Tabela 2. Esses parâmetros são utilizados para converter as velocidades angulares dos rotores em forças de empuxo e torques aerodinâmicos. A adoção de coeficientes b e d segue a formulação usual na modelagem de quadrirotores, também empregada em trabalhos de simulação e controle (Sá, 2012; Lima et al., 2014; Lima, 2015; Iyer; Bansal, 2021).
Tabela 2. Parâmetros dos rotores e hélices utilizados na simulação. Fonte: Adaptado com base na formulação de empuxo e torque de rotores apresentada em (Sá, 2012; Lima, 2015; Beard, 2008).
Parâmetro | Símbolo | Valor | Unidade |
Constante de empuxo | b | 2, 98 × 10−6 | N/(rad/s)2 |
Constante de arrasto/torque | d | 1, 14 × 10−7 | N m/(rad/s)2 |
Inércia do rotor | Jr | 3, 36 × 10−5 | kg m2 |
A dinâmica dos motores foi representada pela equação diferencial simplificada apresentada na Equação 42. Essa representação considera que a velocidade angular do rotor não muda instantaneamente quando uma tensão é aplicada ao motor, mas evolui de acordo com a contribuição da tensão, da força contra eletromotriz e do torque resistente da hélice. Essa consideração é importante para aproximar a simulação do comportamento físico real dos atuadores, conforme discutido por Sá (Sá, 2012) e Silva (Silva, 2015). Os coeficientes adotados são apresentados na Tabela 3.
Tabela 3. Parâmetros da dinâmica simplificada dos motores. Fonte: Parâmetros adotados para simulação com base na estrutura de dinâmica de motores discutida em (Sá, 2012; Silva, 2015; Iyer; Bansal, 2021).
Parâmetro | Símbolo | Valor | Unidade |
Coeficiente quadrático da dinâmica do motor | Am2 | 0, 005 | rad−1 |
Ganho da tensão aplicada | Bm | 900, 0 | rad/(s2V) |
Tensão mínima aplicada | Vmin | 0, 0 | V |
Tensão máxima aplicada | Vmax | 16, 8 | V |
Velocidade angular máxima do rotor | Ωmax | 1300, 0 | rad/s |
A partir dos parâmetros adotados, a velocidade angular necessária para a condição aproximada de voo pairado é dada por
Substituindo os valores utilizados na simulação, obtém-se
A tensão aproximada para manter essa velocidade em regime permanente, considerando a dinâmica simplificada dos motores, é calculada por
Com os parâmetros adotados, resulta em
Esse valor foi utilizado como referência para a geração dos comandos de tensão em torno da condição de sustentação do quadrotor. O peso total do veículo considerado na simulação é
Além dos parâmetros físicos, foram definidos os parâmetros numéricos das simulações. A Tabela 4 apresenta os valores adotados para o passo de integração, tempo final dos ensaios e limites utilizados na geração dos comandos prescritos. O método de integração utilizado foi o Runge–Kutta de quarta ordem, frequentemente empregado em simulações de sistemas dinâmicos não lineares por apresentar boa precisão e implementação direta (Beard, 2008; Iyer; Bansal, 2021; Neto, 2024).
Tabela 4. Parâmetros numéricos utilizados nas simulações.
Parâmetro | Símbolo | Valor | Unidade |
Passo de integração | h | 0, 01 | s |
Tempo final da rota em linha reta | Treta | 14, 0 | s |
Tempo final da trajetória circular | Tcirc | 16, 0 | s |
Inclinação máxima prescrita | ϕmax, θmax | 12, 0 | graus |
Limite inferior de empuxo total | U1,min | 0, 05mg | N |
Limite superior de empuxo total | U1,max | 2, 5mg | N |
Para a análise do comportamento em malha aberta, foram consideradas duas trajetórias geométricas simples. A primeira consiste em uma rota em linha reta no espaço, partindo de (0, 0, 1) m e terminando em (3, 0, 1, 4) m. A segunda consiste em uma trajetória circular horizontal, com raio de 0, 8 m, período de 16 s e altitude constante de 1, 2 m. Essas trajetórias foram utilizadas apenas como referência geométrica para comparação com a resposta do modelo, não havendo realimentação ou controlador de trajetória. A utilização de trajetórias simples permite observar de forma direta o acoplamento entre atitude e deslocamento translacional, além de evidenciar os erros acumulados em malha aberta, como também é discutido em estudos de simulação e navegação de quadrirotores (Muñoz, 2012; Lima, 2015; Okasha; Kralev; Islam, 2022).
Tabela 5. Parâmetros das trajetórias geométricas utilizadas na simulação.
Trajetória | Parâmetro | Valor | Unidade |
Reta 3D | Ponto inicial | (0, 0, 1) | m |
Reta 3D | Ponto final | (3, 0, 1, 4) | m |
Reta 3D | Tempo total | 14, 0 | s |
Círculo | Raio | 0, 8 | m |
Círculo | Altitude | 1, 2 | m |
Círculo | Período | 16, 0 | s |
Para gerar a rota em linha reta de forma suave, foi utilizada uma função polinomial do tipo smoothstep de sétima ordem. Essa escolha evita mudanças abruptas de posição, velocidade e aceleração no início e no fim da trajetória de referência. Já a trajetória circular foi definida por funções senoidais, com posição, velocidade e aceleração compatíveis com o movimento circular desejado.
Por fim, para avaliar a influência de perturbações externas, foi incluído um modelo simplificado de vento aplicado como força externa no referencial inercial. A inclusão de perturbações externas é recorrente em estudos de controle e simulação de quadrirotores, pois permite avaliar a sensibilidade do sistema e a necessidade de estratégias de rejeição de distúrbios (Lima, 2015; Okasha; Kralev; Islam, 2022; Delgado-Reyes et al., 2024). Neste trabalho, o vento foi composto por uma parcela aleatória filtrada e por uma rajada suave aplicada entre 6 s e 10 s. Os principais parâmetros utilizados para essa perturbação são mostrados na Tabela 6.
Tabela 6. Parâmetros do vento utilizado nas simulações com perturbação.
Parâmetro | Símbolo | Valor | Unidade |
Desvio padrão da parcela aleatória | σv | 0, 04 | N |
Constante de filtragem do vento | τv | 0, 7 | s |
Início da rajada | tini | 6, 0 | s |
Fim da rajada | tfim | 10, 0 | s |
Amplitude da rajada em x | Fvx | 0, 18 | N |
Amplitude da rajada em y | Fvy | −0, 08 | N |
Amplitude da rajada em z | Fvz | 0, 03 | N |
Com esses parâmetros, foram realizadas as simulações em malha aberta para os casos sem vento e com vento. Em todos os cenários, as referências geométricas foram utilizadas apenas para gerar comandos prescritos e para comparação posterior com a trajetória simulada. Portanto, os erros de posição calculados não foram realimentados no sistema e não atuaram sobre as tensões aplicadas aos motores.
3.8. Considerações Finais
Neste capítulo, foi apresentada a metodologia adotada para a modelagem e simulação computacional do quadrotor. Inicialmente, foram definidos os referenciais, a configuração física da aeronave e as relações entre as velocidades dos rotores, o empuxo total e os torques responsáveis pelos movimentos de rolagem, arfagem e guinada. Em seguida, o modelo dinâmico foi desenvolvido com base no formalismo de Newton–Euler, considerando o quadrotor como um corpo rígido sujeito à ação do empuxo, da gravidade, dos torques de atitude e dos efeitos giroscópicos dos rotores.
Também foi incluída a dinâmica dos motores, de modo que as velocidades angulares dos rotores deixam de ser impostas diretamente e passam a ser calculadas a partir das tensões aplicadas. Essa inclusão torna a simulação mais próxima do comportamento físico real, pois representa o atraso natural entre o comando elétrico enviado ao motor e a força efetivamente gerada pela hélice. Assim, o modelo passa a integrar duas partes principais: a dinâmica dos atuadores e a dinâmica de corpo rígido do quadrotor.
A implementação computacional foi organizada em espaço de estados, permitindo reunir em uma única estrutura as variáveis de posição, velocidade, orientação, velocidade angular e velocidade dos motores. A solução temporal das equações diferenciais foi realizada pelo método de Runge–Kutta de quarta ordem, aplicado diretamente ao modelo não linear. Dessa forma, cada passo de simulação atualiza simultaneamente os estados dos motores e os estados mecânicos da aeronave.
As simulações propostas são realizadas em malha aberta, sem aplicação de controladores. Portanto, a análise busca observar a resposta natural do sistema diante de entradas prescritas de tensão nos motores. Essa abordagem é adequada ao objetivo deste trabalho, pois permite verificar a coerência física da modelagem, avaliando como variações nas tensões aplicadas aos rotores afetam o empuxo, os torques, os ângulos de atitude, os deslocamentos translacionais e a altitude do quadrotor.
Com isso, a metodologia apresentada fornece a base necessária para a etapa de simulação e análise dos resultados. A partir dela, torna-se possível estudar o comportamento dinâmico do quadrotor de forma mais completa, considerando tanto a resposta dos motores quanto a resposta do corpo rígido da aeronave, sem recorrer a simplificações de controle ou linearização em torno do voo pairado.
4. RESULTADOS
4.1. Rota em Linha Reta
Nesta primeira simulação, foi considerada uma rota de referência em linha reta no espaço tridimensional. A referência foi utilizada apenas como trajetória geométrica desejada, permitindo comparar o caminho planejado com o movimento efetivamente realizado pelo quadrotor.
A Figura 3 apresenta a comparação entre a trajetória de referência, a resposta do modelo sem vento e a resposta do modelo com vento. Observa-se que, no caso sem vento, o quadrotor acompanha qualitativamente o sentido da rota planejada, porém já apresenta um desvio em relação ao ponto final da referência. Esse comportamento é esperado, pois pequenas diferenças entre o movimento desejado e a dinâmica real do sistema não são corrigidas em malha aberta. No caso com vento, o desvio se torna muito mais evidente, principalmente nos eixos x e y, mostrando que a perturbação externa altera a trajetória e provoca afastamento acumulado do caminho planejado.
A evolução temporal das posições é mostrada na Figura 4. Sem vento, a resposta em x permanece próxima da referência durante parte significativa da simulação, mas termina acima do valor desejado. No eixo z, a altitude também acompanha a tendência de subida, com pequeno desvio ao final. Já no eixo y, como a referência é praticamente nula, o deslocamento permanece reduzido no caso sem vento. Com a inclusão do vento, entretanto, a posição em x cresce além da referência, a posição em y passa a apresentar deslocamento negativo significativo e a altitude também se afasta do valor planejado. Esse resultado evidencia que, sem uma malha de controle, a aeronave não possui mecanismo para compensar perturbações externas.
Para quantificar o desvio da trajetória, foi calculado o erro de posição em cada instante de tempo a partir da distância euclidiana entre a posição de referência e a posição simulada:
A partir desse erro, foram avaliados o erro RMS, o erro máximo e o erro final. O erro RMS representa uma medida média quadrática do erro ao longo de toda a simulação, sendo calculado por
A Tabela 7 apresenta os valores obtidos para a rota em linha reta, com e sem vento.
Tabela 7. Erros de trajetória para a rota em linha reta.
Caso | Erro RMS [m] | Erro máximo [m] | Erro final [m] |
Reta 3D sem vento | 0,3041 | 0,8169 | 0,8169 |
Reta 3D com vento | 2,1879 | 5,8096 | 5,8096 |
No caso sem vento, o erro RMS foi de 0, 3041 m, enquanto o erro máximo e o erro final foram de 0, 8169 m. Como o erro máximo coincide com o erro final, observa-se que o desvio cresceu ao longo da trajetória e atingiu seu maior valor ao término da simulação. Mesmo sem perturbações externas, esse resultado mostra que a aplicação de comandos fixos em malha aberta não garante que o quadrotor finalize exatamente na posição desejada.
Com a presença de vento, o erro RMS aumentou para 2, 1879 m, enquanto o erro máximo e o erro final chegaram a 5, 8096 m. Em comparação com o caso sem vento, o erro médio quadrático aumentou aproximadamente sete vezes. Esse crescimento mostra que a perturbação externa provoca um afastamento significativo da trajetória planejada. Além disso, o fato de o maior erro ocorrer no final da simulação indica um comportamento acumulativo: uma vez desviado, o sistema não possui ação corretiva para retornar à referência.
A força de vento aplicada durante a simulação é apresentada na Figura 5. Nota-se a ocorrência de uma rajada principalmente entre 6 s e 10 s, com componente positiva no eixo x, componente negativa no eixo y e uma pequena componente vertical. Essa perturbação é suficiente para alterar a trajetória do quadrotor, pois o sistema em malha aberta não ajusta as tensões dos motores em resposta ao deslocamento causado pelo vento.
Os ângulos máximos de atitude permaneceram pequenos nos dois casos, com valor máximo de arfagem de aproximadamente 0, 7122◦ e valores nulos de rolagem e guinada. Isso indica que, nesta simulação, o aumento do erro está mais associado à ausência de correção de trajetória e à ação da perturbação externa do que a grandes variações angulares. Portanto, o resultado reforça a principal limitação da malha aberta: mesmo com comandos suaves e pequenas inclinações, desvios de posição podem se acumular e se tornar significativos, especialmente quando o quadrotor é submetido a vento.
4.2. Trajetória Circular
Na segunda simulação, foi considerada uma trajetória circular no plano horizontal, com altitude aproximadamente constante. Assim como no caso da reta, a curva de referência foi utilizada apenas como trajetória geométrica desejada, permitindo comparar o comportamento do modelo em malha aberta com o caminho planejado. Como não há realimentação, qualquer desvio acumulado ao longo do movimento não é corrigido durante a simulação.
A Figura 6 apresenta a comparação entre a referência circular, a resposta do modelo sem vento e a resposta com vento. No caso sem vento, observa-se que o quadrotor consegue reproduzir qualitativamente o formato da trajetória, mantendo um caminho próximo ao círculo desejado, embora com pequenas diferenças de amplitude e de posição ao final do percurso. Já no caso com vento, a resposta se afasta de forma significativa da trajetória circular, deixando de preservar a geometria da rota e apresentando um desvio progressivo ao longo do movimento.
A evolução temporal das posições é mostrada na Figura 7. Sem vento, as coordenadas x e y acompanham de forma satisfatória o comportamento senoidal esperado para uma trajetória circular, embora exista um pequeno desvio acumulado ao final da simulação. A coordenada z permanece praticamente constante, como desejado, indicando que a altitude é mantida de forma próxima ao valor de referência nesse cenário. Com a inclusão do vento, entretanto, a trajetória passa a apresentar crescimento excessivo em x, deslocamento negativo acentuado em y e aumento progressivo em z a partir da metade final da simulação. Esse comportamento mostra que o sistema em malha aberta perde rapidamente a capacidade de permanecer próximo da órbita desejada quando submetido a perturbações externas.
Para quantificar o desvio da trajetória, foram novamente utilizados o erro de posição instantâneo, definido na Equação 63, e o erro RMS, definido na Equação 64. A Tabela 8 apresenta os valores obtidos para a trajetória circular, com e sem vento.
Tabela 8. Erros de trajetória para a trajetória circular.
Caso | Erro RMS [m] | Erro máximo [m] | Erro final [m] |
Círculo sem vento | 0,1591 | 0,4162 | 0,4162 |
Círculo com vento | 2,8950 | 7,1459 | 7,1459 |
No caso sem vento, o erro RMS foi de 0, 1591 m, enquanto o erro máximo e o erro final foram de 0, 4162 m. Esses valores são menores do que os obtidos para a trajetória em linha reta sem vento, indicando que, neste cenário específico, a trajetória circular prescrita foi seguida com desvio relativamente pequeno ao longo da simulação. Ainda assim, como o erro máximo coincide com o erro final, observa-se novamente um acúmulo de desvio ao longo do percurso, evidenciando que mesmo na ausência de perturbações o sistema em malha aberta não retorna exatamente ao ponto desejado ao final da trajetória.
Com a presença de vento, o erro RMS aumentou para 2, 8950 m, enquanto o erro máximo e o erro final atingiram 7, 1459 m. Em relação ao caso sem vento, o erro RMS aumentou mais de dezoito vezes, enquanto o erro final se tornou mais de dezessete vezes maior. Esse crescimento expressivo mostra que a trajetória circular é fortemente degradada quando o sistema é submetido a perturbações externas sem qualquer mecanismo de correção. Além disso, como o erro máximo também ocorre no instante final, conclui-se que o desvio não apenas surge, mas continua se acumulando ao longo do movimento.
A força de vento aplicada durante a simulação é apresentada na Figura 8. Assim como no caso da rota em linha reta, observa-se uma rajada predominante entre aproximadamente 6 s e 10 s, com componente positiva no eixo x, componente negativa no eixo y e pequena componente vertical. A atuação dessa perturbação altera progressivamente a órbita do quadrotor, fazendo com que a resposta deixe de descrever uma trajetória aproximadamente fechada e passe a apresentar deriva espacial significativa.
Os ângulos máximos de atitude permaneceram baixos nos dois casos, com valores de aproximadamente 0, 7358◦ para roll, 0, 7448◦ para pitch e 0, 0285◦ para yaw. Assim como observado na rota em linha reta, esses valores indicam que o aumento do erro de trajetória não está associado a grandes inclinações da aeronave, mas sim à ausência de realimentação e à incapacidade do sistema em malha aberta de compensar o vento e corrigir desvios acumulados. Portanto, os resultados da trajetória circular reforçam a mesma conclusão obtida anteriormente: mesmo quando o quadrotor opera com pequenas variações angulares, a falta de controle faz com que perturbações externas provoquem erros de posição significativos e crescentes ao longo do tempo.
5. CONCLUSÕES
Este trabalho apresentou a modelagem dinâmica e a simulação computacional de um VANT do tipo quadrotor, considerando o veículo como um corpo rígido sujeito ao empuxo gerado pelos quatro rotores, aos torques de rolagem, arfagem e guinada, à ação da gravidade e aos efeitos giroscópicos associados à rotação dos rotores. A modelagem foi desenvolvida com base no formalismo de Newton–Euler, permitindo relacionar diretamente as forças e torques aplicados ao veículo com suas acelerações translacionais e rotacionais.
A partir da formulação matemática, foi obtido um modelo dinâmico não linear capaz de representar a evolução temporal das posições, velocidades lineares, ângulos de Euler e velocidades angulares do quadrotor. Além disso, foi incluída a dinâmica dos motores, de modo que as velocidades angulares dos rotores não foram impostas instantaneamente, mas calculadas a partir das tensões aplicadas aos motores. Essa escolha tornou a simulação mais representativa do comportamento físico do sistema, pois considera o atraso natural entre o comando elétrico aplicado e a força efetivamente gerada pelas hélices.
O modelo completo foi organizado em espaço de estados e implementado computacionalmente. A integração temporal das equações diferenciais foi realizada pelo método de Runge–Kutta de quarta ordem, aplicado diretamente ao modelo não linear. Com isso, cada passo de simulação atualizou simultaneamente os estados associados ao corpo rígido e às velocidades dos quatro motores, permitindo analisar o comportamento do quadrotor sem recorrer à linearização em torno do voo pairado.
As simulações foram realizadas em malha aberta, sem realimentação e sem aplicação de controladores. Dessa forma, os sinais de tensão aplicados aos motores foram definidos previamente para cada cenário, e o erro entre a posição desejada e a posição obtida não foi utilizado para corrigir o movimento do veículo. Essa abordagem permitiu observar a resposta natural do modelo e evidenciar a limitação da operação em malha aberta, principalmente em situações nas quais o quadrotor precisa seguir uma trajetória ou rejeitar perturbações externas.
Nos ensaios com trajetória em linha reta, o modelo sem vento apresentou erro RMS de 0, 3041 m e erro final de 0, 8169 m. Já com a inclusão de vento, o erro RMS aumentou para 2, 1879 m e o erro final atingiu 5, 8096 m. Esses resultados mostram que, mesmo quando o sistema acompanha qualitativamente o sentido da rota planejada, pequenos desvios não corrigidos tendem a se acumular ao longo do tempo. A presença de vento intensificou esse efeito, deslocando significativamente o quadrotor em relação à referência.
Na trajetória circular, o comportamento sem vento apresentou erro RMS de 0, 1591 m e erro final de 0, 4162 m, indicando que o modelo conseguiu reproduzir de forma aproximada o formato da rota planejada. Entretanto, com a perturbação de vento, o erro RMS aumentou para 2, 8950 m e o erro final chegou a 7, 1459 m. Esse resultado reforça que a ausência de realimentação impede o sistema de compensar desvios, fazendo com que a trajetória deixe de preservar sua geometria original e passe a apresentar deriva acumulada.
Em ambos os cenários, os ângulos máximos de atitude permaneceram pequenos, o que indica que os erros de trajetória não foram causados por grandes inclinações ou movimentos angulares excessivos. O aumento dos desvios está associado principalmente à própria característica da simulação em malha aberta: uma vez que o quadrotor se afasta da referência, não existe mecanismo de correção para reduzir o erro. Assim, mesmo comandos suaves e pequenas variações angulares podem resultar em erros de posição significativos ao longo do tempo.
Dessa forma, conclui-se que o modelo implementado apresenta comportamento coerente com a física esperada para um quadrotor. A simulação permitiu verificar a influência do empuxo na altitude, o acoplamento entre atitude e deslocamento translacional, a resposta gradual dos motores e a sensibilidade do sistema a perturbações externas. Os resultados também mostram que a simples definição de uma trajetória de referência não é suficiente para garantir seu acompanhamento, sendo necessária uma estratégia de controle em malha fechada para corrigir erros e rejeitar distúrbios.
Portanto, o objetivo principal do trabalho foi atingido, pois foi desenvolvida uma base matemática e computacional capaz de representar e simular a dinâmica não linear de um quadrotor com dinâmica dos motores. Essa base pode ser utilizada em estudos posteriores envolvendo controle de atitude, altitude e posição, identificação de parâmetros, inclusão de modelos mais completos de vento e sensores, além de validação experimental em uma plataforma física.
REFERÊNCIAS BIBLIOGRÁFICAS
BEARD, R. W. Quadrotor Dynamics and Control Rev 0.1. Provo, UT, 2008. Disponível em: https://scholarsarchive.byu.edu/facpub/1325/. Acesso em: 7 ago. 2026.
BELKHEIRI, M.; RABHI, A.; HAJJAJI, A. E.; PéGARD, C. Different linearization control techniques for a quadrotor system. In: 2012 2nd International Conference on Communications, Computing and Control Applications (CCCA). [S.l.]: IEEE, 2012. p. 1–6.
BOUABDALLAH, S. Design and Control of Quadrotors with Application to Autonomous Flying. 155 p. Tese (Tese de Doutorado) — École Polytechnique Fédérale de Lausanne, Lausanne, Switzerland, 2007. Thèse n. 3727. Disponível em: https://infoscience.epfl.ch/entities/publication/4fbc8fc8-1a21-4a08-a24d-98f8302630f8. Acesso em: 7 ago. 2026.
CHOVANCOVá, A.; FICO, T.; CHOVANEC, ; HUBINSKý, P. Mathematical modelling and parameter identification of quadrotor: A survey. Procedia Engineering, v. 96, p. 172–181, 2014. DELGADO-REYES, G.; VALDEZ-MARTíNEZ, J. S.; GUEVARA-LóPEZ, P.; HERNáNDEZ-PéREZ, M. A. Hover flight improvement of a quadrotor unmanned aerial vehicle using PID controllers with an integral effect based on the Riemann–Liouville fractional-order operator: A deterministic approach. Fractal and Fractional, v. 8, n. 11, p. 634, 2024.
DIEBEL, J. Representing Attitude: Euler Angles, Unit Quaternions, and Rotation Vectors. Stanford, CA, 2006. 1–35 p. Technical report.
IYER, A.; BANSAL, H. O. Modelling, simulation, and implementation of PID controller on quadrotors. In: 2021 International Conference on Computer Communication and Informatics (ICCCI). [S.l.]: IEEE, 2021. p. 1–7.
LIMA, G. V. Modelagem dinâmica e controle para navegação de um veículo aéreo não tripulado do tipo quadricóptero. 120 f. Dissertação (Mestrado em Engenharia Elétrica) — Universidade Federal de Uberlândia, Uberlândia, MG, 2015. Disponível em: https://repositorio.ufu.br/handle/123456789/14608. Acesso em: 7 ago. 2026.
LIMA, G. V.; SOUZA, R. M. J. A. d.; MORAIS, A. S. d.; MORAIS, J. S. d. Modelagem dinâmica de um veículo aéreo não tripulado do tipo quadricóptero. In: UNIVERSIDADE FEDERAL DE UBERLâNDIA, Uberlândia, MG. Anais da XII Conferência de Estudos em Engenharia Elétrica (CEEL). 2014. Disponível em: https://www.peteletricaufu.com.br/static/ceel/doc/artigos/artigos2014/ceel2014artigo020r01.pdf. Acesso em: 7 ago. 2026.
MARTINI, S.; VALAVANIS, K. P.; STEFANOVIC, M.; RUTHERFORD, M. J.; RIZZO, A. Correction to the euler lagrange multirotor model with euler angles generalized coordinates. Journal of Intelligent & Robotic Systems, v. 110, p. 17, 2024.
MUñOZ, M. E. P. Modelagem matemática e controle de um quadrimotor. Dissertação (Mestrado em Sistemas Mecatrônicos) — Universidade de Brasília, Brasília, DF, 2012.
NETO, J. B. d. S. Modelagem e análise de dinâmica de voo para VANT quadrirrotor com aplicação de controle PID. 74 f. Monografia (Trabalho de Conclusão de Curso (Graduação em Engenharia de Computação)) — Universidade Federal Rural do Semi-Árido, Pau dos Ferros, RN, 2024.
OKASHA, M.; KRALEV, J.; ISLAM, M. Design and experimental comparison of PID, LQR and MPC stabilizing controllers for parrot mambo mini-drone. Aerospace, v. 9, n. 6, p. 298, 2022.
Sá, R. C. Construção, modelagem dinâmica e controle PID para estabilidade de um veículo aéreo não tripulado do tipo quadrirotor. 67 f. Dissertação (Mestrado em Engenharia de Teleinformática) — Universidade Federal do Ceará, Fortaleza, CE, 2012. Disponível em: http://repositorio.ufc.br/handle/riufc/4097. Acesso em: 7 ago. 2026.
SILVA, M. F. d. Modelagem dinâmica, identificação de parâmetros e controle de um veículo aéreo não tripulado do tipo quadricóptero. 128 f. Dissertação (Mestrado em Sistemas de Comunicação e Automação) — Universidade Federal Rural do Semi-Árido, Mossoró, RN, 2015. Disponível em: https://repositorio.ufersa.edu.br/handle/tede/537. Acesso em: 7 ago. 2026.
1 Universidade Federal de São João del-Rei (UFSJ), São João del-Rei, MG, Brasil. E-mail: syrioigna [clique para visualizar o e-mail]acesse o artigo original para visualizar o e-mail.
2 Universidade Estadual de Campinas (Unicamp), Campinas, SP, Brasil. E-mail: [clique para visualizar o e-mail]acesse o artigo original para visualizar o e-mail.