Maximizing profit in a supply chain by considering advertising and price elasticity of demand

Computers & Industrial Engineering 135 (2019) 265–274

Contents lists available at ScienceDirect

Computers & Industrial Engineering journal homepage: www.elsevier.com/locate/caie

Maximizing profit in a supply chain by considering advertising and price elasticity of demand☆

T

Kamran Kianfar Faculty of Engineering, University of Isfahan, Isfahan 81746-73441, Iran

ARTICLE INFO

ABSTRACT

Keywords: Supply chain management Price elasticity of demand Advertising elasticity Harmony search Mixed-integer quadratically constrained quadratic programming

A three-level supply chain, including manufacturers, distribution centers and customer zones is assessed in this article. The objective is to maximize the gross profit without violating the operational constraints and capacities of production and inventory. Gross profit is achieved by subtracting costs from the final revenue. Different sources of cost related to raw materials, transportation, production and advertising are addressed. The final demand for a product depends on its price and advertising expenditures. Three solution methods of a mathematical model, a harmony search algorithm and a new combined algorithm of both are proposed for the problem. The mathematical model is of type mixed integer quadratically constrained quadratic programming and its characteristics are analyzed. The experimental results reveal that, statistically speaking, the proposed heuristic algorithms converge into optimal solution with gaps of less than 2 percent. A comprehensive sensitivity analysis is run on price and advertising elasticity coefficients, manufacturers and distribution centers' capacity, base demand and unit transportation cost.

1. Introduction All activities beginning from the procurement of raw materials to the finished products delivery, constitute the supply chain. These activities include procurement, production, distribution and sales. The process of managing all these activities is known as supply chain management. The advances made through the non-stop growth in technology, market competitions and rapid growth in economy, force trade organizations to focus on supply chain management. Several factors, like price and advertising, affect the demands in a supply chain. Price elasticity of demand is a measure applied to determine the responsiveness, or elasticity, of the quantity of goods and services in demand against a change in their price of the same. Advertising is one of the main promotional tools used in modern marketing management. Through advertising, manufacturers can influence potential customers’ interest and raise brand awareness and retailers can stimulate the customers' buying behavior. Advertising elasticity coefficient is a tool used in measuring the effect of an increase or a decrease in advertising in a market. This paper studies a three-level supply chain including manufacturers, distribution centers and customer zones. Different types of raw materials and final products are involved in this problem and the objective is to maximize the final profit. The capacities of

☆

manufacturers and distribution centers are limited. Fixed and variable costs of production, transportation costs, raw materials and advertising costs are considered. It is assumed that the final demand depends on the products’ price and advertising expenditures. Three solution methods are proposed for the aforesaid problem: (1) a mixed integer quadratically constrained quadratic programming model, the experimental results of which indicate that it solves the problem up to medium-size instances, (2) a solution method based on harmony search algorithm, one of the most recent meta-heuristic methods inspired by musical orchestra, and (3) a combination of both. The main contribution of this study is threefold: (1) for the first time, the effects of advertising and price elasticity of demand are studied in a three-level supply chain with respect to almost all types of costs, (2) two heuristic methods are developed for the problem and the experimental tests demonstrate their high performance in getting close to near-optimal solutions in reasonable run times. To the knowledge of the authors, combining mathematical models and the harmony search algorithm is a new idea in solving supply chain problems, and (3) the applied sensitivity analysis describes the effects of altering elasticity coefficients, capacity limitations and potential market size in a supply chain with variable demands. Meta-heuristic methods are widely adopted in solving supply chain models. Harmony search is an emerging meta-heuristic algorithm,

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors. E-mail address: [email protected].

https://doi.org/10.1016/j.cie.2019.06.007 Received 17 February 2018; Received in revised form 3 June 2019; Accepted 5 June 2019 Available online 06 June 2019 0360-8352/ © 2019 Elsevier Ltd. All rights reserved.

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

Table 1 Summary of researches on supply chain with variable demand. Authors (year)

supply chain members

Sources of cost

Chen (2014) Chen (2015) Giri and Sharma (2014) Hong et al. (2015) Hull (2005) Kaplan et al. (2011) Karray and Amin (2014) Li et al. (2016) Li et al. (2013) Ma et al. (2013) Mokhlesian and Zegordi (2014) Seifbarghy et al. (2015) SeyedEsfahani et al. (2011) Szmerekovsky and Zhang (2009) Xie and Neyret (2009) Xie and Wei (2009) Yue et al. (2006) Zhang et al. (2013)

Ma, Re Ma, Re, Cu Ma, Re Ma, Re Ma, Cu Su, Ma, Wa, Dc, Re Ma, Re Ma, Re Ma, Re Ma, Re, Cu Ma, Re Ma, Re Ma, Re Ma, Re Ma, Re Ma, Re Ma, Re Ma, Re

Pu, Ls, La Pu, Op, Ad Sa, Ad Ma, Re, Ad Gc Ma, Rm, Ho, Tr Ad Ma Ma, Ret Qu, Mar, Ma Ma, Ho Ma, Or, Qu Ma, Ha, Pu Ad, Ma Ad, Ma, Pu Ad, Ma, Ha, Pu Ad, Ma Ad

This paper

Ma, Dc, Cu

Tr, Rm, Ma, Ad

Constraints

Pc, Be Oc, Mf Gp Pc

Pc, Dc

Demand factors

Pr, Pr, Ad Pr, Pr Pr Pr, Pr, Pr, Pr, Pr Pr, Pr, Pr, Pr, Pr, Pr, Pr,

Ad Ad Ad Ad Gd Ad Qu, Mar Qu Ad Ad Ad Ad Ad Ad

Pr, Ad

Advertising expenditure

solution method

Local

Cooperative

Exact

* * *

Gt Gt Gt Gt Gt Mm Gt Gt Gt Gt Mm Gt Gt Gt Gt Gt Gt Gt

Global

* *

* * *

*

* *

* *

* * *

* * *

*

*

*

Mm

Meta-heuristic

Hga

Hs

