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

De Física Computacional
Ir para navegação Ir para pesquisar
 
(9 revisões intermediárias por 2 usuários não estão sendo mostradas)
Linha 79: 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>
 
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>
</math>


Linha 95: Linha 115:


<math>
<math>
\Delta v_{x_i} = -k v_{x_i} \Delta t + \beta \Delta W_i
\Delta v_{x_i} = -k v_{x_i} \Delta t + \beta \Delta W_{x_i}
</math>
</math>


Linha 103: Linha 123:


<math>
<math>
\Delta v_{y_i} = \left(-g  -k v_{x_i} \right) \Delta t + \beta \Delta W_i
\Delta v_{y_i} = \left(-g  -k v_{x_i} \right) \Delta t + \beta \Delta W_{x_i}
</math>
</math>




Entretanto, o termo <math>B(X,t)</math> depende do processo estocástico e 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:
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 119: Linha 139:
</math>
</math>


Ao subsituir os valores de <math>A(X(t))</math> e de <math>B(X(t))</math> na equação acima, chega-se em:
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.
 


<math>
Com isso, a atenção volta-se para o cálculo de Itô, dado por:
\Delta x_i = v_{x_i} \Delta t
</math>


<math>
<math>
\Delta v_{x_i} = -k v_{x_i} \Delta t + \beta e^{-\frac{y_i}{H}} e^{- \frac{\Delta y_i}{2H}} \Delta W_i
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>


<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:
\Delta y_i = v_{y_i} \Delta t
</math>


<math>
<math>
\Delta v_{y_i} = \left(-g  -k v_{x_i} \right) \Delta t + \beta e^{-\frac{y_i}{H}} e^{-\frac{\Delta y_i}{2H}} \Delta W_i
dX_i(t) = A^{(I)}_i(X(t))dt + \sum_{\alpha=1}^m B_{i,\alpha}(X(t)) \bullet dW_\alpha(t),
</math>
</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>. Agora, é necessário explicitar os <math>\Delta y_i</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.
onde
 
 
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)
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>
</math>


Após a substituição dos valores de <math>A(X(t))</math> e <math>B(X(t))</math>, obtêm-se as seguintes expressões:
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>
<math>
Linha 153: Linha 167:


<math>
<math>
\Delta v_{x_i} = -k v_{x_i} \Delta t + \beta e^{- \frac{y(t)}{H}} \Delta W_i
\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>


Linha 161: Linha 175:


<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_i
\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>
</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 ==
== Resultados ==


Nas figuras abaixo podemos ver as trajetórias de certo objeto para ruído aditivo e ruído multiplicativo considerando o chão em zero.   
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 multiplicativo.gif|400px|Figura 1 - Ruído multiplicativo]][[File:Ruido aditivo.gif|400px|Figura 2 - Ruído aditivo]]
[[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.
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 SCC.png|400px|Figura 3 - Médias para N = 1000]]
[[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.
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.
Linha 180: Linha 200:


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.   
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:Medias N=10000 SCC.png|400px|Figura 4 - Médias para N = 10000]]
[[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 esquerda.       
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:Medias N=1000 CCC.png|400px|Figura 5 - Médias para N = 1000]]
[[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.     
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:Ruidoalto.png|400px|Figura 6 - Médias para ruido alto]]
[[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