Stochastic approach to electric process tomography

(1)

Stochastic approach to electric process tomography

Hanna Katriina Pikkarainen Institute of Mathematics Helsinki University of Technology

Espoo, Finland hanna.pikkarainen@hut.fi Abstract

We consider the process tomography problem of imaging the concentration distribution of a given substance in a fluid moving in a pipeline based on electromagnetic measurements on the surface of the pipe.

We view the problem as a state estimation problem. The concentration distribution is treated as a stochastic process satisfying a stochastic differential equation referred to as the state evolution equation. The measurements are described in terms of an observation equation containing the measurement noise. Our main interest is in the mathematical formulation of the state evolution equation. The time evolution is modelled by a stochastic convection-diffusion equation.

We derive a discrete evolution equation for the concentration distribution by using the stochastic integration theory and the semigroup technique.

Nomenclature

argλ the argument of a complex number λ

A^∗ the adjoint of an operatorA B(E, H) the space of bounded linear

operators fromE intoH B(E) the Borel σ-algebra of a to-

pological space E

C([0, T];E) the space of continuous functions from [0, T] intoE D¯ the closure of a set D

∂D the boundary of a setD D(A) the domain of an operatorA H²(D) the space of functions fromD

into C with square integrable weak derivatives up to order 2 Ker(A) the null space of an operator

A

L¹(0, T;E) the space of Bochner integrable functions from [0, T] into E

L²(D) the space of square integrable functions fromD intoC N(m, Q) the Gaussian measure with

mean m and covariance operator Q

R(λ, A) the resolvent operator of an operator A

Tr(A) the trace of an operatorA U C(D) the space of uniformly con-

tinuous functions fromDinto Rⁿ

U C¹(D) the space of uniformly con- tinuously differentiable functions from DintoRⁿ

ν the outer unit normal vector ρ(A) the resolvent set of an oper-

ator A

(2)

Introduction

We consider the process tomography problem of imaging the concentration distribution of a given substance in a fluid moving in a pipeline based on electromagnetic measurements on the surface of the pipe.

In electrical impedance tomography (EIT), electric currents are applied to electrodes on the surface of an object and the result-

ing voltages are measured using the same electrodes (Figure 1). The conductivity distribution inside the object is reconstructed based on the voltage measurements. The relation between the conductivity and the concentration depends on the process and it is usually non-linear. At least for multiphase mixtures and strong electrolytes, such rela- tions are studied and discussed in the liter- ature.

I U

c(x,t) I

v

Figure 1: EIT in process tomography

In traditional EIT, it is assumed that the object remains stationary during the measurement process. A complete set of measurements, also called a frame, consists of all possible linearly independent injected cur- rent patterns and the corresponding set of voltage measurements. In process tomography, in general we cannot assume that the target remains unaltered during a full set of measurements. Thus conventional re- construction methods cannot be used. The time evolution needs to be modelled prop- erly. We view the problem as a state estimation problem. The concentration distribution is treated as a stochastic process, or a state of the system, that satisfies a stochastic differential equation referred to as the state evolution equation. The measurements are described in terms of an observation equation containing the measurement noise.

Often in a state estimation approach the time variable is assumed to be discrete and the space variable to be finite dimensional.

It is convenient from the practical point

of view. Observations are usually done at discrete times and the compution requires space discretization. Since our interest is in the mathematical formulation of the problem, we assume that the space variable is infinite dimensional. The solution of the state estimation problem is a function valued random variable instead of anRⁿ-valued Gaus- sian distribution.

Our goal is to have a real-time monitoring for a flow in a pipeline. For that reason the computational time has to be minim- ized. Therefore, we use a simple model, the convection-diffusion equation, for the flow.

It is easy to implement and fast to compute.

Since we cannot be sure that other features such as turbulence of the flow do not appear, we use stochastic modelling. Therefore the randomness is due to the lack of informa- tion, not to the intrinsic randomness of the concentration.

The measurements are done in a part of the boundary of the pipe. We get enough in- formation for an accurate computation only

(3)

