Simultaneous and semi-alternating projection algorithms for solving split equality problems.

In this article, we first introduce two simultaneous projection algorithms for solving the split equality problem by using a new choice of the stepsize, and then propose two semi-alternating projection algorithms. The weak convergence of the proposed algorithms is analyzed under standard conditions. As applications, we extend the results to solve the split feasibility problem. Finally, a numerical example is presented to illustrate the efficiency and advantage of the proposed algorithms.

Introduction

Let , and be real Hilbert spaces, let and be two nonempty closed convex sets, and let and be two bounded linear operators.

In this article, we consider the classical split equality problem (SEP), which was first introduced by Moudafi [1]. The SEP can mathematically be formulated as follows: Throughout this paper, assume that SEP (1) is consistent and denote by the solution of SEP (1). Then Γ is closed, convex and nonempty.

The split equality problem (1) is actually an optimization problem with weak coupling in the constraint (see [1] for details) and its interest covers many situations, for instance, in domain decomposition for PDEs, game theory and intensity-modulated radiation therapy (IMRT). In domain decomposition for PDEs, this equals to the variational form of a PDE in a domain that can be decomposed into two non-overlapping subdomains with a common interface (see, e.g., [2]). In decision sciences, this allows to consider agents who interplay only via some components of their decision variables (see, e.g., [3]). In IMRT, this amounts to envisaging a weak coupling between the vector of doses absorbed in all voxels and that of the radiation intensity (see [4] for further details). Attouch [5] pointed out more applications of the SEP in optimal control theory, surface energy and potential games, whose variational form can be seen as a SEP.

Next we present an example, in which a separable optimization problem can be rewritten as a split equality problem.

Example 1.1

Consider the separable optimization problem with and , where and . Assume that f and g are convex and the solution set of problem (2) is nonempty.

Set and . Then the optimization problem (2) equals to the following split equality problem:

A great deal of literature on algorithms for solving SEP has been published, most of which are projection methods [1–3, 6–11]. Based on the classical projection gradient algorithm, Byrne and Moudafi [12] introduced the following algorithm, which is also called the simultaneous iterative method [13]: where , and are the operator (matrix) norms and (or the largest eigenvalues of and ), respectively. To determine stepsize , one needs first to calculate (or estimate) the operator norms and . In general, it is difficult or even impossible. On the other hand, even if we know the norm of A and B, the algorithm (4) method with fixed stepsize may be slow.

In order to deal with this, the authors [9] introduced a self-adaptive projection algorithm, in which the stepsize is computed by using an Armijo search.

Define the function by and the function by The self-adaptive projection algorithm in [9] is defined as follows.

Algorithm 1.1

Given constants , , and . Let and be taken arbitrarily. For  , compute where is chosen to be the largest satisfying

Compute If then set ; otherwise, set .

In fact, Algorithm 1.1 can be seen as an extension of the classical extragradient method first proposed by Korpelevich [14]. Notice that, in Algorithm 1.1, the stepsize of the prediction (5) and that of the correction (7) are equal. Thus these two steps seem to be ‘symmetric’.

Recently, Chuang and Du [15] proposed the following projection algorithm (which is called the hybrid projected Landweber algorithm).

Algorithm 1.2

Given constants and , let and be taken arbitrarily. For  , compute where is chosen via (6) and (8). Compute next iterates and by where and

Note that Algorithm 1.2 with in (9) can be seen as a special case of Tseng’s method [8, 16]. The projections in the second step of Tseng’s method are made onto two nonempty closed convex sets and , other than C and Q. X and Y can be any sets such that the intersections of X and C (and Y and Q) are nonempty, and they may be taken to have simple structures so that the projections onto them are easy to calculate.

Chuang and Du [15] proved the convergence of Algorithm 1.2 and also presented the convergence property of Algorithm 1.2 as follows: where .

The stepsize in Algorithms 1.1 and 1.2 is obtained through the Armijo search (6). In general, the computational cost of a self-adaptive algorithm is large, since one may need to calculate (5) several times to get the stepsize .

To overcome this difficulty, the authors [17] introduced a projection algorithm for which the stepsize does not depend on the operator norms and , and one can directly compute the stepsize instead of using the Armijo search.

