## Abstract

Multiphase flow imaging is a very challenging and critical topic in industrial process tomography. In this article, simulation and experimental results of reconstructing the permittivity profile of multiphase material from data collected in electrical capacitance tomography (ECT) are presented. A multiphase narrowband level set algorithm is developed to reconstruct the interfaces between three- or four-phase permittivity values. The level set algorithm is capable of imaging multiphase permittivity by using one set of ECT measurement data, so-called absolute value ECT reconstruction, and this is tested with high-contrast and low-contrast multiphase data. Simulation and experimental results showed the superiority of this algorithm over classical pixel-based image reconstruction methods. The multiphase level set algorithm and absolute ECT reconstruction are presented for the first time, to the best of our knowledge, in this paper and critically evaluated.

This article is part of the themed issue ‘Supersensing through industrial process tomography’.

## 1. Introduction

The multiphase flows that appear in many industries (i.e. oil and gas, aerospace, chemicals, etc.) are considered complex and unpredictable. Sometimes, a multiphase flow is associated with undesirable deposits and dirt (which are considered as another phase) that disturb the system and cause instrumentation faults. In order to ensure safety and increase profits, the multiphase flow must be adequately controlled [1]. In order to do so, analysis of the multiphase flow is required to determine some of the flow features/properties, such as the flow regime, number of phases, concentration of each phase, etc. Currently, multiphase flows are analysed by means of (i) multiphase flow meters (MPFMs) to get the flow rate of each phase or (ii) process tomography to obtain an image of the cross section of the actual pipe being monitored. Each method has its limitations and is being constantly enhanced and modified. A comprehensive review of the different process tomography techniques was conducted by [2].

Electrical capacitance tomography (ECT) is one method of process tomography used in industrial process monitoring for imaging the cross section of a pipe/vessel. It measures the external capacitance of the enclosed objects to determine the internal permittivity distribution, which is then used to reconstruct an image of the process. ECT was first developed in the late 1980s to image a two-phase flow [3]. Since then, ECT systems have been used successfully in numerous research investigations for industrial multiphase processes, including gas/solid distribution in pneumatic conveyors [4], fluidized beds [5–7], flame combustion [8–10], gas/liquid flows [11], the water/oil/gas separation process [12], water hammer [13], determining the characteristics of the molten metal in the lost foam casting process [14], deposit detection [15] and many others. This paper presents a new application area for ECT as a novel imaging technique and measurement tool for multiphase flows in pipelines. The method presented here demonstrates the ability and potential of ECT imaging for multiphase flow imaging and measurement using absolute ECT data.

Most studies for process ECT employ classical iterative reconstruction techniques, where the inverse problem is solved iteratively by updating the estimated model, thus improving the match between the measured and the calculated capacitance. However, such iterative techniques require stabilization, and thus a least-squares functional is usually augmented by some additional regularization terms. Such regularization may lead to ‘over-smoothing’ the reconstructed images, which makes the characterization of all the different phases in terms of size and contrast a major challenge, especially when one of the phases has low concentration or has relatively higher permittivity. The interaction between the capacitance measurements and the process media is nonlinearly influenced by the inhomogeneity of the flow. Accordingly, applying standard reconstruction that is based on regularization tools usually only results in low-resolution images, especially in high-contrast cases where one phase has higher permittivity value compared with the other phases. This makes classical reconstruction techniques more challenging in terms of correctly distinguishing the other low-permittivity phases, because the high-permittivity phase will dominate and make the problem highly nonlinear. Thus, so far, the ECT system has mainly been used to image such low-contrast non-conductive media, not taking into consideration the more complicated scenarios of multiphase media. However, some conductive media, such as water-continuous flow in oil pipelines, are being investigated recently under different conditions.

In this paper, a multiphase shape reconstruction based on the level set method for ECT is presented, where the multiphase framework for image segmentation developed by [16] and the narrowband level set formulation developed by [17] are modified and employed to tackle this nonlinear problem. Both high- and low-contrast cases of a multiphase reconstruction are studied. In this approach, the permittivity values of all the different phases are approximated and assumed to be known, whereas their size, shape and location are recovered from the data. The aim of our reconstruction is to correctly determine the location of the present phases and their concentration. This is done by estimating the size and contrast of every phase without an emphasis on the exact shape of the phases but rather an estimated one. The multiphase level set shape reconstruction results shown in this paper are some of the first ones using experimental data of ECT. Such information has the potential to be used for MPFMs.

In this paper, §2 introduces the ECT system, §3 describes the image reconstruction method used and developed, whereas §4 shows some simulated and experimental results of the low- and high-contrast cases. Section 5 gives an analysis of the results and explains the superiority of our method over classical pixel-based reconstruction methods.

## 2. Electrical capacitance tomography system

Figure 1 shows a typical 12-electrode ECT system, which comprises mainly three subsystems: the capacitance sensor, the data acquisition unit and the computer unit. The capacitance sensor consists of more than one electrode, an external shield and axial and radial guards as seen in figure 1*a*. Typically, several electrodes are mounted equidistantly on the periphery of the non-conductive object to be imaged (i.e. process pipe or vessel). The adjacent electrodes are separated from each other by a small piece of copper called the axial guard. In order to reduce external electrical noise, eliminate the effects of external grounded objects, and protect the electrodes from damage, a screen is used to cover the electrodes, which can be either conductive or non-conductive. An insulating material is used to fill the gap between the electrodes and the screen. To reduce the standard capacitance between adjacent electrode pairs, the sensor also includes radial guards. The screen and the radial guards are always maintained at earth potential. Different configurations of the radial and axial guards were proposed and studied to improve the system performance.

Numerous measurement protocols can be used for an *N*-electrode system. However, in a typical ECT system, electrodes 1 to *N* are used consecutively as source electrodes, and the capacitance values among all pair combinations are measured. This means that, for one measurement cycle of a 12-electrode ECT sensor, electrode 1 will first act as the source electrode after being excited with a fixed positive voltage, whereas electrodes 2–12 will act like detectors by being held at earth potential. Afterwards, the capacitances between 1–2, 1–3, …, 1–12 are measured concurrently using transducers. Next, electrode 2 will act as the source electrode, and electrodes 3–12 act as detectors, whereas electrode 1 will be inactive to eliminate reciprocity. This process continues until electrode 11 acts as the source electrode and only electrode 12 acts as the detector, whereas electrodes 1–10 are inactive. Commonly, the number of independent capacitance measurements *M* is given by the number of *N*-electrode combinations, taking two at a time, as *M*=*N*(*N*−1)/2.

The data acquisition unit measures and conditions the signals that are acquired from the capacitance sensor. The data acquisition unit works according to the following sequence of events:

Apply a fixed voltage,

*V*_{c}, to the source electrode while keeping the rest of the detectors grounded.Measure the induced currents in the detectors and their capacitive components by the use of phase sensitive detector.

Deduce the capacitance measurements

