Soluções de problemas envolvendo equações
diferenciais em que existem incertezas
Fabio Antonio Dorini
Departamento de Matemática, DAMAT, UTFPR, 80230-901, Curitiba, PR
E-mail: [email protected],
Maria Cristina de Castro Cunha,
Departamento de Matemática Aplicada, IMECC, UNICAMP, 13083-970, Campinas, SP
E-mail: [email protected].
Resumo: Este trabalho objetiva analisar, através de três exemplos, a influência de se considerar
aleatoriedades na solução de equações diferenciais com dados e/ou parâmetros aleatórios. Um
comparativo das médias das soluções das “equações estocásticas” com as soluções das “equações
determinı́sticas simplificadas”, nas quais substituı́mos os parâmetros aleatórios por suas médias,
é apresentado.
Palavras-chave: equações diferenciais estocásticas, variáveis aleatórias, esperança matemática,
função de densidade de probabilidade.
Introdução
A modelagem matemática de fenômenos fı́sicos em meios cujas propriedades não são inteiramente conhecidas pode ser realizada utilizando modelos probabilı́sticos (estocásticos) de tais
propriedades. Os modelos probabilı́sticos compensam a escassez de dados e fornecem margens
de erros para estimativas de laboratório. O escoamento em meios porosos se encaixa neste perfil. Propriedades como a porosidade e a condutividade hidráulica são estimadas, por exemplo,
por meio de amostras do meio poroso. Resultados experimentais sugerem que a distribuição
lognormal representa satisfatoriamente estas propriedades [8]. Quando modelos probabilı́sticos
são acoplados aos parâmetros das equações diferenciais que governam o escoamento de fluidos,
o modelo completo passa a ser descrito por equações diferenciais estocásticas (EDE).
Uma das grandes dificuldades do estudo das EDE está na falta de um cálculo diferencial
estocástico sedimentado. A tarefa de adequar os conceitos básicos das equações diferenciais
determinı́sticas ao caso estocástico tem se mostrado árdua; e nada indica que teremos um cálculo
estocástico operacional num futuro próximo.
Na tentativa de resolver problemas deste tipo, duas vertentes têm sido motivos de estudos
nas últimas décadas. A mais formal segue da interpretação das EDE como equações integrais
num sentido mais geral, inicialmente introduzido por Itô, em meados de 1940, e generalizado por
outros autores. A extensão destas fórmulas às equações diferenciais aparece a partir da década
1960. Por exemplo, em [2] são examinados resultados gerais e recentes da análise funcional
de alguns tipos de EDE sob condições de regularidade exigidas dos dados. Na outra vertente
estão métodos e técnicas matemáticas cujo objetivo é obter informações práticas sobre a função
aleatória (estocástica) que é solução do problema, no sentido probabilı́stico. Estas informações
podem ser a média (esperança matemática), alguns momentos de ordem superior ou, melhor
ainda, a função de densidade de probabilidade (FDP) da solução. O método de Monte Carlo [4]
foi o primeiro nesta direção e se tornou conhecido nas simulações da construção das primeiras
bombas atômicas quando, também, os primeiros computadores estavam sendo usados na computação cientı́fica. A ideia do método é calcular médias, ou outros momentos, usando um grande
180
número de ensaios experimentais. Cada um destes ensaios deve resultar em um problema determinı́stico cuja solução possa ser calculada. Em outra direção, alguns autores procuram equações
diferenciais determinı́sticas cujas soluções podem ser médias, momentos, e até mesmo a FDP.
Este trabalho objetiva analisar, através de três exemplos, a influência de se considerar aleatoriedades na solução de equações diferenciais com dados e/ou parâmetros aleatórios. Um estudo
comparativo das médias das soluções das equações estocásticas com as soluções das equações
determinı́sticas simplificadas, nas quais substituı́mos os parâmetros aleatórios por suas médias,
é apresentado.
Considere um espaço de probabilidade (Ω, F, P) e uma variável aleatória X : Ω → R, com
FDP fX . Sua função de distribuição cumulativa é dada por
Z x
FX (x) = P ({ω ∈ Ω | X(ω) ≤ x}) = P(X ≤ x) =
fX (u) du.
−∞
A média, ou esperança matemática, da variável aleatória X é definida por
Z +∞
E[X] =
xfX (x) dx.
−∞
Nos três problemas/exemplos seguintes, um ou mais parâmetros das equações diferenciais são
variáveis aleatórias. Tendo em vista o interesse em analisar a influência dessas aleatoriedades,
é apresentada, para cada exemplo, um comparativo entre a média da solução estocástica e a
solução da equação simplificada, na qual a variável aleatória X é substituı́da pelo número E[X].
1
Problema de valor de contorno aleatório
O ponto de partida é o seguinte problema de valor de contorno:
µ
¶
dU
d
A
(x) = 1, 0 < x < 1,
−
dx
dx
U (0) = U (1) = 0,
(1)
em que A é uma variável aleatória que assume apenas valores positivos, isto é, P(A ≤ 0) = 0. A
aleatoriedade de A implica que a solução do problema aleatório acima é uma função aleatória.
De fato, para cada realização A(ω), de A, considere o problema determinı́stico:

