# Applications of Fixed-Point and Optimization Methods to the Multiple-Set Split Feasibility Problem

Authors: Yonghong Yao; Rudong Chen; Giuseppe Marino; Yeong Cheng Liou
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

The multiple-set split feasibility problem requires finding a point closest to a
family of closed convex sets in one space such that its image under a linear transformation
will be closest to another family of closed convex sets in the image space.
It can be a model for many inverse problems where constraints are imposed on the
solutions in the domain of a linear operator as well as in the operator’s range. It
generalizes the convex feasibility problem as well as the two-set split feasibility
problem. In this paper, we will review and report some recent results on iterative approaches
to the multiple-set split feasibility problem.


---

## Body

## 1. Introduction

### 1.1. The Multiple-Set Split Feasibility Problem Model

The intensity-modulated radiation therapy (IMRT) has received a great deal of attention recently; for related works, please refer to [1–29]. In intensity modulated radiation therapy, beamlets of radiation with different intensities are transmitted into the body of the patient. Each voxel within the patient will then absorb a certain dose of radiation from each beamlet. The goal of IMRT is to direct a sufficient dosage to those regions requiring the radiation, those that are designated planned target volumes (PTVs), while limiting the dosage received by the other regions, the so-called organs at risk (OAR). The forward problem is to calculate the radiation dose absorbed in the irradiated tissue based on a given distribution of the beamlet intensities. The inverse problem is to find a distribution of beamlet intensities, the radiation intensity map, which will result in a clinically acceptable dose distribution. One important constraint is that the radiation intensity map must be implementable; that is, it is physically possible to produce such an intensity map, given the machine’s design. There will be limits on the change in intensity between two adjacent beamlets, for example.The equivalent uniform dose (EUD) for tumors is the biologically equivalent dose which, if given uniformly, will lead to the same cell kill within the tumor volume as the actual nonuniform dose. Constraints on the EUD received by each voxel of the body are described in dose space, the space of vectors whose entries are the doses received at each voxel. Constraints on the deliverable radiation intensities of the beamlets are best described in intensity space, the space of vectors whose entries are the intensity levels associated with each of the beamlets. The constraints in dose space will be upper bounds on the dosage received by the OAR and lower bounds on the dosage received by the PTV. The constraints in intensity space are limits on the complexity of the intensity map and on the delivery time, and, obviously, that the intensities be nonnegative. Because the constraints operate in two different domains, it is convenient to formulate the problem using these two domains. This leads to a split feasibility problem.The split feasibility problem (SFP) is to find anx in a given closed convex subset C of RJ such that Ax is in a given closed convex subset Q of RI, where A is a given real I by J matrix. Because the constraints are best described in terms of several sets in dose space and several sets in intensity space, the SFP model needs to be expanded into the multiple-set SFP. It is not uncommon to find that, once the various constraints have been specified, there is no intensity map that satisfies them all. In such cases, it is desirable to find an intensity map that comes as close as possible to satisfying all the constraints. One way to do this, as we will see, is to minimize a proximity function.Fori=1,…,I and j=1,…,J, let bi≥0 be the dose absorbed by the ith voxel of the patient’s body, xj≥0 the intensity of the jth beamlet of radiation, and Aij≥0 the dose absorbed at the ith voxel due to a unit intensity of radiation at the jth beamlet. The nonnegative matrix A with entries Aij is the dose influence matrix. Let us assume that we have M constraints in the dose space and N constraints in the intensity space. Let Hm be the set of dose vectors that fulfill the mth dose constraint, and let Xn be the set of beamlet intensity vectors that fulfill the nth intensity constraint.In intensity space, we have the obvious constraints thatxj≥0. In addition, there are implementation constraints; the available treatment machine will impose its own requirements, such as a limit on the difference in intensities between adjacent beamlets. In dosage space, there will be a lower bound on the dosage delivered to those regions designated as planned target volumes (PTV) and an upper bound on the dosage delivered to those regions designated as organs at risk (OAR).Suppose thatSt is either a PTV or an OAR, and suppose that St contains Nt voxels. For each dosage vector b=(b1,…,bI)T, define the equivalent uniform dosage function (EUD function) et(b) by(1.1)et(b)=(1Nt∑i∈St(bi)α)1/α,
where 0<α<1 if St is a PTV, and α>1 if St is an OAR. The function et(b) is convex, for b nonnegative, when St is an OAR and -et(b) is convex, when St is a PTV. The constraints in dosage space take the form(1.2)et(b)≤at,
when St is an OAR, and(1.3)-et(b)≤at,
when St is a PTV. Therefore, we require that b=Ax lie within the intersection of these convex sets. In a summary, we have formulated the constraints in the radiation intensity space RJ and in the dose space RI, respectively, and the two spaces are related by the dose influence matrix A; that is, this problem referred as the multiple-set split feasibility problem (MSSFP) is formulated as follows.(1.4)Findanx∈⋂i=1NXisuchthatAx∈⋂j=1MHj,
which was first investigated by Censor et al. [5]. There are a great deal of literature on the MSSFP; see [5, 7, 8, 18, 19, 22, 23]. In the sequel, there will be involved optimization and variational inequality techniques. For related references, please see [30–42].


### 1.2. Fixed-Point Method

