Lançamento Oblíquo Estocástico: mudanças entre as edições

De Física Computacional
Ir para navegação Ir para pesquisar
 
(34 revisões intermediárias por 3 usuários não estão sendo mostradas)
Linha 17: Linha 17:
</math>
</math>


Para as equações de movimento usamos a segunda lei de Newton nas componentes horizontais  
Para as equações de movimento usamos a segunda lei de Newton nas componentes horizontais e verticais


<math>
<math>
Linha 44: Linha 44:
</math>
</math>


Essas E.D. possuem solução analitica bastando integrar nos limites adequados.
Essas E.D. possuem solução analitica bastando integrar nos limites adequados. A integração numérica por Euler é dirta a partir das relações a cima e nos da uma trajetoria para comparação quando adicionar o termo de ruido:
 
<math>
v_{x,j+1}=v_{x,j}-kv_{x,j}dt
</math>
 
<math>
v_{y,j+1}=v_{y,j}-g dt- kv_{y,j}dt
</math>
 
Utilizando as seguintes condições iniciais <math>x_{0}=0</math> e <math>y_{0}=0</math> podemos atualizar a posição
 
<math>
v_{x}=\frac{dx}{dt}
</math>
 
<math>
x_{j+1}=x_{j}+v_{x,j}dt
</math>
 
<math>
v_{y}=\frac{dy}{dt}
</math>
 
<math>
y_{j+1}=y_{j}+v_{y,j}dt
</math>
 
Para modelar o ruído produzido pela densidade do ar, consideramos dois regimes. No primeiro, assumimos que a densidade é constante ao longo da trajetória, ou seja, <math>\rho=\rho_{0}</math>  e, portanto, trataremos o problema como uma EDE com ruído aditivo, onde <math>B_{(X(t),t)}=\rho_{0}=\beta</math>. Na segunda situação, consideramos que a densidade varia com a altitude do projétil, seguindo a relação <math>\rho=\rho_{0}e^{-\frac{y}{H}}</math> e, nesse caso, abordaremos como uma EDE com ruído multiplicativo, onde <math>B_{(X(t),t)}=\rho_{(X(t))}=\beta e^{-\frac{y}{H}}</math>.


== Equação diferencial estocástica==
== Equação diferencial estocástica==
Linha 51: Linha 79:


<math>
<math>
\frac{dX}{dt} = A(X(t)) + B(X(t))\xi(t).
\frac{dX}{dt} = A(X(t)) + B(X(t))\xi(t),
</math>
</math>


Caso o termo <math>B(X,t)</math> for uma constante, exemplo em que a densidade do ar não muda significativamente, o ruído é dito aditivo. Então, a integração da equação é direta, resultando em:
onde
 
<math>
X = \begin{pmatrix}
x \\ y \\ v_x \\ v_y
\end{pmatrix}, \quad
A = \begin{pmatrix}
v_x \\ v_y \\ -k v_x \\ -g-kv_y
\end{pmatrix}, \quad
B = \begin{pmatrix}
0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 \\
0 & 0 & \beta e^{-\frac{y}{H}} & 0 \\
0 & 0 & 0 & \beta e^{-\frac{y}{H}}
\end{pmatrix}, \quad
\xi = \begin{pmatrix}
0 \\ 0 \\ \xi_x(t) \\ \xi_y(t)
\end{pmatrix}
</math>
 
Caso o termo <math>B(X,t)</math> fosse uma constante, exemplo em que a densidade do ar não muda significativamente, o ruído é dito aditivo. Então, a integração da equação é direta, resultando em:


<math>
<math>
Linha 60: Linha 108:
</math>
</math>


onde <math>\xi(t)dt = dW(t)</math> é o incremento de Wiener, um processo estocástico com distribuição gaussiana de largura <math>\sqrt{t}dt</math>.
onde <math>\xi(t)dt = dW(t)</math> é o incremento de Wiener, um processo estocástico com distribuição gaussiana de largura <math>\sqrt{dt}</math>. Escrevendo a equação em diferenças finitas e substituindo os valores de <math>A(X(t))</math> e de <math>B</math> têm-se:
 
<math>
\Delta x_i = v_{x_i} \Delta t
</math>
 
<math>
\Delta v_{x_i} = -k v_{x_i} \Delta t + \beta \Delta W_{x_i}
</math>
 
<math>
\Delta y_i = v_{y_i} \Delta t
</math>
 
<math>
\Delta v_{y_i} = \left(-g  -k v_{x_i} \right) \Delta t + \beta \Delta W_{x_i}
</math>
 


