Modeling and optimizing an agro-supply chain considering different quality grades and storage systems for fresh products: a Benders decomposition solution approach.

This paper proposes a mathematical model in the context of agro-supply chain management, considering specific characteristics of agro-products to assist purchase, storage, and transportation decisions. In addition, a new method for determining the required quality score of different types of products is proposed based on their loss factors and purchasing costs. The model aims to minimize total cost imposed by purchasing fresh products, opening warehouses, holding inventories, operational activities, and transportation. Two sets of examples, including small and medium-sized problems, are implemented by general algebraic modeling language (GAMS) software to evaluate the model. Then, Benders decomposition (BD) algorithm is applied to tackle the complexity of solving large-sized instances. The results of both GAMS and BD are compared in terms of objective function values and computational time to demonstrate the efficiency of the BD algorithm. Finally, the model is applied in a real case study involving an apple supply chain to obtain managerial insights.

Introduction

Over recent decades, agricultural logistics has drawn the attention of researchers and practitioners due to increasing concerns in economic, healthy feeding, environmental, and political issues. Hence, good management in post-harvest operations, such as packaging, storing, and distributing agro-products, can add significant value to logistics while satisfying customers’ demands and maintaining food security. Therefore, research on agricultural logistics can speed up modernization and competition abilities in agricultural processes (Otchia 2014).

The specific characteristic of agro-products, such as seasonality, perishability, and long procurement time, make it different from other supply chains, therefore, it is necessary to have accurate information about their effects on supply chain management. These characteristics put restrictions on logistic decisions and make deep concerns over resource depletion and the occurrence of agricultural wastes and losses. According to Widodo et al. (2006), 20–60% of agro-products are wasted along supply chains in any country. Therefore, designing efficient inventory and transportation systems for perishable products plays a vital role in cost-saving and freshness-keeping (Coelho and Laporte 2014). From this perspective, cold chains with proper temperature adjustment can significantly reduce losses (Yu and Xiao 2017).

Although efficient agro-supply chain management prevents storage and transportation losses, facility location problem is another issue strongly affecting the security and safety of foods, especially in developing countries (Morton et al. 2008). Given capacity constraints in warehouses, determining warehouse locations and assigning cold storage rooms under controllable temperature conditions play a pivotal role in loss and waste reduction over a supply chain that may yield a complicated problem (Aung and Chang 2014). Forasmuch as the current study deals with agricultural products with a short shelf life, optimal facility location and storage condition decisions can reduce losses due to the dispersed distribution centers, long distances transportation, etc.

This study presents an agro-supply chain based on a case study involving a wide range of decisions such as supplier selection, facility location, inventory control, and transportation planning. This research aims to optimize the shelf life of products via assigning proper cold storage systems that can prevent postharvest losses. In addition, different quality grades are considered for each type of product, whereby those farmers who have gotten the average quality score are selected to satisfy warehouses’ demands. For these purposes, a mixed-integer nonlinear programming (MINLP) problem is formulated considering several operational constraints. Then, the model is linearized and transformed to mixed-integer linear programming (MILP) problem and solved by GAMS/CPLEX solver in small and medium sizes. Moreover, Benders decomposition (BD) algorithm is proposed to tackle the complexity of the problem caused by increasing the number of binary and integer variables in large sized-instances. Then, the proposed model is implemented for a real case study in Iran covering eight months of a year.

The rest of the paper is organized as follows: Sect. 2 gives an overview in the context of agro-supply chains. In Sect. 3, problem definition, notations, assumptions, and the mathematical model are presented. In Sect. 4, BD algorithm is proposed to tackle the complexity of large-sized instances. In Sect. 5, numerical experiments and a case study are performed to demonstrate the practicality of the proposed model. Section 5 draws conclusions and remarks on some future research directions.

Literature review

Economic costs in agro-supply chains have been implied by inefficient management in purchasing, storage, and distribution of agro-products (Ahumada and Villalobos 2009). Although there are copious studies dealing with facility location and transportation problems, they have not been widely discussed for agricultural products. The facility location problem is a critical element of strategic decisions in agro-supply chains, especially in developing countries (Morton et al. 2008). This is because transportation and storage time will become shorter by properly determining facility locations (Morganti and Gonzalez-Feliu 2015). Therefore, dealing with the capacitated facility location problems with applications in storing and distributing agro-products is an essential component to adapt to environmental and competitive aspects (Validi et al. 2014).