Next, we focus on the multiple-set split feasibility problem (MSSFP) which is to find a pointx* such that(1.5)x*∈C=⋂i=1NCi,Ax*∈Q=⋂j=1MQj,
where N,M≥1 are integers, the Ci(i=1,2,…,N) are closed convex subsets of H1, the Qj(j=1,2,…,M) are closed convex subsets of H2, and A:H1→H2 is a bounded linear operator. Assume that MSSFP is consistent; that is, it is solvable, and S denotes its solution set. The case where N=M=1, called split feasibility problem (SFP), was introduced by Censor and Elfving [43], modeling phase retrieval and other image restoration problems, and further studied by many researchers; see, for instance, [2–4, 6, 9–12, 17, 19–21].We useΓ to denote the solution set of the SFP. Let γ>0 and assume that x*∈Γ. Thus, Ax*∈Q1 which implies the equation (I-PQ1)Ax*=0 which in turn implies the equation γA*(I-PQ1)Ax*=0, hence the fixed point equation (I-γA*(I-PQ1)A)x*=x*. Requiring that x*∈C1, we consider the fixed-point equation(1.6)PC1(I-γA*(I-PQ1)A)x*=x*.
We will see that solutions of the fixed point equation (1.6) are exactly solutions of the SFP. The following proposition is due to Byrne [4] and Xu [2].Proposition 1.1.
Givenx*∈H1. Then  x* solves the SFP if and only if  x* solves the fixed point (1.6).This proposition reminds us that (MSSFP) (1.5) is equivalent to a common fixed-point problem of finitely many nonexpansive mappings, as we show below.Decompose MSSFP intoN subproblems (1≤i≤N):(1.7)xi*∈Ci,Axi*∈Q:=⋂j=1MQj.
For each 1≤i≤N, we define a mapping Ti by(1.8)Tix=PCi(I-γi∇f)x=PCi(I-γi∑j=1MβjA*(I-PQj)A)x,
where f is defined by(1.9)f(x)=12∑j=1Mβj∥Ax-PQjAx∥2,
with βj>0 for all 1≤j≤M. Note that the gradient of ∇f is(1.10)∇f(x)=∑j=1MβjA*(I-PQj)Ax,
which is L-Lipschitz continuous with constant(1.11)L=∑j=1Mβj∥A∥2.
It is known that if 0<γi≤2/L, Ti is nonexpansive. Therefore fixed-point algorithms for nonexpansive mappings can be applied to (MSSFP) (1.5).


### 1.3. Optimization Method

Note thatx* solves the MSSFP implies that x* satisfies two properties:(i)
the distance fromx* to each Ci is zero,(ii)
the distance fromAx* to each Qj is also zero.This motivates us to consider the proximity function(1.12)g(x)=12∑i=1Nαi∥x-PCix∥2+12∑j=1Mβj∥Ax-PQjAx∥2,
where {αi} and {βj} are positive real numbers, and PCi and PQj are the metric projections onto Ci and Qj, respectively.Proposition 1.2.
x* is a solution of MSSFP (1.5) if and only if  g(x*)=0.Sinceg(x)≥0 for all x∈H1, a solution of MSSFP (1.5) is a minimizer of g over any closed convex subset, with minimum value of zero. Note that this proximity function is convex and differentiable with gradient(1.13)∇g(x)=∑i=1Nαi(I-PCi)x+∑j=1MβjA*(I-PQj)Ax,
where A* is the adjoint of A. Since the gradient ∇g(x) is L′-Lipschitz continuous with constant(1.14)L′=∑i=1Nαi+∑j=1Mβj∥A∥2,
we can use the gradient-projection method to solve the minimization problem(1.15)min⁡x∈Ω⁡g(x),
where Ω is a closed convex subset of H1 whose intersection with the solution set of MSSFP is nonempty, and get a solution of the so-called constrained multiple-set split feasibility problem (CMSSFP)(1.16)x*∈Ωsuchthatx*solves(1.5).
In this paper, we will review and report the recent progresses on the fixed-point and optimization methods for solving the MSSFP.


## 1.1. The Multiple-Set Split Feasibility Problem Model