Entretanto, quando o termo <math>B(X,t)</math> depende do processo estocástico, o ruído é dito multiplicativo e devemos utilizar um cálculo estocástco especial, por exemplo  o cálculo de Itô e de Stratonovich. No cálculo de Stratonovich, as regras de integração, diferenciação, mudança de variáveis, etc se mantém as mesmas do cálculo usual, e a equação é escrita como:
Entretanto, o termo <math>B(X,t)</math> depende do processo estocástico e o ruído é dito multiplicativo, assim, devemos utilizar um cálculo estocástco especial, por exemplo  o cálculo de Itô e de Stratonovich. No cálculo de Stratonovich, as regras de integração, diferenciação, mudança de variáveis etc, se mantém as mesmas do cálculo usual, e a equação é escrita como:


<math>
<math>
Linha 68: Linha 133:
</math>
</math>


com o detalhe que o argumento de <math>B(X(t))</math> é, na realidade, cálculado na média de X entre os limites do intervalo de integração. Ou seja, ao escrever a equação em diferenças finitas, obtém-se:
com o detalhe que o argumento de <math>B(X(t))</math> é, na realidade, calculado na média de X entre os limites do intervalo de integração. Ou seja, ao escrever a equação em diferenças finitas, obtém-se:


<math>
<math>
Linha 74: Linha 139:
</math>
</math>


Agora, é necessário explicitar <math>\Delta X(t)</math> na equação acima, porém, a dependência de <math>B(X(t))</math> não é linear em relação à altura, o que torna inviável esse procedimento. Com isso, a atenção volta-se para o cálculo de Itô, dado por:
Agora, é necessário explicitar o <math>\Delta X(t)</math> na equação acima, porém a dependência de <math>B(X(t))</math> não é linear em relação à altura, o que torna inviável esse procedimento.  
 
 
Com isso, a atenção volta-se para o cálculo de Itô, dado por:


