Agradecimientos Agradesco

Profesor Eladio Ocaña Y Profesor

Philippe Mahey

Keywords: block optimization problem, block optimization problems, 4, 6 Numerical Example, Uncertainty Environment

We focus in this thesis, on splitting methods which be applied to special optimization or inclusion problems considering its related inclusion problems with an appropriated Lagrangian map.

In a general setting for solving a monotone inclusion problem and obtain splitting algorithms, we develop a generalized proximal point and construct a related map with similar contraction properties as the resolvent map. Our general setting includes popular splitting algorithms.

Also, we show decomposition techniques in order to solve the multi-block version of our model problems, finding adequate formulations of the original problem and then apply a particular algorithm version of our general scheme.

Finally, we apply the splitting method to a large-scale energy production planning problem. 3.1 The separable model with coupling variable (SMCV) . . . . . . . . . 3.2 Separable model with coupling constraints (SMCC) . . . . . . . . . . 3.3 Proximal separation into two sub-blocks . . . . . . . . . . . . . . . . 3.3.1 The separable model with coupling variables . . . . . . . . . . 3.3.2 The separable model with coupling constraints . . . . . . . . . 3.4 Multi-

R ÉSUM É

Dans cette thèse, nous nous concentrons sur les méthodes d'éclatements d'opérateurs monotones (en anglais "splitting") qui sont appliquées à des problèmes d'optimisation convexe ou d'inclusions monotones présentant des structures décomposables. Une des contributions majeures de cette thèse est de considérer des reformulations de type point-selle (Lagrangiennes) de ces problèmes d'inclusion menant vers des algorithmes de type primal-dual. Dans un contexte général pour résoudre un problème d'inclusion monotone et obtenir des algorithmes de splitting, nous développons une méthode de point proximal généralisée et construisons une application associée avec des propriétés de contraction similaire à celle de l'opérateur résolvante classique. Notre configuration générale inclut des algorithmes de splitting connus dans l'abondante littérature qui les applique à des modèles d'apprentissage et de régression parcimonieuse très en vogue ces dix dernières années dans les domaines de l'inférence statistique et du traitement de signal et d'images.

Nous présentons également des techniques de décomposition afin de résoudre la version multi-blocs de certains problèmes structurés en exploitant des reformulations appropriées du problème original et puis en appliquant une version particulière de l'algorithme de notre schéma général.

Finalement, nous appliquons la méthode de splitting généralisée à un problème de planification de la production d'énergie à grande échelle avec une décomposition spatiale.

iii

Introduction

Decomposition techniques have been widely used in Mathematical Programming and Variational Analysis to cope mainly with the large-scale models issued from Decision Systems with many variables and constraints, but also to manage heterogeneous models including discrete choices, uncertainties or even mix of conflictual criteria. The recent explosion in size and complexity of datasets and the increased availability of computational resources has led to what is sometimes called the big data era. The large dimension of these data sets and the often parallel, distributed, or decentralized computational structures used for storing and handling the data, set new requirements on the optimization algorithms that solve these problems. Much effort has gone into developing algorithms that scale favorably with problem dimension and that can exploit structure in the problem as well as the computational environment.

Splitting methods for convex optimization or monotone variational analysis are commonly referred to address the construction of decomposition techniques based on regularization and duality. Indeed, many hard problems can be expressed under the form of a minimization of a sum of terms where each term is given by the composition of a convex function with a linear operator. The main advantage of splitting methods results thus from the fact that they can yield very efficient optimization schemes according to which a solution of the original problem is iteratively computed through solving a sequence of easier subproblems, each one involving only one of the terms constituting the objective function. These algorithms can also handle both differentiable and non smooth terms, the former by use of gradient operators (yielding explicit forward steps) and the latter by use of proximal operators (yielding implicit backward steps), thus giving rise to efficient first-order algorithms.

Since the pioneer works of Martinet [START_REF] Martinet | Régularistion d'inéquations variationnelles par approximations successives[END_REF], Glowinski-Marocco [START_REF] Glowinski | Sur l'approximation par éléments finis d'ordre 1 et la résolution par pénalisation-dualité d'une classe de problèmes de dirichlet[END_REF], Gabay [START_REF] Gabay | Applications of the method of multipliers to variational inequalities[END_REF] and Rockafellar [START_REF] Rockafellar | Monotone operators and the proximal point algorithm[END_REF], many algorithms have been studied for different models. A cornerstone was Lions and Mercier's paper in 1979 [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF] about the Douglas-Rachford's 1 family of splitting methods applied to the following inclusion : Find x ∈ X such that 0 ∈ S(x) + T (x) (V s )

where S and T are maximal monotone operators (typically subdifferential operators of convex functions) on a Hilbert space X. Back to the motone inclusion problem, i.e. to find x ∈ X such that 0 ∈ T (x) where T is maximal monotone, the proximal point method (PPM) constructs the mapping J T λ = (I + λT ) -1 for λ > 0, the resolvent of T with known contractive properties, and transforms the above inclusion into an equivalent fixed-point equation, i.e. x = J T λ (x). (PPM) is thus defined by the following fixed-point iteration :

x t+1 = J T λ (x t ) which corresponds, when T = ∂f , the subdifferential of a convex function f , to the following subproblem :

x t+1 = argmin x f (x) + 1 2λ x -x t 2
That so-called implicit backward step leads to the celebrated Augmented Lagrangian algorithm when f is the dual function associated to the Lagrangian Relaxation of a constrained concave maximization problem with many potential applications (see for instance [START_REF]Lagrangean relaxation[END_REF]).

In general the maximality assumption on T is restrictive, consider for instance the inclusion problem corresponding to the sum of two operators, i.e. inclusion (V). The sum S + T is not necessarily maximal monotone and also its resolvent doesn't necessarily maintain its separability structure. Fortunately Lions and Mercier [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF], solve that disadvantage considering an appropriate operator called after "Douglas-Rachford" operator (cf. [START_REF] Douglas | On the numerical solution of the heat conduction problem in 2 and 3 space variables[END_REF], defined by

G λ = I -J T λ + J S λ [2J T λ -I]
having splitting properties and whose fixed points are closely related with the solution points of problem (V ). Moreover, the maximality assumptions on T and S, ensure that G λ maintains the contractive properties of the single resolvent and the fullness of its domain.

When T = ∂f and S = ∂g, unlike the resolvent map, G λ is not in general the subdifferential of a function [START_REF] Eckstein | Splitting methods for monotone operators with applications to parallel optimization[END_REF], so we need to continue working with monotone operators even in optimization problems in order to obtain splitting methods for optimization problems.

The Douglas-Rachford method leads us to consider in Chapter 2 and 4 the following methodology: given an inclusion problem with separable structure, we construct an appropriated average map with full domain having splitting properties and then apply the fixed point method to this new map.

Gabay [START_REF] Gabay | Applications of the method of multipliers to variational inequalities[END_REF] noticed that the Douglas-Rachford operator is behind the popular Alternate Direction Method of Multipliers (ADMM), considering the dual variational problem associated with the composite model

Minimize f (x) + g(Ax)
where f : IR n → IR and g : IR m → IR are proper lower semi-continuous (lsc, for short) convex functions and A a given (m × n) matrix. The ADMM is an iterative method that consider two sub problems associated with f and g separately at each iteration.

Algorithm (ADMM)      x k+1 = argmin x f (x) + σ 2 Ax -z k + σ -1 y k 2 z k+1 = argmin z g(z) + σ 2 Ax k -z + σ -1 y k 2 y k+1 = y k + σ(Ax k+1 -z k+1 )
where σ is a positive real parameter.

Recently Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF] have presented a unified scheme algorithm for solving the last composite model based on the introduction of additional proximal terms like in Rockafellar's Proximal Method of Multipliers [START_REF] Rockafellar | Augmented lagrangians and applications of the proximal point algorithm in convex programming[END_REF], this algorithm includes a version of a Proximal ADMM and other known algorithms like Chambolle and Pock [START_REF] Chambolle | A first-order primal-dual algorithm for convex programs with applications to imaging[END_REF].

In Chapter 2, we consider an extended model problem coming from an energy production planning problem min (x,y)

[f (x) + g(y) : Ax + By = 0] (P )

where f and g are again proper lsc convex functions, and A and B are matrices of order m × n and m × p, respectively. Its saddle-point formulation in the variational setting is

Find (x, z, ȳ) ∈ IR n × IR p × IR m such that 0 ∈ L(x, z, ȳ) (V L )
where L is the maximal monotone map defined on IR n × IR p × IR m as L(x, z, y)

:=    ∂f (x) ∂g(z) 0    +    0 0 A t 0 0 B t -A -B 0       x z y    .
For solving (V L ), a generalized resolvent of L, defined as J L P = (L + P ) -1 P , is introduced, where P is a symmetric positive semidefinite matrix with special block structure in order to split J L P into a separable structure leaving f and g separated. Then we consider the relaxed fixed-point method applying to J L P , where after changing the variables we obtain the following generalized algorithm Generalized Splitting Scheme (GSS)

zk+1 ∈ argmin g(z) + 1 2 Bz + M -1 u k + Ax k 2 M + 1 2 z -z k 2 V 2 (1)
v k+ 1 2 = γAx k -(γ -1)Bz k + M -1 u k (2) xk+1 ∈ argmin f (x) + 1 2 Ax + v k+ 1 2 + 2γB zk+1 2 M + 1 2 x -x k 2 V 1 (3)
ũk+1 = u k + M (γAx k + (1 -γ)Ax k+1 + B zk+1 ) (4)
(x k+1 , z k+1 , u k+1 ) = ρ(x k+1 , zk+1 , ũk+1 ) + (1 -ρ)(x k , z k , u k ).

(

) 5
where the special cases γ = 0 and γ = 1, include the Shefi-Teboulle algorithm.

To guarantee the convergence of GSS and find an average map behind it, we developed a general setting in order to solve the inclusion problem: finding x ∈ IR n such that 0 ∈ T (x), studying the generalized proximal point associated to the generalized resolvent J T P = (T +P ) -1 P , where P is a symmetric positive semidefinite matrix, and then defining the map G T S as

G T S = S(T + S t S) -1 S t (6)
where S is a matrix satisfying S t S = P . By the monotonicity of T , G T S is 1 2 -average.

In the special case of T = L, we find conditions on the matrix S in order to G L S has full domain, which is used in the proof of the convergence of GSS. Also, considering S = S 3 (defined in Remark 2.4.3), we recover the Douglas-Rachford operator coinciding with λ

1 2 G L S 3 λ -1 2
, showing that its fundamental properties of splitability and 1 2 -average can be deduced from our generalized setting corresponding the Lagrangian map. Also, the general setting help us to find the ergodic and nonergodic rate of convergence of GSS.

Finally, in the last section of Chapter 2, for f i : IR n i → IR (i = 1, • • • , q) and g : IR m → IR proper lsc convex functions, A i and B matrices of order p × n i and p × m, respectively, the following S-Model problem is presented:

inf (x 1 ,••• , xq, z) q i=1 f i (x i ) + g(z) s.t q i=1 A i x i + Bz = 0.
This model is a multiblock generalization of problem (P ). In [START_REF] Dussault | Separable augmented lagrangian algorithm with multidimensional scaling for monotropic programming[END_REF] the authors propose the so called Separable Augmented Lagrangian Algorithm (SALA) in order to solve the S-Model problem by considering g = 0 and B = 0. At each iteration k, the algorithm considers a set of subproblems (one subproblem fo each i = 1, • • • , q) defined by

argmin x f i (x) + 1 2 A i x -z k i + M -1 y k 2 M ( 7
)
where M is a parameter matrix considered symmetric positive definite, solved in parallel processing.

As a generalization, they also consider different parameter matrices depending on each iteration.

Reformulating the S-Model problem in order to apply GSS, we get a new algorithm called "Proximal Multi-block Algorithm (PMA)" where like SALA it considers a set of problems at each iteration k:

argmin x f i (x) + 1 2 A i x -zk+1 i + M -1 i ỹk+1 c 2 M i + 1 2 x -x k i 2 Q i , i = 1, • • • , q,
where M i is symmetric positive definite and Q i symmetric positive semi-definite. Matrix Q i allow us to deal with the not injective case on A i and choosing it appropriately, each sub problem becomes on proximal step of f i .

In Chapter 3, we show some decomposition techniques, which consist in finding an adequate formulation of the original problem in order to apply GSS with particular parameter matrices V 1 , V 2 and M , and consequently find new splitting algorithms.

For every i ∈ {1, • • • , q}, let f i : IR n → IR be proper lsc convex. We consider the separable model with coupling variable (SMCV) defined as

min x q i=1 f i (x) (S)
For this problem we recover algorithm PDA and the algorithm given in [START_REF] Mahey | Proximal decomposition on the graph of a maximal monotone operator[END_REF] by considering a particular parameter. We also show the relationship between these two algorithms.

Adding linear constraints in the model, we obtain the separable model with coupling constraint (SMCC) where A i , i = 1, . . . , q are (p × n i ) matrices

min q i=1 f i (x i ) s.t q i=1 A i x i = 0
Notice that this problem can be seen as a S-Model problem and, since SMCC can be formulated as SMCV, we apply the results obtained for SMCV to SMCC, getting directly two algorithms which can also be recovered from P M A considering special formulation for SMCC.

The precedent algorithms found in this chapter, consider, for each iteration, the proximal step of all family {f i } i=1,••• ,q or separate the family into two sub-family, one consisting of {f i } i=1,••• ,q-1 and the other consisting of f q . Then, the proximal step of all {f i } i=1,••• ,q-1 are found in parallel processing and then, after linear combination of all these values, the proximal step of f q is found.

We show that after special reformulations, we get two splitting algorithms, one for SMCV and the other for SMCC. Each algorithm separates the problem into two sub-block problems, considering the proximal step to one sub-block and then (at a linear combination of the preceding values) the proximal step is found for the other sub-block, both in parallel processing.

In the last part of this chapter, we consider the following multi-block optimization problem

min x=(x 1 ,••• , xq) q i=1 f i (x i ) + g( q i=1 A i x i ) + s(x) (P sc )
where for i ∈ {1, ..., q}, f i : IR n i → IR , g : IR p → IR and s : IR n → IR (n = q i=1 n i ) are proper lsc convex functions, and A i are matrices of order p × n i . This special structure is deduced from the formulation over the Euclidean space of the following stochastic optimization model problem with finite scenarios Ξ min

X∈L E ξ T t=1 g t (X t (ξ), ξ) : s.t X ∈ N and T t=1 B ξ t X t (ξ) = 0, ∀ξ ∈ Ξ (SP )
where L is the linear space of all mapping X from Ξ to IR n := IR n 1 × • • • × IR n T , and N the nonanticipativity subspace of L, and for ξ ∈ Ξ and t = 1, • • • , T , g t (•, ξ) is a proper lsc convex function and B ξ t an m ξ × n t matrix.

As a consequence of obtaining a splitting algorithm for (P sc ), we get a splitting algorithm for problem (SP ).

Coming back to problem (V s ), if one these maps (say T ) is co-coercive, then we can apply the Forward-Backward method or Backward-Forward method, which combine the Backward step of S with the Forward step of T . Notice that the Forward step only needs the value of single value map T unlike the Backward step needs the value of the resolvent map of S which in general is not easy. D. Davis and W. Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] generalized problem (V s ) considering the sum of three maps:

Find x such that 0 ∈ S(x) + T (x) + C(x)

where S and T are again maximal monotone (with T not necessarily single value) and C a co-coercive operator with full domain. They combine separately the Backward steps on S and T , with the Forward step on C, defining the following map

G := I -J T λ + J S λ (2J T λ -I -λC(J T λ ))
which clearly extends the Douglas-Rachford operator G λ and the operators corresponding to the Forward-Backward and Backward-Forward methods.

In Chapter 4, we consider the more general composite monotone inclusion:

Find x such that 0 ∈ S(x) + A t T (Ax) + C(x) (V ar)

where S : IR n -→ -→ IR n and T : IR m -→ -→ IR m are maximal monotone maps, C : IR n → IR n a β-co-coercive with full domain and A an m × n matrix.

For the particular case C = 0, the map G coincides with the Douglas-Rachford operator which in turn is equal to λ 1 2 G L S 3 λ -1 2 as mentioned just after [START_REF] Bot | Admm for monotone operators : convergence analysis and rates[END_REF]. For the general cases (C not necessarily equal to zero) and assuming A injective, we get from [START_REF] Bot | Admm for monotone operators : convergence analysis and rates[END_REF] considering T = L an alternative Lagrangian map and S = D a special matrix, a map G L D which extends the Davis-Yin operator G (coinciding it with λ 1 2 G L D λ -1 2 when A = I) and maintaining similar properties as G, for instance under mild assumptions, G L D is an average map with full domain. Then using the generalized resolvent J L D , we get a new splitting algorithm which converges to a saddle-point of Lagrangian map associated with primal problem (V ar) because, as mentioned previously, G L D is an average map with full domain.

For the general case, where A is not necessarily injective, problem (V ar) is reformulated as

0 ∈ S(x) + A t M 1 2

is injective if and only if

A t M A + V is invertible. So, assuming that condition and applying the former splitting algorithm (for the injective case on A) to (V ar 1 ), we get a splitting algorithm for problem (V ar) in the general setting, which is termed "Generalized splitting algorithm for three operators (GSA3O)":

(GSA3O) zk+1 = (T + M ) -1 (y k + M Ax k ) (8) ỹk+1 = y k + M Ax k -M zk+1 (9)
r k+1 = C((V + A t M A) -1 (V x k + A t M zk+1 )) (10) xk+1
= (S + V + A t M A) -1 V x k + A t M zk+1 -A t ỹk+1 -r k+1 (11) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ). ( 12
)
The structure of problem (V ar) is related to the variational formulation of the minimization of separable convex functions:

Minimize f (x) + g(Ax) + h(x) (13)
where f : IR n → Ī R and g : IR m → Ī R are proper lower semi-continuous convex functions, h : IR n → IR is convex and ( 1β )-Lipschitz-differentiable, and A an m×n matrix.

Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] presents two types of algorithms for solving (13) that we call CA1 and CA2, for simplicity, we have considered with fixed relaxation parameter ρ > 0 and without error term.

Algorithm (CA1)

       xk+1 = (τ ∂f + I n×n ) -1 (x k -τ ∇h(x k ) -τ A t y k ) ỹk+1 = (σ∂g * + I m×m ) -1 (y k + σA(2x k+1 -x k )) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k )
and the other one switch the roles of primal and dual variables:

Algorithm (CA2)        ỹk+1 = (σ∂g * + I m×m ) -1 (y k + σAx k ) xk+1 = (τ ∂f + I n×n ) -1 (x k -τ ∇h(x k ) -τ A t (2ỹ k+1 -y k )) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k )
The main difference between these algorithms is the order of action of the proximal steps. In the same manner, we present two versions of each algorithm proposed in this chapter.

From our general setting, choosing special parameter matrices and different Lagrangian maps, we get different variants of algorithms CA1 and CA2. From the variant of CA1 we recover in particular the recently algorithm proposed by Yang.

Using the same techniques given in Chapter 2, we show the ergodic and nonergodic rates of convergence of the algorithms found in this chapter.

The last model considered in this chapter is concerned with the more general S-Model problem defined as:

inf x=(x 1 ,••• , xq), z q i=1 f i (x i ) + h(x) + g(z) s.t q i=1 A i x i -Bz = 0
where for i = 1, • • • , q, f i : IR n i → IR and g : IR m → IR are proper lsc convex functions, h : IR n → IR is convex and ( 1β )-Lipschitz-differentiable (n = q i=1 n i ), and A i and B are matrices of order p × n i and p × m, respectively.

Rewriting this last problem as a variational inclusion problem having similar structure as problem (V ar), we apply GSA3O getting a new algorithm called "Separable Primal -Dual Variant (SPDV)" which is a generalization of PMA appearing in Chapter 2.

Finally, in Chapter 5, we apply the splitting algorithm developed in the previous chapters in order to solve a model of long-term energy pricing problem. Rewriting this model as problem (P ), we apply algorithm GSS, getting three types of algorithms. Finally, we give some ideas on how to deal with the considered model in the stochastic case where many 'scenarios' with given probabilities need to be included, thus minimizing the total expected cost on the horizon.

9

Chapter 1

Notations, preliminaries and basic results on convex optimization problems 1.1 Notations and basic definitions

Throughout this thesis, we will use the following notations on convex optimization and variational inequality problems . Most of the theoretical material can be found in [START_REF] Rockafellar | Convex Analyis[END_REF][START_REF] Bauschke | Convex Analysis and Monotone Operator Theory in Hilbert Spaces[END_REF].

We will denote by IR the set of extended real numbers IR∪{±∞} = [-∞, +∞]. For a given subset C ⊂ IR n we will denote by cl (C), int (C) and ri (C), the closure, the interior and the relative interior of C, respectively.

For a given set C ⊂ IR n , the orthogonal subspace to C, denoted by C ⊥ , is the linear subspace

C ⊥ = {x * ∈ IR n : x * , x = 0 for all x ∈ C}. For a closed convex set C ⊂ IR n , we denote by Proj C (x) the projection of x ∈ IR n onto C which consists of all ȳ ∈ C satisfying x -ȳ ≤ x -y for all y ∈ C
where • denotes a norm of IR n . Of course, if C = ∅, then ȳ satisfying this inequality is unique if the considered norm is Euclidean (unless otherwise stated, we will use in all the Thesis this type of norm). For F ⊂ IR n , the set Proj C (F ) denotes the collection of all Proj C (x) for x ∈ F .

For a closed convex set C ⊂ IR n , the normal cone of C at a given point x ∈ C is the set denoted by N C (x) and defined as

N C (x) = {x * ∈ IR n : x * , y -x ≤ 0 for all y ∈ C}, assuming N C (x) = ∅ if x / ∈ C. A function f : IR n → IR is said to be convex, if its epigraph epi (f ) := {(x, λ) ∈ IR n × IR : f (x) ≤ λ}
is convex; and concave, if -f is convex. The function f is said to be lower semicontinuous (lsc, for short) at a given point x if for every λ ∈ IR verifying λ < f (x) there exists an open set V containing x such that λ < f (x) for all x ∈ V . This function is lsc if it is lsc at every point of

IR n . Of course, f is lsc if epi (f ) is closed in IR n × IR.
Assuming f convex, it is said to be proper if f (x) > -∞ for all x ∈ IR n and its domain defined as

dom (f ) := {x ∈ IR n : f (x) < +∞} = Proj IR n (epi (f )) is nonempty. Of course, dom (f ) is convex if f is convex.
A function f : IR n → IR is said to be strongly convex (with modulus α > 0) or α-strongly convex if for all x, y ∈ IR n and all t ∈ [0, 1], one has

f ((1 -t)x + ty) ≤ (1 -t)f (x) + tf (y) - α 2 t(1 -t) x -y 2 .
The function f is said to be β-Lipschitz differentiable function (with β > 0) if it is differentiable whose gradient is β-Lipschitz (ie, Lipschitz continuous with constant β). In that case, f is of real value on the whole IR n .

Associated to a function f : IR n → IR, its Fenchel-conjugate is the function

f * defined on IR n as f * (x * ) = sup x∈IR n [ x * , x -f (x)]
and its biconjugate is the function f * * which is the conjugate of f * , ie

f * * (x) = sup x * ∈IR n [ x, x * -f * (x * )].
It is clear that f * and f * * are convex and lsc, and f * * is the greatest lsc convex function bounded from above by f . Moreover, f * * and f coincide if f is proper lsc convex.

The subdifferential of f at a point x ∈ IR n is the set

∂f (x) := {x * ∈ IR n : f (x) + x * , y -x ≤ f (y) for all y ∈ IR n } or equivalently ∂f (x) = {x * ∈ IR n : f (x) + f * (x * ) ≤ x * , x }.
Clearly ∂f (x) = ∅ if x / ∈ dom (f ) or if f is not lsc at x. In general, ∂f (x) is convex and closed, maybe empty. It is nonempty and bounded on int (dom (f )).

Another very important property of subdifferential is its monotonicity. For all (x, x * ), (y, y * ) in the graph of ∂f one has

x * -y * , x -y ≥ 0.
In general, a multivalued map (or simply, map) Γ : IR n -→ -→ IR n is said to be monotone if for all (x, x * ), (y, y * ) in the graph of Γ [is the set consisting of all pair (z, z * ) ∈ IR n × IR n such that z * ∈ Γ(z)] one has

x * -y * , x -y ≥ 0.
It is clear that if Γ is monotone, then its inverse map Γ -1 defined by Γ -1 (x * ) = {x : x * ∈ Γ(x)}, is monotone. So, the monotonicity property can be seen as a property on the graph instead on the map itself.

The map Γ is said to be maximal monotone if for any monotone map Σ : IR n -→ -→ IR n satisfying Γ(x) ⊂ Σ(x) for all x ∈ IR n , one has Γ = Σ. It also follows that Γ is maximal monotone if and only if Γ -1 is maximal monotone.

A very important characterization of the maximality in the monotone sense is given by Minty's theorem [START_REF] Minty | Monotone operators in hilbert spaces[END_REF]. It say that a monotone map Γ : IR n -→ -→ IR n is maximal monotone if and only if the inverse map (I + Γ) -1 , which is single-valued and with full domain. Here I denotes the identity map from IR n into itself.

Analogous to the strongly convexity, a map Γ : IR n -→ -→ IR n is said to be strongly monotone (with modulus α > 0) or α-strongly monotone if Γ -ρI is monotone, i.e. for all (x, x * ), (y, y * ) in the graph of Γ, it holds

x * -y * , x -y ≥ α x -y 2 .
One deduces that f : IR n → IR is strongly convex if and only if its subdifferential ∂f is strongly monotone.

The inverse of a strongly monotone map (with modulus α) is clearly single value and α -1 -Lipschitz.

A map Γ is said to be co-coercive with constant β ( or shortly β-co-coercive) if its inverse Γ -1 is β-strongly monotone. That is, for all (x, x * ), (y, y * ) in the graph of Γ, it holds

x * -y * , x -y ≥ β x * -y * 2 .
One deduce that if Γ is co-coercive with constant β then Γ is β -1 -Lipschitz. When β ≥ 1, the map Γ is nonexpansive. In general, a map Γ is said to be nonexpansive if it is Lipschitz with constant ≤ 1, i.e. if there exists γ ≤ 1 such that for all (x, x * ), (y, y * ) in the graph of Γ, it holds that

x * -y * ≤ γ x -y .

Another important property used in many parts of the thesis is the α-average of a map. A map Γ is said to be α-average if Γ = (1 -α)I + αR where R is a nonexpansive map. A 2 -1 -average map is also called firmly nonexpansive. For example, the resolvent of a maximal monotone map Γ, J Γ := (I +Γ) -1 is firmly nonexpansive (and defined on the whole space).

We finish this section by introducing the following notations that we will use for instance in Chapters 2, 3 and 4. For arbitrary maps T 1 and T 2 and vectors x and y of appropriated dimensions, we denote

T 1 T 2 x y = T 1 (x) T 2 (y) .
Analogously, for two given functions g 1 and g 2 , we denote

(g 1 , g 2 )(z 1 , z 2 ) = g 1 (z 1 ) + g 2 (z 2 )
for all z 1 and z 2 of appropriated dimensions.

The duality scheme

An optimization problem in the mathematical context can be set as

α := inf [f (x) : x ∈ IR n ] (P )
where f : IR n → IR is a given function. Problem (P ) is commonly called primal problem.

In order to develop the duality scheme following Rockafellar's scheme [START_REF] Rockafellar | Convex Analyis[END_REF], we introduce a duality space IR p and a perturbation function ϕ : IR n ×IR p → IR satisfying ϕ(x, 0) = f (x) for all x ∈ IR n . Then the associated perturbed primal problems is defined as

h(u) := inf [ϕ(x, u) : x ∈ IR n ]. (P u ) If ϕ is convex on IR n × IR p then h is convex on IR p ; but if ϕ is lsc it does not neces- sarily imply that h is lsc.
It is clear that h(0) = α, then the duality arises from the idea to find in other way h(0), for this we use the Fenchel-conjugate function h * of h,

h * (u * ) := sup u [ u * , u -h(u)] = ϕ * (0, u * ). The biconjugate of h is h * * (u) = sup u * [ u, u * -ϕ * (0, u * )].
The biconjugate, under some conditions (see [START_REF] Rockafellar | Convex Analyis[END_REF]), allows us to recover the initial function. In general,

β := h * * (0) ≤ h(0) = α.
Then the dual problem is defined as

β = h * * (0) = sup [-ϕ * (0, u * ) : u * ∈ IR p ]. (D)
The primal and dual problem are also related through the Lagrangian function l : IR n × IR p → IR defined as:

l(x, u * ) := inf [ϕ(x, u) -u * , u : u ∈ IR p ].
So, if ϕ is convex proper lsc, the primal and dual problem are respectively:

inf x sup u * l(x, u * ) and sup u * inf x l(x, u * ).
In order to obtain optimal solution of primal and dual problems without duality gap (α = β), the Saddle Point problem arises which consist in finding (x, ȳ) ∈ IR n × IR p such that inf

x l(x, ȳ) = l(x, ȳ) = sup y l(x, y).
Under some regularity condition on ϕ, the primal and dual problems can be respectively formulated as inclusion problems called optimality condition

Find x ∈ IR n such that 0 ∈ π t 1 ∂ϕ(π 1 x) (P opc ) and Find y ∈ IR p such that 0 ∈ π t 2 ∂ϕ * (π 2 y) (D opc )
where