Algorithm 1.3

Choose initial guesses arbitrarily. Assume that the kth iterate , has been constructed and ; then we calculate the th iterate via the formula where the stepsize is chosen in such a way that where . If , then is a solution of SEP (1) and the iterative process stops; otherwise, we set and go onto (12) to evaluate the next iterate .

Note that the choice in (13) of the stepsize is independent of the norms and .

Polyak [18, 19] first introduced the inertial extrapolation algorithms, which were widely studied as an acceleration process. The authors [20] made an inertial modification for Algorithm 1.3 and introduced the following inertial projection methods for SEP.

Algorithm 1.4

Choose initial guesses , arbitrarily. Compute where and the stepsize is chosen in the same way as (13).

They showed the weak convergence of Algorithm 1.4 under some conditions on the inertial parameter .

In fact, Algorithm 1.4 can be seen as a FISTA (fast iterative shrinkage-thresholding algorithm) introduced by Beck and Teboulle [21] to solve the linear inverse problems, if we take the inertial parameter , where and , , and choose a constant stepsize or choose via a backtracking stepsize rule. A shortcoming of the method of Beck and Teboulle is that they could not prove the convergence of the iterative sequence . Chambolle and Dossal [22] improved the choice of the inertial parameter, took , where , and presented the convergence of the iterative sequence .

In this paper, inspired by the work in [17, 23, 24], we introduce two simultaneous projection algorithms by improving the stepsizes and of the second step (7) and (9) in Algorithms 1.1 and 1.2, respectively. We also present two alternating projection algorithms, in which we take an alternating technique in the first step.

The structure of the paper is as follows. In the next section, we present some concepts and lemmas which will be used in the main results. In Section 3, two classes of projection algorithms are provided and their weak convergence is analyzed. In Section 4, we extend the results to the split feasibility problem. In the final section, some numerical results are provided, which show the advantages of the proposed algorithms.

Preliminaries

Let H be a real Hilbert space with the inner product and the induced norm , and let D be a nonempty, closed and convex subset of H. We write to indicate that the sequence converges weakly to x and to indicate that the sequence converges strongly to x. Given a sequence , denote by its weak ω-limit set, that is, any such that there exists a subsequence of which converges weakly to x.

In this paper, an important tool of our work is the projection. Let H be a real Hilbert space and C be a closed convex subset of H. Recall that the projection from H onto C, denoted by , is defined in such a way that, for each , is the unique point in C such that The following two lemmas are useful characterizations of projections.

Lemma 2.1

([25])

Given and . Then if and only if

Lemma 2.2

([25, 26])

For any and , it holds

;

.

Definition 2.1

The normal cone of C at , denoted by , is defined as

Definition 2.2

Let be a point-to-set operator defined on a real Hilbert space H. The operator A is called a maximal monotone operator if A is monotone, i.e., and the graph of A, is not properly contained in the graph of any other monotone operator.

It is clear [27, Theorem 3] that a monotone mapping A is maximal if and only if, for any , if for all , then it follows that .

Lemma 2.3

([26])

Let D be a nonempty, closed and convex subset of a Hilbert space H. Let be a bounded sequence which satisfies the following properties: Then converges weakly to a point in D.

every limit point of lies in D;

exists for every .

Main results

In this section, we present two classes of projection algorithms and establish their weak convergence under standard conditions.

Simultaneous projection algorithms

Let . Define , and let be the adjoint operator of K, then the original problem (1) can be modified as

Note that if the solution set of (15) is nonempty, it equals the following minimization problem: which is a standard (and a simple) problem from the convex optimization point of view. There are many methods for solving the minimization problem (16) such as the classical projection gradient method. Algorithm (4) (also Algorithm 1.3) is the exact projection gradient method when applied to (16).

Inspired by Cai [24] and Dong et al. [17], we propose two new simultaneous projection algorithms by improving the stepsizes in the second step of Algorithms 1.1 and 1.2.

Algorithm 3.1

Given constants , , and , let be taken arbitrarily.

For  , compute where is chosen to be the largest satisfying

Compute next iterates by or where , and If then set ; otherwise, set .

Remark 3.1

