# Smoothing Techniques and Augmented Lagrangian Method for Recourse Problem of Two-Stage Stochastic Linear Programming

Authors: Saeed Ketabchi; Malihe Behboodi-Kahoo
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

The augmented Lagrangian method can be used for solving recourse problems
and obtaining their normal solution in solving two-stage stochastic linear
programming problems. The augmented Lagrangian objective function of
a stochastic linear problem is not twice differentiable which precludes the
use of a Newton method. In this paper, we apply the smoothing techniques
and a fast Newton-Armijo algorithm for solving an unconstrained smooth
reformulation of this problem. Computational results and comparisons are
given to show the effectiveness and speed of the algorithm.


---

## Body

## 1. Introduction

In stochastic programming, some data are random variables with specific possibility distribution [1], which was first introduced by the designer of linear programming problems, Dantzig, in [2].In this paper, we consider the following two-stage stochastic linear program (slp) with recourse which involves the calculation of an expectation over a discrete set of scenarios:(1)minx∈Xf(x)=cTx+ϕ(x),X={x∈ℝn:Ax=b,x≥0},
where
(2)ϕ(x)=E(Q(x,ω))=∑i=1N‍Q(x,ωi)ρ(ωi)
and E shows the expectation of function Q(x,ω) which depend on the random variable ω. The function Q is defined as follows:
(3)Q(x,ω)=miny∈ℝn2{q(ω)Ty∣W(ω)Ty≥h(ω)-T(ω)x},
where A∈ℝm×n, c∈ℝn, and b∈ℝm. Also, in the problem (3) vector of coefficients q(·)∈ℝn2, matrix of coefficients WT(·)∈ℝm2×n2, demand vector h(·)∈ℝm2, and matrix T(·)∈ℝm2×n depend on the random vector ω with support space Ω. The problems (1) and (3) are called master and recourse problems of stochastic programming, respectively.We assume that the problem (3) has a solution for each x∈X and ω∈Ω.In general, the recourse functionϕ(x) is not differentiable everywhere. Therefore, the traditional methods use nonsmooth optimization techniques [3–5]. However, in the last decade, it is proposed smoothing method for recourse function in standard form of recourse problem [6–11]. In this paper, we apply a smooth approximation technique to smooth recourse function that the recourse problem has inequality linear constrained. For more explanation see Section 2. The approximated problem is based on the least two-norm solution of recourse problem. This paper considers the augmented Lagrangian method to obtain least two-norm solution (Section 3). For convenience, Euclidean least two-norm solution of linear programming problem is named normal solution. This effective method contains solving an unconstrained quadratic problem which its objective function is not twice differentiable. To apply a fast Newton method we use the soothing technique and replace plus function by an accurate smooth approximation [12, 13]. In Section 4, the smoothing algorithm and the numerical results are presented. Also, concluding remarks are given in Section 5.We now describe our notation. Leta=[ai] be a vector in ℝn. By a+ we mean a vector in ℝn whose ith entry is 0 if ai<0 and equals ai if ai≥0. By AT we mean the transpose of matrix A, and ∇f(x0) is the gradient of f at x0. For x∈ℝn, ∥x∥ and ∥x∥∞ denote 2-norm and infinity norm, respectively.


## 2. Approximation of Recourse Function