π 1 = I n×n 0 p×p and π 2 = 0 n×n I p×p .
The Saddle Point problem can also be formulated as

Find (x, y) ∈ IR n × IR p such that (0, 0) ∈ (∂ x l) × (∂ y [-l])(x, y) (L opc )
which, in terms of ∂ϕ and assuming ϕ proper lsc convex, it holds that

(z, y) ∈ ∂ϕ(x, u) if only if z ∈ ∂ x l(x, y), u ∈ ∂ y [-l](x, y).
Since ∂ϕ * = (∂ϕ) -1 , then problems (P V ), (D V ) and (L V ) are equivalent to each other in the sense that the mapping intervening in each inclusion problem is the composite of ∂ϕ or ∂ϕ -1 with π 1 or π 2 (and its respective transpose matrix).

Following this construction, [START_REF] Pennanen | A splitting method for composite mappings[END_REF] (see also [START_REF] Ocaña | Un schéma de dualité pour les problèmes d'inéquations variationnelles[END_REF]) has extended this duality scheme to general variational or inclusion problems which can be set as Find x ∈ IR n such that 0 ∈ S(x).

(V )

A perturbation map associated to problem (V ) is a map F satisfying π t 1 F π 1 = S. Then the corresponding dual and lagrangian problems are respectively formulated as Find y ∈ IR p such that 0

∈ π t 2 F -1 (π 2 y) (D V ) and Find (x, y) ∈ IR n × IR p such that (0, 0) ∈ L(x, y) (L V )
where (z, y) ∈ F (x, u) if and only if (z, u) ∈ L(x, y).

Coming back to the duality scheme for an optimization problem, we reformulate the duality and its respective Lagrangian problem through their optimality conditions for some particular classes of optimization problems.

The composite model

A composite model is an optimization problem that can be set as Minimize f (x) + g(Ax)

(P c )
where f : IR n → IR and g : IR p → IR are proper lsc convex functions, and A a p × n matrix. We consider the following perturbation function ϕ :

IR n × IR p → IR defined by ϕ(z, u) = f (x) + g(Ax + u)
and its corresponding dual problem

Minimize f * (-A t y) + g * (y) (D c )
and its Lagrangian function l : IR n × IR p → IR defined by l(x, y) := f (x) -g * (y) + y, Ax .

Under some regularity conditions, the optimal conditions of (P c ) and (D c ) are respectively Find x ∈ IR n such that 0 ∈ ∂f (x) + A t ∂g(Ax)

(P cv ) and Find x ∈ IR n such that 0 ∈ -A∂f * (-A t y) + ∂g * (y). (D cv )
The corresponding Saddle point problem can also be formulated as

Find (x, y) ∈ IR n × IR p such that 0 ∈ ∂f (x) ∂g * (y) + 0 A t -A 0 x y . (L cv ) 1.2.

The separable case

We now consider a more general optimization problem regarding to the previous one min

(x,z) [f (x) + g(z) : Ax + Bz = 0] (P sc )
where f : IR n → IR and g : IR p → IR are again proper lsc convex functions, and A and B are matrices of order m × n and m × p, respectively. It is clear that problem (P sc ) includes the composite model (P c ) by considering B = -I p×p . Conversely, problem (P sc ) can be written as the following composite model: min

(x,z) f (x) + g(z) + δ {0} A B x z .
So, the dual problem is min

y f * (-A t y) + g * (-B t y) (D sc )
and the lagrangian function l sc :

IR n × IR p × IR m → IR is l sc (x, z, y) = f (x) + g(z) + y, Ax + Bz .
The optimality conditions of (P sc ) and (D sc ) and the corresponding saddle point problem are respectively:

Find (x, z) ∈ IR n × IR p st 0 0 ∈ ∂f (x) ∂g(z) + A t B t N {0} (Ax + Bz) (P scv ) and Find x ∈ IR n such that 0 ∈ -A∂f * (-A t y) -B∂g * (-B t y) (D scv ) and Find (x, z, ȳ) ∈ IR n × IR p × IR m such that 0 ∈ L(x, z, ȳ) (L scv )
where L is the maximal monotone map defined on

IR n × IR p × IR m as L(x, z, y) = (∂ x,z l sc )×(∂ y [-l sc ]) =    ∂f (x) ∂g(z) 0    +    0 0 A t 0 0 B t -A -B 0       x z y    . (1.1)

The optimization problem with linear subspace constraints

We note that problem (P sc ) (in particular problem (P c )) can be formulated as

min z∈V ψ(z) (P sbp )
where ψ : IR r → IR is a proper lsc convex function and V a linear subspace of IR r . We present now two equivalent formulations of (P sbp ) like the primal and dual problems by expressing V as Image space and as well as kernel space of suitable matrices.

It is clear that (P sbp ) can be formulated as the composite model:

min ψ(z) + δ V (z) ( P sbp )
whose associated dual problem is

min ψ * (z * ) + δ V ⊥ (z * ) ( D sbp )
which can be formulated as the following optimization problem with linear subspace constraints min

z * ∈V ⊥ ψ * (z * ). (D sbp )
Let K and W be matrices of order p × r and r × n, respectively, such that V = ker(K) = range (W ). Then (P sbp ) and (D sbp ) can also be formulated respectively as the following optimization problems termed primal projection problem

min x∈IR n ψ • W (x) (P V sbp )
and dual projection problem

min u * ∈IR p ψ * • K t (u * ). (D V sbp )
We observe that problems (P sbp ) and (P V sbp ) are defined on different linear spaces and, if W is injective, then problem (P V sbp ) is defined on a linear space of dimension dim(V ) which is less or equal than r, that is, the dimension of the linear space where problem (P sbp ) is defined. Similarly, problems (D sbp ) and (D V sbp ) are also defined on different linear spaces.

As a special case, problems (P c ) and (D c ) can also be formulated respectively as the primal projection and dual projection problems regarding the following optimization problem min

(x,y)∈V ψ(x, y) (Op)
where ψ(x, y) := f (x) + g(y) and

V = range I n×n A = ker -A I p×p .
In the separable case, problem (D sc ) can also be formulated as the dual projection problem related to (Op) by considering V = ker -A -B . Then, by setting R and D matrices of order n × q and p × q, respectively, such that

V = range R D
we obtain two related problems, the primal projection problem

min z∈IR q f (Rz) + g(Dz)
and the optimization problem with linear subspace constrains min

(x * ,y * ) [f * (x * ) + g * (z * ) : R t x * + D t z * = 0].
Remark 1.2.1 Considering IR r = IR n × IR p , the linear subspace V coincides with IR n × {0 p } and the matrices W and K are exactly the projection matrices π 1 and π t 2 , respectively, defined in problems (P opc ) and (D opc ). So,

(ψ • W )(x) = ψ(x, 0) and (ψ * • K t )(u * ) = ψ * (0, u *
) for all x ∈ IR n and u * ∈ IR p and hence problems (P V sbp ) and (P V sbp ) are respectively the primal and dual optimization problems associated to perturbation function ϕ = ψ described in the duality scheme.

The gradient and proximal point methods

The gradient and proximal point methods are apparently the most popular and basic methods to solve an optimization problem or (more generally) an inclusion problem.

Let h : IR r → IR be an ( 1 β )-Lipschitz-differentiable convex function. In order to find a minimizer of h, the gradient method generates, from a given initial point x 0 ∈ IR r , the iterative points defined by:

x k+1 = x k -α∇h(x k ).
The corresponding sequence converges to a minimizer of h if α ∈ ]0, 2 β [ . Regarding the proximal point algorithm in order to find a minimizer of a proper lsc convex function f : IR r → IR, this algorithm constructs the so called Moreau envelope function, an alternative Lipschitz-differentiable convex function having the same minimizers as function f . Then the proximal point algorithm can be recovered from the gradient algorithm applied to the Moreau envelope function.

Given λ > 0, the Moreau envelope function of f is the function f λ defined as

f λ (z) := min x f (x) + 1 2λ x -z 2 . ( * )
One deduce immediately that

• inf f = inf f λ and argmin f = argmin f λ ; and

• f λ is differentiable on IR r and its gradient is

∇f λ (z) = 1 λ (I -(λ∂f + I) -1 )(z)
which is ( 1 λ )-Lipschitz. The set value (λ∂f + I) -1 (z) is singleton whose unique element is the minimizer of problem ( * ).

So, given an initial point z 0 ∈ IR r , the proximal point method generates a sequence defined by

z k+1 = z k -ρλ∇f λ (z k ) = (1 -ρ)z k + ρ(λ∂f + I) -1 (z k ).
This sequence converges to a minimizer of f if ρ ∈ ]0, 2[ .

Application to the dual projection problem

The dual projection problem (D V sbp ), min

u * ∈IR p ψ * • K t (u * )
represents, as we saw, the general formulation of all dual problems described in Subsection 1.2.3. In order to apply the gradient method to this problem, we check under what condition the dual objective function ψ * • K t is Lipschitz-differentiable on the whole IR p . This is the object of the next proposition.

Proposition 1.3.1 With the same notations as before, suppose that ψ is proper lsc α-strongly convex and K a nonzero matrix, then ψ * • K t is differentiable with gradient K∇ψ * K t which is K 2 α -Lipschitz with full domain.

Proof. From the assumptions, ∂ψ * is univalued with full domain. Then the subdiferential of ψ * •K t is K∇φ * K t having full domain. On other hand, for arbitrary points x, y ∈ IR p , we have that

K∇ψ * K t x -K∇ψ * K t y, x -y = ∇ψ * K t x -∇ψ * K t y, K t x -K t y
and since ∂ψ = (∇ψ * ) -1 is α-strongly monotone, we get

K∇ψ * K t x -K∇ψ * K t y, x -y ≥ α ∇ψ * K t x -∇ψ * K t y 2 which implies, if K is nonzero, that K∇ψ * K t x -K∇ψ * K t y, x -y ≥ α K 2 K∇ψ * K t x -K∇ψ * K t y 2 .
The Lipschitz constant of K∇ψ * K t is deduced applying the Cauchy-Schwarz's inequality.

Regarding the proximal point algorithm applied to before dual problem, for a given r × r symmetric positive definite matrix Q, we consider a little more general Moreau envelope function of f denoted by f Q , defined as

f Q (z) := min x f (x) + 1 2 x -z 2 Q . (env)
Similarly to the classical Moreau envelope function f λ , it holds that

• inf f = inf f Q and argmin f = argmin f Q ; • f Q is differentiable on IR r and its gradient is ∇f Q (z) = Q[I -(∂f + Q) -1 Q](z).
The set [(∂f + Q) -1 Q](z) is singleton whose element is the optimal solution of the minimization problem (env);

• Q -1/2 ∇f Q Q -1/2 is 1-Lipschitz on IR r .
The next proposition shows another way to express the Moreau envelope function for the objective dual function.

Proposition 1.3.2 Let ψ : IR r → IR be a proper lsc convex function, K a p × r matrix satisfying Im (K t ) ∩ ri dom (ψ * ) = ∅. For a p × p positive definite matrix M one has

(ψ * • K t ) M -1 (u * ) = -inf x ψ(x) + 1 2 Kx 2 M -u, Kx for all u ∈ IR p . ( * * ) Furthermore, denoting z u := (∂(ψ * • K t ) + M -1 ) -1 M -1 u the minimizer of problem (env) with f = (ψ * •K t ) and Q = M -1 , and x u = (∂ψ+K t M K) -1 K t u a minimizer of problem ( * * ), then z u = u -M Kx u .
Proof. From assumptions, ψ * • K t is proper lsc convex and then

(ψ * • K t ) M -1 = [(ψ * • K t ) M -1 ] * * = [(ψ * • K t ) * + 1 2 • 2 M ] * . Also, (ψ * • K t ) * (v) = inf x [ψ(x) : Kx = v] and hence (ψ * • K t ) M -1 (u) = -inf v inf x [ψ(x) + 1 2 v 2 M -v, u : Kx = v] which implies that (ψ * • K t ) M -1 (u) = -inf x ψ(x) + 1 2 Kx 2 M -u, Kx .
The relationship between the optimal solutions follows from this expression.

On the other hand, since both optimization problems min x f (x) and min

w f Q (w)
have same optimal values and same minimizers, then under the regularization condition given in Proposition 1.3.2, both optimization problems min y∈IR p ψ * • K t (y) and min

v∈IR p -min x ψ(x) + 1 2 Kx 2 M -v, Kx
have also same optimal values and same minimizers. The problem on the right is termed Augmented Dual Problem [START_REF] Rockafellar | Variational Analyis[END_REF] and its objective function is nothing else than ψ + 1 2 K(•) 2

M * • K t and hence the augmented dual problem is the dual projection problem corresponding to problem (P sbp ) with objective function ψ + 1 2 K(•) 2 M . In particular, for a given perturbation function ϕ of f in the duality scheme, the objective function of the dual problem corresponding to ϕ(x, u) + 1 2 u 2 (another perturbation function of f ), results to be Lipschitz-differentiable if ϕ and π 2 satisfy the conditions given in Proposition 1.3.2 considering ψ = ϕ and K = π 2 .

The resolvent map corresponding to the Saddle Point Problem

Rockafellar [START_REF] Rockafellar | Augmented lagrangians and applications of the proximal point algorithm in convex programming[END_REF], considered the proximal point method to a saddle point problem corresponding to a convex optimization problem with inequality constraints, getting the so called "proximal multiplier algorithm". Corresponding to the Saddle Point Problem (L scv ) defined in Section 1.2.2 and associated to the diagonal block symmetric positive definite matrix

P :=    W 1 0 0 0 W 2 0 0 0 M   
where W 1 , W 2 and M are symmetric positive definite matrices of order n × n, p × p and m × m, respectively, we introduce the following resolvent map defined by (L + P ) -1 P.

Since L is maximal monotone, Minty's theorem guarantees that (L + P ) -1 P has full domain and its value at each point (x, z, y) ∈ IR n × IR p × IR m is a singleton whose value is

(L + P ) -1 P (x, z, y) = (η, ν, M -1 (M y + Aη + Bν) : 0 ∈ H 1 (x,z) (η, ν) 0 ∈ H 2 (x,z) (η, ν)
where

H 1 (x,z) (η, ν) := ∂f (η) + W 1 (η -x) + A t y + A t M -1 (Aη + Bν) and H 2 (x,z) (η, ν) := ∂g(ν) + W 2 (ν -z) + B t y + B t M -1 (Aη + Bν).
It is noteworthy that the involved subproblems cannot directly be splitted because of the coupling on their variables η and ν is present and hence the solvability of such subproblems becomes very difficult in practice. In Chapter 2 we present another matrix P avoiding the aforementioned coupling.

The α-average maps

The gradient and proximal point method have common structure in the sense that both methods can be formulated as a relaxed fixed point method [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] for a suitable mapping having the following property defined now.

Definition 1 An operator T is α-average if α ∈ ]0, 1[ and there exists a nonexpansive map N such that T = (1 -α)I + αN .

For instance, the α-average maps involved in the gradient and proximal point methods are respectively β∇h and I -(λ∂f +I) -1 . Baillon-Haddad's theorem shows that if h is convex and ( 1β )-Lipschitz-differentiable, then β∇h is 1 2 -average having full domain. On the other hand, Rockafellar [START_REF] Rockafellar | Convex Analyis[END_REF] (Proposition 12.11) shows that if f is proper lsc convex, then (λ∂f + I) -1 and hence I -(λ∂f + I) -1 are 1 2 -average having full domain.

Another important example of α-average map is given through a maximal monotone map T : IR n -→ -→ IR n and a n × n positive definite matrix M . It is not difficult to show that the resolvent map

M 1 2 (T + M ) -1 M 1 2
which is related with the multidimensional scaling proximal point method [START_REF] Lenoir | A survey of monotone operator splitting methods for the decomposition of convex programs[END_REF] is 1 2 -average having (due Minty) full domain.

In general, an α-average map T satisfies the following important inequality

T x -T y 2 ≤ x -y 2 - 1 -α α (I -T )x -(I -T )y 2
which immediately implies the convergence of the relaxed fixed point algorithm as mentioned in next proposition.

Proposition 1.4.1 Let T : IR n → IR n be an α-average map with full domain and ρ ∈ ]0, α -1 [ . Assume that the fixed point set of T is nonempty. For a given initial point x 0 ∈ IR n consider the following iteration points

x n+1 = ρT (x n ) + (1 -ρ)x n .
Then the corresponding sequence {x n } converges to a fixed point of T .

This important convergence result makes it possible to deal with general monotone inclusion problems by transforming them into fixed point equations corresponding to α-average maps. Moreover, for practical treatments is also important that the corresponding α-average map possesses splitting property.

We give two examples of α-average maps corresponding to inclusion problems for the sum of two and three monotones maps. The first one is due the Douglas & Rachford scheme discussed in [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF] and the second one is due to Davis & Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF].

Douglas-Rachford map

Consider the following inclusion problem for the sum of two maps 0 ∈ S(x) + T (x) where S and T are two maximal monotone maps from IR n into itself. It is well known that S + T is not necessarily maximal monotone which is a condition to apply the proximal point method. Another disadvantage of this method is the absence of splitting structure of its corresponding resolvent map.

Lions and Mercier [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF], reformulate the above inclusion into a fixed-point equation with respect to an appropriated operator called after "Douglas-Rachford" operator, defined by

G λ = I -J T λ + J S λ [2J T λ -I], (1.2)
where J T λ = (λT + I) -1 and J S λ = (λS + I) -1 are the resolvent maps of T and S, respectively. This map is 1 2 -average and, unlike the resolvent map of S + T , it has splitting property and having full domain. This map is behind the popular ADMM algorithm as noticed Eckstein [START_REF] Eckstein | On the douglas-rachford splitting method and the proximal point algorithm for maximal monotone operators[END_REF]. On the other hand, D. O'Connor and L. Vandenberghe [START_REF] O'connor | On the equivalence of the primal-dual hybrid gradient method and douglas-rachford splitting[END_REF] have recently shown that Chambolle-Pock algorithm [START_REF] Chambolle | A first-order primal-dual algorithm for convex programs with applications to imaging[END_REF] is also constructed using this map by considering S and T with special structures in the inclusion problem.

Davis-Yin map

We consider the sum of three maps

0 ∈ S(x) + T (x) + C(x) (V 2 )
where S and T are two maximal monotone maps and C a β-co-coercive function with full domain, all from IR n into itself. It is possible to apply the Douglas-Rachford method considering the sum of S (or T ) and C as a unique map, but in general this procedure does not split S and C (or T and C).

Davis and Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] have recently considered the following map

G := I -J T λ + J S λ [2J T λ -I -λC(J T λ )] (1.3)
having splitting property and defined everywhere of

IR n . It is also α-average for α = 2β 4β-λ , if λ ∈ ]0, 2β[ .

Convergence Study

We recall that x n = O(y n ) means that there exists a positive C such that for all n sufficiently large

x n ≤ C y n .

And x n = o(y n ) means that xn yn converges to 0. Also we say that x n converge linearly to x * if there exists a positive C < 1 such that for all n sufficiently large

x n+1 -x * ≤ C x n -x * .
For example, H. Brezis et P.L. Lions [?] showed that given a monotone map T with at least one zero then the fixed point residual (FPR)

J T λ (x n ) -(x n ) is O( 1 √ k+1
). D. Davis and W. Yin [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] improve this result for any average map with has at least a fixed point getting that its FPR is o(

1 √ k+1 ) and the ergodic FPR is O( 1 k ).
Rockafellar [START_REF] Rockafellar | Monotone operators and the proximal point algorithm[END_REF] showed that if we consider T a strongly monotone map (or more generally T -1 Lipschitz continuous at 0) and T has at least one zero then the proximal point applied to T generates a linearly convergent sequence.

Lions and Mercier [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF] showed that if T in 4.2 is strongly monotone and Lipschitz then the Douglas-Rachford method generates a linearly convergent sequence. Recently Giselsson [START_REF] Giselsson | Tight global linear convergence rate bounds for douglas rachford splitting[END_REF] gave a best upper bound rate as Lions-Mercier and showed linear convergence under other regularity conditions over S or T , proving that in these cases the map G λ is contractive. In the convex case D. Davis and W. Yin [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] showed the ergodic and nonergodic convergence rate of the feasibility and objective function error related to the relaxed Douglas-Rachford method. D. Davis and W. Yin in [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] showed that under regularity assumptions the map G is a contractive map from which the linear convergence is deduced.

In this thesis, we will not focus on the special cases when linear convergence is attained, rather keeping the analysis on global or point wise convergence in the ergodic or non ergodic sense.

Chapter 2

A unified splitting algorithm for composite monotone inclusions Operator splitting methods have been recently concerned with inclusions problems based on composite operators made of the sum of two monotone operators, one of them associated with a linear transformation. We analyze here a general and new splitting method which indeed splits both operator proximal steps, and avoiding costly numerical algebra on the linear operator. The family of algorithms induced by our generalized setting includes known methods like Chambolle-Pock primal-dual algorithm and Shefi-Teboulle Proximal Alternate Direction method of multipliers. The study of the ergodic and non ergodic convergence rates show similar rates with the classical Douglas-Rachford splitting scheme. We end with an application to a multi-block convex optimization model which leads to a generalized Separable Augmented Lagrangian algorithm1 .

Introduction

Composite models involving sums and compositions of linear and monotone operators are very common and still challenging problems like in constrained separable convex optimization or composite variational inequalities. We will consider here composite monotone inclusions of the form (X and Y are Hilbert spaces) :

0 ∈ S(x) + A * T (Ax) (2.1)
where S : X → X and T : Y → Y are maximal monotone operators and A : X → Y is a linear transformation (associated with its adjoint operator A * , which will be denoted by A t when dealing with finite-dimensional spaces).

Most existing monotone operator splitting methods can deal with composite models, for example the Douglas-Rachford family (see [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF]) and its special decomposition versions, the Alternate Direction Method of Multipliers (ADMM) (see [START_REF] Gabay | Applications of the method of multipliers to variational inequalities[END_REF][START_REF] Eckstein | On the douglas-rachford splitting method and the proximal point algorithm for maximal monotone operators[END_REF]) and the Partial Inverse or Proximal Decomposition Algorithm (PDA) (see [START_REF] Spingarn | Applications of the method of partial inverses to convex programming:decomposition[END_REF][START_REF] Mahey | Proximal decomposition on the graph of a maximal monotone operator[END_REF][START_REF] Pennanen | A splitting method for composite mappings[END_REF]).

Lions and Mercier [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF] analyzed the Douglas-Rachford's method (including the limiting case of Peaceman-Rachford splitting, PRS) for the case of the sum of two maximal monotone operators (S + T ), alternating between proximal steps applied to each operator separately. Gabay [START_REF] Gabay | Applications of the method of multipliers to variational inequalities[END_REF] analyzed the case S + A * T A where A is an injective linear transformation (and A * its adjoint), yielding the celebrated Alternative Direction Method of Multipliers (ADMM). Spingarn [START_REF] Spingarn | Partial inverse of a monotone operator[END_REF] studied the case when the operator is the sum of the normal cone of a closed subspace M and a maximal monotone operator T . Later, Pennanen [START_REF] Pennanen | A splitting method for composite mappings[END_REF] showed how to reformulate that model as a monotone inclusion

The first study which explicitly considered an algorithm to solve the composite inclusion which avoids the use of projection (or proximal) steps on the range of A was proposed in [START_REF] Briceno-Arias | A monotone+skew splitting model for composite monotone inclusions in duality[END_REF] (an extension of Spingarn's Partial Inverse to composite models was proposed too in [START_REF] Alghamdi | A primaldual method of partial inverses for composite inclusions[END_REF]). The corresponding algorithms solve the dual problem at the same time, which is defined by :

0 ∈ -AS -1 (-A * y) + T -1 (y)
Many applications surge in the minimization of separable convex functions like :

Minimize f (x) + g(Ax) (2.2)
where f : IR n → IR and g : IR m → IR are proper lower semi-continuous convex functions and A is a given (m × n) matrix. The Dual problem in the sense of Rockafellar-Fenchel theory is :

Minimize f * (-A t y) + g * (y)
where f * is the Fenchel-conjugate of f . Recently Chambolle-Pock [START_REF] Chambolle | A first-order primal-dual algorithm for convex programs with applications to imaging[END_REF] studied model (2.2) and introduced new splitting schemes applied to a Lagrangian formulation of the primal minimization problem. They applied a primal-dual version of (ADMM) to the following saddle-point formulation :

min x max y f (x) -g * (y) + Ax, y
Observe that we could as well define a Lagrangian operator associated with the composite inclusion (2.1) :

L(x, y) = [S(x) + A t y] × [T -1 (y) -Ax] (2.3)
Chambolle and Pock's algorithm relies on two Proximal steps on f and g with an additional extrapolation step (in a similar fashion of Varga's iterative principle [START_REF] Varga | Matrix Iterative Analysis[END_REF]) as summarized below :

     x k+1 = (I + τ ∂f ) -1 (x k -τ A t ȳk ) y k+1 = (I + σ∂g * ) -1 (y k + σAx k+1 ) ȳk+1 = y k+1 + θ(y k+1 -y k )
where (I + τ ∂f ) -1 is the resolvent operator of the subdifferential operator S = ∂f which is known to be defined on the whole space and supposed to be easily computable in a so-called 'backward' proximal step as detailed below.

The difference and presumed advantage of that formulation is the symmetry (considering that x and y can be updated in reverse order) and a potentially decomposable algorithm which depends on three parameters. Their convergence result states that we should choose their values such that στ A 2 < 1.

Observe now that (CPA) can be rewritten using Augmented Lagrangian-like functions by using the Moreau identity (see [START_REF] Moreau | Proximité et dualité dans un espace hilbertien[END_REF]) :

(I + σ∂g * ) -1 (y) + σ(I + σ -1 ∂g) -1 (σ -1 y) = y
Resuming the transformed steps into the following iteration:

Algorithm (CPA)          x k+1 = argmin x f (x) + 1 2τ x -x k + τ A t ȳk 2 z k+1 = argmin z g(z) + σ 2 z -Ax k+1 -σ -1 y k 2 y k+1 = y k + σ(Ax k+1 -z k+1 ) ȳk+1 = y k+1 + θ(y k+1 -y k )
Chambolle and Pock confirmed the expected rate of convergence in O(1/k) and even obtain the accelerated rate of O(1/k 2 ) following the FISTA scheme of Beck and Teboulle [START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF] (thus reaching Nesterov's optimal rates in convex programming [START_REF] Nesterov | Smooth minimization of nonsmooth functions[END_REF]).

In a recent survey, Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF] have presented a unified scheme algorithm for solving model (2.2) based on the introduction of additional proximal terms like in Rockafellar's Proximal Method of Multipliers [START_REF] Rockafellar | Augmented lagrangians and applications of the proximal point algorithm in convex programming[END_REF]. The resulting schemes include a version of a Proximal (ADMM) and other known algorithms like Chambolle-Pock's method (CPA). Indeed, a generic sequential algorithm proposed by Shefi and Teboulle is the following three steps scheme :

Algorithm (STA)      x k+1 = argmin x f (x) + σ 2 Ax -z k + σ -1 y k 2 + 1 2 x -x k 2 M 1 z k+1 = argmin z g(z) + σ 2 Ax k+1 -z + σ -1 y k 2 + 1 2 z -z k 2 M 2 y k+1 = y k + σ(Ax k+1 -z k+1 )
where . M is the norm induced by a symmetric positive definite matrix M , i.e.

x 2 M = x t M x. Algorithm (STA) makes use of alternate minimization steps on the Augmented Lagrangian function associated with the coupling subspace Ax -z = 0. It is noted in [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF] that (CPA) with the choice θ = 1 corresponds exactly to (STA) with M 1 = τ -1 I -σA t A and M 2 = 0 (which implies again that στ A 2 < 1).

Later, Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] extended the model (2.2) and algorithm (CPA) to the case f = F + h where F : IR n → IR is convex and smooth. He relaxed the restriction on the parameters allowing στ A 2 = 1 and also includes the Douglas-Rachford family in the case of A = I (therefore we can say that Chambolle-Pock's method generalized Douglas-Rachford's splitting scheme). Condat showed too that Chambolle-Pock's method is the proximal point method applied to the Lagrangian operator associated with the primal and dual pair of inclusions.

In this chapter we will further extend the algorithms surveyed by Shefi and Teboulle, in order to solve the following convex optimization problem min

(x,z) [f (x) + g(z) : Ax + Bz = 0]. (P )
where f and g are again convex lsc functions and, A and B are two matrices of order m × n and m × p, respectively. It is clear that this problem includes problem (2.2) by considering B = -I p×p .

The primal variational formulation of (P ) is the following

Find (x, z) ∈ IR n × IR p such that 0 0 ∈ ∂f (x) ∂g(z) + A t B t N {0} (Ax + Bz) (2.4) where N C (a) is the normal cone to set C at point a.
The dual variational formulation of (P ) is the following

Find y ∈ IR m such that 0 ∈ -A(∂f ) -1 (-A t y) -B(∂g) -1 (-B t y) (2.5)
In Section 2.2, we propose a generalized proximal point method (GPPM) which was developed implicitly by Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF], where we consider specific assumptions to relax the condition of symmetric positive definiteness of the matrix associated with the resolvent, to authorize matrices which are only symmetric positive semidefinite, maintaining the properties of convergence of the proximal method.

In Section 2.3, we apply GPPM in order to find a zero of the Lagrangian map associated with problem (P ), selecting an appropriate symmetric positive semi definite matrix in order to obtain a Generalized Splitting Scheme (GSS), which includes various known algorithms, for instance both types of algorithms studied by Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF] correspond indeed to particular choices of the parameters in GSS.

In Section 2.4, we define a 1-co-coercive operator G T P related to GPPM, which set of fixed points is related to the zeroes of T . When T is the Lagrangian operator and the matrix P has a special structure as considered in Section 2.3, we show examples where we can get that operator explicitly, in particular we can recover the Douglas-Rachford operator.

In Section 2.5, we investigate the rate of converge of the GSS scheme, in the ergodic and non ergodic sense, analyzing the convergence of the sequences of the optimal values and the constraints violations associated with problem (P ).

Finally, section 2.6 applies the GSS scheme to some general multi-block convex optimization problem with a composite structure. We show the relationship with a separable Augmented Lagrangian algorithm (SALA) introduced in [START_REF] Mahey | Adaptive scaling and convergence rates of a separable augmented lagrangian algorithm[END_REF].

A generalized proximal point method

The classical Proximal Point method is used to solve a monotone inclusion

Find x ∈ IR r such that 0 ∈ T (x) (V )
where T : IR r -→ -→ IR r is a maximal monotone operator. We denote by sol (V ) the solution set of problem (V ). It is closed, convex and may be empty. The iteration exploits the contractive properties of the resolvent operator J T τ = (I + τ T ) -1 to define a sequence given by x k+1 = J T τ (x k ) which converges weakly to a solution of (V ) if it is nonempty.

Following former ideas developed by Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] in the proof of the convergence of a specialized splitting method closely related (CPA), we define the generalized Proximal Point iteration by substituting the classical resolvent by

J T P := (T + P ) -1 P (2.6)
where P is an r × r symmetric positive semidefinite matrix. Since T is monotone, then for any (x, x * ), (x, x * ) ∈ graph (J T P ), one has

x * -x * , P x -P x ≥ P x * -P x * , x * -x * ≥ 0. (2.7)
We deduce immediately the following properties:

• T + P and thereby its inverse (T + P ) -1 are monotone.

• R := P + I r×r -Q is a symmetric positive definite matrix, whenever Q is the orthogonal projection onto the image of P , which implies in particular that Q satisfies QP = P Q = P and Q 2 = Q.

• J T P = J T P Q, where Q is as above. As R is symmetric positive definite, it induces an inner product on IR r , u, v R := Ru, v for all u, v ∈ IR r with its corresponding norm u R := u, u R for all u ∈ IR r .

Hence, from (2.7), for all x, x ∈ dom (QJ T P ) = dom (J T P ),

QJ T P (x) -QJ T P (x), x -x R ≥ QJ T P (x) -QJ T P (x) 2 R ,
which implies that QJ T P is 1-co-coercive wrt R on domain of J T P .

We deduce immediately the following relationship between the solution set of problem (V ) and the fixed points of J T P and QJ T P . Proposition 2.2.1 With the same notations as before, we have

• x ∈ sol (V ) if and only if x is a fixed point of J T P . • v is a fixed point of QJ T P if and only if v = Qx for some x ∈ sol (V ) ∩ J T P (v). Proof.
The first property is directly by definition. The second one follows from the fact that v ∈ QJ T P v if and only if there exists x such that x ∈ J T P (v) satisfying v = Qx. It follows that x ∈ J T P (v) = J T P (Qx) = J T P (x). Using the first equivalence we deduce that x belongs to sol (V ).

Concerning the regularity of J T P , we have • If P is positive definite, then Q = I r×r and R = P . We deduce that J T P = QJ T P and then J T P is 1-co-coercive wrt P on the whole of its domain.

• If P is not positive definite, then J T P may not be single valued. But if it is single valued, then it is continuous on the whole of its domain.

We consider now a relaxed version of the generalized proximal iteration. In connection with the resolvent operator J T P and a real positive parameter ρ, we consider for an arbitrary point x 0 ∈ dom J T P , the sequence {x k } defined by

x k+1 ∈ ρJ T P (x k ) + (1 -ρ)x k . (2.8)
Notice that this sequence is well defined whenever range (ρJ

T P + (1 -ρ)I) ⊆ dom (J T P ).
Concerning the convergence of {x k }, we distinguish the following situations:

• If P is positive definite, then J T P is 1-co-coercive wrt P (hence single valued) with full domain which implies that {x k } converges, for ρ ∈ (0, 2), assuming sol (V ) nonempty. In fact, given x * ∈ sol (V ), the convergence follows from the inequality

x k -x * 2 P ≥ 2 -ρ ρ x k+1 -x k 2 P + x k+1 -x * 2 P .
• In general, since QJ T P is 1-co-coercive wrt R, then for ρ ∈ (0, 2) and assuming that QJ T P has closed domain and nonempty fixed point set (which is equivalently to sol (V ) being nonempty), the sequence {Qx k } is convergent. The convergence of {x k } needs additional assumptions as we show in the following proposition.

Proposition 2.2.2 Let T : IR r -→ -→ IR r be maximal monotone and P be an r × r positive semidefinite matrix. Assuming J T P single valued (which implies that it is continuous) with closed domain and sol (V ) not empty. Then, for ρ ∈ (0, 2), the sequence {x k } converges to some point belonging to sol (V ).

Proof. Since QJ T P is 1-co-coercive wrt R, it is single valued on its domain; and since J T P = J T P Q, then from (2.8) we obtain that

Qx k+1 = ρQJ T P (Qx k ) + (1 -ρ)Qx k . (2.9)
Using again the fact that QJ T P is 1-co-coercive wrt R and, by assumptions with closed domain, ρ ∈ (0, 2) and sol (V ) nonempty, then {Qx k } converges to some point a, which is a fixed point of QJ T P . From Proposition 2.2.1 and the single valuedness assumption, J T P (a) ∈ sol (V ). On the other hand, using the triangular inequality in (2.8) we have

x k+1 -J T P (a) ≤ ρ J T P (Qx k ) -J T P (a) + |1 -ρ| x k -J T P (a) .
Since J T P is continuous, the sequence J T P (Qx k ) -J T P (a) converges to 0. We deduce that {x k } converges to J T P (a). Some examples of specially tailored co-coercive operators will be discussed in Section 2.4.

Generalized splitting algorithms

With the convex minimization problem (P ) defined in Section 2.1, we associate its Lagrangian function defined as

l(x, z, y) = f (x) + g(z) + y, Ax + Bz (2.10)
and then its saddle-point problem in the variational setting

Find (x, z, ȳ) ∈ IR n × IR p × IR m such that 0 ∈ L(x, z, ȳ) (V L )
where L is the maximal monotone map defined on

IR n × IR p × IR m as L(x, z, y) := (∂ x,z l) × (∂ y [-l]) =    ∂f (x) ∂g(z) 0    +    0 0 A t 0 0 B t -A -B 0       x z y    . (2.11)
The map L, as the sum of maximal monotone operators and a skew-symmetric linear operator, satisfies similar inequalities as the subdifferential of a convex-concave bifunction. These inequalities will be used in order to obtain the rate of convergence studied in Section 2.5.

Proposition 2.3.1 For any (d, d * ), ( d, d * ) ∈ graph (L), considering d = (x, z, y) and d = (x, z, ȳ), it holds d -d, d * ≥ l(x, z, ȳ) -l(x, z, y) ≥ d -d, d * . These inequalities are still verified if we consider (d, d * ) ∈ graph (L) and d ∈ dom (f ) × dom (g) × IR m , for the first inequality; and ( d, d * ) ∈ graph (L) and d ∈ dom (f ) × dom (g) × IR m , for the second inequality.
It is well known that, under some regularity conditions, problem (V L ) admits a saddle-point if and only if problem (P ) admits an optimal solution. One instance of such regularity condition is :

There exist x ∈ ri (dom f ) and z ∈ ri (dom g) such that Ax + Bz = 0.

(H)

We now apply to problem (V L ) the relaxed proximal method described in the previous section for a specially tailored matrix P in order to provide a separable structure to the algorithm.

The separable structure on the main step

In this part we describe the main iteration step of the relaxed proximal method given in (2.8) providing a decomposable structure.

We will choose an appropriate symmetric matrix P in order to split (L + P ) -1 or equivalently J L P = (L + P ) -1 P , into a separable structure leaving f and g separated. To that end, given (x, z, ỹ) ∈ IR n × IR p × IR m , we analyze the solution of the following inclusion system: Find (x, z, y) such that

   ∂f (x) ∂g(z) 0    +    0 0 A t 0 0 B t -A -B 0       x z y    +   
      x z y       x z ỹ    .
We introduce now two parameters α, β ∈ IR, and a positive definite matrix M to simplify the third row-block of

P into P 3 = [(1 + α)A (1 + β)B M -1 ]
. So, the last inclusion can be expressed as

y = M ỹ -αM Ax -βM Bz (2.12)
and hence, replacing it in the second block-system, this results in

∂g(z) + (2 + β)B t (M ỹ -αM Ax -βM Bz) + P 21 x + P 22 z z.
So, in order to express this last system eliminating primal variable x, we need to consider

P 21 = α(2 + β)B t M A, obtaining z ∈ (∂g + P 22 -β(2 + β)B t M B) -1 (z -(2 + β)B t M ỹ). (2.13)
Using again (2.12), now in the first block system, we get

∂f (x) + (2 + α)A t (M ỹ -αM Ax -βM Bz) + P 11 x + α(2 + β)B t M Az x which is equivalent to x ∈ (∂f + P 11 -α(2 + α)A t M A) -1 (x -(2 + α)A t M ỹ -2(α -β)A t M Bz). (2.14)
Summarizing the previous sequence in order to get a separable structure, we must first solve system (2.13), then system (2.14) and finally system (2.12). The corresponding matrix P , of order (r × r) with r = n + p + m, is then of the form

P :=    C 1 α(2 + β)A t M B (1 + α)A t α(2 + β)B t M A C 2 (1 + β)B t (1 + α)A (1 + β)B M -1    (2.15)
where C 1 (n × n), C 2 (p × p) are arbitrary symmetric matrices, From the maximality of ∂f and ∂g, the inclusions in (2.13) and (2.14) are indeed equalities if the matrices defined as

W 1 := C 1 -α(2 + α)A t M A and W 2 := C 2 -β(2 + β)B t M B,
are positive definite. In that case (L + P ) -1 is single-valued with full domain and therefore J L P is continuous with full domain. It is clear that P is symmetric. It is positive semidefinite (resp. positive definite) if and only if the matrix

U := C 1 -(1 + α) 2 A t M A (α -β -1)A t M B (α -β -1)B t M A C 2 -(1 + β) 2 B t M B (2.16)
is positive semidefinite (resp. positive definite).

We now list some conditions in order to get a positive semidefinite matrix U :

A1 If C 1 -[(1 + α) 2 + (α -β -1) 2 ]A t M A and C 2 -[(1 + β) 2 + 1]B t M B are positive semidefinite then U is positive semidefinite. A2 If C 1 -[(1 + α) 2 + 1]A t M A and C 2 -[(1 + β) 2 + (α -β -1) 2 ]B t M B are positive semidefinite then U is positive semidefinite. A3 If β ≤ α -1, and C 1 -[(1 + α) 2 + (α -β -1)]A t M A and C 2 -[(1 + β) 2 + (α - β -1)]B t M B are positive semidefinite then U is positive semidefinite. A4 If β = α -1. Then C 1 -(1 + α) 2 A t M A and C 2 -α 2 B t M B are positive semidefinite if only if U is positive semidefinite.
In order to calculate the sequence in (2.8), we first calculate (x k+1 , zk+1 , ỹk+1 ) =

J L P (x k , z k , y k ), which is equal to (x k+1 , zk+1 , ỹk+1 ) = (L + P ) -1    C 1 x k + α(2 + β)A t M Bz k + (1 + α)A t y k α(2 + β)B t M Ax k + C 2 z k + (1 + β)B t y k (1 + α)Ax k + (1 + β)Bz k + M -1 y k   
Then from (2.13), we have that

zk+1 = Jg W 2 z -(β + 2)B t M Ax k (2.17) where z = C 2 z k -(2 + β)(1 + β)B t M Bz k -B t y k and Jg W 2 = (∂g + W 2 ) -1
is the generalized resolvent operator associated with the convex function g.

From (2.14), we have that

xk+1 = Jf W 1 x -2(α -β)A t M B zk+1 (2.18) where x = C 1 x k -(2 + α)(1 + α)A t M Ax k + (α -2β -2)A t M Bz k -A t y k and Jf W 1 = (∂f + W 1 ) -1
is the generalized resolvent operator associated with the convex function f ; and from (2.12), we have that

ỹk+1 = y k + (1 + α)M Ax k + (1 + β)M Bz k -αM Ax k+1 -βM B zk+1 . (2.19)
The sequence in (2.8) is completed with an extrapolation step for a given ρ ∈ (0, 2):

(x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ).
(2.20)

We obtain the following proposition directly from Proposition 2.2.2.

Proposition 2.3.2 Let ρ ∈ (0, 2). Assume that C 1 ∈ IR n×n , C 2 ∈ IR p×p and M ∈ IR m×m are symmetric, with M positive definite; and α, β ∈ IR, such that W 1 and W 2 are positive definite and satisfying one of conditions (A1)-(A4). If sol (V L ) is nonempty, then for an arbitrary (x 0 , z 0 , y 0 ) ∈ IR n+p+m , the sequence (x k , z k , y k ) defined by the sequential update formulas (2.17 → 2.18 → 2.19 → 2.20) converges to some element of sol (V L ).

We will now further reformulate the iteration to show the alternating steps on separable Augmented Lagrangian functions. We introduce the parameter γ = α -β and the matrices defined as

V 1 := W 1 -A t M A and V 2 := W 2 -B t M B.
(2.21)

The conditions (A1) -(A4) become:

A1' If V 1 -(γ -1) 2 A t M A and V 2 -B t M B are positive semidefinite then U is positive semidefinite. A2' If V 1 -A t M A and V 2 -(γ -1) 2 B t M B are positive semidefinite then U is positive semidefinite. A3' If γ ≥ 1. Then V 1 -(γ-1)A t M A and V 2 -(γ-1)B t M B are positive semidefinite then U is positive semidefinite. A4' If γ = 1. Then V 1 and V 2 are positive semidefinite if only if U is positive semidefinite.
We introduce a new primal-dual auxiliary variable

u k := y k + (α -γ + 1)M Ax k + (1 + β)M Bz k ,
to obtain the following updates :

z k+ 1 2 = V 2 z k -B t u k (2.22) zk+1 = J g W 2 [z k+ 1 2 -B t M Ax k ] (2.23)
x k+ 1 2 = V 1 x k -γA t M Ax k + (γ -1)A t M Bz k -A t u k (2.24) xk+1 = J f W 1 [x k+ 1 2 -2γA t M B zk+1 ] (2.25) ũk+1 = u k + γM Ax k + (1 -γ)M Ax k+1 + M B zk+1 (2.26) (x k+1 , z k+1 , u k+1 ) = ρ(x k+1 , zk+1 , ũk+1 ) + (1 -ρ)(x k , z k , u k ) (2.27)
which is equivalent to the following sequential minimization subproblems :

Generalized Splitting Scheme (GSS) zk+1 ∈ argmin g(z) + 1 2 Bz + M -1 u k + Ax k 2 M + 1 2 z -z k 2 V 2
(2.28)

v k+ 1 2 = γAx k -(γ -1)Bz k + M -1 u k (2.29) xk+1 ∈ argmin f (x) + 1 2 Ax + v k+ 1 2 + 2γB zk+1 2 M + 1 2 x -x k 2 V 1 (2.30) ũk+1 = u k + M (γAx k + (1 -γ)Ax k+1 + B zk+1 ) (2.31) (x k+1 , z k+1 , u k+1 ) = ρ(x k+1 , zk+1 , ũk+1 ) + (1 -ρ)(x k , z k , u k ). (2.32)
From Proposition 2.3.2, we obtain the proposition of convergence of (GSS)

Proposition 2.3.3 Let ρ ∈ (0, 2). Assume that V 1 ∈ IR n×n , V 2 ∈ IR p×p and M ∈ IR m×m are symmetric, with M positive definite such that V 1 +A t M A and V 2 +B t M B are positive definite. Let γ ∈ IR such that one of conditions (A1 )-(A4 ) is satisfied. If sol (V L
) is nonempty, then for an arbitrary (x 0 , z 0 , u 0 ) ∈ IR n+p+m , the sequence

(x k , z k , u k ) in (2.28)-(2.32) converges to some element of sol (V L ).
We analyze now the special cases when γ = 0 and γ = 1, which correspond to the two types of algorithms proposed by Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF].

Case γ = 0 From (A1 ), if both matrices V 1 -A t M A and V 2 -B t M B are positive semi-definite then P is a positive semi-definite matrix.
Switching the order (2.28) for (2.30), we get the following algorithm where the primal updates are performed in parallel:

xk+1 ∈ argmin f (x) + 1 2 Ax + Bz k + M -1 u k 2 M + 1 2 x -x k 2 V 1 (2.33) zk+1 ∈ argmin g(z) + 1 2 Ax k + Bz + M -1 u k 2 M + 1 2 z -z k 2 V 2 (2.34) ũk+1 = u k + M (Ax k+1 + B zk+1 ) (2.35) (x k+1 , z k+1 , u k+1 ) = ρ(x k+1 , zk+1 , ũk+1 ) + (1 -ρ)(x k , z k , u k ) (2.

36)

If B = -I p×p , M = cI p×p and ρ = 1 , we obtain the algorithm STA type I proposed by Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF].

Summarizing, from Proposition 2.3.3, we obtain the following proposition of convergence of the sequence defined by (2.33)-(2.36).

Proposition 2.3.4 Let ρ ∈ (0, 2). Assume that V 1 ∈ IR n×n , V 2 ∈ IR p×p and M ∈ IR m×m are symmetric, with M positive definite, such that V 1 +A t M A and V 2 +B t M B are positive definite and V 1 -A t M A and V 2 -B t M B are positive semi-definite. If sol (V L ) is nonempty, then for an arbitrary (x 0 , z 0 , u 0 ) ∈ IR n+p+m , the sequence (x k , z k , u k ) in (2.33)-(2.36) converges to some element of sol (V L ). Case γ = 1
From (A4 ), it holds that V 1 and V 2 are positive semi-definite if only if P is a positive semi-definite matrix. In this case GSS becomes :

zk+1 ∈ argmin g(z) + 1 2 Ax k + Bz + M -1 u k 2 M + 1 2 z -z k 2 V 2 (2.37) ũk+1 = u k + M (Ax k + B zk+1 ) (2.38) xk+1 ∈ argmin f (x) + 1 2 Ax + B zk+1 + M -1 ũk+1 2 M + 1 2 x -x k 2 V 1 (2.39) (x k+1 , z k+1 , u k+1 ) = ρ(x k+1 , zk+1 , ũk+1 ) + (1 -ρ)(x k , z k , u k ) (2.40) If B = -I p×p , M = τ I p×p , V 2 = 0 and V 1 = σ -1 I n×n -τ A t T A such that 1 ≥ στ A 2 ,
then we obtain the over relaxed algorithm proposed by Chambolle-Pock [START_REF] Chambolle | A first-order primal-dual algorithm for convex programs with applications to imaging[END_REF].

Considering ρ = 1 and defining, xk := x k , zk := z k+1 and ūk := u k+1 , then substituting in (2.37)-(2.39) and switching the order, we get the following algorithm

xk+1 ∈ argmin f (x) + 1 2 Ax + B zk + M -1 ūk 2 M + 1 2 x -xk 2 V 1 (2.41) zk+1 ∈ argmin g(z) + 1 2 Ax k+1 + Bz + M -1 ūk 2 M + 1 2 z -zk 2 V 2 (2.42) ūk+1 = ūk + M (Ax k+1 + B zk+1 ) (2.43)
If B = -I p×p and M = cI p×p , we obtain the algorithm STA type II proposed by Shefi and Teboulle [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF], which is called the Proximal Alternating Direction Method (PADM).

Further transformations applied to (2.37)-(2.40) lead us to consider two interesting algorithms. The first of them is obtained by considering V 2 = 0, and considering the auxiliary variables x k+1 , z k , u k , s k to update the relaxed sequences

x k+1 := 1 ρ x k+1 + (1 -1 ρ )x k = xk+1 , z k := 1 ρ z k+1 + (1 -1 ρ )z k = zk+1 , u k := ũk+1 and s k := x k , getting x k+1 ∈ argmin f (x) + 1 2 Ax + B z k + M -1 u k 2 M + 1 2 x -s k 2 V 1
(2.44)

z k+1 ∈ argmin g(z) + 1 2 ρA x k+1 + Bz + M -1 u k + (ρ -1)B z k 2 M
(2.45)

u k+1 = u k + ρM A x k+1 + (ρ -1)M B z k + M B z k+1
(2.46)

s k+1 = ρ x k+1 + (1 -ρ) s k (2.47)
The second interesting algorithm is obtained by considering the auxiliary variables xk , žk , ǔk , šk to update the relaxed sequences xk :

= 1 ρ x k+1 + (1 -1 ρ )x k = xk+1 , žk := 1 ρ z k+1 + (1 -1 ρ )z k = zk+1 , ǔk := ũk+1 and šk := x k , getting žk+1 ∈ argmin g(z) + 1 2 ρAx k + Bz + M -1 ǔk + (ρ -1)B žk 2 M (2.48) ǔk+1 = ǔk + ρM Ax k + (ρ -1)M B žk + M B žk+1 (2.49) šk+1 = ρx k + (1 -ρ)š k (2.50) xk+1 ∈ argmin f (x) + 1 2 Ax + B žk+1 + M -1 ǔk+1 2 M + 1 2 x -šk+1 2 V 1 (2.51)
So, by considering in these two last algorithms B = -I p×p , M = cI p×p and V 1 = 0, the sequences s k and šk becomes unnecessary. Moreover, (2.44)-(2.47) become the generalized ADMM proposed by Eckstein [START_REF] Eckstein | On the douglas-rachford splitting method and the proximal point algorithm for maximal monotone operators[END_REF], and (2.48)-(2.51) become the algorithm 2 considered in [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF].

From Proposition 2.3.3, we obtain the convergence of the sequence (2.37)-(2.40)

Proposition 2.3.5 Let ρ ∈ (0, 2). Assume that V 1 ∈ IR n×n , V 2 ∈ IR p×p and M ∈ IR m×m are symmetric, with V 1 and V 2 positive semi-definite and M positive definite such that V 1 + A t M A and V 2 + B t M B are positive definite. If sol (V L ) is nonempty, then for an arbitrary (x 0 , z 0 , u 0 ) ∈ IR n+p+m , the sequence (x k , z k , u k ) defined in (2.37)-(2.40) converges to some element of sol (V L ).

The co-coercive map associated with GPPM

Lions and Mercier [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF] have transformed an inclusion problem of the sum of two maximal monotone operators (S + T ) into a fixed-point equation with respect to an appropriated operator, the Douglas-Rachford operator, which is 1-co-coercive map and, in order to compute its value at each point of its domain, only local calculations of proximal terms of S and T separately are needed. Eckstein [START_REF] Eckstein | On the douglas-rachford splitting method and the proximal point algorithm for maximal monotone operators[END_REF] later showed the relationship between the splitting algorithm (ADMM) and the fixed-point method applied to a Douglas-Rachford operator, after a suitable linear transformation.

In our general setting, we show in this section that the sequence generated by the generalized proximal point method (GPPM) corresponding to map J T P for arbitrary maximal monote operator T and arbitrary symmetric positive semidefinite matrix P is nothing else but the sequence generated by the fixed point method corresponding to map G T S defined in (2.53), after a linear transformation S (satisfying P = S t S). It leads thus in some sense to a generalization of a Douglas-Rachford operator, keeping the property of 1-co-coercivity.

As pointed out in Subsection 2.3.1, the sequence generated by GPPM for T = L defined in (2.11) and P defined in (2.15) corresponds to the sequence generated by the generalized splitting scheme (GSS) defined in (2.17)-(2.20).

In Section 2.2, we have shown that the sequence generated by GPPM is nothing else but, under the linear transformation Q, the sequence generated by the fixed point method corresponding to the 1-co-coercive wrt R map QJ T P (see (2.9)). But for arbitrary symmetric positive semidefinite matrix P , matrices Q and R are difficult to calculate; when P is symmetric positive definite, then Q = I and R = P . Alternately by considering S such that P = S t S, we define G T S an operator easier to implement than QJ T P and having similar properties, for example, it is 1-co-coercive property wrt the usual norm. In particular, using G T S instead QJ T P , we give an alternative proof of Proposition 2.2.2.

Finally, by considering S = S 3 defined in Remark 2.4.3, one gets G L S 3 = S t 3 (L + S t 3 S 3 ) -1 S 3 which corresponds, under a reparametrization, to the classical Douglas-Rachford operator defined by

M -1 2 S t 3 (L+S t 3 S 3 ) -1 S 3 M 1 2
. In other words, the Douglas-Rachford operator and its fundamental properties of of co-coercivity and splittability will be shown to be a special case of our generalized setting based on the Lagrangian monotone inclusion.

Associated with the r × r symmetric positive semidefinite matrix P introduced in the former section, let consider a q × r matrix S satisfying

P = S t S (2.52)
and then the map G T S : IR q -→ -→ IR q defined as

G T S := S(T + S t S) -1 S t . (2.53)
It follows that

SJ T P = G T S S (2.54)
and hence, from (2.7), we get for all w, w ∈ IR r :

G T S (Sw) -G T S (Sw ), Sw -Sw ≥ G T S (Sw) -G T S (Sw ) 2 .
Since for any s, s ∈ IR q there exist w, w ∈ IR r such that S t Sw = S t s and S t Sw = S t s , we get

G T S (s) -G T S (s ), s -s ≥ G T S (s) -G T S (s ) 2 which means that G T
S is 1-co-coercive with respect to the usual norm.

The following proposition shows in particular that G T S is the Moreau-Yosida regularization of ST -1 S t . This will be used in the examples considered in this Section and in Section 2.6 (Proposition 2.6.1).

Proposition 2.4.1 Let T : IR r -→ -→ IR r be an arbitrary map, S and M two matrices of order q × r and q × q, respectively, with M invertible. For z ∈ IR q the value (ST

-1 S t + M ) -1 M z is nonempty if and only if (T + S t M -1 S) -1 S t z is nonempty. Furthermore, it holds that (ST -1 S t + M ) -1 M z = z -M -1 S(T + S t M -1 S) -1 S t z.
Proof. The proof follows from the two properties:

• x ∈ (ST -1 S t + M ) -1 M z if and only if there exists y ∈ IR m such that S t x ∈
T (y) and z -M -1 Sy = x.

• y * ∈ (T + S t M -1 S) -1 S t z if and only if exists x * ∈ IR r such that S t x * ∈ T (y * ) and z -M -1 Sy * = x * .

Similar to Proposition 2.2.1, we get the relationship between the solution set of problem (V ) and the fixed points of G T S .

Proposition 2.4.2 With the same notations as before, we have

• If z ∈ sol (V ), then Sz is a fixed point of G T S .
• If w is a fixed point of G T S , then w = Sq for some q ∈ sol (V ) ∩ (T + P ) -1 S t w.

We deduce that the set of fixed point of G T S is exactly

S(sol (V )) = {Sw : w ∈ sol (V )}.
Applying S to the sequence {w k } defined in (2.8) and considering the permutation property (2.54), we get:

Sw k+1 = ρG T S (Sw k ) + (1 -ρ)Sw k . (2.55)
This equation gives us another alternative proof of convergence of the sequence {w k } under the same conditions of Proposition 2.2.2. In fact, since G T S is 1-cocoercive and from (2.55), we have that, given w * ∈ sol (V )

Sw k -Sw * 2 - 2 -ρ ρ Sw k+1 -Sw k 2 -Sw k+1 -Sw * 2 ≥ 0. (2.56)
Since rank S t S = rank S t , the domain of G T S is equal to the domain of J T P which is closed, using this fact and from (2.56) we deduce that Sw k converges to some point b, which is a fixed point of G T S . On the other hand, using the triangular inequality and considering w := (T + P ) -1 S t b, we get

w k+1 -w ≤ ρ (T + P ) -1 S t (Sw k ) -w + |1 -ρ| w k -w .
From the continuity of J T P , we deduce the continuity of (T + P ) -1 S t = J T P S + , where S + denotes the Moore-Penrose pseudo-inverse matrix of S. Therefore we deduce that {w k } converges to w.

We now give some explicit expressions of G L S for the Lagrangian operator L and matrix S such that P = S t S, considered in Section 2.3.

Examples of co-coercive operators

C 1 = V 1 + (1 + α) 2 A t M A and C 2 = V 2 + α 2 B t M B,
where V 1 and V 2 are as (2.21) assumed positive semidefinite matrices. In (2.37)-(2.40) matrices V 1 and V 2 are associated with the additional proximal term that will be used in ADMM, which, as we have shown in Subsection 2.3.1 (Case γ = 1), is related to Shefi-Teboulle algorithm type II [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF]. We get :

P =    V 1 + (1 + α) 2 A t M A (1 + α)αA t M B (1 + α)A t (1 + α)αB t M A V 2 + α 2 B t M B αB t (1 + α)A αB M -1    .
The matrix

S 1 =     V 1 2 1 0 0 0 V 1 2 2 0 (1 + α)M 1 2 A αM 1 2 B M -1 2    
satisfies (2.52) and the corresponding map G L S 1 , that applies IR n × IR p × IR m into itself, is defined as

G L S 1 ( x, z, y) =     V 1 2 1 x V 1 2 2 z M 1 2 Ax + M 1 2 Bz + y     where x = (∂f + V 1 + A t M A) -1 (V 1 2 1 x -A t M 1 2 ( y + 2M 1 2 Bz)) z = (∂g + V 2 + B t M B) -1 (V 1 2 2 z -B t M 1 2 y). Note that G L S 1 has full domain if V 1 + A t M A and V 2 + B t M B are assumed positive definite matrices.
Remark 2.4.1 The map G L S 1 is the Douglas-Rachford operator [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF], applied to the two maps

-    V 1 2 1 0 0 I p×p M 1 2 A 0    ∂f 0 -1 - V 1 2 1 0 A t M 1 2 0 I p×p 0 and -    -I n×n 0 0 -V 1 2 2 0 M 1 2 B    0 ∂g -1 - -I n×n 0 0 0 -V 1 2 2 B t M 1 2
The corresponding sum of these two maps is exactly the dual variational map (2.5) associated with the following optimization problem

min (x 1 ,x 2 ,z 1 ,z 2 )∈F (f, 0)(x 1 , x 2 ) + (0, g)(z 1 , z 2 ) where F is the set of all (x 1 , x 2 , z 1 , z 2 ) satisfying    V 1 2 1 0 0 I p×p M 1 2 A 0    x 1 x 2 +    -I n×n 0 0 -V 1 2 2 0 M 1 2 B    z 1 z 2 = 0.
Remark 2.4.2 In the case that V 2 = 0, which corresponds to Chambolle-Pock algorithm as we showed in Subsection 2.3.1, we can restrict the map G L S 1 , and obtain the map D 1 that applies IR n × IR m into itself, where

D 1 (x, u) is V 1 2 1 (∂f + V 1 + A t M A) -1 [V 1 2 1 x -A t M 1 2 (u + 2z)] M 1 2 A(∂f + V 1 + A t M A) -1 [V 1 2 1 x -A t M 1 2 (u + 2z)] + z + u where z = M 1 2 B(∂g + B t M B) -1 B t M 1 2 (-u).
Note if B is injective, then D 1 has full domain.

The map D 1 can be obtained in the form (2.53), considering that when V 2 = 0, the matrix

S 2 = V 1 2 1 0 0 (1 + α)M 1 2 A αM 1 2 B M -1 2 satisfies (2.

52), and we obtain that

D 1 = G L S 2 .
The map D 1 can also be obtained as the Douglas-Rachford operator, applied to the two maps

- V 1 2 1 M 1 2 A (∂f ) -1 -V 1 2 1 A t M 1 2
and

- -I n×n 0 0 M 1 2 B 0 ∂g -1 - -I n×n 0 0 B t M 1 2
The corresponding sum of these two maps is exactly the dual variational map associated with the following optimization problem min

(x,z 1 ,z 2 )∈F f (x) + (0, g)(z 1 , z 2 )
where F is the set of all triples (x,

z 1 , z 2 ) satisfying V 1 2 1 M 1 2 A x + -I n×n 0 0 M 1 2 B z 1 z 2 = 0.
Remark 2.4.3 In the case V 1 = 0 and V 2 = 0, we can restrict the map G L S 1 , and obtain the map D 2 that applies IR m into itself, where

D 2 (x, u) is M 1 2 A(∂f + A t M A) -1 A t M 1 2 [-u -2z] + z + u where z = M 1 2 B(∂g + B t M B) -1 B t M 1 2 (-u).
Note that if A and B are injective, then D 2 has full domain.

The map D 2 can be obtained in the form (2.53), considering that when V 1 = V 2 = 0, the matrix

S 3 = (1 + α)M 1 2 A αM 1 2 B M -1 2 verifies (2.

52), and we obtain that

D 2 = G L S 3 .
The map D 2 can also be obtained as the Douglas-Rachford operator [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF], applied to the two maps

-M 1 2 A(∂f ) -1 (-A t M 1 2 ) and -M 1 2 B(∂g) -1 (-B t M 1 2 ).
The corresponding sum of these two maps is exactly the dual variational map (2.5) associated with the following optimization problem min

(x,y) [f (x) + g(z) : M 1 2 Ax + M 1 2 Bz = 0].
Alternatively we can consider, instead D 2 , the map

D2 := M -1 2 D 2 M 1 2 , i.e D2 (ū) = A(∂f + A t M A) -1 A t M [-ū -2z] + z + ū where z = B(∂g + B t M B) -1 B t M (-ū),
which is co-coercive w.r.t. the metric induced by M .

Example 2.4.2 Let γ = 0 (α = β). We consider in (2.15)

C 1 = (1 + (α + 1) 2 )A t M A + R and C 2 = (1 + (α + 1) 2 )B t M B,
where R is a positive semidefinite matrix. Then V 1 and V 2 in (2.21) are equal to

V 1 = A t M A + R and V 2 = B t M B.
These matrices are associated with the additional proximal term considered in (2.33)-(2.36), which, as we have shown in Subsection 2.3.1 (Case γ = 0), is related to Shefi-Teboulle algorithm type I [START_REF] Shefi | Rate of convergence analysis of decomposition methods based on the proximal method of multipliers for convex minimization[END_REF]. We get :

P =    (1 + (α + 1) 2 )A t M A + R α(2 + α)A t M B (1 + α)A t α(2 + α)B t M A (1 + (α + 1) 2 )B t M B (1 + α)B t (1 + α)A (1 + α)B M -1    .
The matrix

S 4 =    R 1 2 0 0 M 1 2 A -M 1 2 B 0 (1 + α)M 1 2 A (1 + α)M 1 2 B M -1 2    satisfies (2.

52) and hence the value

G L S 4 ( x, z, y) of the corresponding map G L S 4 , that applies IR n × IR m × IR m into itself, is    R 1 2 x M 1 2 Ax -M 1 2 Bz M 1 2 Ax + M 1 2 Bz + y    where x = (∂f + 2A t M A + R) -1 (R 1 2 x + A t M 1 2 ( z -y)) z = (∂g + 2B t M B) -1 B t M 1 2 (-z -y).
Note that G L S 4 has full domain if 2A t M A + R and 2B t M B are assumed positive definite matrices.

Remark 2.4.4 In the case that R = 0, we can restrict the map G L S 4 , and obtain the map D 3 that applies IR m × IR m into itself, where

D 3 ( z, y) is M 1 2 A(∂f + 2A t M A) -1 A t M 1 2 ( z -y) + M 1 2 (∂g + 2B t M B) -1 B t M 1 2 (-z -y) M 1 2 A(∂f + 2A t M A) -1 A t M 1 2 ( z -y) -M 1 2 (∂g + 2B t M B) -1 B t M 1 2 (-z -y) + y
The map D 3 can be obtained as the form (2.53), considering that when V 1 = A t M A and V 2 = B t M B, the matrix

S 5 = M 1 2 A -M 1 2 B 0 (1 + α)M 1 2 A (1 + α)M 1 2 B M -1 2 .
satisfies (2.52), then we obtain that D 3 = G L S 5 .

Rate of Convergence

The global rate of convergence of ADMM and other monotone operator splitting algorithms has motivated many research contributions that we cannot survey here (see [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] for example). We will recover these results for the generalized splitting scheme GSS with no further refinements (like uniform or strong convexity) and will remain in the framework of finite-dimensional spaces (see [START_REF] Attouch | Augmented lagrangian and proximal alternating direction method of multipliers in hilbert spaces. applications to games, pde's and control[END_REF] for similar results in Hilbert spaces). D. Davis and W. Yin [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] have show the ergodic and nonergodic convergence rate of the feasibility and objective function error related to the relaxed PRS and relaxed ADMM, which is a particular case of our general scheme as remarked in Subsection 2.3.1. Similarly, in this Section, without regularity assumption, we show the ergodic and nonergodic convergence rate of the constraint violations (feasibility) and objective function error related to the chain of steps (2.17) → (2.18) → (2.19) → (2.20), defined in Subsection 2.3.1, which is our main sequence associated with primal problem (P ) defined in the first section.

With the same expressions of matrices P and U defined in (2.15) and (2.16), respectively, we get the following identity by using S satisfying P = S t S and explicit expressions of P and U , (x, z, y) 2 P = S(x, z, y

) 2 = (x, z) 2 U + M 1 2 ((1 + α)Ax + (1 + β)Bz) + M -1 2 y 2 .
(2.57) Notice that for γ = 0 (β = α),

(x, z) 2 U = x 2 V 1 -A t M A + z 2 V 2 -B t M B + Ax -Bz 2 M and for γ = 1 (β = α -1), (x, z) 2 U = x 2 V 1 + z 2 V 2 .
(2.58)

Back to the sequence (2.17)-(2.20) and considering w k = (x k , z k , y k ), it holds from definition that

J L P w k = (x k+1 , zk+1 , ỹk+1 ) and w k+1 = ρJ L P w k + (1 -ρ)w k .
(2.59)

The following proposition will be used later in Subsection 2.5.2 in order to estimate an upper bound of the optimal value of problem (P ).

Proposition 2.5.1 With the same notations as before and considering w = (x, z, y) ∈ dom (f ) × dom (g) × IR m , the following inequality holds:

w k -w 2 P - 2 -ρ ρ w k+1 -w k 2 P -w k+1 -w 2 P ≥ 2ρ l(x k+1 , zk+1 , y) -l(x, z, ỹk+1 )
Proof. Let w = (x, z, y) ∈ dom (f ) × dom (g) × IR m . Since P (w k -J L P w k ) ∈ L(J L P w k ), then using Proposition 2.3.1, it holds that

J L P w k -w, P (w k -J L P w k ) ≥ l(x k+1 , zk+1 , y) -l(x, z, ỹk+1 ). ( 2

.60)

On the other hand, from the symmetry of P , it holds 2ρ

J L P w k -w, P (w k -J L P w k ) = w k -w 2 P - 2 -ρ ρ w k+1 -w k 2 P -w k+1 -w 2 P
So, replacing this last expression in (2.60), we get the desired inequality.

In particular, from the inequality of last proposition, we get

w k -w 2 P -w k+1 -w 2 P ≥ 2ρ l(x k+1 , zk+1 , y) -l(x, z, ỹk+1 ) . ( 2

.61)

This inequality will be used in Theorem 2.5.1 for approximating the optimal value of problem (P ).

We note that Proposition 2.5.1 is a general version of the inequality given in Proposition 2 of [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] which is obtained by considering A = I = -B, M = γ -1 I and P as in Remark 2.4.3, w = (x, x, 0) (which implies M -1 2 S 3 w = x), z = M -1 2 S 3 z k , and

z + = (T P RS ) λ (z) = (M -1 2 G L S 3 M 1 2 ) 2λ (M -1 2 S 3 z k ) = M -1 2 S 3 w k+1 .
Similarly, Proposition 2.5.1 is also a general version of the one given in Proposition 11 of [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF] by considering M = γI and P as in Remark 2.4.3; (x * , z * , ȳ * ) and z * fixed points of G L S 3 and (T P RS ) λ = (M

1 2 G L S 3 M -1 2 ) 2λ , respectively; w k satisfying M 1 2 S 3 w k = z k and w = (x * , z * , 0) such that M 1 2 S 3 w = M 1 2 S 3 (x * , z * , ȳ * ) -ȳ * = z * -w * where w * = J γ(-B)(∂g) -1 (-B t ) (z * ).

Bounding the fixed-point residual

The fixed-point residual of operator ρG T S + (1 -ρ)I q×q is the sequence with general term (ρG

T S + (1 -ρ)I q×q )Sw k -Sw k 2
which, from (2.55), is equal to

Sw k+1 -Sw k 2 .
Since ρ ∈ (0, 2), then ρG T S +(1-ρ)I q×q is non expansive and hence { Sw k+1 -Sw k } is non increasing. Summing over k = 0, • • • , N -1 in (2.56), we get

Sw k -Sw k-1 2 ≤ ρ (2 -ρ)k Sw 0 -Sw * 2 .
(2.62)

On the other hand, using the Jensen's inequality, we get

Sw k -Sw 0 2 ≤ 2 Sw k -Sw * 2 + 2 Sw 0 -Sw * 2 ≤ 4 Sw 0 -Sw * 2
and hence

1 N N k=1 (Sw k -Sw k-1 ) 2 = 1 N 2 Sw N -Sw 0 2 ≤ 4 N 2 Sw 0 -Sw * 2 .
(2.63)

Notice that upper bounds (2.62) and (2.63) can also be deduced respectively from Theorem 1 "Notes on Theorem 1" and Theorem 2 developed in D. Davis and W. Yin [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF].

Bounding the saddle-point gap

We consider the following ergodic sequences defined as: for N ≥ 1,

x N := 1 N N k=1 xk , z N := 1 N N k=1
zk and

y N := 1 N N k=1
ỹk .

Theorem 2.5.1 With the same notations as before, we get the following rate of convergence:

• Ergodic Convergence: for any w = (x, z, y) ∈ dom (f ) × dom (g) × IR m l(x k , z k , y) -l(x, z, y k ) ≤ 1 2ρk
Sw 0 -Sw 2 .

(2.64)

• Nonergodic Convergence: for any w * = (x * , z * , y * ) ∈ sol (V L ) l(x k+1 , zk+1 , y * ) -l(x * , z * , ỹk+1 ) ≤ 1 + |1 -ρ| ρ ρ(2 -ρ)(k + 1) Sw 0 -Sw * 2 . (2.65)
Proof. Summing (2.61) over k = 0, • • • , N -1, and applying the Jensen's inequality to the convex functions l(•, •, y) -l(x, z, •) for arbitrary fixed elements x ∈ dom (f ), z ∈ dom (g) and y ∈ IR m , where l is the lagrangian function defined in (2.10) of Section 2.3, we deduce the desired ergodic convergence.

Given w * ∈ sol (V L ) and considering w = w * in (2.60), we get

G L S Sw k -Sw * , Sw k -G L S Sw k ≥ l(x k+1 , zk+1 , y * ) -l(x * , z * , ỹk+1 ) ≥ 0
and hence, from the Cauchy-Schwarz inequality and (2.55), we obtain

1 ρ G L S Sw k -Sw * Sw k+1 -Sw k ≥ l(x k+1 , zk+1 , y * ) -l(x * , z * , ỹk+1 ). ( 2

.66)

On other hand, from (2.55) and since { Sw k+1 -Sw * } is non increasing, we get

G L S Sw k -Sw * = 1 ρ (Sw k+1 -Sw * )+(1- 1 ρ )(Sw k -Sw * ) ≤ 1 + |1 -ρ| ρ Sw 0 -Sw * .
So, replacing this last expression and inequality (2.62) in expression (2.66), we deduce the desired nonergodic convergence.

Bounding the constraint violation

We consider, for N ≥ 1,

x N := 1 N N k=1 x k-1 and z N := 1 N N k=1 z k-1 .
We get the following result Theorem 2.5.2 With the same notations as before, for any w * ∈ sol (V L ), we get the following rate of convergence:

• Ergodic Convergence:

(x k -x k , z k -z k ) 2 U + Ax k + Bz k 2 M ≤ 4 ρ 2 k 2 Sw 0 -Sw * 2 .
• Nonergodic Convergence:

(x k -x k-1 , zk -z k-1 ) 2 U + Ax k + B zk 2 M ≤ 1 (2 -ρ)ρk Sw 0 -Sw * 2 .
Proof. From (2.59) we have

w k -w k-1 = ρ(x k -x k-1 , zk -z k-1 , ỹk -y k-1
) and hence, from (2.19), we get

w k -w k-1 = ρ(x k -x k-1 , zk -z k-1 , M [(1+α)Ax k-1 +(1+β)Bz k-1 -αAx k -βB zk ]). (2.67) Summing over k = 1, • • • , N , we obtain 1 N N k=1 (w k -w k-1 ) = ρ(x N -x N , z N -z N , M [(1+α)A x N +(1+β)B z N -αAx N -βBz N ]).
Then from (2.57), we get

1 ρ 2 1 N N k=1 (w k -w k-1 ) 2 P = (x N -x N , z N -z N ) 2 U + Ax N + Bz N 2 M
and hence, given w * ∈ sol (V L ), we deduce from (2.63) the ergodic rate of convergence for constraint violations.

Using (2.67), from (2.57), we get

1 ρ 2 w k -w k-1 2 P = (x k -x k-1 , zk -z k-1 ) 2 U + Ax k + B zk 2 M
and hence, from (2.62), we deduce the nonergodic rate of convergence for constraint violations.

We note that the particular case γ = 1, V 1 = 0 and V 2 = 0, which implies that U = 0, the two terms (

x k -x k , z k -z k ) 2 U and (x k -x k-1 , zk -z k-1 ) 2
U of inequalities in Theorem 2.5.2 are null and hence we recover the Theorem 15 of [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF].

Application to the decomposition of multiblock optimization problems

To conclude our study, we consider the application of the generalized scheme GSS to the decomposition of some block structured convex optimization problems. For i ∈ {1, ..., q}, let f i : IR n i → IR and g : IR m → IR are proper lsc convex functions, A i and B matrices of order p × n i and p × m, respectively. We consider the following S-Model problem: inf

(x 1 ,••• , xq, z) q i=1 f i (x i ) + g(z) s.t q i=1 A i x i + Bz = 0.
This problem has been analyzed by many authors (see [START_REF] Lenoir | A survey of monotone operator splitting methods for the decomposition of convex programs[END_REF] for instance). We rewrite it into two different forms, (B 1 ) and (B 2 ), but with the same structure as (BP ) defined below, then we rewrite (BP ) as problem (P ) also defined below. Finally, we apply the algorithm (2.37)- (2.40) to this last problem.

The S-Model problem is equivalent to inf

(x 1 ,••• , xq, z) q i=1 f i (x i ) + g(z) + δ {0} ( q i=1 A i x i + Bz). (B 1 )
In this formulation the function g can be viewed as a function f i . The associated dual problem of (B 1 ) is inf

y * q i=1 (f * i • A t i )y * + (g * • B t )y * . (Ds)
Now, by considering n = q i=1 n i and f : IR n → IR defined as f (x) := q i=1 f i (x i ), the problem (Ds) can be written as

inf y * f * •    A t 1 . . . A t q    y * + (g * • B t )y * .
This is a composite problem whose associated dual problem is inf

(x 1 ,••• , xq) q i=1 f i (x i ) + (g * • B t ) * • (- q i=1 A i x i ). (B 2 )
We observe that in this last problem we reduce the number of variables considered in the S-Model problem and the function g acts now as regularization function.

Using the same notations as before, we define a problem having the same structures as problems (B 1 ) and (B 2 ):

V P = inf (x 1 ,••• , xq) q i=1 f i (x i ) + (g * • B t ) * • ( q i=1 A i x i ).
(BP )

In order to apply the splitting algorithm to problem (BP), we reformulate it to an appropriate optimization problem. To do it, consider

K := I p×p • • • I p×p ∈ IR p×pq and A :=    A 1 . . . A q    ∈ IR pq×n .
So, problem (BP ) can be formulated as

inf x∈IR n , z∈IR pq f (x) + (g * • B t ) * • Kz : Ax -z = 0 . (P )
Notice that this last formulation problem has a good separable structure.

We apply to problem (P ) the algorithm (2.37)-(2.40) developed in Subsection 2.3.1 (Case γ = 1). We assume that g verifies the following identity

∂[(g * • B t ) * • K] = K t (B(∂g) -1 B t ) -1 K.
The saddle-point problem of (P ) is

Find (x, z, ȳ) ∈ IR n × IR pq × IR pq such that 0 ∈ L(x, z, ȳ) (VL)
where L is the maximal monotone map defined on

IR n × IR pq × IR pq as L(x, z, y) :=    ∂f (x) K t (B(∂g) -1 B t ) -1 Kz 0    +    0 0 A t 0 0 I -A I 0       x z y    .
For i ∈ {1, ..., q}, let M i be an p × p symmetric positive definite matrix and Q i be an n i × n i symmetric positive semi-definite matrix.

In order to take advantage of the separability of f , we take V 1 = diag([Q 1 , ..., Q q ]) and M = diag([M 1 , ..., M q ]), and we consider V 2 = 0 in order to calculate z k+1 using alone the resolvent of ∂g. So, the related algorithm (2.37)-(2.40) takes the following structure:

zk+1 = (K t (B(∂g) -1 B t ) -1 K + M ) -1 (M Ax k + y k ) (2.68) ỹk+1 = y k + M (Ax k -zk+1 ) (2.69) xk+1 = (∂f + A t M A + V 1 ) -1 (V 1 x k + A t M zk+1 -A t ỹk+1 ) (2.70) (x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ) (2.71)
Because of the diagonal structure of expression (2.70) the calculation of xk+1 is realized in parallel: for i ∈ {1, • • • , q},

xk+1 i = (∂f i + A t i M i A i + Q i ) -1 (Q i x k i + A t i M i zk+1 i -A t i ỹk+1 i ).
Now, in order to calculate zk+1 , the following identity is relevant:

Proposition 2.6.1 With the same notations as before, the following identity holds

(K t (B(∂g) -1 B t ) -1 K + M ) -1 M = I -M -1 K t Σ(I -B(∂g + B t ΣB) -1 B t Σ)K
where Σ is a p × p matrix defined by

Σ := (KM -1 K t ) -1 = q i=1 M -1 i -1
.

Proof. From Proposition 2.4.1, we have

(K t (B(∂g) -1 B t ) -1 K + M ) -1 M = I -M -1 K t (B(∂g) -1 B t + KM -1 K t ) -1 K
and hence by combining it with the following identity

(B(∂g) -1 B t + Σ -1 ) -1 Σ -1 = I -ΣB(∂g + B t ΣB) -1 B t
obtained also from Proposition 2.4.1, we get the desired identity.

So, using the identity of this last proposition, we can obtain an equivalent expression of ỹk+1 in (2.69) but with a more tractable expression for computational purpose :

ỹk+1 = K t Σ(I -B(∂g + B t ΣB) -1 B t Σ)K(Ax k + M -1 y k ).
(2.72)

It follow in particular that ỹk+1 ∈ range K t and, by considering y k ∈ range K t in (2.71), we have that y k+1 ∈ range K t and hence all the block components of ỹk+1 (similarly of y k+1 ) are equal. We denote by ỹk+1 c (resp y k+1 c ) such a block component of ỹk+1 (resp y k+1 ). Then,

ỹk+1 c = Σ(I -B(∂g + B t ΣB) -1 B t Σ)K(Ax k + M -1 K t y k c ) (2.73)
By denoting

ζ k+1 := (∂g + B t ΣB) -1 B t (Σ q j=1 (A j x k j ) + y k c )
we obtain, from (2.73),

ỹk+1 c = y k c + Σ( q j=1 (A j x k j ) -Bζ k+1 ). ( 2

.74)

On the other hand, from (2.69), we get

zk+1 = Ax k + M -1 K t (y k c -ỹk+1 c )
which combining with (2.74), we deduce that for i ∈ {1, • • • , q},

zk+1 i = A i x k i -M -1 i Σ q j=1 (A j x k j ) -Bζ k+1 .
Therefore we obtain the following algorithm, called "Proximal Multi-block Algorithm".

Proximal Multi-block Algorithm (PMA)

For i ∈ {1, • • • , q} set Q i ∈ IR n i ×n i symmetric positive semi-definite, M i ∈ IR p×p symmetric positive definite. Set Σ = q i=1 M -1 i -1 .
Then for an arbitrary (x 0 , z 0 , y 0

c ) ∈ IR n × IR pq × IR p Step 1. Find ζ k+1 such that ζ k+1 = argmin g(w) + 1 2 Bw - q j=1 (A j x k j ) -Σ -1 y k c 2 Σ
Step 2. Find zk+1 For all i ∈ {1, ..., q} do Find zk+1 i such that

zk+1 i = A i x k i -M -1 i Σ q j=1 (A j x k j ) -Bζ k+1 .
Step 3. Find ỹk+1

c such that ỹk+1 c = y k c + Σ q j=1 (A j x k j ) -Bζ k+1 .
end for

Step 4. Find xk+1 For all i ∈ {1, ..., q} do Find xk+1

i such that

xk+1 i = argmin f i (x i ) + 1 2 A i x i -zk+1 i + M -1 i ỹk+1 c 2 M i + 1 2 x i -x k i 2 Q i end for
Step 5. Find (x k+1 , z k+1 , y k+1 c

) (x k+1 , z k+1 , y k+1 c ) = ρ(x k+1 , zk+1 , ỹk+1 c ) + (1 -ρ)(x k , z k , y k c ).
The next proposition gives conditions in order to guarantee the convergence of PMA. The proof is a direct consequence of Proposition 2.3.5. Proposition 2.6.2 Let ρ ∈ (0, 2). For i ∈ {1, ..., q}, assume that Q i ∈ IR n i ×n i and M i ∈ IR p×p are symmetric, with Q i positive semi-definite and M i positive definite such that Q i + A t i M i A i is positive definite. If sol (VL) is nonempty, then for an arbitrary (x 0 , z 0 , y 0 c ) ∈ IR n × IR pq × IR p , the sequence (x k , z k , K t y k c ) generated by (PMA) converges to some element of sol (VL).

The Separable Augmented Lagrangian Algorithm (SALA) with multidimensional scaling has been proposed in [START_REF] Dussault | Separable augmented lagrangian algorithm with multidimensional scaling for monotropic programming[END_REF] to solve a special case of the S-Model where g = 0 and B = 0 . This algorithm can be recovered if instead of applying the algorithm (2.37)-(2.40) to problem (P ), we consider the algorithm (2.41)-(2.43) with V 1 = V 2 = 0. Therefore SALA is a particular version of (PMA). The advantages of (PMA) are twofold: 1) the inclusion of the relaxing term ρ ∈ (0, 2), which enables the accelaration of the algorithm, and 2) the additional proximal term x i -x k i 2 Q i considered in the subproblems of Step 4, which improves the strong convexity of the proximal subproblem when we choose an adequate matrix Q i . More specifically, considering σ i and τ i positive numbers holding σ i τ i A i 2 ≤ 1 and choosing M i and Q i matrices defined as

M i = σ i I p×p and Q i = τ -1 i I n i ×n i -σ i A t i A i , 56 the conditions about matrices Q i , M i and Q i + A t i M i A i in
i = argmin f i (x i ) + 1 2τ i x i -x k i -τ i [σ i A t i zk+1 i -σ i A t i A i x k i -A t i ȳk+1 ] 2
which has an explicit solution in some particular cases, for instance f i (x i ) = x i 1 .

57

Chapter 3

Decomposition techniques

In the last section of previous chapter we have developed a general algorithm termed Proximal Multi-block Algorithm (PMA) in order to solve a S-Model problem as in Section 2.6. This algorithm uses the proximal step of all the family {f i } i=1,••• ,q through a parallel processing for each iteration.

The special case considering g ≡ 0 and B ≡ 0 in the previous S-model, which will be termed separable model with coupling constraint (SMCC) is considered in Section 3.2. We propose an alternative way to use the proximal step of each f i , separating the problem into two sub-block problems and considering the proximal step to one sub-block and then (at a linear combination of the preceding values) the proximal step is found for the other sub-block, both in parallel processing.

For this purpose, we first study in Section 3.1 the separable model with coupling variable (SMCV) using the duality scheme developed in Chapter 1, finding an adequate formulation for that problem, allowing recovery two know algorithms with multi-scaling parameters, and their relationship. We also show the numerical behavior of these two algorithms. Since SMCC can be formulated as a SMCV, we apply the results obtained for SMCV to SMCC, getting another way to recover SALA and DSALA, the last also contained in P M A.

In Section 3.3, we get two splitting algorithms, one for SMCV and the other for SMCC. Each algorithm separates the problem into two arbitrary sub-block problems considering the proximal step of all functions corresponding to each sub-block in parallel processing similar to the one obtained in the precedent section.

All found algorithms in this chapter are consequences of finding adequate formulations of the original problem and then applying popular algorithms. In this vein, in the last section we study an especial block optimization problem, that will be used for solving a stochastic optimization model problem, reformulating it and applying the multi-scaling ADMM (considering V 1 = 0 and V 2 = 0 in (2.41)-(2.43)).

In [START_REF] Chen | The direct extension of admm for multi block convex minimization problems is not necessarily convergent[END_REF] the authors show that "The direct extension of ADMM is not necessarily convergent". The authors in [START_REF] Wang | Solving multiple block separable convex minimization problems using two-block alternating direction method of multipliers[END_REF] showed an alternative way to deal with a multi-block optimization problem, reformulating it into two-block problem in order to apply the ADMM.

The main difference of the corresponding algorithm found by these authors regarding the one proposed in Section 3.2 of this chapter is the step order of subproblems intervening in the algorithm.

On the other hand, in [START_REF] Eckstein | A simplified form of block-iterative operator splitting and an asynchronous algorithm resembling the multi-block alternating direction method of multipliers[END_REF] the author gives another way to solve problem SMCC using projective splitting methods.

The separable model with coupling variable (SMCV)

In this section we present the sum problem (or consensus problem as termed in [START_REF] Boyd | Distributed optimization and statistical learning with the alternating direction method of multipliers[END_REF]). In this book the authors give important references on recent advances and applications regarding the SMCV.

For every i ∈ {1, • • • , q}, let f i : IR n → IR be proper lsc convex. We consider the following SMCV inf

x∈IR n q i=1 f i (x) (S)
which by defining f : IR nq → Ī R as f (z) := q i=1 f i (z i ) for z = (z 1 , • • • , z q ) with z i ∈ IR n , and W the nq × n matrix defined by W = (I n×n • • • I n×n ) t , that problem can be set as

inf x∈IR n f • W (x) (P V s )
which is termed Primal projected problem regarding to the following optimization problem with linear subspace constraint

inf z q i=1 f i (z i ) (3.1) s.t z ∈ V (3.2)
where

V := range (W ) = {z = (z 1 , • • • , z q ) ∈ IR nq : z 1 = • • • = z q }.
Clearly, this problem is also equivalent to inf z∈IR nq f (z) + δ V (z).

( P s )

Notice that if we apply the Douglas-Rachford's method (DRM) to problem ( P s ), we obtain the so called [START_REF] Lenoir | A survey of monotone operator splitting methods for the decomposition of convex programs[END_REF] "Proximal Decomposition algorithm (PDA)".

On the other hand, by making V = ker(K) where K is the n(q -1) × nq matrix defined by K = I n(q-1)×n(q-1) -D with D being the n(q -1) × n matrix defined by D = (I n×n • • • I n×n ) t , we obtain the following problem so called Dual projected problem:

inf u * ∈IR n(q-1) f * • K t u * . (D V s )
Now, by considering h : IR n(q-1) → IR defined as h(y)

:= q-1 i=1 f i (y i ) for y = (y 1 , • • • , y q ) with y i ∈ IR n , problem (D V s ) can be formulated as inf u * ∈IR n(q-1) h * (u * ) + f * q • -D t u * . ( DV s )
Since D t D = (q -1)I, applying the DRM to problem ( DV s ), we obtain a suitable decomposition algorithm called "Dual Proximal Decomposition algorithm (DPDA)".

By considering

V = ker(K) in problem (3.1)-(3.2), it is transformed in inf (u,s) q-1 i=1 f i (u i ) + f q (s) s.t u i -s = 0, i = 1, • • • , q -1
which is exactly the dual problem of (D V s ). In [START_REF] Boyd | Distributed optimization and statistical learning with the alternating direction method of multipliers[END_REF] the authors termed this problem "Global variable consensus with regularization", where f q represents the regularization function. In order to solve such a problem, they apply ADMM, which is equivalent to apply the Douglas-Rachford's method to problem ( DV s ).

Summarizing the previous discussions, PDA is DRM applied to problem (P s ), and DPDA is DRM applied to problem (D s ).

We now describe these two algorithms:

Proximal Decomposition algorithm (PDA)

For every i ∈ {1, • • • , q -1} set M i ∈ IR n×n symmetric positive definite. Then for an arbitrary (x 0 , y 0 , u 0 ) ∈ IR nq × IR n × IR nq ,

Step 1. Find x k+1 For all i ∈ {1, ..., q} do find x k+1 i the optimal value of problem inf

x f i (x) + 1 2 x -y k + M -1 i u k i 2 M i end for
Step 2. Find y k+1 such that

y k+1 = q i=1 M i -1 q j=1 (M i x k+1 i + u k i )
Step 3. Find u k+1 For all i ∈ {1, ..., q} do find u k+1 i such that

u k+1 i = u k i + M i (x k+1 i -y k+1 )
end for Dual Proximal Decomposition algorithm (DPDA)

For every i ∈ {1, • • • , q -1} set M i ∈ IR n×n symmetric positive definite matrix.
Then for an arbitrary (x 0 , y 0 , u 0 ) ∈ IR n(q-1) × IR n × IR n(q-1)

Step 1. Find x k+1 For all i ∈ {1, ..., q -1} do find x k+1 i the optimal value of the problem inf

x f i (x) + 1 2 x -y k + M -1 i u k i 2 M i end for
Step 2. Find y k+1 the optimal value of the problem

inf x    f q (x) + 1 2 x - q-1 i=1 M i -1 q-1 j=1 (M i x k+1 i + u k i ) 2 q-1 i=1 M i   
Step 3. Find u k+1 For all i ∈ {1, ..., q -1} do find u k+1 i such that

u k+1 i = u k i + M i (x k+1 i -y k+1 )
end for Remark 3.1.1 If f q = 0, then DPDA becomes PDA applied to problem

inf x∈IR n q-1 i=1 f i (x).
So, apply PDA to problem (S) is equivalent to apply DPDA to the following problem

inf x∈IR n q i=1 f i (x) + 0(x).
Example 3.1.1 (Numerical illustration for the sum of three operators) We apply the methods PDA and DPDA in order to solve the following problem

min x=(x 1 ,x 2 )∈IR 2 1 0 0 2 x, x + 8 0 0 1 x, x + 5 0 0 7 x, x
which evidently x * = 0 is the unique solution.

We consider, for i ∈ {1, 2, 3}, M i = λI in the algorithm of (P DA); and for j ∈ {1, • • • , q -1}, we set M j = λI in the algoritm of (DP DA).

The graph in (fig 3.1) describes the relationship of parameter λ for PDA and DPDA versus the necessary number of iterations in order to get an approximation of the optimal value with an error ( x k 1 -x * ) less than 10 -8 , 10 -25 and 10 -60 respectively.

Separable model with coupling constraints (SMCC)

For i ∈ {1, • • • , q}, let f i : IR n i → IR be proper lsc convex function and A i be a p × n i matrix. We consider the SMCC which can be expressed as

inf x q i=1 f i (x i ) s.t q i=1 A i x i = 0.
This is a minimization problem over a linear subspace whose dual projection problem is inf

u * ∈IR p q i=1 (f * i • A t i )u * (D V sm )
which is a SMCV as described in the previous section. So, if we apply the algorithms PDA and DPDA to this problem, we obtain respectively the Separable Augmenting Lagrangian Algorithm (SALA) and the Dual Separable Augmenting Lagrangian Algorithm (DSALA). In order to apply PDA or DPDA we will use the equivalent expression of the Moreau envelope function to each composite function f i • A t i as showed in Proposition 2.4.1.

We now describe the algorithms SALA and DSALA.

Separable Augmenting Lagrangian Algorithm (SALA)

For every i ∈ {1, • • • , q}, set M i ∈ IR p×p symmetric positive definite. Set Σ = ( q i=1 M -1 i ) -1 (similarly to the given in Proposition 2.6.1). Then for an arbitrary (x 0 , y 0 , u 0 ) ∈ IR

q i=1 n i × IR p × IR pq , Step 1. Find x k+1
For all i ∈ {1, ..., q} do find x k+1 i the optimal value of the problem inf

x f i (x) + 1 2 A i x -M -1 i y k + u k i 2 M i end for
Step 2. Find y k+1 such that

y k+1 = y k -Σ( q j=1 (A i x k+1 i ))
Step 3. Find u k+1 For all i ∈ {1, ..., q} do find u k+1 i such that

u k+1 i = -A i x k+1 i + M -1 i Σ( q i=1 (A i x k+1 i ))
end for

Dual Separable Augmenting Lagrangian Algorithm (DSALA)

For i ∈ {1, • • • , q -1}, let M i ∈ IR p×p symmetric positive definite. Set Σ = ( q-1 i=1 M -1 i ) -1 .
Then for an arbitrary (x 0 , y 0 , u 0 ) ∈ IR q-1 i=1 n i × IR p × IR p(q-1)

Step 1. Find x k+1 For all i ∈ {1, ..., q -1} do find x k+1 i the optimal value of the problem inf

x f i (x) + 1 2 A i x -M -1 i y k + u k i 2 M i end for
Step 2. Calculate w k+1 the optimal value of the problem

inf x f q (x) + 1 2 A q x -Σ -1 y k + q-1 j=1 (A i x k+1 i ) 2 Σ
Step 3. Find y k+1 such that

y k+1 = y k -Σ( q-1 j=1 (A i x k+1 i ) + A q w k+1 )
Step 4. Find u k+1 For all i ∈ {1, ..., q -1} do find u k+1 i such that

u k+1 i = -A i x k+1 i + M -1 i Σ( q-1 i=1 (A i x k+1 i ) + A q w k+1 )
end for Remark 3.2.1

• If A q = 0, then DSALA becomes SALA applied to problem

inf x q-1 i=1 f i (x i ) s.t q-1 i=1 A i x i = 0.
• For i ∈ {1, ..., q}, let f i : IR n i → IR and g : IR p → IR proper lsc convex functions, and A i be p × n i matrix. We consider the Multi-block problem

min x q i=1 f i (x i ) + g( q i=1
A i x i ) also termed in [START_REF] Boyd | Distributed optimization and statistical learning with the alternating direction method of multipliers[END_REF] as Sharing problem, where f i is called local cost function and g the shared objective.

It is clear that Multi-block problem includes SMCC problem by considering g = δ {0} . Conversely, we can formulate the Multi-block problem as a SMCC by making x q+1 = q i=1 A i x i :

inf x q i=1 f i (x i ) + g(x q+1 ) s.t q i=1 A i x i -x q+1 = 0.
So, by applying the SALA or DSALA algorithms to this last problem, we get a splitting algorithm for Multi-block problem.

Example 3.2.1 (Numerical illustration for the sum of operators) We consider the algorithms SALA and DSALA in order to solve the following constrained problem min

x 1 ,x 2 ,x 3 1 0 0 2 x 1 , x 1 + 8 0 0 1 x 2 , x 2 + 5 0 0 7 x 3 , x 3 subject to 1 -1 0 2 x 1 + 3 4 2 1 x 2 + 1 3 -1 5 x 3 = 0 whose solution is x * 1 = x * 2 = x * 3 = 0.

The separable model with coupling variables

Let r and s be two positive numbers such that r + s = q and r ≥ s. We separate the block of functions into two sub-blocks: the first one consisting of r functions {f i } i=1,••• ,r and the second one consisting of s functions {f i+r } i=1,••• ,s , and the matrices coupling these two sub-blocks: B 1 and B 2 of order n(q -1)×nr and n(q -1)×ns, respectively, defined as

B 1 =         I n×n C . . . C I n(r-s)×n(r-s)         and B 2 =         -C . . . -C -I n×n -U         with C = I n×n I n×n t
and U a n(r-s)×n matrix such that

U t = I n×n • • • I n×n . Since V = ker B 1 B 2 , then problem (3.1)-(3.2) can be equivalently formu- lated as inf (x,z) r i=1 f i (x i ) + s j=1 f j+r (z j ) (3.3) s.t B 1 x + B 2 z = 0. (3.4)
On the other hand, the diagonal structure of B t 1 B 1 and B t 2 B 2 allows to solve in parallel processing the proximal step of each sub-block when the algorithm (2.41)-(2.43) with V 1 = 0 and V 2 = 0, is applies to problem (3.3)- (3.4). Explicitly B t 1 B 1 and B t 2 B 2 are respectively equal to:

        I n×n 2I n×n . . . 2I n×n I n(r-s)×n(r-s)         and       2I n×n . . . 2I n×n ηI n×n      
where η = r -s + 1.

So, applying the algorithm (2.41)-(2.43) with V 1 = 0, V 2 = 0 and M = λ (which is a slight variant of ADMM) to problem (3.3)-(3.4) by considering f (x) = r i=1 f i (x i ) and g(z) = s j=1 f j+r (z j ), we get the following algorithm:

x k+1 = (∂f + λB t 1 B 1 ) -1 (λB t 1 B 2 z k + B t 1 u k ) z k+1 = (∂g + λB t 2 B 2 ) -1 (λB t 2 B 1 x k+1 + B t 2 u k ) u k+1 = u k + λ(B 1 x k+1 + B 2 z k+1 )
Step 3. Find u k+1 For all i ∈ {1, ..., q -1} do if

1 ≤ i ≤ 2s -2 then a i = x k+1 (i-(i mod 2))/2+1 , b i = -z k+1 (i+(i mod 2))/2 end if If i ≥ 2s -1 then a i = x k+1 i-s+1 , b i = -z k+1 s end if Calculated u k+1 i such that u k+1 i = u k i + λ(a k+1 i + b k+1 i ) end for

The separable model with coupling constraints

Since SMCC can be formulated as a SMCV (D V sm ), we can get from the previous algorithms given for the last model an splitting algorithm for the SMCC which separates it into two sub-block problems: the first one corresponding to {f i } i=1,••• ,r and the second one corresponding to {f r+j } j=1,••• ,s . Then, the proximal step of all functions corresponding the first sub-block is obtained by parallel processing, and then, using these values, the proximal step of the other functions is also obtained in parallel processing.

So, the SMCC can also be formulated as inf

(x,z) r i=1 (f * i • A t i )(x i ) + s j=1 (f * r+j • A t r+j )(z j ) (3.6) s.t B 1 x + B 2 z = 0 (3.7)
where B 1 and B 2 are the matrices corresponding to expression (3.4).

Similarly to algorithm (2sb-PDA), applying the algorithm (2.41)-(2.43) with V 1 = 0, V 2 = 0 and M = λ to this new formulation problem we get the following algorithm called (2sb-SALA).

Two Sub-blocks Separable Augmenting Lagrangian Algorithm (2sb-SALA)

Set the finite sequence {α i } i=1,••• ,q defined in (3.5) and λ a positive real number.

Then for an arbitrary (x 0 , z 0 , u 0 ) ∈ IR rn × IR sn × IR n(q-1)

Step 1. Find x k+1 For all i ∈ {1, ..., r}

do if i = 1 then ξ i = -z k i , µ i = u k i , end if if 2 ≤ i ≤ s then ξ i = -z k i-1 -z k i , µ i = u k 2(i-1) + u k 2i-1 , end if if i ≥ s + 1 then ξ i = -z k s , µ i = u k i+s-1 , end if find xk+1 i the optimal value of the problem inf x f i (x) + 1 2λα i A i x + λξ i + µ i 2 calculated x k+1 = -α -1 i ξ i -(λα i ) -1 µ i -(λα i ) -1 A i xk+1 i end for
Step 2. Find z k+1 For all j ∈ {1, ..., s} do if j ≤ s -

1 then ζ j = -x k+1 i -x k+1 j+1 , ν j = -u k 2j-1 -u k 2j , end if if j = s then ζ i = -r i=s x k+1 i , ν i = -r i=s u k i+s-1 , end if find zk+1 i the optimal value of the problem inf z f r+j (z) + 1 2λα r+j A r+j z + λζ j + ν j 2 calculated x k+1 = -α -1 j+r ζ j -(λα r+j ) -1 ν j -(λα r+j ) -1 A r+j zk+1 j end for Step 3. Find u k+1 For all i ∈ {1, ..., q -1} do if 1 ≤ i ≤ 2s -2 then a i = x k+1 (i-(i mod 2))/2+1 , b i = -z k+1 (i+(i mod 2))/2 end if If i ≥ 2s -1 then a i = x k+1 i-s+1 , b i = -z k+1 s end if Calculated u k+1 i such that u k+1 i = u k i + λ(a k+1 i + b k+1 i ) end for

Multi-block optimization problem

We consider the following problem

min x=(x 1 ,••• , xq) q i=1 f i (x i ) + g( q i=1 A i x i ) + s(x) (P sc )
where for i ∈ {1, ..., q}, f i : IR n i → IR , g : IR p → IR and s : IR n → IR (n = q i=1 n i ) are proper lsc convex functions, and A i are matrices of order p × n i .

By considering f : IR n → IR defined as f (x) = q i=1 f i (x i ) and matrices K and A of order p × pq and pq × n, respectively, defined by

K := I p×p • • • I p×p and A :=    A 1 . . . A q    ,
then (P sc ) can be set as min

x f (x) + (g • K, s) A I n×n x.
For every i ∈ {1, ..., q}, set M i an p × p symmetric positive definite matrix and Q i an n i × n i symmetric positive definite matrix. We define the blocks diagonal matrices

M = diag (M 1 , • • • , M q ) and Q = diag (Q 1 , • • • , Q q ).
We apply the algorithm (2.41)-(2.43) with V 1 = 0 and V 2 = 0, which is equivalent to the scaling ADMM. So, in order to take advantage of the separability of ∂f we consider the matrix of scaling defined as M = diag ( M , Q).

x k+1 = (∂f + Q + A t M A) -1 (A t M z k 1 + Qz k 2 -A t u k 1 -u k 2 ) (3.8) z k+1 1 = (K t ∂gK + M ) -1 ( M Ax k+1 + u k 1 ) (3.9) z k+1 2 = (∂s + Q) -1 (Qx k+1 + u k 2 ) (3.10) u k+1 1 = u k 1 + M (Ax k+1 -z k+1 1 ) (3.11) u k+1 2 = u k 2 + Q(x k+1 -z k+1 2 ) (3.12)
Then using similar techniques described in Section 6 of Chapter 2, we get the following algorithm called "Separable Multi-block for sum of three blocks function (SMS3BF)" (SMS3BF)

For every i ∈ {1, • • • , q}, set M i and Q i symmetric positive definite matrices of order p × p and n i × n i , respectively, and Σ = ( q i=1 M -1 i ) -1 . Then for arbitrary (x 0 , y 0 , z 0 , v0

, u 0 ) ∈ IR n × IR pq × IR n × IR p × IR n , Step 1. Find x k+1
For all i ∈ {1, ..., q} do find x k+1 i such that

x k+1 i = (∂f i + Q i + A t i M i A i ) -1 (A t i M i y k i + Q i z k i -A t i vk -u k i )
end for

Step 2. Calculate w k+1

w k+1 = (∂g + Σ) -1 (Σ q j=1 (A j x k+1 j ) -vk )
Step 3. Find y k+1 For all i ∈ {1, ..., q} do find y k+1 i such that

y k+1 i = A i x k+1 i -M -1 i Σ( q j=1 (A j x k+1 j ) + w k+1 ).
end for

Step 4. Calculate s k+1 z k+1 = (∂s + Q) -1 (Qx k+1 + u k )

Step 5. Calculate vk+1 vk+1 = vk + Σ( q j=1 (A j x k+1 j ) + w k+1 ).

Step 6. find u k+1 For all i ∈ {1, ..., q} do find u k+1 i such that

u k+1 i = u k i + Q i (x k+1 i -z k+1 i )
end for Remark 3.4.1 If s ≡ 0 in problem (P sc ), then this problem becomes a multi-block problem

min x q i=1 f i (x i ) + g q i=1 A i x i
which is treated in Section 6 of Chapter 2, considering B = -I. In that case, the k-th iteration in Step 4 satisfies

Qz k+1 = Qx k+1 + u k
and hence from Step 6, we get u k+1 i = 0 for all i ∈ {1, ..., q}. So, the (k + 1)-th iteration in Step 1, becomes

x k+2 i = (T i +Q i +A t i D i A i ) -1 (A t i D i y k+1 i +Q i x k+1 i -A t i vk+1
), z k+2 = x k+2 and u k+2 = 0 implying that the variables z and u in the algorithm turn out obsolete. The resulting algorithm becomes equivalent to algorithm (P M A).

Aplication to a stochastic problem

We consider a stochastic problem with finite scenarios, which can be reformulated as a problem defined over a Euclidean linear space, having the same structure of problem (P sc ) considered at the beginning of this section. Then we apply (SMS3BF) in order to obtain a splitting algorithm for such a stochastic problem. Let consider a finite set Ξ of scenarios and a corresponding positive probability function p. We also consider the linear space L consisting of all mapping from Ξ to

IR n := IR n 1 × • • • × IR n T ,
and the following inner product related

X, Y L = ξ∈Ξ p(ξ) X(ξ), Y (ξ) (3.13)
Set E := card(Ξ) and G := {1, • • • , T }. We consider the linear space

IR nE := IR n 1 E × • • • × IR n T E
, and the following related inner product

((x ξ t ) ξ∈Ξ ) t∈G , ((y ξ t ) ξ∈Ξ ) t∈G IR nE = ξ∈Ξ p(ξ) (x ξ t ) t∈G , (y ξ t ) t∈G .
There is a relationship between IR nE and L through the following isomorphic mapping W :

(IR nE , •, • IR nE ) → (L, •, • L ) such that for x := ((x ξ t ) ξ∈Ξ ) t∈G , the value W (x) ∈ L satisfies W (x)(ξ) = (x ξ 1 , • • • , x ξ T ) for all ξ ∈ Ξ.
For each t ∈ G, we consider A t , a partition of Ξ such that A t+1 is a refinement of A t . The nonanticipativity subspace of L is defined as (see (1.3) in [START_REF] Rockafellar | Scenarios and policy aggregation in optimization under uncertainty[END_REF])

N := {X ∈ L : X t is constant on each A ∈ A t for t ∈ G}.
For t ∈ G and ξ ∈ Ξ, we consider B ξ t an m ξ × n t matrix and then C(ξ

) := ker [ B ξ 1 • • • B ξ T ] .
The following stochastic optimization problem is considered

min X∈L E ξ T t=1 g t (X t (ξ), ξ) : s.t X ∈ N and T t=1 B ξ t X t (ξ) = 0, ∀ξ ∈ Ξ . (SP )
Reformulating this problem in the Euclidean linear space IR nE , we get min

x=((x ξ t ) ξ∈Ξ ) t∈G ξ∈Ξ p(ξ) T t=1 g t (x ξ t , ξ) + ξ∈Ξ δ {0m ξ } ( T t=1 B ξ t x ξ t ) + δ W -1 N (x)
This last problem can be reformulated having the same structure of (P sc ) considered at the beginning of this section:

min x=((x ξ t ) ξ∈Ξ ) t∈G T t=1 ξ∈Ξ p(ξ)g t (x ξ t , ξ) + δ ξ∈Ξ {0m ξ } ( T t=1 B t (x ξ t ) ξ∈Ξ ) + δ W -1 N (x) (3.14)
where

B t :=    B ξ 1 t . . . B ξ E t    .
Since algorithm SMS3BF applied to last problem uses the proximal step of δ W -1 N through the isomorphic mapping W , we get an equivalent expression for the general resolvent of ∂δ W -1 N with respect to matrix

Q := λdiag (Q 1 , • • • ,Q T )
, where Q t is defined as

Q t :=    p(ξ 1 )I nt×nt . . . p(ξ E )I nt×nt    .
Proposition 3.4.1 Given x = ((x ξ t ) ξ∈Ξ ) t∈G ∈ IR nE and M defined as before. Then the resolvent value z := (∂δ W -1 N +Q) -1 Q(x) can be calculated as follows: for every t ∈ G and every A ∈ A t ,

z ξ t = 1 ξ ∈A p(ξ ) ξ ∈A p(ξ )x ξ t for all ξ ∈ A. Proof. Set Y ∈ L and P Y = Proj N Y . Then Y -P Y, η L = 0 for all η ∈ N and hence Q(W -1 Y -W -1 P Y ), W -1 η = 0 for all η ∈ N . Since W -1 P Y ∈ W -1 N , we get (∂δ W -1 N + Q) -1 Q(W -1 Y ) = W -1 P Y.
So, using the equivalent expression of P Y given by Rockafellar and Wets [START_REF] Rockafellar | Scenarios and policy aggregation in optimization under uncertainty[END_REF], we deduce the result.

We now apply algorithm SMS3BF considering for every t ∈ G, the parameter matrices

Q t = Q t and M t = diag (p(ξ 1 ) M ξ 1 t , • • • , p(ξ E ) M ξ E t
), Then we obtain the following algorithm, called "Time Scenarios Decomposition (TSD)" (TSD)

For every t ∈ G, ξ ∈ Ξ, set M ξ t a symmetric positive definite matrix. For every ξ ∈ Ξ, set also Σ ξ = ( T t=1 ( M ξ t ) -1 ) -1 and p = ξ∈Ξ m ξ . Then for an arbitrary (x 0 , y 0 , z 0 , v 0 , u

0 ) ∈ IR nE × IR pT × IR nE × IR p × IR nE , Step 1. find x k+1
For all t ∈ {1, ..., T } and ξ ∈ Ξ do find (x k+1 ) ξ t such that

(x k+1 ) ξ t = g t (•, ξ) + λI nt×nt + (B ξ t ) t M ξ t B ξ t -1 (η) where η = (B ξ t ) t M ξ t B ξ t (y k ) ξ t + (z k ) ξ t -p(ξ) -1 B ξ t (v k ) ξ -p(ξ) -1 (u k ) ξ t end for
Step 2. Calculate y k+1 For all t ∈ {1, ..., T } and ξ ∈ Ξ do

(y k+1 ) ξ t = B ξ t (x k+1 ) ξ t -( M ξ t ) -1 Σ ξ T i=1 B ξ i (x k+1 ) ξ i .
end for

Step 3. Calculate z k+1 For all t ∈ {1, ..., T } do

set A ∈ A t (z k+1 ) ξ t = ξ ∈A p(ξ )[(x k+1 ) ξ t + p(ξ) -1 (u k ) ξ t ] / ξ ∈A p(ξ ), ∀ξ ∈ A.
Step 4.

Calculate vk+1

For all ξ ∈ Ξ do

(v k+1 ) ξ = (v k ) ξ + p(ξ)Σ ξ T i=1 B ξ i (x k+1 ) ξ i .
end for

Step 5. Calculate u k+1 For all t ∈ {1, ..., T } and ξ ∈ Ξ do

(u k+1 ) ξ t = (u k ) ξ t + p(ξ)[(x k+1 ) ξ t -(z k+1 ) ξ t ] end for
Douglas-Rachford map is recovered from (2.53) considering the associated lagrangian map and a special matrix (see Remark 2.4.3). On the other hand, the case C = 0 and under mild assumption, we construct an average map with similar properties as the Davis-Yin map (recovering it when A = I) using the same definition (2.53), but considering a variant of the lagrangian map associated to (V ar) and a special matrix. Then, we construct a generalized resolvent map defined in (2.6) deducing from it new splitting algorithms in order to solve problem (V ar).

The structure of problem (V ar) is related to the variational formulation of the minimization of separable convex functions:

Minimize f (x) + g(Ax) + h(x) (4.1)
where f : IR n → Ī R and g : IR m → Ī R are proper lower semi-continuous convex functions, h : IR n → IR is convex and ( 1β )-Lipschitz-differentiable, and A an m×n matrix.

Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] considers problem (4.1) where under some regularity conditions it can be equivalently written as the following inclusion problem

0 0 ∈ ∂f (x) + ∇h(x) ∂g * (y) + 0 A t -A 0 x y . (V )
To solve such problem, the author applies the Forward-Backward method, getting two algorithms described just below where, for simplicity, we fixed the relaxation parameter ρ > 0 and without error term, termed for us CA1 and CA2:

Algorithm (CA1)        xk+1 = (τ ∂f + I n×n ) -1 (x k -τ ∇h(x k ) -τ A t y k ) ỹk+1 = (σ∂g * + I m×m ) -1 (y k + σA(2x k+1 -x k )) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ) and Algorithm (CA2)        ỹk+1 = (σ∂g * + I m×m ) -1 (y k + σAx k ) xk+1 = (τ ∂f + I n×n ) -1 (x k -τ ∇h(x k ) -τ A t (2ỹ k+1 -y k )) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k )
The main difference between CA1 and CA2 is the parameter matrices chosen in the Forward-Backward method. Choosing special parameter matrices and Lagrangian maps in our general setting, we get different variants of algorithms CA1 and CA2. For example, a variant of CA1 is the recently algorithm YA proposed by O'connor [START_REF] O'connor | On the equivalence of the primal-dual hybrid gradient method and douglas-rachford splitting[END_REF], a slight variant of algorithm PD3O proposed by Yan [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF], in order to solve model (4.1): Algorithm (YA)

x k+1 = (τ ∂f + I n×n ) -1 (x k -τ ∇h(x k ) -τ A t y k ) y k+1 = (σ∂g * + I m×m ) -1 (y k + σA(2x k+1 -x k + τ ∇h(x k ) -τ ∇h(x k+1 )))
Note that CA1 and YA differ on their second update and on the choose of parameter ρ, considering ρ = 1 in the second algorithm. Through numerical experiment, Yan [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF] noticed that YA has more advantages than CA1 (considering ρ = 1).

Similarly, our variant algorithm of CA2 has the same advantage like YA respect to CA1, ie. has large range of acceptable parameters ensuring convergence and better numerical result.

In Chapter 2 we have proposed some splitting algorithms for the following separable optimization problem min

(x,z) [f (x) + g(z) : Ax + Bz = 0] (P 0 )
where f and g are convex and A and B two matrices of appropriated dimensions.

In practice (see for instance [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF]) f and g have the form (k + h) where h is convex differentiable and k not necessarily differentiable but with a tractable proximal step. So, we need to propose an appropriate algorithm such that instead of finding the value of (∂k + ∇h + Q) -1 at a given point, where Q is a symmetric positive definite matrix, the algorithm must use the proximal step of ∂k and the evaluation of ∇h, separately, in oder to obtain a splitting structure. So, we will assume that problem (P 0 ) has the following form: min

(x,z) [f (x) + h 1 (x) + g(z) + h 2 (z) : Ax + Bz = 0]. (P )
where f and g are convex lsc functions, h i (i = 1, 2) is convex and ( 1 β i )-Lipschitzdifferentiable, and A and B two matrices of order m×n and m×p, respectively. It is clear that this problem includes problem (4.1) by considering B = -I p×p and h 2 = 0.

In Section 4.1, we analyze the case A injective (m ≥ n), constructing an average map and a related splitting algorithm for solving (V ar). Then, in order to obtain another algorithm switching the proximal steps of S and T with respect to algorithm (4.4)-(4.7) found in this section, we construct an appropriated average map, getting also ergodic and nonergodic convergence.