<math>
<math>
dX(t) = \left[A(X(t)) + \frac{1}{2} \frac{\partial B}{\partial X}B(X(t))\right] dt + B(X(t)) \bullet dW(t)
dX(t) = \left[A(X(t)) + \frac{1}{2} \frac{\partial B}{\partial X}B(X(t))\right] dt + B(X(t)) \bullet dW(t)
</math>
</math>
No caso em que <math>B</math> não depende de <math>X(t)</math>, a derivada que aparece na equação é zero e é recuperado o caso de ruído aditivo já mencionado. Aqui, vale ressaltar que <math>B(X(t))</math> é uma matriz e está sendo derivada com relação ao vetor <math>X(t)</math>, então o cálculo para as componentes de <math>X_i(t)</math> são:
<math>
dX_i(t) = A^{(I)}_i(X(t))dt + \sum_{\alpha=1}^m B_{i,\alpha}(X(t)) \bullet dW_\alpha(t),
</math>
onde
<math>
A^{(I)}_i(X(t)) = A_i(X(t)) + \frac{1}{2} \sum_{j=1}^n \sum_{\alpha=1}^m \frac{\partial B_{i,\alpha}}{\partial X_j} B_{j,\alpha}(X(t)).
</math>
Após a substituição dos valores de <math>A_i(X(t))</math> e <math>B_i(X(t))</math>, obtêm-se as seguintes expressões:
<math>
\Delta x_i = v_{x_i} \Delta t
</math>
<math>
\Delta v_{x_i} = \left[ - k v_{x_i} -\frac{\beta^2}{2H} e^{-\frac{2}{H} y_i} \right] \Delta t + \beta e^{- \frac{y(t)}{H}} \Delta W_{x_i}
</math>
<math>
\Delta y_i = v_{y_i} \Delta t
</math>
<math>
\Delta v_{y_i} = \left[ -g - k v_{x_i} -\frac{\beta^2}{2H} e^{-\frac{2}{H} y_i} \right] \Delta t + \beta e^{- \frac{y_i}{H}} \Delta W_{y_i}
</math>
É importante ressaltar que os incrementos de Wiener que aparecem nas equações para velocidade são independentes entre si, isto é, deverá ser gerado duas variáveis aleatórias com distribuição gaussiana de largura <math>\sqrt{\Delta t}</math>. Além disso, um detalhe interessante é que a equação para <math>\Delta v_x</math> é semelhante à equação para <math>\Delta v_y</math> em que <math>g = 0</math>.
== Resultados ==
Aqui temos alguns gráficos comparativos, nesses dois primeiros comparamos os dois modelos para a mesma distribuição de ruídos.
[[File:multeuler.gif|400px|Figura 1 - Ruído multiplicativo]][[File:Aditivo_euler.gif|400px|Figura 2 - Ruído aditivo]]
Nas figuras abaixo podemos ver as trajetórias de certo objeto para ruído aditivo e ruído multiplicativo para diferentes distribuições. 
[[File:Ruido_multiplicativo2.gif|400px|Figura 1 - Ruído multiplicativo]][[File:Ruido_aditivo2.gif|400px|Figura 2 - Ruído aditivo]]
Esperávamos que as média para ruído aditivo e multiplicativos fossem iguais, dado que o ruído não interfere no cálculo da média pois a probabilidade de puxar em qualquer direção é a mesma, porém como pode ser visto no gráfico abaixo há divergência para o ruído multiplicativo.
[[File:N=1000_SCC2.png|400px|Figura 3 - Médias para N = 1000]]
Tanto a divergência da médias como o aumento da dispersão da trajetórias podem ser explicados pelo modelo que estamos usando na EDE-ito <math>(\beta e^{- \frac{y_i}{H}})</math>, nesse modelo conforme temos y variando para valores negativos a exponencial explode consequentemente temos aumento na intensidade do ruído, esse aumento pode ser observado no gráfico abaixo.
[[File:Exp.png|400px|Figura 7 - exp]]
Temos que para poucos passos o ganho do ruído em qual direção causara grande divergência, ou seja, para termos vermos a convergência para média do ruído multiplicativo precisamos incrementar o número de trajetórias, a figura a seguir foi feita para um aumento de dez vezes no número de passos comparado ao da figura anterior. 
[[File:N=10000_SCC2.png|400px|Figura 4 - Médias para N = 10000]]
Nos gráficos abaixo temos as trajetórias para o ruído multiplicativo, as médias para os dois métodos e a solução por euler levando em conta o solo, para as trajetória é visível a dispersão em torno da solução de euler que acontece por causa do ruído, porém como comentado anteriormente essa dispersão não influência na trajetória média, da mesma forma que as duas primeiras figuras conforme a diminuição das velocidades pelo termo dissipativo temos maior influência do ruído no trajeto que pode até mesmo ocasionar em um retorno  do objeto na direção da posição inicial porém a média continua somente para a direita.     
[[File:N=1000_CCC2.png|400px|Figura 5 - Médias para N = 1000]]
Agora para um maior valor de ruído podemos observar o gráfico a seguir, onde notamos que a média, na região de contato com o solo, fica acima da solução de euler, isso ocorre pela condição de contorno que estamos utilizando, quando a posição da partícula chega em zero invertemos o sinal da velocidade mas ao fazermos isso alteramos também o ruído que posteriormente faria a partícula tender para baixo porém agora, com a troca de sinal, estará contribuído para cima assim temos uma deslocamento da média para cima. Neste caso temos que o deslocamento para o aditivo é maior do que o do multiplicativo logo o ruído do aditivo é mais intenso.   
[[File:Ruidoalto2.png|400px|Figura 6 - Médias para ruido alto]]

Edição atual tal como às 17h59min de 19 de agosto de 2024

O lançamento oblíquo é um clásico problema da mecanica, onde um projétil e lançado com velocidade inicial em uma direção que faz um angulo com a horizontal. Nesse casso iremos considerar o movimento com arrasto em que a força de resistencia do ar é oposta ao movimento e proporcional a velocidade

O objetivo é descrever o movimento do projétil em duas dimensões, levando em consideração tanto a força da gravidade quanto a resistência do ar então adicionar um termo de ruido no sistema e analizar seu comportamento.

Equações de Movimento

Começamos decompondo as velocidades nas direções x e y

Para as equações de movimento usamos a segunda lei de Newton nas componentes horizontais e verticais

onde:

  • é a massa do projétil
  • é a aceleração da gravidade
  • é o coeficiente de arrasto

Introduzindo e podemos reescrever as equações como diferenciais de primeira ordem.

Essas E.D. possuem solução analitica bastando integrar nos limites adequados. A integração numérica por Euler é dirta a partir das relações a cima e nos da uma trajetoria para comparação quando adicionar o termo de ruido:

Utilizando as seguintes condições iniciais e podemos atualizar a posição

Para modelar o ruído produzido pela densidade do ar, consideramos dois regimes. No primeiro, assumimos que a densidade é constante ao longo da trajetória, ou seja, e, portanto, trataremos o problema como uma EDE com ruído aditivo, onde . Na segunda situação, consideramos que a densidade varia com a altitude do projétil, seguindo a relação e, nesse caso, abordaremos como uma EDE com ruído multiplicativo, onde .