*C*from the given potential difference*V*_{c}and the measured current*Q*according to*C*=*Q*/*V*_{c}.Convert the capacitance measurements into voltage signals using a capacitance-to-voltage converter.

Condition (i.e. multiplex, amplify, filter, etc.) and digitize the voltage signals.

Send the conditioned and digitized voltage signals to a computer for further processing and display.

The computer unit receives the voltage signals from the data acquisition unit and interprets them again into capacitance measurements in order to use them for the image reconstruction algorithm to determine the permittivity distribution (i.e. the image).

## 3. Image reconstruction

### (a) The forward problem

The ECT forward problem needs to be solved before solving the ECT inverse problem. The ECT mathematical system can be modelled by Poisson's equation, as it can be treated as an electrostatic field problem. Assuming no free charge within the field, Poisson's equation can be simplified to Laplace's equation,
3.1where *ε*_{0} is the free space permittivity, and *ε*(*x*,*y*) and *u*(*x*,*y*) are the two-dimensional permittivity distribution and the electrical potential distribution, respectively. Laplace's equation is solved with the following Dirichlet boundary conditions:
3.2Here *V* _{c} is a fixed voltage applied to the source electrode, *Γ*_{i} is the spatial locations of electrode *i*, and *Γ*_{guards} and *Γ*_{screen} are the spatial locations of the guards and screen, respectively. In other words, (3.2) states that the boundary condition is a fixed voltage *V* _{c} applied to the source electrode while the detectors, guards and screen are kept grounded. The measured capacitance can be calculated by
3.3where *V* _{ij} is the voltage difference between the detector *j* and the source electrode *i* (i.e. *V* _{i}−*V* _{j}), *Γ*_{j} is the surface of the detector *j*, is a unit vector normal to *Γ*_{j} and *u* is the electrical potential distribution. The finite element method (FEM) is used to calculate the capacitance and the sensitivity distribution [15]. This study presents a nonlinear multiphase level set algorithm, which will work with one set of capacitance data (i.e. absolute capacitance data). In this case, the background data of free space are used for a simple scaling calibration of the forward model.

### (b) Level set-based inverse problem

Osher & Sethian [18] first developed the level set method in the 1980s for modelling the front propagation of surfaces under curvature. Since then, it has been used in many disciplines, such as image processing, fluid mechanics, optimal design, computer graphics, computational geometry and inverse problems. Santosa [19] was the first to develop two computational approaches for inverse problems based on the level set method in which the desired unknown is a region in or . The first approach results in a Hamilton–Jacobi equation, which is a nonlinear time-dependent partial differential equation for the level set function whose evolution minimizes the residual in the data fit. Whereas the second approach is based on an optimization method (i.e. Gauss–Newton) that generates a sequence of level set functions that reduces the residual. For electrical and electromagnetic tomography in particular, some initial studies were conducted using different types of level set techniques for shape reconstruction [16,20–22].

In the level set method, the constructed shapes are given as the zero level set of a one- or higher-dimensional function called the ‘level set function’. By updating this function, the constructed shapes can move accordingly, and thus the topological changes can be performed automatically. Hence, instead of finding *ε*_{int}, the level set technique is able to represent boundaries of different shapes by the zero level set of a predefined level set function *ψ*. For a two-dimensional problem, suppose *D* is the characteristic set of interest. Then the boundary of the inclusion ∂*D* is the interface between the two materials and is given by the zero level set of the function *ψ* as
3.4Thus, the image parameter at each point, *ε*(*r*), can be represented in terms of a level set function as
3.5Accordingly, the capacitance data can be described as a function of a level set function *C*=*F*(*ψ*(*ε*(*r*))).

Equation (3.5) represents the mapping that assigns to a given level set function the corresponding permittivity distribution. This can be denoted as *ε*=*φ*(*ψ*(*r*)). Representing *ε*(*r*) through a level set function *ψ*(*r*) instead of representing it in terms of *D* makes the image reconstruction process simpler, as no *a priori* assumption about the topology (i.e. spatial location) or the nature (i.e. shape and size) of *D* is needed. On the other hand, owing to the nonlinear dependence of *ε* on *ψ*, this approach will lead to a nonlinear inverse problem even if the forward problem was linearized. Note that, in this inverse problem statement, *ε*(*r*) is not the unknown; instead, it only performs as an intermediate parameter to arrive at the real unknown, which is the level set function *ψ*(*r*).

In order to facilitate the dependence of the forward map on small changes in the shape, we look at the variation of *ε* caused by a variation in *ψ*. Let *r* be a point on the boundary of *D* which is represented by (3.4). Suppose *ψ*(*r*) and *r* are perturbed by a small variation, giving *ψ*′(*r*)=*ψ*(*r*)+*δψ*(*r*) and *r*′=*r*+*δr*, respectively. This will result in changing the region *D* to a new region denoted by *D*′. Accordingly, *ε* will vary to *ε*′=*ε*+*δε*. In figure 2, note that *r*, which was on the boundary of *D* and on *ψ*(*r*)=0, now belongs to *ε*_{int} of *ψ*(*r*) and the boundary of *D*′(*ψ*′(*r*)=0) according to (3.5). Figure 2 shows how updating the level set function will change the status of each point *r* to belong to either *ε*_{int} or *ε*_{ext} of the new level set *ψ*′(*r*).

Equation (3.5) is suitable for two-phase reconstruction; however, it can be extended to cover multiple phases, as several authors have suggested. One of the earliest and most common methods for multiphase image reconstruction is based on the colour level set technique, which was first proposed by Vese & Chan [16] for image segmentation. In this technique, the evolution of up to 2^{p} different phases can be described using *p* level set functions. However, in this study, a modification of the colour level set method suggested by [23] was employed for the ECT reconstruction. In this approach, *m* different phases are described by using *p*=*m*−1 level set functions. Therefore, three level set functions are used to represent four phases (figure 3), where each of the four phases is allowed to have any complex topology. Accordingly, (3.5) can be extended to read
3.6

This level set function *ψ* is given as an *a priori* assumption at first and then updated in such a way that moves the boundaries accordingly. Such updates are chosen in a way that reduces a given cost function. One of the first employed techniques for updating moving boundaries is based on a time evolution approach leading to a Hamilton–Jacobi equation. However, another approach for shape reconstruction using a level set is the use of optimization methods as an alternative to shape evolution, which was first suggested by Santosa [19]. To find the update of the level set, different optimization techniques can be employed. This study employs a modified optimization approach based on the Gauss–Newton method to derive the updates for the level set function.

In the optimization approach, a sequence of surfaces *ψ*_{k}(*r*) such that *D*_{k}→*D*, where ∂*D*_{k}={*r*:*ψ*(*r*)=0}, is generated. To solve the inverse problem in a least-squares sense, the minimizer of the objective functional *J*(*ε*) is found:
3.7Here *α*>0 is the regularization parameter that can be tuned with different mathematical techniques and Reg(*ψ*(*ε*(*r*))) is the regularization term, which is an arbitrary functional that takes different forms depending on the regularization method applied, and the argument in the generalized equation can be set according to the requirements of the solution. Different regularization methods have been used to enhance the level set method, such as curve shortening [16,17,23]. The regularization term has much less effect on the problem when using a level set method, which makes the problem better posed. Using a modified Gauss–Newton approach for minimizing the nonlinear least-squares functional *J*(*ψ*) in equation (3.7), we arrive at the following equation for the level set evolution:
3.8Here *B*=*SK* and *K*=*cχ*. In (3.8) the variables have the following meaning:

—

*L*is the regularization matrix. In this paper, we considered*L*^{T}*L*=*I*; alternatively other regularization matrices can also be used, such as Laplacian.—

*ξ*is the step size (relaxation) parameter, which is the magnitude of change in shape in every update. The step size can stay constant for all iterations or can be changed with the iterations. For this study, the step size was chosen by trial and error and was held constant for all iterations.—

*B*is the narrowband Jacobian matrix, which represents the sensitivity values of only those few elements surrounding ∂*D*of all different regions.—

*S*is the sensitivity matrix for all elements of the FEM mesh.—

*K*is the discretized form of*ε*=*φ*(*ψ*(*r*)), which includes the indicator function*χ*of a small narrow band of half-width an FEM element centred at ∂*D*and a constant*c*as a normalization factor.

As can be noticed from equation (3.8), the inverse problem of the interfaces between two different materials is solved using a narrowband method where only the surrounding pixels are inverted. This technique tends to improve the computational cost and the condition number of the discrete inverse problem by decreasing the number of unknowns compared with the traditional pixel-based image reconstruction.

Both the regularization parameter *α* and the step size parameter *ξ* need to be tuned in order to arrive at their optimal choice, which depends on the mesh density, the permittivity contrast and the initial guess of the permittivity values. In general, the larger *α* is chosen, the smoother the updates will be and the harder for the scheme to perform splitting of the topological changes. In this study, *α* was manually tuned for each case. However, other techniques can be applied to make the reconstruction process more automatic.

The initial guess of the level set functions can be any kind of function; however, the most common choice is a distance function *ψ*^{k}=−dist(∂*D*). In this paper, we choose three separated circular inclusions to be the initial guesses for all three level set functions as
3.9where (*X*_{0},*Y* _{0}) are the Cartesian coordinates of the centre and *ρ* is the radius of the circular inclusion. Figure 4 shows a flowchart of the level set ECT image reconstruction algorithm developed for this study.

## 4. Evaluation of multiphase monitoring using level set method

To demonstrate the power of this algorithm for multiphase monitoring, simulated and experimental results of low-contrast and high-contrast cases are shown and discussed. For the low-contrast case, more than two different phases, which all have low permittivity values, are imaged and analysed. This case can be related to real scenarios such as when a pipeline is carrying *ε*_{1} = gas, *ε*_{2}= oil and *ε*_{3}= low-permittivity deposits. On the other hand, for the high-contrast case, a combination of low-permittivity and high-permittivity phases is imaged and analysed. This case can be related to real scenarios such as when a pipeline is carrying *ε*_{1}= gas, *ε*_{2}= oil, *ε*_{3}= deposits and *ε*_{4}= water. This case is relatively harder to image using traditional pixel-based methods, as the inverse problem becomes highly nonlinear when a high-permittivity component is introduced, as it dominates the region. FEM mesh models were designed and used for the inverse problem, where 1367 three-node non-homogeneous triangular elements (i.e. pixels) corresponding to 739 nodes were generated to create the mesh.

### (a) Simulation results

For the sake of comparison, the traditional Tikhonov image reconstruction simulation results for both cases are shown and discussed to illustrate the advantage of the shape-based reconstruction methods over the classical pixel-based methods.

#### (i) Case 1: low-contrast simulation

A four-phase simulation was conducted, and the results for the low-contrast case for both the level set method and the standard Tikhonov image reconstruction method are shown in figure 5 along with the true shape for the sake of comparison. The permittivity values used in the simulation were *ε*_{1}=4,*ε*_{2}=3,*ε*_{3}=2 and *ε*_{4}=1. In order to construct the level set image, 250 iterations were used, and the norms of the error between the simulated and calculated capacitances for all iterations are shown in figure 6. The level set reconstructed image shown in figure 5*b* is based on the minimum residual norm, which in this case appeared in iteration 155. This shows that the best result might not occur in the latest iteration, rather somewhere in the iteration process. Accordingly, the algorithm might just terminate at the minimum residual iteration to reduce the computation cost.

Figure 7 shows the histograms of the permittivity distribution per pixel for the true shape, the level set reconstruction and the standard Tikhonov image reconstruction, respectively. The obvious similarity between the true shape and the level set reconstruction histograms demonstrates the superiority of the level set method over the traditional Tikhonov reconstruction. The histogram of the level set method shown in figure 7*b* shows only four phases, and each pixel strictly belongs to one and only one of the four possible phases. On the other hand, the histogram of the standard Tikhonov image reconstruction method shown in figure 7*c* shows a large range of permittivity values, which makes such methods harder to threshold, and hence information about the individual phases is harder to obtain.

#### (ii) Case 2: high-contrast simulation

A four-phase simulation was conducted, and the results for the high-contrast case for both the level set method and the standard Tikhonov image reconstruction method are shown in figure 8 along with the true shape for the sake of comparison. The permittivity values used for the four-phase simulation were *ε*_{1}=80,*ε*_{2}=3,*ε*_{3}=2 and *ε*_{4} = 1. In order to construct the level set image, 250 iterations were used, and the norms of the error between the simulated and calculated capacitances for all iterations are shown in figure 9. The level set reconstructed image shown in figure 8*b* is based on the minimum residual norm, which in this case appeared in iteration 172. Figure 10 shows the histograms of the permittivity distribution per pixel for the true shape, the level set reconstruction and the Tikhonov image reconstruction, respectively.

From figure 8*a*,*b*, it can be seen that the high-permittivity phase of water is dominating the region and thus it is hard to see with the naked eye the lower-permittivity phases. However, all phases do exist in the data, and this is evident by the histograms of the permittivity distribution per pixel shown in figure 10*a*,*b*. In addition, the number of reconstructed pixels per phase was extracted to calculate the phase volume fraction and it is shown in table 1.

#### (iii) Simulated phase area fraction

By definition, the phase area fraction is the cross-sectional area locally occupied by the phase of a multiphase region, divided by the cross-sectional area of the entire region (e.g. pipe) at the same local point. This can be translated digitally by dividing the number of pixels per phase over the total number of pixels. Consequently, the area fractions of all phases were calculated according to the equation below and are given in table 1: 4.1

Calculating the phase area fraction demonstrates the ability to extract information using ECT by the level set method. The phase fraction is a key parameter for determining other important parameters such as density, viscosity, average velocity, pressure drop and heat transfer.

### (b) Experimental results

For the experimental evaluation, a 12-electrode sensor was used, which was built in the University of Bath Engineering Tomography Laboratory (ETL). The sensor has an inner diameter of 151 mm, whereas the data acquisition unit used is a commercial PTL300E-TP-G capacitance measurement unit from Process Tomography Limited (figure 1*b*). The unit can collect sets of capacitance data at 100 frames s^{−1} with an effective resolution of 0.1 fF and measurement noise level better than 0.07 fF. Matlab was used for data collection, image reconstruction and processing. Four different experimental models were prepared in the laboratory in order to examine the ability of the algorithm to detect and distinguish different phases.

#### (i) Case 3: low-contrast experiment

The experimental results for the low-contrast case for the level set method are shown in figure 11 along with the experimental set-ups for the sake of comparison. Two experimental set-ups were investigated:

— Set-up 1 is a three-phase set-up where two wooden rods of permittivity

*ε*_{1}=2.3 and a plastic rod of permittivity*ε*_{2}=3 were reconstructed against an air background of permittivity*ε*_{3}=1.— Set-up 2 is a four-phase set-up where a plastic rod of permittivity

*ε*_{1}=3 and a wooden rod of permittivity*ε*_{2}=2.3 were divided by a book of permittivity*ε*_{3}=3.85 and reconstructed against an air background of permittivity*ε*_{4}=1. This set-up is more challenging, because all phases are close in permittivity values and are closely located to each other.

Such simple experimental set-ups are sufficient for a static study as long as they provide different phases (i.e. different permittivity values) and can be considered as a valid first step towards studying multiphase process fluids.

In this case, 250 iterations were used to construct the level set image, and the norms of the error between the measured and calculated capacitances for all iterations are shown in figure 12. The level set reconstructed images shown in figure 11 are based on the minimum residual norm. Standard Tikhonov reconstruction was not produced for the experimental data, because absolute capacitance data were used to generate the level set images where only one measured dataset was used without the need for a reference dataset. Whereas standard Tikhonov reconstruction needs a reference dataset to produce good results. Hence, using standard Tikhonov reconstruction when absolute data are used could not generate meaningful images.

Figure 13 shows the histograms of the permittivity distribution per pixel. Comparison between the true shape and the level set reconstruction histograms is hard to achieve for experimental set-ups; however, these histograms demonstrate the superiority of the level set method over the traditional Tikhonov reconstruction. This can be seen in the constant permittivity values of the level set method (each constant represents one phase).

#### (ii) Case 4: high-contrast experiment

The experimental results for the high-contrast case for the level set method are shown in figure 14 along with the experimental set-ups for the sake of comparison. Two set-ups were investigated:

— Set-up 3 is a three-phase set-up where water of permittivity

*ε*_{1}=80 (ignoring the thin plastic bottle) and two wooden rods of permittivity*ε*_{2}=2.3 were reconstructed against an air background of permittivity*ε*_{3}=1.— Set-up 4 is a four-phase set-up where water of permittivity

*ε*_{1}=80, a plastic rod of permittivity*ε*_{2}=3 and two wooden rods of permittivity*ε*_{3}=2.3 were reconstructed against an air background of permittivity*ε*_{4}=1. This set-up is more challenging, because wood and plastic are close in permittivity values and are closely located to each other.

In this case, 250 iterations were used to construct the level set images, and the norms of the error between the measured and calculated capacitance for all iterations are shown in figure 15. The level set reconstructed images shown in figure 14 are based on the minimum residual norm. Because water with its high permittivity is dominating the area, it is hard to distinguish the lower-permittivity phases from the reconstructed images. However, the number of pixels per phase can be extracted as seen in the histogram (figure 16) of the permittivity distribution per pixel. Accordingly, the phase area fraction can be extracted and analysed as in table 2.

## 5. Discussion

One objective of this study is to demonstrate the ability of the level set method to detect and quantify all phases despite any previous knowledge. From the results, it can be seen that this objective can be better achieved by the level set method than by traditional regularization-based image reconstruction. This can be observed from the reconstructed level set images, which show patterns similar to the true patterns shown in the simulated and experimental results. Another advantage of the level set shape reconstruction method over traditional image reconstruction methods is the ability to accurately estimate the phase volume fraction, as illustrated in tables 1 and 2. This is because, unlike traditional image reconstruction methods, the level set method does not require any *a priori* information except for the permittivity values of each phase.

The results shown in this paper are designed in order to study the potential of this algorithm for multiphase imaging. The algorithm shows promising results; however, it is also open to enhancement. The performance of this method and hence the accuracy can be further investigated and enhanced by refining the following features:

— Generally, the permittivity values may be known and supplied to the algorithm as

*a priori*assumption, but methods can be extended where both the permittivity values and the topology are determined simultaneously [24].— A pre-step of determining the best initial guess of the level set functions can improve the speed of convergence.

— The level set model can be improved by modifying the number/type of used level set functions and/or the mapping

*ε*=*φ*(*ψ*(*r*)) shown in (3.6).— The tuning parameters (i.e. the step size, the regularization parameter) and the regularization matrix can be chosen differently using other methods to arrive at the best choices.

— In order to adapt this method to online monitoring, the tuning parameters can be automated.

— Applying a less dense mesh or stopping the iteration process if the residual has reached a specified set point can enhance the speed of convergence.

— Applying a denser mesh or allowing further iterations may enhance the resolution of the image.

— This technique can also be extended to three-dimensional image reconstruction.

— For a typical PC like the one used for this study, the nonlinear level set algorithm presented here takes 34 s to run the experimental data for a mesh of 1367 elements (the mesh used for this study), and 65 s for a mesh of 2064 elements. Such computation times are not ideal for online measurements, as they will not allow real-time measurements of the process. Hence, other optimization techniques can be applied to arrive at the best results in terms of image resolution and computation time. One important optimization would be to implement this algorithm in high-performance computing such as GPU.

In addition to the fact that the level set approach does promise to make an ill-posed problem better behaved, it offers other advantages such as the following:

— the level set method allows one to perform numerical computations of a geometrical object containing curves and surfaces on a fixed Cartesian grid without the need to parametrize such objects;

— the level set method has the ability to follow shapes with changing topology (i.e. split in to multiregions, develops holes or vice versa);

— the level set method does not require any

*a priori*assumptions about the topology (i.e. spatial location) or the nature (i.e. shape and size) of the inclusion;— the level set method is well suited for high-contrast image reconstruction with sharp boundaries, because the high-contrast inclusion is incorporated explicitly in the level set model without it dominating the entire region and overshadowing the other low-contrast inclusions, unlike traditional regularization-based techniques, which have the limitation of oversmoothing the reconstructed image especially when a high-contrast inclusion is dominating the region.

Because of such advantages, the level set method has a great potential for shape reconstruction applications and thus for multiphase image reconstruction.

## 6. Conclusion

The results in this paper demonstrate the potential and feasibility of the proposed level set technique for detecting, locating and characterizing all the different phases possibly appearing in a multiphase process. Simulation and experimental results demonstrate that a multiphase level set-based algorithm is capable of imaging three- or four-phase permittivity contrasts. Such promising results can also be extended to image more than four phases by adding more level set functions to the model accordingly. The nonlinear nature of the proposed algorithm as well as prior knowledge about the contrast permittivity values enabled us to recover those permittivity profiles by using one set of capacitance data. This is critically important in real-life applications where absolute ECT data can be used to recover multiphase flow data without requiring reference data. The hope is that the results of this paper will help in better understanding and further extending the application of ECT for multiphase flow imaging.

## Authors' contributions

E.A.H. carried out the experiments, wrote the multiphase code and drafted the manuscript. M.S. conceived the idea, supervised E.A.H. and revised the manuscript critically. M.S. and E.A.H. both tested and evaluated the real data and simulated results and gave final approval for publication.

## Competing interests

The authors declare that they have no competing interests.

## Funding

We received no funding for this study.

## Footnotes

One contribution of 10 to a theme issue ‘Super-sensing through industrial process tomography’.

- Accepted March 2, 2016.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.