The intensity-modulated radiation therapy (IMRT) has received a great deal of attention recently; for related works, please refer to [1–29]. In intensity modulated radiation therapy, beamlets of radiation with different intensities are transmitted into the body of the patient. Each voxel within the patient will then absorb a certain dose of radiation from each beamlet. The goal of IMRT is to direct a sufficient dosage to those regions requiring the radiation, those that are designated planned target volumes (PTVs), while limiting the dosage received by the other regions, the so-called organs at risk (OAR). The forward problem is to calculate the radiation dose absorbed in the irradiated tissue based on a given distribution of the beamlet intensities. The inverse problem is to find a distribution of beamlet intensities, the radiation intensity map, which will result in a clinically acceptable dose distribution. One important constraint is that the radiation intensity map must be implementable; that is, it is physically possible to produce such an intensity map, given the machine’s design. There will be limits on the change in intensity between two adjacent beamlets, for example.The equivalent uniform dose (EUD) for tumors is the biologically equivalent dose which, if given uniformly, will lead to the same cell kill within the tumor volume as the actual nonuniform dose. Constraints on the EUD received by each voxel of the body are described in dose space, the space of vectors whose entries are the doses received at each voxel. Constraints on the deliverable radiation intensities of the beamlets are best described in intensity space, the space of vectors whose entries are the intensity levels associated with each of the beamlets. The constraints in dose space will be upper bounds on the dosage received by the OAR and lower bounds on the dosage received by the PTV. The constraints in intensity space are limits on the complexity of the intensity map and on the delivery time, and, obviously, that the intensities be nonnegative. Because the constraints operate in two different domains, it is convenient to formulate the problem using these two domains. This leads to a split feasibility problem.The split feasibility problem (SFP) is to find anx in a given closed convex subset C of RJ such that Ax is in a given closed convex subset Q of RI, where A is a given real I by J matrix. Because the constraints are best described in terms of several sets in dose space and several sets in intensity space, the SFP model needs to be expanded into the multiple-set SFP. It is not uncommon to find that, once the various constraints have been specified, there is no intensity map that satisfies them all. In such cases, it is desirable to find an intensity map that comes as close as possible to satisfying all the constraints. One way to do this, as we will see, is to minimize a proximity function.Fori=1,…,I and j=1,…,J, let bi≥0 be the dose absorbed by the ith voxel of the patient’s body, xj≥0 the intensity of the jth beamlet of radiation, and Aij≥0 the dose absorbed at the ith voxel due to a unit intensity of radiation at the jth beamlet. The nonnegative matrix A with entries Aij is the dose influence matrix. Let us assume that we have M constraints in the dose space and N constraints in the intensity space. Let Hm be the set of dose vectors that fulfill the mth dose constraint, and let Xn be the set of beamlet intensity vectors that fulfill the nth intensity constraint.In intensity space, we have the obvious constraints thatxj≥0. In addition, there are implementation constraints; the available treatment machine will impose its own requirements, such as a limit on the difference in intensities between adjacent beamlets. In dosage space, there will be a lower bound on the dosage delivered to those regions designated as planned target volumes (PTV) and an upper bound on the dosage delivered to those regions designated as organs at risk (OAR).Suppose thatSt is either a PTV or an OAR, and suppose that St contains Nt voxels. For each dosage vector b=(b1,…,bI)T, define the equivalent uniform dosage function (EUD function) et(b) by(1.1)et(b)=(1Nt∑i∈St(bi)α)1/α,
where 0<α<1 if St is a PTV, and α>1 if St is an OAR. The function et(b) is convex, for b nonnegative, when St is an OAR and -et(b) is convex, when St is a PTV. The constraints in dosage space take the form(1.2)et(b)≤at,
when St is an OAR, and(1.3)-et(b)≤at,
when St is a PTV. Therefore, we require that b=Ax lie within the intersection of these convex sets. In a summary, we have formulated the constraints in the radiation intensity space RJ and in the dose space RI, respectively, and the two spaces are related by the dose influence matrix A; that is, this problem referred as the multiple-set split feasibility problem (MSSFP) is formulated as follows.(1.4)Findanx∈⋂i=1NXisuchthatAx∈⋂j=1MHj,
which was first investigated by Censor et al. [5]. There are a great deal of literature on the MSSFP; see [5, 7, 8, 18, 19, 22, 23]. In the sequel, there will be involved optimization and variational inequality techniques. For related references, please see [30–42].


## 1.2. Fixed-Point Method

Next, we focus on the multiple-set split feasibility problem (MSSFP) which is to find a pointx* such that(1.5)x*∈C=⋂i=1NCi,Ax*∈Q=⋂j=1MQj,
where N,M≥1 are integers, the Ci(i=1,2,…,N) are closed convex subsets of H1, the Qj(j=1,2,…,M) are closed convex subsets of H2, and A:H1→H2 is a bounded linear operator. Assume that MSSFP is consistent; that is, it is solvable, and S denotes its solution set. The case where N=M=1, called split feasibility problem (SFP), was introduced by Censor and Elfving [43], modeling phase retrieval and other image restoration problems, and further studied by many researchers; see, for instance, [2–4, 6, 9–12, 17, 19–21].We useΓ to denote the solution set of the SFP. Let γ>0 and assume that x*∈Γ. Thus, Ax*∈Q1 which implies the equation (I-PQ1)Ax*=0 which in turn implies the equation γA*(I-PQ1)Ax*=0, hence the fixed point equation (I-γA*(I-PQ1)A)x*=x*. Requiring that x*∈C1, we consider the fixed-point equation(1.6)PC1(I-γA*(I-PQ1)A)x*=x*.
We will see that solutions of the fixed point equation (1.6) are exactly solutions of the SFP. The following proposition is due to Byrne [4] and Xu [2].Proposition 1.1.
Givenx*∈H1. Then  x* solves the SFP if and only if  x* solves the fixed point (1.6).This proposition reminds us that (MSSFP) (1.5) is equivalent to a common fixed-point problem of finitely many nonexpansive mappings, as we show below.Decompose MSSFP intoN subproblems (1≤i≤N):(1.7)xi*∈Ci,Axi*∈Q:=⋂j=1MQj.
For each 1≤i≤N, we define a mapping Ti by(1.8)Tix=PCi(I-γi∇f)x=PCi(I-γi∑j=1MβjA*(I-PQj)A)x,
where f is defined by(1.9)f(x)=12∑j=1Mβj∥Ax-PQjAx∥2,
with βj>0 for all 1≤j≤M. Note that the gradient of ∇f is(1.10)∇f(x)=∑j=1MβjA*(I-PQj)Ax,
which is L-Lipschitz continuous with constant(1.11)L=∑j=1Mβj∥A∥2.
It is known that if 0<γi≤2/L, Ti is nonexpansive. Therefore fixed-point algorithms for nonexpansive mappings can be applied to (MSSFP) (1.5).


## 1.3. Optimization Method

Note thatx* solves the MSSFP implies that x* satisfies two properties:(i)
the distance fromx* to each Ci is zero,(ii)
the distance fromAx* to each Qj is also zero.This motivates us to consider the proximity function(1.12)g(x)=12∑i=1Nαi∥x-PCix∥2+12∑j=1Mβj∥Ax-PQjAx∥2,
where {αi} and {βj} are positive real numbers, and PCi and PQj are the metric projections onto Ci and Qj, respectively.Proposition 1.2.
x* is a solution of MSSFP (1.5) if and only if  g(x*)=0.Sinceg(x)≥0 for all x∈H1, a solution of MSSFP (1.5) is a minimizer of g over any closed convex subset, with minimum value of zero. Note that this proximity function is convex and differentiable with gradient(1.13)∇g(x)=∑i=1Nαi(I-PCi)x+∑j=1MβjA*(I-PQj)Ax,
where A* is the adjoint of A. Since the gradient ∇g(x) is L′-Lipschitz continuous with constant(1.14)L′=∑i=1Nαi+∑j=1Mβj∥A∥2,
we can use the gradient-projection method to solve the minimization problem(1.15)min⁡x∈Ω⁡g(x),
where Ω is a closed convex subset of H1 whose intersection with the solution set of MSSFP is nonempty, and get a solution of the so-called constrained multiple-set split feasibility problem (CMSSFP)(1.16)x*∈Ωsuchthatx*solves(1.5).
In this paper, we will review and report the recent progresses on the fixed-point and optimization methods for solving the MSSFP.


