KALMAN AND PARTICLE FILTERS Tutorial 10 H. R. B. Orlande, M. J. - - PowerPoint PPT Presentation

kalman and particle filters
SMART_READER_LITE
LIVE PREVIEW

KALMAN AND PARTICLE FILTERS Tutorial 10 H. R. B. Orlande, M. J. - - PowerPoint PPT Presentation

KALMAN AND PARTICLE FILTERS Tutorial 10 H. R. B. Orlande, M. J. Colao, G. S. Dulikravich, F. L. V. Vianna, W. B. da Silva, H. M. da Fonseca and O. Fudym Department of Mechanical Engineering Escola Politcnica/COPPE Federal University of Rio de


slide-1
SLIDE 1

KALMAN AND PARTICLE FILTERS

  • H. R. B. Orlande, M. J. Colaço, G. S. Dulikravich,
  • F. L. V. Vianna, W. B. da Silva, H. M. da Fonseca and O. Fudym

Department of Mechanical Engineering Escola Politécnica/COPPE Federal University of Rio de Janeiro, UFRJ Rio de Janeiro, RJ, Brazil helcio@mecanica.ufrj.br

1

Tutorial 10

slide-2
SLIDE 2

SUMMARY

  • State estimation problems
  • Kalman filter
  • Particle filter
  • Examples: Lumped system and 1D Heat conduction
  • Applications
  • Conclusions
slide-3
SLIDE 3

STATE ESTIMATION PROBLEM State Evolution Model: Observation Model:

( , )

k k k k

 z h x n

x

n

R  x Subscript k = 1, 2, …, denotes an instant tk in a dynamic problem = state variables to be estimated

v

n

R  v

z

n

R  z

n

n

R  n = state noise = measurements = measurement noise xk = fk (xk-1, uk-1, vk-1) u  R np = input variable

slide-4
SLIDE 4

STATE ESTIMATION PROBLEM Definition: The state estimation problem aims at

  • btaining information about xk based on the state

evolution model and on the measurements given by the observation model.

slide-5
SLIDE 5

The solution of the inverse problem within the Bayesian framework is recast in the form of statistical inference from the posterior probability density, which is the model for the conditional probability distribution of the unknown parameters given the

  • measurements. The measurement model incorporating the related

uncertainties is called the likelihood, that is, the conditional probability of the measurements given the unknown parameters. The model for the unknowns that reflects all the uncertainty of the parameters without the information conveyed by the measurements, is called the prior model. BAYESIAN FRAMEWORK

slide-6
SLIDE 6

BAYESIAN FRAMEWORK The formal mechanism to combine the new information (measurements) with the previously available information (prior) is known as the Bayes’ theorem:

( ) ( ) ( ) ( ) ( )

posterior

       x z x x x z z

where posterior(x) is the posterior probability density, (x) is the prior density, (z|x) is the likelihood function and (z) is the marginal probability density of the measurements, which plays the role of a normalizing constant.

slide-7
SLIDE 7

The evolution-observation model is based on the following assumptions : (i) The sequence

k

x

for k = 1, 2, …, is a Markovian process, that is,

1 1 1

( , , , ) ( )

k k k k

 

 

 x x x x x x

(ii) The sequence

k

z

for k = 1, 2, …, is a Markovian process with respect to the history of

k

x

, that is,

1

( , , , ) ( )

k k k k

   z x x x z x

(iii) The sequence

k

x

depends on the past observations only through its own history, that is,

1 1: 1 1

( , ) ( )

k k k k k

 

  

 x x z x x

State Evolution Model: Observation Model:

( , )

k k k k

 z h x n

xk = fk (xk-1, uk-1, vk-1) STATE ESTIMATION PROBLEM

slide-8
SLIDE 8

State Evolution Model: Observation Model:

( , )

k k k k

 z h x n

xk = fk (xk-1, uk-1, vk-1)

  • 1. The prediction problem, concerned with the determination of

1: 1

( )

k k

x z ;

  • 2. The filtering problem, concerned with the determination of

1:

( )

k k

 x z ;

  • 3. The fixed-lag smoothing problem, concerned the determination of

1:

( )

k k p

x z

, where

1 p  is the fixed lag;

  • 4. The whole-domain smoothing problem, concerned with the determination of

1:

( )

k K

 x z , where

1:

{ , 1, , }

K i i

K   z z

is the complete sequence of measurements.

Different problems can be considered: STATE ESTIMATION PROBLEM

slide-9
SLIDE 9

FILTERING PROBLEM By assuming that is available, the posterior probability density is then obtained with Bayesian filters in two steps: prediction and update ( ) ( )    x z x

1:

( )

k k

 x z

slide-10
SLIDE 10

(x0) Prediction (x1) Update (x1 |x0) (z1 |x1) (x1| z1) Prediction (x2 | z1) Update (x2 |x1) (z2 |x2) (x2| z1:2) (x0) Prediction (x1) Update (x1 |x0) (z1 |x1) (x1| z1) Prediction (x2 | z1) Update (x2 |x1) (z2 |x2) (x2| z1:2)

slide-11
SLIDE 11

THE KALMAN FILTER

  • Evolution and observation models are linear.
  • Noises in such models are additive and Gaussian, with

known means and covariances.

  • Optimal solution if these hypotheses hold.

k k k k

  z H x n

State Evolution Model: Observation Model:

  • F and H are known matrices for the linear evolutions of the state x and of the
  • bservation z, respectively.
  • G is matrix that determines how the control u affects the state x.
  • Vector s is assumed to be a known input .
  • Noises v and n have zero means and covariance matrices Q and R, respectively.

x𝑙

− = F𝑙x𝑙−1 + G𝑙u𝑙−1 + s𝑙−1+v𝑙−1

slide-12
SLIDE 12

THE KALMAN FILTER

Prediction:

1 k k k  

 x F x

1 T k k k k k  

  P F P F Q

Update:

 

1 T T k k k k k k k   

  K P H H P H R ( )

k k k k k k  

   x x K z H x

 

k k k k 

 P I - K H P

x𝑙

− = F𝑙x𝑙−1 + G𝑙u𝑙−1 + s𝑙−1+v𝑙−1

slide-13
SLIDE 13

THE PARTICLE FILTER

  • Monte-Carlo techniques are the most general and robust for

non-linear and/or non-Gaussian distributions.

  • The key idea is to represent the required posterior density

function by a set of random samples (particles) with associated weights, and to compute the estimates based on these samples and weights.

  • Introduced in the 50’s, but no much used until recently because
  • f limited computational resources.
  • Particles degenerated very fast in early implementations, i.e.,

most of the particles would have negligible weight. The resampling step has a fundamental role in the advancement of the particle filter.

slide-14
SLIDE 14

Sampling Importance Resampling (SIR) Algorithm

(Ristic, B., Arulampalam, S., Gordon, N., 2004, Beyond the Kalman Filter, Artech House, Boston)

slide-15
SLIDE 15

Sampling Importance Resampling (SIR) Algorithm

(Ristic, B., Arulampalam, S., Gordon, N., 2004, Beyond the Kalman Filter, Artech House, Boston)

slide-16
SLIDE 16

Sampling Importance Resampling (SIR) Algorithm

(Ristic, B., Arulampalam, S., Gordon, N., 2004, Beyond the Kalman Filter, Artech House, Boston)

Although the resampling step reduces the effects of degeneracy, it introduces other practical problems:

  • Limitation in the parallelization.
  • Particles that have high weights are statistically selected many

times: Loss of diversity, known as sample impoverishment, specially if the evolution model errors are small.

slide-17
SLIDE 17

Sampling Importance Resampling (SIR) Algorithm

(Ristic, B., Arulampalam, S., Gordon, N., 2004, Beyond the Kalman Filter, Artech House, Boston)

Step 1 For 1, , i N  draw new particles x

i k from the prior

density

 

1

x xi

k k

and then use the likelihood density to calculate the correspondent weights

 

z x

i i k k k

w   . Step 2 Calculate the total weight

1 N i w k i

T w

  and then normalize the particle weights, that is, for 1, , i N  let

1 i i k w k

w T w

 Step 3 Resample the particles as follows : Construct the cumulative sum of weights (CSW) by computing

1 i i i k

c c w

  for 1, , i N  , with c  . Let 1 i  and draw a starting point

1

u from the uniform

distribution

1

0, U N      For 1, , j N  Move along the CSW by making

 

1 1

1

j

u u N j

   While

j i

u c  make 1 i i   . Assign sample

j i k k

