SISTEMAS LINEARES ESTOCÁSTICOS COM ATRASO GOVERNADO POR
PROCESSO DE BERNOULLI
Alessandro N. Vargas∗, João B. R. do Val†
∗
Universidade Tecnológica Federal do Paraná, UTFPR,
Av. Alberto Carazzai 1640, 86300-000 Cornelio Procópio-PR, Brasil;
∗
Basque Center for Applied Mathematics, BCAM,
Alameda de Mazarredo 14, E-48009 Bilbao, Vizcaya, España
(de 1o. janeiro à 28 de fevereiro de 2011).
†
Universidade Estadual de Campinas (UNICAMP)
Depto. de Telemática-FEEC, 13081-970, Campinas, SP, Brasil
Emails: [email protected], [email protected]
Resumo— Este artigo apresenta resultados que caracterizam a dinâmica do valor esperado do estado de
sistemas lineares estocásticos sujeitos a atraso. O atraso no sistema depende de um processo aleatório de Bernoulli,
com uma determinada probabilidade de ocorrência de sucessos e falhas, com aplicações imediatas em sistemas
sujeitos a falhas nos canais de transmissão. Apresentamos uma aplicação real, realizando o controle de nı́vel
num tanque, que emprega transmissão de dados sob uma rede sujeita a perda de pacotes, ilustrando assim a
potencialidade do uso dos resultados desenvolvidos.
Palavras-chave—
sistemas estocásticos, sistemas com atraso, processos sujeitos a falhas
Abstract— This paper presents results that characterize the dynamics of the expected value operator of the
system state for linear stochastic systems subject to delay. The delay in the system depends on a random
Bernoulli process with a fixed probability of ocurrence of successes and failures, specially adapted for systems
subject to failures on the communication links. We also present a real application that implements the control
level in a tank, considering the data transmitted over a network subject to packet dropouts, thus illustrating the
usefulness of the derived approach.
stochastic systems, systems with delay, processes subject to failures
Keywords—
1
Introdução
Sistemas sujeitos a atraso e falhas nas linhas de comunicação são muito comuns na prática, e muito
esforço têm-se empreendido nos últimos anos no
estudo e desenvolvimento de novos resultados para
essa classe de sistemas, veja por exemplo os trabalhos (Mahmoud, 2010), (Verriest and Michiels, 2009), (Hespanha and Naghshtabrizi, 2007),
(Xia et al., 2009), (Moayedi et al., 2011) como
uma pequena amostra. O resultado deste trabalho
contribui nessa direção por considerar um caso especial de sistemas lineares com atraso dependendo
de uma representação particular de sucessos e falhas na comunicação e governado por um processo
de Bernoulli.
Vamos agora apresentar o sistema sob estudo. Considere um espaço de probabilidade fixado (Ω, F, P ), e associado a este considere o sistema
x(k + 1) = Ax(k) + Bx(k − δk ) + Dw(k),
(1)
∀k ≥ 0, com x0 ∈ Rn e δ0 = 0, no qual xk e
wk , k ≥ 0 são processos tomando valores respectivamente em Rr e Rq . A entrada de ruı́do {wk }
Pesquisa financiada em parte pela FAPESP Processo 03/06736-7 e CNPq Processos 471557/2009-9 e
304856/2007-0, e em parte pela Fundación Carolina (Madrid, Espanha) – Programa “Movilidad de profesores e investigadores Brasil-España. C.2010”.
forma um processo iid com média nula e matrix
de covariância igual a identidade para cada k ≥ 0.
Associado a (1), consideramos um processo aleatório usual de Bernoulli com probabilidade de sucessos e falhas dada por p e 1−p, respectivamente,
e tomando valores “0” e “1”.
Os instantes de tempo
0 < k 1 < k2 < . . . < ki < . . .
representam o tempo no qual a informação é recebida com sucesso pelo sistema (1). De fato, o
processo de contagem das variáveis aleatórias de
Bernoulli formam um processo binomial que se encarrega da contagem de fracassos correspondentes
ao atraso, e quando ocorre o sucesso o instante
ki registra essa ocorrencia e reinicia o processo de
contagem fazendo δki = 0, veja a Fig. 1. Note
que o processo de atraso satisfaz a regra
δk = k − ki ,
k = ki , . . . , ki+1 − 1.
(2)
A principal contribuição deste trabalho é caracterizar a expressão da dinâmica do valor esperado de x(k). Em outras palavras, conseguimos
determinar a função ϕk : Rn × · · · × Rn 7→ Rn
para cada k > 0 tal que
E[x(k)] = ϕk (E[x(k − 1)], E[x(k − 2)], . . . , x0 ),
veja a expressão em (17) como comparação. Cabe
ressaltar que a função ϕk , para ser bem definida,
ou ainda de modo equivalente a
8