## 2. Some Concepts and Tools

AssumeH is a Hilbert space and C is a nonempty closed convex subset of H. The (nearest point or metric) projection, denoted PC, from H onto C assigns for each x∈H the unique point PCx∈C in such a way that(2.1)‖x-PCx‖=inf⁡{‖x-y‖:y∈C}.Proposition 2.1.
Basic properties of projections are(i)
〈x-PCx,y-PCx〉≤0 for all x∈H and y∈C;(ii)
∥x-PCx∥2≤∥x-y∥2-∥y-PCx∥2  for all x∈H and y∈C;(iii)
〈x-y,PCx-PCy〉≥∥PCx-PCy∥2  for all x,y∈H, and equality holds if and only if x-y=PCx-PCy. In particular, PC is nonexpansive; that is,
(2.2)‖PCx-PCy‖≤‖x-y‖,
for all x,y∈H;(iv)
ifC is a closed subspace of H, then PC is the orthogonal projection from H onto C:
(2.3)x-PCx⊥C,or〈x-PCx,y〉=0∀x∈H,y∈C.Definition 2.2.
The operator(2.4)Qλ=(1-λ)I+λPC
is called a relaxed projection, where λ∈(0,2) and I is the identity operator on H.
A mappingR:H→H is said to be an averaged mapping if R can be written as an average of the identity I and a nonexpansive mapping T:
(2.5)R=(1-α)I+αT,
where α is a number in (0,1) and T:H→H is nonexpansive.Proposition2.1(iii) is equivalent to saying that the operator S=2PC-I is nonexpansive. Indeed, we have(2.6)∥Sx-Sy∥2=∥2(PCx-PCy)-(x-y)∥2=4∥PCx-PCy∥2-4〈PCx-PCy,x-y〉+∥x-y∥2≤∥x-y∥2.
Consequently, a projection can be written as the mean average of a nonexpansive mapping and the identity:(2.7)PC=I+S2.
Thus projections are averaged maps with α=1/2. Also relaxed projections are averaged.Proposition 2.3.
LetT:H→H be a nonexpansive mapping and  R=(1-α)I+αT an averaged map for some  α∈(0,1). Assume  T has a bounded orbit. Then, one has the following. (1)
R is asymptotically regular; that is,
(2.8)lim⁡n→∞‖Rn+1x-Rnx‖=0,
for all x∈H.(2)
For anyx∈H, the sequence {Rnx} converges weakly to a fixed point of T.Definition 2.4.
LetA be an operator with domain D(A) and range R(A) in H.(i)
A is monotone if for all x,y∈D(A),
(2.9)〈Ax-Ay,x-y〉≥0.(ii)
Given a numberν>0. A is said to be ν-inverse strongly monotone (ν-ism) (or cocoercive) if(2.10)〈Ax-Ay,x-y〉≥ν∥Ax-Ay∥2,x,y∈H.
It is easily seen that a projectionPC is a 1-ism.Proposition 2.5.
GivenT:H→H, let  V=I-T be the complement of  T. Given also  S:H→H, then one has the following.(i)
T is nonexpansive if and only if V is 1/2-ism.(ii)
IfS is ν-ism, then, for γ>0, γS is ν/γ-ism.(iii)
S is averaged if and only if the complement I-S is ν-ism for some ν>1/2.
The next proposition includes the basic properties of averaged mappings.Proposition 2.6.
Given operatorsS,T,V:H→H, then one has the following. (i)
IfS=(1-α)T+αV for some α∈(0,1) and if T is averaged and V is nonexpansive, then S is averaged.(ii)
S is firmly nonexpansive if and only if the complement I-S is firmly nonexpansive. If S is firmly nonexpansive, then S is averaged.(iii)
IfS=(1-α)T+αV for some α∈(0,1), T is firmly nonexpansive and V is nonexpansive, then S is averaged.(iv)
IfS and T are both averaged, then the product (composite) ST is averaged.(v)
IfS and T are both averaged and if S and T have a common fixed point, then(2.11)Fix⁡(S)⋂Fix⁡(T)=Fix⁡(ST).Proposition 2.7.
Consider the variational inequality problem(VI).
(2.12)Findapointx∈Csuchthat〈Ax*,x-x*〉≥0,∀x∈C,
where C is a closed convex subset of a Hilbert space H and A is a monotone operator on H. Assume that VI (2.12) has a solution and A is ν-ism. Then for 0<γ<2ν, the sequence {xn} generated by the algorithm
(2.13)xn+1=PC(xn-γAxn),n≥0,
converges weakly to a solution of the VI (2.12).An immediate consequence of Proposition2.7 is the convergence of the gradient-projection algorithm.Proposition 2.8.
Letf:H→R be a continuously differentiable function such that the gradient  ∇f is Lipschitz continuous:
(2.14)‖∇f(x)-∇f(y)‖≤L‖x-y‖,x,y∈H.
Assume that the minimization problem
(2.15)min⁡x∈C⁡f(x)
is consistent, where C is a closed convex subset of H. Then, for 0<γ<2/L, the sequence {xn} generated by the gradient-projection algorithm
(2.16)xn+1=PC(xn-γ∇f(xn))
converges weakly to a solution of (2.15).