x x  Assign sample

1 j k

w N  

  • Weights are easily evaluated and

importance density easily sampled.

  • Sampling of the importance

density is independent of the measurements at that time. The filter can be sensitive to outliers.

  • Resampling is applied every

iteration, which can result in fast loss of diversity of the particles.

slide-18
SLIDE 18

EXAMPLE: Lumped System

( ) ( ) ( ) d t mq t m t dt h     for t > 0    for t = 0 where ( ) ( ) t T t T 

  T T 

 

h m c L  

slide-19
SLIDE 19

EXAMPLE: Lumped System Two illustrative cases are examined: (i) Heat Flux q(t) = q0 constant and deterministically known; (ii) Heat Flux q(t) = q0 f(t) with unknown time variation.

  • Plate is made of aluminum ( = 2707 kgm-3, c = 896 Jkg-1K-1),

with thickness L = 0.03 m, q0 = 8000 Wm-2, =20 oC, h = 50 Wm-2K-1 and T0 = 50 oC.

  • Measurements of the transient temperature of the slab are

assumed available. These measurements contain additive, uncorrelated, Gaussian errors, with zero mean and a constant standard deviation sz.

  • The errors in the state evolution model are also supposed to be

additive, uncorrelated, Gaussian, with zero mean and a constant standard deviation s.

slide-20
SLIDE 20

(i) Heat Flux q(t) = q0 constant and deterministically known

The analytical solution for this problem is given by: ( ) (1 )

mt mt

q t e e h  

 

   (A.3) The only state variable in this case is the temperature ( )

k k

t    since the applied heat flux q0 is constant and deterministically known, as the other parameters appearing in the formulation. By using a forward finite-differences approximation for the time derivative in equation (A.1.a), we

  • btain:

1

(1 )

k k

mq m t t h         (A.4) Therefore, the state and observation models given by equations (3.a,b) are obtained with:

[ ]

k k

  x [(1 )]

k

m t    F

k

q m t h         s [1]

k 

H

2

[ ]

k 

s  Q

2

[ ]

k z

s  R

(A5.a-f)

k k k k

  z H x n x𝑙

− = F𝑙x𝑙−1 + G𝑙u𝑙−1 + s𝑙−1+v𝑙−1

slide-21
SLIDE 21

(ii) Heat Flux q(t) = q0 f(t) with unknown time variation

The analytical solution for this problem is given by: ( ) ( )

t mt mt t

mq t e e f t dt h  

  

             

(A.6) In this case, the state variables are given by the temperature ( )

k k

t    and the function that gives the time variation of the applied heat flux, that is, ( )

k k

f t f  .As in the case examined above, the applied heat flux q0 is constant and deterministically known, as the other parameters appearing in the formulation. By using a forward finite-differences approximation for the time derivative in equation (A.1.a), we obtain the equation for the evolution of the state variable ( )

k k

t    :

1 1

(1 )

k k k

mq m t t f h   

           (A.7) A random walk model is used for the state variable ( )

k k

f t f  , which is given in the form:

1 1 k k k

f f 

 

  (A.8) where k-1 is Gaussian with zero mean and constant standard deviation srw.

slide-22
SLIDE 22

(ii) Heat Flux q(t) = q0 f(t) with unknown time variation

k k k k

  z H x n x𝑙

− = F𝑙x𝑙−1 + G𝑙u𝑙−1 + s𝑙−1+v𝑙−1

Therefore, the state and observation models given by equations (3.a,b) are obtained with:

k k k

f         x (1 ) 1

k

mq m t t h               F

k

       s (A.9.a-c) [1 0]

k 

H

2 2 k rw 

s s          Q

2

[ ]

k z

s  R

(A.9.d-f)

slide-23
SLIDE 23

Linear Heat Conduction Problem

2 2

1 T T t x      

in 0 < x < L, for t > 0

T 

at x=0, for t > 0

* T T 

at x=L, for t > 0

* T T 

for t=0, in 0 < x < L

Explicit finite-differences:

1 k k  

 T FT S

1 N

T T            Τ (1 2 ) (1 2 ) (1 2 ) (1 2 ) r r r r r r r r r r                      F * rT              S

EXAMPLE: 1D Heat conduction