Here, related literature with emphasis on three core parts of agro-supply chains, purchasing, facility location, and transportation decisions, has been reviewed to pinpoint the research gap. In this context, Paam et al. (2016) and Ganesh Kumar et al. (2017) have comprehensively reviewed agro-supply chain models. Then, Banasik et al. (2017) analyzed tradeoffs between economic and environmental benefits of an agro-food closed-loop supply chain with location-allocation decisions and production planning. They attempted to valorize agro-wastes by applying recycling concepts. Ghezavati et al. (2017) proposed a mathematical model for a distribution planning system considering the freshness and ripeness of products. They proposed BD algorithm to solve their problem in a reasonable time. Orjuela-Castro et al. (2017) presented a MILP model for a multi-product facility location problem. They aim to minimize the costs of transportation, installation, and fruit losses and show how intermediate logistics platforms affect fruit losses. Soto-Silva et al. (2017) formulated a two-stage agro-supply chain model for a case study in Chile. The first stage contains particular decisions on supplier selection and transportation planning and the second stage associates with establishing different types of cold chambers to store products and transportation planning. Allaoui et al. (2018) proposed a two-stage programming model for a sustainable agro-food supply chain. In the first stage, a hybrid analytic hierarchy process with an ordered weighted averaging method was used for assigning weights to sustainable indicators. In the second stage, a multi-objective mathematical model was proposed with respect to the results in the first stage. Paam et al. (2018) proposed an agro-supply chain model to optimize the shelf lives of apples with different quality grades in cold storage rooms. Liu et al. (2018) studied optimal purchase decisions and inventory retrieval policies for a single-product finite period agro-supply chain to maximize total profit under fluctuating prices. Jonkman et al. (2019) investigated the impacts of location-allocation, inventory, and transportation decisions on an agro-food supply chain under uncertain conditions. Liu et al. (2019) developed a transportation model for fresh agricultural products considering the impacts of ethephon concentration, time, and temperature on the maturity of soluble solids. They used ant colony algorithm to solve their model. Onggo et al. (2019) assessed stochastic demand in a multi-period inventory routing problem and attempted to minimize inventory, transportation, and shortage costs. They used a local metaheuristic algorithm to tackle the complexity of their model. Paam et al. (2019) proposed a new model for inventory management optimization in an apple supply chain. They analyzed the effects of two scenarios in storage rooms on fruit losses. Yu et al. (2019) investigated the impacts of horizontal and vertical integration on the quantity/quality losses and retailer profits.

The prior work on agro-supply chains is summarized in Table 1 to reveal existing gaps in the literature. The majority of the research has investigated multi-level supply chains under several constraints, without sufficient consideration on specific characteristics of agro-products, farmers, and distribution processes. However, a survey of practical experiences illustrates significant impacts of these properties on spoilage rate, in turn, on the amount of loss and system cost. Therefore, this paper takes several actions to fill the gaps, which are vividly explained as follows:

Table 1

A literature review on agro-supply chain models

Researchers	Products		SLD	Operational constraints					Solution procedures
	single	multiple		Inventory balance	Transportation decision	Facility location	Quality scores	CFD
Banasik et al. (2017)	–		–		–		–	–	ε-constraint
Ghezavati et al. (2017)		–	–		–		–	–	Solver-BD
Orjuela-Castro et al. (2017)	–		–	–	–		–	–	Meta-heuristic
Soto-Silva et al. (2017)	–						–		Solver
Allaoui et al. (2018)	–		–			–	–	–	Solver
P Paam et al. (2018)		–			–	–	–		Solver
Jonkman et al. (2019)	–		–		–		–	–	Solver
Onggo et al. (2019)		–	–			–	–	–	Simulation-meta-heuristic
P. Paam et al. (2019)		–			–	–	–		Solver
Current study	–								Solver-BD

SLD decision on shelf life, CFD decision on cold storage facilities

This study is an extension of the work done by Soto-Silva et al. (2017) in a way to be applicable for an apple supply chain in Iran. In their study, unlike many agro-supply chain models, three types of refrigeration technologies for maintaining different types of agro-products have been proposed that can make a privileged contribution in this context. In spite of their primary step towards sustainability of agro-products, it lacks an in-depth discussion on adjustment of storage conditions and its effects on inventory decisions. However, practical accessibility gave us valuable insights to design a mathematical model for adjusting storage conditions according to costs imposed on the system and its effects on the durability of agro-products. In addition, the storage of agro-products for a long-time mainly depends on their initial quality. Therefore, this study aims to manage the quality of agro-products by giving scores to producers and accordingly selecting them. To the best of the authors’ knowledge, there is no previous effort on quality management despite its importance in handling the shelf life of agro-products. Overall, the main contributions of the paper are as follows:

A MINLP model has been presented to jointly optimize facility location, transportation, inventory decisions in an agro-supply chain.

A new method is proposed to assign weights to suppliers based on their loss factor. Accordingly, decisions on supplier selection are made.

The current research seeks to optimize the shelf life of products through adjusting appropriate storage conditions.

The proposed MINLP method is linearized and transformed to a MILP model.

BD algorithm has been utilized to tackle the complexity of the proposed model and solve it in a reasonable time.

The proposed model is implemented in a real case study in Iran to demonstrate its applicability and efficiency.

Problem definition

This study is an extension of the proposed model by Soto-Silva et al. (2017), who addressed supplier selection, transportation, and facility location decisions in an agro supply chain. However, their study lacks in-depth discussions on quality management and storage and inventory planning. Therefore, we formulate a multi-period MINLP model under some operational constraints to intertwine supplier selection, transportation, storage operation, inventory planning, and facility location decisions in an agro-supply chain. Afterwards, we reformulate it to a MILP model and propose BD algorithm to solve it. Figure 1 is meant to show the process of this paper.

Fig. 1

Overview of the paper

Figure 2 is a perspective of the proposed agro-supply chain, including producers, warehouses, and a wholesaler. Fresh products should be purchased from the producers who gain the required score based on their quality grades. This approach satisfies the required quality of products to be kept for a long time in cold chambers, which is more profitable for the system. The products are shipped to the warehouses by the private fleet. The reason is that producers, either small or medium sizes, are often located in rural and sparsely populated areas. The long distances between orchards and facilities and busy harvest seasons require the producers to provide flexible operating by their private fleet. The processing operations should carry before the storage of products, by which products are sorted based on their varieties, sizes, etc. Then, each type of products regardless of their sizes, color, etc., are kept in a cold storage room, and the temperature conditions should be adjusted to optimize the shelf life of products. Afterward, the products are delivered to the wholesaler to satisfy customers’ demands. Herein, transportation decisions are made by warehouses. In total, the proposed model aims to minimize the total cost of the agro-supply chain model while satisfying the constraints.