µ
¶
dU
d

A(ω)
(x, ω) = 1, 0 < x < 1,
−
dx
dx

U (0, ω) = U (1, ω) = 0,
cuja solução é obtida integrando duas vezes os lados da equação e impondo as condições de
contorno:
x − x2
, 0 ≤ x ≤ 1.
U (x, ω) =
2 A(ω)
Observe que U (x, ω) > 0 para 0 < x < 1 e ω ∈ Ω. A fim de calcular a FDP, fU (q; x), q > 0,
da solução calcula-se, inicialmente, a função de distribuição cumulativa
¶
¶
¶
µ
µ
µ
x − x2
x − x2
x − x2
≤q =1−P A≤
= 1 − FA
.
FU (q; x) = P (U (x) ≤ q) = P
2A
2q
2q
Derivando FU (q; x) com relação à variável q obtém-se
µ
¶
µ
¶
x − x2
ξ(x)
ξ(x)
x − x2
d
f
=
f
,
fU (q; x) = FU (q; x) =
A
A
dq
2q 2
2q
q2
q
(2)
onde ξ(x) = (x − x2 )/2.
181
1.1
Exemplo
Considere a variável aleatória A, uniformemente distribuı́da no intervalo [1, 2], isto é, fA (q) = 1
quando 1 < q < 2 e é nula fora deste intervalo. De (2) segue que

ξ(x)
 ξ(x)
,
< q < ξ(x),
fU (q; x) =
q2
2