The general case (A not necessarily injective) is analyzed in Section 4.2 by reformulating (V ar) as a problem preserving his original structure but corresponding to an injective matrix. Then, applying the results of Section 4.1, we get two general algorithms in this setting.

So, taking special parameter matrices intervening in the two previous general algorithms, we get in Section 4.3 two new algorithms (Alg1) and (Alg2) closely related with Condat's splitting algorithms (CA1) and (CA2). We show that (Alg1) and (Alg2) can also be obtained from Davis-Yin algorithm [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF]. The rate of convergence of the two new algorithms is also analyzed.

In Section 4.4, we study problem (P ) by reformulating it in order to apply the same procedure of Section 4.1 corresponding to injective matrix.

In Section 4.5, we apply the algorithm developed in Section 4.2 to a general multi-block convex optimization problem.

Finally, in Section 4.6 a numerical example is given.

Matrix A injective

Problem (V ar) can be set as a composite monotone inclusions by considering the sum of S and C as a unique map, which, similarly to the one given in [START_REF] Oré | A unified splitting algorithm for composite monotone inclusions[END_REF], is equivalent to the primal problem

Find (x, z) ∈ IR n ×IR m such that 0 0 ∈ S(x) + C(x) T (z) + A t -I N {0} (Ax-z)
where, as usual, N K (a) is the normal cone to set K at point a. This problem is in turn equivalent to the following saddle-point inclusion problem

Find (x, z, ȳ) ∈ (IR n × IR m ) × IR m such that 0 ∈ L(x, z, ȳ) (V L )
where L is the map defined on (

IR n × IR m ) × IR m as L(x, z, y) :=    S(x) + C(x) T (z) 0    +    0 0 A t 0 0 -I -A I 0       x z y    .
From [START_REF] Oré | A unified splitting algorithm for composite monotone inclusions[END_REF] we know that the lagrangian map L allows us to find splitting algorithms alternating between general proximal steps applied to S+C and T , separately. In order to also split the general proximal map of S + C into the general proximal map of S and the evaluation of C, we consider an alternative Lagrangian map which considers S and C defined with different variables.

Take M an arbitrary m × m positive definite matrix. Since A is of full-rank, then A t M A is invertible and hence, from he third row-block of expression L in (V L ), we get x = (A t M A) -1 A t M z. One deduces that the solution set of problem (V L ), denoted by sol (V L ), coincides with the solution set of the following inclusion problem

Find (x, z, ȳ) ∈ (IR n × IR m ) × IR m such that 0 ∈ L(x, z, ȳ) (V L )
where L is the map defined on (

IR n × IR m ) × IR m as L(x, z, y) :=    S(x) + C((A t M A) -1 A t M z) T (z) 0    +    0 0 A t 0 0 -I -A I 0       x z y    .
The maps S and T being defined on different variables, we will get an additional separation of these maps by introducing an appropriate regularization map similarly to the one given in Section 4 of [START_REF] Oré | A unified splitting algorithm for composite monotone inclusions[END_REF].

The average map G L S : an appropriate regularization map

First recall that the generalized resolvent operator associated with a maximal monotone map T and a positive semi-definite linear map P is defined by J T P := (T +P ) -1 P which is not necessarily defined on the whole space. Similarly, the corresponding map G T S := S(T + P ) -1 S t where S satisfies P = S t S, is not necessarily defined on the whole space but has contractive properties under additional assumptions (see Section 2.4 of Chapter 2).

Back to the Lagrangian map L defined above, we take a symmetric positive definite matrix M of order m × m and the matrices

P =    A t M A 0 A t 0 0 0 A 0 M -1    and D = M 1 2 A 0 M -1 2 ,
which satisfy P = D t D. We define the generalized resolvent map

J L P := ( L + P ) -1 P
and then the map

G L D := D( L + D t D) -1 D t .
It is clear that the set of fixed points of J L P and G L D are respectively sol (V L ) and

D(sol (V L )) := {M 1 2 Ax + M -1 2 ȳ : -A t ȳ ∈ S(x) + C(x), ȳ ∈ T (Ax)}.
81 Also, by simple calculations we get

G L D = I -JT M + JS A t M A [2 JT M -I -C • JT M ] (4.2)
where

JT M = M 1 2 (T + M ) -1 M 1 2 , JS A t M A = M 1 2 A(S + A t M A) -1 A t M 1 2 , and C = M 1 2 A(A t M A) -1 C(A t M A) -1 A t M 1 2 .
Notice that C = 0 implies L monotone and thereby the co-coercivity of G L D as shown in expression just before of Proposition 2.4.1. If C = 0, the map G L D is still co-coercive by making an additional condition on matrix M as shows the following proposition :

Proposition 4.1.1 Let S : IR n -→
-→ IR n and T : IR m -→ -→ IR m be maximal monotone maps, C : IR n → IR n a β-co-coercive function with full domain, and A an m × n injective matrix. Assume that

(A t M A) -1 ∈ ]0, 2β[ . Then, G L D is α-average with full domain, where α := 2β 4β-(A t M A) -1 < 1.
Proof. The fullness of the domain of G L D is deduced from the maximality of S and T and the fullness of the domain of C.

From expression (4.2) the α-average of map G L D is deduced from Proposition 2.1 of Davis-Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] taking into account that JT M and JS A t M A are both 1-co-coercive, and

1 (A t M A) -1 C is β-co-coercive.
We note that the particular case A = I and M = γI in M

-1 2 G L D M 1
2 correspond to the average map associated to the sum of three maps defined in [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF].

Remark 4.1.1 Using Proposition 2.4.1 the map G L D can also be deduced from the average map associated to the sum of three maps defined in [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] applying to the following equivalent problem of (V ar)

0 ∈ (M 1 2 AS -1 A t M 1 2 ) -1 (y)+M -1 2 T M -1 2 (y)+M 1 2 A(A t M A) -1 C(A t M A) -1 A t M 1 2 (y).
This equivalent problem is deducted as follows, first we note that problem (V ar) is equal to the composite problem

0 ∈ S(x) + A t M 1 2 M -1 2 T M -1 2 + M 1 2 A(A t M A) -1 C(A t M A) -1 A t M 1 2 M 1 2 A(x)
then the dual problem is equal to

0 ∈ -M 1 2 AS -1 (-A t M 1 2 y)+ M -1 2 T M -1 2 + M 1 2 A(A t M A) -1 C(A t M A) -1 A t M 1 2 -1 (y)
finally taking the dual again considering the last problem as sum of two map problem, we obtained the desirable equivalent problem.

Constructing the splitting algorithm

In connection with the resolvent operator J L P and a real positive parameter ρ, we consider for an arbitrary point w 0 ∈ dom J L P , the sequence {w k } defined by

w k+1 ∈ ρJ L P (w k ) + (1 -ρ)w k . (4.3)
Denoting w k := (x k , z k , y k ) and wk+1 := (x k+1 , zk+1 , ỹk+1 ) = J L P w k , we get from

(4.3) that zk+1 = (T + M ) -1 (y k + M Ax k ) (4.4) ỹk+1 = y k + M Ax k -M zk+1 (4.5)
r k+1 = C((A t M A) -1 A t M zk+1 ) (4.6) xk+1 = (S + A t M A) -1 A t M zk+1 -A t ỹk+1 -r k+1 (4.7) (x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ) (4.8)
and

Dw k+1 = ρG L D ( Dw k ) + (1 -ρ) Dw k . (4.9)
It is worth to mention that { Dw k } is the sequence generated by the fixed point method corresponding to operator ρG

L D +(1-ρ)I which is ρα-average if ρ ∈ ]0, α -1 [ because G L
D is α-average (see Proposition 4.1.1), and hence the sequence { Dw k } converges if sol (V ar) is nonempty.

The next proposition concerns the convergence of w k ; its proof is similar to the given in Section 2.4 of Chapter 2. Proposition 4.1.2 With the same hypothesis given in Proposition 4.1.1, let ρ ∈ ]0, α -1 [ and assume that sol (V ar) is nonempty. Then for an arbitrary (x 0 , z 0 , y 0 ) ∈ (IR n × IR p ) × IR m , the sequence (x k , z k , y k ) defined by the sequential update formulae (4.4) -(4.8), converges to some element of sol (V L ).

Proof. From the comments given just before this proposition, { Dw k } converges to some b ∈ sol (V L ), which is a fixed point to G L D . Now, from (4.3) and considering w := ( L + P ) -1 D t b, we get by using the triangular inequality

w k+1 -w ≤ ρ ( L + P ) -1 D t ( Dw k ) -w + |1 -ρ| w k -w .
On the other hand, the map J L P being continuous, ( L + P ) -1 D t = J T P D + is also continuous, where D + denotes the Moore-Penrose pseudo-inverse matrix of D, and hence the convergence of {w k } to w is deduced. Remark 4.1.2 We note that the sequence z k in (4.4)-(4.8) is only used in the final step of the algorithm, so we can discard it and consider the final step as

(x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ).
Applying (4.4)-(4.8) to the optimization problem (4.1), we get the following sequence:

zk+1 = argmin g(z) + 1 2 z -M -1 y k -Ax k 2 M (4.10) ỹk+1 = y k + M Ax k -M zk+1 (4.11) rk+1 = A(A t M A) -1 ∇h((A t M A) -1 A t M zk+1 ) (4.12) xk+1 = argmin f (x) + 1 2 Ax -zk+1 + M -1 ỹk+1 + rk+1 2 M (4.13) (x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ). (4.14)

Switching the proximal step

Applying the forward-backward method to a lagrangian inclusion problem (V), Condat [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] obtains two splitting algorithms CA1 and CA2 corresponding to two appropriated parameter matrices. The main difference between these algorithms is the order of action of the proximal steps.

In the same manner, we present an algorithm switching the order of action of the proximal steps regarding algorithm (4.4)-(4.8). For this purpose, it is not only necessary to find an appropriate matrix but also consider another alternative Lagrangian map.

For a given m × m positive definite matrix M (similarly to the given at the beginning of this section), we consider the map L defined on (

IR n × IR m ) × IR m as L(x, z, y) :=    S(x) T (z) + M A(A t M A) -1 C(x) 0    +    0 0 A t 0 0 -I -A I 0       x z y   
whose set of their zeroes is exactly

{(x * ,M A(A t M A) -1 C(x * ) + y * , z * ) : (x * , y * , z * ) ∈ sol (V L )}.
Note that when C ≡ 0, this new Lagrangian coincides with Lagrangian L corresponding to problem (V L ) defined at the beginning of this section.

Analogously to the one given in Subsection 4.1.1, we define J L P and G L D corresponding to matrices

P =    0 0 0 0 M -I 0 -I M -1    and D = 0 M 1 2 -M -1 2 , which satisfies P = D t D. It follows that G L D = I -JS A t M A + JT M [2 JS A t M A -I -C • JS A t M A ] (4.15)
where JT M , JS A t M A and C are defined just after of expression (4.2).

We note that the main difference between G L D and G L D (the last one defined in (4.2)) is the switching position of JT M and JS A t M A in their expressions. So, under the same conditions giving in Proposition 4.1.1, we get that G L D is also α-average. The fixed point iteration method applied to J L P generates the following sequences

xk+1 = (S + A t M A) -1 (A t M z k -A t y k ) (4.16)
ỹk+1 = y k + M Ax k+1 -M z k (4.17)

r k+1 = M A(A t M A) -1 C(x k+1 ) (4.18) zk+1 = (T + M ) -1 M Ax k+1 + ỹk+1 -r k+1 (4.19) (x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ). (4.20)
We note that the main difference between this algorithm and the one of (4.4)-(4.8) is the order position of the (general) proximal step of S and T . Remark 4.1.3 Briceño [START_REF] Briceno-Arias | Forward douglas-rachford splitting and forward partial inverse method for solving monotone inclusions[END_REF] has analyzed the particular case A = I and T = N V where V is a linear sub-space of a Hilbert space H. He proposed two alternative methods where the first one was obtained through the composition of two special average maps, and the other one through the forward-backward method applied to the sum problem corresponding to the partial inverse of map S with respect to V and a special co-coercive map. The considered model is

Find x ∈ H such that 0 ∈ m i=1 S i (x) + C(x) (S v )
where for i ∈ {1, • • • , m}, S i is maximal monotone and C co-coercive, all defined on H. The aforementioned algorithms were applied to an appropriate reformulation of this model. Considering H = IR n , an alternative reformulation of problem

(S v ) is Find x ∈ IR n such that 0 ∈ S 1 (x) + D t T (Dx) + C(x) (RS v )
where

T := (S 2 , • • • , S m )

Rate of Convergence

This part is dealing with the rate of convergence of algorithm (4.10)-(4.14) . In this direction, the next proposition gives an upper bound estimation of the saddle-point gap of optimization problem (4.1) defined in the introduction of this chapter.

Proposition 4.1.3 With the same notations as before and considering w = (x, z, y) ∈ dom (f ) × dom (g) × IR m , the following inequality holds:

w k -w 2 P -γ w k+1 -w k 2 P -w k+1 -w 2 P ≥ 2ρ l(x k+1 , zk+1 , y) -l(x, z, ỹk+1 )
where

γ = 1 ρ 2 -ρ -(A t M A) -1
β and l(x, z, y) = f (x) + h(x) + g(z) + y, Ax + Bz .

Proof. Since wk+1 = J L P w k , one has L( wk+1 ) ∈ P (w k -wk+1 ) and hence

L ( wk+1 ) ∈ P (w k -wk+1 ) -    ∇h( z k+1 ) 0 0   
where L is L without the term C, and

z k+1 = (A t M A) -1 A t M zk+1 . Note that L = (∂ x,z l ) × (∂ y [-l ]), where l (x, z, y) = f (x) + g(z) + y, Ax -z .
From Prop. 3 given in [START_REF] Oré | A unified splitting algorithm for composite monotone inclusions[END_REF] and denoting w := (x, z, y) ∈ dom (l), we get wk+1 -w, P (w k -wk+1 ) -xk+1 -x, ∇h( z k+1 ) ≥ l (x k+1 , zk+1 , y)-l (x, z, ỹk+1 ). (

On the other hand, since h is convex and β -1 -Lipschitz differentiable, we have

1 2β z k+1 -xk+1 2 + ∇h( z k+1 ), xk+1 -z k+1 ≥ h(x k+1 ) -h( z k+1 ) and h(x) ≥ h( z k+1 ) + ∇h( z k+1 ), x -z k+1 .
Then, summing the three last inequalities we get wk+1 -w, P (w

k -wk+1 ) + 1 2β z k+1 -xk+1 2 ≥ l(x k+1 , zk+1 , y) -l(x, z , ỹk+1 ). (4.22)
We now find an appropriate upper bound for 1 2β z k+1 -xk+1 2 . From the definition of G L D we have

D wk+1 = G L D ( Dw k ) = Dw k -M 1 2 zk+1 + M 1 2 Ax k+1
and hence by using w k+1 = ρ wk+1 + (1 -ρ)w k we get

1 ρ D(w k+1 -w k ) = D( wk+1 -w k ) = M 1 2 Ax k+1 -M 1 2 zk+1 . (4.23)
Since

z k+1 -xk+1 = (A t M A) -1 A t M 1 2 [M 1 2 zk+1 -M 1 2 Ax k+1 ]
we get, from (4.23), the desirable appropriate upper bound

1 2β z k+1 -xk+1 2 ≤ (A t M A) -1 2β 1 ρ D(w k+1 -w k ) 2 = (A t M A) -1 2βρ 2 w k+1 -w k 2 P .
On the other hand, from the symmetry of matrix P , it holds 2ρ wk+1 -w, P (w k -wk+1 ) = w k -w 2 P -

2 -ρ ρ w k+1 -w k 2 P -w k+1 -w 2 P .
So, replacing the two last expressions into (4.49), we get the desired inequality.

Upper bound of fixed-point residual

Set M and ρ satisfying the hypothesis of Proposition 4.1.2, i.e,

α -1 = 2 - (A t M A) -1 2β > ρ > 0.
With this condition ρG L D + (1 -ρ)I is ρα-average and hence from (4.9), we have

Dw k -Dw * 2 -θ Dw k+1 -Dw k 2 -Dw k+1 -Dw * 2 ≥ 0 (4.24) where θ = 1 ρ 2 -ρ -(A t M A) -1 2β
> 0 and w * ∈ sol (V L ).

Then using similar argument given in Section 2.5.1 of Chapter 2, we get that

Dw k -Dw k-1 2 ≤ 1 kθ Dw 0 -Dw * 2 (4.25) and 1 N N k=1 ( Dw k -Dw k-1 ) 2 ≤ 4 N 2 Dw 0 -Dw * 2 . (4.26)
These two relations can also be deduced respectively from Theorem 1 "Notes on Theorem 1" and Theorem 2 given in [START_REF] Davis | Convergence rate analysis of several splitting schemes[END_REF].

Bounding the saddle-point gap

We consider the following ergodic sequences defined for N ≥ 1 as

x N := 1 N N k=1 xk , z N := 1 N N k=1
zk and

y N := 1 N N k=1
ỹk .

We get the following result

(A t M A) -1 β ≥ ρ > 0.
The following rate of convergence are deduced:

• Ergodic Convergence: for any w = (x, z, y)

∈ dom (f ) × dom (g) × IR m , l(x k , z k , y) -l(x, z, y k ) ≤ 1 2ρk
Dw 0 -Dw 2 . (4.27)

• Nonergodic Convergence: for any w * = (x * , z * , y * ) ∈ sol (V L ),

l(x k+1 , zk+1 , y * ) -l(x * , z * , ỹk+1 ) ≤ α 1 √ k + 1 + α 2 k + 1 Dw 0 -Dw * 2 (4.28)
where

α 1 = 1 + |1 -ρ| ρ 2 √ θ and α 2 = (A t M A) -1 2βρ 2 θ .

Linear convergence

In the case when A = I (sum of three maps problem), D. Davis and Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] got a linear convergence under additional regularity condition (co-coercive, strong montone or lipschitz properties) over the monotone maps. In the case when matrix A is injective, since the map G L D can be obtained from Davis-Yin map associated to a sum of three maps problem (see remark 4.1.1), and noting that if S γ-strong monotone then (M

1 2 AS -1 A t M 1 2 ) -1 is γ A t M A -strong monotone and if T is θ-co- coercive then M -1 2 T M -1 2 is θ M -co-coercive.
We deduce the linear convergence of algorithm (4.4)-(4.8) and (4.16)-(4.20), if we consider the additional hypothesis: S strong monotone and T is co-coercive.

The general case on matrix A

We now consider problem (V ar) without assuming matrix A injective. In order to cope with this rank deficiency, we reformulate the problem as

0 ∈ S(x) + A t M 1 2 V 1 2 M -1 2 T M -1 2 0 M 1 2 Ax V 1 2 x + C(x) (V ar 1 )
where M and V are two symmetric matrices of order m × m and n × n, respectively, with V positive semi-definite and M positive definite.

It is clear that matrix M

1 2 A V 1 2
is injective if and only if A t M A + V is invertible. So, applying the algorithm described in Section 4.1 for matrix M = I we get a splitting algorithm for problem (V ar) in the general setting.

It is important to note that formulation (V ar 1 ) is motivated by the optimization problem defined in (4.1). Indeed, using the same notations given at the end of the first section of Chapter 1, problem (4.1) can be formulated as min

(x,z 1 ,z 2 )∈F f (x) + h(x) + (g, 0)(z 1 , z 2 ) (•)
where F is the set of all triples (x, z 1 , z 2 ) satisfying

M 1 2 A V 1 2 x + -M 1 2 0 0 -I z 1 z 2 = 0.
It is clear that the optimal solution set of problem (•) consists of all (x, Ax, V Notice that problem (•) has the same structure as problem (P 0 ). Indeed, by taking f 1 (x) = f (x) + h(x) and f 2 = (g, 0) and the matrices

B 1 = M 1 2 A V 1 2 and B 2 = -M 1 2 0 0 -I ,
the dual variational formulation of problem (•) consists in finding a zero of the sum of two composite maps consisting of

- M 1 2 A V 1 2 (∂f + ∇h) -1 -A t M 1 2 V 1 2
and

K := - -M 1 2 0 0 -I ∂g 0 -1 - -M 1 2 0 0 -I .
Then the dual problem of this sum problem is

0 ∈ ∂f (x) + ∇h(x) + A t M 1 2 V 1 2 K -1 M 1 2 Ax V 1 2 x .
Since K -1 has the following expression

K -1 = M -1 2 0 0 I ∂g 0 M -1 2 0 0 I = M -1 2 ∂gM -1 2 0 ,
then the last inclusion problem can be set as

0 ∈ ∂f (x) + A t M 1 2 V 1 2 M -1 2 ∂g(x)M -1 2 0 M 1 2 Ax V 1 2 x + ∇h(x).

The main algorithm for non injective operators

Like to map L and matrix D corresponding to problem (V ar), we denote respectively by L and D the map and matrix corresponding to problem (V ar 1 ). It is clear that

L is defined on (IR n × IR m × IR n ) × (IR m × IR n ) and its value is L (x, z 1 , z 2 , y 1 , y 2 ) =        S(x) + C((A t M A + V ) -1 (A t M 1 2 z 1 + V 1 2 z 2 )) M -1 2 T M -1 2 (z 1 ) 0 z 2 0 y 1 0 y 2        +        0 0 0 A t M 1 2 V 1 2 0 0 0 -I 0 0 0 0 0 -I -M 1 2 A I 0 0 0 -V 1 2 0 I 0 0               x z 1 z 2 y 1 y 2        .
The set of zeroes of L is

{(x, M 1 2 z, V 1 2 x, M -1 2 ȳ, 0 y 2 ) : (0 x , 0 z , 0 y ) ∈ L(x, z, ȳ)},
where 0 w denotes the zero vector of the w-space.

On the other hand, the corresponding map defined in (4.2) associated with L is denoted by G L D which applies IR m × IR n into itself and whose value at (u, x) is

  u -JT u + M 1 2 A JS V 1 2 x + A t M 1 2 (2 JT u -u) -C(V 1 2 x + A t M 1 2 JT u) V 1 2 JS V 1 2 x + A t M 1 2 (2 JT u -u) -C(V 1 2 x + A t M 1 2 JT u)   where JT = M 1 2 (T + M ) -1 M 1 2 , JS = (S + V + A t M A) -1 and C = C(A t M A + V ) -1 . Set r = (V + A t M A) -1 . From Proposition 4.1.1, G L D is 2β 4β-r -average with full domain, if r ≤ 2β.
Applying the algorithm described in (4.4)-(4.8), we get that

zk+1 1 = M 1 2 (T + M ) -1 M 1 2 (y k 1 + M 1 2 Ax k ) zk+1 2 = y k 2 + V 1 2 x k ỹk+1 1 = y k 1 + M 1 2 Ax k -zk+1 1 ỹk+1 2 = y k 2 + V 1 2 x k -zk+1 2 r k+1 = C((A t M A + V ) -1 (A t M 1 2 zk+1 1 + V 1 2 zk+1 2 )) xk+1 = (S + A t M A + V ) -1 A t M 1 2 zk+1 1 + V 1 2 zk+1 2 -A t M 1 2 ỹk+1 1 -V 1 2 ỹk+1 2 -r k+1 (x k+1 , z k+1 , y k+1 , r k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 , rk+1 ) + (1 -ρ)(x k , z k , y k , r k ).
Defining the new variables η k = σAx k + y k -σz k -σAr k and ηk = σAx k + ỹkσz k -σAr k , we obtain the following algorithm Algorithm (Alg1)

             xk+1 = (τ S + I n×n ) -1 (x k -τ A t η k -r k ) rk+1 = τ C(x k+1 ) ηk+1 = (σT -1 + I m×m ) -1 (η k + σA(2x k+1 -x k ) + σAr k -σAr k+1 ) (x k+1 , η k+1 , r k+1 ) = ρ(x k+1 , ηk+1 , rk+1 ) + (1 -ρ)(x k , η k , r k )
This algorithm is indeed a variant of CA1 for solving optimization problem (4.1).

In particular when ρ = 1, we recover YA algorithm [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF] which is described in the introduction of this chapter.

On the other hand, applying (4.30)-(4.34) with M = σI m×m and V = τ -1 I n×n -σA t A, we get the following sequence:

zk+1 = (T + σI m×m ) -1 (y k + σAx k ) ỹk+1 = y k + σAx k -σz k+1 r k+1 = C((τ -1 I n×n -σA t A)τ x k + τ σA t zk+1 )) xk+1 = (S + τ -1 I n×n ) -1 (τ -1 I n×n -σA t A)x k + σA t zk+1 -A t ỹk+1 -r k+1 (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ).
Then using that (σT -1 + I m×m ) -1 = I -σ(T + σI m×m ) -1 , and eliminating the term zk+1 , we obtain the following algorithm:

Algorithm (Alg2)              ỹk+1 = (σT -1 + I m×m ) -1 (y k + σAx k ) r k+1 = τ C(x k -τ A t (ỹ k+1 -y k )) xk+1 = (τ S + I n×n ) -1 (x k -τ A t (2ỹ k+1 -y k ) -r k+1 ) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ).
So, by considering σ, τ and ρ positive parameters such that στ A 2 ≤ 1, τ < 2β and ρ < 4β-τ 2β , then applying Proposition 4.2.1, we deduce the convergence of sequence (x k , Ax k , y k ) to an optimal solution of the lagrangian problem corresponding to problem (4.1).

With respect to the convex problem (4.1), algorithm Alg2 is a variant of CA2 by changing ∇h(x k ) by ∇h(x k -τ A t (ỹ k+1 -y k )).

Relationship with the Condat's method

When C = 0, algorithms Alg1 and Alg2 are exactly CA1 and CA2 respectively. Otherwise, when C = 0, they are different. We consider the lagrangian function defined as l (x, y) = f (x) -g * (y) + Ax, y

and its corresponding maximal monotone map L defined on IR n × IR m as

L (x, y) := (∂ x l ) × (∂ y [-l ]) = ∂f (x) ∂g * (y) + 0 A t -A 0 x y . (4.43)
The next inequalities which are immediately deduced by definition will be used later in Proposition 4.3.1. Notice that algorithms CA1 and CA2 generate the sequences w k i = (x k , y k ) and wk i = (x k , ỹk ), for i = 1, 2, respectively, which satisfy the following inclusion

L ( wk+1 i ) + P i ( wk+1 i -w k i ) - ∇h(x k ) 0
where P 1 and P 2 are matrices defined as

P 1 = 1 τ I -A t -A 1 σ I and P 2 = 1 τ I A t A 1 σ I . (4.44)
As showed in the proof of Theorem 2 regarding the relaxed primal-dual algorithm given by Chambolle-Pock [START_REF] Chambolle | A first-order primal-dual algorithm for convex programs with applications to imaging[END_REF], it holds that

2ρ[L(x k+1 , y) -L(x, ỹk+1 )] ≤ w k -w 2 P i -wk+1 -w 2 P i - 2 -ρ ρ w k+1 -w k 2 U i
where U 1 and U 2 are two matrices defined as

U 1 = ( 1 τ -1 β(2-ρ) )I -A t -A 1 σ I and U 2 = ( 1 τ -1 β(2-ρ) )I A t A 1 σ I ,
and L defined as

L(x, y) = f (x) + h(x) -g * (y) + Ax, y . (4.45)
The last inequality is fundamental to deduce the ergodic convergence of CA1 and CA2 as showed in the aforementioned Chambolle-Pock's paper.

On the other hand, algorithm Alg1 generates the sequences νk = (x k , ηk ), ξ k = (x k -r k , η k ) and ξk = (x k -rk , ηk ) satisfying

L (ν k+1 ) + P 1 ( ξk+1 -ξ k ) - ∇h(x k+1 ) 0 . (4.46)
Similarly, algorithm Alg2 generates the sequences w k = (x k , y k ) and wk+1 = (x k+1 , ỹk+1 ) satisfying

L ( wk+1 ) + P 2 ( wk+1 -w k ) - ∇h(x k -τ A t (ỹ k+1 -y k )) 0 . (4.47)
From (4.47) and (4.46), we deduce the following fundamental upper bound of the saddle-point gap. This will be used later in the analysis of the ergodic and nonergodic rates of convergence. Proposition 4.3.1 With the same notations as before; let us assume that στ A 2 ≤ 1, then for any w = (x, y) ∈ dom (f ) × dom (g * ), the following inequalities hold:

• For the sequences generated by algorithm Alg1, it holds that

ξ k -w 2 P 1 -λ ξ k+1 -ξ k 2 P 1 -ξ k+1 -w 2 P 1 ≥ 2ρ[L(x k+1 , y) -L(x, ηk+1 )]
where w = w -(τ ∇h(x), 0)

• For the sequences generated by algorithm Alg2, it holds that

w k -w 2 P 2 -λ w k+1 -w k 2 P 2 -w k+1 -w 2 P 2 ≥ 2ρ[L(x k+1 , y) -L(x, ỹk+1 )] where λ = 1 ρ 2 -ρ -τ β .
Proof. From (4.47), and applying Proposition 4.3.1, we get νk+1 -w, P 1 (ξ k -ξk+1 ) -xk+1 -x, ∇h(x k+1 ) ≥ l (x k+1 , y) -l (x, ηk+1 ).

using that νk = ξk + (r k , 0) in the last inequality, we have ξk+1 -w, P 1 (ξ k -ξk+1 ) -u k+1 -x, ∇h(x k+1 ) ≥ l (x k+1 , y) -l (x, ηk+1 ).

