Loading [MathJax]/jax/output/SVG/jax.js
Research article

A polynomial local optimality condition for the concave piecewise linear network flow problem

  • Received: 12 October 2020 Accepted: 04 December 2020 Published: 09 December 2020
  • MSC : 90B06, 90C26, 90C46

  • This paper studies the local optimality condition for the widely applied concave piecewise linear network flow problem (CPLNFP). Traditionally, for CPLNFP the complexity of checking the local optimality condition is exponentially related to the number of active arcs (i.e., arcs in which the flow is at the breakpoints). When the scale of CPLNFP is large, even local optimality is unverifiable and the corresponding local algorithms are inefficient. To overcome this shortcoming, a new local optimality condition is given. This new condition is based on the concavity and piecewise linearity of the cost function and makes full use of the network structure. It is proven that the complexity of the new condition is polynomial. Therefore, using the new condition to verify the local optimality is far superior to using the traditional condition, especially when there are many active arcs.

    Citation: Zhibin Nie, Shuning Wang, Xiaolin Huang. A polynomial local optimality condition for the concave piecewise linear network flow problem[J]. AIMS Mathematics, 2021, 6(3): 2094-2113. doi: 10.3934/math.2021128

    Related Papers:

    [1] Woocheol Choi, Young-Pil Choi . A sharp error analysis for the DG method of optimal control problems. AIMS Mathematics, 2022, 7(5): 9117-9155. doi: 10.3934/math.2022506
    [2] Anouar Ben Mabrouk, Abdulaziz Alanazi, Zaid Bassfar, Dalal Alanazi . New hybrid model for nonlinear systems via Takagi-Sugeno fuzzy approach. AIMS Mathematics, 2024, 9(9): 23197-23220. doi: 10.3934/math.20241128
    [3] Aly R. Seadawy, Bayan Alsaedi . Contraction of variational principle and optical soliton solutions for two models of nonlinear Schrödinger equation with polynomial law nonlinearity. AIMS Mathematics, 2024, 9(3): 6336-6367. doi: 10.3934/math.2024309
    [4] Erfeng Xu, Wenxing Xiao, Yonggang Chen . Local stabilization for a hyperchaotic finance system via time-delayed feedback based on discrete-time observations. AIMS Mathematics, 2023, 8(9): 20510-20529. doi: 10.3934/math.20231045
    [5] Zahra Pirouzeh, Mohammad Hadi Noori Skandari, Kamele Nassiri Pirbazari, Stanford Shateyi . A pseudo-spectral approach for optimal control problems of variable-order fractional integro-differential equations. AIMS Mathematics, 2024, 9(9): 23692-23710. doi: 10.3934/math.20241151
    [6] Jie Liu, Zhaojie Zhou . Finite element approximation of time fractional optimal control problem with integral state constraint. AIMS Mathematics, 2021, 6(1): 979-997. doi: 10.3934/math.2021059
    [7] Hsien-Chung Wu . Numerical method for solving the continuous-time linear programming problems with time-dependent matrices and piecewise continuous functions. AIMS Mathematics, 2020, 5(6): 5572-5627. doi: 10.3934/math.2020358
    [8] Shima Soleimani Manesh, Mansour Saraj, Mahmood Alizadeh, Maryam Momeni . On robust weakly $ \varepsilon $-efficient solutions for multi-objective fractional programming problems under data uncertainty. AIMS Mathematics, 2022, 7(2): 2331-2347. doi: 10.3934/math.2022132
    [9] Mei Li, Wanqiang Shen . Integral method from even to odd order for trigonometric B-spline basis. AIMS Mathematics, 2024, 9(12): 36470-36492. doi: 10.3934/math.20241729
    [10] Şuayip Toprakseven, Seza Dinibutun . Error estimations of a weak Galerkin finite element method for a linear system of $ \ell\geq 2 $ coupled singularly perturbed reaction-diffusion equations in the energy and balanced norms. AIMS Mathematics, 2023, 8(7): 15427-15465. doi: 10.3934/math.2023788
  • This paper studies the local optimality condition for the widely applied concave piecewise linear network flow problem (CPLNFP). Traditionally, for CPLNFP the complexity of checking the local optimality condition is exponentially related to the number of active arcs (i.e., arcs in which the flow is at the breakpoints). When the scale of CPLNFP is large, even local optimality is unverifiable and the corresponding local algorithms are inefficient. To overcome this shortcoming, a new local optimality condition is given. This new condition is based on the concavity and piecewise linearity of the cost function and makes full use of the network structure. It is proven that the complexity of the new condition is polynomial. Therefore, using the new condition to verify the local optimality is far superior to using the traditional condition, especially when there are many active arcs.



    The network flow problem (NFP) is a fundamental and important task in network optimization. NFP has been widely applied in transportation, communication, construction project [1], manufacturing [2], network design [3], water resource management [4], and other fields. NFP is to minimize a cost function under linear constraints. Thus, its complexity strongly depends on the category of the cost function. Using a linear function as the cost is the simplest case. But in practice, the cost function is often concave due to the marginal utility. As described in [5], the concave cost is or can be well approximated by a concave piecewise linear function. Moreover, NFP with general nonlinear cost functions can be transformed into NFP with concave cost functions on an expanded network [6], so it is of great significance to investigate the concave piecewise linear network flow problem (CPLNFP).

    Let G=(N,A) be a directed network with n nodes and m arcs, where N and A are the sets of nodes and arcs, respectively. Each arc (i,j)A is associated with a flow xij, a capacity uij, and a cost function fij(xij). For the i-th node iN, the supply or the demand flow is denoted by bi. If bi>0, then node i is a supply node; if bi<0, then node i is a demand node; otherwise, node i is a transshipment node. Due to the balance of network flow, the total supply flow must be equal to the total demand flow, which means iNbi=0. With those notations, a CPLNFP defined on network G can be formulated as follows:

    min (i,j)Afij(xij)s.t.{j|(i,j)A}xij{j|(j,i)A}xji=bi,iN,0xijuij,(i,j)A, (1.1)

    where fij(xij) is a continuous concave piecewise linear function that can be expressed as

    fij(xij)={c1ijxij+0,      xij[0,λ1ij],c2ijxij+d2ij,    xij[λ1ij,λ2ij],           csijxij+dsij,    xij[λs1ij,λsij],           csijijxij+dsijij,  xij[λsij1ij,uij]. (1.2)

    In (1.2), fij(xij) has sij segments, and the interval [0,uij] is accordingly divided into sij subintervals with breakpoints λ1ij,λ2ij,,λsij1ij. To make notations consistent, we set d1ij=0, λ0ij=0 and λsijij=uij. Then, in each subinterval [λs1ij,λsij], the objective function fij(xij) becomes a linear function with slope csij and intercept dsij. Due to the continuity and concavity of fij(xij), we have csijλsij+dsij=cs+1ijλsij+ds+1ij,s=1,2,,sij1, and c1ij>c2ij>>csijij.

    CPLNFP has been widely applied in many areas, including facilities location [7] and inventory management [8]. Since CPLNFP is nonconvex, the global optimality cannot be guaranteed, and researchers focus on its local optimality, which also serves as the basis of many algorithms [9,10]. The classical local algorithms include the dynamic slope-scaling procedure (DSSP) [11] and the dynamic cost updating procedure (DCUP) [12]. DSSP uses a linear function to approximate the concave cost function based on the iterative solution. DCUP relaxes the CPLNFP as a bilinear programming problem that can be transformed into a series of linear NFPs.

    As a concave piecewise linear optimization over a polyhedron, the locally optimal solution of the CPLNFP can be achieved at some vertices of the feasible domain. Algebraically, the vertex refers to the basic feasible solution of the constraint system in (1.1). In many local algorithms, such as the DCUP and DSSP, the vertex is obtained by solving a linear NFP and then its local optimality is checked.

    By definition, the feasible domain of a piecewise linear function can be divided into many regions in which the function stays linear. Each feasible point, including the vertex, belongs to one or more regions. Moreover, for each point, we can certainly find a neighbourhood whose intersection with the feasible domain is contained in the union of regions that this point belongs to. Figure 1 illustrates these direct but important properties with a two-dimensional concave piecewise linear function. For example, it can be seen that the feasible domain can be divided into four regions and v1 belongs to two regions: Ω1 and Ω2.

    Figure 1.  A piecewise linear function g(y1,y2) defined on a two-dimensional polyhedron that can be divided into four regions.

    From the above discussing, for concave piecewise linear minimization over a polyhedron, a vertex that is locally optimal is equivalent means that it is optimal in all the regions that it belongs to. This local optimality condition is geometrically intuitive: If we can find a region in which this vertex is not optimal, there must be a feasible descent direction in this region, so this vertex is not locally optimal. For example, in Figure 1, v4 is a locally optimal point if and only if v4 is an optimal point in Ω3 and Ω4. For more detailed discussions about local optimality for piecewise linear programming, see [13], [14] and [15]. Following the idea of checking optimality for all regions, we come to the traditional local optimality condition for the CPLNFP.

    A vertex v in the feasible domain is locally optimal for the CPLNFP if and only if v is optimal in all the regions that v belongs to.

    Consider a vertex v in the feasible domain of the CPLNFP. Based on the characteristics of the cost function, the regions that v belongs to can be obtained by restricting the flow associated with arc (i,j) to a subinterval of [0,uij] that vij belongs to. Let m1 be the number of active arcs (i.e., arcs in which the flow is at the breakpoints), and then clearly the number of regions that v belongs to is 2m1 since the flow in active arcs belongs to two subintervals at the same time. In each region, the concave piecewise linear cost function is a linear function, from which a linear NFP is defined with the constraints of the CPLNFP; whether v is optimal in this region is equivalent to whether v is optimal for this linear NFP. Therefore, to verify the local optimality, one needs to check the optimality of v for 2m1 linear NFPs, which is obviously inefficient when the problem dimension, the number of segments of the cost function, and thus the number of active arcs increase. In [12], the number of active arcs (m1) is assumed to be zero and based on this assumption, DCUP was designed to obtain a local optimum. However, this assumption exclude many real situations. The interested readers could find a numerical experiment in Appendix A, showing the estimated m1 for different CPLNFPs. In Figure 2, we plot the results: m1 is not always zero and becomes large and unacceptable as the problem dimension increases, from which it follows that using the traditional local optimality condition is inefficient.

    Figure 2.  The number of linear NFPs that need to be verified is 2m1, and this figure plots the increasing trend of m1 with respect to the problem dimension. For example, when m=8000,sij=10,(i,j)A, m1 is approximately equal to 42 and 2m1 is more than 4×1012, which is certainly unacceptable.

    In this paper, we provide a new local optimality condition, which is based on the network structure and the concavity of fij(xij),(i,j)A. The complexity of implementing the proposed condition is polynomial with the network parameters. After briefly introducing the background, the proposed condition, the proof, and the strategy of implementation will be given in section 3.

    To better represent the parameters associated with the arcs in A, we first define a tuple collection Θ. Each tuple in Θ has m elements, having one-to-one correspondence with the arcs in A. For any tuple γΘ, the element corresponding to arc (i,j)A is a real number and denoted as γij. In other words, the tuple γ can be expressed as

    γ=(γij,(i,j)A).

    Based on Θ, we can easily represent the parameters associated with the arcs in A. For example, we use xΘ to denote flows in all the arcs in A, and express the feasible domain of the CPLNFP as below,

    Ω={xΘ|{j|(i,j)A}xij{j|(j,i)A}xji=bi,iN          0xijuij,(i,j)A}.

    By this way, the CPLNFP can be succinctly written as

     min  f(x)=(i,j)Afij(xij)s.t.  xΩ. (2.1)

    As discussed, the feasible domain can be divided into a number of regions obtained by restricting xij,(i,j)A to one of the subintervals in [0,uij]. Evidently, the cost function stays linear in each region. Before introducing the new local optimality condition, we give the following definitions to describe the related regions and linear NFPs for a feasible point ˆxΩ.

    Definition 1. For any ˆxΩ and (i,j)A, we define Ψij(ˆxij) as follows:

    Ψij(ˆxij)={1},if ˆxij=0,Ψij(ˆxij)={s},if λs1ij<ˆxij<λsij,Ψij(ˆxij)={s,s+1},if ˆxij=λsij,ssij,Ψij(ˆxij)={sij},if ˆxij=uij.

    Ψij(ˆxij) denotes the index(es) of the subinterval(s) that ˆxij belongs to. If ˆxij is at the breakpoints, then ˆxij belongs to two subintervals, and Ψij(ˆxij) has two elements. For example, if ˆxij=λ2ij, then ˆxij belongs to both the second and the third subinterval, and thus Ψij(ˆxij)={2,3}.

    Definition 2. For any ˆxΩ, the arc set A can be partitioned into two subsets as follows:

    {(i,j)A1,if  (i,j)A,ˆxijBij,(i,j)A2,if  (i,j)A,ˆxijBij, (2.2)

    where Bij={λ1ij,λ2ij,,λsij1ij} consists of all breakpoints in [0,uij]. For any (i,j)A1, Ψij(ˆxij) has two elements. We call arcs in A1 active arcs and use m1 to represent the number of active arcs.

    As discussed, we can randomly select an element from Ψij(ˆxij),(i,j)A and restrict the flow xij to the corresponding subinterval to obtain a region that ˆx belongs to. Considering all possible selections, a set collection Φ(ˆx) can be used to indirectly index all the regions that ˆx belongs to as follows.

    Φ(ˆx)={zΘ|zijΨij(ˆxij),(i,j)A}.

    In above definition, "zΘ" is just used to show that any tuple φΦ(ˆx) has m elements that have a one-to-one correspondence with the arcs in A. In fact, φij is an integer in [1,sij] since φijΨij(ˆxij). The region associated with φ is determined by selecting the φij-th subinterval for any arc (i,j)A and can be denoted as

    Ωφ={xΩ|λφij1ijxijλφijij,(i,j)A}.

    Clearly, the number of regions that ˆx belongs to is 2m1 because Ψij(ˆxij) has two elements for any (i,j)A1. Since the feasible domain Ω is constrained by equality constraints, Ωφ might have only one point ˆx. In this case, when checking the local optimality, we do not need to check the optimality of ˆx in this region because ˆx is the only point in this region and is naturally optimal. However, in practice, we do not know in advance which region has only one point, so we have to treat all the regions equally and check the optimality of ˆx in this region. When Ωφ has only one point, we can easily know that ˆx is optimal in this region, and the process of checking the optimality of ˆx is equal to verifying that Ωφ has a unique point. Therefore, we only need to check the optimality of ˆx in all the regions that ˆx belongs to, and consequently, the region with only one point becomes irrelevant.

    In Ωφ, the cost function is linear and can be expressed as

    fφ(x)=(i,j)A(cφijijxij+dφijij).

    Due to the concavity of f(x), we can derive the following property that will be used in the proof of Theorem 1.

    Property 1. For any ˆxΩ, we have f(ˆx)fφ(ˆx).

    Proof. This property can be obtained directly: f(ˆx) can be expressed as f(ˆx)=(i,j)Amin{c1ijˆxij+d1ij,c2ijˆxij+d2ij,,csijijˆxij+dsijij}.

    Definition 3. The linear NFP associated with φ is defined as

     P(φ):  min  fφ(x)s.t.   xΩ. (2.3)

    Similarly, the linear NFPs associated with ˆxΩ are defined as P(φ),φΦ(ˆx). The number of linear NFPs associated with ˆx is 2m1.

    Now we use a CPLNFP instance to illustrate the above concepts, including Ψij(ˆxij), Φ(ˆx), φ and P(φ).

    Example 1. Consider a CPLNFP defined on a network with 6 nodes and 8 arcs, as shown in Figure 3. Following the formulation (1.1), this CPLNFP is formulated as

    Figure 3.  A practical CPLNFP.
    min f12(x12)+f13(x13)+f23(x23)+f24(x24)+f34(x34)+f35(x35)+f46(x46)+f56(x56)s.t.  x12+x13=7;x24+x23x12=0x34+x35x13x23=0;x46x23x34=0x56x35=0;x46x56=70x126;0x134;0x2330x244;0x345;0x3580x462;0x567

    where

    f12(x12)={5x12,       x12[0,3]3x12+6,  x12[3,6];  f13(x13)={6x13,       x13[0,2]4x13+4,  x13[2,4];
    f23(x23)={3x23,       x23[0,1]2x23+1,  x23[1,3]; f24(x24)={6x24,       x24[0,2]3x24+6,  x24[2,4];
    f34(x34)={4x34,       x34[0,3]3x34+3,  x34[3,5]; f35(x35)={4x35,        x35[0,5]2x35+10,  x35[5,8];
    f46(x46)={6x46,       x46[0,1]3x46+3,  x46[1,2]; f56(x56)={6x56,        x56[0,5]4x56+10,  x56[5,7].

    The corresponding tuple collection Θ can be expressed as

    Θ={(γ12,γ13,γ23,γ24,γ34,γ35,γ46,γ56)},

    and the feasible domain can be expressed as Ω={xΘ|x12+x13=6,x24+x23x12=0,x34+x35x13x23=0,x46x23x34=0,x56x35=0,x46+x56=6,0x126,0x134,0x233,0x244,0x345,0x358,0x462,0x567}.

    For a feasible solution ˆx=(3,4,1,2,0,5,2,5), we know Ψ12(ˆx12)={1,2}, Ψ13(ˆx13)={2}, Ψ23(ˆx23)={1,2}, Ψ24(ˆx24)={1,2}, Ψ34(ˆx34)={1}, Ψ35(ˆx35)={1,2}, Ψ46(ˆx46)={2}, Ψ56(ˆx56)={1,2}. Obviously, there are 5 active arcs: (1,2), (2,3), (2,4), (3,5) and (5,6). By definition, we have Φ(ˆx)={(z12,z13,z23,z24,z34,z35,z46,z56)|z12{1,2},z13{2},z23{1,2},z24{1,2},z34{1},z35{1,2},z46{2},z56{1,2}}, which consists of 25=32 tuples. For φ=(1,2,1,2,1,2,2,2)Φ(ˆx), the linear NFP associated with φ can be expressed as

    P(φ):min 5x12+4x13+4+3x23+3x24+6+4x34+2x35+10+3x46+3+4x56+10s.t.  xΩ

    Let V represent the set of vertices of Ω, or in other words, the set of the basic feasible solutions of the constraint system in (1.1). Based on Definition 3, we can derive the following local optimality condition.

    Theorem 1. A vertex vV is locally optimal for the CPLNFP if and only if v is optimal for all the linear NFPs associated with it.

    Proof. If there is a P(φ),φΦ(v) for which v is not optimal, then there must be a feasible descent direction dΘ satisfying v+εdΩ and fφ(v+εd)<fφ(v) for a sufficiently small positive number ε. Since fφ(v)=f(v), we will have f(v+εd)fφ(v+εd)<f(v), which contradicts that v is locally optimal for the CPLNFP.

    For a sufficiently small δRm, if v+δΩ, then v+δ must belong to the union of regions that v belongs to. Therefore, there must be a region Ωφ,φΦ(v) containing v+δ, from which it follows that f(v+δ)=fφ(v+δ). Since v is optimal for all the linear NFPs associated with it, we have f(v)=fφ(v)fφ(v+δ)=f(v+δ), which means v is locally optimal for the CPLNFP.

    In the worst case, checking the local optimality of v needs to consider the optimality of v for 2m1 linear NFPs, when one relies on Theorem 1. It is noted that the feasible domain of P(φ) is the same as that of the CPLNFP, not Ωφ. Thus, one does not need to worry about the change of the feasible domain, but only the change of the objective function when checking the optimality of v for these linear NFPs.

    The network simplex algorithm (NSA) (for detail see [16], [17] and[18]) is a classical and efficient algorithm to solve the linear NFP. When solving P(φ) by the NSA, the vertex vV as the basic feasible solution obtained in each iteration is associated with the following parameters:

    Basis: Mathematically, the basis is a triple (T,L,U), which is a partitioning of the arc set A as follows:

    {The arcs in T constitute a spanning tree of network G;(i,j)L, if  (i,j)T,vij=0;(i,j)U, if  (i,j)T,vij=uij. (2.4)

    The iteration of the NSA is essentially the iteration of the basis. Therefore, the basis is the core of iteration and others parameters are related to the vertex through the basis.

    For ease of expression, we set ¯T=LU. An arc in T is called a basis arc and an arc in ¯T is called a nonbasis arc. By definition, we have A1¯T=ϕ. Similar to the partitioning of A, T can also be partitioned into T1 and T2 as follows:

    {(i,j)T1, if (i,j)T,(i,j)A1;(i,j)T2, if (i,j)T,(i,j)A1. (2.5)

    By definition, we have T1=A1.

    Price: The price pi for node i is defined as

    pi=pr(h,t)Path(i)+Tcφhtht+(h,t)Path(i)Tcφhtht, (2.6)

    where pr is the price for the root node of the spanning tree, and Path(i)+T and Path(i)T are the arc sets consisting of the forward arcs and backward arcs in the path of the spanning tree from the root node to node i, respectively. Generally, the root node of the spanning tree and its price pr are specified in advance and will not change during iteration. In this paper, pr is set as zero.

    Reduced cost: The reduced cost rij for arc (i,j)A is defined as

    rij=cφijij+pjpi. (2.7)

    Note that, for any (i,j)T, the reduced cost rij is always zero by this definition.

    Optimality condition: A vertex vV is optimal for P(φ) if the reduced costs associated with v satisfy the following condition:

    {rij0,(i,j)L,rij0,(i,j)U. (2.8)

    In this section, we first derive the detailed local optimality condition by using the NSA to determine the optimality of the vertex for linear NFPs and then prove that the derived condition can be implemented with polynomial complexity in network parameters.

    Let vφ be the optimal vertex to P(φ). For vertex v, we first randomly select a φΦ(v) and then verify the optimality of v for P(φ). If v is not optimal for P(φ), we can immediately determine that v is not locally optimal for the CPLNFP. Therefore, we only need to focus on the case that v=vφ and propose a local optimality condition for vφ.

    Traditionally, we need to verify the optimality of vφ for the remaining 2m11 linear NFPs associated with it to determine the local optimality. The optimality of vφ for P(φ) can be verified through solving P(φ) by the NSA. In this paper, it is assumed that vφ is nondegenerate and is thus associated with only one spanning tree. Here, the statement that vφ is nondegenerate means that the basic feasible solution corresponding to vφ is nondegenerate. We say that vφ is nondegenerate from a geometrically intuitive perspective. Similarly to the discussion on simplex algorithm, dealing with degeneracy is skilful and requires complicated representations. In this paper, we focus on nondegenerate case and leave the discussion on degeneracy for further study. When solving the remaining linear NFPs associated with vφ, we can use vφ as the initial solution, and thus the parameters associated with vφ can be utilized. For example, the basis associated with vφ will not change when solving the linear NFPs associated with vφ.

    For the sake of distinction, in the following discussion, we use rφij and rφij to represent the reduced cost associated with vφ in solving P(φ) and P(φ), respectively. Based on the optimality condition (2.8) for the linear NFP, the local optimality condition for vφ can be expressed as follows.

    Theorem 2. The vertex vφ is locally optimal for CPLNFP if and only if rφij,(i,j)¯T satisfy condition (2.8) for any φΦ(vφ).

    Traditionally, for each of the linear NFPs associated with vφ, we need to calculate the reduced costs to determine the optimality of vφ. In other words, we have to calculate the reduced cost for each nonbasis arc 2m1 times.

    In this paper, we consider the optimality of vφ for the linear NFPs associated with it from the point of view of each nonbasis arc rather than each linear NFP. According to Theorem 2, for each (i,j)¯T, we actually only need to check whether the extremum of the reduced cost when solving the linear NFPs associated with vφ satisfies condition (2.8). Here, the extremum refers to the minimum for (i,j)L and the maximum for (i,j)U. If the extremum of the reduced cost for arc (i,j) satisfies (2.8), we can conclude that the reduced cost for arc (i,j) always satisfies (2.8) when solving all the linear NFPs associated with vφ. Based on this idea, we actually only need to calculate the extremum of the reduced cost for each nonbasis arc, which means we only need to calculate the reduced cost for each nonbasis arc once. To better represent the extremum, we define ψ(ij) as follows.

    ψ(ij)={argminφ{rφij|φΦ(vφ)},if (i,j)L;argmaxφ{rφij|φΦ(vφ)},if (i,j)U. (3.1)

    With this definition, rψ(ij)ij represents the extremum of the reduced cost for arc (i,j) when solving the linear NFPs associated with vφ and P(ψ(ij)) is the linear NFP that makes the reduced cost reach the extremum. Notice that although ψ(ij) is defined on an exponential-sized set, but we do not need actually calculate it in that way. In the next subsection, we will provide a strategy to obtain ψ(ij) and rψ(ij)ij in complexity O(n).

    Theorem 3. The local optimality condition shown in Theorem 2 can be checked with polynomial complexity.

    Proof. Since rψ(ij)ij can be obtained with complexity O(n) for any nonbasis arc and the number of nonbasis arcs is mn, the total complexity for checking the condition shown in Theorem 2 is O(n(mn)), which is polynomial in the network parameters.

    One key procedure of checking local optimality is to obtain ψ(ij) and rψ(ij)ij for any nonbasis arc (i,j) with complexity O(n), which is represented as the following theorem.

    Theorem 4. For vertex vφ, ψ(ij) and rψ(ij)ij can be obtained with complexity O(n) for any (i,j)¯T.

    Proof. This is a constructive proof, in which a strategy to obtain ψ(ij) for any nonbasis arc (i,j) with complexity O(n) is given. This strategy is based on the difference in the reduced cost between solving P(φ) and the remaining linear NFPs associated with vφ.

    For any of the remaining linear NFPs: P(φ) (φΦ(vφ), φφ), we need to calculate rφij to check the optimality of vφ. Since rφij has been obtained when solving P(φ), we only need to consider the difference between rφij and rφij. Moreover, for any nobasis arc (i,j), vφij equals 0 or uij, which means that vφij belongs to only one subinterval (i.e., the first subinterval or the last subinterval), and thus we have φij=φij and cφijij=cφijij. Therefore, based on the definition of the reduced cost, the difference between rφij and rφij can be expressed as

    Δφφrij=ΔφφpjΔφφpi,(i,j)¯T, (3.2)

    where Δφφpi and Δφφpj are the increments of pi and pj when the linear NFP changes from P(φ) to P(φ), respectively.

    By definition, the price is related to the arcs in the spanning tree. So we only need to consider the effect of basis arcs on Δφφpi and Δφφpj. For ease of expression, we use Δφhtφhtrij, Δφhtφhtpi and Δφhtφhtpj to represent the effects of arc (h,t)T on Δφφrij, Δφφpi and Δφφpj, respectively. Then, we have

    Δφφrij=(h,t)TΔφhtφhtrij=(h,t)T(ΔφhtφhtpjΔφhtφhtpi) (3.3)

    Moreover, for any basis arc (h,t), if (h,t)T1, we have φht=φht according to the definition of T1, which means arc (h,t) has no effect on Δφφpi and Δφφpj. Therefore, when calculating Δφφrij, we need to focus on the arcs in T1.

    Definition 4. For any (i,j)¯T, the unique cycle formed by arc (i,j) and the spanning tree is called the derived cycle and is denoted as CTij.

    By this definition, CTij can be expressed as

    CTij=(i,j)Path(i)TPath(j)T.

    where Path(i)T and Path(j)T are the arc sets consisting of the arcs in the path of the spanning tree from the root node to node i and j, respectively.

    Lemma 1. For any (h,t)T1 and (i,j)¯T, if (h,t)CTij, then arc (h,t) has no effect on Δφφrij.

    Proof. If (h,t)CTij, then there are three cases for any (h,t)T1.

    Case 1: (h,t)Path(i)+T, (h,t)Path(j)+T;

    Case 2: (h,t)Path(i)T, (h,t)Path(j)T;

    Case 3: (h,t)Path(i)+TPath(i)T, (h,t)Path(j)+TPath(j)T.

    Clearly, without considering other arcs, we have Δφhtφhtpi=Δφhtφhtpj for all three cases based on the definition of price in (2.6). Therefore, we have Δφhtφhtrij=0 which means that arc (h,t) has no effect on Δφφrij.

    Therefore, only (h,t)T1CTij contributes to Δφφrij. Let GTij=T1CTij; we have

    Δφφrij=(h,t)GTijΔφhtφhtrij=(h,t)GTij(ΔφhtφhtpjΔφhtφhtpi) (3.4)

    According to the location and direction, the arcs in GTij can be divided into the following four subsets:

    GT1ij={(h,t)GTij | (h,t) belongs to the path of the spanning tree from the root node to node j, and it has the reverse direction as (i,j) in the cycle CTij}.

    GT2ij={(h,t)GTij | (h,t) belongs to the path of the spanning tree from the root node to node j, and it has the same direction as (i,j) in the cycle CTij}.

    GT3ij={(h,t)GTij | (h,t) belongs to the path of the spanning tree from the root node to node i, and it has the same direction as (i,j) in the cycle CTij}.

    GT4ij={(h,t)GTij | (h,t) belongs to the path of the spanning tree from the root node to node i, and it has the reverse direction as (i,j) in the cycle CTij}.

    Now, for any arc (h,t) in these four subsets, we use an example to illustrate the effect of arc (h,t) on Δφφrij.

    Example 2. Examples of these four subsets are shown in Figure 4. In the derived cycle CTij, it is assumed that (e,f),(g,f),(b,c),(d,c)T1; then, we have (e,f)GT1ij, (g,f)GT2ij, (b,c)GT3ij, and (d,c)GT4ij.

    Figure 4.  Four cases of arc (h,t)GTij.

    For any (h,t)GTij, if vφhtht=λkht, then φht{k,k+1}. Regardless of the influence of other arcs, Δφhtφhtrij can be calculated as follows.

    (1) Case 1: As shown in Figure 4(a), (h,t)=(e,f)GT1ij. Based on (2.6), we have Δφhtφhtpi=0 and Δφhtφhtpj=cφijhtcφijht. Therefore, according to (3.2), Δφhtφhtrij can be calculated by Δφhtφhtrij=cφijhtcφijht.

    (2) Case 2: As shown in Figure 4(b), (h,t)=(g,f)GT2ij. Based on (2.6), we have Δφhtφhtpi=0 and Δφhtφhtpj=cφijhtcφijht. Therefore, according to (3.2), Δφhtφhtrij can be calculated by Δφhtφhtrij=cφijhtcφijht.

    (3) Case 3: As shown in Figure 4(c), (h,t)=(b,c)GT3ij. Based on (2.6), we have Δφhtφhtpi=cφijhtcφijht and Δφhtφhtpj=0. Therefore, according to (3.2), Δφhtφhtrij can be calculated by Δφhtφhtrij=cφijhtcφijht.

    (4) Case 4: As shown in Figure 4(d), (h,t)=(d,c)GT4ij. Based on (2.6), we have Δφhtφhtpi=cφijhtcφijht and Δφhtφhtpj=0. Therefore, according to (3.2), Δφhtφhtrij can be calculated by Δφhtφhtrij=cφijhtcφijht.

    This example illustrates the effect of a single arc on Δφφrij. Considering the effect of all the arcs in GTij, we can calculate Δφφrij by

    Δφφrij=(h,t)GT1ij(cφhthtcφhtht)+(h,t)GT2ij(cφhthtcφhtht)+(h,t)GT3ij(cφijhtcφhtht)+(h,t)GT4ij(cφhthtcφhtht). (3.5)

    Based on the concavity of fij(xij), we have ckij>ck+1ij. For any arc (i,j)L, since ψ(ij)=argminφ{rφij|φΦ(vφ)}, the value of ψ(ij)ht for any (h,t)A can be obtained as follows.

    ψ(ij)ht={k,if (h,t)GT1ijGT4ij;k+1,if (h,t)GT2ijGT3ij;φht,else. (3.6)

    Accordingly, the difference between rψ(ij)ij and rφij can be calculated as

    Δφψ(ij)rij=(h,t)GT1ijGT4ij(cφhthtckht)+(h,t)GT2ijGT3ij(ck+1htcφhtht). (3.7)

    Similarly, for any arc (i,j)U, since ψ(ij) = argmaxφ{rφij|φΦ(vφ)}, the value of ψ(ij)ht for any (h,t)A can be obtained as follows.

    ψ(ij)ht={k+1,if (h,t)GT1ijGT4ij;k,if (h,t)GT2ijGT3ij;φht,else. (3.8)

    Accordingly, the difference between rψ(ij)ij and rφij can be calculated as

    Δφψ(ij)rij=(h,t)GT1ijGT4ij(cφhthtck+1ht)+(h,t)GT2ijGT3ij(ckhtcφhtht). (3.9)

    In summary, for any nonbasis arc (i,j), if we have already obtained CTij, then it is clear that the computational complexity of obtaining ψ(ij) and Δφψ(ij)rij by (3.6) and (3.7) (or (3.8) and (3.9)) is O(n). Based on the NSA, the computational complexity of searching for a derived cycle CTij is also O(n) through backtracking the spanning tree from node i and node j to the root node. Then, rψ(ij)ij can be obtained by rψ(ij)ij=rφij+Δφψ(ij)rij, which requires only an addition operation. Therefore, the total complexity of obtaining ψ(ij) and Δrψ(ij)ij is still O(n).

    At this point, the constructive proof of Theorem 4 is completed.

    Using the strategy in the proof of Theorem 4 to obtain ψ(ij), the procedure for checking the local optimality of a vertex v can be expressed as Algorithm 1.

    Algorithm 1 Check the local optimality of v.
    Require: Feasible domain Ω; fij(xij),(i,j)A, v, flag=1;
    1: Select an index tuple φΦ(v);
    2: Solve the linear NFP P(φ) and obtain rφij;
    3: if v is not optimal for P(φ) then
    4:    flag:=0;
    5: else
    6:   for (i,j)¯T do
    7:      Obtain ψ(ij) by (3.6) or (3.8);
    8:      Calculate Δφψ(ij)rij by (3.7) or (3.9);
    9:      Calculate rψ(ij)ij by rψ(ij)ij=rφij+Δφψ(ij)rij;
    10:     if rψ(ij)ij satisfies condition (2.8) then
    11:        Go to the next arc in ¯T.
    12:     else
    13:       flag:=0;
    14:        End the current for-loop.
    15:     end if
    16:   end for
    17: end if
    18: if flag == 1 then
    19:    v is locally optimal for CPLNFP;
    20: else
    21:    v is not locally optimal for CPLNFP;
    22: end if

    Theorem 5. The complexity of checking the local optimality of v by Algorithm 1 is polynomial in the network parameters.

    Proof. If v is not optimal for P(φ), we determine that v is not locally optimal for CPLNFP by verifying the optimality of v for only one linear NFP. In this case, the computational complexity of verifying the optimality of v for P(φ) is O(mn) according to the NSA. If v is optimal for P(φ) (i.e., v=vφ), we check the optimality of v(vφ) by the strategy stated in the proof of Theorem 4 in complexity O(n(mn)). Therefore, the total complexity of checking the local optimality of v by Algorithm 1 is still O(n(mn)), which is polynomial in the network parameters.

    Traditionally, as mentioned earlier, for any vertex v, we need to verify the optimality of v for 2m1 linear NFPs P(φ),φΦ(v). The total computational complexity is O(2m1(mn)) since the complexity of checking v for a linear NFP is O(mn). It is clear that for large-scale CPLNFPs with 2m1>>n, Algorithm 1 is far superior to the traditional implementation of the local optimality condition.

    Example 3. For the CPLNFP instance in Example 1, we can easily obtain that ˆx=(3,4,1,2,0,5,2,5) is a nondegenerate vertex. Moreover, ˆx is also an optimal solution to P(φ) for φ=(1,2,1,2,1,2,2,2), which means that vφ=ˆx=(3,4,1,2,0,5,2,5). When solving P(φ) by the NSA, the basis associated with vφ is

    {T={(1,2),(2,3),(2,4),(3,5),(5,6)},L={(3,4)},U={(1,3),(4,6)}. (3.10)

    Figure 5 shows the basis, where the basis arcs are represented by solid lines and nonbasis arcs are represented by dashed lines. By definition, we have ¯T=LU={(1,3),(3,4),(4,6)} and T1={(1,2),(2,3),(2,4),(3,5),(5,6)}. Let node 1 be the root of T. The prices associated with vφ are p1=0,p2=5,p3=8,p4=8,p5=10,p6=14 when solving P(φ) by the NSA. Accordingly, the reduced costs for nonbasis arcs are rφ13=4+(8)0=4,rφ34=4+(8)(8)=4,rφ46=3+(14)(8)=3. It can be seen that rφ13, rφ34 and rφ46 satisfy the optimality condition (2.8), which in turn verifies the optimality of ˆx for P(φ).

    Figure 5.  The basis associated with vφ.

    For arc (1,3)U, the derived cycle is CT13={(1,2),(2,3),(1,3)} and the intersection of CT13 and T1 is GT13={(1,2),(2,3)}. As discussed in section 3.2, GT13 can be divided into four subsets: GT113={(1,2),(2,3)}, GT213=ϕ, GT313=ϕ and GT413=ϕ. According to (3.6) and (3.7), we have ψ(13)=(2,2,2,2,1,2,2,2) and Δφψ(13)r13=(53)+(32)=3. Then, rψ(13)13 can be calculated as rψ(13)13=rφ13+Δφψ(13)r13=4+3=1<0. Clearly, rψ(13)13 still satisfies condition (2.8).

    Similarly, for arc (3,4)L and arc (4,6)U, we can obtain ψ(34)=(1,2,2,1,1,2,2,2),Δφψ(34)r34=(23)+(36)=4,rψ(34)34=44=0 and ψ(46)=(1,2,1,2,1,1,2,1),Δφψ(46)r46=(33)+(42)+(64)=4,rψ(46)46=3+4=1>0. Obviously, rψ(46)46 does not satisfy condition (2.8). Therefore, vφ is not locally optimal for the CPLNFP instance in Example 1. In practice, when rψ(46)46 is obtained, one can determine that vφ is not locally optimal and obtain a better solution by solving P(ψ(46)). In the process of checking the local optimality of vφ, we calculate the reduced cost for each nonbasis arc only once.

    Traditionally, we need to check the optimality of vφ for 32 linear NFPs. Moreover, we need to calculate the reduced cost for each nonbasis arc every time we solve one of these 32 linear NFPs. For example, when solving P(φ) where φ=(1,2,2,2,1,1,2,2)Φ(ˆx), we need to recalculate prices as: p1=0,p2=5,p3=7,p4=8,p5=11,p6=15, and then check the optimality of vφ by calculating the reduced costs for nonbasis arcs as rφ13=4+(7)0=3,rφ34=4+(8)(7)=3,rφ46=3+(15)(8)=4. Therefore, if we want to check the optimality of vφ for 32 linear NFPs, we need to calculate the reduced cost for each nonbasis arc 32 times.

    The motivation of our study comes from the demand for the local optimal solution to the CPLNFP. The efficient implementation of the local optimality condition is the key to verify and obtain a local optimal solution. In this paper, we propose a polynomial local optimality condition for the nondegenerate vertex based on the mathematical properties of linear NFPs associated with this vertex. The proposed local optimality condition makes use of the network characteristics of the CPLNFP and is suitable for large-scale CPLNFPs.

    This project was supported by the National Natural Science Foundation of China (61473165, 61977046).

    All authors declare no conflicts of interest in this paper.

    This appendix estimates the growth trend of the number of active arcs (m1) in the CPLNFPs defined on benchmark networks.

    Test problems used to estimate the growth trend of m1 are generated randomly. Each test problem consists of two parts: the network and the concave piecewise linear cost function. The network is generated randomly with parameters n and m by the benchmark network generator NETGEN [19]. The total supply flow is equal to the total demand flow and is set as 10×n. The capacity uij for arc (i,j) is generated along with G, as stated in [19]. Note that parameters n and m are correlated; thus, we use n/m to denote these two parameters. Without loss of generality, it is assumed that sij remains the same over all arcs and is expressed as ˆs. The concave piecewise linear cost function fij(xij) for each arc (i,j)A is generated according to the following three steps.

    Step 1: The breakpoints (λ1ij,λ2ij,,λˆs1ij) are randomly generated in (0, uij), which satisfies λ1ij<λ2ij<<λˆs1ij.

    Step 2: The slopes (c1ij,c2ij,,cˆsij) are randomly generated in [cmin, cmax], which satisfies c1ij>c2ij>>cˆsij. In this paper, we set cmin=5 and cmax=5ˆs.

    Step 3: The intercepts (d1ij,d2ij,,dˆsij) are generated iteratively based on the continuity of fij(xij) as follows:

    d1ij=0;

    For s=2 to ˆs do

    * dsij=cs1ijλs1ij+ds1ijcsijλs1ij.

    Considering that various practical problems require various network dimension parameters n/m and different piecewise linear cost function parameters ˆs, we set ten different values for parameter n/m and four different values for parameter ˆs, which can be expressed as as Ξ = {10/20 20/40 40/100 60/400 80/600 100/1000 200/2000 300/3000 400/6000 500/8000} and S={3,5,7,10}, respectively. Problems with the same n/m and ˆs form a problem set, which is denoted as Pn,m,ˆs. In each problem set Pn,m,ˆs, there are 20 randomly generated test problems.

    In this paper, the problem group that consists of all test problems is denoted as P. According to parameter ˆs, P is divided into four problem subgroups P1, P2, P3 and P4.

    P1={Pn,m,3|n/mΞ},P2={Pn,m,5|n/mΞ},P3={Pn,m,7|n/mΞ},P4={Pn,m,10|n/mΞ}.

    Clearly, there are 200 test problems in Pk,k=1,2,3,4 and a total of 800 test problems in P.

    To study the relationship between m1 and the parameters n, m and ˆs, we estimate the value of m1 in problem groups P1, P2, P3 and P4. The test problems in each problem set Pn,m,ˆs have the same parameters n, m and ˆs; thus, they are estimated as a whole. The average of the estimated values of m1 for all 20 problems is used to estimate the overall level of m1 for Pn,m,ˆs. For a test problem, each vertex in the feasible region corresponds to a value of m1. Since the number of vertices is combinatorial, it is not realistic to find values of m1 for all vertices. Here, we propose a reasonable sampling method to estimate the value of m1.

    For ease of expression, we use the concept of P(φ) that is a linear NFP defined from one of the regions in the feasible domain. The formal definition of P(φ) will be given in (2.3). The linear NFPs defined from all the possible regions in the feasible domain can be expressed as

    ˆP={P(φ)|φˆΦ},

    where ˆΦ is defined as

    ˆΦ={φΘ|φij{1,2,,ˆs},(i,j)A}.

    The definition of Θ can be seen in section 2.1.

    For a vertex v in the feasible domain, if v is not optimal for P(φ),φˆP, then we can determine v is not locally optimal for the CPLNFP by solving any one of these linear NFPs. Therefore, we do not have to explore the value of m1 for this vertex. We should focus on vertices that are at least optimal for one linear NFP in ˆP. It is appropriate to sample vertices from the optimal solutions to problems in ˆP.

    For a test problem pP with the feasible region Ω, the process of sampling ns vertices is shown as follows.

    For k=1:ns do

    - Randomly select φ satisfying φˆΦ;

    - Solve P(φ) by the network simplex algorithm ([16]) and obtain the optimal vertex vφ;

    - Calculate the value of m1 associated with vφ;

    - γk=m1;

    - ˆΦ:=ˆΦφ;

    End

    The average of γ1,γ2,,γns is used to estimate the value of m1 in this problem.

    ηp=γ1+γ2++γnsns. (A.1)

    In our experiments, we set ns=n. The overall level of m1 for Pn,m,ˆs is estimated by

    ζ=pPn,m,ˆsηp20. (A.2)

    It can be seen in Table 1, as the scale of the problem increases, the value of ζ gradually increases and the value of 2ζ becomes very large. For example, we have 2ζ1.71×1047>>500 for test problems in P500,8000,10. Moreover, for the same problem scale, the value of ζ also increases as the value of ˆs increases. The specific growth trend of ζ can be seen in Figure 2 in Introduction.

    Table 1.  The overall level of the value of m1.
    Problemn1020406080100200300400500
    Groupm204010040060010002000300060008000
    P1(ˆs=3)ζ0.150.340.801.241.782.253.855.408.1210.50
    P2(ˆs=5)ζ0.320.982.032.423.554.357.7511.2014.9818.31
    P3(ˆs=7)ζ0.601.203.123.965.566.4012.4217.0521.8828.91
    P4(ˆs=10)ζ0.872.143.945.947.229.9318.1526.8433.9141.58

     | Show Table
    DownLoad: CSV


    [1] J. Xu, Y. Tu, Z. Zeng, A nonlinear multiobjective bilevel model for minimum cost network flow problem in a large-scale construction project, Math. Probl. Eng., 2012 (2012), 115–123.
    [2] V. Prahalad, K. Mathur, R. H. Ballou, An efficient generalized network-simplex-based algorithm for manufacturing network flows, J. Comb. Optim., 15 (2008), 315–341. doi: 10.1007/s10878-007-9080-6
    [3] F. S. Salman, R. Ravi, J. N. Hooker, Solving the Capacitated Local Access Network Design Problem, Informs J. Comput., 20 (2008), 243–254.
    [4] C. D'Ambrosioa, A. Lodi, S. Wiese, et al., Mathematical programming techniques in water network optimization, Eur. J. Oper. Res., 243 (2015), 774–788.
    [5] T. L. Magnanti, D. Stratila, Separable concave optimization approximately equals piecewise linear optimization, In: Integer Programming and Combinatorial Optimization, Springer, Berlin, 2004.
    [6] B. W. Lamar, A method for solving network flow problems with general non-linear arc costs, In: Network Optimization Problems: Algorithms Applications and Complexity, World Scientific, Singapore, 1993.
    [7] A. Saif, S. Elhedhli, A Lagrangian heuristic for concave cost facility location problems: the plant location and technology acquisition problem, Optim. Lett., 10 (2016), 1087–1100. doi: 10.1007/s11590-016-0998-4
    [8] Y. Perlman, I. Levner, Perishable Inventory Management in Healthcare, Journal of Service Science and Management, 7 (2014), 11–17. doi: 10.4236/jssm.2014.71002
    [9] S. Yan, Y. L. Shih, W. T. Lee, A particle swarm optimization-based hybrid algorithm for minimum concave cost network flow problems, J. Global Optim., 49 (2011), 539–559. doi: 10.1007/s10898-010-9548-2
    [10] Z. Xu, K. Liu, X. Xi, S. Wang, Method of hill tunneling via weighted simplex centroid for continuous piecewise linear programming, 2015 Conference on Decision and Control, 24 (2015), 301–316.
    [11] D. Kim, P. M. Pardalos, Dynamic slope scaling and trust interval techniques for solving concave piecewise linear network flow problems, Networks, 35 (2000), 216–222. doi: 10.1002/(SICI)1097-0037(200005)35:3<216::AID-NET5>3.0.CO;2-E
    [12] A. Nahapetyan, P. M. Pardalos, A bilinear relaxation based algorithm for Concave Piecewise Linear Networks Flow Problems, J. Ind. Manag. Optim., 3 (2007), 71–85.
    [13] A. R. Conn, M. Mongeau, Discontinuous piecewise linear optimization, Math. Program., 80 (1998), 315–380.
    [14] L. T. H. An, P. D. Tao, The DC (Difference of Convex Functions) Programming and DCA Revisited with DC Models of Real World Nonconvex Optimization Problems, Ann. Oper. Res., 133 (2005), 23–46. doi: 10.1007/s10479-004-5022-1
    [15] X. Huang, J. Xu, X. Mu, S. Wang, The hill detouring method for minimizing hinging hyperplanes functions, Comput. Oper. Res., 39 (2012), 1763–1770. doi: 10.1016/j.cor.2011.10.017
    [16] G. B. Dantzig, Linear Programming and Extensions, Journal of Industrial and Management Optimization, Princeton University Press, NJ, 1963.
    [17] P. Kovacs, Minimum-cost flow algorithms: an experimental evaluation, Optim. Methods Softw., 30 (2015), 94–127. doi: 10.1080/10556788.2014.895828
    [18] M. Holzhauser, S. O. Krumke, C. Thielen, A Network Simplex Method for the Budget-Constrained Minimum Cost Flow Problem, Eur. J. Oper. Res., 259 (2017), 864–872. doi: 10.1016/j.ejor.2016.11.024
    [19] D. Klingman, A. Napier, J. Stutz, NETGEN: A program for generating large scale capacitated assignment, transportation, and minimum cost flow network problems, Manage. Sci., 20 (1973), 814–821.
  • Reader Comments
  • © 2021 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(1444) PDF downloads(24) Cited by(0)

Figures and Tables

Figures(5)  /  Tables(1)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog