Articles in Advance, pp. 1–19 ISSN 0041-1655 (print) — ISSN 1526-5447 (online) http://dx.doi.org/10.1287/trsc.2015.0586 © 2015 INFORMS Small is Beautiful: A Framework for Evaluating and Optimizing Live-Cube Compact Storage Systems Nima Zaerpour Faculty of Economics and Business Administration, VU University Amsterdam, 1081 HV Amsterdam, The Netherlands,

[email protected]

Yugang Yu School of Management, University of Science and Technology of China, Hefei 230026, China,

[email protected]

René B. M. de Koster Rotterdam School of Management, Erasmus University, 3000 DR Rotterdam, The Netherlands,

[email protected]

W arehouses occupy much space and land, which has become increasingly scarce in many parts of Europe, Asia, and the United States, particularly close to areas where demand is generated, such as large cities. This paper studies live-cube compact storage systems that may solve this space shortage problem as they do not require travel aisles. Each stored unit load is accessible individually and can be moved in x and y directions by a shuttle as long as an empty location is available, comparable to the well-known 15-puzzle in which 15 numbered tiles slide within a 4 × 4 grid. When multiple empty locations are available on a level, the shuttles can cooperate to create a virtual aisle for fast retrieval of a desired unit load. A lift moves the unit loads across different levels in z direction. Such storage systems are increasingly used in different service sectors like car parking, warehousing, and container handling, but so far they have hardly been studied. For live-cube systems, many research questions still have to be answered, including cycle time calculations, cost comparisons, and energy requirements. In this paper, we first derive simple to use closed-form formulas for expected retrieval time of an arbitrary unit load and validate the quality of these formulas by comparing them with a real application. Second, we propose and solve a mixed-integer nonlinear model to optimize system dimensions by minimizing the retrieval time. We obtain closed-form expressions for minimum retrieval time that are simple to apply in practice. Third, we compare the investment, operational costs, and energy consumption of live-cube systems with traditional systems based on a real application. Keywords: logistics; material handling; service sectors; live-cube compact storage system; response time History: Received: February 2013; revision received: July 2014; accepted: September 2014. Published online in Articles in Advance. 1. Introduction than traditional low-bay warehouses) because they Warehouses form vital nodes in any supply chain, require transport aisles between any two racks. These decoupling supply from demand, providing an as- travel aisles may consume about 35% of all stor- sortment rapidly accessible for customers; allowing age space and contribute to high building costs (Gue grouping of flows and products, thereby reducing 2006). Space scarcity has forced many companies to transport costs; and creating manufacturing post- look for even more compact storage systems. The dis- ponement options, thereby reducing inventory costs. advantage of many of these systems, e.g., deep-lane With increased trade, in many regions the number storage, is that stored loads are no longer individually of warehouses and the space used for warehous- accessible without reshuffling loads stored in front ing has grown rapidly. Unfortunately, warehouses of it. This increases response times of these systems, occupy much land and infrastructure, and land has making them primarily fit for companies with large become increasingly scarce in many parts of Europe, stocks per product so that reshuffling can be avoided. Asia, and the United States, particularly in areas Typical applications can be found in manufacturing with major customer concentrations, such as large environments or in refrigerated storage. cities. It is possible to construct high-bay warehouses, In this paper, we propose a system that may solve which are usually automated (so-called AS/R sys- this space problem, with stored loads that are indi- tems), but these still consume much space (albeit less vidually accessible, while avoiding reshuffling. We 1 Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 2 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS have been inspired by recently developed automated load’s storage level. In case of a single empty location, car park systems, where each car is stored on a moving a requested unit load to the depot is com- shuttle that all can move autonomously and simul- parable to solving a “15-puzzle,” which consists of taneously in both x and y directions and vertical 15 numbered tiles in a random order in a 4 × 4 frame transport is carried out by lifts (an example can be (Slocum and Sonneveld 2006). We assume the I/O found in the last part of the movie MI4). In fact, point is located at the lower left corner of the system many system manufacturers are currently developing (at the ground level) and we furthermore assume a such so-called “live-cube” storage systems in differ- random storage policy is used, which implies all stor- ent variants with a growing number of implementa- age locations have an equal chance of serving a stor- tions in automated parking systems, warehousing and age job. Random storage requires the least data since cross-docking, and container handling (e.g., Hyundai no product information is used in determining the Elevator Co. Ltd. 2013; Automotion Parking Systems storage location and it is broadly studied in the litera- 2013; Eweco 2014; Wöhr 2013; OTDH 2013; Navstors ture (e.g., Hausman, Schwarz, and Graves 1976; Bozer 2013; Navpak 2013; EZ-Indus 2013). However, many and White 1984; Goetschalckx and Ratliff 1990; Lee research questions still have to be answered, includ- and Elsayed 2005; de Koster, Le-Duc, and Yu 2008). ing cycle time calculations for given sizes and storage In addition, random storage is frequently used as a strategies; in particular, cost comparisons have to be benchmark in comparison for system response time made. In addition, buildings and operations have to with other storage policies (e.g., Hausman, Schwarz, become more sustainable, so energy consumption is and Graves 1976; Lee and Elsayed 2005). also a point of consideration. Our first contribution is the derivation of expected We address these issues in this paper: retrieval time formulas for any system configuration. 1. What is the performance of a live-cube storage Building and land restrictions might force managers system of a given size and configuration? to evaluate the performance of nonoptimal system 2. What is the optimal configuration of a live-cube storage system that minimizes the response time, configurations. We obtain the expected retrieval time given a required storage capacity? formulas for an arbitrary retrieval corresponding 3. Is a live-cube system better (in terms of invest- to all possible (16 in total) system configurations ment and operational costs and environment) than a through a comprehensive decomposition procedure. traditional system? What is the impact of system con- To simplify the problem analysis, we assume the coor- figuration on the investment and operational costs of dinates of the system to be continuous. The assump- such systems? tion of continuous dimensions is commonly made in The main components of the live-cube system are other papers (e.g., Hausman, Schwarz, and Graves multiple levels of storage grids, shuttles, a lift, and a 1976; Yu and de Koster 2009). In addition, we con- depot (or input/output (I/O) point) as shown in Fig- sider a system with multiple empty locations on each ure 1. Shuttles can move in x and y directions while level. The shuttles can then cooperate to create a vir- carrying a unit load. Each unit load holds only one tual aisle for fast retrieval of a requested unit load. product. All storage locations and unit loads are stan- To evaluate the quality of these closed-form formu- dard and have the same size. Shuttles simultaneously las (obtained using a continuous-space approxima- move unit loads into the empty locations (at least tion), a Monte Carlo simulation is performed based one empty location needed on each level) to maneu- on a real (discrete) application. Our second contri- ver the requested unit load to the lift at the unit bution is deriving closed-form expressions for the y x Unit load Shuttle z Lift z y x Figure 1 (Color online) A Live-Cube Compact Storage System with a Lift Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 3 optimal dimensions of a live-cube compact storage and conveyors take care of the movement of unit system minimizing the expected retrieval time. To loads in depth direction. They optimally solve the obtain these dimensions, we propose a mixed-integer problem by decomposing it into three cases. nonlinear model. We split the model into several solv- Compared with deep-lane compact storage sys- able submodels such that the optimal solution for tems (conveyor-based and satellite-based), live-cube the model can be obtained. Our third contribution systems ensure individual access to each unit load, is to compare live-cube systems with traditional sys- reducing the response times to customers. Gue and tems in terms of investment and operational costs, Kim (2007) are the first to study a single-level live- energy consumption, and carbon dioxide (CO2 ) emis- cube system (which they call a “puzzle-based” stor- sions based on a real application. In addition, we age system) at an operational level. For systems with investigate the impact of a live-cube system configu- a single empty location, they develop a method to ration on the investment and operational costs of such maneuver unit loads to the I/O point minimizing the systems. retrieval time of an arbitrary unit load. In addition, for Traditional automated material handling systems systems with multiple empty locations, they propose have been studied extensively in the literature over a heuristic retrieval method yielding short retrieval the past 30 years. For a general review of the design times. To our knowledge, our paper is the first study and control of automated material handling sys- that investigates the impact of system configuration tems, we refer to Johnson and Brandeau (1996) and on response time, investment, and operational costs Roodbergen and Vis (2009). Only a few papers in of multilevel live-cube systems. Our study differs the literature study compact storage systems (e.g., from Gue and Kim’s (2007) work as we focus on the Stadtler 1996; Sari, Saygin, and Ghouali 2005; Gue and optimal design of multilevel live-cube systems which Kim 2007; de Koster, Le-Duc, and Yu 2008; Yu and minimizes the system response time, using a simple de Koster 2009; Zaerpour, Yu, and de Koster 2015). method to maneuver unit loads to the I/O point. Stadtler (1996) and Zaerpour, Yu, and de Koster (2015) The remainder of this paper is organized as follows. study compact storage systems at an operational deci- Section 2 describes the problem and gives the nota- sion level aiming at shorter response times. Stadtler tions and model assumptions. Section 3 derives the (1996) introduces a five-module storage and retrieval expected retrieval time for multilevel live-cube sys- assignment model for a deep-lane compact storage tems. Section 4 proposes the mathematical model and system where a lift takes care of vertical movement obtains the optimal system dimensions for multilevel and a satellite machine moves the unit load inside live-cube systems. Section 5 discusses further exten- the storage lane for deep storage or retrieval. To sions. Section 6 discusses a real case study. Section 7 access an individual unit load stored in depth, block- concludes the paper. ing unit loads stored in front of the unit load need to be removed first. Zaerpour, Yu, and de Koster (2015) also study a deep-lane compact storage sys- 2. Problem Description and tem where a crane is responsible for movements of Assumptions unit loads in horizontal and vertical directions and Consider a system with length l (in travel time units), a satellite linked to the crane moves the unit loads width or depth w (in travel time units), and height in deep storage lanes. They propose a new storage h (in travel time units). For the sake of convenience assignment policy in order to minimize total retrieval and without loss of generality, we suppose that the time of customer requests. Sari, Saygin, and Ghouali travel time in length direction of the system is not (2005), de Koster, Le-Duc, and Yu (2008), and Yu and less than the travel time in the depth direction; i.e., de Koster (2009) aim at shorter response times by l ≥ w (see also Bozer and White 1984 and Eynan and improving the physical design of compact storage sys- Rosenblatt 1994). tems. Sari, Saygin, and Ghouali (2005) study a deep- The retrieval time of a unit load at a given level lane compact storage system where a lift takes care depends on its (x1 y) position and the distribution of of vertical movement of unit loads and conveyors are empty locations on its trajectory to the lift. In case only in charge of movement of unit loads in deep stor- one empty location is available, Gue and Kim (2007) age lanes. To access a unit load, other unit loads also show it equals 4x + 2y − 8 (if x > y). However, in prac- need to be moved by the rotating conveyor. They tice utilizations of not more than 95% are common in derive closed-form expressions for travel time of unit systems containing 1,000 unit loads or more (de Koster loads for any system design and validate their results 1996). This implies that at one level sufficient empty using simulation. de Koster, Le-Duc, and Yu (2008) locations are commonly available to create a virtual and Yu and de Koster (2009) study the optimal design aisle to maneuver a requested unit load to the lift at the of a deep-lane compact storage system where a crane same level without interfering with other unit loads moves unit loads in horizontal and vertical directions (a sufficient condition for this is that the number of Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 4 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS empty locations on a level equals at least the maxi- maxx + y z where x + y is from Equation (1) and z is mum of the number of rows and number of columns; the time needed for the lift to go from the I/O point for more details see Online Supplement 1 (available as to level z. supplemental material at http://dx.doi.org/10.1287/ 2. The time needed for the lift to return to the I/O trsc.2015.0586)). This assumption does not prohibit point from level z. a potentially high utilization of the system while it Thus, the retrieval time of a unit load located at reduces the retrieval time of any unit load significantly (x y z to the I/O point for a multilevel system can compared to the results of Gue and Kim (2007). Fig- be calculated by ure 2 shows an example of a situation where a car needs to be moved out of a single-level live-cube park- tx y z = maxx + y z + z (2) ing system. Block movement is possible, i.e., moving multiple unit loads simultaneously, if the movement Equation (2) gives an estimate for the retrieval time path is a straight line. As a result, a virtual aisle can of a given unit load. To evaluate the performance of be created in a negligible constant time by simulta- the system, the expected retrieval time of an arbitrary neously moving interfering cars. The desired car can unit load must be derived. This expression can then then be moved to the I/O point without any further be used to optimize the dimensions of the system, interference. leading to minimum expected retrieval time. Proposition 1. With an increasing size of a live-cube system, the constant time required to create a virtual aisle for a requested unit load can be neglected. 3. Derivation of the Expected Proof. See Online Supplement 1. Retrieval Time To calculate the expected retrieval time of an arbi- Therefore, for a system where the I/O point is at the trary unit load, the probability density function of the lower left corner, the retrieval time of a unit load at retrieval time must be obtained. This is done by split- location (x y) of a single-level live-cube storage sys- ting the system into several cases and then obtain- tem, tx y, can be calculated by ing the probability density function for each case individually. tx y = x + y (1) The retrieval time (t) of a unit load in a multilevel sys- In Equation (1) we assume the coordinates of the tem at position (x y z) is calculated by maxx + y z system to be continuous. In a multilevel live-cube sys- + z. The set x y z maxx + y z + z ≤ t includes tem, (x y z) is a location (in travel time units) where all of the locations with coordinates (x y z) where z refers to the level coordinate (in travel time units) maxx + y z + z ≤ t. Let T = maxX + Y Z + Z, in a vertical direction. When idle, the lift waits at the where X, Y , and Z are random variables represent- ground level (this assumption will be relaxed in §5.2). ing the length, depth, and height coordinates of the The retrieval time of a unit load consists of the fol- load to be retrieved, respectively. The cumulative dis- lowing two components: tribution function of T is denoted by F t and the 1. The time needed to bring the unit load to the lift. probability density function of T is denoted by f t. Since the movements in x and y directions are inde- The expected retrieval time of an arbitrary unit load pendent of movement in z direction, this time equals in a multilevel live-cube storage system with a given Retrieval request A solution by creating a virtual aisle y I/O I/O x Figure 2 (Color online) A Virtual Aisle in a Single-Level Live-Cube Parking System Note. The requested car is shown by a different color. Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 5 capacity can be calculated as referring to a specific configuration of the system. maxw+l h+h There are four cases of system configuration: ET = tf t dt • Case A. h ≤ w. t=0 • Case B. w ≤ h ≤ l. where t represents the retrieval time for retrieval loca- • Case C. l ≤ h ≤ l + w. tion x y z and f t represents the probability den- • Case D. l + w ≤ h. sity function of retrieval time, t0 ≤ t ≤ maxw +l h+ These cases can again be split into several subcases. h. The cumulative distribution function of retrieval As an illustration, in §3.1 the calculation of ET for time, F t, can be calculated as one subcase of case A (subcase A.1) is explained. Sec- tion 3.2 gives the formulas of ET for all other cases. F t = P T ≤ t = P maxX + Y Z + Z ≤ t 3.1. Computation of Expected Retrieval Time of = P X + Y + Z ≤ t ∧ 2Z ≤ t an Arbitrary Unit Load for Case A The two conditions, X +Y +Z ≤ t and 2Z ≤ t are not In this section, we give the detailed procedure of cal- independent and we have to consider both conditions culating ET for one of the subcases of case A (sub- simultaneously, which increases the complexity of cal- case A.1: l ≤ 2h). Given the configuration of subcase culations. The two following inequalities form the A.1, Figures 5(a)–5(g) show how the region T ≤ t shape that includes all of the locations with retrieval changes with increasing retrieval time t. time less than or equal to t: x + y + z ≤ t and 2z ≤ t where x y z ≥ 0. Figure 3(a) illustrates the shape that Definition 1 (Critical Retrieval Time). The crit- includes all of the locations with retrieval time less ical retrieval time is a time instant at which the region than or equal to t. For any value of retrieval time t, T ≤ t changes shape, resulting in a different volume the probability that the random variable T is less than formula. or equal to t is given by Equation (3) It can be shown that for case A, the region T ≤ t changes at the following values for t: F t = P T ≤ t volume of the region T ≤ t in the system w l 2h w + h l + h w + l l + w + h (4) = (3) volume of the system For subcase A.1, we can arrange the critical retrieval In Figure 3(a), two planes, x + y + z = t and 2z = t, are times as follows: contour planes where the retrieval time of locations located on these two planes equals t. w ≤ l ≤ 2h ≤ h + w ≤ h + l ≤ l + w ≤ l + w + h (5) However, the region in Figure 3(a) is also restricted by the system dimensions (Figure 3(b)). Therefore, Equation (5) gives seven different possible positions the shape in Figure 3(a) can be truncated in differ- of the retrieval time t 0 ≤ t < w; w ≤ t < l l + w ≤ ent ways depending on the system dimensions and t ≤ l + w + h, each corresponding to one of the figures retrieval time t. Figure 4 illustrates all possible shapes in Figure 5. for different system dimensions and retrieval times t For Figure 5(a), we have 0 ≤ t < w. Because t < w, that can be obtained through a tedious decomposition the shape does not touch the y = w plane of the cube. process. In addition, since w ≤ l, the shape also does not touch It turns out that the calculation of ET can be clas- the x = l plane of the cube. Furthermore, z = t/2 repre- sified into four different complementary cases, each sents the top surface of the shape. However, t < w and also w < 2h (see Equation (5)) and therefore t < 2h. (a) (0, t/2, t/2) (b) (0, w, h) Hence, the plane z = t/2 of the shape does not touch (0, 0, t/2) (0, 0, h) the plane z = h of the cube. (l, w, h) z z For Figure 5(b), w ≤ t < l. Because w ≤ t the shape (t/2, 0, t/2) touches the y = w plane of the cube. The point t − (l, 0, h) (0, t, 0) w w 0 in the y = w plane has a retrieval time equal (0, w, 0) to tt − w + w, which is greater than w. The explana- I/O tion of the shape of the polytope in x and z directions I/O y y (l, w, 0) is similar to the previous shape. In Figure 5(c), l ≤ t < 2h. Because l ≤ t, the shape x (t, 0, 0) touches the x = l plane of the cube. The retrieval x (l, 0, 0) time of the corner point l t − l 0 of the shape has Figure 3 (Color online) (a) The Region Formed by x + y + z ≤ t and retrieval time tt − l + l, which is greater than l. The 2z ≤ t; (b) Cubic Shape of a Multilevel System with its explanation of the shape in y and z directions is sim- Corner Points ilar to the previous shape. Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 6 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS (b) t/2 ≥ h and w ≤ t ≤ w + h (a) t/2 ≥ h and t ≤ w (c) t/2 ≥ h and l ≤ t ≤ w + h (d) t/2 ≥ h and w + h ≤ t ≤ l and t ≤ l (e) t/2 ≥ h and w + h ≤ t (f) t/2 ≥ h and (g) t/2 ≥ h and (h) t/2 ≥ h and and l ≤ t ≤ min{l + w, l + h} l+w≤t≤l+h l+h≤t≤l+w max{l + h, l + w} ≤ t ≤ l + w + h (j) t/2 ≤ h and w ≤ t ≤ 2w (i) t/2 ≤ h and t ≤ w (k) t/2 ≤ h and l ≤ t ≤ 2w (l) t/2 ≤ h and 2w ≤ t ≤ l and t ≤ l (m) t/2 ≤ h and 2w ≤ t (n) t/2 ≤ h and (o) t/2 ≤ h and (p) t/2 ≤ h and and l ≤ t ≤ l + w l + w ≤ t ≤ 2l 2l ≤ t ≤ 2(l + w) 2(l + w) ≤ t ≤ 2h Figure 4 (Color online) All Possible Shapes of the Region T ≤ t in the System The other shapes of Figure 5 can be explained Figure 5(a) is given by similarly. t/2 t−z We now derive the volumes of the first three shapes t − y − z dy dz = 7t 3 /48 z=0 y=0 in Figure 5. To derive the volumes of the first three shapes, we consider the base of the shape to be in the In Figure 5(b) y and z are bounded as follows: 0 ≤ y ≤ yz plane. minwt −z, 0 ≤ z ≤ t/2. The volume is given by For Figure 5(a), y and z are bounded: 0 ≤ y ≤ t − z t/2 minw t−z and 0 ≤ z ≤ t/2. The volume of the shape can be t − y − z dy dz z=0 y=0 obtained by integrating the function minl t − y − z over the mentioned boundaries. In this case minl t − Figure 5(b) shows the line y = t −z intersects y − z always equals t − y − z, and so the volume of y = w plane at point 0wt −w. For any z ≤ t −w, Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 7 (a) 0 ≤ t < w (b) w ≤ t < l (c) l ≤ t < 2h (0, w, t – w) (t – w, w, 0) (l, 0, t – l) (l, t – l, 0) (d) 2h ≤ t < h + w (e) h + w ≤ t < h + l (f) h + l ≤ t < l + w (0, t – h, h) (t – w – h, w, h) (t – h, 0, h) (l, t – l – h, h) (g) l + w ≤ t ≤ l + w + h (l, w, t – l – w) Figure 5 (Color online) Region of T ≤ t for Subcase A.1 l ≤ 2h with Increasing t minwt −z = w and for any z > t −w, minwt −z = the removed pyramid equals t − l3 /6. The resulting t −z. Thus, the volume of Figure 5(b), where w ≤ t < l, volume of Figure 5(c) is given by is given by t−w w t−w w t/2 t−z t − y − z dy dz z=0 y=0 t − y − z dy dz + t − y − z dy dz z=0 y=0 z=t−w y=0 t/2 t−z + t − y − z dy dz − 16 t − l3 = 1 48 −t 3 + 24tt − ww + 8w 3 z=t−w y=0 = 1 48 −t 3 − 8t − l3 + 24tt − ww + 8w 3 In Figure 5(c), y and z are bounded as follows: 0 ≤ y ≤ minw t −z, 0 ≤ z ≤ t/2. Figure 5(c) is a truncated The volumes of the other figures can be derived version of Figure 5(b); the plane x = l has removed similarly. The formulas representing the volumes of a pyramid shaped top in x direction. As a result, Figures 4(a)–4(p) are given in Online Supplement 2. the volume (third shape) = volume (second shape) − The cumulative distribution function for each volume (removed pyramid). boundary in Figure 5 can now be derived by dividing The height of the removed pyramid equals t − l, its volume by the volume of the cube: Vcube = wlh. By and the base of the pyramid is a right triangle with taking the derivative, the probability density function both legs equaling t − l. Consequently, the volume of is obtained. Corresponding to Figures 5(a)–5(c), the Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 8 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS w+h cumulative and density functions are given by Z 1 1 + t − 4−l +t52 + 4−3h2 +6ht −3t 2 5 7t 3 /48 2h 2 6 Ft<w 4t5 = 1 tw 1 1 wlh + + 4t −w5w dt dF 21t 2 /48 2 2 wlh ft<w 4t5 = t<w = 1 Z l+h 1 dt wlh 1 + t − 4−l +t52 +hw dt 1/484−t 3 +24t4t −w5w +8w 3 5 w+h 2 wlh Fw≤t<l 4t5 = 1 Z l+w 1 wlh + t − 4l −t +w52 dFw≤t<l −3t 2 +24tw +244t −w5w l+h 2 fw≤t<l 4t5 = = 1 dt 48hlw 1 1 + 4h+l −t +w52 dt 1/484−t 3 −84t −l53 +24t4t −w5w +8w 3 5 2 wlh Fl≤t<2h 4t5 = 1 wlh Z l+w+h 4h+l −t +w52 + t dt dFl≤t<2h l+w 2wlh fl≤t<2h 4t5 = dt h3 +12hlw +12lw4l +w5 −244l −t52 −3t 2 +24tw +244t −w5w = 0 = 0 24lw 48hlw The probability density function for other bound- This expression is only valid for subcase A.1 (l ≤ 2h). aries in Figure 5 can be calculated similarly. This However, case A can be divided into five different sub- results in cases because of five different possible arrangements  of the critical retrieval times (see Equation (4)). 421t 2 /485/4wlh5 if 0 ≤ t < w1 • Subcase A.1. l ≤ 2h.     • Subcase A.2. w ≤ 2h and 2h < l ≤ w + h.   −3t 2 + 24tw + 244t − w5w  if w ≤ t < l1      48hlw • Subcase A.3. w ≤ 2h and l > w + h.    • Subcase A.4. w > 2h and 2h < l ≤ w + h. −244l − t52 − 3t 2 + 24tw + 244t − w5w      • Subcase A.5. w > 2h and l > w + h. 48hlw    The expected retrieval time 4E6T 75 can be calculated if l ≤ t < 2h1   in a similar way for subcases A.2–A.5. The results are     summarized in Proposition 2.  −3h2 − 34l − t52 + 6ht − 3t 2 + 3tw + 34t − w5w      6hlw  f 4t5 = Proposition 2. The expected retrieval times (E6T 7) for if 2h ≤ t < w + h1    all subcases of case A (h ≤ w) can be formulated as Equa- tion (7)  −41/254l − t52 + hw      if w + h ≤ t < l + h1 hlw   E6TA 7 = E6TA1 7 = E6TA2 7 = E6TA3 7 = E6TA4 7 = E6TA5 7    h + 24l − t + w5   h3 + 12hlw + 12lw4l + w5    if l + h ≤ t < l + w1 2lw = if h ≤ w0 (7)      24lw 4h + l − t + w52      if l + w ≤ t ≤ l + w + h1 Proof. See Online Supplement 3. 2hlw       = 0 otherwise.  3.2. Expected Retrieval Time of an Arbitrary Unit (6) Load for Cases B, C, and D The expected retrieval time E6T ] for cases B, C, and D From this the total expected retrieval time of an can be derived in a fashion similar to case A. Proposi- arbitrary unit load for subcase A.1 can be calcu- tions 3–5 summarize the results. All of the proofs can lated as be found in Online Supplement 4. Z w 21t 3 1 E6TA1 7 = · dt Proposition 3. The expected retrieval time (E6T 7) of 0 48 wlh Z l t 2 tw 1 an arbitrary unit load for case B (w ≤ h ≤ l5 is given by 1 + t − + + 4t −w5w dt Equation (8) w 16 2 2 wlh Z 2h t 2 1 4h3 + 6h2 42l − w5 − w 3 + 4h43l2 + 3lw + w 2 5 + t − − 4−l +t52 E6TB 7 = l 16 2 24hl tw 1 1 if w ≤ h ≤ l0 (8) + + 4t −w5w dt 2 2 wlh Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 9 Proposition 4. The expected retrieval time (E6T 7) of the optimal solution, the optimal values of the length an arbitrary unit load for case C (l ≤ h ≤ l + w) is given and width of the system are equal (w ∗ = l∗ 5. by Equation (9) Proposition 6. In the optimal solution of case A, the E6TC 7 optimal values of the length and width of the system are h4 +l4 +6h2 4l −w52 +w 4 −4h3 4l +w5−4h4l +w53 equal (w ∗ = l∗ ). =− 24hlw Proof. See Online Supplement 5. if l ≤ h ≤ l +w0 (9) Because w ∗ = l∗ , w can be replaced with l in the Proposition 5. The expected retrieval time (E6T 7) of model. Furthermore, since wlh = V , we can replace l an arbitrary unit load for case D (l + w ≤ h5 is given by with 4V /h51/2 . Thus, the objective function can be Equation (10) written as a function of h 12h2 + 2l2 + 3lw + 2w 2 h4h3 + 12V + 244V /h53/2 5 E6TD 7 = if l + w ≤ h0 (10) E6TA 4h57 = 0 12h 24V 4. General Model and Optimization The optimal value of h, minimizing the objective The results obtained in §3 can be used to optimize function can be obtained by the first order condition the system dimension by minimizing the expected √ retrieval time. To optimize the dimensions of a live- dE6TA 7 h4 + 3hV − 3V V = = 00 cube system, we propose the following mixed-integer dh 6hV general model (Model MG): The following equations give the optimal closed- Model MG form expressions for l, w, and h that satisfy the con- min X ui E6Ti 71 (11) dition 0 < h ≤ w: i∈8A1B1C1D9 √ 1/3 h∗ 4V 5 = −2 + 12 411 − 3 135 subject to lwh = V 1 (12) √ 1/3 1/3 1/3 + 12 411 + 3 135 V l −w ≥ 01 (13) X = 00874V 1/3 1 (21) ui = 11 (14) √ 1/3 −2 + 12 411 − 3 135 i∈8A1B1C1D9 w ∗ 4V 5 = l∗ 4V 5 = √ 1/3 1/3 −1/2 1/3 uA 4w −h5 ≥ 01 (15) + 12 411 + 3 135 V uB 4h−w5 ≥ 01 (16) = 10069V 1/3 0 (22) uB 4l −h5 ≥ 01 (17) The second derivative test proves that the objective uC 4h−l5 ≥ 01 (18) function is a convex function and therefore the opti- mal dimensions give the minimum objective value uC 4l +w −h5 ≥ 01 (19) p uD 4h−l −w5 ≥ 01 (20) d 2 E6TA 7 h2 3 V /h = + ≥ 00 Decision variables: l > 01w > 01h > 01 and ui ∈ 80119 dh2 2V 4h2 for i ∈ 8A1B1C1D90 4.2. Optimizing the Dimensions in Cases B, C, and D Equation (11) minimizes the expected retrieval The following equations give the optimal closed-form time E[T ]. Constraint (12) makes sure that the given expressions for l, w, and h for case B (w ≤ h ≤ l5 as capacity is achieved. Constraint (13) ensures the length a function of the volume of the system. The detailed is at least equal to the width of the system. Con- calculation is given in Online Supplement 6. straint (14) guarantees exactly one of the cases is con- sidered in the objective function. Constraints (15)–(20) √ h∗ 4V 5 = w ∗ 4V 5 = 4−3 + 1551/3 V 1/3 = 00955V 1/3 1 take care of the feasibility of the solutions of each case. To solve Model MG, we first solve each case indi- V 1/3 l∗ 4V 5 = √ = 10094V 1/3 0 vidually and obtain the optimal solution in each 4−3 + 1552/3 case. Next, by comparing the optimal solutions obtained from the four cases, the optimal solution of For case C (l ≤ h ≤ l + w5, the optimal closed-form Model MG can be found. expressions for l, w, and h as a function of the volume of the system are 4.1. Optimizing the Dimensions in Case A (h ≤ w) For this case the following proposition shows that in h∗ 4V 5 = w ∗ 4V 5 = l∗ 4V 5 = V 1/3 0 Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 10 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS For case D (l + w ≤ h the optimal closed-form expressions for l, w, and h as a function of the volume 15 of the system are V 1/3 w ∗ V = l∗ V = = 0793V 1/3 E[T *] 10 Case D 21/3 Case C h∗ V = 22/3 V 1/3 = 1587V 1/3 Case B 5 Case A 4.3. Optimal Solution of Model MG By substituting the results of §§4.1 and 4.2 in the objective function, the optimal values of ETA , ETB , 200 400 600 800 1,000 ETC , and ETD as a function of volume of the sys- V tem, V , can be obtained √ Figure 6 (Color online) Optimal ET for Cases A, B, C, and D versus ETA∗ = 38 12 68+22/3 9947 −429 131/3 System Volume √ 1/3 +22/3 9947 +429 131/3 V 1/3 = 1531V 1/3 (23) 5. Further Discussion and Extension In this section, we investigate how the expected 2/3 √ 1/3 ETB∗ = 3 4 132 +5 15 V 1/3 = 1538V 1/3 (24) retrieval time changes by departing from optimal dimensions (see §5.1). Section 5.2 discusses extensions 37 1/3 of our model. ETC∗ = V = 1542V 1/3 (25) 24 1/3 ETD∗ = 24 55 1 V = 1819V 1/3 (26) 5.1. Further Discussion 2 Section 3 derives four closed-form expressions for the Figure 6 illustrates the minimum ET of four cases optimal expected retrieval time for a live-cube system for varying system volumes. by dividing the solution space into four regions. Fig- As can be seen from Equations (23)–(26), the opti- ure 7(a) illustrates which case has to be considered mal solution of case A gives the minimum E[T ] for for any volume, height, and shape factor (w/l). As the Model MG. Therefore, for any given volume of the shape factor decreases, case B becomes the dominant system, the following dimensions give the minimum feasible region. Moreover, by increasing the height of expected retrieval time (see Equations (21) and (22)). the system, the feasible region moves from case A to √ h∗ V = −2+ 12 11−3 131/3 case B to case C and finally to case D. For any given √ volume v, a projection of Figure 7(a) on the plane 1/3 + 12 11+3 131/3 V 1/3 V = v gives all possible values of height and shape √ factor for each case. For example, Figure 7(b) shows w ∗ V = l∗ V = −2+ 12 11−3 131/3 the projection of the three-dimensional plot on the √ 1/3 −1/2 1/3 plane V = 106 . For V = 106 , A∗ (h = 8744, w/l = 1), + 12 11+3 131/3 V B∗ (95.57, 0.8728), C∗ (100 1), and D∗ (15874 1) are and u∗A = 1 u∗B = u∗C = u∗D = 0 the optimal solutions of cases A, B, C, and D with (a) (b) A* C* D* 1.0 Case D 1.0 B* Case C 0.8 0.6 Case C Case D w/l 0.5 w/l Case A Case A 0.4 Case B 0 200 0.2 Case B 150 500,000 100 h 0 50 0 50 100 150 200 V 1 × 106 0 h Figure 7 (Color online) (a) Feasible Regions: Cases A, B, C, and D for Varying Volume, Height, and Shape Factor; (b) Projection of 3D Plot on the Plane V = 106 Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 11 the objective values 153.097, 153.789, 154.167, and Figure 8 shows the gap for different values of h and 181.889, respectively. The arrows in Figure 7(b) illus- V . Figure 8(a) illustrates Equation (27) for a relatively trate that starting from a random feasible point, the large range of V (V ∈ 401 106 75. Figure 8(b) presents the global optimum (A∗ ) can always be obtained by con- contour plot of Equation (27) for different values of h tinuously modifying the parameters into the direction and V . Figure 8 illustrates that for any given volume, of the arrows. the gap increases with an increasing height of the sys- Table 1 summarizes the obtained results from tem. However, for smaller values of the volume, the §§3 and 4, including the closed-form formulas for the increase is more dramatic. expected retrieval time and optimal dimensions for We also investigate how variation of shape factor each case individually. (w/l5 affects the expected retrieval time of a live-cube As can be seen from Table 1, the optimal dimen- system. Since it was assumed that w ≤ l, the shape fac- sions in case A lead to the minimum expected tor can obtain any value in 401 17. Equation (28) deter- retrieval time. However, the optimal dimensions in mines the gap between the expected retrieval time of cases B and C also result in a very close to optimal varying w/l and optimal expected retrieval time given (less than 1% gap) expected retrieval time. This means the volume is fixed (V = 106 ) if the company is willing to have a taller system with a smaller footprint, the optimal dimensions of cases B E6T 4w1 l1 h57 − E6T ∗ 7 Gapw/l = or C can be used, which give a gap of less than 1% E6T ∗ 7 from the optimal expected retrieval time. A cubic in  time system (l = w = h5 appears to be not optimal but  4E6TA 4w1 l1 h57 − E6T ∗ 75/E6T ∗ 7 if h ≤ 44w/l5V 51/3 1   has a gap of only 1% from the optimal solution.    4E6TB 4w1 l1 h57 − E6T ∗ 75/E6T ∗ 7  To check how the configuration of a live-cube stor-    1/3 1/3 (28)  age system affects the expected retrieval time, we  if 44w/l5V 5 < h ≤ 44l/w5V 5 1   have varied h for a fixed volume of the system. The = 4E6TC 4w1 l1 h57 − E6T 75/E6T 7 ∗ ∗ gap is determined by Equation (27). The right-hand  1/3  if 44l/w5V 5 < h   side of Equation (27) should be evaluated at w = l. In     ≤ 444l/w5 + 4w/l5 + 25V 51/3 1 Equation (27), boundaries are calculated according to  4E6TD 4w1 l1 h57 − E6T ∗ 75/E6T ∗ 7    the conditions of the cases. If the condition of more   if h > 444l/w5 + 4w/l5 + 25V 51/3 0  than one case is satisfied, the case with a minimum objective value is considered. Figure 9 illustrates the three-dimensional and con- E6T 4h57 − E6T 7 ∗ tour plots of Equation (28) for varying values of shape Gaph = factor (w/l) and h. Figure 9(a) demonstrates that for a E6T ∗ 7  decreasing shape factor (w/l) for any given h, the gap ∗ ∗ 4E6TA 4h57 − E6T 75/E6T 7  increases. The increase is more dramatic for smaller 1/3 if h ≤ V 1  values of h (see Figure 9(b)).     4E6TC 4h57 − E6T ∗ 75/E6T ∗ 7  = (27)  if V 1/3 < h ≤ 22/3 V 1/3 1  5.2. Extensions  4E6TD 4h57 − E6T ∗ 75/E6T ∗ 7 In this section, we propose some extensions of our prob-     if h > 22/3 V 1/3 0 lem and discuss the possible solution approaches.  Table 1 Optimal Retrieval Times and System Dimensions for Each Case Case E6T 7 formula Optimal dimensions Aa h 3 + 12hlw + 12lw 4l + w 5 h ∗ 4V 5 = 0087446V 1/3 24lw w ∗ 4V 5 = l ∗ 4V 5 = 1006937V 1/3 E6TA∗ 7 = 1053097V 1/3 B 4h 3 + 6h 2 42l − w 5 − w 3 + 4h43l 2 + 3lw + w 2 5 h ∗ 4V 5 = w ∗ 4V 5 = 0095573V 1/3 24hl l ∗ 4V 5 = 1009479V 1/3 E6TB∗ 7 = 1053789V 1/3 C h 4 + l 4 + 6h2 4l − w 52 + w 4 − 4h 3 4l + w 5 − 4h4l + w 53 h ∗ 4V 5 = w ∗ 4V 5 = l ∗ 4V 5 = V 1/3 − 24hlw E6TC∗ 7 = 1054167V 1/3 D 12h 2 + 2l 2 + 3lw + 2w 2 w ∗ 4V 5 = l ∗ 4V 5 = 00793701V 1/3 12h h ∗ 4V 5 = 1058740V 1/3 E6TD∗ 7 = 1081889V 1/3 a The overall optimal solution is obtained in case A according to §4.3. Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 12 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS (a) (b) 1,000,000 2,000 4,000 6,000 15,000 10,000 1,000,000 500,000 8,000 %Gap V 5,000 0 0 500,000 10,000 V 12,000 5,000 14,000 h 10,000 0 0 0 5,000 10,000 h Figure 8 (Color online) (a) 3D Plot; (b) Contour Plot of the Gap (%) Between ET h and ET ∗ for Varying h and V ∈ 0 106 Minimizing the maximum response time. The objective By comparing the two above objective values, it is of Model MG is to minimize the expected retrieval obvious that the optimal dimensions that minimize time. However, it is also possible to minimize the max- the maximum response time of a live-cube system are imum response time of a live-cube compact storage l∗ V = w ∗ V = h∗ V = V 1/3 , that is, a cubic in time system. The maximum response time of the live-cube system. Moreover, the optimal solution of Model MG system is the retrieval time of the unit load located at (l∗ V = w ∗ V = 1069V 1/3 , h∗ V = 0874V 1/3 ) has a the location farthest from the I/O point. Therefore, the gap of only 1% from this optimal solution. In conclu- maximum response time is maxl + w h + h. sion, the optimal solution of Model MG has a gap of The corresponding mathematical model can be for- only 1% from the optimal solution of Model MRT and mulated as follows (Model MRT): vice versa. Location of the lift. In our problem, the lift is as- min maxl + w h + h sumed to be located at the left corner of the system subject to lwh = V and when idle it waits at the ground floor. However, decision variables l > 0 w > 0 and h > 0 by changing the location of the lift, we can calculate the retrieval time differently. For example, if the lift is The model can be solved optimally by decompos- located at the middle of the front face of the system, ing it into two submodels and using Lagrange multi- then the travel time in x direction of any arbitrary unit pliers for each submodel. load to the lift at the same level can be reduced. To • If l + w > h, then the objective function becomes solve such a problem, the travel time formula and the l + w + h, and therefore the optimal dimension sizes expected retrieval time should be derived and then and objective value for this case are l∗ V = w ∗ V = the same model as Model MG can be constructed. h∗ V = V 1/3 and 3V 1/3 , respectively. However, it is also possible to consider the system • If l + w ≤ h, then the objective function becomes as two similar subsystems where the lift is located at 2h, and therefore the objective function is minimized the corner of each subsystem and then results of our when h = l + w. Thus the optimal dimension sizes and model can still be applied. The volume of each sub- objective value are l∗ V = w ∗ V = V /21/3 , h∗ V = system will be the half volume of the system. Thus, 2V /21/3 , and 4V /21/3 , respectively. the optimal dimensions and expected retrieval time of (a) (b) 1.0 80 60 40 20 100 150 1.0 0.5 %Gap 100 w/l 50 120 0 140 0.5 w/l 160 180 50 h 0 100 0 50 100 h Figure 9 (Color online) (a) 3D Plot; (b) Contour Plot of the Gap (%) Between ET w l h and ET ∗ for Varying w /l, h, and V = 106 Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 13 Table 2 Comparison of the Expected Retrieval Times Obtained from the Closed-Form Formula and Simulation for the Base Example a Corresponding case E6TD 7 from derivation (sec) E6TS 7 from simulation (sec) %Gap Footprint (m2 ) Base example A 25.66 26002 ± 0039 1.44 624 a Gap = 4—E6TS 7 − E6TD 7—/E6TD 75 × 100%. a live-cube system with a lift at the middle position First, we evaluate the quality of closed-form for- are h∗ 4V 5 = 0069406V 1/3 , w ∗ 4V 5 = 0084876V 1/3 , l∗ 4V 5 = mulas using a Monte Carlo simulation based on this 1069753V 1/3 , and E6T ∗ 7 = 1021513V 1/3 . Results show real (discrete) system. Next, we investigate the impact that a lift at the middle position leads to a 20% reduc- of the number of empty locations on system perfor- tion in expected retrieval time. mance. Also, we compare the investment and opera- In general, our studied problem can provide the tional costs of this live-cube system with a traditional basic analysis for any other live-cube compact stor- system of the same capacity. Furthermore, we investi- age system, including live-cube systems with multiple gate the impact of live-cube system dimensions on the lifts. For a system with two lifts, optimal locations of investment cost of such systems. Finally, we compare the lifts can improve the performance of a live-cube energy consumption and CO2 emissions of live-cube system. Such a problem appears to be complicated systems with traditional systems. because two unit loads can be retrieved simultane- ously. Thus, moving a shuttle to create an empty loca- 6.1. Evaluation of Closed-Form Expressions tion for one retrieving unit load might interfere with To obtain the closed-form formulas, a continuous sys- retrieving the other unit load. However, it is still pos- tem is assumed, but in reality the systems are discrete. sible to use our results because a system with two lifts Therefore, we test the accuracy of the closed-form for- can be decomposed into subsystems each with one lift mulas for estimating the response time of a discrete at the corner. It turns out that the optimal locations system. As it can be expected that the continuous- of two lifts depend on the system configuration. Lifts located at two opposite sides of the system minimize space approximation is better for large-scale systems, the expected retrieval time if l ≤ 2w. In case l > 2w, we also include small-size systems. both lifts should be located at one side of the system, First, we present the results for the live-cube sys- each with equal distance to the middle and corner of tem illustrated in the appendix (base example). Then, the system (for the proof, see Online Supplement 7). we perform a sensitivity analysis to test the quality of closed-form formulas for different sources of varia- tion in the input parameters. All instances have been 6. A Case Study simulated using MATLAB software on a portable com- In this section, we consider a real live-cube sys- puter Intel Pentium M 1.86 GHz with 512 MB RAM tem from an international company that produces and Windows XP professional. The results for the base live-cube parking systems. Although the system can example are shown in Table 2. The expected retrieval be designed in any desired size by the customer, we focus on a medium-sized live-cube parking sys- time for each instance is the average of 1,000 randomly tem (a total capacity of 192 cars) in an urban set- generated retrieval requests. The confidence intervals ting (urban setting is defined as anywhere the land are given at a 95% confidence level. price is above $800 per square meter; Sollohub 2010). In the subsequent tables, we vary system capacity Input parameters describing this live-cube system are (Table 3) and height tiers (Table 4) whereas the other shown in the appendix. Compared with conventional parameters stay the same. In Table 3, when varying multistory car parks, live-cube parking systems take the capacity, the number of tiers in x and y direc- 40%–50% less footprint area for the same number of tions are equal; i.e., l/w = 1. In Table 4, when varying parking slots. height tiers, the system capacity is constant. Table 3 Results of Sensitivity Analysis When Varying System Capacity System capacity Corresponding E6TD 7 from E6TD 7 from (no. of locations) case derivation (sec) simulation (sec) %Gap 100 A 19030 19063 ± 0036 1072 500 A 40069 41014 ± 0064 1011 1,000 A 45072 45039 ± 0066 0071 2,000 B 56088 57019 ± 0075 0057 4,000 D 86073 87010 ± 1040 0043 Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 14 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS Table 4 Results of Sensitivity Analysis When Varying Height Tiers Height Corresponding E6TD 7 from E6TS 7 from (no. of tiers) case derivation (sec) simulation (sec) %Gap 2 A 42009 41060 1015 4 A 32004 32000 0013 6 A 27081 27097 0058 12 C 27011 27047 1031 24 D 46039 46069 0063 Table 5 Impact of Number of Empty Locations on Average Retrieval Time in a Live-Cube Parking System No. of empty Capacity Storage E[T ] Trade-off locations Land (m2 ) Volume (m3 ) (No. of cars) densitya (%) (sec)b Base example max8—I—1 —J—9 624.0 14,601.6 192 (24 cars per level) 80 33.48 Space vs. retrieval time 50% 998.4 23,362.5 192 (24 cars per level) 50 31.75 —I— + —J— 748.8 17,521.9 192 (24 cars per level) 66 31.42 a Storage capacity/total number of locations. b The values are upper bounds for the average retrieval time. Tables 2 through 4 show that the closed-form formu- with a larger system having more maneuvering space. las give a very precise approximation of the expected We do this by considering the live-cube system as retrieval times of real discrete systems, with errors less explained in the appendix. We investigate the trade- than 2%. In addition, the sensitivity analysis proves off between the system size and retrieval time for a that by varying the system parameters, the closed- fixed storage capacity. Table 5 illustrates the results. form formulas can still give very good approximations As Table 5 shows, in the base example, with few of the expected retrieval times of real discrete systems empty locations per level, virtual aisles can be created, even for extreme cases such as very small-size systems resulting in an average retrieval time of 33.5 seconds. (e.g., 100 storage locations). Moreover, we observe that To create more maneuvering space, the system size (in the gap between the results of the closed-form formula terms of footprint and volume) needs to be increased. and simulation decreases with an increase in capacity This results in a larger system with lower storage of the system (see Table 3). This is because the conti- density and higher investment and operational costs, nuity assumption becomes more accurate if the size of whereas the reduction in retrieval time is slight. The the system divided by the size of a storage location is reason is that the longer retrieval time due to a larger sufficiently large. system offsets the effect of the elimination of reposi- tioning move. This comparison shows that a live-cube 6.2. Impact of Number of Empty Locations system with max8—I—1 —J —9 empty locations per level can To obtain the closed-form formulas, we have achieve short response times to the customers while neglected the constant time required to create a vir- the available land is utilized efficiently. tual aisle for a requested unit load. However, depend- 6.3. Cost Analysis ing on the number of empty locations, a live-cube In this section, we first compare the investment and system might require zero, one, or two repositioning operational costs of a live-cube parking system with moves to create a virtual aisle (before the requested a traditional multistory car park of the same capacity. unit load is retrieved). For a live-cube system with According to the results, for a live-cube parking sys- 50% storage capacity utilization at each level (e.g., tem, the savings obtained from the smaller land area cars located on white fields of a chessboard), the might offset higher costs of construction compared to virtual aisle can be created while the requested car a traditional multistory car park. We then investigate is moving toward the lift location (no repositioning the impact of the system configuration on the invest- move). When the number of empty locations equals ment cost and response time of live-cube parking sys- —I— + —J — (number of rows plus number of columns), at tems. The results show that a system designer might most one repositioning move is required, and when prefer to slightly deviate from the optimal design the number of empty locations equals max8—I—1 —J —9, at because a significant reduction in investment can be most two repositioning moves are required to create achieved. a virtual aisle (see Online Supplement 1). Therefore, Cost comparison with traditional systems. To compare for the same storage capacity, a smaller live-cube sys- the total investment and operational costs of a live- tem with fewer empty locations needs to be compared cube system with a traditional system, we introduce Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 15 two parameters: construction cost ratio, , and foot- print ratio, . The construction cost ratio, , shows how costly a live-cube system is compared to a tra- ditional system in terms of technology and building. This ratio can be obtained as follows: ct 10,000 5 = (29) Cinc($) 0 ct 4 –10,000 ct where ct and are the investment in technology and 1.0 3 building per storage location for a traditional system and a live-cube system, respectively. The construction 1.5 2 cost for a live-cube system is generally higher than 1 the one for a traditional system; i.e.,  > 1. The foot- 2.0 print ratio, , compares the required footprint for each Figure 10 (Color online) Incremental Investment Cost cinc Graph for storage location in traditional and live-cube systems. Different Values of and Footprint ratio can be obtained as follows:  The parameters related to traditional systems such = (30)  as cost of land per storage location, cl ; average square meters of land required for each storage location, ; where  and  represent the average square meters and the cost of building and technology per storage of land required for each storage location for a tradi- location, ct , can simply be estimated because such tional system and a live-cube system, respectively. In systems are commonly used in practice. Therefore, a traditional system, each storage location normally in order to estimate the investment cost of a newly needs a larger footprint than in a live-cube system; designed live-cube system, it only suffices to estimate i.e.,  > 1. the land savings () and the additional construction The incremental investment cost for a storage loca- cost per location () of the live-cube system. tion, cinc , can now be formulated as follows (i.e., the We now consider one real medium-size multistory average investment cost of a traditional system minus car park described in Monahan (2012). The related the average investment cost of a live-cube system per parameters are ct = 16000 USD/slot,  = 297 m2 / storage location) slot, cl = 800 USD/m2 . Figure 10 illustrates the incre- mental investment cost function (cinc ) for different val- 1 cinc = cl  1 − + ct 1 −  (31) ues of  and . The dark colored area of the graph  includes scenarios where cinc is positive (i.e., the live- where cl is the cost of land per square meter. How- cube system is cheaper) and the yellow-colored area ever, in order to make a fair comparison, the opera- of the graph includes scenarios where cinc is negative tional costs (e.g., labor, maintenance, energy, etc.) also (i.e., traditional multistory car park is cheaper). need to be considered. The net present value (NPV) of For a live-cube parking system with the same the incremental operational cost, sinc , over the accept- capacity, the required footprint can be reduced able payback period,  max can be calculated as follows approximately 50% (i.e.,  = 2; see Table A.1 in the (incremental operational cost is the annual average appendix). As can be seen from Figure 10, up to a operational cost of the traditional system minus the 70% increase in the average construction cost per loca- annual average operational cost of a live-cube system tion ( ≤ 17) still leads to a live-cube system that is per storage location) cheaper than a traditional car park. Although live- cube parking systems might impose a higher annual  max sinc maintenance cost, they will significantly reduce CO2 NPVi  max = (32) 1 + it emissions and fuel consumption (see §6.3). t=0 Impact of system configuration on investment cost. where i is the discount rate and t is the time period. A live-cube system with optimal configuration min- To compare a live-cube system with a traditional sys- imizing the retrieval time does not necessarily need tem, the NPV of the incremental operational cost the minimum investment. In this section, we compare and the incremental investment cost (cinc ) have to be different live-cube system configurations in terms of considered. If cinc + NPV > 0, the live-cube system is retrieval time and (investment) cost. Table 6 compares cheaper and if cinc + NPV < 0, the traditional system systems with different configurations with the sys- is cheaper in terms of total cost (investment and oper- tem in the base example in terms of retrieval time ational). and footprint. To make a fair comparison in Table 6, Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 16 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS Table 6 Impact of System Configuration on Retrieval Time A system designer might prefer to design a system and Footprint with 12 height tiers instead of the system in the base Height (no. of tiers) a GapE6T 7 (%) b GapFootprint (%) example because the decrease in cost is more favor- able than the slight increase in retrieval time. 2 64 300 4 25 100 6.4. Impacts on Environment 6 8 33 The use of the latest technologies in a live-cube park- 12 6 −33 24 81 −67 ing system makes it a relatively “green” system. For a vertical movements, live-cube parking systems oper- 4GapE6T 7 = 4E6TGiven Height 7 − E6TBase Example 75/E6TBase Example 75 × 100%. b ate with electrically powered lifts and for horizontal GapFootprint = 44FootprintGiven Height −FootprintBase Example 5/FootprintBase Example 5 ×100%. transport, shuttles are used. A combination of lifts and shuttles leads to significantly reduced fossil fuel systems with a different number of height tiers are and energy consumption and CO2 emissions com- assumed to have the same capacity as the system in pared to traditional multistory car parks where the the base example. driver needs to drive the car inside the system. As Table 6 shows, the systems with fewer height Table 7 compares the average energy consumption tiers (i.e., two, four, and six tiers) than in the base and CO2 emissions (for storing and retrieving a car) example (eight tiers) both have longer retrieval times of the given live-cube system in the base example and and require a larger footprint area. For example, a sys- a traditional multistory car park of the same capacity tem with two height tiers is on average 64% slower to per level (a total of 192 cars) and for different types respond to customers and requires three times more of power plants generating the energy needed for land compared to the system in the base example. operation (lighting, ventilation, moving the cars). The The systems with more height tiers than in the base CO2 emissions of the live-cube parking systems and example (e.g., 12, 24) still have longer retrieval times traditional multistory car parks are calculated using but require a smaller footprint area. The material han- life cycle analysis. Parameters describing both stor- dling cost of such live-cube systems are mainly the age systems and detailed calculations are given in the cost of shuttles and lift. Because the storage capacity appendix. of all instances in Table 6 is the same as the capacity of As Table 7 shows, a live-cube parking system sig- the system in the base example, the material handling nificantly reduces the energy consumption and CO2 cost for all instances is approximately identical. In emissions compared to the traditional multistory car addition, systems with almost the same height have park. This savings is even more significant for CO2 approximately equal building costs (excluding cost of emissions if electricity is provided by a fossil-fuel land) in an urban setting (de Koster 1996). Therefore, power plant. In a live-cube parking system, the lifts for in Table 6, the cost of the system with 12 height tiers vertical movements require powerful engines to hoist is about 33% lower than the cost of the system in heavy cars. However, the energy consumption and the base example because it takes up 33% smaller CO2 emissions of such systems are significantly lower land compared to the system in the base example. than traditional systems for two main reasons: (1) the On the other hand, the average retrieval time of the average total travel time of a live-cube system is much system with 12 height tiers is 6% longer than the aver- shorter than a traditional system and (2) because in age retrieval time of the system in the base example. live-cube parking systems car engines remain off and Table 7 Energy Consumption and CO2 Emissions of Live-Cube Parking System and Traditional Multistory Car Park Generated by Fossil-fuel power plant Nuclear power plant Biomass-fuel power plant Live-cube Traditional Live-cube Traditional Live-cube Traditional Parking type parking car park parking car park parking car park Average CO2 272 4,369 2 200 1 184 emission (gram/car) Average S/R energy 0034 4094 0034 4094 0034 4094 consumption (kWh/car) Average lighting energy 0000 0025 0000 0025 0000 0025 consumption (kWh/car) Average ventilation energy 0000 5000 0000 5000 0000 5000 consumption (kWh/car) Total (kWh/car) 0034 10019 0034 10019 0034 10019 Note. Input data retrieved from Hyundai Elevator Co. Ltd. (2013) and Multi-Storey Car Park (2013). Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 17 the whole storage and retrieval process is performed as class-based storage policy and compare the results automatically, there is no need for exhaust ventilation with the results obtained here in this study. It is also and lighting. interesting to study the impact on makespan of opti- mally sequencing a group of retrievals. We know from AS/RS literature that savings of 20% to 70% can be 7. Conclusion and Future Research achieved compared to first-come-first-served sequenc- To make decisions on storage and retrieval solutions, ing (Han et al. 1987; Yu and de Koster 2012). As shut- managers have to trade off multiple criteria, including tles can move unit loads to the lift on different levels investment and operational costs, performance (aver- simultaneously, improvements might even be larger age and maximum response time), environmental for live-cube compact storage systems. Although we impact, land and space availability, and availability have studied live-cube compact storage systems with of labor. Compact, live-cube storage systems bring one and two lifts, results for live-cube compact stor- new solutions that often appear to outperform tra- age systems with multiple lifts may also prove worth- ditional solutions in space requirements, investment, while investigating. response time, and emissions. In general, operations in buildings with a larger footprint are not only more expensive than in smaller buildings but, in view of Supplemental Material Supplemental material to this paper is available at http://dx longer travel distances, are often also less responsive .doi.org/10.1287/trsc.2015.0586. and environmentally friendly. The retrieval time in the live-cube parking system is only 30% of that in Acknowledgments a traditional car park. However, admittedly, waiting In this research, the first author was partly supported by for a car to be retrieved might be perceived differ- Dinalog, in the project Cargo Driven Intermodal Transport. ently than walking in the garage and driving out. Still, The second author was supported by the Thousand Young research on these systems is just at its infancy and Scholar Program grants, the National Science Foundation of much work still needs to be done. As a first step, China (NSFC) for Distinguished Youth Scholars [71225002] in this study, we propose an evaluative framework and NSFC [Grants 71110107024 and 71121061]. to help storage facility managers in their decision making process. This framework includes response Appendix. Parameters Describing a Traditional time analysis of different live-cube system configura- Multistory Car Park and a Live-Cube tions, system configuration optimization, investment Parking System and operational cost comparisons with traditional Table A.1 presents the parameters describing a multistory systems, and energy consumption and CO2 emission car park and a live-cube parking system (base example). Table A.2 gives the parameters describing the live-cube calculations. parking system with a varying number of empty locations. Several research questions regarding the live-cube We consider a typical multistory car park with an effi- compact storage system remain open. Although we cient layout without external ramps (for the layout, see consider a live-cube parking system in our case study, Multi-Storey Car Park 2013). In addition, to make a fair a stronger case should still be made for other compet- comparison, we assume both alternatives have the same ing systems. It is possible to study the live-cube com- capacity per level (24 cars). The building dimensions can pact storage system with other storage policies such be calculated based on parking slot dimensions, system Table A.1 Parameters Describing a Traditional Multistory Car Park and a Live-Cube Parking System Average Parking Capacity Parking slot Building Engine Speed travel time system (no. of cars) dimensions (mm) dimensions (m) power (kW) (m/min) per S/R (sec) Live-cube parking 192a (24 cars per level) width = 3,250 width = 19.5 1.5 (per level for Shuttle speed in 66 x-movement), x-direction (vx 5 = 45, length = 6,400 length = 32.0 3 (per level for Shuttle speed in y -movement), y -direction (vy 5 = 100, height = 2,925 height = 23.4 22 (per lift for Lift speed 4vz 5 = 90 z-movement) Traditional car park 192 (24 cars per level) width = 3,500 width = 25.0 100 (average engine Driving speed = 300 178 power of a car), length = 6,500 length = 51.0 30 (ventilation per level), height = 3,300 height = 26.4 1.5 (lighting per level) a Capacity of the system = total number of locations (length × width × height = 6 × 5 × 8 = 240)—empty locations required to create virtual aisles (required empty location per level × height in number of levels = 6 × 8). Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage 18 Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS Table A.2 Parameters Describing a Live-Cube Parking System for Varying Number of Empty Locations Building Building No. of empty Capacity dimensions Building dimensions locations (no. of cars) (no. of tiers) dimensions (m) (sec) max8—I—1 —J—9 192 (24 cars per level) width = 5 width = 19.5 width = 15.36 length = 6 length = 32.0 length = 21.66 height = 8 height = 23.4 height = 23.4 50% 192 (24 cars per level) width = 6 width = 26.0 width = 19.20 length = 8 length = 38.4 length = 30.30 height = 8 height = 23.4 height = 13.65 —I— + —J— 192 (24 cars per level) width = 6 width = 19.5 width = 19.20 length = 6 length = 38.4 length = 21.66 height = 8 height = 23.4 height = 13.65 Note. The building dimensions are based on Table A.1 and have been chosen to meet the number of empty locations constraint whereas the length and width should be approximately equal, to minimize retrieval time. capacity, and aisle width in traditional car parks. Conse- de Koster MBM (1996) Ontwerp van een eenvoudig magazijn. quently, the average storage and retrieval time of a car Duijker JP, de Koster MBM, Ploos van Amstel MJ, eds. Prak- can be obtained based on building dimensions, shuttles, tijkboek Magazijnen/Distributiecentra (Kluwer, Deventer, The Netherlands), 3.10A.01–3.10A.20. lift speeds (in live-cube systems), and maximum allowed de Koster MBM, Le-Duc T, Yu Y (2008) Optimal storage rack design driving speed (in car parks). By using average storage and for a 3-dimensional compact AS/RS. Internat. J. Production Res. retrieval times (in hours) and engine power (in kW), the 46(6):1495–1514. average energy consumption (in kWh) for each alterna- Eweco (2014) Space parking optimization technology (SPOT). tive can be calculated (considering the average energy con- Accessed February 1, 2014, https://www.youtube.com/ watch?v=NozbyQJAnlQ. sumption to create virtual aisles in a live-cube system). In Eynan A, Rosenblatt MJ (1994) Establishing zones in single- addition, for a multistory car park, lighting and ventilation command class-based rectangular AS/RS. IIE Trans. 26(1): energy consumption have to be considered. To obtain the 38–46. average lighting and ventilation energy consumption per EZ-Indus (2013) UCW container storage system. Accessed Febru- ary 1, 2013, http://www.ezindus.com. car, we assume a 16-hour shift and (on average) a four- Goetschalckx M, Ratliff HD (1990) Shared storage policies based hour duration of stay per car. This means each level stores on the duration stay of unit loads. Management Sci. 36(9): and retrieves 96 cars per shift (at each level, 24 cars per 1120–1132. four-hour period). During a shift, 24 (1.5 kW × 16 h) and Gue KR (2006) Very high density storage systems. IIE Trans. 38(1): 480 kWh are consumed for lighting and ventilating each 79–90. Gue KR, Kim BS (2007) Puzzle-based storage systems. Naval Res. level, respectively. Consequently, the average lighting and Logist. 54(5):556–567. ventilation energy consumption per stored car can be calcu- Han MH, McGinnis LF, Shieh JS, White JA (1987) On sequencing lated. Input data for a live-cube system and a multistory car retrievals in an automated storage/retrieval system. IIE Trans. park are retrieved from Hyundai Elevator Co. Ltd. (2013) 19(1):56–66. and Multi-Storey Car Park (2013), respectively. The CO2 Hausman WH, Schwarz LB, Graves SC (1976) Optimal storage assignment in automatic warehousing systems. Management emissions of both alternatives can be obtained by using life Sci. 22(6):629–638. cycle analysis for different types of power plants. For one Hyundai Elevator Co. Ltd. (2013) Hyundai cart type auto park- kWh energy consumption provided by fossil fuel, nuclear, ing system. Accessed February 1, 2013, http:// hyundaielevator and biomass fuel power plants, 80, 6, and 3 grams of CO2 .co.kr/eng/parking/car/automobile.jsp. are emitted to the environment, respectively (see Manicore Johnson ME, Brandeau ML (1996) Stochastic modeling for auto- mated material handling system design and control. Trans- 2013). Therefore, based on the energy consumption in the portation Sci. 30(4):330–350. live-cube parking and multistory car park, the average CO2 Lee MK, Elsayed EA (2005) Optimization of warehouse storage emissions per stored car can be calculated. capacity under a dedicated storage policy. Internat. J. Production Res. 43(9):1785–1805. Manicore (2013) A tool for companies and office activities: The ”Carbon Inventory” of ADEME. Accessed February 1, References 2013, http://www.manicore.com/anglais/missions_a/carbon Automotion Parking Systems (2013) Park, swipe, leave systems. _inventory.html. Monahan D (2012) Man vs. machine: Is robotic parking right for Accessed February 1, 2013, http://www.automotionparking. your project? Internat. Parking Institute 34–37. Accessed July 5, com/index.php. 2012, http://www.walkerparking.com/wp-content/uploads/ Bartholdi JJ, Hackman ST (2011) Warehouse and distribution sci- 2013/05/Monahan-Robotic-Parking-Pkg-Prof-Sept-2012.pdf. ence: Release 0.95. School of Industrial and Systems Engi- Multi-Storey Car Park (2013) Multi storey car parks. Accessed neering, Georgia Institute of Technology, Atlanta, http://www February 1, 2013, http://www.multi-storey-car-parks.com. .warehouse-science.com. Navpak (2013) Shipboard package handling system. Accessed Bozer YA, White JA (1984) Travel-time models for automated stor- February 1, 2013, http://www.agilesystems.com/navpak age/retrieval systems. IIE Trans. 16(4):329–338. %20frame.htm. Zaerpour, Yu, and de Koster: Evaluating and Optimizing Live-Cube Compact Storage Transportation Science, Articles in Advance, pp. 1–19, © 2015 INFORMS 19 Navstors (2013) Naval automatic stowage and retrieval system. Stadtler H (1996) An operational planning concept for deep lane Accessed February 1, 2013, http://www.agilesystems.com/ storage systems. Production Oper. Management 5(3):266–282. navstors_frame.htm. Wöhr (2013) Wöhr parksafe. Accessed February 1, 2013, http:// OTDH (2013) Magic black box. Accessed February 1, 2013, http:// www.woehr.de/en/project/items/liverpool-parksafe-583.html. www.odth.be/nl/solutions-nl. Yu Y, de Koster MBM (2009) Optimal zone boundaries for two- Roodbergen KJ, Vis IFA (2009) A survey of literature on automated class-based compact three-dimensional automated storage and storage and retrieval systems. Eur. J. Oper. Res. 194(2):343–362. retrieval systems. IIE Trans. 41(3):194–208. Sari Z, Saygin C, Ghouali N (2005) Travel-time models for flow-rack Yu Y, de Koster MBM (2012) Sequencing heuristics for storing and automated storage and retrieval systems. Internat. J. Advanced retrieving unit loads in 3D compact automated warehousing Manufacturing Tech. 25(9–10):979–987. systems. IIE Trans. 44(2):69–87. Slocum J, Sonneveld D (2006) The 15 Puzzle Book (Slocum Puzzle Zaerpour N, Yu Y, de Koster MBM (2015) Storing fresh produce Foundation, Beverly Hills, CA). for fast retrieval in an automated compact cross-dock system. Sollohub D (2010) The architecture of parking. J. Architectural Ed. Production Oper. Management, ePub ahead of print January 25, 64(1):161–162. http://dx.doi.org/10.1111/poms.12321. View publication stats