Supply chain members | Cu: customer | Dc: distribution center | Ma: manufacturer | Re: retailer | Su: supplier | Wa: warehouse|. Source of cost | Ad: advertising | Gc: general cost | Ha: handling | Ho: holding | La: local advertising | Ls: lost sale | Ma: manufacturing | Mar: marketing | Op: operating | Or: ordering | Pu: purchasing | Qu: quality | Re: remanufacturing | Ret: retail | Rm: raw material | Sa: sales | Tr: transportation |. Constraints | Be: bullwhip effect | Dc: distribution centers capacity | Gp: green production | Mf: material flow | Oc: operational constraints | Pc: production capacity |. Demand factors | Ad: advertising | Gd: green degree | Mar: marketing | Pr: price | Qu: quality |. Solution method | Gt: game theory | Hga: hybrid genetic algorithm | Hs: harmony search | Mm: mathematical modelling |.

supply chain performance and organizations are analyzed. The combination of green marketing and supply chain decisions is a challenging topic for researchers (see Liu, Kasturiratne, and Moizer (2012) and Brindley and Oxborrow (2014)). Some researchers have applied the elasticity concept in designing supply chains. Giri and Sharma (2014) study a supply chain with one manufacturer and two competing retailers based on advertising-dependent demand. Hull (2005) develops a model describing the performance of supply chains based on their elasticity of supply and demand. Seifbarghy, Nouhi, and Mahmoudi (2015) study a two-level supply chain including a manufacturer and a retailer where the customer demand is subject to price and quality of final product. Kaplan, Türkay, Karasözen, and Biegler (2011) apply the price elasticity of demand concept and examine the price-centric behavior of customers. In an article by Li, Wang, and Yan (2013), pricing, ordering and advertising coordination in a supply chain are assessed where retailers face stochastic demand depending on retailer's price and advertising expenditures. Ma, Wang, and Shang (2013) study quality and marketing effort-dependent demands in a two-stage supply chain. Table 1 summarizes some studies about supply chain problems with variable demand. From this table, it is obvious that our paper is the first study that considers variable demands (price and advertising-dependent) in a three-level supply chain with some different sources of cost. Also, the proposed solution method is novel in this area of research. The rest of this article is organized as follows: the problem description is explained in Section 2. A mixed integer quadratically constrained quadratic program (MIQCQP) is provided in Section 3. Some theorems about near-optimal solutions of the problem and two new heuristic algorithms named HS and HS-QCP are developed in Section 4. In Section 5, the performance of the proposed methods is examined through numerical experiments and some sensitivity analysis will check the effect of changing the main parameters of the supply chain on costs and profits. Finally, the concluding remarks as well as some directions for future research are given in Section 6.

applied in handling many optimization problems in recent years. Banyai (2011) develops a harmony search approach for optimizing processes in green supply chain. This approach includes revenue, warehousing, transportation, recycling and disposal costs. Alaei and Khoshalhan (2015) study a one-buyer, multi-vendor, coordination model in a multi-objective centralized supply chain. They propose four meta-heuristic algorithms: (1) particle swarm optimization, (2) scatter search, (3) population based harmony search, and (4) harmony search based on cultural algorithm. A multiple-buyer multiple-vendor supply chain with stochastic demands is studied by Taleizadeh, Akhavan Niaki, and Barzinpour (2011), where a harmony search algorithm is applied to solve the integer nonlinear programming model of the problem. Another study on harmony search application in supply chain management is run by Wong and Guo (2010), where the sales forecasting in fashion retail supply chains is studied. Recently, the concepts of marketing and supply chains are combined in many studies. Chen (2015) studies the effects of pricing and cooperative advertising policies in a two-echelon dual-channel supply chain. The objective is to assist managers to understand the interaction between the upstream and downstream entities of a dual channel structure. Aust and Buscher (2014) review some cooperative advertising models in supply chain management similar to the article presented by Karray and Amin (2014). Other studies on application of advertising in supply chains are run by Hong, Xu, Du, and Wang (2015), Kozlenkova, Hult, Lund, Mena, and Kekec (2015), Zhang, Gou, Liang, and Huang (2013), Chen (2014) and Aust (2015). Mokhlesian and Zegordi (2014) apply a nonlinear multidivisional bi-level programming model for the problem of coordinating inventory and pricing decisions in a competitive supply chain. The study by Guillén, Badell, Espuña, and Puigjaner (2006) involves the financial decisions in an integrated planning of chemical supply chains. Pricing policies are studied by Li, Zhu, Jiang, and Li (2016) in a competitive dual-channel green supply chain in both centralized and decentralized cases. Aligning marketing strategies in supply chains is assessed by Green, Whitten, and Inman (2012), where the correlations among marketing strategy alignment, 266

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

2. Problem description

in this formulation are:

• Final demand is deterministic, but it is a function of selling price and advertising expenditures. • All customer demands must be met. • The time horizon consists of a single period. • Inventory costs are not included in the formulation. • Each plant can manufacture every product type. • No product is spoiled during the transportation process or during storage in distribution centers. • The distances between plants, distribution centers and customers