7
δk
6
x(ki + s) = As +
5
4
s−1
X
j=0
+

Aj B  x(ki )
s−1
X
As−j Dw(ki + j)
(4)
j=0
3
2
sempre que s ∈ {1, . . . , ki+1 − ki }. Visando simplificar a notação, defina
1
0
k0
k1 k2
k3
k4k5k6
s
M (s) = A +
s−1
X
Aj B,
s = 1, 2, . . . ,
(5)
j=0
Figura 1: Uma realização para o processo aleatório
de atraso {δk } com reinı́cio nos instantes de chegada
(sucesso) ki , i ≥ 0. O sucesso ou falha de Bernoulli
define o tempo de ocorrencia de cada ki .
requer que o raio espectral da matriz pA seja estritamente menor que um. Além disso, as expressões
encontradas nesse trabalho possibilitam a avaliação numérica de E[x(k)], e como subproduto é
possı́vel avaliar a estabilidade ou estabilidade assintótica de E[x(k)] quando k tende a infinito.
O resultado deste trabalho possui interesse
particular para o controle de sistemas via rede
com perdas de pacotes, pois a representação do
atraso que adotamos é geral o suficiente e pode ser
usada para incorporar falhas nas linhas de comunicação. Para ilustrar essa propriedade, apresentamos nesse artigo uma aplicação real que controla
o nı́vel de água num tanque através de uma rede
emulada sujeita a perda de pacotes de informação,
veja a Seção 3.
Este artigo está organizado da seguinte maneira. A Seção 2 apresenta a notação necessária, definições, e resultado principal. Finalmente,
apresentamos na Seção 3 uma aplicação num processo de nı́vel num tanque.
2
ν(s, k) =
s−1
X
As−j Dw(k + j),
k, s = 1, 2, . . . ,
j=0
(6)
com ν(0, k) = 0 para cada k ≥ 0. Observe desde
(4)-(6) que
x(ki + s) = M (s)x(ki ) + ν(s, ki ).
(7)
O processo {δk } é governado, por hipótese,
pelo processo de Bernoulli iid, e por isto o tempo
entre chegadas (definindo o tempo entre a ocorrência de dois sucessos sequenciais) é definido por
si = ki+1 − ki ,
∀i ≥ 0,
(8)
que por sua vez segue a distribuição geométrica
(Bertsekas and Tsitsiklis, 2002, Ch. 5).
Pr(si = n) = p(1 − p)n−1 ,
∀i ≥ 0, ∀n ≥ 1. (9)
Combinando (4)-(8), temos
x(ki+1 ) = M (si )x(ki ) + ν(si , ki ),
∀i ≥ 0. (10)
Comentário 1 Afirmamos que
E[ν(si , ki )] = 0,
i = 0, 1, . . . .
(11)
De fato, através da hipótese que o processo {si }
é iid e independente do processo de ruı́do aditivo
{w(ki )}, podemos escrever
Notação e resultado principal
Os números reais e naturais são denotados respectivamente por R e N. O superscript ′ indica
o transposto da matriz; k · k representa a norma
Euclideana em Rn . Seja σ(·) o operador de raio
espectral. A função indicadora do conjunto C é
representada por 11C .
Vamos agora relembrar o sistema (1). A identidade em (2) nos habilita a reescrever (1) de modo
equivalente a
k ∈ {ki , . . . , ki+1 − 1}
⇒ x(k + 1) = Ax(k) + Bx(ki ) + Dw(k),
com M (0) = I e
(3)

E
sX
i −1
j=0

Asi −j Dw(ki + j)

= E
sX
i −1
j=0

