In this work, we propose a total least-squares formulation for recovering multiple medium coefficients for a class of inverse medium problems that are governed by the forward model of the form:
where is a bilinear operator on , is the state variable, while represents one or multiple unknown coefficients in the model that are to be recovered, under some measurement data of . Here , and are three Hilbert spaces, and is an observation map from to .
In many applications, it is often required to recover multiple coefficients simultaneously. For instance, the diffusive optical tomography (DOT) aims at recovering the diffusion and absorption coefficients and from the governing equation [1, 2]:
using the Cauchy data collected at the boundary of :
using some electric or magnetic measurement data or .
The inverse reconstruction of multiple medium coefficients is generally much more technical and difficult than the single coefficient case. We shall proposed a total least-squares formulation with appropriate regularization to transform the inverse problem into an optimization problem. The total least-squares philosophy is not uncommon. One conventional approach for inverse medium problems is to seek an optimal parameter from a feasible set such that it minimizes an output least-squares functional of the form
where solves the forward model (1.1) when is given, is a regularization term and is a regularization parameter. We refer the readers to [2, 5, 6] for more details about this traditional approach. A relaxed variational method of the least-squares form was proposed and studied in [7, 8] for the impedance computed tomography. The least-squares functional consists of a residual term associated with the governing equation while the measurement data is enforced at the feasible set. Different from the aforementioned approaches, we shall follow some basic principle of a total least-squares approach from  and treat the governing equation and the data fitting separately, along with a regularization. That is, we look for an optimal parameter from and a state variable from together such that they minimize an extended functional of the form
This functional combines the residual of model equation, data fitting and constraints on parameters in the least-squares criterion. The combination in (1.5) results in a regularization effect to treat the model equation and allows a robustness as a quasi-reversibility method . Compared with the conventional approaches, the domain of is much more regular and the semi-norm defined by the formulation is much stronger. More precisely, the total least-squares formulation aims to find a solution pair simultaneously in a smooth class and is less sensitive to the noise and uncertainty in the inverse model. Another important feature of this formulation is that the functional is quadratical and convex with respect to each variable and , if the regularization is chosen to be quadratical and convex, while the traditional one in (1.4) is highly nonlinear and non-convex in general. This special feature facilitates us naturally to minimize the functional effectively by the alternating direction iterative (ADI) method [11, 12] so that only two quadratical and convex suboptimizations are required in terms of the variables and respectively at each iteration.
In addition to the functional (1.5) that uses the residual of the forward model (1.1), we will also address another least-squares functional that makes use of the equivalent first-order system of the forward model (1.1) and replaces the first term in (1.5) by the residuals of the corresponding first-order system. Using the first-order system has been the fundamental idea in modern least-squares methods in solving second-order PDEs [13, 14, 15]. The advantages of using first-order formulations are much more significant to the numerical solutions of inverse problems, especially when we aim at simultaneously reconstructing multiple coefficients as we do in this work. First, the multiple coefficients appear in separated first-order equations, hence are naturally decoupled. This would greatly reduce the nonlinearity and enhance the convexity of the resulting optimization systems. Second, the first-order formulation relaxes the regularity requirement of the solutions in the resulting analysis.
A crucial step to an effective reconstruction of multiple coefficients is to seek some reasonable initial approximations to capture some general (possibly rather inaccurate) geometric and physical profiles of all the unknown multiple coefficients. This is a rather technical and difficult task in numerical solutions. For this purpose, we shall propose to adopt the direct sampling-type method (DSM) that we have been developing in recent years (cf. [17, 18, 26, 27]). Using the index functions provided by DSM, we shall determine a computational domain that is often much smaller than the original physical domain, then the restricted index functions on the computational domain serve as the initial guesses of the unknown coefficients. In this work, we will apply a newly developed DSM , where two groups of probing and index functions are constructed to identify and decouple the multiple inhomogeneous inclusions of different physical nature, which is different from the classical DSMs targeting the inhomogeneous inclusions of one single physical nature. As we shall see, DSMs turn out to be very effective and fast solvers to provide some reasonable initial approximations.
The rest of the paper is structured as follows. In Section 2, we justify the well-posedness of the least-squares formulation for the general inverse medium problems. In Section 3, we propose an alternating direction iterative method for solving the minimization problem and prove the convergence of the ADI method. We illustrate in Section 4 how this total least-squares method applies to a concrete inverse problem, by taking the DOT problem as a benchmark problem. We present numerical results in Section 5 for a couple of different types of inhomogeneous coefficients for the DOT problem to demonstrate the stability and effectiveness of this proposed method. Throughout the paper, , and denote generic constants which may differ at each occurrence.
2 Well-posedness of the least-squares formulation for inverse medium problem
Recall that to solve the inverse medium problems modeled by (1.1), we propose the following least-squares formulation
This section is devoted to the well-posedness of the total least-squares formulation (2.1), namely, the existence of a solution to (2.1) and the conditional stability of the reconstruction with respect to the measurement. To provide a rigorous justification of the well-posedneness, we present several assumptions on the least-squares formulation, which are minimal for the proof. We will verify these checkable assumptions in Section 4 for a concrete example of such inverse medium problems.
Let us first introduce several notations. For simplicity, for a given (resp. ), we will write (resp. ) as
We denote the subdifferential of the regularization term at by , and denote the inner products of the Hilbert spaces , and by , and respectively.
2.1 Existence of a minimizer
We present the following assumptions on the regularization term and operators and in the forward model:
The regularization term is strictly convex and weakly lower semicontinuous. Furthermore, is also coercive , i.e., .
This assumption implies that the level set defines a bounded set in .
Given a constant , for in the level set , , where with a specific side constraint that is not in the least-squares formulation (2.1), is a closed linear operator and is uniformly coercive, i.e., the graph norm satisfies uniformly in for some constant , and thus has a unique solution in .
Under Assumption 2, we can define the inverse operator , which is uniformly bounded by the coercivity of . We also need the following assumption on the sequentially closedness of operators and .
The operators and are weakly sequentially closed, i.e., if a sequence converges to weakly in , then the sequence converges to weakly in and the sequence converges to weakly in .
Then we can verify the existence of the minimizers to the least-squares formulation (2.1).
Since and are nonempty, there exists a minimizing sequence in such that
By Assumptions 1–2, is a coercive functional and the graph norm is uniformly coercive, thus it follows (1.5) that the sequence is uniformly bounded. Then there exists a subsequence of , still denoted by , and some such that converges to weakly in and converges to weakly in . As and are weakly sequentially closed by Assumption 3, there hold
From the weak lower semicontinuity of the norms and , we have
Together with the lower semicontinuity of the regularization term , we can deduce that
Hence it follows (2.3) that is indeed a minimizer of the functional in . ∎
2.2 Conditional stability
In this subsection, we present some conditional stability estimates of the total least-squares formulation (2.1) for the general inverse medium problems. First we introduce two pairs and that satisfy
Letting be the unique minimizer of (2.1) in a neighborhood of , we study the approximation error to illustrate the stability of the least-squares formulation (2.1) with respect to the measurement and also the term in the governing equation (1.1). Denote the residual of the governing equation by . As is the local minimizer of functional in (1.5), we have . Therefore, by the definition of we have the inequality
which directly leads to the following observation on :
If is coercive, (2.7) provides a rough estimate of the reconstruction with respect to the data noise.
We can further derive an estimate of the approximation error , under the following assumption on the operator .
There exists a norm on such that for any ,
Assumption 4 holds when belongs to a finite rank subspace of and for all non-zero . We can now deduce the following result of the approximation error .
Using the bilinear property of , one can rewrite the difference as
It follows Assumption 4 that there exists a norm such that
where is a constant, which completes the proof. ∎
The rest of this section is devoted to verifying the consistency of the least-squares formulation (2.1) as the noise level of measurement goes to zero, which is an essential property of a regularization scheme. If we choose an appropriate regularization parameter according to the noise level of the data, we can deduce the convergence result of the reconstructed coefficients associated with the regularization parameter . More precisely, given a set of exact data , we consider a parametric family such that . In the rest of this section, we denote the functional in (1.5) with and by , and the minimizer of by . Then we justify the consistency of the least-squares formulation (2.1) by proving the convergence of the sequence of minimizers to the minimum norm solution  of the system (2.5) as .
As is the minimizer of , there holds
By definition, , and thus is uniformly bounded. Following the similar argument in the proof of Theorem 1, there exist a subsequence of , still denoted as , and some such that converges to weakly in . By Assumption 3, we have
From (2.13) one can also derive that
Therefore, as , will converge to satisfying
Recall that one has an estimate of from (2.14) that
which leads to
as . Using the lower semicontinuity of functional , one obtain that
3 ADI method and convergence analysis
An important feature of the least-squares formulation (2.1) is that the functional is quadratical and convex with respect to each variable and . This unique feature facilitates us naturally to minimize the functional effectively by the alternating direction iterative (ADI) method [11, 12] so that only two quadratical and convex suboptimizations of one variable are required at each iteration. We shall carry out the convergence analysis of the ADI method in this section for general inverse medium problems.
Alternating direction iterative method for the minimization of (1.5).
Given an initial pair , find a sequence of pairs for as below:
Given , find by solving
Given , find by solving
We shall establish the convergence of the sequence generated by the ADI method, under Assumptions 1–3 on the least-squares formulation (2.1). For this purpose, we would like to introduce the Bregman distance  with respect to ,
Using the optimality condition satisfied by the minimizer of (3.1), we deduce
Similarly, from the optimality condition satisfied by the minimizer of (3.2), we obtain
Taking , and in (3.3), we can derive that
The following equality hold for the first term in (3.6):
which implies that is bounded. Then we can conclude using Assumption 2 that forms a Cauchy sequence and thus converges to some . Since is uniformly bounded for all , we can derive that converges to some from the strict convexity of . As the sequence is monotone decreasing, there exists a limit . Following the similar argument in the proof of Theorem 1, we conclude that by Assumption 3. This completes the proof of the convergence. ∎
4 Diffusive optical tomography
We take the diffusive optical tomography (DOT) as an example to illustrate the total least-squares approach for the concrete inverse medium problem in this section. We will introduce the mixed regularization term, present the least-squares formulation of the first-order system of DOT, and then verify the assumptions in Section 2 for the proposed formulation. We shall use the standard notations for Sobolev spaces. The objective of the DOT problem is to determine the unknown diffusion and absorption coefficients and simultaneously in a Lipschitz domain () from the model equation
with a pair of Cauchy data on the boundary , i.e.,
Throughout this section, we shall use the notation in the total least-squares formulation, where denotes the Neumann to source map. To define the Neumann to source map, we first introduce the boundary restriction mapping on , i.e., denotes the boundary value of . Then we use to denote the trace operator , which is formally defined to be the unique linear continuous extension of as an operator from onto . Using Riesz representation theorem, there exists a function in , denoted by , such that for any ,
which will be denoted by in the least-squares formulation.
4.1 Mixed regularization
In this subsection, we present the mixed regularization term for the DOT problem. As the regularization term in the least-squares formulation (2.1) shall encode the priori information, e.g., sparsity, continuity, lower or upper bound and other properties of the unknown coefficients, it is essential to choose an appropriate regularization term for a concrete inverse problem to ensure satisfactory reconstructions. In this work, we introduce a mixed – regularization term for a coefficient :
In this revision we change the semi-norm of into the full norm in the regularization to ensure the strict convexity of . where is given by
and and are the lower and upper bounds of the coefficient respectively. The first term is the regularization term, the second term corresponds to the regularization, and the third term enforces the reconstruction to meet the constrains of the coefficient.
In practice, the regularization enhances the reconstruction and helps find geometrically sharp inclusions, but might cause spiky results. The regularization generates the reconstructions with overall clear structures, while the retrieved background may be blurry. Compared with other more conventional regularization methods, this mixed regularization technique in (4.2) combines two penalty terms and effectively promotes multiple features of the target solution, that is, it enhances the sparsity of the solution while it still preserves the overall structure at the same time. The scalar parameters and need to be chosen carefully to compromise the two regularization terms.
4.2 First-order formulation of DOT
As we have emphasized in the Introduction, it may have some advantages to make use of the residuals of the first-order system of the model equation (4.1), instead of the residual of the original equation in the formulation (2.1), when we aim at recovering two unknown coefficients and simultaneously. Similarly to the formulation (2.1), we now have , and the operator is given by
where we have introduced an auxiliary vector fluxand each entry of is of a first-order form such that two coefficients are separated naturally. Clearly, is still bilinear with respect to the state variables and coefficients . Using the first-order system, we can then come to the following total least-squares functional:
where is the trace operator, and and are the corresponding mixed regularization terms of and defined as in (4.2), , are lower and upper bounds of , and , are lower and upper bounds of . We shall minimize (4.4) over , that is, , .
We will apply the ADI method to solve the least-squares formulation of and :
Given an initial guess , we find a sequence for as below:
Given , , find , by solving
Given , , find , by solving
4.3 Well-posedness of the least-squares formulation for DOT
Recall that we have proved the well-posedness of the least-squares formulation in Section 2 for general inverse medium problems. This part is devoted to the verification of assumptions in Section 2 for the formulation (4.5) to ensure its well-posedness.