A supply chain is the higher abstract level of the actions taking place from conversion of raw materials and delivery of final products to the end costumers. Within this abstract structure, there exist three processes of sourcing, manufacturing and delivery of final products. A three-level supply chain including manufacturing plants, distribution centers and customer zones are addressed in this article. In this problem, a variety of products are involved which must be transported from plants to distribution centers and then to the final customers. There are limitations on plants' production capacity and distribution centers' holding capacity. Each type of product needs a specific time to be manufactured in a plant and occupies a fixed volume of space in distribution centers; consequently, the available production time in each plant and storage capacity in distribution centers are limited. The objective here is to maximize the gross profit in the given supply chain without violating the operational constraints and capacities of production and inventory. Gross profit is defined as subtracting all costs from the final revenue. Revenue is obtained from multiplying unit price of products by their final demand volume. Here, it is assumed that the final demand is not fixed and it depends on the products' price and advertising expenditure. Different types of costs related to raw materials, transportation, production and advertising are addressed in this study. Transportation costs include the costs from plants to distribution centers and the costs from distribution centers to customer zones. The production costs are divided into two groups of fixed and variable costs. There exist two types of advertising, global and local, in the supply chain, each with its own expenditures. A company's global advertising is planned to encourage potential consumers to become accustomed to its brand and to develop brand knowledge and preference. Retailers resort to local advertising to encourage customer's buying behavior. Price elasticity of demand is a concept that determines the correlation between price and demand of a product (Lysons & Gillingham, 2003). Most products, except for luxury items, have a negative price elasticity value (Nicholson & College, 2012); therefore, an increase in product's price leads to a decrease in product's demand. The demand decrease rate depends on the product's price elasticity value. Advertising elasticity is a measure of an advertising campaign's effectiveness in generating new sales. Good advertising will result in a positive effect on demand. Simon and Arndt (1980) review more than 100 articles and deduce that the shape of advertising-sales response function can perfectly be represented as diminishing returns. Similarly to many models in the literature, the expected demand function in this article is modeled in a multiplicative form, (SeyedEsfahani, Biazaran, & Gharakhani, 2011; Szmerekovsky & Zhang, 2009; Xie & Neyret, 2009; Xie & Wei, 2009; Yue, Austin, Wang, & Huang, 2006). The customer demand d (p , a, a') in Eq. (1) is taken from Xie and Wei (2009) and SeyedEsfahani et al. (2011); and it depends on selling price p , retailer’s local advertising expenditure a and manufacturer’s global advertising expenditure a'

d (p , a, a') = g (p)·h (a, a') = (D0

. p)·(K a + K ' a' )

zone are known; that is, the transportation unit costs are predetermined.

3. Mathematical model A mathematical model, named Model1, is presented here. This model maximizes the final profit in a three-level supply chain where the variable demand is obtained from Eq. (1). The indices, parameters and variables of this model are defined as follows: Indices i = 1, 2, , I j = 1, 2, , J k = 1, 2, , K m = 1, 2, , M r = 1, 2, , R Parameters CTm,i, j

CTm,j, k RMm,r CRr CFm,i CVm, i Capi Capj

Timem,i Massm Dm0 , k

Amount of time for producing one unit of product m in plant i Inventory capacity of distribution centers occupied by one unit of product m Reference demand of customer k for product m

K'm Variables pm, k

Selling price of product m to customer zone k

dm, k ym, i, j ym, j, k

xm, i z m, i am, k a'm