Asi −j D E[w(ki + j)].
Essa identidade e a hipótese que E[wk ] ≡ 0 produzem o resultado em (11).
O próximo resultado define uma matriz determinı́stica que será útil mais adiante.
Lema 1 Se σ(pA) < 1, então é válida
E[M (si )] = H(p),
i = 0, 1, . . .
(12)
Aplicando o operador de valor esperado em
ambos os lados de (10), e considerando (11) e (16),
obtemos
E[x(ki+1 )] = E[M (si )x(ki )] + E[ν(si , ki )]
no qual q = 1 − p e
= E[M (si )]E[x(ki )].
H(p) = pA[I + qA(I − qA)−1 ]
+ [I + qA(I − qA)−1 ]B.
(13)
Prova: Pela definição,
E[M (si )] =
∞
X
O resultado segue imediatamente do Lema 1.
✷
Estamos agora em condições de apresentar o
principal resultado deste artigo.
Teorema 3 Se σ(pA) < 1, então para cada k =
0, 1, . . ., temos que a fórmula em (17) é válida.
M (s) Pr(si = s).
s=1
Prova: Considere um instante de tempo k fixado.
Então pode ocorrer i = 0, 1, . . . , k sucessos até o
instante k. Considerando que ki denota o instante
do i-ésimo sucesso, podemos escrever a igualdade
Usando (5) e (9) na identidade acima, obtemos
E[M (si )] =
∞
X
M (s)q s−1 p
s=1
=p
∞
X
s=0

q s As+1 +
s
X
j=0

Aj B  .
E[x(k)] =
(14)
s=1 j=0
(15)
O segundo termo dentro dos parenteses em
(15) é igual a
qA2 (I − qA)−1 ,
e o terceiro é igual a
E[x(k)11ki ≤k<ki+1 ]
= E[x(k)|ki ≤ k < ki+1 ] Pr{ki ≤ k < ki+1 }.
(19)
A probabilidade de ocorrência do evento {ki ≤
k < ki+1 } é equivalente a ocorrerem i sucessos em
k tentativas, que por sua vez representa uma v.a.
com distribuição Binomial dada por (Bertsekas
and Tsitsiklis, 2002, Sec. 5.1)
k i
Pr{ki ≤ k < ki+1 } =
p (1 − p)k−i . (20)
i
Pr{k0 ≤ k < k1 } = (1 − p)k ,
Estes fatos nos habilitam a concluir que
+ qB + qA(I − qA)−1 B
= H(p),
x(k) = M (k)x0 + ν(k, 0).
✷
(22)
Podemos concluir de (19)-(22) com i = 0 que
O próximo resultado mostra que a matriz
H(p) é suficiente para caracterizar a dinâmica do
valor esperado do estado no caso especı́fico do
tempo entre chegadas.
E[x(k)11k0 ≤k<k1 ]
= E[x(k)|k0 ≤ k < k1 ] Pr{k0 ≤ k < k1 }
= E[M (k)x0 |k0 ≤ k < k1 ](1 − p)k
= (1 − p)k M (k)x0 .
Lema 2 Se σ(pA) < 1, então
(23)
Esta última igualdade, combinada com (18)–(20),
nos leva a concluir que
∀i ≥ 0,
no qual H(p) satisfaz (13).
E[x(k)] = (1 − p)k M (k)x0
Prova: Da hipótese de que si e x(ki ) são independentes para cada i ≥ 0, temos que
+
k
X
i=1
E[M (si )x(ki )] = E[M (si )]E[x(ki )].
(21)
o que significa que nos instantes 1, 2, . . . , k ocorrem somente fracassos. Isto nos permite obter de
(7) que (relembre que k0 = 0)
E[M (si )] = p(A + B) + pqA2 (I − qA)−1
que por sua vez prova o resultado.
(18)
No caso particular quando i = 0, temos de
(20) que
q
I + A(I − qA)−1 B.
p
E[x(ki+1 )] = H(p)E[x(ki )],
E[x(k)11ki ≤k<ki+1 ].
i=0
Note que
Fazendo um simples reajuste nos termos de (14),
é possı́vel mostrar que E[M (si )] é idêntico a


