Solution of the Capacitated Lot-Sizing Problem with Remanufacturing (CLSPR) in a General Way with the Help of Simulation and Relaxation

1. Introduction

2. Problem Description

• The demand for new products and remanufactured products are separate and backlog is not allowed.
• The manufacturing capacity is sufficient to meet the demands in each period, in particularly we have: the capacity can satisfy the demands for new products and remanufactured products simultaneously, i.e.

  $${\displaystyle \sum _{j=1}^t(d_j^n+d_j^r)\le \sum _{j=1}^t{c}_j,{\displaystyle \forall t=1,\cdots ,T}}$$

  (1)
• Initial and end inventory stocks are zero, ¹ i.e.

  $$\begin{array}{ccc}\hfill {\displaystyle i_0^n=0,}& {\displaystyle i_0^r=0,}\hfill & \hfill {\displaystyle i_0^{re}=0}\\ \hfill {\displaystyle i_T^n=0,}& {\displaystyle i_T^r=0}\hfill & \hfill \end{array}$$

  (2)
  The demand will be fully satisfied if the final inventories are zero.
• The quantity of returned products can satisfy the demand for remanufactured products i.e.

  $${\displaystyle \sum _{j=1}^t{d}_j^r\le \sum _{j=1}^t{r}_j^{re}{\displaystyle \forall t=1,\cdots ,T}.}$$

  (3)
• In economic terms, inventory holding cost of returned products is less than that of remanufactured products.

  $${\displaystyle \sum _{j=t}^T{h}_j^{re}\le \sum _{j=t}^T{h}_j^r,{\displaystyle \forall t=1,\cdots ,T}}$$

  (4)
  This hypothesis can be found in [39] too.

3. Model A (Relaxation)

4. Model B (Simulation)

4.1. Low Discrepancy Sequences