As mentioned the objective function of (1) is nondifferentiable. This disadvantage property occurs on the recourse function. In this section, there is an attempt to approximate it to a differentiable function.Using dual of the problem (3), function Q(x,ω) can be written as follows:
(4)Q(x,ω)=maxz∈ℝm2(h(ω)-T(ω)x)Tzs.t.W(ω)z=q(ω),z≥0.
Unlike the linear recourse function, the quadratic recours function is differentiable. Thus in this paper, the approximation is based on the following quadratic problem with helpful properties:
(5)Qϵ(x,ω)=maxz∈ℝm2(h(ω)-T(ω)x)Tz-ϵ2∥z∥2s.t.W(ω)z=q(ω),z≥0.
The next theorem shows that, for the sufficiently small ϵ>0, the solution of this problem is the normal solution of the problem (4).Theorem 1.
For functionsQ(x,ω) and Qϵ(x,ω) introduced in (4) and (5), the following can be presented:(a)
∃ϵ->0 such that, for each ϵ∈(0,ϵ-], the solution for the problem (5) is the normal solution for the problem (4).(b)
For eachϵ>0, function Qϵ(x,ω) is differentiable with respect to x.(c)
The gradient of functionQϵ(x,ω) at point x is
(6)∇Qϵ(x,ω)=-TT(ω)zϵ*(x,ω)in whichzϵ*(x,ω) is the solution of the problem (5).Proof.
To prove (a), refer to [14, 15].
Also, (b) and (c) can be easily proved considering that functionQϵ(x,ω) is the conjugate of function
(7)p(z)={ϵ2∥z∥2,z∈Z,∞,z∉Z,
where
(8)Z={z∈ℝm2:W(ω)z=q(ω),z≥0},
and Theorems (26-3) and (23-5) in [16].Using the approximated recourse functionQϵ(x,ω), we can define a differentiable approximation function to the objective function of (1):
(9)fϵ(x)=cTx+∑i=1N‍Qϵ(x,ωi)ρ(ωi).
By (6), the gradient of above function exists and is obtained by
(10)∇fϵ(x)=c+∑i=1N‍∇Qϵ(x,ωi)ρ(ωi)=c-∑i=1N‍TT(ωi)zϵ*(x,ωi)ρ(ωi).
This approximation has paved the way to use the optimization algorithm for master problem (1) in which the objective function is substituted by fϵ(x)(11)minx∈Xfϵ(x).
In [7], it is considered slp problem with inequality constrained in master problem and equality constrained in recourse problem. Also, in Theorem 2.3 of [7], it is shown that a solution of the approximated problem is a good approximation to a solution of master problem. Here we can express a similar theorem for the problem (1) by using the similar technique in the proof of Theorem 2.3 in [7].Theorem 2.
Consider the problem (1). Then, for any x∈X, there exists an ϵ-(x)>0 such that for any ϵ∈(0,ϵ-(x)](12)|f(x)-fϵ(x)|≤ϵ2M,
where M is defined as follows:
(13)M=maxi=1,2,…,N∥zϵ*(x,ωi)∥2.
Let x* be a solution of (1) and xϵ* a solution of (11). Then, there exists an ϵ->0 such that for any 0<ϵ≤ϵ-(14)max{f(xϵ*)-f(x*),fϵ(x*)-fϵ(xϵ*)}≤ϵ2M.
Further, one assumes that f or fϵ are strongly convex on X with modulus μ>0. Then,
(15)∥x*-xϵ*∥≤Mϵμ.According to Theorem1, it can be found that for obtaining the gradient of function fϵ(x) in each iteration, we need the normal solution of N linear programming problems (4). In this paper, the augmented Lagrangian method [17] is used for this purpose.


## 3. Smooth Approximation and Augmented Lagrangian Method

In the augmented Lagrangian method, the unconstrained maximization problem is solved which gives the project of a point on the solution set of the problem (4).Assume thatz^ is an arbitrary vector. Consider the problem of finding the least 2-norm projection z^* of z^ on the solution set Z* of the problem (4)
(16)12∥z^*-z^∥2=minz∈Z*12∥z-z^∥2,Z*={z∈ℝm2:W(ω)z=q(ω),ξTz=Q(x,ω),z≥0}.
In this problem, vector x and random variable ω are constants; therefore, for simplicity, this is assumed to be ξ=h(ω)-T(ω)x, and function Q^(ξ) is defined in a way that Q^(ξ)=Q(x,ω).Considering that the objective function of the problem (16) is strictly convex, its solution is unique. Let us introduce the Lagrangian function for the problem (16) as follow:
(17)L(z,p,β,z^,ξ,ω)=12∥z-z^∥2+pT(W(ω)z-q(ω))+β(ξTz-Q^(ξ)),
where p∈ℝn2 and β∈ℝ are Lagrangian multipliers and ξ, z^ are constant values. Therefore, the dual problem of (16) becomes
(18)maxβ∈ℝmaxp∈ℝn2minz∈ℝ+m2L(z,p,β,z^,ξ,ω).
By solving the inner minimization of the problem (18), duality of the problem (16) is obtained:
(19)maxβ∈ℝmaxp∈ℝn2L^(p,β,z^,ξ),
where duality function is
(20)L^(p,β,z^,ξ)=qT(ω)p-12∥(z^+WT(ω)p+βξ)+∥2+βQ^(ξ)+12∥z^∥2.
The following theorem states that if β is sufficiently large, solving the inner maximization of (19) gives the solution of the problem (16).Theorem 3 (see [17]).
Consider the following maximization problem(21)maxp∈ℝn2S(p,β,z^,ξ,ω)
in which β, z^, and ξ are constants, and function S(p,β,z^,ξ) is introduced as follows:
(22)S(p,β,z^,ξ,ω)=qT(ω)p-12∥(z^+WT(ω)p+βξ)+∥2.
Also, assume that the set Z* is nonempty, and the rank of submatrix Wl of W corresponding to nonzero components of z^* is n2. In such a case, there is β* which for all β≥β*, z^*=(z^+WTp(β)+βξ)+ is the unique and exact solution for the problem (16), where p(β) is the point obtained from solving the problem (21).Also, in special conditions, the solution for the problem (3) can be also obtained and the following theorem expresses this issue.Theorem 4 (see [17]).
Assume that the solution setZ* is nonempty. For each β>0 and z^∈Z*, y*=p(β)/β is one exact solution for the linear programming problem (3), where p(β) is the solution for the problem (21).According to the theorems mentioned above, augmented Lagrangian method presents the following iteration process for solving the problem (16):
(23)pk+1∈argmaxp∈ℝn2{qT(ω)p-12∥(zk+W(ω)Tp+βξ)+∥2},zk+1=(zk+WT(ω)pk+1+βξ)+,
where z0 is an arbitrary vector and here we can use zero vector as initial vector for obtaining normal solution of the problem (4).We note that the problem (23) is a concave problem and its objective function is piecewise quadratic and is not twice differentiable. Applying the smoothing techniques [18, 19] and replacing x+ by a smooth approximation, we transform this problem to a twice continuously differentiable problem.Chen and Mangasarian [19] introduced a family of smoothing functions, which is built as follows. Let ρ:R→[0,∞) be a piecewise continuous density function satisfying
(24)∫‍-∞+∞ρ(s)ds=1,∫‍-∞+∞|s|ρ(s)ds<∞.
It is obvious that the derivative of plus function is step function, that is, (x)+=∫‍-∞xδ(t)dt, where the step function δ(x) is defined 1 if x>0 and equals 0 if x≤0. Therefore, a smoothing approximation function of the plus function is defined by
(25)φ(x,α)=∫‍-∞xψ(t,α)dt,
where ψ(x,α) is smoothing approximation function of step function and is defined as
(26)ψ(x,α)=∫‍-∞xαρ(αt)dt.
By choosing
(27)ρ(s)=e-s(1+e-s)2,
specific cases of these approaches are obtained as follows:
(28)ψ(x,α)=11+e-αx≈δ(x),φ(x,α)=x+1αlog(1+e-αx)≈(x)+.
The function φ with a smoothing parameter α is used here to replace the plus function of (22) to obtain a smooth reformulation of function (22):
(29)S^(p,β,z^,ξ,β,ω,α):=qT(ω)p-∥φ(z^+WT(ω)p+βξ,α)∥2.
Therefore, we have the following iterative process instead of (23) and (28):
(30)pk+1∈argmaxp∈ℝn2{qT(ω)p-∥φ(zk+WT(ω)p+βξ,α)∥2},zk+1=φ(zk+WT(ω)pk+1+βξ,α).
It can be shown that as the smoothing parameter α approaches infinity any solution of smooth problem (29) approaches the solution of the equivalent problem (22) (see [19]).We begin with a simple lemma that bounds the square difference between the plus functionx+ and its smooth approximation φ(x,α).Lemma 5 (see [13]).
Forx∈ℝ and |x|<ϖ(31)φ2(x,α)-(x+)2≤(log(2)α)2+2ϖαlog(2),
where φ(x,α) is the φ function of (28) with smoothing parameter α>0.Theorem 6.
Consider the problems (21) and
(32)maxp∈ℝn2S^(p,β,z^,ξ,ω,α).
Then, for any p∈ℝn2 and α>0(33)|S(p,β,z^,ξ,ω)-S^(p,β,z^,ξ,ω,α)|≤(log(2)α)2+2Mm2log(2)α,
where M is defined as follows:
(34)M=max1≤i≤m2|z^i+WiT(ω)p+βξi|.
Let p* be a solution of (21) and pα* a solution of (32). Then
(35)max{S^S(p*,β,z^,ξ,ω)-S(pα*,β,z^,ξ,ω),S^(pα*,β,z^,ξ,ω,α)-S^(p*,β,z^,ξ,ω,α)}≤(log(2)α)2+2Mm2log(2)α.
Further, one assumes that WT is a full rank matrix. Then,
(36)∥p*-pα*∥≤μ(2log(2)α)2+8Mm2μlog(2)α.Proof.
For anyα>0 and p∈ℝn2(37)φ(z^i+WiT(ω)p+βξi,α)≥(z^i+WiT(ω)p+βξi,α)+.
Hence
(38)|S(p,β,z^,ξ,ω)-S^(p,β,z^,ξ,ω,α)|=S(p,β,z^,ξ,ω)-S^(p,β,z^,ξ,ω,α)=∥φ(z^+WT(ω)p+βξ,α)∥2-∥(z^+WT(ω)p+βξ,α)+∥2=∑i=1m2‍((z^i+WiT(ω)p+βξi,α)+2φ2(z^i+WiT(ω)p+βξi,α)-(z^i+WiT(ω)p+βξi,α)+2).
By using Lemma 5, we get that
(39)|S(p,β,z^,ξ,ω)-S^(p,β,z^,ξ,ω,α)|≤∑i=1m2‍((log(2)α)2+2|z^i+WiT(ω)p+βξi,α|log(2)α)≤(log(2)α)2+2log(2)α×∑i=1m2‍max1≤i≤m2|z^i+WiT(ω)p+βξi|=(log(2)α)2+2Mm2log(2)α.
From above inequality, we have
(40)S^(p,β,z^,ξ,ω,α)≤S(p,β,z^,ξ,ω)≤S^(p,β,z^,ξ,ω,α)+(log(2)α)2+2Mm2log(2)α.
Therefore
(41)S(p*,β,z^,ξ,ω)-S(pα*,β,z^,ξ,ω)≤S^(p*,β,z^,ξ,ω,α)-S^(pα*,β,z^,ξ,ω,α)+(log(2)α)2+2Mm2log(2)α≤(log(2)α)2+2Mm2log(2)α,S^(pα*,β,z^,ξ,ω,α)-S^(p*,β,z^,ξ,ω,α)≤S(pα*,β,z^,ξ,ω)-S(p*,β,z^,ξ,ω)+(log(2)α)2+2Mm2log(2)α≤(log(2)α)2+2Mm2log(2)α.
Suppose that WT is full rank. Then the Hessian of S^ is negative definite, and S^ is strongly concave on bounded sets. By the definition of strong concavity, for any γ∈(0,1),
(42)S^(γpα*+(1-γ)p*,β,z^,ξ,ω,α)-γS^(pα*,β,z^,ξ,ω,α)-(1-γ)S^(p*,β,z^,ξ,ω,α)≥12μγ(1-γ)∥pα*-p*∥2.
Let γ=1/2, then
(43)18μ∥pα*-p*∥2≤S^(12(pα*+p*),β,z^,ξ,ω,α)-12(S^(pα*,β,z^,ξ,ω,α)+S^(p*,β,z^,ξ,ω,α))≤S^(pα*,β,z^,ξ,ω,α)-12(S^(pα*,β,z^,ξ,ω,α)+S^(p*,β,z^,ξ,ω,α))≤12(S^(pα*,β,z^,ξ,ω,α)-S^(p*,β,z^,ξ,ω,α))≤12(log(2)α)2+Mm2log(2)α.Considering the advantage of the twice differentiability of the objective function of the problem (32) allows us to use a quadratically convergent Newton algorithm with an Armijo stepsize [20] that makes the algorithm globally convergent.


## 4. Numerical Results and Algorithm

In each iteration of the process (30), one concave, quadratic, unconstrained maximization problem is solved. For solving it, the fast Newton method can be used.In the algorithm, the Hessian matrix may be singular, thus we use a modified Newton. The direction in each iteration for solving (30) is obtained through the following relation:
(44)ds=-(∇p2S^(p,β,z^,ξ,ω,α)-δIn2)-1(∇pS^(p,β,z^,ξ,ω,α)),ps+1=ps+λsds,
where δ is a small positive number, In2 is the identity matrix of order n2, and λs is the suitable step length that Armijo algorithm is used for determining it (see Algorithm 1).Algorithm 1:Newton method with the Armijo rule.
Choose az0∈Rm2, α>1,ι>1, ϵ>0 be error tolerance and δ is a small positive
number.
i
=
0;
While  ∥zi-zi-1∥∞≥ϵ
Choose ap0∈Rn2 and set k=0.
While  ∥∇pS^(pk,β,zi,ξ,ω,α)∥∞≥ϵ
Choose  λk=  max{s,sσ,sσ2,…}such that
S
^
(
p
k
,
β
,
z
i
,
ξ
,
ω
,
α
)
-
S
^
(
p
k
+
λ
k
d
k
,
β
,
z
i
,
ξ
,
ω
,
α
)
≥
-
λ
i
μ
∇
(
S
^
(
p
k
,
β
,
z
i
,
ξ
,
ω
,
α
)
)
T
d
k,
where,
d
k
=
-
(
∇
2
S
^
(
p
k
,
β
,
z
i
,
ξ
,
ω
,
α
,
)
-
δ
I
n
2
)
-
1
∇
S
^
(
p
k
,
β
,
z
i
,
ξ
,
ω
,
α
)
,
s
>
0  be a constant,
σ
∈
(
0,1
)  andμ∈(0,1).
Put
p
k
+
1
=
p
k
+
λ
k
d
k, k=k+1  andα=ια.
end
Set
p
i
+
1
=
p
k
+
1
,
z
i
+
1
=
φ
(
z
i
+
W
T
(
ω
)
p
i
+
1
+
β
ξ
,
α
)  andi=i+1.
endThe proposed algorithm was applied to solve some recourse problems. Table1 compares this algorithm with CPLEX v. 12.1 solver for quadratic convex programming problems (5). As is evident from Table 1, most of recourse problems could be solved more successful by the algorithm which is based on smooth augmented Lagrangian Newton method (SALN) than CPLEX package (for illustration see the problems 21–25 in Table 1). This algorithm gives us high accuracy and the solution with minimum norm in suitable time (see last column of Table 1). Also, we can find that CPLEX is better than the algorithm proposed for some recourse problems in which the matrices are approximately square (Ex. line 5–12).Table 1
Comparative between smooth augmented Lagrangian Newton method (SALN) and CPLEX solver.
N. P
Recourse problemn2×m2×d
Solver
∥
W
z
-
q
∥
∞
|
Q
^
(
ξ
)
-
ξ
T
z
|
∥
z
∥
Time (second)
P1
50
×
50
×
0.68
SALN
9.9135
e
-
011
7.2760
e
-
012
41.0362
0.3292
CPLEX
1.0717
e
-
007
5.2589
e
-
007
41.0362
0.3182
P2
100
×
105
×
0.4
SALN
1.8622
e
-
010
4.3074
e
-
009
57.1793
0.1275
CPLEX
6.5591
e
-
008
1.6105
e
-
006
57.1826
0.1585
P3
150
×
150
×
0.5
SALN
3.1559
e
-
010
5.6361
e
-
009
74.9098
0.6593
CPLEX
5.9572
e
-
010
1.2014
e
-
009
74.9098
0.1605
P4
200
×
200
×
0.5
SALN
5.3819
e
-
010
1.0506
e
-
008
85.6646
0.2530
CPLEX
1.1734
e
-
007
1.3964
e
-
006
85.6646
0.1820
P5
300
×
300
×
0.5
SALN
9.4178
e
-
010
2.7951
e
-
008
102.4325
2.1356
CPLEX
4.1638
e
-
010
4.4456
e
-
009
102.4325
0.1830
P6
350
×
350
×
0.5
SALN
1.2787
e
-
009
2.6226
e
-
008
110.4189
3.2102
CPLEX
7.7398
e
-
010
3.4452
e
-
009
110.4189
0.2116
P7
450
×
500
×
0.05
SALN
1.6564
e
-
010
4.7094
e
-
009
124.1204
0.7807
CPLEX
9.0949
e
-
013
1.6371
e
-
011
124.1205
0.2606
P8
500
×
550
×
0.04
SALN
1.0425
e
-
010
1.0241
e
-
009
134.3999
0.7567
CPLEX
6.1618
e
-
011
1.9081
e
-
009
134.3999
0.2660
P9
700
×
800
×
0.6
SALN
4.6139
e
-
009
2.3908
e
-
007
153.9782
7.1364
CPLEX
9.1022
e
-
009
1.3768
e
-
007
153.9906
0.8435
P10
900
×
1100
×
0.4
SALN
5.3396
e
-
009
1.5207
e
-
007
178.1151
7.9383
CPLEX
4.5020
e
-
011
2.0373
e
-
010
178.1163
1.3643
P11
1500
×
2000
×
0.1
SALN
1.1289
e
-
008
2.0696
e
-
007
231.5284
13.2493
CPLEX
7.5886
e
-
011
1.1059
e
-
009
231.5286
2.1343
P12
1000
×
2000
×
0.01
SALN
3.6398
e
-
010
7.9162
e
-
009
198.4905
5.2620
CPLEX
6.8212
e
-
013
1.1642
e
-
010
198.4922
0.6752
P13
1000
×
5000
×
0.001
SALN
2.7853
e
-
010
1.9281
e
-
010
190.0141
0.2709
CPLEX
4.8203
e
-
010
1.5425
e
-
009
190.0142
0.2162
P14
1000
×
10000
×
0.001
SALN
1.1221
e
-
010
1.4988
e
-
009
212.2416
0.4867
CPLEX
7.9094
e
-
009
2.6691
e
-
007
212.3453
0.3353
P15
1000
×
1
e
5
×
0.001
SALN
9.7702
e
-
010
2.3283
e
-
009
231.7930
3.7511
CPLEX
2.7285
e
-
012
1.8044
e
-
009
231.8763
1.2472
P16
1000
×
1
e
6
×
0.0002
SALN
9.9432
e
-
013
2.0464
e
-
012
121.3937
9.0948
CPLEX
9.9098
e
-
006
1.1089
e
-
004
121.3940
3.7848
P17
100
×
1
e
6
×
0.001
SALN
3.9563
e
-
010
5.6607
e
-
009
73.8493
6.8537
CPLEX
2.0082
e
-
003
1.1777
e
-
002
73.9412
2.5582
P18
10
×
1
e
4
×
0.001
SALN
2.2737
e
-
013
7.9581
e
-
013
19.5735
0.0166
CPLEX
2.2737
e
-
013
1.1369
e
-
013
19.5739
0.1386
P19
10
×
1
e
6
×
0.001
SALN
1.2478
e
-
009
1.0710
e
-
008
18.8192
5.6399
CPLEX
5.9615
e
-
004
6.4811
e
-
005
18.8863
2.5623
P20
10
×
5
e
6
×
0.01
SALN
2.0425
e
-
008
1.1816
e
-
008
20.8470
28.7339
CPLEX
4.9966
e
+
004
3.3500
e
+
005
0.0000
1.9482
P21
100
×
1
e
6
×
0.01
SALN
3.8654
e
-
012
7.7875
e
-
012
42.0698
7.8895
CPLEX
1.7994
e
-
004
6.2712
e
-
005
42.0931
8.8324
P22
100
×
1
e
4
×
0.1
SALN
1.1084
e
-
012
2.1600
e
-
012
39.4563
0.1440
CPLEX
1.2518
e
-
005
1.3174
e
-
005
39.4563
0.3099
P23
100
×
1
e
5
×
0.05
SALN
1.7053
e
-
012
6.3121
e
-
012
43.5944
1.1379
CPLEX
3.9827
e
-
006
2.8309
e
-
006
43.5944
1.5729
P24
1000
×
5
e
4
×
0.1
SALN
1.3074
e
-
012
2.4102
e
-
011
117.4693
15.1065
CPLEX
6.7455
e
-
007
6.7379
e
-
006
117.4699
20.3967
P25
1000
×
1
e
5
×
0.08
SALN
2.5011
e
-
012
1.0516
e
-
011
116.6964
23.2540
CPLEX
9.4298
e
-
006
5.5861
e
-
004
116.6964
33.4319The test generator generates recourse problems. These problems are generated using the MATLAB code show in Algorithm2.Algorithm 2:
%
Sgen: Generate random solvable recourse problems:
%
Input: m,n,d(ensity); Output: W,q,
ξ
;
m=input(’Entern2:’)
n=input(’Enterm2:’)
d=input(’Enter d:’)
pl=inline(’(abs(x)+x)/2’)
W=sprand(
n
2
,
m
2
,d);W=100*(W-0.5*spones(W));
z=sparse(10*pl(rand(
m
2
,1)));
q=W*z;
y=spdiags((sign(pl(rand(
n
2
,1)-rand(
n
2
,1)))),0,
n
2
,
n
2
)
*5*((rand(
n
2
,1)-rand(
n
2
,1)));
ξ=W’*y-10*spdiags((ones(m2,1)-sign(z)),0,m2,m2)*ones(m2,1));
format short e; nnz(W)
/
prod(size(W))The algorithm considered for solving several recourse problems was run on a computer with 2.5 dual-core CPU and 4 GB memory in MATLAB 7.8 programming environment. Also, in the generated problems, recourse matrixW is the Sparse matrix (n2×m2) with the density d. The constants β and δ in the above algorithm in (44) were selected 1 and 10-8, respectively.In Table1, the second column indicates the size and density of matrix W, the forth column indicates the feasibility of the primal problem (4), and the next column indicates the error norm function of this problem (the MATLAB code of this paper is available from the authors upon request).


## 5. Conclusion

In this paper, a smooth reformulation process, based on augmented Lagrangian algorithm, was proposed for obtaining the normal solution of recourse problem of a stochastic linear programming. This smooth iterative process allows us to use a quadratically convergent Newton algorithm, which accelerates obtaining the normal solution.Table1 shows that the proposed algorithm has appropriate speed in most of the problems. This result, specifically, can be observed in recourse problems with the matrix of coefficients in which the number of constraints is noticeably more than the number of variables. The more challenging is solving the problems which their coefficient matrix is square (the numbers of constraints and variables get closer to each other), and more time is needed by the algorithm for solving the problem.



---

*Source: 735916-2013-06-06.xml*