The function g (p) is the impact of the selling price on the final demand. Here, it is assumed that g (p) decreases with respect to price in a linear sense. Parameter D 0 is the potential demand in market when products are at no cost to customers, and is the price elasticity coefficient. The function h (a, a') is the advertising expenditures impact on demand, known as the sale response function. Parameters K and K' are positive constants reflecting the local and global advertising effect on increasing demand. The demand function in Eq. (1) is an increasing and concave function with respect to local and global advertising expenditures with a feature that additional advertising expenditure generates continuous diminishing returns. The other assumptions applied

Transportation cost for product m from plant i to distribution center j Transportation cost for product m from distribution center j to customer zone k Amount of raw material r for one unit of product m Unit cost for raw material r Fixed production cost for product m in plant i Variable production cost for one unit of product m in plant i Production (time) capacity of plant i Inventory capacity of distribution center j

Coefficient of price elasticity of demand for product m Coefficient of local advertising for product m and customer zone k Coefficient of global advertising for product m

m

Km, k

(1)

Plants Distribution centers Customer zones Products Raw materials

Final demand of product m by customer zone k Amount of product m transported from plant i to distribution center j Amount of product m transported from distribution center j to customer zone k Amount of product m manufactured in plant i A binary variable showing if product m is manufactured in plant i Local advertising expenditures for product m and customer zone k Global advertising expenditures for product m

The objective function indicates the final profit in supply chain, which is equal to difference between the total revenue and the sum of all costs. The revenue from a product is obtained by multiplying its final demand by the selling price. The other terms in Eq. (2) represent transportation, raw materials, fixed and variable production costs and advertising costs, respectively. The flow constraints for plants, distribution centers and customer zones are expressed in Eqs. (3)–(5). The plant production capacities and inventory capacities of distribution centers are implemented through constraints (6) and (7). Based on inequality (8), the setup costs must be paid for a plant in order to be able to meet some 267

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

products’ demand. The final demand of products based on their base demand, selling price and the assigned advertising expenditures, are calculated through Eq. (9). The last relation shows the types and domains of the decision variables.

Max

(dm, k · pm, k )

m, k

(CTm, i, j· ym, i, j ) +

m , i, j

x m, i ·RMm, r · CRr m, r

m , j, k

am, k +

m, i

0

m, j

(4)

i

dm, k =

ym, j, k

m, k

Capi

i

(6)

m

(Massm· ym, i, j )

Capj

j

Capi z m, i Timem, i

dm, k =

(Dm0 , k

0

0

0 0 0 0

0 0 0 0

0 0 0 0

0 0 0 0

0 0 0 0

0 0 0 0

m, i

Q12 =

(8)

m · pm, k )(K m, k

am, k + K 'm a'm )

2

m, k

0

z m, i

i , j, k, m {0, 1}

Most of the optimization solvers handle quadratic models more efficiently compared with general nonlinear models. In order to convert the above model into quadratic form, the terms am, k and a'm are

substituted by variables am, k and a'm . Hence, Eqs. (2) and (9) are replaced by Eqs. (11) and (12), while other constraints in the model remain unchanged. The model extracted from Eqs. (3)–(8) and (10)–(12) is named Model1 of mixed integer quadratically constrained type.

m, k

(dm, k · pm, k )

m , i, j

(CTm, i, j· ym, i, j ) +

x m, i ·RMm, r · CRr m, r

i

m, k

dm, k = (Dm0 , k

m , j, k

(CTm, j, k· ym, j, k )

(CFm, i·z m, i) m, i

(CVm, i·x m, i ) m, i

2

am2 , k +

a'm

(11)

m

m · pm, k )(K m, k ·am, k

+ K 'm ·a'm)

m, k

(12)

In order to do some mathematical optimization analysis on Model1, consider a general quadratically constrained quadratic programming (QCQP) of the form

min xT Q0 x + 2q0T x +

0

s. t . xT Qp x + 2qpT x +

p

0

0 0 0 0

0 0 0 0

p = 1, 2,

,P

RMm, r . CRr CFm, i

,

0

0

K m, k

2

2

K m, k 0 0 0 0 0

2

K 'm 0 0 0 0 0

K 'm

0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

1 0 1 Dm, k Km 0 0 0 0 0 2 2

T

,

=0

It is seen that none of the matrices Q0 and Q12 are positive semidefinite and hence, Model1 is not convex. In some cases, semi-definite programming (SDP) and second order cone programming (SOCP) relaxations provide exact optimal solutions for non-convex quadratic optimization problems. For example, Kim and Kojima (2003) adopt the matrix notation Mp: =[ p qpT ; qp Qp](0 p P ) based on the model n× n QCQP (13). They define a family of symmetric matrices A to be almost off-diagonal nonpositive if there exists a sign vector { 1, + 1} n such that [Ap ]i, j i j 0 (1 i < j n, 0 p P ) . They show that if all matrices Qp (0 p P ) in a QCQP model are almost off-diagonal nonpositive, then the optimal value of SDP relaxation is equal to the optimal value of QCQP (13). In theory, the optimal solution for Model1 remains unchanged when all the equality constraints are transformed to “less than or equal to” = [+ 1, + 1, 1, 1, 1, 1, 1, 1] inequalities. The vector sign T ; q12 Q12] converts the matrices M0: =[ 0 q0T ; q0 Q0] and M12: =[ 12 q12 into the almost off-diagonal nonpositive form. However, regarding the constraints (3)–(8) in Model1, there is no sign vector that makes all the corresponding matrices almost off-diagonal nonpositive and hence, semi-definite programming and second order cone programming relaxations are not applicable in our case. Based on the fact that Model1 is non-convex and has some binary variables, some heuristic algorithms are developed in the next section.

(10)

i, m

0

1 0 Dm, k Km, k 2

q12 =

(9)

0

0 0 0 0 0

12

ym, i, j , ym, j, k , dm, k , xm, i , am, k , am

Max

0 0 0 0

=0

(7)

m, i

x m, i

0 0 0 0

0

(5)

j

(Timem, i ·x m, i )

1 I 2 m.k

where Im is the identity matrix of size m with ones on the main diagonal and zeroes elsewhere. All the constraints in Model1 except Eq. (12) are linear and can easily be represented by quadratic terms of the QCQP (13). Regarding Eq. (12) as the only quadratic constraint in Model1, the related QCQP matrices are

(3)

ym, i, j

0 1 I 2 m.k

r

(2)

j

ym, j, k =

0

T

s.t.

k

0 0 0 0 0 0 0 0

q0 = 0 0 0 0 CTm, i, j CTm, j, k CVm, i +

m

ym, i, j = xm, i

0 0

m, i

a'm

m, k

Q0 =

(CVm, i·x m, i )

m, i

0 0

0

(CTm, j, k· ym, j, k )

(CFm, i·z m, i)

i

Im . k 0 0 Im

4. Heuristic algorithms Here, some theorems are developed regarding near optimal solutions for the problem at hand and later, the results of these theorems are applied in designing the algorithms. For this purpose, the fixed production cost is removed from the list of the costs, and the capacity constraints for the distribution centers and plants are eliminated from definition of the problem.

(13)

n × n are symmetric matrices, q n and where Qp for p p 0 p P . If we represent Model1 in this form, the vector of decision variables will be x = [am, k a'm pm, k dm, k ym, i, j ym, j, k xm, i z m, i ]T and regarding the objective function Eq. (11) we have,

268

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

Fig. 1. The outline for Harmony Search algorithm.

Let parameter TCm, k represent the sum of raw materials costs and average of variable production costs and transportation costs assigned to one unit product m sold in customer zone k . This parameter can be estimated as follows:

TCm, k =

i

(RMm, r . CRr ) +

CVm, i I

r

+

i, j

CTm, i, j

I. J

+

j

then, we get: Revm =0 pm, k [Dm0 ,k

CTm, j, k

pm, k =

(14)

J

Based on the definition in Section 2, the effect of advertising expenditures on the demand of a product is Km, k am, k + K 'm a'm which will be presented as Am, k in the rest of this article. The newly introduced parameter

m, k

=

0 Dm ,k

TCm, k

2 m

2

2

dm, k = (Dm0 , k

m

m · TCm, k ][K m, k

m · pm, k ) Am, k

Dm0 , k

=

m · pm, k )(Km, k

am, k + K 'm a'm )

am, k ]

am, k + K 'm a'm ] = 0

(16)

m, k

m ·TCm, k

2

2

Am, k =

m

m, k Am, k

Theorem 4.2. Without considering the fixed production costs and the production and inventory capacities, the near optimal values for global and local advertising expenditures are as follows. Proof. we have:

Proof. The total revenue generated by selling product m in all customer zones is expressed in Eq. (15); where, TCm, k is the total costs except advertising. Consequently, the local and global advertising costs are subtracted from the total revenue, in a separate manner. TCm, k )(Dm0 , k

2

+

TCm, k + = TCm,k + 2

(17)

Theorem 4.1. Without considering the fixed production costs and the production and inventory capacities, the near optimal values for products’ pm, k = TCm, k + m, k price and demand consist of: and dm, k = m m, k Am, k .

[(pm, k

m · pm, k

will be applied in facilitating the cal-

culations and simplifying the results of the given theorems.

Revm =

2

Dm0 ,k

Revm =0 am, k am, k =

a'm

m · K m, k

2

Dm0 , k 2

m

TCm, k 2

2 2

=

m · K m, k

2

2 m, k

(18)

k

(15)

269

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

calculated in customer zone k . There exists the possibility that these demand values may violate the capacity of the distribution centers or plants and this algorithm should adjust the demands accordingly (2nd and 3rd loops in Algorithm 1). If the demand dm, k has decreased due to the capacity limitation, then the algorithm increases the price, pm, k , according to Eq. (1). The next two loops select the distribution centers sort and plants based on Jmsort , j and Im, i orders to which the maximum possible demand of each product is assigned. The final revenue is yield by subtracting the sum of transportation, raw materials, advertising and production costs from the total revenue of selling products.

Revm =0 a'm (pm, k

a'm =

2

2

m · K 'm

2

m · pm, k ) K 'm

2

k

=

TCm, k )(Dm0 , k

m, k k

(19)

4.1. Harmony search algorithm

Algorithm1: Algorithm HS

In computer science and operations research fields, harmony search (HS) is a phenomenon-mimicking algorithm inspired by the improvisation process of musicians proposed by Geem, Kim, and Loganathan (2001). Harmony search is based on the assumption that the purpose of music is to achieve a perfect state of harmony analogous to finding the optimality in an optimization process. A musician always intends to produce a piece of music with perfect harmony; while, an optimal solution to an optimization problem is the best solution regarding the given objectives and limited by some constraints. Both processes tend to produce the best or optimum. The similarities between two processes can be applied in developing new algorithms which would learn from one another. The HS algorithm is a vivid example for transforming a qualitative improvisation process into some quantitative rules through idealization. An outline for harmony search algorithm is presented in Fig. 1. In the rest of this paper, the expression Unif [a , b] denotes a random number from continuous uniform distribution between a and b. Also, the expression UnifInt [a, b] represents a random number from discrete uniform distribution between a and b . The initial population contains the randomly generated vectors x h where h = 1, 2, , H and H denotes the number of harmonies in the harmony memory. At each iteration of the algorithm, a new vector x new is generated where its eth element is denoted by xenew . With the probability of harmony memory consideration rate (HMCR), the element xenew in the new vector x new is chosen from harmonies x h in the current harmony memory; otherwise, it is generated between the lower and upper bounds for the eth element (le and ue ) in a random manner. The obtained value for xenew is further examined to determine whether it needs some pitch adjustment or not. Pitch adjustment procedure is responsible for local changes in new harmonies. With the probability of pitch adjustment rate (PAR), the elements of the new vector x new are tuned. The band width (BW) factor is to control the local search around the variable xenew . After the new harmony vector is generated, if its objective function value is better than the objective value of the worst harmony vector x worst , it will replace the worst harmony vector. In this proposed algorithm, the size of decision vector x h is M × (2K + I + J + 1) . Each vector will be divided into the five matrixes of:

Step1:

For all products m and customer zones k do

Step2:

End for If j Capj <

dm, k

max{0, (Dm0 , k m, k

m ·pm, k )(Km, k

am,k + K 'm a'm )}

(dm,k ·Massm) then

For all products m and customer zones k do

dm, k pm, k

Step3:

dm, k · 1 m

j Capj m, k (dm, k·Massm )

Dm0 , k

dm, k Km, k am, k + K 'm a'm

End for End if If

i

Capi <

dm,k ·

m, k

i Timem, i I

then

For all products m and customer zones k do

dm, k

dm, k ·

i Capi i Timem, i I

m, k dm, k .

pm, k

Step4:

1 m

Dm0 , k

dm, k Km, k am, k + K 'm a'm

End for End if For all products m and customer zones k do

Select a distribution center j' based on Jmsort , j order

Assign the maximum possible demand of product m and customer zone k to center j' Update the remaining capacity for distribution center j' and remaining demanddm, k Update the inventory invm, j' of product m in distribution center j' Step5:

End for For all products m and distribution centers j do Select a plant i' based on Imsort , i order

Assign the maximum possible inventory of product m in distribution center j to planti' Update the production capacity of plant i' and remaining inventoryinvm, j

Step6:

Update the production amount of product m in planti' End for Calculate the final revenue based on obtained costs and income

The lower and upper bounds, le and ue , of the elements in decision vector x h are determined empirically. According to some numerical tests, the best intervals for the first three blocks of the vector x h consist a'm pm, k [0.8pm, k , 1.4pm, k ], am, k [0.2am, k , 1.2am, k ] of and [0.2a'm , 1.2a'm]. The parameters pm, k , am, k and a'm are calculated with respect to Theorems 4.1 and 4.2. The elements in the last two blocks of x h are chosen from the interval [1, 1000] in an arbitrary manner and the algorithm will sort them based on their priority order of importance. Based on the empirical tests, the best harmony memories’ size in this problem is H = 40 . By generating 15 new harmonies at every iteration, the algorithm keeps a good balance between run time and quality of solutions. The HMCR and PAR values are set to 0.8 and 0.2, respectively. The initial value for band width of variable x h is 0.1(ue le ), while these values are reduced by a factor of 0.995 at each iteration. The stopping occurs at 150 consecutive iterations without improving the objective function.

1. size M × K , describing the assignment of price values, pm, k 2. size M × K , describing the assignment of local advertising expenditures, am, k 3. size M × 1, describing the assignment of global advertising expenditures, am' 4. size M × J , describing the order of selecting distribution centers for providing each product, Jmsort ,j 5. size M × I , describing the order of selecting plants for providing each product, Imsort ,i sort The orders of selecting distribution centers, Jmsort , j , and plants, Im, i , for providing products’ demands are accomplished by sorting the elements of the last two matrixes. An important feature of this algorithm is the manner by which a decision vector x h is decoded into a feasible solution and the related final revenue is calculated. First, the demand dm, k of a product m is

270

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

4.2. Algorithm HS-QCP

impact on the results of this algorithm; hence, it is eliminated in Algorithm 3. The parameters of Algorithm 3 are set empirically and through some Max test problems. The selected values for parameters p0,Min m, k , p0, m, k and pchange are 0.8, 1.2 and 0.05, respectively. The results indicate that HMCR = 0.8 is an appropriate value. As to harmony memories rates, the harmony memories size in this algorithm is set to 25. Here, one new solution is generated at every iteration. The stopping is subject to 150 iterations without improving the best objective value found so far.

A heuristic algorithm is developed through combining the harmony search algorithm and the mathematical model from Section 2. At each iteration of algorithm HS-QCP, the variables pm, i and z m, i are fixed through the harmony search procedure and then, the relaxed mathematical model is solved. This mathematical model is of quadratic programming type since the objective function is quadratic and all the constraints are linear and the unfixed decision variables are continuous. It is obvious that the objective function here is convex and the GAMS/ CPLEX solver can easily find optimal solutions. The algorithm HS-QCP is composed of two phases; one generates the initial solutions and the other is about improving iterations. In the first phase, (Algorithm 2), some initial solutions are generated and are added to the harmony memory. In Step1 of Algorithm 2, parameter Usedcap is calculated which represents the minimum used capacity between plants and distribution centers occupied by base demands Dmo , k . For each harmony h in Step2, variables pm, k and z m, i are randomly generated and are fixed in Model1. Then, the other variables for each solution h are obtained via solving Model1. Min Max , pinit ]· pm, k in a The prices pm, k are selected from interval [pinit random manner, where pm, k is taken from Eq. (16). Based on numerical Min = 0.8 and experiments, the values of these parameters are set at pinit Max pinit = 1.2 . The term

0 Dm ,k m

Algorithm 3: Improving iterations for Algorithm HS-QCP Step1:

zmAve ,i

Step2:

Step2:

m, k

i Timem, i . D 0 m, k I i Capi

,

min

0 Dm ,k m

pitMin1, m,k + (1

pchange

pitMax , m, k

pitMax 1, m, k

pchange

Else pm, k

pitMin1,m, k ) Max (pit 1, m,k 1)

min

0 Dm ,k m

Max , Unif [pitMin , m, k , pit , m, k ] pm, k

randomly chosen from a solutionin HM

End if End for For all products m and plants i do If Unif [0, 1] > HMCR then z m, i UnifInt [0, 1] Else if zmAve , i > Unif [0, 1] then

z m, i 1 Else z m, i 0 End if End for Solve Model1 with fixed variables p andz If current profit is better than the worst profit in HM then Replace the worst solution in HM with the current solution End if Repeat the main loop until the stopping criteria is met End for

0 m, k (Massm . Dm, k ) j Capj

For all harmonies h in harmony memory do For all products m do For all customer zones k do

pm, k

pitMin , m, k

pm, k

Algorithm 2: Initial solution for Algorithm HS-QCP

min

H

If Unif [0, 1] > HMCR then

in Algorithm 2 is applied for preventing pm, k

UsedCap

h zm ,i

For all iterations it do For all products m and customer zones k do

from taking negative values. The binary variables z i, m are generated in a manner that Usedcap percent of them are equal to 1 and the remainder are equal to 0. In the innermost IF-block in Algorithm 2, a random number from interval [0,1] is generated and if this number is less than or equal to Usedcap then we set z i, m = 1 and otherwise, we set z i, m = 0 .

Step1:

h

5. Numerical experiments

Min Max , Unif [pinit , pinit ]. pm, k

To examine the efficiency of the proposed solution methods, several numerical experiments are run. Four different groups of test problems are designed based on Table 2. The columns of this table respectively indicate the number of plants, distribution centers, customer zones, products and raw materials. All parameters are generated in a random manner and the intervals from which the parameters are generated are tabulated in Table 3. The signs Unif and UnifInt represent continuous and discrete uniform distributions, respectively. The mathematical model Model1 is coded in GAMS v24.1 and the solver COIN-OR Bonmin is applied for solving it. Algorithm HS-QCP is solved in GAMS v24.1 environment and CPLEX v12 is applied for solving the relaxed mathematical model at each iteration. The algorithm HS is solved in Visual Studio C++ 2010 environment. All the methods are run on a personal computer equipped with 3.1 GHz Pentium IV processor and 4 GB of RAM. The trend of convergence for algorithms HS and HS-QCP executed

End for For all plants i do If Unif [0, 1] Usedcap then

z m, i 1 Else z m, i 0 End if End for End for Solve Model1 with fixed variables p andz Save the solution in harmony memory End for

The improving iterations are described in Algorithm 3. In Step1 of this algorithm, parameter z mAve , i is calculated representing the average of z m, i in all harmonies h in the harmony memory. At each iteration of Step2, the variables pm, k and z m, i are randomly generated and are fixed in Model1. Then, Model1 is solved and if its profit is better than the worst profit in HM then, the worst solution in HM is replaced by the current solution. At every iteration of Algorithm 3, the values for variables pm, k and z m, i are generated with respect to HMCR probability in a random manner and are selected from one solution in HM with probability 1HMCR. In order to improve the convergence of this algorithm, the intervals for generating pm, k are gradually tightened as they approach pm, k . The speed of tightening procedure is controlled by parameter pchange which is set at 0.1, experimentally. Results obtained from the experiments reveal that pitch adjustment procedure has no significant

Table 2 Characteristics of four groups of test problems.

Group Group Group Group

271

G1 G2 G3 G4

I

J

K

M

R

5 5 10 15

8 12 20 30

15 20 35 60

10 15 25 40

6 12 20 25

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

Table 3 Parameters' values used in creating test problems. CTm,i, j

UnifInt [10, 100]

Timem,i

Unif [1, 7]

CTm,j, k

UnifInt [50, 150]

Massm

Unif [1, 4]

RMm,r

Unif [1, 5]

Dm0 , k

CRr

UnifInt [8, 20]

CFm,i CVm, i Capi

UnifInt [10000000, 15000000] UnifInt [200, 500] Unif [0.8, 1.2]·

UnifInt [35000, 50000]

Unif [30, 45]

m

m

0 k Dm, k

Km, k

Unif [0.001, 0.003] Unif [0.0003, 0.0005]

K'm Capj

i Timem, i I

Unif [0.8, 1.2]·

m

I

0 k Dm, k Massm J

Fig. 2. Convergence of algorithms HS and HS-QCP.

on 10 replications of problem instances in group G3 is shown in Fig. 2. It is evident that both proposed algorithms benefit from good convergence to near optimal solutions. The results of solving Model1 from Section 3 and algorithms HS and HS-QCP on different groups of test problems are tabulated in Table 4. Ten replications are executed for each algorithm and each problem size. The data regarding the average of these replications are tabulated in Table 4. Model1 reaches optimal solution for group G1 problems in 1039 s in average. In case of groups G2 and G3 this model only reaches the initial relaxed continuous solutions in a predefined 3600 s time limit. The relaxed continuous solutions are applied as upper bounds in calculating optimality gaps of algorithms HS and HS-QCP in groups G2

and G3. As to group G4 test problems, Model1 did not reach the relaxed solutions within the predefined time limit (1 h); therefore, with respect to the better solutions obtained in previous groups, algorithm HS-QCP is considered as the reference and the optimality gaps for algorithm HS are calculated after their comparison with results obtained from HSQCP. In group G1, both HS and HS-QCP generate optimality gaps of less than 2%. Optimality gaps for groups G2 and G3 are less than 7% compared to the upper bound from Model1, revealing the efficiency of both HS and HS-QCP in approaching near optimal solutions even in large-scale problems. The last column in Table 4 reveals that the solutions provided through HS-QCP are slightly better than that of the HS.

Table 4 Results of implementing Model1 and algorithms HS and HS-QCP on test problems. Group

No. equations

No. variables

Relaxed sol.

Final sol.

Gap (%)

Binary

Continuous

Time (s)

No. iter

Time (s)

No. iter

G1

Model 1 HS-QCP HS

500 500 —

50 — —

2116 1807 —

59.29 — —

1802 — —

1039 65.1 99.99

9678 297.9 511.8

— 1.13 1.39

G2

Model 1 HS-QCP HS

954 954 —

75 — —

5496 4882 —

69.64 — —

1256 — —

— 128.4 156.05

— 428.1 434.4

— 5.05 5.71

G3

Model 1 HS-QCP HS

2787 2787 —

250 — —

25,406 23,632 —

457 — —

3154 — —

— 392.7 444.06

— 469.3 470.7

— 6.11 7.85

G4

Model 1 HS-QCP HS

7252 7252 —

600 — —

97,846 93,007 —

— — —

— — —

— 1513.2 667.98

— 540.4 279.1

— — 4.34

272

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

Fig. 5. Sensitivity analysis on the unit transportation costs.

Fig. 3. Sensitivity analysis on the capacity of production and distribution centers.

5.1. Sensitivity analysis A comprehensive sensitivity analysis is run on the main parameters in the supply chain under study. Test problems in group G3 and algorithm HS-QCP are selected for implementing the sensitivity analysis. The capacities of plants and distribution centers are first multiplied by coefficients 0.25, 0.5, 1.5 and 2 and then, the changes on amounts of profits, costs, demands and prices in the supply chain are examined. By increasing the capacity values, obviously the satisfied demands and the final profits would increase (see Fig. 3). The changes in base demand values (D 0 ) are bar-charted in Fig. 4. The more the markets’ potential, the more the demand for products and consequently, an increase in the final costs and profit in the supply chain. As observed in Fig. 4, selling prices can be increased when the potential demand expands in a market, thus, higher revenue and more profit. As observed in Fig. 5, any increase in unit transportation costs will lead to a decrease in demands, other related costs, final revenues and profits. When parameter CT gains 1.5 times of its initial value, the total transportation costs increase; however, transportation costs decrease at CT = 2 because of extreme reduction in demand. The correlation between selling prices and parameter CT is not significant. Price elasticity coefficient ( ) in Fig. 6 has a reverse correlation with demand and the final price of products. It is because when customers become very sensitive to price changes, their demands decrease in drastic slope in response to price increase. The fixed transportation costs are not significantly affected by the changes in parameters. It means almost a fixed number of vehicles are needed while the parameters alter. The sensitivity analysis for advertising elasticity coefficient is bar-charted in Fig. 7. An increase in advertising elasticity does not have a significant effect on price, while a strict increase is evident as to demand. The final profits and revenues change consistent with changes in advertising elasticity coefficients. In general, it can be claimed that bigger market size (parameter D 0 )

Fig. 6. Sensitivity analysis on price elasticity coefficient.

Fig. 7. Sensitivity analysis on advertising elasticity coefficient.

will yield higher profits. Price elasticity has a reverse effect on selling strength in the market. Bigger advertising elasticity coefficients can make the market more prosperous; hence, an improvement in profitability in the supply chain. 6. Conclusion This article studied a three-level supply chain considering two of the most important factors affecting the demands in a supply chain, the price and advertising elasticity. Considering price and advertising-dependent demands next to different types of cost sources makes the problem close to the real world problem. A mathematical model of type mixed integer quadratically constrained is presented and a harmony search algorithm is developed to solve this problem. The third proposed method combines the principles of harmony search and mathematical modeling to develop a new heuristic algorithm. In order to examine the quality of these proposed methods, four groups of test problems of different size are generated. The mathematical model is only able to solve the small-size test instances in one hour time limit while the proposed two algorithms have the optimality gaps

Fig. 4. Sensitivity analysis on the base demand values. 273

Computers & Industrial Engineering 135 (2019) 265–274

K. Kianfar

less than 2% compared to upper bound from the mathematical model in small and medium-size instances. The optimality gap from upper bound is less than 7% even in large-scale instances, indicating the efficiency of these proposed algorithms. The sensitivity analyses reveal that the final profit increases through an increase in the capacity of plants and distribution centers. The market potential demand has a positive effect on the final demand and the total profit in the supply chain. Any increase in unit transportation costs leads to a decrease in demand, other costs and the final profit. There exists a reverse correlation between price elasticity and final price of products and advertising elasticity coefficients do not change the prices but increase the final demand. Many possible future research paths may be defined for the current study. The base demand values or capacity of manufacturers and distribution centers can be considered as stochastic or fuzzy values. Considering the elasticity of other factors affecting the demand in a supply chain like price of substitute products is another intriguing issue to be studied.

167, 12–22. Hull, B. (2005). The role of elasticity in supply chain performance. International Journal of Production Economics, 98(3), 301–314. Kaplan, U., Türkay, M., Karasözen, B., & Biegler, L. T. (2011). Optimization of supply chain systems with price elasticity of demand. INFORMS Journal on Computing, 23(4), 557–568. Karray, S., & Amin, S. H. (2014). Cooperative advertising in a supply chain with retail competition. International Journal of Production Research, 53(1), 88–105. Kim, S., & Kojima, M. (2003). Exact solutions of some nonconvex quadratic optimization problems via SDP and SOCP relaxations. Computational Optimization and Applications, 26, 143–154. Kozlenkova, I. V., Hult, G. T. M., Lund, D. J., Mena, J. A., & Kekec, P. (2015). Review the role of marketing channels in supply chain management. Journal of Retailing, 91(4), 586–609. Li, Liying, Wang, Yong, & Yan, Xiaoming (2013). Coordinating a supply chain with price and advertisement dependent stochastic demand. The Scientific World Journal, 2013, 1–12. https://doi.org/10.1155/2013/315676. Li, B., Zhu, M., Jiang, Y., & Li, Z. (2016). Pricing policies of a competitive dual-channel green supply chain. Journal of Cleaner Production, 112(3), 2029–2042. Liu, S., Kasturiratne, D., & Moizer, J. (2012). A hub-and-spoke model for multi-dimensional integration of green marketing and sustainable supply chain management. Industrial Marketing Management, 41(4), 581–588. Lysons, K., & Gillingham, M. (2003). Purchasing and supply chain management: Financial. Times Prentice Hall. Ma, P., Wang, H., & Shang, J. (2013). Supply chain channel strategies with quality and marketing effort-dependent demand. International Journal of Production Economics, 144(2), 572–581. Mokhlesian, M., & Zegordi, S. H. (2014). Application of multidivisional bi-level programming to coordinate pricing and inventory decisions in a multiproduct competitive supply chain. The International Journal of Advanced Manufacturing Technology, 71(9), 1975–1989. Nicholson, W., & College, A. (2012). Microeconomic theory: Basic principles and extensions. Dartmouth College: Cengage Learning. Seifbarghy, M., Nouhi, K., & Mahmoudi, A. (2015). Contract design in a supply chain considering price and quality dependent demand with customer segmentation. International Journal of Production Economics, 167, 108–118. SeyedEsfahani, M. M., Biazaran, M., & Gharakhani, M. (2011). A game theoretic approach to coordinate pricing and vertical co-op advertising in manufacturer–retailer supply chains. European Journal of Operational Research, 211(2), 263–273. Simon, J. L., & Arndt, J. (1980). The shape of the advertising response function. Journal of Advertising Research, 20(4), 11–28. Szmerekovsky, J. G., & Zhang, J. (2009). Pricing and two-tier advertising with one manufacturer and one retailer. European Journal of Operational Research, 192(3), 904–917. Taleizadeh, A. A., Akhavan Niaki, S. T., & Barzinpour, F. (2011). Multiple-buyer multiplevendor multi-product multi-constraint supply chain problem with stochastic demand and variable lead-time: A harmony search algorithm. Applied Mathematics and Computation, 217(22), 9234–9253. Wong, W. K., & Guo, Z. X. (2010). A hybrid intelligent model for medium-term sales forecasting in fashion retail supply chains using extreme learning machine and harmony search algorithm. International Journal of Production Economics, 128(2), 614–624. Xie, J., & Neyret, A. (2009). Co-op advertising and pricing models in manufacturer-retailer supply chains. Computers & Industrial Engineering, 56(4), 1375–1385. Xie, J., & Wei, J. C. (2009). Coordinating advertising and pricing in a manufacturer–retailer channel. European Journal of Operational Research, 197(2), 785–791. Yue, J., Austin, J., Wang, M.-C., & Huang, Z. (2006). Coordination of cooperative advertising in a two-level supply chain when manufacturer offers discount. European Journal of Operational Research, 168(1), 65–85. Zhang, J., Gou, Q., Liang, L., & Huang, Z. (2013). Supply chain coordination through cooperative advertising with reference price effect. Omega, 41(2), 345–353.

Declaration of Competing Interest None. References Alaei, S., & Khoshalhan, F. (2015). A hybrid cultural-harmony algorithm for multi-objective supply chain coordination. Scientia Iranica, 22(3), 1227–1241. Aust, G. (2015). Vertical cooperative advertising in supply chain management. Switzerland: Springer International Publishing. Aust, G., & Buscher, U. (2014). Cooperative advertising models in supply chain management: A review. European Journal of Operational Research, 234(1), 1–14. Banyai, T. (2011). Optimisation of a multi-product green supply chain model with harmony search. In B. Katalinic (Ed.). DAAAM International Scientific (pp. 15–30). Vienna, Austria: DAAAM International. Brindley, C., & Oxborrow, L. (2014). Aligning the sustainable supply chain to green marketing needs: A case study. Industrial Marketing Management, 43(1), 45–55. Chen, T. H. (2014). On the impact of cooperative advertising and pricing for a manufacturer-retailer supply chain. Journal of Industrial and Production Engineering, 31(7), 417–424. Chen, T. (2015). Effects of the pricing and cooperative advertising policies in a twoechelon dual-channel supply chain. Computers & Industrial Engineering, 87, 250–259. Geem, Z. W., Kim, J. H., & Loganathan, G. V. (2001). A new heuristic optimization algorithm: Harmony search. Simulation, 76(2), 60–68. Giri, B., & Sharma, S. (2014). Manufacturer's pricing strategies in cooperative and noncooperative advertising supply chain under retail competition. International Journal of Industrial Engineering Computations, 5(3), 475–496. Green, K. W., Jr, Whitten, D., & Inman, R. A. (2012). Aligning marketing strategies throughout the supply chain to enhance performance. Industrial Marketing Management, 41(6), 1008–1018. Guillén, G., Badell, M., Espuña, A., & Puigjaner, L. (2006). Simultaneous optimization of process operations and financial decisions to enhance the integrated planning/ scheduling of chemical supply chains. Computers & Chemical Engineering, 30(3), 421–436. Hong, X., Xu, L., Du, P., & Wang, W. (2015). Joint advertising, pricing and collection decisions in a closed-loop supply chain. International Journal of Production Economics,

274