Fig. 2

Structure of the proposed agro-supply chain

In order to formulate the mathematical model, the following assumptions are considered. In addition, indices, parameters, and decision variables are presented in Table 2.

Table 2

Indices, parameters, and decision variables

Indices

i: Index of products (i∈I) t,g: Index of periods (t,g∈T) c: Index of temperature conditions (c∈C) w: Index of potential warehouses (w∈W) k: Index of vehicles (k∈K) j: Index of producers (j∈J) s: Index of quality grades (s∈S) v: Index of cold chambers (v∈V)

Parameters

CT_jw: Operational costs to transport products from producer j to warehouse w P_w: Operational costs in warehouse w HW_c: Holding cost of products under temperature condition c CFT_k: Fixed cost of transporting products by vehicle type k CW_w: Fixed cost of opening a warehouse at location w W_w: Available capacity of a cold chamber in warehouse w FW_w: Number of available cold chambers in a warehouse if it is opened at location w CV_k: Transportation capacity of vehicle type k CQ_ijs: Purchase cost of product i with quality grade s supplied by producer j T_ic: Shelf life of product i under temperature condition c D_it: Quantity of demand for product i in period t LF_ijs: Loss factor of product i with quality grade s purchased from producer j CF_i: Supply capacity of product type i WEP_ijs: Quality score of product i with quality grade s supplied by producer j a_it: Binary parameter; 1 if product i is available in period t, 0 otherwise α_i: Minimum acceptable quality score for product i

Variables

QF_ijst: Quantity of product i with quality grade s purchased from producer j in period t LOSS_ijs: Quantity of product loss i with quality grade s purchased from producer j QFW_ijwvct: Quantity of product i purchased from producer j by warehouse w to be kept in cold chamber v under temperature condition c in period t QW_iwvcktg: Quantity of stored product i in warehouse w in cold chamber v under temperature condition c in period g to deliver to the wholesaler by vehicle type k in period t I_iwvctg: Product i stored in warehouse w in cold chamber v under temperature condition c in period g and remained until period t TE_iwvct: Storage time for product i in warehouse w in cold chamber v under temperature condition c in time period t y_kt: Number of vehicles k for transporting products from warehouses to the wholesaler in period t A_ijst: Binary variable, 1 if product i with quality grade s is supplied by producer j in period t, 0 otherwise AC_iwvct: Binary variable; 1 if cold chamber v is turned on in warehouse w under temperature condition c in period t, 0 otherwise AW_w: Binary variable; 1 if a warehouse is opened at location w, 0 otherwise

: Index of products

: Index of periods

: Index of temperature conditions

: Index of potential warehouses

: Index of vehicles

: Index of producers

: Index of quality grades )

: Index of cold chambers

Operational costs to transport products from producer to warehouse

: Operational costs in warehouse

: Holding cost of products under temperature condition

: Fixed cost of transporting products by vehicle type

: Fixed cost of opening a warehouse at location

: Available capacity of a cold chamber in warehouse

: Number of available cold chambers in a warehouse if it is opened at location w

: Transportation capacity of vehicle type

: Purchase cost of product with quality grade supplied by producer

: Shelf life of product under temperature condition

: Quantity of demand for product in period

: Loss factor of product with quality grade purchased from producer

: Supply capacity of product type i

: Quality score of product with quality grade supplied by producer

: Binary parameter; 1 if product is available in period , 0 otherwise

: Minimum acceptable quality score for product i

: Quantity of product with quality grade purchased from producer in period

: Quantity of product loss with quality grade s purchased from producer

Quantity of product purchased from producer by warehouse to be kept in cold chamber under temperature condition c in period

: Quantity of stored product in warehouse in cold chamber under temperature condition in period to deliver to the wholesaler by vehicle type in period

: Product stored in warehouse in cold chamber under temperature condition in period and remained until period t

Storage time for product in warehouse in cold chamber under temperature condition in time period

: Number of vehicles for transporting products from warehouses to the wholesaler in period

: Binary variable, 1 if product with quality grade is supplied by producer in period , 0 otherwise

: Binary variable; 1 if cold chamber is turned on in warehouse under temperature condition in period , 0 otherwise

: Binary variable; 1 if a warehouse is opened at location , 0 otherwise

Assumptions

Shortages are not allowed.

Three different quality grades are considered for products in orchards.

Each of the quality grades for each variety of products imposes different amounts of losses.

Each variety of products is available for purchase at a specific period. Therefore, they are purchased with sufficient quality to be stored in cold chambers in order to meet customers’ demands for a while.

Decisions on the type of temperature conditions should be made for each cold chamber in warehouses.

Mixed-integer nonlinear model

In this section, a MINLP model is presented as follows:

The objective function (1) minimizes total cost consisting of six components, purchasing cost, operational costs including packing, laborers, and transportation costs (from producers to warehouses), operational costs in warehouses, fixed cost of opening warehouses, fixed cost of transportation (from warehouses to the wholesaler), inventory holding costs. Constraints (2) ensure that the purchased products are of satisfactory quality with reasonable prices. Constraints (3) ensure that only one type of product should be stored in each cold chamber under a specific temperature condition. Constraints (4) indicate the effects of temperature conditions on the shelf life of products. Constraints (5) ensure that products with the same quality grade are purchased from a producer. Constraints (6) show the relation between variables and . Constraints (7) ensure that the amount of product in warehouses should not exceed the amount of purchased product. Equations (8) and (9) show the inventory level of products in warehouses where and , respectively. Constraints (10) state the limited storage capacity of cold chambers in warehouses. Constraints (11) are the demand response. Constraints (12) represent the limitation in the number of cold chambers in a warehouse. Constraints (13) are capacity constraints of vehicles for transporting products from warehouses to the wholesaler. Equations (14) calculate the losses after purchasing products from producers. Constraints (15) indicate the available capacity of producers to supply product . Constraints (16) define types of decision variables which are non-negative continuous, binary, or integer variables.