s
∞ X
∞
X
X
q s Aj B  .
(qA)s +
p (A + B) + A
s=1
k
X
(16)
E[x(k)|ki ≤ k < ki+1 ]
k i
p (1 − p)k−i .
i
(24)
k
E[x(k)] = (1 − p) M (k)x0 +
k
X
i=1
k−1
X
t−1 i
M (k − t)E[x(t)]
p (1 − p)t−i
i
−
1
t=i
! k−1 i
k i
k−i
i
p (1 − p)
p (1 − p)k−i .
+ H(p) x0
i−1
i
Considere agora i > 0. Perceba que o evento
{ki ≤ k < ki+1 , ki = t} implica que no instante
ki = t ocorre um sucesso e nos instantes subsequentes t + 1, . . . , k ocorrem fracassos. Temos então de (7) que
x(k) = M (k − t)x(t) + ν(k − t, t).
(25)
(17)
h(k)
Network
Podemos agora empregar (25) para escrever a
identidade
E[x(k)|ki ≤ k < ki+1 , ki = t]
= E[M (k − t)x(t) + ν(k − t)]
= M (k − t)E[x(t)].
(26)
Figura 2: Diagrama esquemático representando o
No caso particular em que ki = k em (26), temos
do Lema 2 que
processo de nı́vel de água e as linhas de comunicação correspondentes entre o processo e o computador,
de acordo com a aplicação da Seção 3.
E[x(k)|ki ≤ k < ki+1 , ki = k]
= E[x(ki )] = H(p)i x(k0 ) = H(p)i x0 .
(27)
Todavia, devido a igualdade condicional
E[x(k)|ki ≤ k < ki+1 ]
=
k
X
E[x(k)|ki ≤ k < ki+1 , ki = t] Pr{ki = t}
t=i
podemos recuperar as expressões em (26) e (27)
para escrever
E[x(k)|ki ≤ k < ki+1 ]
=
k−1
X
M (k − t)E[x(t)] Pr{ki = t}
t=i
+ H(p)i x0 Pr{ki = k}.
(28)
Note que Pr{ki = t}, que representa a probabilidade de ocorrencia do i-ésimo sucesso no instante t, obedece a distribuição de Pascal dada por
(Bertsekas and Tsitsiklis, 2002, Sec. 5.1)
t−1 i
Pr{ki = t} =
p (1 − p)t−i , t = i, . . . , k.
i−1
(29)
Finalmente, juntando (24), (28) e (29), obtemos imediatamente o resultado em (17).
✷
3
Aplicação: Processo de nı́vel de água
numa rede sujeita a perda de pacotes
É bastante conhecido que perda de pacotes ocorrem frequentemente em rede de dados, e em tais
situações é extremamente desejável verificar o
comportamento do sistema, e como via de regra a estabilidade do sistema global é um dos
principais objetivos a serem alcançados no projeto (Hespanha and Naghshtabrizi, 2007), (Zhang
and Yu, 2007), (Sun and Qin, 2011), (Ishido
et al., 2011), (Moayedi et al., 2011).
De modo a verificar o comportamento de
longo prazo do sistema sujeito a perda de pacotes
na comunicação, e como subproduto a sua estabilidade em termos do valor esperado do estado,
apresentamos uma aplicação real de controle de nı́vel num tanque no qual as linhas de comunicação
da realimentação dependem de linhas de transmissão que estão sujeitas a perda de pacotes, veja Fig.
2 para uma representação pictorial. O processo de
tentativas de Bernoulli emula uma rede através da
implementação do esquema de falhas e sucessos no
transito das informações entre o computador e o
sistema de nı́vel.
O aparato experimental usado é baseado no
Módulo 2325 - Controle de Nı́vel e Temperatura fabricado pela empresa Datapool Eletrônica Ltda,
Brasil, associado com um cartão de aquisição de
dados da National Instruments USB-6008 para
implementar a linha fı́sica com o computador, veja
a Fig. 3. A rede é emulada no software Matlab,
que implementa as linhas de realimentação, observando que estas estão sujeitas à variável de Bernoulli com probabilidade de sucessos fixada em p.
O computador é também usado para implementar
fisicamente um controlador remoto phase-lag. O
p = 0.01
p = 0.3
p = 0.7
p = 0.99
kE[x(k)]k
300
200
100
0
0
20
40
60
k
Figura 3: Aparato experimental utilizado no procedimento descrito na Seção 3.
perı́odo de amostragem usado é de 35 milisegundos.
O sistema linear estocástico representando o
processo de nı́vel de água no tanque é dado por
x(k + 1) = Ax(k) + By(k) + Dω(k),
y(k) = Cx(k − δk ), k ≥ 0,
Figura 4: Valor da norma do valor esperado de x(k),
para cada k ≥ 0, considerando quatro distintos valores
para a probabilidade p, de acordo com o modelo do
processo de nı́vel de água da Seção 3.
Nı́vel de água medido no tanque
12.8
p = 0.99
p = 0.01
11.2
(30)
h(k)
no qual {y(k)} representa a saı́da do sistema, e os
demais elementos coincidem com (1). Os parâmetros do sistema são
1.996439
−0.996454
A=
,
0.999996 −3.518042 · 10−6
1.878949
B = 10−3 ×
,
1.876749
C = 10−3 × 1.876762 1.874540 , D = 0.05I.
O nı́vel de água dentro do tanque, dado em centı́metros, satisfaz a expressão
h(k) = 4y(k) + 11.64,
9.6
8
0
4
k(×104 )
8
12
Figura 5: Nı́vel de água medido no tanque, dado em
centı́metros, para duas realizações distintas do processo de Bernoulli com duas probabilidades de sucesso
opostas, de acordo com a aplicação da Seção 3.
∀k ≥ 0.
Avaliando numericamente a fórmula (17) para
o modelo em (30), percebemos que para cada 0 <
p < 1 o sistema (30) é assintoticamente estável na
média (Arnold, 1974, p.188), (Vargas and do Val,
2010), ou seja,
lim kE[x(k)]k = 0,
k→∞
veja a Fig. 4 para a representação gráfica que
considera quatro probabilidades distintas.
Este comportamento de estabilidade obtido
numericamente também é confirmado no experimento, veja a Fig. 5. A figura apresenta a resposta real do processo considerando transmissão
de dados pela rede com probabilidades de sucesso alta e baixa, cujos valores são p = 0.99
e p = 0.01 respectivamente. Percebemos claramente que o controle no nı́vel de água no tanque
é mais pobre no segundo caso quando comparado
com o primeito, e este fato é razoável pois uma pequena probabilidade de sucessos representa uma
alta probabilidade de falhas que torna o ativamento do atuador da bomba de controle menos
eficiente.
Em resumo, mesmo no cenário real de falhas
de transmissão de dados em sistemas de controle
via rede, o equipamento mantém um comportamento estável conforme indicado pela avaliação
numérica obtida através da expressão em (17).
4
Comentários finais
Este trabalho desenvolve a expressão da dinâmica
do valor esperado do estado do sistema (1), ou de
modo equivalente, o resultado do Teorema 3 provê
a expressão para calcular E[x(k)]; essa expressão
é função de
Mahmoud, M. S. (2010). Switched Time-Delay
Systems: Stability and Control, Springer.
E[x(k − 1)], E[x(k − 2)], . . . , x0 .
Moayedi, M., Foo, Y. and Soh, Y. (2011). Networked LQG control over unreliable channels,
Int. J. Robust Nonlinear Control . doi:
10.1002/rnc.1822.
Como o sistema estocástico considera o atraso
governado por um processo de Bernoulli, as aplicações para sistemas de controle via rede sujeitos
a perda de pacotes são evidentes e uma aplicação
real é apresentada para ilustrar este fato.
Referências
Arnold, L. (1974). Stochastic differential equations: Theory and applications, WileyInterscience.
Bertsekas, D. P. and Tsitsiklis, J. N. (2002). Introduction to Probability, Athena Scientific; 1st
ed.
Hespanha, J. P. and Naghshtabrizi, Y. X. (2007).
A survey of recent results in networked control systems, Proc. IEEE Special Issue on
Technology of Networked Control Systems,
pp. 138–162.
Ishido, Y., Takaba, K. and Quevedo, D. E. (2011).
Stability analysis of networked control systems subject to packet-dropouts and finitelevel quantization, Systems Control Lett.
60(5): 325–332.
Sun, Y. and Qin, S. (2011). Stability of networked
control systems with packet dropout: an average dwell time approach, IET Control Theory Appl. 5(1): 47–53.
Vargas, A. N. and do Val, J. B. R. (2010). Average cost and stability of time-varying linear systems, IEEE Trans. Automat. Control
55: 714–720.
Verriest, E. I. and Michiels, W. (2009). Stability analysis of systems with stochastically
varying delays, Systems Control Lett. 58(1011): 783–791.
Xia, Y., Fu, M. and Shi, P. (2009). Analysis and
Synthesis of Dynamical Systems with TimeDelays, Springer. (Lecture Notes in Control
and Information Sciences).
Zhang, W.-A. and Yu, L. (2007). Output feedback stabilization of networked control systems with packet dropouts, IEEE Trans. Automat. Control 52(6): 1705–1710.