0, caso contrário.
Observe, para fins de comparação, que a esperança matemática de U (x) é dada por
Z +∞
Z ξ(x)
¯
¢
ξ(x)
ln (2) ¡
¯ξ(x)
x − x2 .
E[U (x)] =
qfU (q; x)dq =
q 2 dq = ξ(x) ln (q) ¯ξ(x)/2 =
q
2
0
ξ(x)/2
(3)
Por outro lado, a solução da equação simplificada (usando E[A] = 3/2 no lugar de A em (1)) é
Ũ (x) =
¢
1¡
x − x2 .
3
(4)
Como ln (2)/2 ≈ 0.3466 as funções dadas em (3) e (4) não diferem muito para este caso, como
mostra a Figura 1(à esquerda) - na qual é apresentada, também, a média gerada pelo método
de Monte Carlo, com 100 000 realizações da variável aleatória A. Maiores detalhes sobre como
se obter a média pelo método de Monte Carlo podem ser encontrados em [4] ou [8].
0.1
0.25
0.08
0.2
0.06
0.15
0.04
0.1
0.02
0.05
0
0
0.2
0.4
0.6
0.8
1
x
0
0
0.2
0.4
0.6
0.8
1
x
Figura 1: Comparativo entre E[U (x)] (linha sólida), Ũ (x) (linha pontilhada), e a média da solução
obtida pelo método de Monte Carlo (cı́rculos), 0 ≤ x ≤ 1.
1.2
Exemplo
Considere agora o caso em que A é lognormal, A = exp (ϕ), ϕ ∼ N [0, 1]. A FDP de A é dada
por
µ
¶
[ln (q)]2
1
.
fA (q) = √ exp −
2
q 2π
Tendo em vista (2), obtém-se
µ
µ
¶
¶
[ln (ξ(x)/q)]2
ξ(x)
ξ(x)
1
ξ(x)
exp −
= 2 √
=
fU (q; x) = 2 fA
q
q
q
2
2π[ξ(x)/q]
¶
µ
[ln (q) − ln (ξ(x))]2
1
√
,
(5)
exp −
=
2
q 2π
ou seja, U (x) é também lognormal, U (x) = exp (ψ), ψ ∼ N [ln (ξ(x)), 1]. Usando (5) é possı́vel
calcular todos os momentos estatı́sticos de U (x). Neste caso, denotando o número de Euler, ou
neperiano, por “e”, a esperança matemática de U (x) é dada por
¶
µ
√
¢
¡
¢
√
e¡
1
= e ξ(x) =
x − x2 ≈ 0.8244 x − x2 ,
E[U (x)] = exp ln (ξ(x)) +
2
2
182
√
e a solução da equação simplificada (considerando E[A] = e no lugar de A em (1)) é
¢
¡
¢
1 ¡
Ũ (x) = √ x − x2 ≈ 0.3033 x − x2 .
2 e
Na Figura 1(à direita) é esboçado E[U (x)] e Ũ (x), bem como a média da solução gerada pelo
método de Monte Carlo, com 100 000 realizações de A. Maiores detalhes sobre momentos estatı́sticos de uma variável aleatória lognormal podem ser encontrados em [8].
2
Equação da advecção linear aleatória
Considere, agora, a equação da advecção linear aleatória
∂
∂
U (x, t) + V
U (x, t) = 0,
∂t
∂x
U (x, 0) = U0 (x),
t > 0,
x ∈ R,
(6)
onde a velocidade, V , é uma variável aleatória e a condição inicial, U0 (x), é uma função aleatória.
Em [7] é utilizado o método das caracterı́sticas [5] e a lei de probabilidade total [6] para mostrar
que a FDP da solução de (6) no ponto (x, t), fQ (q; x, t), é dada por
µ
¶
Z ∞
x − x0
1
fV
fU0 (q; x0 ) dx0 =
fU (q; x, t) =
(7)
t
−∞ t
Z ∞
=
fV (v) fU0 (q; x − vt) dv = EV [fU0 (q; x − V t)] .
−∞
Além disso, no caso em que U0 (x) = u(x) é uma função determinı́stica, (7) pode ser apresentada
como
Z
+∞
fU (q; x, t) =
−∞
fV (v) δ (q − u(x − vt)) dv,
(8)
onde δ é a distribuição delta de Dirac.
2.1
Exemplo
Considere o caso em que V ∼ N (µ, σ 2 ) e a condição inicial u(x) é a função de Heaviside,
u(x) = 1, se x > 0, e u(x) = 0, se x < 0. Esta condição inicial tem sido usada por vários autores
para estudar zonas de mistura em concentrações de substâncias. De (8), segue que
Z x/t
Z +∞
fU (q; x, t) = δ(q − 1)
fV (v)dv + δ(q)
fV (v)dv,
−∞
x/t
isto é, U (x, t) é a variável aleatória de Bernoulli com
¸
·
Z x/t
Z x/t
(ξ − µ)2
1
dξ =
exp −
P(U (x, t) = 1) =
fV (v)dv = √
2σ 2
2πσ −∞
−∞
µ
¶
Z (x−µt)/(√2σt)
x − µt
1
1
−τ 2
=√
,
e dτ = 1 − erfc √
2
π −∞
2σt
onde erfc(x) é a função erro complementar. Além disso,
1
E[U (x, t)] = 1 . P(U (x, t) = 1) + 0 . P(U (x, t) = 0) = 1 − erfc
2
µ
x − µt
√
2σt
¶
.
A Figura 2 ilustra a média, E[U (x, t)], da solução de (6) em t = 0.5 e t = 0.1, −4 ≤ x ≤ 4, com
V ∼ N (0.5, 1). Também, apresenta a solução da equação simplificada Ũ (x, t) = u(x − E[V ]t) =
u(x − 0.5t), e a solução produzida pelo método de Monte Carlo, com 10 000 realizações de V . A
difusão observada na Figura 2, causada pela aleatoriedade de V, é comprovada pela Proposição
2.1 e pelo Colorário 2.1, apresentados a seguir:
183
1.5
1.5
1
1
0.5
0.5
0
0
−0.5
−4
−2
0
x
2
4
−0.5
−4
−2
0
x
2
4
Figura 2: Comparativo entre E[U (x, t)] (linha sólida), Ũ (x, t) (linha pontilhada), e a média da solução
obtida pelo método de Monte Carlo (cı́rculos); t = 0.5 (à esquerda) e t = 1.0 (à direita).
Proposição 2.1. No caso em que V ∼ N (µ, σ 2 ), a FDP de U (x, t), fU (q; x, t), satisfaz à
seguinte equação de advecção-difusão
∂
∂2
∂
fU (q; x, t) + µ
fU (q; x, t) = σ 2 t 2 fU (q; x, t).
∂t
∂x
∂x
(9)
Demonstração. Derivando (7) em x e t observa-se que (9) é válida se a função h(x, t) =
(1/t) fV ((x − x0 )/t) for tal que
∂
∂2
∂
h(x, t) + µ
h(x, t) = σ 2 t 2 h(x, t).
∂t
∂x
∂x
Como V ∼ N (µ, σ 2 ), segue que sua FDP é dada por
¸
·
1
(x − µ)2
,
fV (x) = √
exp −
2σ 2
2πσ
cujas derivadas são
d
fV (x) = −fV (x)
dx
·
x−µ
σ2
¸
e
·
¸
d2
fV (x)
(x − µ)2
fV (x) =
−1 +
.
dx2
σ2
σ2
Deste modo,
¸
·
∂
1
(x − x0 )(x − x0 − µt)
,
h(x, t) = 2 h(x, t) −t +
∂t
t
σ2t
¸
·
∂
1
x − x0 − µt
e
h(x, t) = − 2 h(x, t)
∂x
t
σ2
¸
·
∂2
1
(x − x0 − µt)2
.
h(x, t) = 2 3 h(x, t) −t +
∂x2
σ t
σ2t
Portanto,
¸
·
∂
∂
1
(x − x0 )(x − x0 − µt) µt (x − x0 − µt)
−
=
h(x, t) + µ
h(x, t) = 2 h(x, t) −t +
∂t
∂x
t
σ2 t
σ2 t
¸
·
2
1
(x − x0 − µt)2
2 ∂
=
σ
t
= 2 h(x, t) −t +
h(x, t).
t
σ2t
∂x2
Corolário 2.1. A média, E[U (x, t)], da solução de (6), com V ∼ N (µ, σ 2 ), satisfaz à equação
de advecção-difusão (9), isto é,
∂
∂
∂2
E[U (x, t)] + µ
E[U (x, t)] = σ 2 t 2 E[U (x, t)].
∂t
∂x
∂x
184
3
Problema Burgers-Riemann aleátório
Considere o caso não-linear obtido pelo problema de Riemann aleatório para a equação de
Burgers (introduzida por J. M. Burgers em [1])
∂
1 ∂ 2
U (x, t) +
U (x, t) = 0, t > 0,
∂t
½ 2 ∂x
UL , se x < 0,
U (x, 0) =
UR , se x > 0,
x ∈ R,
(10)
em que os estados iniciais, UL e UR , são variáveis aleatórias. Observe que aqui a aleatoriedade
aparece apenas na condição inicial, transformando o problema como um todo numa equação
diferencial aleatória. Para uma simples realização UL (ω) e UR (ω), de UL e UR , respectivamente,
tem-se o problema de Burgers-Riemann determinı́stico:
∂
1 ∂ 2
u(x, t, ω) +
u (x, t, ω) = 0, t > 0,
∂t
½ 2 ∂x
UL (ω), se x < 0,
u(x, 0, ω) =
UR (ω), se x > 0.
x ∈ R,
(11)
Soluções fisicamente corretas de (11), isto é, soluções de entropia, são ondas de rarefação ou
ondas de choque [5]. Em [3] é mostrado que a solução de (10), em (x, t), é a função aleatória
U (x, t) = UL XR− + βXR0 + UR XR+ ,
onde β = x/t, e XR− , XR0 e XR+ são as funções caracterı́sticas dos conjuntos mutuamente
exclusivos definidos na Figura 3. Observe que a função caracterı́stica de R0 , por exemplo, é
definida por: XR0 = 1 se (UL , UR ) ∈ R0 , e XR0 = 0 se (UL , UR ) 6∈ R0 .
UR
R0 (β)
P
β
β
R− (β)
UL
2β
R+ (β)
Figura 3: Ilustração das regiões de integração.
O m-ésimo momento estatı́stico de U (x, t) é dado por
Z Z
E[U m (x, t)] =
um
L fUL UR (uL , uR )duL duR +
−
R
Z Z
Z Z
m
+β
fUL UR (uL , uR )duL duR +
R0
R+
um
R fUL UR (uL , uR )duL duR .
Para a obtenção dos momentos estatı́sticos da solução de (10), usando a expressão acima,
deve-se calcular três integrais duplas para cada valor de β. Em alguns casos é possı́vel calcular
tais integrais exatamente, como ilustrado no exemplo a seguir.
185
3.1
Exemplo
Considere o caso em que UL e UR são variáveis aleatórias independentes e uniformemente distribuı́das no intervalo [−a, a]. Alguns cálculos mostram que a média da solução de (10) é dada
por
(
β
− 2 (sgn(β)β − a)2 , if − a ≤ β ≤ a,
E[U (x, t)] =
4a
0,
caso contrário.
com β = x/t, sgn(β) = −1 se β < 0, sgn(β) = 0 se β = 0, e sgn(β) = 1 se β > 0. Na Figura
4 é ilustrada a média, E[U (x, t)], da solução de (10) em t = 0.4 e t = 0.8, −1 ≤ x ≤ 1. UL e
UR são independentes e uniformemente distribuı́das em [−1, 1]. Convém ressaltar que a solução
da equação simplificada de (10), com E[UL ] = E[UR ] = 0 no lugar de UL e UR , é identicamente
nula.
0.06
0.06
0.04
0.04
0.02
0.02
0
0
−0.02
−0.02
−0.04
−0.04
−0.06
−1
−0.5
0
0.5
1
−0.06
−1
−0.5
0
0.5
1
Figura 4: Ilustração de E[U (x, t)] em t = 0.4 (à esquerda) e t = 0.8 (à direita), −1 ≤ x ≤ 1.
Agradecimentos
Os autores agradecem o apoio financeiro do CNPq (Processo 472992/2008-2) e ao Prof. Saulo
Pomponet Oliveira (DMAT-UFPR), pelas interessantes discussões (sobre o Problema 1) .
Referências
[1] J. M. Burgers, A mathematical model illustrating the theory of turbulance. Ad. Appl. Mech.
1 (1948) 171–179.
[2] P. L. Chow, “Stochastic Partial Differential Equations”, Chapman & Hall/CRC, New York,
2007.
[3] M. C. C. Cunha, F. A. Dorini, Statistical moments of the solution of the random BurgersRiemann problem. Mathematics and Computer in Simulation 79 (2009) 1440–1451.
[4] G. S. Fishman, “Monte Carlo: concepts, algorithms and applications”, Springer-Verlag,
New York, 1996.
[5] R. J. LeVeque, “Numerical Methods for Conservation Laws”, Birkhäuser, Berlin, 1992.
[6] A. Papoulis, “Probability, Random Variables, and Stochastic Processes”, McGraw-Hill,
Inc., New York, 1984.
[7] L. T. Santos, F. A. Dorini, M. C. C. Cunha, The probability density function to the random
linear transport equation. Applied Mathematics and Computation 216 (2010) 1524–1530.
[8] D. Zhang, “Stochastic Methods for Flow in Porous Media - Coping with Uncertainties”,
Academic Press, San Diego, 2002.
186