slide-24
SLIDE 24
  • Concrete with thermal diffusivity  = 4.9x10-7 m2/s
  • Standard-deviation for the measurement errors = 2 oC
  • Final time = 250 seconds
  • Measurements available in the region every 1 second
  • L = 0.1 m
  • N = 50 internal nodes

EXAMPLE: 1D Heat conduction

slide-25
SLIDE 25

50 100 150 200 250 0.05 0.1

  • 20

20 40 60 80 Time, s Exact Temperature x, m T, C 50 100 150 200 250 0.05 0.1

  • 20

20 40 60 80 Time, s Measured Temperature x, m T, C

slide-26
SLIDE 26

50 100 150 200 250 15 20 25 30 35 40 45 50 55 60 Time, s T(x=0.01 m), C Exact Measured Predicted 99% Confidence Interval

Kalman Filter (Standard-deviation for the evolution model errors of 1 oC)

slide-27
SLIDE 27

Kalman Filter (Standard-deviation for the evolution model errors of 0.5 oC)

50 100 150 200 250 15 20 25 30 35 40 45 50 55 60 Time, s T(x=0.01 m), C Exact Measured Predicted 99% Confidence Interval

slide-28
SLIDE 28

Kalman Filter (Standard-deviation for the evolution model errors of 5 oC)

50 100 150 200 250 10 20 30 40 50 60 70 Time, s T(x=0.01 m), C Exact Measured Predicted 99% Confidence Interval

slide-29
SLIDE 29

Particle Filter (Standard-deviation for the evolution model errors of 1 oC)

50 100 150 200 250 15 20 25 30 35 40 45 50 55 60 Time, s T(x=0.01 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

CPU Time Kalman filter: 0.5 s Particle filter: 6.9 s

slide-30
SLIDE 30

Particle Filter (evolution model errors with uniform distribution in [-1,1] oC)

50 100 150 200 250 15 20 25 30 35 40 45 50 55 60 Time, s T(x=0.01 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

slide-31
SLIDE 31

Nonlinear Heat Conduction Problem

( ) ( ) T T C T K T t x x              in 0 < x < L, for t > 0 T x    at x=0, for t > 0 ( ) * T k T q x    at x=L, for t > 0

* T T 

for t=0, in 0 < x < L

 

3

1 2 T A

C T A A e  

 

3

1 2 T B

k T B B e  

Graphite

20 40 60 80 100 120 140 160 180 200 500 1000 1500 2000 2500 1x10

6

2x10

6

3x10

6

4x10

6

5x10

6

6x10

6

Capacidade térmica Capacidade Térmica Volumetrica (J/m

3.°C)

Temperatura (°C) Condutividade térmica Condutividade térmica (W/m.°C)

Thermal conductivity Heat capacity Temperature (oC) Thermal Conductivity (Wm-1oC-1) Volumetric Heat Capacity (Jm-3oC-1)

20 40 60 80 100 120 140 160 180 200 500 1000 1500 2000 2500 1x10

6

2x10

6

3x10

6

4x10

6

5x10

6

6x10

6

Capacidade térmica Capacidade Térmica Volumetrica (J/m

3.°C)

Temperatura (°C) Condutividade térmica Condutividade térmica (W/m.°C)

Thermal conductivity Heat capacity Temperature (oC) Thermal Conductivity (Wm-1oC-1) Volumetric Heat Capacity (Jm-3oC-1)

EXAMPLE: 1D Heat conduction

slide-32
SLIDE 32
  • T* = 20 oC
  • q* = 105 W/m2
  • L = 0.01 m
  • 50 finite-volumes
  • Final time = 90 s
  • Time step = 1 s
  • Errors in the state evolution and observation models were

supposed to be additive, Gaussian, uncorrelated, with zero mean and constant standard-deviations

  • Standard-deviation for the state evolution model = 5 oC
  • Standard-deviation for the observation model = 10 oC

EXAMPLE: 1D Heat conduction

slide-33
SLIDE 33

20 40 60 80 0.005 0.01

  • 100

100 200 300 400 500 600 Time, s Exact Temperature x, m T, C 20 40 60 80 0.005 0.01

  • 100

100 200 300 400 500 600 Time, s Measured Temperature x, m T, C

slide-34
SLIDE 34

10 20 30 40 50 60 70 80 90

  • 100

100 200 300 400 500 600 Time, s T(x=0.0002 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

slide-35
SLIDE 35

10 20 30 40 50 60 70 80 90

  • 100

100 200 300 400 500 600 Time, s T(x=0.001 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

slide-36
SLIDE 36

10 20 30 40 50 60 70 80 90

  • 300
  • 200
  • 100

100 200 300 400 500 600 700 Time, s T(x=0.002 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

slide-37
SLIDE 37

10 20 30 40 50 60 70 80 90

  • 200
  • 100

100 200 300 400 500 600 Time, s T(x=0.004 m), C Particle Filter with Re-sampling - Number of Particles = 50 Exact Measured Predicted 99% Confidence Interval

slide-38
SLIDE 38

APPLICATIONS Estimation of position-dependent transient heat source

38 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

 

2 2 2 2

( , , ) T T T h g x y t C k k T T t e e x y

          

at 0 < x < L, 0 < y < L, for t > 0 T x    at x = 0 and x = L, for t > 0 T y    at y = 0 and y = L, for t > 0 T T  for t = 0, at 0 < x < L and 0 < y < L (4)

slide-39
SLIDE 39

Estimation of position-dependent transient heat source

39 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

APPLICATIONS

slide-40
SLIDE 40

Estimation of position-dependent transient heat source

40 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

1 1 k k k k k T  

  T = T J P v

1 k k k P 

 P = IP v

1 1 2 2

( ) ( ) ( )

k k k k k k k M M

L t T T t t L t T T L t T T t

  

                        J

APPLICATIONS

slide-41
SLIDE 41

Estimation of position-dependent transient heat source

41 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

1 1 k k k k k T  

  T = T J P v

1 k k k P 

 P = IP v

1 2 k k k k M

T T T                T

1 1 1 2 2 2 k k k k k k k k M k M k M

a H G a H G a H G                                                                  P

APPLICATIONS

slide-42
SLIDE 42

Estimation of position-dependent transient heat source

42 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

APPLICATIONS

slide-43
SLIDE 43

Estimation of position-dependent transient heat source

43 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

APPLICATIONS

slide-44
SLIDE 44

Estimation of position-dependent transient heat source

44 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

APPLICATIONS

slide-45
SLIDE 45

Estimation of position-dependent transient heat source

45 H.Massard, F. sepulveda, H.R.B. Orlande and O. Fudym, 2010, Kalman Filtering for thermal diffusivity and transient source term mapping from infrared images, Inverse Problems, Design and Optimization Symposium João Pessoa, Paraíba, Brazil, August 23-27, 2010

APPLICATIONS

slide-46
SLIDE 46

Estimation of the location of the solidification front and the intensity of a line heat sink

46 Silva, W. B. Orlande, H. R. B.; Colaço, M. J., Fudym, O., 2011, Application Of Bayesian Filters To A One- Dimensional Solidification Problem, 21st Brazilian Congress of Mechanical Engineering, Natal, RN, Brazil.

APPLICATIONS

Solid Liquid S(t) Interface

T

Line Heat Sink T (r,t)

slide-47
SLIDE 47

47

APPLICATIONS

The mathematical formulation for the solid phase is given as

   

, , 1 1 ( )

s s s

T r t T r t r in r S t and t r r r t                 

(1.a) while the liquid phase is described as

   

, , 1 1 ( )

l l l

T r t T r t r in S t r and t r r r t                  

(1.b)

 

,

l i

T r t T in r and t   

(1.c)

 

,

l i

T r t T in t and r   

(1.d) At the interface between liquid and solid phases, the following conditions must be satisfied

   

, , ( )

s l m

T r t T r t T in r S t and t    

(1.e)

   

, , ( ) ( )

s l s l

T r t T r t S t k k L in r S t and t r r t            (1.f)

slide-48
SLIDE 48

48

APPLICATIONS

An analytical solution of this problem can be obtained for this physical problem and it is given by [8]:

 

 

2 2

, ( ) 4 4

s s m i i s s

Q r T r t T E E r S t k t                      

(2.a)

   

2 2

, ( ) 4

i m l i i s s i l

T T r T r t T E S t r t E                              

(2.b) where the eigenvalues λ and the solidification front S (t) are given by

 

2 2

2 2

4

s l

l i m s s s i l

k T T Q e e L k E

   

      

 

           (3.a) ( ) 2

s

S t t   

(3.b) In the above equations

i

T is the uniform initial temperature,

m

T is the melting temperature of the material, L is the

latent heat of solidification of the material,  is the density,

s

k and

l

k are the thermal conductivities of the solid and

liquid phases, respectively,

s

 and

l

 are the thermal diffusivities of the solid and liquid phases, respectively, and

s

T and

l

T are temperatures of the solid and liquid phases, respectively.

slide-49
SLIDE 49

Estimation of the location of the solidification front and the intensity of a line heat sink

49 Silva, W. B. Orlande, H. R. B.; Colaço, M. J., Fudym, O., 2011, Application Of Bayesian Filters To A One- Dimensional Solidification Problem, 21st Brazilian Congress of Mechanical Engineering, Natal, RN, Brazil.

APPLICATIONS

The physical problem defined by Eqs. (1.a-f) was solved analytically, where we used the following data, corresponding to solidifying water:

25

i

T C  

,

m

T C   ,

2

0.00118

s

m s  

,

2

0.000146

l

m s   , 2.22

s

w k m c  

, 0.61

l

w k m c   ,

3

997.1kg m  

,

J 80 kg L 

. The line heat sink was supposed to have a constant value equals to W Q = 50 m . In this work, the measurements (for the observation model) were obtained at r=0.01 m. The simulated noisy measurements were uncorrelated, additive, Gaussian, with zero mean and constant standard deviation equal to 5% of the maximum temperature. Figures 3.a,b show the transient measurements obtained after applying such constant line heat sink, with and without errors, respectively.

slide-50
SLIDE 50

APPLICATIONS Estimation of the location of the solidification front and the intensity of a line heat sink

Silva, W. B. Orlande, H. R. B.; Colaço, M. J., Fudym, O., 2011, Application Of Bayesian Filters To A One- Dimensional Solidification Problem, 21st Brazilian Congress of Mechanical Engineering, Natal, RN, Brazil.

slide-51
SLIDE 51

Auxiliary Sampling Importance Resampling (ASIR) Algorithm

(Ristic, B., Arulampalam, S., Gordon, N., 2004, Beyond the Kalman Filter, Artech House, Boston)

Step 1 For i=1,...,N draw new particles xk

i from the prior density

(xk|xi

k-1) and then calculate some characterization of xk,

given xi

k-1, as for example the mean i k=E[xk|xi k-1]. Then

use the likelihood density to calculate the correspondent weights wi

k=(zk|i k)wi k-1

Step 2 Calculate the total weight t=i wi

k and then normalize the

particle weights, that is, for i=1,...,N let wi

k = t-1 wi k

Step 3 Resample the particles as follows : Construct the cumulative sum of weights (CSW) by computing ci=ci-1+wi

k for i=1,...,N, with c0=0

Let i=1and draw a starting point u1 from the uniform distribution U[0,N-1] For j=1,...,N Move along the CSW by making uj=u1+N-1(j-1) While uj>ci make i=i+1 Assign sample xj

k=xi k

Assign sample wj

k=N-1

Assign parent ij=i Step 4 For j=1,...,N draw particles xk

j from the prior density

(xk|xij

k-1), using the parent ij, and then use the likelihood

density to calculate the correspondent weights wj

k=(zk|xj k) / (zk|ij k)

Step 5 Calculate the total weight t=j wj

k and then normalize the

particle weights, that is, for j=1,...,N let wj

k = t-1 wj k

  • The advantage of ASIR over SIR is that

it naturally generates points from the sample at k-1, which, conditioned on the current measurement, are most likely to be close to the true state.

  • The resampling is based on some point

estimate i

k that characterize (xk|xi k-1),

which can be the mean i

k=E[(xk|xi k-1)]

  • r simply a sample of (xk|xi

k-1). If the

state evolution model noise is small, (xk|xi

k-1) is generally well characterized

by i

k, so that the weights wi k are more

even and the ASIR algorithm is less sensitive to outliers than the SIR

  • algorithm. On the other hand, if the state

evolution model noise is large, the single point estimate i

k in the state space may

not characterize well (xk|xi

k-1) and the

ASIR algorithm may not be as effective as the SIR algorithm.

slide-52
SLIDE 52

APPLICATIONS Estimation of the location of the solidification front and the intensity of a line heat sink

Silva, W. B. Orlande, H. R. B.; Colaço, M. J., Fudym, O., 2011, Application Of Bayesian Filters To A One- Dimensional Solidification Problem, 21st Brazilian Congress of Mechanical Engineering, Natal, RN, Brazil.

Bayesian filter Number of Particles (NP) Time RMS error for the solidification front (m) RMS error for the line heat sink intensity (W/m) SIR 100 0.008 min. 9x10-3 1.55 SIR 1000 0.997 min. 2x10-3 1.78 SIR 5000 11.047 min. 1x10-4 0.34 ASIR 100 0.161 min. 7.9x10-5 0.15

slide-53
SLIDE 53

APPLICATIONS Estimation of the location of the solidification front and the intensity of a line heat sink

Silva, W. B. Orlande, H. R. B.; Colaço, M. J., Fudym, O., 2011, Application Of Bayesian Filters To A One- Dimensional Solidification Problem, 21st Brazilian Congress of Mechanical Engineering, Natal, RN, Brazil.

0.00 0.01 0.02 0.03 0.04 0.05 0.06

50 100 150 200 250 300 350 400

TIME (seconds)

LOCATION OF THE SOLIDIFICATION FRONTIER (m)

ESTIMATED FRONTIER ASIR FILTER NP=100 REAL FRONTIER

48 49 50 51 52 53 54 50 100 150 200 250 300 350 400

TIME (seconds) LINE HEAT SINK (W/m) ESTIMATED FLUX ASIR FILTER NP=100 REAL FLUX

slide-54
SLIDE 54

APPLICATIONS Estimation of Unknown Heat Flux in Natural Convection

Colaço, M. J., Orlande, H. R. B.; Silva, W. B., Dulikravich, G., 2011, Application Of A Bayesian Filter To Estimate Unknown Heat Fluxes In A Natural Convection Problem, ASME 2011 International Design Engineering Technical Conferences & Computers and Information in Engineering Conference, IDETC/CIE 2011, August 29-31, Washington, DC.

PHYSICAL PROBLEM The physical problem under picture in this paper involves the transient laminar natural convection

  • f a fluid inside a two-dimensional square cavity.

The fluid is initially at rest and at the uniform temperature Tc. At time zero, the bottom and top surfaces are subjected to time-dependent heat fluxes q1(t) and q2(t), respectively. The left and right surfaces are subjected to constant temperatures Tc and Th, respectively. The fluid properties are assumed constant, except for the density in the buoyancy term, where we consider Boussinesq’s approximation valid.

slide-55
SLIDE 55

APPLICATIONS Estimation of Unknown Heat Flux in Natural Convection

Colaço, M. J., Orlande, H. R. B.; Silva, W. B., Dulikravich, G., 2011, Application Of A Bayesian Filter To Estimate Unknown Heat Fluxes In A Natural Convection Problem, ASME 2011 International Design Engineering Technical Conferences & Computers and Information in Engineering Conference, IDETC/CIE 2011, August 29-31, Washington, DC.

slide-56
SLIDE 56

APPLICATIONS Estimation of Unknown Heat Flux in Natural Convection

Colaço, M. J., Orlande, H. R. B.; Silva, W. B., Dulikravich, G., 2011, Application Of A Bayesian Filter To Estimate Unknown Heat Fluxes In A Natural Convection Problem, ASME 2011 International Design Engineering Technical Conferences & Computers and Information in Engineering Conference, IDETC/CIE 2011, August 29-31, Washington, DC.

Exact Estimated

slide-57
SLIDE 57

APPLICATIONS Estimation of Unknown Heat Flux in Natural Convection

Colaço, M. J., Orlande, H. R. B.; Silva, W. B., Dulikravich, G., 2011, Application Of A Bayesian Filter To Estimate Unknown Heat Fluxes In A Natural Convection Problem, ASME 2011 International Design Engineering Technical Conferences & Computers and Information in Engineering Conference, IDETC/CIE 2011, August 29-31, Washington, DC.

Exact Estimated

slide-58
SLIDE 58

APPLICATIONS Estimation and Control of the Temperature Field in Oil Pipelines

Vianna, F., Orlande, H. R. B., Dulikravich, G., 2010, Optimal Heating Control To Prevent Solid Deposits In Pipelines, V European Conference On Computational Fluid Dynamics - ECCOMAS CFD 2010

X-tree Flowline

(1) (2)

Measurement Point (1) (2) R=0 R=1

The physical problem examined in this work is based on a critical operational condition involving a pipeline shutdown situation, where the produced fluid is assumed stagnant. An optimal control approach was used to drive the predicted temperatures above a reference level.

 

 Q Heat flux on boundary surface

slide-59
SLIDE 59

APPLICATIONS The dimensionless mathematical formulation for this one- dimensional unsteady heat diffusion problem is given by The fluid was considered as homogeneous, isotropic and with constant thermophysical properties.

     

R R R R R R                , 1 , ,

2 2

1, R    

     

     Q R Bi R R     , , 1, R   

1 ,  R  1, R    

is the dimensionless temperature distribution into the medium.

  , R

slide-60
SLIDE 60
  • Control strategy is in accordance with the optimum control theory for linear

problems.

  • The aim of the associated optimal control problem is to find the control inputs uk

(heat flux on boundary surface) that minimizes the difference between the fluid temperature field and a desired profile . Where u* and x* refer to the steady values of the control input and state variables

APPLICATIONS

slide-61
SLIDE 61

In terms of the linear quadratic regulator problem, the optimal values of the control input are obtained by minimizing the following quadratic cost functional: where the weighting matrices Q and R are symmetric positive definite. The solution to the optimal control problem is the state feedback control law: where the discrete-time state feedback gain K is of the form

APPLICATIONS

slide-62
SLIDE 62

Matrix S is the steady state solution to the discrete-time Riccati equation: Thus, the control input can be calculated from the control law as However, when state variables are not directly available for control, an observer (KALMAN FILTER OR PARTICLE FILTER ) was built to estimate the state variables from the input and output variables of the system.

APPLICATIONS

slide-63
SLIDE 63

Simulated measurements with standard deviation of 3 oC in the

  • bservation model error

APPLICATIONS

slide-64
SLIDE 64

Evolution model with standard deviation of 0.01 oC KALMAN FILTER

slide-65
SLIDE 65

Temperature Predicted with the Kalman Filter and the Exact Temperature

Evolution model with standard deviation of 0.01 oC

 = 0.69

slide-66
SLIDE 66

Temperature Predicted with the Kalman Filter and the Exact Temperature

Evolution model with standard deviation of 0.01 oC

 = 1.1

slide-67
SLIDE 67

Evolution model with standard deviation of 1 OC KALMAN FILTER

slide-68
SLIDE 68

Temperature Predicted with the Kalman Filter and the Exact Temperature

Evolution model with standard deviation of 1 oC

 = 0.61

slide-69
SLIDE 69

Temperature Predicted with the Kalman Filter and the Exact Temperature

Evolution model with standard deviation of 1 oC

 = 0.8

slide-70
SLIDE 70

Number of Particles: 200 Evolution model with standard deviation of 1 oC PARTICLE FILTER

slide-71
SLIDE 71

Temperature Predicted with the Particle Filter and the Exact Temperature

Evolution model with standard deviation of 1 oC

 = 0.69

slide-72
SLIDE 72

Evolution model with standard deviation of 1 oC

 = 1.2 Temperature Predicted with the Particle Filter and the Exact Temperature

slide-73
SLIDE 73

Simulated measurements with standard deviation of 3 oC in the

  • bservation model error

Two-dimensional case with Particle Filter Observer

Standard deviation of the evolution model error of 1oC. Standard deviation of the evolution model error of 3oC.

(1) (2) (3)

slide-74
SLIDE 74

CONCLUSIONS

  • Kalman filter provides optimal solutions for linear-Gaussian

evolution-observation models.

  • Particle filter is the most general and robust technique for non-

linear models and/or non-Gaussian distributions.

  • ASIR algorithm is faster and requires less particles than the SIR

algorithm.

slide-75
SLIDE 75

ACKNOWLEDGEMENTS

  • Professors Denis Maillet, Philippe Le Masson and Yann

Favennec.

  • Conselho Nacional de Desenvolvimento Científico e

Tecnológico (CNPq).

  • Coordenação de Aperfeiçoamento de Pessoal de Nível

Superior (CAPES).

  • Fundação Carlos Chagas Filho de Amparo à Pesquisa do

Estado do Rio de Janeiro (FAPERJ).

  • METTI sponsors.