where u k+1 = rk+1 + x k -r k + τ A t (η k+1 -η k ). Then, since h is convex and β -1 -Lipschitz-differentiable, we have

∇h(x k+1 ), xk+1 -x - β 2 ∇h(x k+1 ) -∇h(x) 2 ≥ h(x k+1 ) -h(x)
and from the properties of norm we have

1 2β xk+1 -u k+1 2 ≥ xk+1 -u k+1 , ∇h(x k+1 ) -∇h(x) - β 2 ∇h(x k+1 ) -∇h(x) 2 .
Then, summing the three last inequalities we get ξk+1 -w , P 1 (ξ k -ξk+1 ) + 1 2β xk+1 -u k+1 2 ≥ L(x k+1 , y) -L(x, ηk+1 ). (4.48) where w = w -(τ ∇h(x), 0). Now we find an appropriate upper bound for 1 2β xk+1 -u k+1 2 . From the expression of P 1 we have (x, y) 2

P 1 = τ -1 x -τ A t y 2 + y 2 σ -1 I-τ AA t
Then since στ A 2 ≤ 1, we have that σ -1 I -τ AA t is positive definite matrix, so we get

xk+1 -u k+1 2 = xk+1 -rk+1 -x k + r k -τ A t (η k+1 -η k ) 2 ≤ τ ξk+1 -ξ k 2 P 1 = τ ρ 2 ξ k+1 -ξ k 2 P 1
On other hand, from the symmetry of matrix P 1 , it holds

2ρ ξk+1 -w , P 1 (ξ k -ξk+1 ) = ξ k -w 2 P 1 - 2 -ρ ρ ξ k+1 -ξ k 2 P 1 -ξ k+1 -w 2 P 1 .
So, replacing the two last expressions into (4.48), we get the desired inequality of the first item. Now we proof of second item, from (4.46), and applying Proposition 4.3.1, we get wk+1 -w, P 2 (w k -wk+1 ) -xk+1 -x, ∇h( z k+1 ) ≥ l (x k+1 , y) -l (x, ỹk+1 ).

where z k+1 = x k -τ A t (ỹ k+1 -y k ). Then since h is convex diferentiable with ∇h β -1 -Lipschitz, we have

1 2β z k+1 -xk+1 2 + ∇h( z k+1 ), xk+1 -z k+1 ≥ h(x k+1 ) -h( z k+1 ) and h(x) ≥ h( z k+1 ) + ∇h( z k+1 ), x -z k+1 .
Then, summing the three last inequalities we get wk+1 -w, P 2 (w k -wk+1 ) + 1 2β z k+1 -xk+1 2 ≥ L(x k+1 , y) -L(x, ỹk+1 ). (4.49)

Then using the same techniques as the first item, we get the desired inequality of the second item.

Relationship with the Davis-Yin's method

Recently D. O'Connor and L. Vandenberghe [START_REF] O'connor | On the equivalence of the primal-dual hybrid gradient method and douglas-rachford splitting[END_REF], noticed that algorithm YA (described early in the introduction of this chapter) can be deduced from the algorithm developed by D. Davis and W. Yin [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] by means a reformulation of the sum of three special operators. Now we show that algorithm Alg1 and Alg2 are also obtained from Davis-Yin's algorithm considering these adhoc reformulations.

Fixing A = I in problem (V ar) and applying algorithms Alg1 and Alg2 with σ = τ -1 , we obtain ỹk+1 = (τ -1 T -1 + I m×m ) -1 (y k + τ -1 x k ) (4.50)

r k+1 = τ C(x k -τ (ỹ k+1 -y k )) (4.51) xk+1 = (τ S + I n×n ) -1 (x k -τ (2ỹ k+1 -y k ) -r k+1 ) (4.52) (x k+1 , y k+1 ) = ρ(x k+1 , ỹk+1 ) + (1 -ρ)(x k , y k ). (4.53)
and the switching algorithm

xk+1 = (τ S + I n×n ) -1 (x k -τ η k -r k ) (4.54) rk+1 = τ C(x k+1 ) (4.55) ηk+1 = (τ -1 T -1 + I m×m ) -1 (η k + τ -1 (2x k+1 -x k ) + τ -1 r k -τ -1 rk+1 )(4.56) (x k+1 , η k+1 , r k+1 ) = ρ(x k+1 , ηk+1 , rk+1 ) + (1 -ρ)(x k , η k , r k ). (4.57)
Notice that [START_REF] O'connor | On the equivalence of the primal-dual hybrid gradient method and douglas-rachford splitting[END_REF] these two algorithms can also be obtained directly from the Davis-Yin's algorithm. The following rate of convergence are deduced:

• Ergodic Convergence: for any (x, y) ∈ dom (f ) × dom (g * ) and i = 1, 2,

L(ζ k i , y) -L(x, ν k i ) ≤ 1 2ρk (ω 0 i , ν 0 i ) -(x -µ i , y) 2 P i (4.58)
where µ i = (2 -i)τ ∇h(x).

• Nonergodic Convergence: for any (x * , y * ) ∈ sol (V ) and i = 1, 2,

L( ζk+1 i , y) -L(x, νk+1 i ) ≤ α 1 √ k + 1 + α 2 k + 1 (ω 0 i , ν 0 i ) -(x * -µ * i , y * ) 2 P i (4.59) where µ * i = (2 -i)τ ∇h(x * ), α 1 = 1 + |1 -ρ| ρ 2 √ θ and α 2 = (A t M A) -1 2βρ 2 θ .
Respect to the rate of constraint violations we have For any (x * , y * ) ∈ sol (V ) and i = 1, 2, the following rate of convergence are obtained by setting u * i = (x * -(2 -i)τ C(x * ), y * ):

• Ergodic Convergence:

ωk i -ω k i 2 V + σ Aω k i -A ω k i + (-1) i σ -1 (ν k i -ν k i ) 2 ≤ 4 k 2 ρ 2 (ω 0 i , ν 0 i ) -u * i 2 P i .
• Nonergodic Convergence:

ωk i -ω k-1 i 2 V +σ Aω k i -Aω k-1 i +(-1) i σ -1 (ν k i -ν k-1 i ) 2 ≤ 1 αk (ω 0 i , ν 0 i )-u * i 2 P i .
where α = ρ(2 -ρ -τ 2β ) and Remark 4.3.1 Considering the sequence (x k , z k , s k ) generated by PD3O [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF], the sequence (x k , y k ) = (x k , s k+1 ) is generated by YA (or equivalently by Alg1 with ρ = 1) using τ = γ and σ = δ. Moreover, the sequence r k in YA verifies z k+1 = x k -r k -τ A t y k . One deduces that

ξ k -(x -τ ∇h(x), s) 2 P 1 = τ -1 x -∇h(x) -τ A t s -z k+1 2 + s -s k+1 2 δ -1 I-γAA t = τ -1 (z, s) -(z k+1 , s k+1 ) 2 I,M
where (a, b) 2 I,M = a 2 + b 2 M with M = γ δ (I -γδAA t ) defined in [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF]. The last equality relations and the upper bound given in Proposiion 4.3.1, allow us recover the upper bound [START_REF] Minty | Monotone operators in hilbert spaces[END_REF] given in Theorem 2 of [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF] related to the aforementioned sequence

(x k , z k , s k ).
Notice that Theorem 2 of [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF] is exactly the ergodic convergence of Theorem 4.3.1.

General separable optimization problem

Following the same scheme described in Section 4.2 regarding the optimization problem (4.1), we reformulate problem (P ) keeping the same structure of problem (V ar), where the involved matrix is injective and then we apply the algorithm developed in Section 4.1.

Set M, V 1 and V 2 symmetric matrices of order m×m, n×n and p×p, respectively, with V 1 positive semi-definite and M and V 2 positive definite. Then problem (P ) can be formulated as min

(x 1 ,x 2 ,z 1 ,z 2 )∈F (f + h 1 , h 2 )(x 1 , x 2 ) + (g, 0)(z 1 , z 2 ),
where F denotes the set of all (x 1 , x 2 , z 1 , z 2 ) satisfying

    M 1 2 A 0 V 1 2 1 0 0 V 1 2 2     x 1 x 2 +    M 1 2 B 0 0 -I n -V 1 2 2 0    z 1 z 2 = 0.
The dual problem of its variational formulation consists in finding a zero of the sum of the two composition maps given by

-     M 1 2 A 0 V 1 2 1 0 0 V 1 2 2     ∂f + ∇h 1 ∇h 2 -1 - A t M 1 2 V 1 2 1 0 0 0 V 1 2 2 and -    M 1 2 B 0 0 -I n -V 1 2 2 0    ∂g 0 -1 - B t M 1 2 0 -V 1 2 2 0 -I n 0 .
Then the dual formulation of this sum problem is 0 ∈ ∂f 0

x 1 x 2 + ∇h 1 ∇h 2 x 1 x 2 + A t M 1 2 V 1 2 1 0 0 0 V 1 2 2 G     M 1 2 A 0 V 1 2 1 0 0 V 1 2 2     x 1 x 2 , ( 4
.60) where the map G is the inverse of the composite map defined by

   -M 1 2 B 0 0 I n V 1 2 2 0    ∂g 0 -1 -B t M 1 2 0 V 1 2 2 0 I n 0
which is clearly monotone.

The next proposition gives an explicit expression of the resolvent of G and thereby its maximal monotonicity by Minty's theorem.

(G + I) -1 (x, y, z) = (-M 1 2 Bu, y, V 1 2 2 u),
where I denotes the identity map of order m + n + p, and

u = (∂g + V 2 + B t M B) -1 (-B t M 1 2 x + V 1 2 2 z). Proof. Since (G + I) -1 = I -(G -1 + I) -1 , then using Proposition 2.4.1 of Chapter 2 we obtain (G + I) -1 =    -M 1 2 B 0 0 I n V 1 2 2 0    ∂g + B t M B + V 2 I n -1 -B t M 1 2 0 V 1 2 2 0 I n 0
from which the desired equality is deduced.

Observe that problems (4.60) and (V ar) have same structure where the involved matrix in the first one is injective and whose corresponding maps verify the properties required in Proposition 4.1.1. So, the corresponding algorithm described by the sequential update formulae (4.4) -(4.8) converges to a solution of its corresponding saddle-point problem, whose solution set is

{(x, z, M 1 2 Ax, V 1 2 1 x, V -1 2 2 z, M -1 2 ȳ, 0 y 2 , -V -1 2
2 ∇h 2 (z)) : (0 x , 0 z , 0 y ) ∈ L(x, z, ȳ)}, where L is the classical Lagrangian map corresponding to problem (P ), which is defined as

L(x, z, y) :=    ∂f (x) + ∇h 1 (x) ∂g(z) + ∇h 2 (z) 0    +    0 0 A t 0 0 B t -A -B 0       x z y    .
Proposition 4.4.2 With the same notations as before, we assume that V 1 ∈ IR n×n , V 2 ∈ IR p×p and M ∈ IR m×m are symmetric with V 2 and M positive definite such that

V 1 + A t M A is positive definite and satisfying that max( V -1 2 , (V 1 + A t M A) -1 ) ∈ ]0, 2 min(β 1 , β 2 )[. If sol (V L)
is nonempty, then for arbitrary points (x 0 , z 0 , y 0 ) ∈ (IR n × IR p ) × IR m and r 0 2 ∈ IR p , the sequence (x k , z k , y k ) generated by (4.61) -(4.65) converges to some element of sol (V L).

Application to the decomposition of multiblock optimization problems

In this section we extend the algorithm (PMA) described in Chapter 2 in order to solve the more general S-Model defined below. This extension uses a similar reformulation as described in the mentioned chapter and has similar structure as problem (V ar) described in Section 4.2.

Our S-Model problem is as follows

inf x=(x 1 ,••• , xq), z q i=1 f i (x i ) + h(x) + g(z) s.t q i=1 A i x i -Bz = 0,
where f i : IR n i → IR (i ∈ {1, ..., q}) and g : IR m → IR are proper lsc convex functions, h : IR n → IR is convex and ( 1 β )-Lipschitz-differentiable (n = q i=1 n i ), and A i and B are matrices of order p × n i and p × m, respectively.

It is clear that this problem is equivalent to inf

(x=(x 1 ,••• , xq), z) q i=1 f i (x i ) + g(z) + δ {0} Bz - q i=1 A i x i + h(x) (B 1 )
or again

V P = inf x=(x 1 ,••• , xq) q i=1 f i (x i ) + (g * • B t ) * q i=1 A i x i + h(x). (B 2 )
We note that (B 2 ) has the same structure as the optimization problem (4.1) given in the introduction of this chapter. If we apply Alg1 or Alg2, we obtain an algorithm with separable structure but unfortunately with two important disadvantages: the necessity to know the norm of A in order to choose the parameters for the convergence result, and also the necessity to know all values of (τ ∂f i + I) -1 at arbitrary points (parameter τ beings equal for all i).

So, we consider a reformulation of (B 2 ) allowing us to choose in an independently manner parameters for each block, where we only need to know the norm of each A i separatly.

By considering f : IR n → IR defined as f (x) = q i=1 f i (x i ), and matrices K and A of order p × pq and pq × n, respectively, defined by

K := I p×p • • • I p×p and A :=    A 1 . . . A q    , problem (B 2 ) can be set as inf x=(x 1 ,••• , xq) f (x) + [(g * • B t ) * • K](Ax) + h(x). (P S-mod )
Note that this new problem has also the same structure as problem (4.1) but with a good separable structure since f and A have separable structure for blocks. Then, we apply to this last problem algorithm GSA3O developed in Section 4.2. Alternatively, we can apply algorithm (4.36) -(4.40), but for simplicity we will not do it in this work.

Regarding function g, we assume that

∂[(g * • B t ) * • K] = K t (B(∂g) -1 B t ) -1 K.
The saddle-point problem of (P S-mod ) is

Find (x, z, ȳ) ∈ IR n × IR pq × IR pq such that 0 ∈ L(x, z, ȳ) (VL)
where L is the maximal monotone map defined on

IR n × IR pq × IR pq as L(x, z, y) :=    ∂f (x) + ∇h(x) K t (B(∂g) -1 B t ) -1 Kz 0    +    0 0 A t 0 0 -I -A I 0       x z y    .
For i ∈ {1, • • • , q}, let M i be an p × p symmetric positive definite matrix and Q i be an n i × n i symmetric positive semi-definite matrix.

In order to take advantage of the separability of f , we consider the diagonal matrices

V = diag (Q 1 , • • • , Q q ) and M = diag (M 1 , • • • , M q ).
So, the related algorithm GSA3O has the following structure:

zk+1 = (K t (B(∂g) -1 B t ) -1 K + M ) -1 (M Ax k + y k ) (4.66) ỹk+1 = y k + M (Ax k -zk+1 ) (4.67) r k+1 = ∇h((A t M A + V ) -1 (V x k + A t M zk+1 )) (4.68) xk+1 = (∂f + A t M A + V ) -1 (V x k + A t M zk+1 -A t ỹk+1 -r k+1 ) (4.69) (x k+1 , z k+1 , y k+1 ) = ρ(x k+1 , zk+1 , ỹk+1 ) + (1 -ρ)(x k , z k , y k ).
(4.70)

We finally get the following algorithm Separable Primal-Dual Variant (SPDV)

For i ∈ {1, • • • , q} set Q i ∈ IR n i ×n i symmetric positive semi-definite, M i ∈ IR p×p symmetric positive definite. Set Σ = q i=1 M -1 i -1 . Then for an arbitrary (x 0 , z 0 , y 0 c ) ∈ IR n × IR pq × IR p
Step 1. Find ζ k+1 such that

ζ k+1 = argmin g(w) + 1 2 Bw - q j=1 (A j x k j ) -Σ -1 y k c 2 Σ
.

Step 2. Find zk+1 For all i ∈ {1, ..., q} do Find zk+1

i such that

zk+1 i = A i x k i -M -1 i Σ q j=1 (A j x k j ) -Bζ k+1 .
Step 3. ) such that

r k+1 = ∇h((A t M A + Q) -1 (Qx k + A t M zk+1 ))
where

A = diag([A 1 , ..., A q ]), Q = diag([Q 1 , ..., Q q ]) and M = diag([M 1 , ..., M q ])
Step 5. Find xk+1 For all i ∈ {1, ..., q} do Find xk+1 i such that

xk+1 i = argmin f i (x i ) + 1 2 A i x i -zk+1 i + M -1 i ỹk+1 c + A i rk+1 i 2 M i + 1 2 x i -x k i + rk+1 i 2 Q i where rk+1 i = (A t i M i A i + Q i ) -1 r k+1 i . end for Step 6. Find (x k+1 , z k+1 , y k+1 c ) (x k+1 , z k+1 , y k+1 c ) = ρ(x k+1 , zk+1 , ỹk+1 c ) + (1 -ρ)(x k , z k , y k c ).
The following converge result is deduced.

Proposition 4.5.1 Assume that Q i and M i are symmetric matrices of order n i ×n i and p × p, respectively, with Q i positive semi-definite and M i positive definite such that Q i +A t i M i A i is positive definite and satisfying

(Q i +A t i M i A i ) -1 ∈ ]0, 2β[ . Let ρ ∈ ]0, α -1 [ where α := 2β 4β-max i { (V i +A t i M i A i ) -1 } q i=1
. If sol (VL) is nonempty, then for an arbitrary (x 0 , z 0 , y 0 c ) ∈ IR n × IR pq × IR p , the sequence (x k , z k , K t y k c ) generated by (SPDV) converges to some element of sol (VL).

M i = σ i I p×p and Q i = τ -1 i I n i ×n i -σ i A t i A i ,
the subproblem in Step 5 of the Algorithm (PMA) becomes

xk+1 i = argmin f i (x i ) + 1 2τ i x i -x k i -τ i [σ i A t i zk+1 i -σ i A t i A i x k i -A t i ȳk+1 c -r k+1 ] 2 .
If in addition the positive parameters σ i and τ i are chosen satisfying σ i τ i A i 2 ≤ 1 and τ i < 2β, then the conditions on matrices Q i , M i and Q i +A t i M i A i in Proposition 4.5.1 are immediately verified and thereby the sequence (x k , z k , K t y k c ) generated by (SPDV) converges to some element of sol (VL) if nonempty.

Numerical Example

We consider the problem (commonly referred as fused lasso) with the least squares loss as in [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF] min

x 1 2 Qx -b 2 2 + µ 1 x 1 + µ 2 Ax 1 (4.71) where Q ∈ IR n×p , b ∈ IR n and A =      -1 1 -1 1 • • • • • • -1 1      ∈ IR (p-1)×p .
We consider n = 200 , p = 4000, µ 1 = 20 and µ 2 = 200. Moreover, matrix Q and vector b quoted from [START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF].

We compare the algorithms CA1, CA2 with their variants algorithms Alg1 and Alg2 considering ρ = 1 as relaxing parameter. Notice that in this case Alg1 coincides with Y A. We choose as parameters

τ = γ 1 A 2 and σ = 1 8τ
and we implement the last algorithms considering three variant for γ as γ 1 = 1, γ 2 = 1.5 and γ 3 = 1.99. For CA1 and CA2 we only use γ 1 because for γ 2 and γ 3 the convergence is not guaranteed.

We observe that all algorithms have the same behavior for γ 1 , but for γ 2 and γ 3 algorithms Y A and Alg2 (= New) have more velocity and maintain the same behavior. 111 the local demands (d zτ ) are random processes and we consider additional terms i zτ that are random input of the storage. The distribution of these random processes is supposed to be known and generally approximated by a finite set of historical scenarios.

Since we are working with random variables affecting the dynamic equations defining the state and control decisions, we need to consider the nonanticipativity constraints that rule the sequence of decisions when the successive realizations of the random values are revealed at each stage. More specifically, given the following constraints :

• The demand equation is given by

p zτ + u zτ + e∈z + f eτ - e∈z - f eτ = d zτ
where z + (resp. z -) is the set of outgoing (res. ingoing) arcs incident to zone z.

• The storage x zτ obeys the dynamics x z,τ +1 = x zτ -u zτ + i zτ where i zτ is a random input in the storage.

• The quantity transported through line f eτ satisfies the capacity constraints as

f eτ ∈ F eτ
• The variables x zt , u zt and p zt satisfy some constraints x zτ ∈ X zτ , u zτ ∈ U zτ , p zτ ∈ P zτ • The control variable f eτ , u zτ and p zτ should then satisfy the nonanticipativity equations p zτ , u zτ , f eτ F τ .

i.e. p zτ , u zτ , f eτ are F τ -measurable where F τ is a σ-field defined as

F τ = σ({(d zτ , i zτ ) : z ∈ Z, τ ∈ [0, t]})
Summarizing, we consider the following multistage stochastic problem

p zτ + u zτ + e∈z + f eτ - e∈z - f eτ = d zτ , ∀z ∈ Z, τ ∈ [0, T -1] (5.2) x z,τ +1 = x zτ -u zτ + i zτ , ∀z ∈ Z, τ ∈ [0, T -1] (5.3) x zτ ∈ X zτ , u zτ ∈ U zτ , p zτ ∈ P zτ ∀z ∈ Z, τ ∈ [0, T -1] (5.4) f eτ ∈ F eτ ∀e ∈ E, τ ∈ [0, T -1] (5.5) x z,0 = xz0 , ∀z ∈ Z (5.6) p zτ , u zτ F τ ∀z ∈ Z, τ ∈ [0, T -1] (5.7) f eτ F τ ∀e ∈ E, τ ∈ [0, T -1] (5.8)

Solution of a deterministic formulation

In the first part of this chapter, we study the deterministic case of the last model. The modeling of cost functions c z and l eτ is borrowed from [START_REF] Dallagi | A stochastic zonal decomposition algorithm for network energy management problems[END_REF].

Cost on the final state The hydroelectric production cost is negligible in the considered model. On the other hand, we add a cost on the final state x zT → Ψ(x zT ) to penalize the excess of water reserves, defined as

Ψ(x zT ) = c f in z max{0, x z0 -x zT }.

Thermic production

The thermic production cost is a piecewise affine and convex function of the production levels p zτ . It will be defined with a given number of stages j = 1, . . . , Q zτ , each one associated with a given slope c j zτ valid in the interval [P j zτ , P j+1 zτ ]. We need obviously 0 ≤ c 1 zτ ≤ c 2 zτ ≤ • • • ≤ c Qzτ zτ to obtain an increasing convex function of p zτ . The cost function is thus defined by

g zτ (p zτ ) =                c 1 zτ p zτ if 0 ≤ p zτ ≤ P 1 zτ c 1 zτ P 1 zτ + c 2 zτ (p zτ -P 1 zτ ) if P 1 zτ ≤ p zτ ≤ P 2 zτ . . . . . . Qzτ -1 j=1 c j zτ P j zτ + c Qzτ zτ (p zτ -P Qzτ -1 zτ ) if P Qzτ -1 zτ ≤ p zτ ≤ P Qzτ zτ +∞ else
where we defined P j zτ = P j zτ -P j-1 zτ with P -1 zτ = 0.

The cost function g zτ is reformulated as

g zτ (p zτ ) =      inf θ Qzτ j=1 c j zτ θ j zτ s.t.
Qzτ j=1 θ j zτ = p zτ 0 ≤ θ j zτ ≤ P j zτ where θ j zτ represents the production of stage j.

Interzonal transfer costs For an arc e = (z, z ) ∈ E interconnecting two zones z and z , the flow transfer during period τ is the variable f eτ which is bounded by 0 ≤ f eτ ≤ κ eτ . The transfer cost is linear and given by l eτ (f eτ ) = c inter e f eτ .

Failure cost

The failure quantity corresponds to the part of demand not satisfied during period τ in zone z. It will be denoted by η zτ . It is penalized by the cost c f ail >> c j zτ .

Using former formulations of the cost functions, the problem becomes min We rewrite this problem in the context of problem (P ) defined in Chapter 2, in order to apply GSS algorithm also described in that chapter.

(θ,u,f,η,v) T -1 τ =0 z∈Z Qzτ j=1 c j zτ θ j zτ + c f ail η zτ + e∈E l eτ f eτ + z∈Z c f in z v z (5.9) satisfying u zτ + Qzτ j=1 θ j zτ - e∈z + f eτ + η zτ = d zτ - e∈z - f eτ , ∀z ∈ Z, τ ∈ [0, T -1](5.10) x z,τ +1 = x zτ -u zτ + i zτ , ∀z ∈ Z, τ ∈ [0, T -1] (5.11) X min zτ ≤ x zτ ≤ X max zτ , ∀z ∈ Z, τ ∈ [0, T -1] (5.12) 0 ≤ u zτ ≤ U max
Setting q zτ = ((θ j zτ ) j∈Qzτ , u zτ , η zτ , x z,τ +1 ) ∈ IR Qzτ +3 and f τ = (f eτ ) e∈E , relation (5.10) becomes

A zτ q zτ -B z f τ = d zτ , ∀z ∈ Z, ∀τ ∈ [0, T -1]
where B z is the row z of the incidence node-arc matrix for graph G, and the matrix A zτ defined as A zτ = 1 1×Qzτ 1 1 0

Considering q z = (q zτ ) τ ∈[0,T -1] , the objective cost function (5.9) can be rewritten as

z∈Z k z (q z , v z ) + e∈E T -1 τ =0 l eτ f eτ ,
where the zonal cost of production K z is equal to

k z (q z , v z ) = T -1 τ =0 Qzτ j=1 c j zτ θ j zτ + c f ail η zτ + c f in z v z .
We introduce the set of zonal constraints associated to the production C z which are the constraints (5.11)-(5.16) except (5.15) and the set of interzonal transfer constraint C eτ which is the constraint (5.15). Then the planning problem (5.9)-(5.16) now reads min (q, v, f ) z∈Z k z (q z , v z ) + In summary, we rewrite (5.17 Since this problem has the same structure as problem (P ), we apply GSS considering different parameters. In fact, we apply slight modifications of ADMM, Chambolle-Pock, and Spingarn algorithms, which correspond to the GSS algorithms with different parametrized matrices. Notice that each resultant sub-problem is a classical proximal step and can be solvable in parallel procesing with respect to z and eτ indices respectively. Now, in order to get a more explicit form the last algorithm, we find explicit expressions of the norms of A t and Bt (for this purpose we assume in the original model (5.9)-(5.16) that the graph is complete), and also explicit expressions of A t (Aw k + Bf k -d + λ -1 y k ) and B t (Aw k+1 + Bf k -d + λ -1 y k ). Since that in the original model (5.9)-(5.16) the graph is complete, then given n zones, the n × n(n -1) matrix B holds that BB * = 2nI n×n -2(1 n×n ).

ADMM applied to the dynamic model

(5.30)

The following proposition shows some properties of matrix with this structure.

Proposition 5.2.1 Set x, y ∈ IR. We consider the n × n matrix P = x(I n×n ) + y(1 n×n ).

Then,

• P has x and x + ny as unique eigenvalues,

• P 2 = max{|x|, |x + ny|}

• if x / ∈ {0, -ny}, then P is invertible and (s k+1 , g k+1 ) = P V +a ((w k+1 , f k+1 ) + λ -1 (d k , r k ))

P -1 = 1 x (I n×n
d k+1 = d k + λ(w k+1 -s k+1 ) r k+1 = r k + λ(f k+1 -g k+1 )
In order to get a more explicit expression of this algorithm, we need to find an adequate manner to express the projection on the affine space V + a . For these we assume that Q On the other hand, since V is the kernel of matrix DA B , one deduces that the projection P V (w, f ) is equal to (w -A t D t (DAA t D t + B B t ) -1 (DAw + Bf ), f -B t (DAA t D t + B B t ) -1 (DAw + Bf )).

To finish, we calculate the inverse of DAA t D t + B B t . Since when decomposition is the final objective to cope with the curse of dimensionality (see [START_REF] Rockafellar | Scenarios and policy aggregation in optimization under uncertainty[END_REF][START_REF] Shapiro | Stochastic Programming[END_REF][START_REF] Barty | Decomposition of large scale stochastic optimal control problems[END_REF]).

Coming back to the stochastic model problem (5.1)-(5.8), we rewrite it in the context of problem (P ).

Set Notice that Aw and Bf in (5.33) are random vectors because w and f are so. They are defined by Aw(ξ) = A(w(ξ)) and Bf (ξ) = B(f (ξ)) for all ξ ∈ Ξ.

We apply algorithm (2.41)-(2.43) with V 1 = 0, V 2 = 0 and M = λI (which is a slight variant of ADMM) to last problem, assuming the random variable space of finite dimension (finite scenarios) with inner product induced by the expectation function, getting the following algorithm.

Stochastic Application 1(SA1)

Step 1. For all z ∈ Z do where

Q = diag (Q z 1 , • • • , Q zn ), with Q z i = diag (Q z i 0 , • • • .Q z i (T -1) ).
Therefore, since (5.42) is a coupling linear constraint, similar to our original model (5.32)-(5.36), we get the following algorithm for solving (5.

f eτ F τ ∀e ∈ E, ∀ τ ∈ [0, T -1]
Step 3. Dual update

y k+1 = y k + λ(Q Aw k+1 + Q Bf k+1 -Qd)
Choosing U zτ equal to (d zτ , i zτ ), we have that Q zτ B z f k t and Q zτ y k zτ are not noise, on the contrary are function of (d zτ , i zτ ). Therefore we can apply DP for solve the sub-problem of Step 1.

In a future work we will try to apply the algorithm TSD, developed in Subsection 3.4.1 of Chapter 3, to sub-problem of Step 1 of algorithm SA1.

Conclusion

The contributions of this thesis are disseminated in the 5 chapters with different relative importance. Chapters 2 and 4 contain the main basic algorithms for two or more operators and they are nearly self-contained. In Chapter 2, the main point is the generalized splitting scheme which includes most of the classical primal-dual splitting methods as particular cases associated with the choices of the blocks of the matrix involved in a generalized proximal point method (a method constructed in this thesis) applied to the saddle-point inclusion problem.

The general setting shows us the relationship between the splitting algorithms and the fixed point algorithms corresponding to special average maps. This general setting also gives us a common point of view of the splitting and convergence properties of the classical primal-dual splitting methods deduced from this general setting, allowing us to improve them by adding for example multi-scaling parameters and a relaxed parameter.

The separable models for multi-block constrained optimization are studied in Chapter 3 and many new decomposition algorithms are derived with block separable augmented Lagrangian subproblems. One of these algorithms (SMS3BF) is applied to a stochastic model defined as (SP ), splitting the nonanticipative constrains and the linear temporal coupling constraints.

In Chapter 4, a Lipschitz-differentiable function or its corresponding co-coercive map is added to the models proposed in Chapter 2. So we get extended version of this algorithm, where we add a Forward step correspond to that function in there formulation. Notice that, under mild assumptions, these extended algorithms inherit the properties of the generalized splitting schemes of Chapter 2.

Finally, in Chapter 5, some of the new algorithms considered in the thesis are applied to an applicative model, the stochastic multistage production planning problem with a limited set of numerical experiments based on randomly generated data sets. More is to be done to further validate the proposed algorithms considering the different coupling inherent to the model. An important open question is concerned with the tuning of the numerous parameters which influence these splitting methods. Even if the theoretical convergence rate analysis presented in Chapters 2 and 4 are not surprising, either in the ergodic or non ergodic sense, the numerical behavior of the splitting methods is still very sensitive to the choice of the scaling parameters, as already observed in the literature (see [START_REF] Eckstein | Splitting methods for monotone operators with applications to parallel optimization[END_REF][START_REF] Lenoir | A survey of monotone operator splitting methods for the decomposition of convex programs[END_REF]).

1

 1 Let γ = 1 (β = α -1), We consider in(2.15)

Figure 3 . 1 :

 31 Figure 3.1: iteration vs parameter

  and D is an n(m -1) × n matrix defined by D = (I n×n • • • I n×n ) t . Notice that problem (RS v ) has the same structure as model (V ar) considered at the beginning of this chapter, when D is injective. So we can apply algorithm (4.4)-(4.8) or his switched version (4.16)-(4.20), getting in both cases splitting algorithms by considering M = λI. These splitting algorithms combine proximal steps on each S i with the forward step on C, because D t D = (m -1)I n×n .

Theorem 4 . 1 . 1

 411 With the same notations as before. Set M and ρ satisfying 2 -

1 2

 1 x), where x is an optimal solution of problem (4.1).

Lemma 4 . 3 . 1

 431 For any (d, d * ), ( d, d * ) ∈ graph (L ), considering d = (x, y) and d = (x, ȳ), it holds d -d, d * ≥ l (x, ȳ) -l (x, y) ≥ d -d, d * . These inequalities are still verified if we consider (d, d * ) ∈ graph (L ) and d ∈ dom (f ) × dom (g * ), for the first inequality; and ( d, d * ) ∈ graph (L ) and d ∈ dom (f ) × dom (g * ), for the second inequality.

Theorem 4 . 3 . 1

 431 With the same notations as before. Set σ, τ and ρ satisfying 2 -τ β ≥ ρ > 0, and 1 ≥ στ A 2 .

Theorem 4 . 3 . 2

 432 With the same notations as before. Set σ, τ and ρ satisfying 2 -τ 2β ≥ ρ > 0, and 1 ≥ στ A 2 .

Proposition 4 . 4 . 1

 441 With the same notations as before, for given (x, y, z) ∈ IR m × IR n × IR p , it holds that

4 . 1 ,

 41 Find r k+1 = (r k+1 • • • , r k+1 q

Remark 4 . 5 . 1

 451 Choosing matrices

z

  δ h , ∀z ∈ Z, τ ∈ [0, T -1] (5.13) 0 ≤ θ j zτ ≤ P j zτ δ h , ∀z ∈ Z, τ ∈ [0, T -1] (5.14) 0 ≤ f eτ ≤ κ eτ δ h , ∀e ∈ E, τ ∈ [0, T -1] (5.15) v z ≥ 0, v z ≥ (x z0 -x zT ), ∀z ∈ Z (5.16)

A

  zτ q zτ -B z f τ = d zτ , ∀(z, τ ) ∈ Z × [0, T -1] (5.18) (q z , v z ) ∈ C z , ∀z ∈ Z,(5.19)f eτ ∈ C eτ , ∀(e, τ ) ∈ E × [0, T -1]. (5.20) Considering f = (f τ ) τ ∈[0,T -1] and d z = (d zτ ) τ ∈[0,T -1] , relation ((5.18)) becomes A z q z -Bz f = d z , ∀z ∈ Z (5.21)whereA z = diag (A z0 , • • • , A z(T -1) ) and Bz = diag (B z , • • • , B z ). Then considering w = (q z , v z )z∈Z and d = (d z ) z∈Z , we have that (5.21) becomes Aw + Bf = d where A = diag ([A z 1 0 T ×1 ], • • • , [A zn 0 T ×1 ]) and B =

  w z ∈ C z , ∀z ∈ Z,(5.24)f eτ ∈ C eτ , ∀e ∈ E, ∀τ ∈ [0, T -1]. (5.25)

- 1 τl eτ f eτ + λ 2 - 1 τA 1 Step 1 . 2 2 f- 1 τl eτ f eτ + 1 2r 2 f -f k + r 2 λ 2 y

 12111221222 We apply the algorithm (2.41)-(2.43) with V 1 = 0, V 2 = 0 and M = λI (which is a slight variant of ADMM) to the problem (5.22)-(5.25) obtaining the algorithmw k+1 ∈ argmin wz∈Cz z∈Z k z (w z ) + λ 2 Aw + Bf k -d + λ -1 y k 2 (5.26) f k+1 ∈ argmin feτ ∈Ceτ e∈E T =0 Aw k+1 + Bf -d + λ -1 y k 2 (5.27)y k+1 = y k + λ(Aw k+1 + Bf k+1 -d) (5.28) Since Aw + Bf k -d + λ -1 y k 2 = z∈Z T =0 zτ q zτ -B z f τ -d zτ -λ -1 y k zτ 2the corresponding minimization problem in (5.26) can be solvable in parallel processing with respect to z indices. Similarly, the corresponding minimization problem in (5.27) can also be solvable in parallel processing with respect to τ indices.Application Zonal subproblemsFor all z ∈ Z do Find q k+1 z = ((θ j zτ ) k+1 j∈Qzτ , u k+1 zτ , η k+1 zτ , x k+1 z,τ +1 ) τ ∈[0,T -1] and v k+1 z θ j zτ + c f ail η zτ + λ Qzτ j=1 θ j zτ + u zτ + η zτ + B z f k τ -d zτ + λ -1 y k zτ 2 +c f in z v zwe obtain a variant of Chambolle-Pock algorithm. This variant is applied to problem (5.22)-(5.25) obtaining the following algorithmw k+1 ∈ argmin wz∈Cz z∈Z k z (w z ) + 1 2r 1 w -w k + r 1 λA t (Aw k + Bf k -d + λ -1 y k ) k+1 ∈ argmin feτ ∈Ceτe∈E T =0 Bt (Aw k+1 + Bf k -d + λ -1 y k ) k+1 = y k + λ(Aw k+1 + Bf k+1 -d).

For

  g = ((g zτ ) τ ∈[0,T -1] ) z∈Z , the following expressions hold (A t g) z = qt z 0 and ( Bt g) τ = z∈Z g zτ B t zwhere qz = (q zτ ) τ ∈[0,T -1] , with qzτ = g zτ 1 Qzτ ×1 1 1 0 .

2 Step 1 .For all e ∈ E do Find f k+1 eτ solution of min feτ l eτ f eτ + 1 2r 2 f eτ -f k eτ + r 2 λa e 2 s- 1 τ 2 w -s k + λ -1 d k 2 f- 1 τl eτ f eτ + λ 2 f

 21212212 ) -y x(x + ny)(1 n×n ).Since matrix B is a row permutation of the nT ×n(n-1)T matrix diag (B, • • • , B), then B 2 = B 2 . Therefore, from(5.30) and the last proposition, we conclude that (5.29) and (5.31) and choosing the parameters r 1 , r 2 ∈ IR satisfying1 ≥ r 1 λ2n and 1 ≥ r 2 λ(2 + max zτ {Q zτ })we get from the previous algorithm, the desired more explicit algorithm:Application Zonal subproblems For all z ∈ Z do calculated For all τ ∈ {0, • • • , T -1} do b τ := ( Qzτ j=1 (θ j zτ ) k+1 + u k+1 zτ + η k+1 zτ + B z f k+1 τ -d zτ + λ -1 y k zτ ) 1 Qzτ ×1 1 1 0 end for Set b = (b) τ ∈[0,T -1] Find q k+1 z = ((θ j zτ ) k+1 j∈Qzτ , u k+1 zτ , η k+1 zτ , x k+1 z,τ +1 ) τ ∈[0,T -1] and v k+1 z zτ + c f ail η zτ + c f in z v z + 1 2r 1 q z -q k z + r 1 λb 2 τ +1 = x zτ -u zτ + i zτ , ∀ τ ∈ [0, T -1] X min zτ ≤ x zτ ≤ X max zτ , ∀ τ ∈ [0, T -1] 0 ≤ u zτ ≤ U max z δ h , ∀ τ ∈ [0, T -1] 0 ≤ θ j zτ ≤ P j zτ δ h , ∀τ ∈ [0, T -1] v z ≥ 0, v z ≥ (x z0 -x zT ), end for Step 2. Network subproblem For all τ ∈ {0, • • • , T -1} do calculated zτ ) k+1 + u k+1 zτ + η k+1 zτ + B z f k τ -d zτ + λ -1 y k zτ A * z .t 0 ≤ f eτ ≤ κ eτ δ h . end for end for Step 3. Dual update For all (z, τ ) ∈ Z × [0, T -1] do y k+1 zτ = y k zτ + λ( Qzτ j=1 (θ j zτ ) k+1 + u k+1 zτ + η k+1 zτ + B z f k+1 τ -d zτ )end for5.2.3 PDA applied to the dynamic modelWe reformulate the problem (5.22)-(5.25) and apply Spingarn method. Seth(w, f ) = z k z (q z , v z ) + e∈E T =0 l eτ f eτ , V = {(w, f ) : Aw + Bf = 0}and also a such that (A B)a = d. Then, problem (5.22)-(5.25) can be set as min(w,f ) h(w, f ) (w, f ) ∈ V + a w z ∈ C z , ∀z ∈ Z, f eτ ∈ C eτ , ∀e ∈ E, τ ∈ [0, T -1].122 Using the Spingarn's algorithm to that problem, we get w k+1 ∈ argmin wz∈Cz z k z (w z ) + λ k+1 ∈ argmin feτ ∈Ceτ e∈E T =0 -g k + λ -1 r k 2

  z 1 τ = • • • = Q znτ . Since P V +a = a -P V (a ) + P V , it suffices to determinate the projection over V . Set B = diag (B, • • • , B) a permutation matrix D such that D B = B, then V = {(w, f ) : DAw + Bf = 0} and DAA t D t = diag ( Ã0 Ãt 0 , • • • , Ã(T-1) Ãt (T -1) ),where, for τ = 0, • • • , T -1, Ãτ := diag (A z 1 τ , • • • , A znτ ).

3 Step 1 .=

 31 DAA t D t + B B t = diag ( Ã0 Ã * 0 + BB * , • • • , ( Ã(T-1) Ãt (T -1) , +BB t ),it suffices to determinate the inverse of each diagonal block Ãτ Ãtτ + BB t (τ = 0, • • • , T -1). From definition, Ãτ Ãt τ + BB t = diag (Q z 1 τ + 2, • • • , Q znτ + 2) + 2nI n×n -2(1 n×n )and, by consideringQ τ = Q z 1 τ = • • • = Q znτ , we deduce, from Proposition 5.2.1, that ( Ãτ Ãt τ + BB t ) -1 = 1 Q τ + 2n + 2 I n×n + 2 (Q τ + 2n + 2)(Q τ + 2)1 n×n 123 and hence the desired explicit expression of the projection is deduced. In particular,a -P V (a ) = (A t D t b, Bt D t b) where b = (b τ ) τ ∈{0,••• ,T -1} is defined as b τ := 1 Q τ + 2n + 2 d τ + 2 z∈Z d zτ (Q τ + 2n + 2)(Q τ + 2) 1 1×n .Therefore, the last previous algorithm can be set in the following contextApplication Zonal subproblems Find q k+1 z = (((θ j zτ ) k+1 j∈Qzτ , u k+1 zτ , η k+1 zτ )) τ ∈[0,T -1] and p k+1 z zτ + c f ail η zτ + c f in z v z + λ 2 (q z , p z ) -s k z + λ -1 d k z 2 s.t x z,τ +1 = x zτ -u zτ + i zτ , ∀ τ ∈ [0, T -1] X min zτ ≤ x zτ ≤ X max zτ , ∀ τ ∈ [0, T -1] 0 ≤ u zτ ≤ U max z δ h , ∀ τ ∈ [0, T -1] 0 ≤ θ j zτ ≤ P j zτ δ h , ∀τ ∈ [0, T -1] v z ≥ 0, v z ≥ (x z0 -x zT ),Step2. Network subproblem For all (e, τ ) ∈ E × {0, • • • , T -1} do Find f k+1 eτ solution of min feτ l eτ f eτ + λ 2 f eτ -g k eτ + λ -1 r k eτ 2 s.t f eτ ∈ C eτ .

1 τ

 1 w zτ = (p zτ , u zτ , x zτ ), w z = (w zτ ) τ ∈[0,T -1] , w = (w zτ ) z∈Z , d = ((d zτ ) τ ∈[0,T -1] ) z∈Z and f = ((f eτ ) e∈E ) τ ∈[0,T -1] z (p zτ ) + ψ z (x zT ) + T -=0 e∈E l e (f eτ )(5.32)Aw + Bf = d (5.33) w z ∈ C z , ∀z ∈ Z,(5.34)f eτ ∈ C eτ , ∀e ∈ E, ∀τ ∈ [0, T -1].(5.35)w zτ F τ , f eτ F τ ∀z ∈ Z, ∀e ∈ E, ∀τ ∈ [0, T -1](5.36)considering matrices A and B defined asA = diag (A z 1 , • • • , A zn ) and B = A z = diag ([1 1], • • • , [11]) and B z = diag (B z , • • • , B z ).

For 1 τ

 1 each τ ∈ [0, T -1] and z ∈ Z, we consider a random variable U zτ F τ , then we consider the approximate version of our problem (5.32)-[START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF]z (p zτ ) + ψ z (x zT ) + T -=0 e∈E l e (f eτ )(5.37)IE[p zτ + u zτ + B z f τ |U zτ ] = IE[d zτ |U zτ ], ∀z ∈ Z, ∀τ ∈ [0, T -1](5.38) w z ∈ C z , ∀z ∈ Z,(5.39)f eτ ∈ C eτ , ∀e ∈ E, ∀τ ∈ [0, T -1].(5.40)w zτ F τ , f eτ F τ , ∀z ∈ Z, ∀e ∈ E, ∀τ ∈ [0, T -1]. (5.41) Set matrix Q zτ (τ = 0, • • • , T -1 and z ∈ Z = {z 1 , • • • , z n }) satisfying Q zτ h = IE[h|U zτ ],then restriction (5.38) become Q Aw + Q Bf = Qd (5.42)

2 Qλ 2 z∈ZQy k zτ 2 f

 222 zτ , u zτ , Q zτ B z f k t , Q zτ y k zτ , d zτ ) + ψ z (x zT ) x z,τ +1 = x zτ -u zτ + i zτ , ∀ τ ∈ [0, T -1] x zτ ∈ X zτ , u zτ ∈ U zτ , p zτ ∈ P zτ , ∀ τ ∈ [0, T -1] p zτ , u zτ F τ , ∀ τ ∈ [0, T -1]where the function L is defined asL(p, u, f , ỹ, d) = c z (p) + Q zτ (p + u -d) + f ), ỹ + λ zτ (p + u -d) + f 2 zτ (p k+1 zτ + u k+1 zτ + B z f τ -d zτ ) + λ -1 eτ ∈ F eτ ∀e ∈ E, ∀ τ ∈ [0, T -1]

  Proposition 2.6.2 are verified and hence the subproblem in Step 4 of the Algorithm (PMA) becomes xk+1

  We now calculate the norms of A t and Bt . To calculate the norm of A t , note thatAA t = diag (A z 1 A t z 1 , • • • , A zn A t zn ) and A z A t z = diag (Q z 1 0 + 2, • • • , Q zn(T -1) + 2). Then,

	A 2 2 = 2 + max

zτ {Q zτ }.

(5.29)

This chapter corresponds to the paper[START_REF] Oré | A unified splitting algorithm for composite monotone inclusions[END_REF] submitted to Journal of Convex Analysis

We consider, for i ∈ {1, • • • , q}, M i = λI in SALA; and for j ∈ {1, • • • , q -1, M j = λI in DSALA.

The graph in (fig 3 .2) describes the relationship of parameter λ versus the necessary number of iterations in order to get an approximation of the optimal value with an error ( x k 1 -x * 1 ) less than 10 -6 , 10 -25 and 10 -50 respectively.

Proximal separation into two sub-blocks

The algorithm PDA described in the first section of this chapter considers the proximal step of all family {f i } i=1,••• ,q in parallel processing. Unlike this algorithm, DPDA separates the family into two sub-families or sub-blocks: one consisting of q -1 functions, f 1 , • • • , f q-1 , and the other consisting of f q . Then, the proximal step of each function of sub-block {f i } i=1,••• ,q-1 is found in parallel processing and then, at a linear combination of all these values, the proximal step of f q is found.

B. He and X.Yuan [START_REF] He | Alternating direction method of multipliers for linear programming[END_REF] considered a linear programming model and a corresponding algorithm separating the problem into two adequate sub-blocks. They then determinate the proximal step of all functions corresponding to one of these sub-block in parallel processing and, using these values, determinate the proximal step for all functions corresponding to the second sub-block in parallel processing too.

We show in this section that this procedure can also be applied for general setting by separating a given problem into two arbitrary sub-blocks. In practice, is more adequate to separate into two sub-blocks taking into account the difficulty to determinate their proximal steps. Since B t 1 B 1 is diagonal, the calculation of xk+1 can be realized in parallel processing: for i ∈ {1, • • • , r},

Similarly, since B t 2 B 2 is diagonal, the calculation of xk+1 is also realized in parallel processing:

where {α i } i=1,••• ,q is defined as

Then we get the following algorithm.

Two Sub-blocks Proximal Decomposition algorithm (2sb-PDA)

Set the finite sequence {α i } i=1,••• ,q previously defined and λ a positive number. Then for an arbitrary (x 0 , z 0 , u 0 ) ∈ IR rn × IR sn × IR n(q-1) ,

Step 1. Find x k+1 For all i ∈ {1, ..., r}

Chapter 4

A new splitting algorithm for inclusion problems mixing a composite monotone plus a co-coercive operator

In this chapter we consider the following composite monotone inclusion:

where S : IR n -→ -→ IR n and T : IR m -→ -→ IR m are maximal monotone maps, C : IR n → IR n is β-co-coercive with full domain and A an m × n matrix. Regarding this formulation, D. Davis and W. Yin's [START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF] analyze the particular case A = I (thus m = n). They reformulate (V ar) as a fixed-point equation with respect to an appropriate average map with similar properties as the Douglas-Rachford map considered for the sum of two maximal monotone maps (C = 0). In fact this new map recovers the Douglas-Rachford map when C = 0. Monotone inclusion problems with three operators have gained a recent increase of interest in the community of splitting methods. It is motivated by many inverse problems in different fields like data analysis, image processing or machine learning when parcimony is a challenge for very large data sets. Primal-dual splitting methods were proposed in the literature as extensions of the classical splitting schemes for two operators, mainly the Douglas-Rachford family [START_REF] Lions | Splitting algorithms for the sum of two nonlinear operators[END_REF], which lead to decomposing the corresponding proximal steps for each operator separately (see [START_REF] Briceno-Arias | A monotone+skew splitting model for composite monotone inclusions in duality[END_REF][START_REF] Combettes | Primal-dual splitting algorithm for solving inclusions with mixtures of composite, lipschitzian and paralle-sum type monotone operators[END_REF][START_REF] Vũ | A splitting algorithm for dual monotone inclusions involving cocoercive operators[END_REF][START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF][START_REF] Komodakis | Playing with duality : an overview of recent primal-dual approaches for solving large-scale optimization problems[END_REF][START_REF] Davis | A three-operator splitting scheme and its optimization applications[END_REF][START_REF] Latafat | Asymmetric forward-backward-adjoint splitting for solving monotone inclusions involving three operators[END_REF][START_REF] Bot | Admm for monotone operators : convergence analysis and rates[END_REF][START_REF] Yan | A new primal-dual algorithm for minimizing the sum of three functions with a linear operator[END_REF][START_REF] Liang | Local linear convergence analysis of primaldual splitting methods[END_REF]).

The case C = 0, was studied in Chapter 2, where we have shown that the Proof. From Proposition 4.1.3 we get

and applying the Jensen's inequality to the convex functions l(•, •, y) -l(x, z, •) for arbitrary fixed element (x, z, y) ∈ dom (f ) × dom (g) × IR m , where l is the Lagrangian function defined in Proposition 4.1.3, the desired ergodic convergence is deduced.

Regarding the nonergodic convergence, for w * ∈ sol (V L ) and considering w = w * in (4.3.1), we get

and hence, from the Cauchy-Schwarz inequality and (4.9), we obtain

.29) On other hand, from (4.9) and since { Dw k+1 -Dw * } is non increasing, we get

So, replacing this last expression and inequality (4.25) in expression (4.29), we deduce the desired nonergodic convergence.

Constraint violations

Following the same arguments given in the proof of Theorem 2.5.2 of Chapter 2, we get the following result.

For any w * ∈ sol (V L ), the following rate of convergence are obtained:

• Ergodic Convergence:

• Nonergodic Convergence:

Notice that ỹk+1 2 = 0, and hence, by considering y 0 2 = 0, we get y k 2 = 0, which implies in particular zk+1

we deduce from the previous sequences our main algorithm termed "Generalized splitting algorithm for three operators (GSA3O)":

We finally deduce the following proposition directly from Proposition 4.1.2.

Proposition 4.2.1 Assume that V ∈ IR n×n and M ∈ IR m×m are symmetric, with V positive semi-definite and M positive definite, such that V + A t M A is positive definite and satisfying that

Similarly to expression (4.9), the sequence {ζ k := (x k , y k )} generated by algorithm (GSA3O), satisfies the following relation

where

Switching the proximal step

In a similar way described for the injective context, we get now an algorithm where the order position of the proximal steps corresponding to maps S and T in algorithm GSA3O are switched. In order to do that, like to map L and matrix D corresponding to problem (V ar), we denote respectively by L and D the map and matrix for problem (V ar 1 ). Then, the map

where

So, by defining

and assuming for arbitrary points x 0 , r 0 , y 0 2 and z

Hence, by denoting

we deduce from the previous chain of sequences our desired new algorithm:

Notice that the sequence {ξ k = (x k -r k , z k , y k )} corresponding to the sequence (4.36)-(4.40), satisfies the following relation

where

A variant of primal-dual Condat's algorithms

In this section, considering special parameter matrices M and V corresponding to algorithms (4.36)-(4.40) and (4.30)-(4.34), we deduce two algorithms which can be seen as variants of the primal-dual Condat's algorithms CA1 and CA2 [START_REF] Condat | A primal-dual splitting method for convex optimization involving lipschitzian, proximable and linear composite terms[END_REF] and their respective construction relationship. We also show the relationship of these new algorithms with Davis-Yin's algorithm, and also their ergodic and nonergodic rate of convergence. Applying (4.36)-(4.40) considering M = σI m×m and V = τ -1 I n×n -σA t A, we get the following sequence:

Now, in order to recover algorithms Alg1 and Alg2 from algorithms (4.50)-(4.53) and (4.54)-(4.57), respectively, we consider the following inclusion problem corresponding with the sum of three operators, as defined in [START_REF] O'connor | On the equivalence of the primal-dual hybrid gradient method and douglas-rachford splitting[END_REF] for optimization problems :

where Ṽ = (τ σ)

Then the aforementioned two algorithms applied to problem (V ar 2 ) are exactly Alg1 and Alg2.

The fact that problem (V ar 2 ) can be deduced from problem (V ar 1 ) follows from the following steps: We first apply the dual formulation to (V ar), which consists in finding y ∈ IR m such that 0 ∈ T -1 (y) -A(S + C) -1 (-A t y).

Then reformulate it as (V ar 1 ) considering M = I, resulting

Finally, the dual formulation of the last inclusion problem considering Ṽ = -V is exactly (V ar 2 ). Conversely, using the same previous arguments, we can show that problem (V ar 1 ) is deduced from problem (V ar 2 ).

Rate of convergence

By Proposition 4.1.1, G L D and G L D are 2β 4β-τ -average maps, if τ < 2β, and their corresponding fixed point sets are respectively

Now, corresponding to algorithms Alg1 and Alg2, we consider, for i = 1, 2, the following sequences:

and the ergodic sequences, for N ≥ 1,

and, associated to matrix P i defined in (4.44), we consider the norm (x, y) 2

where i = 1 corresponds to algorithm Alg1, and i = 2 to Alg2.

Using the upper bound of the saddle-point gap and the fixed-point residual developed in Sections 2.5.2 and 2.5.3 of Chapter 2 the following rate of converge are deduced in the two next results: sequences is reduced to, for k ≥ 1,

). (4.65)

Notice that the update sequences (4.61) -(4.65) is also satisfied for k = 0, if we consider y 0 2 = 0, and y 0 3 and

for arbitrary r 0 2 and z 0 3 .

From Proposition 4.1.2, we get the following convergence result of the sequences described by (4.61) -(4.65) Chapter 5

Application to stochastic problems

In this chapter, we will consider a large-scale production planning problem with a multiple separable structure which is favorable to the use of the splitting techniques which constitute the heart of the present thesis. Motivated by a long-term energy production planning problem, we analyze here a stochastic optimal control problem where three levels of coupling structure are present, namely:

• the coupling of the dynamic equations w.r.t. time intervals;

• the coupling of the scenario tree w.r.t. the so-called non-anticipativity constraints, i.e. which force the decisions at some period t to be the same for scenarios with identical past history before t;

• the spatial coupling which interconnects the local subsystems

Most approaches in the literature to treat stochastic multistage optimization problems use a scenario tree where the non-anticipativity constraints at each node of the tree are dualized to allow the temporal decomposition as if the model was deterministic (see [START_REF] Shapiro | Stochastic Programming[END_REF] for instance).

The stochastic optimization model

We study the model problem presented in [START_REF] Mahey | Decomposition methods for a spatial model for long-term energy pricing problems[END_REF]. Consider a set of agents Z (geographical zones, markets) with interconnections between them given by a graph (E ⊂ Z × Z). Given a finite period time {0, • • • , T -1}, for each agent (z ∈ Z) there are a production (p zτ ), demand (d zτ ), storage (x zτ ) and interchange (f eτ ) of a commodity (electricity, gas). The objective is to minimize the cost associated with the production, interchange, usage of storage of a commodity, in order to satisfy the demand. The usage of the storage (u zτ ) and the storage (x zτ ) evolve in time satisfying a dynamic equation. Finally, uncertainty affects the following data:

end for

Chambolle-Pock applied to the dynamic model

Considering in algorithm (2.41)-(2.43) the positive parameters r 1 , r 2 such that 1 ≥ r 1 λ B 2 and 1 ≥ r 2 λ Ā 2 , and the parameter matrices V 1 , V 2 and M defined as

then calculate

and finally calculate

Step 4. Dual update Calculate r k+1 = ((r k+1 ) e∈E ) τ and d k+1 = ( q k+1 , p k+1 ) where q k+1 zτ := (( θ j zτ ) k+1 j∈Qzτ , u k+1 zτ , η k+1 zτ ) and pk+1

Uncertainty Environment

The general case of stochastic production planning models has been studied by many authors and we will not detail the different discussions which are behind these models

Step 3. Dual update y k+1 = y k + λ( Aw k+1 + Bf k+1 -d).

Notice that the sub-problem of Step 2, can be solvable by the progressing hedging algorithm proposed by Rockafellar and Wets [START_REF] Rockafellar | Scenarios and policy aggregation in optimization under uncertainty[END_REF] because the only restriction coupling e and τ is the nonanticipativity constraint.

The sub-problem of Step 1 is a stochastic optimal control (SOC) which has less variables than original problem. But the white noise assumption over the random variable (d zτ , i zτ ) τ ∈[0,T -1] is not enough in order to apply dynamic programming, because we have two another families of random variables (f Following the same ideas presented in [START_REF] Mahey | Decomposition methods for a spatial model for long-term energy pricing problems[END_REF], we now reformulate problem (5.32)-(5.36) considering information relaxation in order to obtain a variant of sub-problem corresponding to Step 1 of previous algorithm (SA1) where we can apply DP for solve it.