from a segment of the pipe. It would be nat- ural to choose the domain of the model to be the segment of the pipe. If the domain is re- stricted to be a segment of the pipe, we have to use some boundary conditions in the in- put and the output end of the segment. The choice of boundary conditions has an effect on the solution. The most commonly used boundary conditions do not represent the actual circumstances in the pipe. Therefore, we do not do the domain restriction. We assume that the pipe is infinitely long. With the assumption we derive the state evolution model. The concentration distribution, which we are actually interested in, is the restriction of the solution to the evolution model to a segment of the pipe.

This problem has been considered in the articles [Sep01a, Sep01b, Sep01c] and in the proceeding papers [Sep01d, Sep01e, Sep02, Sep03, Ruu03]. Those articles and proceeding papers concentrate on the numerical im- plementation of the problem. Our main interest is in the mathematical formulation of the state evolution equation and the estimation problem in general. We refer to those articles and references in them con- cerning the observation equation. The rigor- ous knowledge of the stochastic nature of the state evolution equation is essential for solv- ing the electric process tomography problem.

Mathematical preliminaries Analytic semigroup

LetE be a separable Banach space. A family {U(t)}t≥0 of bounded linear operators from E intoE is called asemigroup if

(i) U(t)U(s) =U(t+s) for alls, t≥0 and (ii) U(0) =I.

A semigroup U(t) is strongly continuous, if for all x∈E the functiont7→U(t)x is continuous in the interval [0,∞). It is analytic, if the function t7→U(t) can be extended to

be an analytic function from a sector {λ∈C:λ6= 0, |argλ| ≤β}

with some β ∈ (0, π) to the space of all bounded linear operators fromE intoE.

The linear operatorA defined by D(A) :=

½

x∈E :∃ lim

t→0⁺

U(t)x−x t

¾

and for allx∈ D(A) Ax:= lim

t→0⁺

U(t)x−x t

is the infinitesimal generator of the semig- roupU(t).

Definition 1. A linear operatorA:D(A)⊆ E → E is sectorial, if there exist constants ω∈R, θ∈(π/2, π) and M >0 such that

(i) S_θ,ω:={λ∈C:λ6=ω, |arg(λ−ω)|<

θ} ⊂ρ(A) and

(ii) kR(λ, A)k_B(E) ≤ _|λ−ω|^M for all λ ∈ Sθ,ω.

Theorem 2. [Lun95, Section 2.1] Let A : D(A)⊆E →E be a sectorial operator. The operator family {U(t)}t≥0 defined by

U(0)x:=x and

U(t)x:= 1 2πi

Z

ω+γr,η

e^tλR(λ, A)x dλ for all t > 0 and x ∈ E, where r > 0, η ∈ (π/2, θ) and γr,η is a curve {λ ∈ C :

|argλ| = η, |λ| ≥ r} ∪ {λ ∈ C : |argλ| ≤ η, |λ|=r} oriented counterclockwise, is the analytic semigroup generated by the operator A. The semigroup{U(t)}_t≥0is strongly continuous if and only if the domain D(A) is dense in E.

Analytic semigroups can be used to solve initial value problems.

(4)

Theorem 3. [Lun95, Section 2.1] Let A : D(A) ⊆ E → E be a sectorial operator.

If U(t) is the analytic semigroup generated by the operator A and u₀ ∈ D(A), then the solution to the initial value problem