Let . Then we have (see Section 4.4.1 in [28]) It is easy to see Define the function by and the function by

By setting and , Algorithm 3.1 can be rewritten as follows:

For  , compute where is chosen to be the largest satisfying

Compute next iterates and by or where , and If then set ; otherwise, set .

For convenience, we call the projection algorithms which use update forms (19) (or (28)) and (20) (or (29)) Algorithm 3.1(I) and Algorithm 3.1(II), respectively.

Remark 3.2

For Algorithm 3.1, we can get the following conclusions:

The only difference between Algorithm 3.1(II) and Algorithm 1.1 is that they use different stepsizes in the definitions of and .

There are two differences between Algorithm 3.1(I) and Algorithm 1.2. Firstly, the stepsize in (10) of Algorithm 3.1(I) is larger than that in (36) of Algorithm 1.2. Secondly, there are no projections on the second step (19).

Remark 3.3

By the definitions of in (21), the projection equation (17) can be written as So, from Lemma 2.1 we have

Lemma 3.1

The search rule (18) is well defined. Besides , where

Proof

Obviously, . In the latter case, we know that must violate inequality (18), that is, On the other hand, we have Consequently, we get (34). □

Lemma 3.2

Let and be generated by Algorithm 3.1, and let and be given by (21) and (36), respectively. Then we have

By the Cauchy-Schwarz inequality, we have

By using , we have So, we get (37). □

Lemma 3.3

Let and be generated by Algorithm 3.1, and let be given by (21). Then, for all , we have

Take arbitrarily , that is, , and . By setting in (33), we get which implies that It is easy to show that So we have which implies (39). □

Theorem 3.1

Let be generated by Algorithm 3.1(I). If Γ is nonempty, then we have and converges weakly to a solution of SEP (1).

Let , that is, , and . Then, from (39), we have which yields (41). Since , (42) implies that the sequence is decreasing and thus converges. Moreover, is bounded. This implies that From the definition of , Lemmas 3.1 and 3.2, we have which implies and Using (43), we get and By the boundedness of K, we get

Let , then there exists a subsequence of which converges weakly to ẑ. From (44), it follows that the subsequence also converges weakly to ẑ. We will show that ẑ is a solution of SEP (1). The weak convergence of to and lower semicontinuity of the squared norm imply that that is, . By noting that the equality in (17) can be rewritten as and that the graph of the maximal monotone operator is weakly-strongly closed, and by passing to the limit in the last inclusions, we obtain, from (44) and (45), that Hence . Now we can apply Lemma 2.3 to to get that the full sequence converges weakly to a point in Γ. This completes the proof. □

Remark 3.4

By using Remark 3.1, the contraction inequality (41) can be rewritten as follows: It is obvious that the in (31) is larger than that in (10). Let . Comparing (46) and (11), we claim that Algorithm 3.1(I) has a better contraction property than Algorithm 1.2.

Theorem 3.2

Let be generated by Algorithm 3.1(II). Assume that Γ is nonempty. Then we have Furthermore, converges weakly to a solution of SEP (1).

Let , that is, , and . Using Lemma 2.2(ii), we have By setting in (33), we get It holds Substituting (50) in the right-hand side of (49) and using , we obtain From (40), we get So, adding (51) and (52) and using the definition of , we obtain Adding (48) and (53), we obtain (47). Employing arguments which are similar to those used in the proof of Theorem 3.1, we obtain that the whole sequence weakly converges to a solution of SEP (1), which completes proof. □

Remark 3.5

Comparing inequalities (41) and (47), we conclude that Algorithm 3.1(II) seems to have a better contraction property than Algorithm 3.1(I) since is closer to than when is the same.

Semi-alternating projection algorithms

Inspired by Algorithm 2.2 in [15] and based on Algorithm 3.1, we present two semi-alternating projection algorithms, whose name comes from an alternating technique taken in the first step.

Algorithm 3.2

Given constants , , and , let and be taken arbitrarily.

For  , compute where is chosen to be the largest satisfying Compute next iterates and by or where , and If then set ; otherwise, set .

For convenience, we call the projection algorithms which use update forms (56) and (57) Algorithm 3.2(I) and Algorithm 3.2(II), respectively.