Equação diferencial estocástica

A equação diferencial estocástica a ser resolvida é:

onde

Caso o termo fosse uma constante, exemplo em que a densidade do ar não muda significativamente, o ruído é dito aditivo. Então, a integração da equação é direta, resultando em:

onde é o incremento de Wiener, um processo estocástico com distribuição gaussiana de largura . Escrevendo a equação em diferenças finitas e substituindo os valores de e de têm-se:


Entretanto, o termo depende do processo estocástico e o ruído é dito multiplicativo, assim, devemos utilizar um cálculo estocástco especial, por exemplo o cálculo de Itô e de Stratonovich. No cálculo de Stratonovich, as regras de integração, diferenciação, mudança de variáveis etc, se mantém as mesmas do cálculo usual, e a equação é escrita como:

com o detalhe que o argumento de é, na realidade, calculado na média de X entre os limites do intervalo de integração. Ou seja, ao escrever a equação em diferenças finitas, obtém-se:

Agora, é necessário explicitar o na equação acima, porém a dependência de não é linear em relação à altura, o que torna inviável esse procedimento.


Com isso, a atenção volta-se para o cálculo de Itô, dado por:

No caso em que não depende de , a derivada que aparece na equação é zero e é recuperado o caso de ruído aditivo já mencionado. Aqui, vale ressaltar que é uma matriz e está sendo derivada com relação ao vetor , então o cálculo para as componentes de são:

onde

Após a substituição dos valores de e , obtêm-se as seguintes expressões:

É importante ressaltar que os incrementos de Wiener que aparecem nas equações para velocidade são independentes entre si, isto é, deverá ser gerado duas variáveis aleatórias com distribuição gaussiana de largura . Além disso, um detalhe interessante é que a equação para é semelhante à equação para em que .

Resultados

Aqui temos alguns gráficos comparativos, nesses dois primeiros comparamos os dois modelos para a mesma distribuição de ruídos. Figura 1 - Ruído multiplicativoFigura 2 - Ruído aditivo


Nas figuras abaixo podemos ver as trajetórias de certo objeto para ruído aditivo e ruído multiplicativo para diferentes distribuições.

Figura 1 - Ruído multiplicativoFigura 2 - Ruído aditivo


Esperávamos que as média para ruído aditivo e multiplicativos fossem iguais, dado que o ruído não interfere no cálculo da média pois a probabilidade de puxar em qualquer direção é a mesma, porém como pode ser visto no gráfico abaixo há divergência para o ruído multiplicativo.

Figura 3 - Médias para N = 1000

Tanto a divergência da médias como o aumento da dispersão da trajetórias podem ser explicados pelo modelo que estamos usando na EDE-ito , nesse modelo conforme temos y variando para valores negativos a exponencial explode consequentemente temos aumento na intensidade do ruído, esse aumento pode ser observado no gráfico abaixo.

Figura 7 - exp

Temos que para poucos passos o ganho do ruído em qual direção causara grande divergência, ou seja, para termos vermos a convergência para média do ruído multiplicativo precisamos incrementar o número de trajetórias, a figura a seguir foi feita para um aumento de dez vezes no número de passos comparado ao da figura anterior. Figura 4 - Médias para N = 10000

Nos gráficos abaixo temos as trajetórias para o ruído multiplicativo, as médias para os dois métodos e a solução por euler levando em conta o solo, para as trajetória é visível a dispersão em torno da solução de euler que acontece por causa do ruído, porém como comentado anteriormente essa dispersão não influência na trajetória média, da mesma forma que as duas primeiras figuras conforme a diminuição das velocidades pelo termo dissipativo temos maior influência do ruído no trajeto que pode até mesmo ocasionar em um retorno do objeto na direção da posição inicial porém a média continua somente para a direita.

Figura 5 - Médias para N = 1000


Agora para um maior valor de ruído podemos observar o gráfico a seguir, onde notamos que a média, na região de contato com o solo, fica acima da solução de euler, isso ocorre pela condição de contorno que estamos utilizando, quando a posição da partícula chega em zero invertemos o sinal da velocidade mas ao fazermos isso alteramos também o ruído que posteriormente faria a partícula tender para baixo porém agora, com a troca de sinal, estará contribuído para cima assim temos uma deslocamento da média para cima. Neste caso temos que o deslocamento para o aditivo é maior do que o do multiplicativo logo o ruído do aditivo é mais intenso. Figura 6 - Médias para ruido alto