The generation of random numbers with a computer is not possible Knuth [41]. As John von Neumann said: Any one who considers arithmetical methods of producing random digits is, of course, in a state of sin, [42]. An excellent overview of the methods of generating pseudo-random numbers are available [e.g. [41,43,44].

In this section we explain the number-theoretical concept of discrepancy. Then, we introduce the Halton sequence which is probably the easiest low-discrepancy ² number generation method to describe.

Definition 1. 

Let $\{z_1,\cdots ,z_N\}$ a sequence of real numbers with $0<z_i<1,i=1,\cdots ,N$. The discrepancy $D_N$ for the sequence is defined as

$$D_N=\underset{l}{sup}\left|S_N\left(l\right)-N\left|l\right|\right|$$

(34)

where l is any subinterval $[a,b)\subseteq [0,1],$$\left|l\right|=b-a,$ and $S_N\left(l\right)$ denotes the number of elements of the sequence, that belongs to the interval $l.$

A measure for how a sequence of real numbers $\{z_1,\cdots ,z_N\}$ , $a<z_i<b,i=1,\cdots ,N$ is equidistributed on an interval $[a,b]$ is the discrepancy $D_N$. Low-discrepancy sequences, also known asquasirandomsequences, are numbers that are better equidistributed in a given volume than pseudo-random numbers.

Remark 3. 

A sequence $\{z_1,\cdots ,z_N\},0<z_i<1,i=1,\cdots ,N$ of real numbers is said equidistributed on the interval $[0,1]$ if $D_N=o\left(N\right),N\to \infty ,$ [35].

The ( $\mathbb{QRN}$ ) of Halton, Sobol and Niederreiter have a low discrepancy $D_N=O(ln\left(N\right)/N)$. While pseudorandon sequences have a discrepancy $D_N=O(1/\sqrt{N}),$ [45]. Figure 2 uses two-dimensional projection of a pseudorandom sequence and of a low-discrepancy (Halton) sequence to demonstrate the fundamental difference between the two classes of sequences.

The desirable properties of a sequence of this ( $\mathbb{QRN}$ ) may be summarized as follows see [37]:

• the least period length should be sufficiently large,
• it should have littie intrinsic structure (such as lattice structure),
• it should have good statistical properties,
• the algorithm generating the sequence should be reasonably efficient.

It’s easy to generate sequences of Halton with the following Algorithm 1.

The following Halton sequences of Table 5 are constructed according to Algorithm 1 that uses a prime number as its base.

Remark 4. 

To generate the n-th Halton point in a sequence consider the base $b-$ary expansion of a $n={\sum }_{i=0}^{\infty }{a}_i{b}^i$ where the $b-$ary coefficients $a_i\in \{0,\cdots ,b-1\}.$ Then the n-th Halton point is $H\left(n\right)={\sum }_{i=0}^{\infty }{a}_i{b}^{-i-1}.$ It’s easy to build a Halton-sequence with the following observation: If $a_0<b-1$ then $H(n+1)=H\left(n\right)+1/b$ else if $a_0=b-1$ then $H(n+1)=H\left(n\right)-(1-b^k-b^{k+1})$ where $k=min\{i\ge 0:a_i\ne b-1\}$ (details see [38]). This method is very efficient and will be used in this paper.

4.2. Notation

We use the following notation ${\displaystyle \forall t=1,\cdots ,T}$

Name	Meaning
${\displaystyle d_{\left[t\right]}^n}$ =	${\displaystyle \sum _{j=1}^t{d}_j^n}$	${\displaystyle d_{\left[0\right]}^n}$ =	0
${\displaystyle d_{\left[t\right]}^r}$ =	${\displaystyle \sum _{j=1}^t{d}_j^r}$	${\displaystyle d_{\left[0\right]}^r}$ =	0
${\displaystyle r_{\left[t\right]}^{re}}$=	${\displaystyle \sum _{j=1}^t{r}_j^{re}}$	${\displaystyle r_{\left[0\right]}^{re}}$=	0
${\displaystyle x_{\left[t\right]}^n}$ =	${\displaystyle \sum _{j=1}^t{x}_j^n}$	${\displaystyle x_{\left[0\right]}^n}$ =	0
${\displaystyle x_{\left[t\right]}^r}$ =	${\displaystyle \sum _{j=1}^t{x}_j^r}$	${\displaystyle x_{\left[0\right]}^r}$ =	0

The notation $x\in {[a,b]}_p$ means $x=(b-a)z+a,a\le b,0\le z\le 1$ and the number z is simulated with the following distribution

$$\begin{array}{cc}z=\hfill & \left(\begin{array}{ccc}{z}_h& 0& 1\\ p& q& q\end{array}\right)\hfill \\ \hfill & 0\le p\le 1,q={\displaystyle \frac{1}{2}}(1-p),z_{h}isaHaltonnumber\hfill \end{array}$$

(35)

We generate a pseudorandom number $0\le g\le 1$ to decide, that values take z, see Algorithm 2.

In Model B (Section 4) we will investigate the class of problems with the condition

$$d_t^n+d_t^r=d_t\le c_t,\forall t$$

(36)

If this condition is not satisfied, the problem is easily solved with Model A (Section 3).

Algorithm 2:Simulation of z with distribution (35)

4.5. Simulation of the Objective Function

Let R be a matrix with Halton’s $\mathbb{QRN}$

$$R=\left(\begin{array}{ccc}R(1,1)& \cdots & R(1,T)\\ R(2,1)& \cdots & R(2,T)\end{array}\right)$$

We simulate the production starting in period t=1. The calculation of the objective function is carried out with Algorithm 3 using proposition 1, 2 and 3.

Algorithm 3:Calculation of the objective function

Further information on the complexity-theoretic approach to randomness can be found in [45,46] and [47].

Remark 7. 

 

• The generation of the matrix $R[N\times T]$ requires $O(N\times T)$ operations.
• Die evaluation of the function $\phi (x^n,x^r)$ with Algorithm 3 requires $O\left(T\right)$ operations

Then, the computational complexity of the simulation with N points and T periods is $O(N\times T).$

5. Numerical Experiments

5.3. Model B

[30] have shown that several families of CLSP are NP-hard. For the construction of the NP-hard instances (2200 instances) we follow the findings of [31]. They use the following notation $Nr/\alpha /\beta /\gamma /\sigma$, where $Nr,\alpha ,\beta ,\gamma$ and σ specify respectively the number of items, a special structure for the setup costs, the holding costs, production costs, and capacities. In this paper [31] show that the following class 2/C/G/A/C is NP-hard. For this reason we have created 2200 instances, where the set-up costs per product and capacities per instance are constant. The holding costs do not necessarily follow a specified pattern, the production costs can be chosen arbitrarily. The maximum calculation time for Gurobi is 600 seconds. The heuristic operates according to the number of simulations, which gradually increases with the number of periods. The largest group T300 uses 130 seconds. The parameters for generating data sets (see Table 8) use the following notation $x\in [a;b]\iff x=(b-a)\theta +a,a\le b,$ where $0\le \theta \le 1$ is a random number, that means the values are uniformly distributed on the interval $[a,b].$

The important assumption in model B is: the capacities for each TI is constant, the set-up costs in each TI and for each product are constant. These parameters vary between a minimum and a maximum depending on the T parameter.

The capacities for each TI vary according to the parameter T. For example if $T=15$ the capacities are between 600 and 800. If $T=300$ the capacities vary between 3000 and 5500. Analogously the other parameters. The remaining parameters according to Table 8.

For the simulation we used following parameters:

$$N_T=2^{15}T$$

$N_T$ is the Halton’s numbers used for $T\in \{15,30,\cdots ,300\}$ periods.

All instances for $T=15,30\cdots ,300$ have the same schema Step 1 until Step 5. We used $\tau =4$ (see Table 9).

Algorithm 4:Heuristic: blind search

5.3.1. Results of Model B

We generated 200 random instances for each problem class (PC) and with a Box-plot we compared the feasible solutions $G_{Ti}(x^n,x^r)$ of the problem (16)-(24 ) found by Gurobi 9.0.3 with the solution $S_{Ti}(x^n,x^r)$ of the simulation presented in this paper and clearly see the similarity of the results found (see Figure A4 and Figure A5).

With the following notation we present the results (see Table A3).

T	$\in \{15,30,\cdots ,300\}$
m	$=200$
$S_{Ti}(x^n,x^r)$	$i=1,\cdots ,m$
$G_{Ti}(x^n,x^r)$	$i=1,\cdots ,m$
${\mu }_T$	$={\displaystyle \frac{{\sum }_{i=1}^m{S}_{Ti}(x^n,x^r)}{m}}$
${\mu }_T^{*}$	$={\displaystyle \frac{{\sum }_{i=1}^m{G}_{Ti}(x^n,x^r)}{m}}$
$CPU\_S_T$	$={\displaystyle \frac{{\sum }_{i=1}^{m}CP{U}_{Ti}}{m}}$ with Simulation
$CP{U}_T^{*}$	$={\displaystyle \frac{{\sum }_{i=1}^{m}CP{U}_{Ti}^{*}}{m}}$ with Gurobi

In Table A3 the relative error for Total costs (Tc) and CPU-time for every problem class was calculated as

$$RelativeAvg.Error\left(Tc\right)={\displaystyle \frac{{\mu }_T-{\mu }_T^{*}}{{\mu }_T^{*}}}$$

$$RelativeAvg.Error\left(CPU\right)={\displaystyle \frac{CPU\_S_T-CP{U}_T^{*}}{CP{U}_T^{*}}}$$

The graphical comparison is shown in Figure A3.

This is exactly how we calculated the relative errors in inventory cost and setup cost. Attached in Appendix B are the results of the average total cost, average CPU time, average Inventory costs. Figure A3 (average CPU time) clearly shows that the Gurobi admissible solution for the PC from T150 to T300 are not optimal and that the CPU of the simulation is much faster.

The simulation in average determines the setup costs and return inventory costs always higher than Gurobi. On the other hand, the inventory costs of Gurobi are higher than the simulation. What can we say about the quality of the solution of the problems? The simulation could not give a better solution than Gurobi’s solution. Gurobi solved the (PC) problems up to T120 in an optimal way. The problems from T150 to T300 were not solved optimally by Gurobi, since it would take too much time due to the NP-hard category of the problem. The simulation found feasible solutions much faster than Gurobi. Here is the advantage of the simulation, the simplicity of its implementation and the speed in finding an acceptable solution.

6. Conclusions and Outlook

References

1. Manne, A.S. Programming of economic lot sizes. Management science 1958, 4, 115–135. [Google Scholar] [CrossRef]
2. Wagner, H.M.; Whitin, T.M. Dynamic version of the economic lot size model. Management science 1958, 5, 89–96. [Google Scholar] [CrossRef]
3. Beltrán, J.L.; Krass, D. Dynamic lot sizing with returning items and disposals. IIe transactions 2002, 34, 437–448. [Google Scholar] [CrossRef]
4. Elmaghraby, S.E. The Economic Lot Scheduling Problem (ELSP): Review and Extensions. Management Science 1978, 24, 587–598. [Google Scholar] [CrossRef]
5. Brahimi, N.; Dauzere-Peres, S.; Najid, N.M.; Nordli, A. Single item lot sizing problems. European Journal of Operational Research 2006, 168, 1–16. [Google Scholar] [CrossRef]
6. Buschkühl, L.; Sahling, F.; Helber, S.; Tempelmeier, H. Dynamic capacitated lot-sizing problems: a classification and review of solution approaches. Or Spectrum 2010, 32, 231–261. [Google Scholar] [CrossRef]
7. Thierry, M.; Salomon, M.; Van Nunen, J.; Van Wassenhove, L. Strategic issues in product recovery management. California management review 1995, 37, 114–136. [Google Scholar] [CrossRef]
8. Amin, S.H.; Zhang, G.; Eldali, M. A review of closed-loop supply chain models. Journal of Data, Information and Management 2020, 2, 279–307. [CrossRef]
9. Agrawal, S.; Singh, R.K.; Murtaza, Q. A literature review and perspectives in reverse logistics. Resources, Conservation and Recycling 2015, 97, 76–92. [CrossRef]
10. Ayres, R.; Ferrer, G.; Van Leynseele, T. Eco-efficiency, asset recovery and remanufacturing. European Management Journal 1997, 15, 557–574. [Google Scholar] [CrossRef]
11. Inderfurth, K. Optimal policies in hybrid manufacturing/remanufacturing systems with product substitution. International Journal of Production Economics 2004, 90, 325–343. [Google Scholar] [CrossRef]
12. Tang, O.; Teunter, R. Economic lot scheduling problem with returns. Production and Operations Management 2006, 15, 488–497. [Google Scholar] [CrossRef]
13. Teunter, R.; Kaparis, K.; Tang, O. Multi-product economic lot scheduling problem with separate production lines for manufacturing and remanufacturing. European journal of operational research 2008, 191, 1241–1253. [Google Scholar] [CrossRef]
14. Sahling, F. A Column-Generation Approach for a Short-Term Production Planning Problem in Closed-Loop Supply Chains. BuR- Business Research 2013, 6(1), 55–75. [CrossRef]
15. Cunha, J.O.; Konstantaras, I.; Melo, R.A.; Sifaleras, A. On multi-item economic lot-sizing with remanufacturing and uncapacitated production. Applied Mathematical Modelling 2017, 50, 772–780. [Google Scholar] [CrossRef]
16. Sifaleras, A.; Konstantaras, I. Variable neighborhood descent heuristic for solving reverse logistics multi-item dynamic lot-sizing problems. Electronic Notes in Discrete Mathematics 2015, 47, 69–76. [Google Scholar] [CrossRef]
17. Cunha, J.O.; Kramer, H.H.; Melo, R.A. Effective matheuristics for the multi-item capacitated lot-sizing problem with remanufacturing. Computers & Operations Research 2019, 104, 149–158. [CrossRef]
18. Helmrich, R.; Jans, M.; van den Heuvel, W.; Wagelmans, A. Economic lot-sizing with remanufacturing: Complexity and efficient formulations. IISE Transactions 2014, 46(1), 67–86. [Google Scholar] [CrossRef]
19. Kilic, O.; van den Heuvel, W. Economic lot sizing with remanufacturing: Structural properties and polynomial-time heuristics.IISE Transactions 2019,51:12, 1318–1331. [CrossRef]
20. Richter, K.; Weber, J. The reverse Wagner/Whitin model with variable manufacturing and remanufacturing cost. International Journal of Production Economics 2001, 71, 447–456. [Google Scholar] [CrossRef]
21. Richter, K.; Sombrutzki, M. Remanufacturing planing for the reverse Wagner/Whitin models. Journal of Operational Research 2000, 121, 304–315. [Google Scholar] [CrossRef]
22. Teunter, R.; Bayindir, Z.; van den Heuvel, W. Dynamic lot sizing with product returns and remanufacturing. Int. Journal of Production Research 2006, pp. 4377–4400. [CrossRef]
23. Schulz, T. A new Silver-Meal basic heuristic for the single-item dynamic lot sizing problem with returns and remanufacturing. International Journal of Production Research 2011, p. 2519–2533. [CrossRef]
24. Cunha, J.; Melo, R. A computational comparison of formulations for the economic lot-sizing with remanufacturing. Computers & Industrial Engineering 2016, 92, 72–81. [CrossRef]
25. Piñeyro, P.; Viera, O. The economic lot-sizing problem with remanufacturing and heterogeneous returns: formulations, analysis and algorithms. International Journal of Production Research 2022, 60, 3521–3533. [Google Scholar] [CrossRef]
26. Lahmar, H.; Dahane, M.; Mouss, K.; Haoues, M. Multi-objective production planning of new and remanufactured products in hybrid production system. IFAC-PapersOnLine 2022, 55, 275–280. [CrossRef]
27. Zouadi, T.; Yalaoui, A.; Reghioui, M. Lot sizing and pricing problem in a recovery system with returns and one-way substitution option: Novel cost benefit evaluation based approaches. IFAC-PapersOnLine 2019, 52, 36–41. [CrossRef]
28. Şenyiğit, E.; Erol, R. New lot sizing heuristics for demand and price uncertainties with service-level constraint. International Journal of Production Research 2010, 48, 21–44. [Google Scholar] [CrossRef]
29. Zhang, Z.H.; Jiang, H.; Pan, X. A Lagrangian relaxation based approach for the capacitated lot sizing problem in closed-loop supply chain. International Journal of Production Economics 2012, 140, 249–255. [Google Scholar] [CrossRef]
30. Florian, M.; Lenstra, J.K.; Rinnooy Kan, A. Deterministic production planning: Algorithms and complexity. Management science 1980, 26, 669–679. [Google Scholar] [CrossRef]
31. Bitran, G.R.; Yanasse, H.H. Computational complexity of the capacitated lot size problem. Management Science 1982, 28, 1174–1186. [Google Scholar] [CrossRef]
32. Dixon, P.S. Multi-Item Lot-Sizing with Limited Capacity. dissertation, University of Waterloo, Ontario, 1979.
33. Dziuba, D.; Almeder, C. New construction heuristic for capacitated lot sizing problems. European Journal of Operational Research 2023. [Google Scholar] [CrossRef]
34. Maes, J.; McClain, J.; Van Wassenhove, L. Multilevel capacitated lotsizing complexity and LP-based heuristics. Journal of Operational Research 1991, 53, 131–148. [Google Scholar] [CrossRef]
35. Sobol’, I.M. Calculation of improper integrals using uniformly distributed sequences. Soviet Math. Dokl. 1973, 14, 734–738. [Google Scholar]
36. Halton, J.H. On the efficiency of certain quasi-random sequences of points in evaluating multi-dimensional integrals. Numerische Mathematik 1960, 2, 84–90. [Google Scholar] [CrossRef]
37. Niederreiter, H. Quasi-Monte Carlo methods and pseudo-random numbers. Bull. Am. Math. Soc. 1978, 84, 957–1041. [Google Scholar] [CrossRef]
38. Klinger, B. Numerical integration of singular integrands using low-discrepancy sequences. Computing 1997, 59, 223–236. [Google Scholar] [CrossRef]
39. Van den Heuvel, W.; Wagelmans, A.P. Four equivalent lot-sizing models. Operations Research Letters 2008, 36, 465–470. [Google Scholar] [CrossRef]
40. Maes, J.; Van Wassenhove, L. A simple heuristic for the multi item single level capacitated lotsizing problem. Operations research letters 1986, 4, 265–273. [Google Scholar] [CrossRef]
41. Knuth, D.E. Art of computer programming: Seminumerical algorithms; Vol. 2, Addison-Wesley Professional, 2014.
42. von Neumann, J. Various techniques used in connection with random digits. Monte Carlo Method, Appl. Math. Series 1951, 12, 36–38.
43. Niederreiter, H. Random number generation and quasi-Monte Carlo methods; SIAM, 1992.
44. Tezuka, S. Uniform random numbers: Theory and practice; Vol. 315, Springer Science & Business Media, 2012.
45. Chaitin, G.J. Information, Randomness and Incompleteness. Papers on Algorithmic Information Theory. 1987, p. 236. , p. 236.
46. Kolmogorov, A.N.; Uspenskii, V.A. Algorithms and randomness. Theory of Probability & Its Applications 1988, 32, 389–412. [CrossRef]
47. Schnorr, C. Zufälligkeit und Wahrscheinlichkeit. Lecture Notes in Math. 1971, 218, 109. [Google Scholar]