Remark 3.6

By the definitions of and in (58), the projection equation (54) can be written as So, from Lemma 2.1 we have

Lemma 3.4

The search rule (55) is well defined. Besides , where

Obviously, . In the latter case, we know that must violate inequality (55), that is, On the other hand, we have So, we get (62). □

Lemma 3.5

Let and be generated by Algorithm 3.2 , and let , and be given by (58) and (59), respectively. Then we have

By the Cauchy-Schwarz inequality, we have It holds So, we obtain From the definition of F and G, we have Since by (66), we get Combining (65) and (67), we complete the proof. □

Lemma 3.6

Let and be generated by Algorithm 3.2, and let and be given by (58). Then, for all , we have

Take arbitrarily , that is, , and . By setting in (61), we get which implies that By the definition of F and G, we have and So, we complete the proof. □

Theorem 3.3

Let be generated by Algorithm 3.2(I). If Γ is nonempty, then we have and converges weakly to a solution of SEP (1).

Let , that is, , , and . Then we have Similarly, we get Adding the above inequalities and using Lemma 3.6, we have which yields (69). Since , (70) implies that the sequence is decreasing and thus converges. Moreover, and are bounded. This implies that From the definition of , Lemmas 3.4 and 3.5, we have which implies and Using (71), we get and Hence, we get

Let , then there exist two subsequences and of and which converge weakly to x̂ and ŷ, respectively. From (72), it follows that and also converge weakly to x̂ and ŷ, respectively. We will show that is a solution of SEP (1). The weak convergence of to and the lower semicontinuity of the squared norm imply that that is, . By noting that the two equalities in (54) can be rewritten as and that the graphs of the maximal monotone operators and are weakly-strongly closed, and by passing to the limit in the last inclusions, we obtain, from (72) and (73), that Hence .

Now we can apply Lemma 2.3 to to get that the full sequence converges weakly to a point in Γ. This completes the proof. □

Remark 3.7

Employing arguments which are similar to those used in Remark 3.4, comparing (69) and (2.48) in [15], we conclude that Algorithm 3.2(I) has a better contraction property than the hybrid alternating CQ-algorithm in [15].

Theorem 3.4

Let be generated by Algorithm 3.2(II). If Γ is nonempty, then we have and converges weakly to a solution of SEP (1).

Let , that is, , , and . Using Lemma 2.2(ii), we have Similarly, we get Adding the above inequalities, we obtain By setting in (61), we get It holds Similarly, we get Substituting (77) and (78) in the right-hand side of (76) and using and , we obtain From (68), we have So, adding (79) and (80) and using the definition of , we obtain Adding (75) and (81), we obtain (74). Employing arguments which are similar to those used in the proof of Theorem 3.3, we obtain that the whole sequence weakly converges to a solution of SEP (1), which completes proof. □

Applications

The split feasibility problem (SFP) formulated as follows: was originally introduced in Censor and Elfving [29]. The SFP 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 has a variety of specific applications in real world, such as medical care, image reconstruction and signal processing (see [30-33] for details).

In fact, the SEP is equivalent to the SFP. Firstly, observe that the equality in (1) equals to So, if we define the linear and bounded operator , then the SEP becomes a special case of the SFP with the operator L (e.g., see [34, 35]).

On the other hand, if and , then the split equality problem (1) reduces to the split feasibility problem.

Based on this equivalence, we can construct iterative algorithms for the SEP by using the algorithms for the SFP if the operator is easily computed. We also can extend the algorithms for the SEP to the SFP.

Next, we present an algorithm for the SFP based on Algorithm 3.2.

Define the function by and the function by

Algorithm 4.1

Given constants , , and , let and be taken arbitrarily.

For  , compute where is chosen to be the largest satisfying Compute next iterates and by or where , and If then set ; otherwise, set .

Using Theorem 3.3, we get the convergence of Algorithm 4.1.

Theorem 4.1

Let be generated by Algorithm 4.1. If the set of solutions of the SFP is nonempty, then converges weakly to a solution of SFP (82).

Remark 4.1

Similarly, it is easy to extend Algorithm 3.1 and Theorems 3.1 and 3.2 to the SFP. Here we omit it.