## 3. Iterative Methods

In this section, we will review and report the iterative methods for solving MSSFP (1.5) in the literature.It is not hard to see that the solution setSi of the subproblem (1.7) coincides with Fix(Ti), and the solution set S of MSSFP (1.5) coincides with the common fixed-point set of the mappings Ti. Further, we have (see [9, 18])(3.1)S=⋂i=1NFix⁡(Ti)=Fix⁡(TN…T2T1).
By using the fact (3.1), we obtain the corresponding algorithms and the convergence theorems for the MSSFP.Algorithm 3.1.
The Picard iterations are(3.2)xn+1=TN…T1xn=PCN(I-γ∑j=1MβjA*(I-PQj)A)…PC1(I-γ∑j=1MβjA*(I-PQj)A)xn,n≥0.Theorem 3.2 (see [8]).
Assume that the MSSFP (1.5) is consistent. Let  {xn} be the sequence generated by the Algorithm 3.1, where  0<γ<2/L with  L given by (1.11). Then  {xn} converges weakly to a solution of the MSSFP (1.5).Algorithm 3.3.
Parallel iterations are(3.3)xn+1=∑i=1NλiTi=∑i=1NλiPCi(I-γ∑j=1MβjA*(I-PQj)A)xn,n≥0,where λi>0 for all i such that ∑i=1Nλi=1, and 0<γ<2/L with L given by (1.11).Theorem 3.4 (see [8]).
Assume that the MSSFP (1.5) is consistent. Then the sequence  {xn} generated by the Algorithm 3.3 converges weakly to a solution of the MSSFP (1.5).Algorithm 3.5.
Cyclic iterations are(3.4)xn+1=T[n+1]xn=PC[n+1](I-γ∑j=1MβjA*(I-PQj)A)xn,n≥0,
where T[n]=TnmodN with the mod function taking values in {1,2,…,N}.Theorem 3.6 (see [8]).
Assume that the MSSFP (1.5) is consistent. Let  {xn} be the sequence generated by the Algorithm 3.5, where  0<γ<2/L with  L given by (1.11). Then  {xn} converges weakly to a solution of the MSSFP (1.5).Note that the MSSFP (1.5) can be viewed as a special case of the convex feasibility problem of finding x* such that(3.5)x*∈⋂i=1pCi.
In fact, (1.5) can be rewritten as(3.6)x*∈⋂i=1N+MCi,
where CN+i:={x∈H1:A-1x∈Qj},1≤j≤M.However, the methodologies for studying the MSSFP (1.5) are actually different from those for the convex feasibility problem in order to avoid usage of the inverse A-1. In other words, the methods for solving the convex feasibility problem may not apply to solve the MSSFP (1.5) straightforwardly without involving the inverse A-1. The CQ algorithm of Byrne [1] is such an example where only the operator A (not the inverse A-1) is relevant.Since every closed convex subset of a Hilbert space is the fixed point set of its associating projection, the convex feasibility problem becomes a special case of the common fixed-point problem of finding a pointx* with the property(3.7)x*∈⋂i=1MFix⁡(Ti).
Similarly, the MSSFP (1.5) becomes a special case of the split common fixed-point problem [19] of finding a point x* with the property(3.8)x*∈⋂i=1NFix⁡(Ui),Ax*∈⋂j=1MFix⁡(Tj),
where Ui:H1→H1(i=1,2,…,N) and Tj:H2→H2(j=1,2,…,M) are nonlinear operators. By using these facts, recently, Wang and Xu [17] presented another cyclic iteration as follows.Algorithm 3.7 (cyclic iterations).
Take an initial guessx0∈H1, choose γ∈(0,2/L) and define a sequence {xn} by the iterative procedure:
(3.9)xn+1=PC[n](xn+γA*(PQ[n]-I)Axn),n≥0.Theorem 3.8 (see [17]).
The sequence{xn}, generated by Algorithm 3.7, converges weakly to a solution of MSSFP (1.5) whenever its solution set is nonempty.Since MSSFP (1.5) is equivalent to the minimization problem (1.15), we have the following gradient-projection algorithm.Algorithm 3.9.
Gradient-projection algorithmis(3.10)xn+1=PΩ(xn-γ∇g(xn))=PΩ(xn-γ(∑i=1Nαi(I-PCi)xn+∑j=1MβjA*(I-PQj)Axn)),n≥0.
Censor et al. [5] proved in finite-dimensional Hilbert spaces that Algorithm 3.9 converges to a solution of the MSSFP (1.5) in the consistent case. Below is a version of this convergence in infinite-dimensional Hilbert spaces.Theorem 3.10 (see [8]).
Assume that0<γ<2/L′, where  L′ is given by (1.14). The sequence  {xn} generated by the Algorithm 3.9 weakly converges to a point  z which is a solution of the MSSFP (1.5) in the consistent case and a minimizer of the function  p over  Ω in the inconsistent case.Consequently, Lopez et al. [18] considered a variant version of Algorithm 3.9 to solve (1.16).Algorithm 3.11.
Gradient-projection algorithm is(3.11)xn+1=PΩ(xn-γn∇g(xn))=PΩ(xn-γn(∑i=1Nαi(I-PCi)xn+∑j=1MβjA*(I-PQj)Axn)),n≥0.Theorem 3.12 (see [18]).
Assume that0<liminfn→∞γn≤limsupn→∞γn<2/L′, where  L′ is given by (1.14). The sequence  {xn} generated by the Algorithm 3.11 weakly converges to a solution of (1.16).Remark 3.13.
It is obvious that Theorem3.12 contains Theorem 3.10 as a special case.Perturbation Techniques
Consider the consistent (1.16) and denote by S′ its nonempty solution set. As pointed in the previous, the projection PC, where C is a closed convex subset of H, may bring difficulties in computing it, unless C has a simple form (e.g., a closed ball or a half-space). Therefore some perturbed methods in order to avoid this inconvenience are presented.
We can use subdifferentials when{Ci}, {Qj}, and Ω are level sets of convex functionals. Consider
(3.12)Ci={x∈H1:ci(x)≤0},Qj={y∈H2:qj(y)≤0},Ω={x∈H1:ω(x)≤0},
where ci,ω:H1→R and qj:H2→R are convex functionals. We iteratively define a sequence {xn} as follows.Algorithm 3.14.
The initialx0∈H1 is arbitrary; once xn has been defined, we define the (n+1)th iterate xn+1 by
(3.13)xn+1=PΩn(xn-γn(∑i=1Nαi(I-PCin)xn+∑j=1MβjA*(I-PQjn)Axn)),n≥0,
where
(3.14)Ωn={x∈H1:ω(xn)+〈ξn,x-xn〉≤0},(ξn∈∂ω(xn)),Cin={x∈H1:ci(xn)+〈ξin,x-xn〉≤0},(ξin∈∂ci(xn)),Qjn={y∈H2:qj(Axn)+〈ηj′y-Axn〉≤0},(ηjn∈∂qj(Axn)).Theorem 3.15 (see [18]).
Assume that each of the functions{ci}i=1N,  ω, and  {qj}j=1M satisfies the property: it is bounded on every bounded subset of  H1 and  H2, respectively. (Note that this condition is automatically satisfied in a finite-dimensional Hilbert space.) Then the sequence  {xn} generated by Algorithm 3.14 converges weakly to a solution of (1.16), provided that the sequence  {γn} satisfies
(3.15)0<liminf⁡n→∞γn≤limsup⁡n→∞γn<2L′,
where the constant L′ is given by (1.14).Now consider general perturbation techniques in the direction of the approaches studied in [20–22, 44]. These techniques consist on taking approximate sets which involve the ρ-distance between two closed convex sets A and B of a Hilbert space:(3.16)dρ(A,B)=sup⁡{‖PAx-PBx‖:x∈H1,‖x‖≤ρ}.
Let {Ωn}, {Cin}, and {Qjn} be closed convex sets which are viewed as perturbations for the closed convex sets Ω, {Ci}, and {Qj}, respectively. Define function gn by(3.17)gn(x)=12∑i=1Nαi∥x-PCinx∥2+12∑j=1Mβj∥Ax-PQjnAx∥2.
The gradient ∇gn of gn is(3.18)∇gn(x)=∑i=1Nαi(I-PCin)x+∑j=1MβjA*(I-PQjn)Ax.
It is clear that ∇gn is Lipschitz continuous with the Lipschitz constant L′ given by (1.14).Algorithm 3.16.
Let an initial guessx0∈H1 be given, and let {xn} be generated by the  Krasnosel’skii-Mann iterative algorithm:
(3.19)xn+1=(1-tn)xn+tnPΩn(I-γ∇gn)xn=(1-tn)xn+tnPΩn(xn-γ(∑i=1Nαi(I-PCin)xn+∑j=1MβjA*(I-PQjn)Axn)),n≥0.In [8], Xu proved the following result.Theorem 3.17 (see [8]).
Assume that the following conditions are satisfied.(i)
0<γ<2/L′.(ii)
∑n=1∞tn(1-tn)=∞.(iii)
For eachρ>0, 1≤i≤N, and 1≤j≤M, there hold ∑n=1∞tndρ(Ωn,Ω)<∞, ∑n=0∞tndρ(Cin,Ci)<∞, and ∑n=0∞tndρ(Qjn,Qj)<∞.Then the sequence {xn} generated by Algorithm 3.16 converges weakly to a solution of MSSFP (1.5).Lopez et al. [18] further obtained a general result by relaxing condition (ii).Theorem 3.18 (see [18]).
Assume that the following conditions are satisfied.(i)
0<γ<2/L′.(ii)
tn∈[0,4/(2+γL′)] for all n (note that tn may be larger than one since 0<γ<2/L′) and
(3.20)∑n=0∞tn(42+γL′-tn)=∞.(iii)
For eachρ>0, 1≤i≤N, and 1≤j≤M, there hold ∑n=1∞tndρ(Ωn,Ω)<∞,∑n=0∞tndρ(Cin,Ci)<∞, and ∑n=0∞tndρ(Qjn,Qj)<∞.
Then the sequence{xn} generated by Algorithm 3.16 converges weakly to a solution of (1.16).Corollary 3.19.
Assume that the following conditions are satisfied.(i)
0<γ<2/L′.(ii)
tn∈[0,4/(2+γL′)] for all n (note that tn may be larger than one since 0<γ<2/L′) and(3.21)∑n=0∞tn(42+γL′-tn)=∞.
Then the sequence {xn} generated by
(3.22)xn+1=(1-tn)xn+tnPΩ(xn-γ(∑i=1Nαi(I-PCi)xn+∑j=1MβjA*(I-PQj)Axn)),n≥0,
converges weakly to a solution of the MSSFP (1.5).Note that all above algorithms only have weak convergence. Next, we will consider some algorithms with strong convergence.Algorithm 3.20.
The Halpern iterations are(3.23)xn+1=αnu+(1-αn)Tnxn=αnu+(1-αn)PC[n+1](I-γ∑j=1MβjA*(I-PQj)A)xn,n≥0.Theorem 3.21.
Assume that the MSSFP (1.5) is consistent,  0<γ<2/L with  L given by (1.11), and  {αn} satisfies the conditions (for instance,  αn=1/n for all  n≥1) (C1)
limn→∞αn=0,(C2)
∑n=0∞αn=∞,(C3)
∑n=0∞|αn+1-αn|<∞ or limn→∞(αn+1/αn)=1.
Then the sequence{xn} generated by the Algorithm 3.20 converges strongly to a solution of the MSSFP (1.5) that is closest to u from the solution set of the MSSFP (1.5).Next, we consider a perturbation algorithm which has strong convergence.Algorithm 3.22.
Given an initial guessx0∈H1, let {xn} be generated by the perturbed iterative algorithm
(3.24)xn+1=γnu+(1-γn)PΩn(xn-γ(∑i=1Nαi(I-PCin)xn+∑j=1MβjA*(I-PQjn)Axn)),n≥0.Theorem 3.23 (see [18]).
Assume that the following conditions are satisfied.(i)
0<γ<2/L′.(ii)
limn→∞tn=0 and ∑n=0∞tn=∞.(iii)
For eachρ>0, 1≤i≤N, and 1≤j≤M, there hold ∑n=1∞tndρ(Ωn,Ω)<∞, ∑n=0∞tndρ(Cin,Ci)<∞, and ∑n=0∞tndρ(Qjn,Qj)<∞.
Then the sequence{xn} generated by Algorithm 3.22 converges in norm to the solution of (1.16) which is nearest to u.Corollary 3.24.
Assume that the following conditions are satisfied.(i)
0<γ<2/L′.(ii)
limn→∞tn=0 and ∑n=0∞tn=∞.
Then the sequence{xn} generated by
(3.25)xn+1=tnu+tnPΩ(xn-γ(∑i=1Nαi(I-PCi)xn+∑j=1MβjA*(I-PQj)Axn)),n≥0,
converges in norm to a solution of the MSSFP (1.5).Regularized Methods
Consider the following regularization:(3.26)gα(x):=g(x)+12α∥x∥2=12∑i=1Nαi∥x-PCix∥2+12∑j=1Mβj∥Ax-PQjAx∥2+12α∥x∥2,
where α>0 is the regularization parameter. We can compute the gradient ∇gα of gα as
(3.27)∇gα=∑i=1Nαi(I-PCi)+∑j=1MβjA*(I-PQj)A+αI.
It is easily see that ∇gα is Lα-Lipschitz continuous with constant
(3.28)Lα=∑i=1Nαi+∑j=1Mβj∥A∥2+α.
It is known that ∇gα is strongly monotone.
Consider the following regularized minimization problem(3.29)min⁡x∈Ωgα(x),
which has a unique solution denoted by xα.Theorem 3.25.
The stronglimα→0xα exists and equals  x̃, the minimum-norm solution of (1.16).Algorithm 3.26.
Given an initial pointx0∈Ω.Define a sequence {xn} by the iterative algorithm
(3.30)xn+1=PΩ(I-γn∇gαn)xn=PΩ((I-αnγn)xn-γn∑i=1Nαi(I-PCi)xn-γn∑j=1MβjA*(I-PQj)Axn),n≥0.Theorem 3.27 (see [18]).
Assume the sequences{αn} and  {γn} satisfy the conditions: (i)
0<γn<αn/Lαn2 for all (large enough) n;(ii)
αn→0;(iii)
∑n=1∞αnγn=∞;(iv)
(|γn-γn-1|+|αnγn-αn-1γn-1|)/(αnγn)2→0.
Then the sequence{xn} generated by Algorithm 3.26 strongly converges to the minimum norm solution of (1.16).Self-Adaptive Methods
Consider the following constrained minimization problem:(3.31)min⁡{g(x),x∈Ω},
where g(x) is defined as in (1.12) and Ω⊂RN is the same auxiliary simple nonempty closed convex set as in (1.16). This optimization problem is proposed by Censor et al. [5] for solving the constrained MSSFP (1.5) in the finite-dimensional Hilbert spaces. We know that a point x*∈Ω is a stationary point of problem (3.31) if it satisfies
(3.32)〈∇g(x*),x-x*〉≥0,∀x∈Ω.
Thus, from Proposition 2.8, we can use a gradient projection algorithm below to solve the MSSFP which was developed by Censor et al. ([5, 24]):
(3.33)xn+1=PΩ(xn-γ∇g(xn)),
where
(3.34)γ∈(0,2L).
Note that the above method of Censor et al. is the application of the projection method of Goldstein [45] and Levitin and Polyak [46] to the variational inequality problem (3.32), which is among the simplest numerical methods for solving variational inequality problems. Nevertheless, the efficiency of this projection method depends greatly on the choice of the parameter γ. If one chooses a small s to ensure that it satisfies the condition (3.34) such that it guarantees the convergence of the iterative sequence, the recursion leads to slow speed of convergence. On the other hand, if one chooses a large step size to improve the speed of convergence, the generated sequence may not converge. In real applications for solving variational inequality problems, the Lipschitz constant may be difficult to estimate, even if the underlying mapping is linear, the case such as the MSSFP.
To overcome the difficulty in estimating the Lipschitz constant, He et al. [47] developed a self-adaptive method for solving variational inequality problems, where the constant step size γ in the original Goldstein-Levitin-Polyak method is replaced by a sequence of parameters {γn} and {γn} is selected self-adaptively. The numerical results reported in He et al. [47] have shown that the self-adaptive strategy is valid and robust for solving variational inequality problems. The efficiency of their modified algorithm is not affected by the initial choice of the parameter; that is, for any given initial choice γ0, the algorithm can adjust it and finally find a “suitable” one. Thus, there is no need to pay much attention to the choice of the step size as that of the original Goldstein-Levitin-Polyak method. Moreover, the computational burden at each iteration is not much larger than that of the original Goldstein-Levitin-Polyak method. Later, their method is extended to a more flexible self-adaptive rule by Han and Sun [25].
Motivated by the self-adaptive strategy, Zhang et al. [23] proposed the following method for solving the MSSFP by using variable step sizes, instead of the fixed step sizes as in Censor et al. [5, 24].Algorithm 3.28.
(S1) Given a nonnegative sequenceτnwith∑n=0∞τn<∞,δ∈(0,1),μ∈(0,1),ρ∈(0,1),ϵ>0,β0>0,and arbitrary initial point x0, set γ0=β0 and n=0.(S2)
Find the smallest nonnegative integerln such thatβn+1=μlkγk and(3.35)xn+1=PΩ(xn-βn+1∇g(xn)),
which satisfies
(3.36)βn+1∥∇g(xn)-∇g(xn+1)∥2≤(2-δ)〈xn-xn+1,∇g(xn)-∇g(xn+1)〉.(S3)
If(3.37)βn+1∥∇g(xn)-∇g(xn+1)∥2≤ρ〈xn-xn+1,∇g(xn)-∇g(xn+1)〉,
then set γn+1=(1+τn+1)βn+1; otherwise, set γn+1=βn+1.(S4)
If∥e(xn,βn)∥≤ϵ, stop; otherwise, setn:=n+1 and go to (S2).Theorem 3.29 (see [23]).
The proposed Algorithm3.28 is globally convergent.Remark 3.30.
This new method is a modification of the projection method proposed by Goldstein [45] and Levitin and Polyak [46], where the constant step size β in their original method is replaced by an automatically selected one, βk, per iteration. This is very important, since it helps us avoid the difficult task of selecting a “suitable” step size.The following self-adaptive projection method was introduced by Zhao and Yang [7], which was adopted by using the Armijo-like searches to solve the MSSFP.Algorithm 3.31.
Given constantsβ>0,σ∈(0,1),γ∈(0,1), let x0 be arbitrary. For n=0,1,…,calculate
(3.38)xn+1=PΩ(xn-τn∇g(xn)),
where τn=βγln and ln is the smallest nonnegative integer l such that
(3.39)g(PΩ(xn-βγl∇g(xn)))≤g(xn)-σ〈∇g(xn),xn-PΩ(xn-βγl∇g(xn))〉.Algorithm3.31 need not to estimate the Lipschitz constant of ∇g or compute the largest eigenvalue of the matrix ATA, and the step-size τn is chosen so that the objective function g(x) has a sufficient decrease. It is in fact a special case of the standard gradient projection method with the Armijo-like search for solving the constrained optimization problem (3.31).The following convergence result for the gradient projection method with the Armijo-like searches solving the generalized convex optimization problem (3.31) ensures the convergence of Algorithm 3.31.Theorem 3.32.
Letg∈CΩ1 be pseudoconvex and  {xn} be an infinite sequence generated by the gradient projection method with Armijo-like searches. Then, the following conclusions hold: (1)
limn→∞g(xn)=inf{g(x),x∈Ω};(2)
Ω*, which denotes the set of the optimal solutions to (3.31), is nonempty if and only if there exists at least one limit point of {xn}. In this case, {xn} converges to a solution of (3.31).However, we find that, in each iteration step of Algorithm3.31, it costs a large amount of work to compute the orthogonal projections PCi and PQj. In what follows, we consider the case that the projections are not easily calculated, and we consider a relaxed self-adaptive projection method for solving the MSSFP. In detail, the MSSFP and the convex sets Ci and Qj in this part should satisfy the following assumptions.(1)
The solution set of the constrained MSSFP is nonempty.(2)
The setsCi, i=1,2,…,t, are given by
(3.40)Ci={x∈RN∣ci(x)≤0},
where ci:RN→R are convex functions. The sets Qj, j=1,2,…,r are given by
(3.41)Qj={y∈RM∣qj(y)≤0},
where qj:RM→R are convex functions.(3)
For anyx∈RN, at least one subgradient ξ∈∂ci(x) can be calculated, where ∂ci(x) is a generalized gradient, called subdifferential of ci(x) at x, and it is defined as follows:(3.42)∂ci(x)={ξi∈RN∣ci(z)≥ci(x)+〈ξi,z-x〉∀z∈RN}.
For any y∈RM, at least one subgradient ηj∈∂qj(y) can be calculated, where ∂qj(y) is a generalized gradient, called subdifferential of qj(y) at y and is defined as follows:(3.43)∂qj(y)={ηj∈RM∣qj(u)≥qj(y)+〈ηj,u-y〉∀u∈RM}.
In the kth iteration, let(3.44)Cin={x∈RN∣ci(xn)+〈ξin,x-xn〉≤0},
where ξin is an element in ∂ci(xn):(3.45)Qjn={y∈RM∣qj(Axn)+〈ηjn,y-Axn〉≤0},
where ηjn is an element in ∂qj(Axn).Define(3.46)gn(x):=12∑i=1tαi∥x-PCin(x)∥2+12∑j=1rβj∥Ax-PQjn(Ax)∥2.
Obviously,(3.47)∇gn(x)=∑i=1tαi(I-PCin)x+∑j=1rβjAT(Ax-PQjn(Ax)).Algorithm 3.33.
Givenγ>0,ρ∈(0,1),μ∈(0,1) let x0 be arbitrary. For n=0,1,2,…,compute
(3.48)x̅n=PΩ(xn-τn∇gn(xn)),
where τn=γρln and ln is the smallest nonnegative integer l such that
(3.49)‖∇gn(xn)-∇gn(x̅n)‖≤μ∥xn-x̅n∥τn.
Set
(3.50)xn+1=PΩ(xn-τn∇gn(x̅n)).Theorem 3.34 (see [7]).
The sequence{xn} generated by Algorithm 3.33 converges to a solution of the MSSFP.



---

*Source: 927530-2012-04-24.xml*