Model linearization and reformulation

In the proposed model, the upper limit of index is considered as a decision variable making the model nonlinear. In these cases, if–then decisions should be made to remove spoiled products from warehouses. Indeed, if the storage time exceeds the shelf life of products, the products should be removed from warehouses (If . Conversely, if the storage time is less than the shelf life of products, Eq. (8) are used to calculate the inventory level of products (If ). In order to handle the nonlinear formulation and transform the MINLP model to MILP, a linearization technique is derived from Chen et al. (2010).

Consider the following nonlinear equivalents where and b’ are the right-hand side of inequalities and and are the functions of variable . In this case, if is less than , then . The linearization of the if–then constraint is applied as follows:where BIG and y are a positive big number and binary variable, respectively. Accordingly, the following equations are presented to linearize the if–then constraints in this paper.

Solution procedure

The complexity of MILP problems rises by increasing the number of binary and integer variables, consequently, has a strong influence upon the computational time (Keshavarz-Ghorbani and Pasandideh 2021). Therefore, it is necessary to apply efficient solution procedures to overcome the computational burden. BD algorithm is known as one of the powerful algorithms to tackle the complexity of NP-hard problems. BD algorithm uses the concept of "divide and conquers" to reformulate the original problem into a master problem (MP) and dual subproblem (DSP) in a limited number of iterations. Nowadays, advantages of the BD algorithm, such as a limited number of iterations, less computational time, and guaranteed convergence for large-scale problems, make it suitable for many optimization problems (e.g., Pishvaee et al. 2010; Shaw et al. 2016; Ghezavati et al. 2017; Kayvanfar et al. 2018; Babazadeh et al. 2019; Azadi et al. 2019; Mardan et al. 2019; Naderi et al. 2019; Alkaabneh et al. 2020; Gong and Zhang 2021).

The primary steps of the BD algorithm work iteratively by outer linearization and relaxation. In other words, the original problem is decomposed into an MP and subproblem (SP). The MP is a pure integer problem and includes integer and binary variables. The MP is solved and puts its solutions among the SP. In this step, the SP turns into DSP to make the feasible solution space independent from integer solutions. Then the DSP is solved, and its optimal bounded solutions/ unbounded solutions are added as feasibility/optimality cutting planes into the MP. This procedure reduces the computational time compared to the original problem. The convergence of BD comes from closing upper and lower bounds obtained from solving the DSP and MP, respectively. Figure 3 illustrates the steps of the iterative BD algorithm.

Fig. 3

The flowchart of BD

Benders master problem

To implement BD for the proposed model, the MP is called to minimize a part of the objective function associated with binary and integer variables with respect to constraints (3)–(5), and (12). That is a relaxation of the original problem to obtain lower bounds.

Benders subproblem

The SP includes a part of the objective function associated with continuous variables subject to constraints (2), (6)–(11), (13)–(15) that can be seen in Table 3.

Table 3

Equations related to the SP and its corresponding dual variables

Equations	Dual variables
MinZ=∑_i∈I∑_j∈J∑_s∈S∑_t∈Ta_itCQ_ijsQF_ijst+∑_i∈I∑_w∈W∑_v∈V∑_c∈C∑_k∈K∑_t∈T∑_g∈TE_icwtP_wQW_iwvcktg+∑_i∈I∑_j∈J∑_w∈W∑_v∈V∑_c∈C∑_t∈TCT_jwQFW_ijwvct+∑_i∈I∑_w∈W∑_v∈V∑_c∈C∑_t∈T∑_g∈TE_icwtHW_cI_iwvctg	–
QF_ijst<BIGA_ijst	β1_ijst
∑_w∑_v∑_cQFW_ijwvct≤∑_sa_itQF_ijst-LOSS_ijst	β2_ijt
I_iwvctg=I_iwvct-1g-∑_kQW_iwvcktg	β3_iwvctg
I_iwvctg=∑_jQFW_ijwvct-∑_kQW_iwvcktg	β4_iwvctg
∑_jQFW_ijwvct+∑_gI_iwvct-1g≤W_wAC_iwvct	β5_iwvct
∑_w∑_v∑_c∑_k∑_gQW_iwvcktg≥D_it	β6_it
∑_i∑_w∑_v∑_c∑_gQW_iwcktg≤CV_ky_kt	β7_kt
∑_j∑_sa_it(QF_ijst-LOSS_ijst)≤CF_i	β8_it
∑_j∑_sa_itWEP_ijsQF_ijst≥α_i∑_j∑_sa_itQF_ijst	β9_it

Dual subproblem

The DSP makes the feasible solution space independent from integer solutions and reaches the same objective function value as the SP. Therefore, the DSP is used instead of the SP to obtain upper bounds for the original problem. In this paper, the DSP is formulated as follows (Table 4):

Table 4

Dual variables associated with the sub-problem

Dual variables
β1_ijst≥0∀i∈I,j∈J,s∈S,t∈T;β2_ijt≥0∀i∈I,j∈J,t∈T;β3_iwvctg≥0∀i∈I,w∈W,v∈V,c∈C,t∈T,g∈T;β4_iwvctg∀i∈I,w∈W,v∈V,c∈C,t∈T,g∈T; β5_iwvct≥0∀i∈I,w∈W,v∈V,c∈C,t∈T;β6_it≥0∀i∈I;t∈T;β7_kt≥0∀k∈K,t∈T;β8_it≥0∀i∈I,t∈T;
β9_it≥0∀i∈I,t∈T;

The DSP and its corresponding variables are presented in Table 5.

Table 5

The DSP and its corresponding variables

DSP	Variables
Maxδ=∑_i∑_j∑_s∑_t-BIGA¯_ijstβ1_ijst+∑_i∑_w∑_v∑_c∑_t-W_wAC¯_iwvctβ5_iwvct+∑_i∑_tD_itβ6_it+∑_t∑_k-CV_ky¯_ktβ8_kt+∑_i∑_t-CF_iβ9_it	–
a_it(-β1_ijst+(1-LF_ijs)β2_ijt-(1-LF_ijs)β8_it+(WEP_ijs-α_i)β9_it)≤CQ_ijs	QF_ijst
-β2_ijt+β4_iwvctg-β5_iwvct≤CT_jw	QFW_ijwvct
-β3_iwvctg+β6_it-β7_kt≤P_w	QW_iwvcktgt>g; t-g≤TE¯_iwvcg
-β4_iwvctg+β6_it-β7_kt≤P_w	QW_iwvcktgt=g
-β3_iwvctg+β3_iwvct+1g-β5_iwvct+1≤HW_c	I_iwvctgt>g; t-g≤TE¯_iwvcg
β3_iwvct+1g+β4_iwvctg-β5_iwvct+1≤HW_c	t=g

After solving the DSP, an optimal or feasible cut is added to the MP according to Eq. (24) and an extreme ray of DSP. BD algorithm terminates whenever the stop criterion .where is the feasible solution space of the MP, is a small non-zero number. In addition, UB and LB are the upper and lower bounds obtained from solving DSP and MP, respectively. The pseudocode of BD is given as follows:

Experimental results

This section is presented in three parts. In the first part, numerical test problems, including small, medium, and large-sized instances, are executed for assessing the efficiency and applicability of the proposed model and solution approach. These test problems are randomly generated in Excel and presented in Tables 6 and 7. In the second part, the model is implemented in a case study in Iran to derive valuable insights. Finally, a sensitivity analysis is carried out to help decision-makers determining key parameters which should be carefully tracked in the presented model.

Table 6

Data generation functions

Parameters	Distribution functions
CT_jw	∼Uniform(1,2)
P_w	∼Uniform(1,2)
HW_c	∼Uniform(1,2)
CFT_k	∼Uniform(1,2)
FW_w	∼Uniform(10,20)
W_w	∼Uniform(50,100)
CW_w	∼Uniform(1,2)
CV_k	∼Uniform(5,20)
D_it	∼Uniform(10,100)
CF_i	∼Uniform(50,200)
α_i	∼Uniform(0.6,0.8)

Table 7

Data generation functions for products with different quality grades

s	s1	s2	s3
CQ_ijs	[900,1000)	[1000,1500)	[1500,2000]
LF_ijs	[0.2,0.28]	[0.1,0.2)	[0.01,0.1)

Validation of BD

The presented MILP model and BD algorithm are coded in GAMS 27.2 optimization software and solved by a CPLEX solver. Both CPLEX solver and BD algorithm are executed on a computer with Intel@ Core™ i7-CPU 2.20 GHz and RAM 8.00 GB. Different sizes are presented in Table 8, and numerical results are presented in Table 9, including optimal values of objective functions, the lower and upper bounds of BD, computational time, and gaps between upper and lower bounds of BD and CPLEX solver. BD algorithm continues once it reaches termination criteria which are as follows:

Table 8

Generated small, medium, and large-sized instances

Instances	I	J	W	C	S	V	K
P1	3	3	3	3	3	5	3
P2	4	3	3	3	3	5	3
P3	4	4	4	3	3	5	3
P4	3	8	5	3	3	5	3
P5	5	8	5	3	3	5	3
P6	8	8	5	3	3	5	3
P7	10	8	5	3	3	5	3
P8	10	10	8	3	3	5	3
P9	10	15	10	3	3	5	3
P10	12	15	12	3	3	5	3
P11	12	15	15	3	3	5	3
P12	15	15	20	3	3	10	3
P13	15	20	15	3	3	10	3
P14	20	20	10	3	3	10	3
P15	25	20	10	3	3	10	3
Min	3	3	3	3	3	5	3
Max	25	20	20	3	3	10	3

Table 9

Obtained results from solving small, medium, and large scales instances with GAMS and BD

Size	Instances	CPLEX solver		Benders decomposition algorithm					PRE (%)
		Objective function	CM (s)	LB	UB	CM (s)	Iter	GAP(%)	LB	UB
Small scale	P1	1,692,661	18	1,692,661	1,692,661	31	9	0	0.000	0.000
	P2	1,937,866	25	1,937,512	1,937,893	17	12	0.020	0.018	0.001
	P3	2,014,646	40	2,014,071	2,014,690	20	12	0.031	0.029	0.002
	P4	1,513,379	23	1,509,382	1,513,478	24	11	0.271	0.264	0.007
	P5	2,578,430	35	2,571,032	2,578,488	30	16	0.290	0.287	0.002
Medium-scale	P6	4,488,737	9757	4,476,869	4,489,851	63	19	0.290	0.264	0.000
	P7	5,775,525	28,795	5,773,605	5,787,357	97	19	0.238	0.033	0.205
	P8	5,649,849	> 28,795	5,635,241	5,650,601	495	11	0.273	0.259	0.013
	P9	5,773,698	> 28,795	5,758,896	5,774,913	255	19	0.278	0.256	0.021
	P10	6,652,182	> 28,795	6,652,159	6,653,543	461	22	0.021	0.000	0.020
	P11	6,794,693	> 28,795	6,778,839	6,797,817	621	17	0.280	0.233	0.046
Large scale	P12	Out of memory		8,342,697	8,367,630	2199	26	0.299	–	–
	P13			10,164,640	10,188,150	2555	29	0.231	–	–
	P14			11,570,161	11,600,005	3654	21	0.258	–	–
	P15			13,682,170	13,692,152	3861	34	0.073	–	–

LB lower bound, UB upper bound, s seconds, CM computational time, GAP relative gap between upper and lower bounds, PRE percentage relative gaps between lower or upper bound and optimal values, Iter number of iterations

The number of iterations reaches 200

A relative gap between lower and upper bounds become less than

In addition, the CPLEX solver terminates once the absolute and relative gaps become less than 0.001.

The results in Table 9 point out that the CPLEX solver reaches optimal solutions in less than one minute, whereas the computational time of BD is more than the CPLEX solver in several small sized-instances. As increasing the number of variables in medium sizes, the lower and upper bounds of BD converge to optimal solutions less than an hour, which is remarkably faster than the CPLEX solver. It should be noted, for large sized-instances, the CPLEX solver cannot find optimal solutions and releases an "OUT OF MEMORY" error. In these cases, BD can be utilized to cope with complexity and reach optimal solutions in a reasonable time. Figure 4 demonstrates the changes of computational time for BD and CPLEX solver, in which BD outperforms CPLEX solver, especially in medium and large-sized instances. According to the results, the computational time of the CPLEX solver increases whatever the size of the problem increases. This affirms the efficiency of BD in solving the model in a reasonable time.

Fig. 4

Comparison of the computational time between GAMS and BD solvers

Besides, the average gap between the upper bounds of BD and the CPLEX solver is about 0.01684. It turns out BD converges to optimal solutions with a near-zero gap. In order to verify the efficiency of BD and its convergence, three instances in small, medium, and large sizes are derived from GAMS/CPLEX solver and presented as follows:

Figure 5 shows that BD provides optimality convergence in less than 17 iterations for a small-sized test (no. p5).

Fig. 5

The convergence of lower and upper bounds of BD for a small sized-instance (no. p5)

Figure 6 shows that BD provides optimality convergence in less than 12 iterations for a medium-sized test (no. p8).

Fig. 6

The convergence of lower and upper bounds of BD for a medium sized-instance (no. p8)

Figure 7 demonstrates that BD provides optimality convergence in less than 23 iterations for a large-sized test problem (no. p14).

Fig. 7

The convergence of lower and upper bounds of BD for a large sized-instance (no. p14)

A case study

This section devotes a discussion on an apple supply chain belonging to a domestic distributor in Iran to demonstrate the application of the proposed model. The company decides to establish warehouses at four potential locations, namely Tehran, Fars, Isfahan, and Razavi Khorasan, and activate cold chambers under proper temperature conditions. They often purchase different varieties of apples from producers who truck products by private vehicle fleet and located in Tehran, Alborz, Isfahan, West Azerbaijan, East Azerbaijan, and Semnan provinces. The company purchases good storage apples based on quality characters (i.e., texture, cosmetic appearance) to remain for a while in warehouses and satisfy the wholesaler’s demands in different months of a year. Tables 10, 11, 12, 13, 14, 15 and 16 represent input data related to cold storage construction projects (www.damoon-co.com) and purchasing, transportation, and storage planning for eight months. The average demand is obtained based on historical data when the company faces shortages, and also demands should be met by any variety of apples. In addition, the purchasing cost of apples with different quality grades are available in fruit and vegetable fields organization in Iran, and information about temperature conditions in cold chambers are available in (https://pdfsara.ir/container-storage-services-atmospheric-plan/).

Table 10

Input data related to the varieties of apples

Apples	T_ic			CF_i×10^3Ton	a_it
	c1	c2	c3		t1	t2	t3	t4	t5	t6	t7	t8
Red delicious	1	3	8	3.5	–	–	–	–	1	1	–	–
Golden delicious	1	3	8	3.5	–	–	–	–	1	1	–	–
Golab	1	0	0	4	–	–	1	–	–	–	–	–
Jonathan	1	3	8	3	–	–	–	–	1	1	–	–
Granny Smith	1	3	8	3	–	–	–	–	1	1	–	–
Fuji	1	3	8	3	–	–	–	–	1	1	–	–
Qabaleh	1	3	8	0.7	–	–	–	–	1	1	–	–
Baljei	1	3	8	4	1	–	–	–	–	–	–	–
Soltani	1	3	0	3	–	1	–	–	–	–	–	–
Mashhadi	1	3	0	3.5	–	1	1	–	–	–	–	–
Gardeshiri	1	0	0	3	–	–	–	1		–	–	–
Asadkuhi	1	3	0	3	1	–	–	–	–	–	–	–

Table 11

Operational costs including transportation, packing, and sorting apples

CT_jw($)	Warehouses
producers	Tehran	Fars	Isfahan	Razavi Khorasan
Tehran	12.12	21.76	14.08	21.20
Alborz	12.12	21.76	14.08	21.20
Isfahan	14.36	14.34	12.10	26.01
West Azerbaijan	19.34	29.26	29.18	16.06
East Azerbaijan	6.51	2.48	19.31	26.20
Semnan	21.2	27.14	26.01	5.04

Table 12

Input data related to the vehicles

	k1	k2	k3	k4
CV_k(Ton)	5	10	15	25
CFT_k($)	18.92	25.36	33.2	47.52

Table 13

Input data related to warehouses

	Tehran	Fars	Isfahan	Razavi Khorasan
CW_w(104 $)	33.66	11.00	18.33	29.33
W_w(Ton)	1000	500	500	1000
FW_w	10	6	10	8
P_w(10–4 $))	3.31	0.99	1.66	2.64

Table 14

Holding cost of apples under the different temperature conditions

	c1	c2	c3
HW_c($)	29.4	50.4	159.92

Table 15

Loss factors for different types of apples with different quality grades and minimum acceptable quality score based on expert’s opinion

LF_ijs(%)	s1						s2						s3						α_i
	j1	j2	j3	j4	j5	j6	j1	j2	j3	j4	j5	j6	j1	j2	j3	j4	j5	j6
Red delicious	22	22	21	21	21	21	14	12	12	14	17	13	9	9	7	9	6	5	0.634
Golden Delicious	20	24	25	23	22	20	13	19	17	10	14	11	6	6	7	9	8	6	0.769
Golab	23	21	20	24	24	22	15	13	17	13	10	11	7	5	6	7	5	5	0.710
Jonathan	22	24	26	24	20	21	15	11	13	11	12	11	5	8	7	8	8	5	0.660
Granny Smith	21	24	25	22	25	23	14	12	15	17	19	18	7	6	7	7	5	6	0.658
Fuji	21	21	23	22	23	21	14	16	10	13	11	15	6	7	7	6	6	7	0.645
Qabaleh	24	20	24	26	24	22	16	14	10	15	17	15	9	5	9	6	8	5	0.670
Baljei	21	22	21	20	23	20	14	13	15	11	12	11	5	6	8	7	8	8	0.771
Soltani	24	22	23	26	24	22	15	18	17	16	15	18	8	9	8	6	5	7	0.613
Mashhadi	23	22	26	25	25	21	13	15	19	12	10	16	6	9	6	9	9	7	0.700
Gardeshiri	24	22	23	25	25	23	14	19	11	16	14	13	7	8	5	5	8	6	0.800
Asadkuhi	23	23	23	26	24	24	13	15	16	15	14	15	5	7	7	6	7	5	0.716

Table 16

Demand for different types of apples

D_it×10^3(Ton)	t1	t2	t3	t4	t5	t6	t7	t8
Red delicious	–	–	–	–	1.000	0.900	1.560	1.400
Golden Delicious	–	–	–	–	0.650	0.780	1.000	1.080
Golab	–	–	2.530	1.438	–	–	–	–
Jonathan	–	–	–	–	1.350	1.246	1.321	1.350
Granny Smith	–	–	–	–	1.630	1.630	2.430	2.150
Fuji	–	–	–	–	1.450	1.525	2.350	2.096
Qabaleh	–	–	–	–	0.554	0.324	0.268	–
Baljei	1.500	1.000	1230	1.551	1.470	1.580	–	–
Soltani	–	1.450	1375	–	–	–	–	–
Mashhadi	–	1.250	1250	1.235	–	–	–	–
Gardeshiri	–	–	–	1.542	2.235	–	–	–
Asadkuhi	1.400	1.245	1.115	–	–	–	–	–
D_t×10^3(Ton)	2.900	4.945	7.500	5.766	10.0339	7.985	8.929	8.076

The quality score of fruits supplied by different producers () is calculated based on the minimum and maximum values of loss factor for a product with quality grade s as follows:where is the minimum loss factor for product with quality grade s and is the maximum loss factor for product with quality grade . Equation (26) illustrates the quality score of producer supplying product with quality grade , which equals 1.

Numerical results and discussion on a case study

The model has 184,572 discrete variables, 330,049 continuous variables, and 642,748 equations for the case study and is implemented by both CPLEX solver and BD. Herein, the efficiency of the proposed model and related solutions have been discussed. The lower and upper bounds of BD become 25,289,297 and 25,289,302, respectively. In addition, the optimal solution is obtained by the CPLEX solver, which is 25,289,302. The results indicate that BD converges to the optimal solution in a reasonable time. Figure 8 shows the convergence of BD for the case study.

Fig. 8

The convergence of the lower and upper bounds of BD for the case study

First, producers who gain good scores are selected to supply products. The results show a of apples are of first quality grade, are of second grade, and are of third quality grade. In this respect, the total cost of purchasing apples with different quality grades is about US$ 996,000.

According to the results, it is beneficial to open warehouses in Tehran, Fars, and Isfahan provinces at a total cost of US$ 6,600,000. In this case, the operational costs in warehouses become US$ 147. It is worth mentioning that the presented model is able to show the temperature condition in each cold chamber. Three temperature conditions are utilized for maintaining apples with different shelf life. Accordingly, 73.554% 19.008%, and 7.438% of products are kept under the temperature levels , , and , respectively at a total cost of US$ 17,520,000. In addition, operational costs, including laborer, packing, and transportation costs from producers to warehouses are about US$ 67,422. Transportation cost between warehouses to the wholesaler is about US$ 107,349 for 2262 trips, of which 1, 1, 3, and 2256 are performed by vehicles types 1, 2, 3, and 4, respectively.

Sensitivity analysis

Sensitivity analysis is performed on the case study to determine the behavior of the model by changing the main parameters. Clearly, there are many reasons affecting the selling price of products and may lead to fluctuation in demands. Changes in the selling price of products may be caused by specific economic conditions (e.g., economic impacts of COVID-19 outbreak), irreversible damage during production (e.g., frost, floods), etc. The customers’ demand may be triggered by these factors or changes in taste and expectation, thereby becoming lesser or higher. From this point of view, it is vital to assess the model behavior during changes in parameters and . In addition, this paper gives four vehicle types for fruits transportation and two cold storage rooms for 500 or 1000 Tons of apples that may change. The containers of different capacities can easily fit onto different forms of trucks and provide other options for transportation planning. Besides, redesign of cold rooms of less capacity may affect the feasibility of the problem. Finally, two remaining parameters are inventory holding costs and operational costs, which may be changed by embracing new technologies, identifying inefficiencies and unused services, etc. Therefore, this section provides an overview of the changes and their impacts to revise the plan after the occurrence. The results are summarized in Table 17.

Table 17

The results of sensitivity analysis

Parameters	Changes (%)	Objective function values ($)	Percentage changes (%)
HW_c	50	34,048,377	34.63
	25	29,668,828	17.31
	− 25	20,909,811	− 17.31
	− 50	16,530,329	− 34.63
W_w	50	23,458,018	− 7.24
	25	24,281,249	− 3.98
	− 25	27,137,910	7.31
	− 50	Infeasible	−
D_it	50	Infeasible	−
	25	Infeasible	−
	− 25	15,945,152	− 36.94
	− 50	10,358,746	− 59.03
CW_w	50	28,589,282	13.04
	25	26,939,292	6.52
	− 25	23,639,312	− 6.52
	− 50	21,989,421	− 13.04
CQ_ijs	50	25,785,923	1.96
	25	25,537,897	0.98
	− 25	25,040,161	− 0.98
	− 50	24,790,887	− 1.97
CV_k	50	25,253,548	− 0.14
	25	25,267,852	− 0.08
	− 25	25,325,095	0.14
	− 50	25,396,595	0.42

The results in Table 17 draw up the sensitivity of the model to parameters and , where an increase of 25% in parameter makes the model infeasible, and a decrease of 50% in parameter reduce the total cost by 59%. Indeed, it is necessary to adjust capacity to match it with the changes in demands. We suppose that the available capacity of warehouses or producers is not sufficient in response to increased demands, therefore, the model becomes infeasible. In this case, the capacity of either warehouses or producers may step up to meet the demand. In addition, it is expected to purchase fewer products and avoid overstocks when demand decreases. This may cause a decrease in purchasing costs, operational costs, and inventory holding costs. Besides, a 50% increase in leads to a 34.63% increase in the total cost. This indicates the importance of temperature condition adjustment and its impacts on the total cost.

For further investigation, the effects of the changes in parameters , , on operational costs, inventory holding costs, and purchasing costs are presented as follows:

Figure 9 demonstrates the effects of the changes on operational costs.

Fig. 9

Impacts of the changes on operational costs

The changes in parameter have insignificant opposite impacts on the operational costs.

The changes in parameter have the greatest impacts on the operational costs, and a 25% decrease in demand leads to an 18% decrease in operational costs.

Figure 10 shows the impacts of the changes on inventory holding costs.

Fig. 10

Impacts of the changes on inventory holding costs

A decrease in parameter significantly reduces inventory holding costs. Indeed, a decrease of 25% in leads to a 25% decrease in inventory holding costs.

An increase in parameter significantly affects the inventory holding costs, and a 25% decrease in leads to a 36% decrease in inventory holding costs.

Figure 11 shows the impacts of the changes on purchase costs.

Fig. 11

Impacts of the changes on purchase costs

Changes in parameter have the greatest impacts on the total purchase cost. Indeed, a 25% decrease in leads to a 27% decrease in total purchase cost.

Total purchase cost increases by increasing the amount of parameter .

Conclusion and future directions

Nowadays, an appropriate configuration of agro-supply chains becomes all the more critical issues owing to specific characteristics of agro-products. Hence, in this paper, a multi-product multi-period MINLP model has been developed for an agro-supply chain contributing to purchase, facility location, and transportation decisions. Forasmuch as adjustment of temperature condition affects the quality of products and consequently their shelf lives, this study determines proper temperature conditions based on the varieties of products. This model aims to minimize the total cost of the supply chain considering several operational constraints.

The primary step to solve the proposed model is to linearize the MINLP model and reformulate it into a MILP problem. Then, three sets of numerical examples, small, medium, and large sizes, have been presented to evaluate the performance of the mathematical model. The problem was coded in GAMS software and executed for small and medium test problems. Whatever the size of the problem increases, the computational time increases and the problem becomes difficult to solve. Therefore, BD algorithm has been implemented to tackle the complexity of the problem. Then GAMS/CPLEX solver and BD were compared in terms of objective function values and computational time. The results depict the convergence of BD to optimal solutions in a reasonable time.

Afterward, a precise assessment of the model was carried out by implementing a real case study. The results turn out the optimal purchase decisions of required quality score to be kept for a long time and number of cold chambers under proper temperature conditions. In addition, a sensitivity analysis was performed on the main parameters to predict the effects of the variations.

For future research, first, the model can be revisited for conditions under uncertain demands, thereby proposing efficient scenarios to overcome their associated shortage. Another possible extension of the model is to analyze purchase decisions under discount. This paper assumes only the products of the required quality for long time storage should be purchased. However, a new framework is promising to accommodate the distribution of products with short and long shelf lives that need precise inventory management.

Abbreviation	Description
MINLP	Mixed-integer nonlinear programming
BD	Benders decomposition
MILP	Mixed-integer linear programming
MP	Master problem
SP	Subproblem
DSP	Dual sub-problem
SLD	Decision on shelf life
CFD	Decision on cold storage facilities