Numerical examples

In this section, we use the numerical example in [8] to demonstrate the efficiency and advantage of Algorithms 3.1 and 3.2 by comparing them with Algorithms 1.1, 1.2 and 1.3.

We denote the vector with all elements 0 by , and the vector with all elements 1 by in what follows. In the numerical results listed in the following table, ‘Iter.’ and ‘Sec.’ denote the number of iterations and the cpu time in seconds, respectively. For Algorithms 1.1, 1.2, 3.1 and 3.2, ‘InIt.’ denotes the number of total iterations of finding suitable .

Example 5.1

The SEP with , , , , where , and are all generated uniformly randomly.

In the implementations, we take as the stopping criterion. Take the initial value , .

We make comparison of Algorithms 1.1, 1.2, 1.3, 3.1, 3.2 and FISTA with different J, N, M, and report the results in Tables 1, 2, 3 and Figure 1. We take the stepsize via a backtracking stepsize rule. For comparison, we tried to choose best parameters through numerical experiments. We take , , , and in Algorithms 1.1, 1.2, 3.1 and 3.2. And we take in Algorithm 1.3. So the stepsize is chosen in such a way that We take , and for FISTA with backtracking (see [21]). For comparison, the same random values are taken in each test for different algorithms.

Table 1

Computational results for Example with

J		50	100	150	200	250
Algorithm 1.1	Iter.	3263	90,378	297,135	65,795	31,655
	Inlt.	13,748	377,483	864,321	172,925	82,979
	Sec.	2.188	36.672	110.563	25.500	13.781
Algorithm 1.2	Iter.	8732	194,940	434,539	82,993	43,689
	Inlt.	46,182	1,069,758	1,327,537	225,913	141,813
	Sec.	6.797	106.234	189.078	38.094	23.422
Algorithm 1.3	Iter.	336	2012	4302	1327	676
	Sec.	0.063	0.406	1.125	0.406	0.219
FISTA	Iter.	1389	2580	3787	3260	2491
	Inlt.	1397	2589	3796	3270	2501
	Sec.	0.391	0.453	0.734	0.656	0.563
Algorithm 3.1(I)	Iter.	199	793	4342	482	718
	Inlt.	254	913	5206	1145	1148
	Sec.	0.094	0.156	0.953	0.156	0.188
Algorithm 3.1(II)	Iter.	152	764	2176	580	697
	Inlt.	184	860	2232	604	769
	Sec.	0.078	0.219	0.406	0.125	0.125
Algorithm 3.2(I)	Iter.	81	2208	5704	1022	362
	Inlt.	178	4680	7562	2261	846
	Sec.	0.063	0.469	1.891	0.250	0.094
Algorithm 3.2(II)	Iter.	65	952	2311	629	263
	Inlt.	80	976	2375	653	288
	Sec.	0.031	0.188	0.516	0.156	0.063

Table 2

Computational results for Example with

J		50	100	150	200	250
Algorithm 1.1	Iter.	5639	14,609	36,702	895,632	364,304
	Inlt.	29,524	43,313	107,748	2,740,752	1,179,620
	Sec.	2.813	6.297	19.125	551.250	259.656
Algorithm 1.2	Iter.	9175	33,306	109,165	2,481,066	566,203
	Inlt.	28,810	139,247	389,839	10,985,694	3,467,821
	Sec.	3.906	19.922	69.797	2079.344	649.234
Algorithm 1.3	Iter.	2559	9995	40,713	1,535,172	353,573
	Sec.	0.531	3.063	16.563	793.875	219.938
FISTA	Iter.	2158	3078	3092	17,010	9264
	Inlt.	2167	3088	3102	17,021	9275
	Sec.	0.375	0.656	0.797	4.922	2.984
Algorithm 3.1(I)	Iter.	123	186	1069	26,742	8307
	Inlt.	131	690	1779	32,790	15,519
	Sec.	0.031	0.063	0.359	7.750	5.141
Algorithm 3.1(II)	Iter.	136	171	726	17,575	3765
	Inlt.	160	187	807	18,007	3813
	Sec.	0.063	0.125	0.188	5.047	1.922