(u⁰(t) =Au(t), u(0) =u₀ is

u(t) =U(t)u₀ for all t >0.

Theorem 4. [Paz83, Sections 4.2–4.3] Let A : D(A) ⊆E → E be a sectorial operator with dense domain. Let U(t) be the analytic semigroup generated by the operatorA, u0∈E and f ∈L¹(0, T;E). Then the non- homogeneous initial value problem

(u⁰(t) =Au(t) +f(t), u(0) =u₀

has a unique weak solution given by the formula

u(t) =U(t)u₀+ Z t

0

U(t−s)f(s)ds for all0≤t≤T. The weak solution belongs to the space C([0, T];E).

Stochastic analysis in infinite dimensions

Let (Ω,F,P) be a probability space and (E,G) a measurable space. A function X : Ω → E such that the set {ω ∈ Ω : X(ω) ∈ A} belongs toF for each A∈ G is called a random variable from (Ω,F,P) into (E,G).

A random variable is calledsimpleif it takes only a finite number of values. A random variable X is Gaussian, if its distribution L(X) is a Gaussian measure on (E,G).

Let E be a separable Banach space. An E-valued random variable X is said to be Bochner integrable if

Z

Ω

kX(ω)kP(dω)<∞.

If X is Bochner integrable, the integral R

ΩX dP can be defined and is denoted by E(X). If

Z

Ω

kX(ω)k^pP(dω)<∞

for somep≥1, it is said that X belongs to the spaceL^p(Ω,F,P;E).

IfHis a separable Hilbert space andX, Y ∈ L²(Ω,F,P;H), the covariance operator of X and the correlation operator of X and Y are defined by

Cov(X)h:=E[(h, X−EX)(X−EX)]

and

Cor(X, Y)h:=E[(h, Y −EY)(X−EX)]

for allh∈H.

Proposition 5. [PZ92, Proposition 1.10]

Let E be a separable Banach space, X a Bochner integrable E-valued random variable defined on (Ω,F,P) and G a σ-algebra contained in F. There exists a unique, up to a set ofP-probability zero, Bochner integrable E-valued random variable Z, measurable with respect toG, such that

Z

A

X dP= Z

A

Z dP for allA∈ G.

The random variable Z is denoted by E(X|G) and called the conditional expect- ation of X given G. We use the notation E(X|Y) := E(X|σ(Y)) where σ(Y) is the σ-algebra generated by the random variable Y.

Let {Fi}i∈I be a family of sub-σ-algebras of F. These σ-algebras are said to be independent, if for every finite subset J ⊂ I and every family{Ai}i∈J such thatAi∈ Fi, i∈J,

P Ã

\

i∈J

Ai

!

=Y

i∈J

P(Ai).

(5)

Random variables {Xi}_i∈I areindependent, if the σ-algebras {σ(Xi)}i∈I are independent.

A family X = {X(t)}_t∈[0,T_] of E-valued random variables defined on Ω is called a stochastic process. We set X(t, ω) :=

X(t)(ω) for all t ∈ [0, T] and ω ∈ Ω.

Functions X(·, ω) are called the trajectories ofX. A stochastic process iscontinuous if its trajectories are continuous P-almost surely. A stochastic process X is said to be Gaussian, if for any n ∈ N and for all t₁, t₂, . . . , tn ∈ [0, T] the Eⁿ-valued random variable (X(t₁), X(t₂), . . . , X(tn)) is Gaus- sian.

Let E be a separable Hilbert space and Q a bounded linear symmetric non-negative trace class operator with KerQ={0}.

Definition 6. A E-valued stochastic process W(t), t ∈ [0, T], is called a Q-Wiener process if

(i) W(0) = 0, (ii) W is continuous,

(iii) W has independent increments, i.e., W(u) −W(t) and W(s) −W(r) are independent if 0 ≤r < s ≤t < u≤T and

(iv) L(W(t)−W(s)) =N(0,(t−s)Q) for 0≤s≤t≤T.

Let T > 0 be fixed. A family of σ-algebras {Ft}_t∈[0,T]is called afiltration, ifFs⊆ Ft⊆ F for all s, t ∈ [0, T], s ≤ t. We denote by F_t⁺ the intersection of allFs where t < s≤T, i.e.,

F_t⁺ := \

t<s≤T

Fs.

Definition 7. The filtration {Ft}_t∈[0,T_] is said to be normal if

(i) F0 contains all A ∈ F such that P(A) = 0,

(ii) Ft=F_t⁺ for all t∈[0, T].

W is a Q-Wiener process with respect to a filtration{Ft}_t∈[0,T_] if

(i) W(t) isFt-measurable,

(ii) W(t+h)−W(t) is independent of Ft

for allh≥0 andt, t+h∈[0, T].

We denote by PT the σ-algebra generated by sets of the form

((s, t]×F, 0≤s < t≤T, F ∈ Fs, {0} ×F, F ∈ F0.

A measurable function from ([0, T]×Ω,PT) into (E,B(E)) is called apredictable process.

A B(E, H)-valued process Φ(t), t ∈ [0, T], taking only a finite number of values is said to be elementary if there exist a sequence 0 = t0 < t1 < . . . < tk = T and a sequence Φ₀,Φ₁, . . . ,Φ_k−1 of B(E, H)-valued simple random variables such that Φm is Ftm-measurable and

Φ(t) = Φ_m fort∈(t_m, t_m+1] for all m = 0,1, . . . , k −1. We define the stochastic integral for elementary processes Φ by the formula

Z t

0

Φ(s)dW(s) :=

k−1

X

m=0

Φm(W(t_m+1∧t)−W(tm∧t)) where s∧t = min(s, t). The definition of the stochastic integral can be extended to all B(E, H)-valued predictable processes Φ such that

E Z T

0

Tr [Φ(s)QΦ^∗(s)] ds <∞.

Let{Ft}_t∈[0,T_]be a normal filtration andW aQ-Wiener process with respect to{Ft}.

Theorem 8. [PZ92, Chapter 5]Let U(t) be an analytic semigroup generated by an operator A with dense domain D(A) ⊂ H, f aH-valued predictable process with Bochner

(6)

integrable trajectories on[0, T],B a bounded linear operator fromE intoH andX₀ aF0- measurable random variable. If

Z T

0

Tr [U(s)BQB^∗U^∗(s)] ds <∞, then the stochastic differential equation (dX(t) = [AX(t) +f(t)] +BdW(t),

X(0) =X₀

has exactly one weak solution given by the formula

X(t) =U(t)X0+ Z t

0

U(t−s)f(s)ds+

+ Z t

0

U(t−s)B dW(s) for all 0≤t≤T.

State estimation

Let D ⊂ Rⁿ be a domain that corresponds to the object of interest. We denote by X = X(t, x), x ∈ D, a distributed parameter describing the state of the object – the unknown distribution of a physical target – at time t ≥ 0. We assume that we have a model for the time evolution of the parameter X. We assume that instead of a deterministic functionX is a stochastic process satisfying a stochastic differential equation. This allows us to incorporate phenom- ena such as modelling uncertainties into the model.

Let Y =Y(t) denote a quantity that is dir- ectly observable at times t ∈ I, I = {tk : tk < tk+1, k∈N}. We assume that the de- pendence ofY upon the stateXis known up to observation noise and modelling errors.

The state estimation system consists of a pair of equations

dX(t) =F(t, X, R) +dW(t), (1) Y(t) =G(t, X, S), t∈I. (2) Equation (1) is called the state evolution equation and is to be interpreted as a

stochastic differential equation in which the function F is the evolution model function andR=R(t) andW =W(t) are stochastic processes. The processW is called thestate noise. Equation (2) is called theobservation equation. The function Gis the observation model function andS =S(t) is a stochastic process, theobservation noise.

The state estimation problem can be formu- lated as follows: Estimate the state X satisfying an evolution equation of the type (1) based on the observed values ofY(t), when t is in a given subset of I. Estimators of the stateXare calculated by taking conditional expectations with respect to the measurements. The most commonly used estimators are the predictor

E(X(tk)|Y(ti), i= 1, . . . , k−1), (3) the filter

E(X(tk)|Y(ti), i= 1, . . . , k) (4) and the smoother

E(X(tk)|Y(ti), i∈I). (5) Predictor (3) is based on the history at the previous time step, Filter (4) on the cur- rent history and Smoother (5) on the whole measurement set.

We consider the special case in which observations are obtained by EIT measurements and in which the physical target can be described by the convection-diffusion equation.

Mathematical formulation of the state evolution equation

We examine a concentration distribution in a fluid moving in a pipe with a velocity distribution defined by the laminar flow equation by doing electric measurements at the boundary of the pipe. Let κ =κ(x) be the diffusion coefficient andv=v(x) be the velocity of the flow. The diffusion coefficient and the velocity distribution are assumed to

(7)

be known and stationary. In addition, the incompresessibility condition ∇ ·v = 0 is valid in Dand the flow is tangential at the pipe walls, i.e.,ν·v= 0 at∂D. We assume that the concentration distributionC(t) is a stochastic process satisfying the stochastic differential equation

dC(t) = [LC(t) +f(t)]dt+BdW(t) (6) for everyt >0 with the initial value

C(0) =C₀. (7) The operator L is the deterministic convection-diffusion operator

L:D(L)→L²(D)

c7→ ∇ ·(κ∇c)−v· ∇c with the domain

D(L) =

½

c∈H²(D) : ∂c

∂ν



∂D

= 0

¾

where D = {x = (x₁, x⁰) ∈ Rⁿ : |x⁰| < R}

is an infinitely long pipe with R > 0. The boundary condition at the boundary of the pipe is included in the domain of the operator L. We assume that there is no diffusion through the pipe walls. We model with f(t) a possible control of the system. Since the control term is known, if the state of the system is known, we may assume that f(t) is anL²(D)-valued predictable process.

The term BdW(t) is a source term repres- enting possible modelling errors, where B is a bounded linear operator from L²(D) to itself and W(t) is an L²(D)-valued Wiener process.

We use the semigroup technique to solve the stochastic convection-diffusion equation (6) with the initial value (7).

Theorem 9. [Pik, Chapter 5] The operator L is sectorial, if the diffusion coefficient κ is positive and bounded from below, κ(x) ≥δ > 0 for all x ∈ D, and the diffu-¯ sion coefficient and the velocity of the flow satisfy the conditions

(κ: ¯D→R, κ∈U C¹( ¯D), v: ¯D→Rⁿ, v∈U C( ¯D).

Hence under these requirements L gener- ates a strongly continuous analytic semig- roupU(t).

Consequently, if c₀ ∈ L²(D) and f ∈ L¹(0, T;L²(D)), by Theorem 4 the deterministic version

(c⁰(t) =Lc(t) +f(t) c(0) =c₀

of the stochastic convection-diffusion equation (6) with the initial value (7) has a unique weak solution

c(t) =U(t)c0+ Z t

0

U(t−s)f(s)ds for allt∈[0, T].

By Theorem 8 the weak solution to the stochastic convection-diffusion equation (6) with the initial value (7) has the similar form as the deterministic one.

Theorem 10. [Pik, Chapter 5] Let {Ft}, t ∈ [0, T], be a normal filtration and W a Wiener process with respect to the filtration. If f has Bochner integrable trajectories on [0, T] and C₀ is F0-measurable, then the stochastic convection-diffusion equation (6) with the initial value (7) has exactly one weak solution given by the formula

C(t) =U(t)C₀+ Z t

0

U(t−s)f(s)ds+

+ Z t

0

U(t−s)B dW(s) for allt∈[0, T].

Discrete evolution model without control

We assume that there is no control in our system, i.e.,f ≡0. Then the solution of the stochastic convection-diffusion equation (6) with the initial value (7) is

C(t) =U(t)C₀+ Z t

0

U(t−s)B dW(s)

(8)

for all t∈[0, T]. If the initial value C₀ is a Gaussian random variable, then the concentration distributionC is a Gaussian process.

Since C is a continuous process, it can also be interpreted as a C([0, T];L²(D))-valued random variable. Hence the realizations of the concentration distribution belong to the same space as the solutions of the deterministic case.

We assume that the measurements are done in a discrete set of timest_k. We use the nota- tionC_k =C(t_k) and ∆_k=t_k+1−tk. Then a discrete evolution model for the concentration distribution is

C_k+1 =U(∆_k)C_k+W_k (8) where

Wk = Z tk+1

t_k

U(t_k+1−s)B dW(s).

The term W_k can be seen as a state noise.

The state noise Wk is a Gaussian random variable with mean 0 and covariance operator given by the integral

Z tk+1

t_k

U(t_k+1−s)BQB^∗U^∗(t_k+1−s)ds.

Furthermore, the state noises at different time steps are uncorrelated, i.e.,

Cor(Wk, Wl) = 0 for all k6=l.

The discrete evolution model (8) combined with an observation model enables us to cal- culate some estimator for the concentration distribution.

The discrete evolution model used in the numerical implementations in the papers men- tioned before is of the same form as Equa- tion (8), but there is an extra term coming from the boundary condition used in the in- put end of the computational domain. The matrix appearing in the place of the operator U(∆_k) and the state noise of the finite dimensional evolution equation depend on the discretization scheme. To develop

discretization invariant estimation methods, [Las02], it is important to study continuous stochastic models. This is one of the main motivations of this work. Knowledge of the properties of the operatorU(∆_k) and the state noiseW_k is essential for the imple- mentation of the problem. Further results will be published in [Pik].

References

[Las02] Sari Lasanen. Discretizations of Generalized Random Variables with Applications to Inverse Prob- lems. PhD thesis, University of Oulu, Oulu, Finland, 2002.

[Lun95] Alessandra Lunardi. Analytic Semigroups and Optimal Reg- ularity in Parabolic Problems, volume 16 ofProgress in Nonlinear Differential Equations and Their Applications. Birkh¨auser, Basel, Switzerland, 1995. ISBN 3-7643- 5172-1.

[Paz83] Amnon Pazy. Semigroups of Lin- ear Operators and Applications to Partial Differential Equations, volume 44 of Applied Mathem- atical Sciences. Springer-Verlag, New York, U.S.A., 1983. ISBN 0- 387-90845-5.

[Pik] Hanna Katriina Pikkarainen. A Mathematical Model for Electrical Impedance Process Tomography.

PhD thesis, Helsinki University of Technology, Espoo, Finland. (in preparation).

[PZ92] Giuseppe Da Prato and Jerzy Za- bczyk.Stochastic Equations in In- finite Dimensions, volume 44 of Encyclopedia of Mathematics and its Applications. Cambridge Uni- versity Press, Cambridge, United Kingdom, 1992. ISBN 0-521- 38529-6.

(9)

[Ruu03] A.R. Ruuskanen, A. Sepp¨anen, S. Duncan, E. Somersalo, and J.P.

Kaipio. Optimal control in process tomography. InProceedings of 3rd World Congress on Industrial Pro- cess Tomography, pages 245–251.

The Virtual Centre for Industrial Process Tomography, 2003. ISBN 0-853-16240-9.

[Sep01a] A. Sepp¨anen, M. Vauhkonen, P.J.

Vauhkonen, E. Somersalo, and J.P. Kaipio. State estimation with fluid dynamical evolution models in process tomography – an ap- plication with impedance tomography. Inverse Problems, 17:467–

484, 2001.

[Sep01b] A. Sepp¨anen, M. Vauhkonen, E. Somersalo, and J.P. Kaipio.

State space models in process tomography – approximation of state noise covariance. Inverse Problems in Engineering, 9:561–

585, 2001.

[Sep01c] A. Sepp¨anen, M. Vauhkonen, P.J.

Vauhkonen, E. Somersalo, and J.P. Kaipio. Fluid dynamical models and state estimation in process tomography: Effect due to inaccuracies in flow fields. Journal of Electronic Imaging, 10(3):630–

640, 2001.

[Sep01d] A. Sepp¨anen, M. Vauhkonen, E. Somersalo, and J.P. Kaipio. Ef-

fects of inaccuracies in fluid dynamical models in state estimation of process tomography. In Pro- cess Imaging for Automatic Con- trol, volume 4188 ofProceedings of SPIE, pages 69–80, 2001.

[Sep01e] A. Sepp¨anen, M. Vauhkonen, P.J.

Vauhkonen, E. Somersalo, and J.P. Kaipio. State estimation in three dimensional impedance imaging – use of fluid dynamical evolution models. In Proceedings of 2nd World Congress on Industrial Process Tomography, pages 198–

207. The Virtual Centre for Indus- trial Process Tomography, 2001.

[Sep02] A. Sepp¨anen, M. Vauhkonen, E. Somersalo, and J.P. Kaipio. In- ference of velocity fields based on tomographic measurements in process industry. In4th International Conference on Inverse Problems in Engineering, 2002.

[Sep03] A. Sepp¨anen, L. Heikkinen, T. Sa- volainen, E. Somersalo, and J.P.

Kaipio. An experimental eval- uation of state estimation with fluid dynamical models in process tomography. InProceedings of 3rd World Congress on Industrial Pro- cess Tomography, pages 541–546.

The Virtual Centre for Industrial Process Tomography, 2003. ISBN 0-853-16240-9.