Algorithm 3.2(I)	Iter.	83	808	477	27,199	10,584
	Inlt.	182	1713	1140	31,301	12,084
	Sec.	0.063	0.125	0.147	13.078	5.563
Algorithm 3.2(II)	Iter.	43	235	297	15,515	7331
	Inlt.	66	251	322	15,839	7520
	Sec.	0.006	0.063	0.094	7.094	3.750

Table 3

Computational results for Example with

J		50	100	150	200	250
Algorithm 1.1	Iter.	3477	6828	15,749	48,555	255,440
	Inlt.	11,943	21,084	60,677	336,456	1,326,061
	Sec.	1.688	3.844	12.672	67.453	334.891
Algorithm 1.2	Iter.	10,464	21,319	32,055	122,743	483,468
	Inlt.	30,876	86,575	162,387	420,341	2,185,371
	Sec.	5.281	16.063	32.156	111.891	600.625
Algorithm 1.3	Iter.	648	6647	16,810	44,817	132,873
	Sec.	0.188	2.781	9.500	32.734	118.203
FISTA	Iter.	2346	2931	4040	3804	6977
	Inlt.	2355	2941	4051	3815	6989
	Sec.	0.500	0.750	1.250	1.344	3.141
Algorithm 3.1(I)	Iter.	109	236	343	814	2077
	Inlt.	151	278	415	975	2518
	Sec.	0.031	0.094	0.188	0.344	7.078
Algorithm 3.1(II)	Iter.	168	188	262	756	1106
	Inlt.	180	200	268	792	1142
	Sec.	0.063	0.109	0.125	0.218	0.438
Algorithm 3.2(I)	Iter.	117	168	818	1718	1725
	Inlt.	128	222	986	2000	2408
	Sec.	0.063	0.063	0.438	1.063	1.281
Algorithm 3.2(II)	Iter.	82	83	373	582	1240
	Inlt.	98	108	400	644	1285
	Sec.	0.031	0.078	0.199	0.406	0.563

Figure 1

Numbers of projections with .

The numerical results are listed in Tables 1, 2, 3 and Figures 1-6, from which we can get some conclusions:

Figure 6

Numbers of matrix-vector evaluations with .

Algorithm 1.2 behaves worst, and Algorithm 1.1 is superior to it, while inferior to Algorithms 1.3, 3.1 and 3.2.

The numbers of projections and matrix-vector evaluations that Algorithms 1.3, 3.1 and 3.2 need are close when are small. However, the numbers of projections and matrix-vector evaluations that Algorithms 3.1 and 3.2 need are less than those of Algorithm 1.3 as M, N become bigger.

In Figures 1, 3 and 5, the number of projections of Algorithm 3.1(I) and (II) (or Algorithm 3.2) is close although the iteration number of Algorithm 3.1(II) is less than that of Algorithm 3.1(I). The reason is that two projections are needed in Algorithm 3.1(II) while one projection is needed in Algorithm 3.1(I) per each iteration.

Figure 3

Numbers of projections with .

Figure 5

Numbers of projections with .

In Tables 1, 2 and 3, Algorithm 3.1(II) (or Algorithm 3.2(II)) has better performance than Algorithm 3.1(I) (or Algorithm 3.2(I)), maybe because the projections onto C and Q are very simple.

From Figures 1-5, it is observed that there exist peak values for Algorithms 1.3, 3.1 and 3.2, while FISTA has no peak values and is better than Algorithms 1.3, 3.1 and 3.2 near the peak values for some cases. However, for the other values of M, N, J, Algorithms 3.1 and 3.2 behave better than FISTA.

Conclusion

In this article we introduce two simultaneous projection algorithms and two semi-alternating projection algorithms to solve the SEP. We present larger stepsizes in (31) and (59) than those in Algorithms 2.1 and 2.2 in [15], which leads to a better contraction property and faster convergence speed of Algorithms 3.1 and 3.2. The weak convergence for the proposed methods is proved under standard conditions.

A numerical experiment is provided to illustrate that, except for FISTA, Algorithms 1.3, 3.1 and 3.2 have peak values. It is thus natural to combine our methods with inertial effects. This is one of our future research topics.