1 Introduction

A fundamental problem in open-pit mine planning is the determination of the ultimate pit (UP), which consists of finding the contour of the mine that maximizes the difference between profits obtained from minerals minus extraction costs. The ultimate pit is very relevant in mine planning because it allows for a rough estimation of the value of a mine, and permits the determination of the region of interest on a mineral deposit. The problem is formulated in a very simplistic way, representing the mine as a set of uniform blocks (block model) with a predefined economic value and taking into account only structural requirements, such as the maximum slope angle on the pit walls. A first formulation of the UP problem can be traced back to 1965: the seminal paper (Lerchs and Grossman 1965) was one of the first works to propose an algorithm to solve this problem, and their methodology has been widely used in the mining industry. Since then, many authors have used different implementations and techniques to solve this problem in a variety of settings (Francois-Bongarcon and Guibal 1984; Caccetta 2007; Coléou 1988; Dagdelen 1986; Whittle Programming Pty Ltd 1998).

The UP problem is relevant in practice mainly for two reasons. First, it provides a rough estimate of the total ore and metal tonnage that can be potentially extracted, defining the region of interest of the deposit that deserves further analysis. Second, the possibility of solving UP problems efficiently allows one to generate so-called nested pits (Lerchs and Grossman 1965): by solving the problem for different prices—using “revenue factors”—of the mineral, it is possible to obtain a sequence of pits such that higher factors generate larger pits that contain the smaller ones. The main purpose of the nested pits technique is to derive phases (or pushbacks) used to guide the sequence of extraction of the mine for later production scheduling considering more operational constraints. In fact, nested pits have been the basis of most mine planning software during the last 30 years.

A crucial element of the UP problem is the economic value assigned to each block, which in general considers—at least—the grade, the price and the cost of extracting and processing a block to recover the minerals. The vast majority of works concerning the UP problem assume that these parameters of the problem are known. However, one of the most relevant sources of uncertainty comes from the estimation of the grades of minerals—commonly referred to as ore grade—in each block. Geostatisticians have developed several methodologies to study this problem, considering the strong spatial correlation that characterizes ore grade distributions in a deposit. The most common technique is to extract samples from the deposit via drill holes, to infer the distribution of the ore grades based on those samples and to generate scenarios using techniques such as kriging (Cressie 1990). However, most of this rich information about the deposit and its ore grades is not used on further steps of the mine planning process, and only one single grade value per block—in general, the average grade of each block—is utilized for solving the UP problem.

There are few works that incorporate the distribution of ore grades (or at least, a finite set of scenarios for the ore grades) into a mine planning problem. In Leite and Dimitrakopoulos (2007) and Consuegra and Dimitrakopoulos (2010), the authors show how the net present value of a mine can be improved by using techniques based in stochastic programming. In Marcotte and Caron (2013), expected profits are considered in the objective function, and the authors conclude that larger profits can be obtained by incorporating uncertainty into the UP formulation with respect to the deterministic case and that the relative gains of the stochastic approach increase with the treatment costs. In Espinoza et al. (2013), the authors move one step further and study the risk-averse UP problem, replacing the expected value by risk measures such as the conditional value-at-risk (CVaR) (Rockafellar and Uryasev 1997), which has been widely used in applications in finance and energy.

The idea behind the nested pits technique can be interpreted as a way of incorporating risk aversion into the mine planning problem. The smaller pits are those that are profitable, even if the ore grades (or metal prices) are multiplied by a small revenue factor. However, it is known that this simplistic approach does not provide the expected protection under bad scenarios (Espinoza et al. 2013). We propose a detailed study of the risk-averse UP problem, considering ore grades as an uncertain parameter of each block. We generate a set of nested pits by varying the risk level of the problem instead of considering revenue factors. First, we propose some properties—risk nestedness and additive consistency—that we believe a risk measure should have in the context of the UP problem and show that the only risk measure that satisfies those properties is the entropic risk measure. In particular, the CVaR does not satisfy any of these properties and is probably not the best choice to capture risk aversion in the mine planning context. Second, we derive conditions under which the entropic risk measure, which is a convex risk measure primarily used with utility functions over capital investments (Föllmer and Schied 2002), generates nested pits by varying the risk-aversion level of the decision maker.

We apply our methodology to a small case study of a self-constructed instance to illustrate the gains of our approach. We also verify our proposed methodology in a real-world mine, showing how ore grade uncertainty can drastically change the solutions obtained, while the common nested-pit approach is unaffected by the stochastic nature of the ore grade.

The structure of this paper is as follows. Section 2 provides an overall view of the UP problem, and we present a two-stage stochastic formulation and summarize the concept of nested pits. In Sect. 3, we define risk in the UP problem and propose desirable properties that risk measures should have. In Sect. 4, we show the results of our methodology using a small instance to visualize the differences between ours and the nested-pit methodology and using a real mine to show the potential of its use in real-life instances. Section 5 concludes the work and points out avenues for future research.

2 Ultimate pit problem and nested pits

2.1 Ultimate pit problem

2.1.1 Block model and its economic values

A block model is a representation of the deposit, constructed by a three-dimensional array of units—usually of equal size—called blocks (Sevim and Lei 1998). Each block has a set of parameters including—at least—a tonnage, representing the amount of material in each block and an ore grade \(g_b\), which determines the proportion of mineral in the block. We denote by B the set of blocks of the mine. From the position of each block, we can also estimate a set of precedences for each block, which is the set of blocks that must be extracted before a given block is removed from the mine, to avoid geotechnical difficulties. In practice, these precedences are computed to maintain a maximum slope angle on the walls of the pit (Khalokakaie et al. 2000). We denote by \(P\subseteq B\times B\) the set of all precedences; that is, if \((b,b') \in P\), then to order to extract block b, we must extract block \(b'\). Figure 1 presents an example of this representation.

Fig. 1
figure 1

Example of a block model of a mine (Sevim and Lei 1998)

The economic value of each block depends on its destination, which in turn depends, among other factors, on the metal grades of the block. In the simplest case, we assume that there is only one grade and two possible destinations for each block after extracting it: either to process it to recover the mineral or to drop it in a waste dump. Depending on this destination, we say that a block is either ore or waste.

Let \(c^e_b\) be the cost of extracting block b and \(c^p_b\) be the additional cost of processing block b after being extracted. In practice, these costs are usually a fixed value per ton of material, but we can consider more complex formulas, for example, by including contaminants or considering the hardness of the block. As a convention, \(c^e := \{c^e_b\}_{b \in B}\) and \(c^p := \{c^p_b\}_{b \in B}\) will be defined as the vector of all extracting costs and processing costs for every block \(b \in B\), respectively.

When processing a block, a fraction of its contained metal is recovered, which generates an income. We denote this income by \(r_b\cdot g_b\), where \(g_b >0\) is the ore grade of block b, and \(r_b >0 \) includes other factors needed to compute the profit. In the simplest case, \(r_b\) is the product of the tonnage of block b, the metal recovery rate of the process, and the price per unit of the recovered metal. Since we are interested in the study of ore grade uncertainty, we assume that \(r_b\) has a deterministic value for every block \(b\in B\).

With these parameters, the economic value \(v_b\) of an extracted block \(b\in B\) is given by the following formula:

$$\begin{aligned} v_b = {\left\{ \begin{array}{ll} r_b\cdot g_b - c_b^p - c_b^e & \text { if } b \text { is extracted and processed (ore),} \\ -c_b^e & \text { if } b \text { is extracted but not processed (waste).} \end{array}\right. } \end{aligned}$$
(2.1)

For the ultimate pit problem, the destination of each block can be decided a priori: we only process a block if its income pays off the processing costs (\(r_b\cdot g_b - c_b^p>0\)). In the case that \(r_b\) and \(c^p_b\) have a fixed value for all blocks in B, the threshold grade \(g_\text {cut}=c^p/r\) is known as the cut-off grade, and blocks are classified as either ore or waste depending if their grade \(g_b\) is greater than \(g_\text {cut}\).

2.1.2 UP with deterministic ore grades

The UP problem consists of selecting the set of blocks to extract and process that which maximizes the total value of the mine.

Using our previous notation, we can formulate the UP problem as the following mixed-integer optimization problem:

$$\begin{aligned} {\begin{array}{rl} \text {UP}(g,r) = \min \limits _{x^e,x^p}& \quad \sum \limits _{b \in B} (c_b^p - r_b g_b)x_b^p + c_b^e x_b^e\\ s.t.&\quad x_{b}^e \le x_{b'}^e \ \forall (b,b') \in P, \\ &\quad x_b^p \le x_b^e \ \forall b\in B, \\ &\quad x^e_b,x^p_b \in \{0,1\} \ \forall b \in B, \end{array}} \end{aligned}$$
(2.2)

where the decision of extracting and processing each block \(b \in B\) is represented by the binary variables \(x_b^e\) and \(x_b^p\), respectively, and \(g := \{g_b\}_{b \in B}\) represents the vector of ore grades for every block \(b \in B\). We remark that we wrote the problem as a minimization of cost, instead of maximization of profit (the negative of the objective function in (2.2)), to be consistent with most of the risk literature. The first set of constraints represents the extraction precedences for every block, and the second set of constraints conditions the processing of a block to its extraction.

As mentioned before, variables \(x^p\) can be eliminated from the problem because the destination of a block can be decided a priori, so we can consider the extraction decision \(x^e_b\) only and replace its objective coefficient by \(-v_b\) (see (2.1)). Nonetheless, since this simplification is not always valid for the stochastic version of the UP problem, we will use formulation (2.2) as a starting point in the next section.

2.1.3 UP with uncertainty on ore grades

The first step in developing a risk-averse model is to define a way in which to handle the ore grade uncertainty of the mine. For each \(b\in B\), let \({\tilde{g}}_b\) be the random variable that represents the ore grade of block b. Note that the decision of whether a block is ore or waste depends on the realization of \({\tilde{g}}_b\), so we cannot decide a priori its economic value as in the deterministic case. Hence, it is natural to introduce these random ore grades in the UP problem as a two-stage stochastic optimization problem. In the first stage, a decision must be made before the uncertainty is revealed, and a second-stage decision is made once the uncertainty is known, as a recourse action. In the context of the UP problem, we define the extraction decision as a first-stage variable and the processing decision as a second-stage variable.

Therefore, the UP problem under uncertainty can be formulated as follows:

$$\begin{aligned} {\begin{array}{rl} \text {UP}_U = \min \limits _{x^e}& \quad \sum \limits _{b \in B} c_b^e x_b^e + \rho _{\alpha }\left( Q(x^e,{\tilde{g}}, r)\right) \\ s.t.&\quad x_b^e \le x_{b'}^e \ \forall (b,b') \in P,\\ &\quad x_b^e \in \{0,1\} \ \forall b \in B, \end{array}} \end{aligned}$$
(2.3)

where

$$\begin{aligned} {\begin{array}{rl} Q(x^e, {\tilde{g}}, r) = \min \limits _{x^p}& \quad \sum \limits _{b \in B} (c_b^p - r_b {\tilde{g}}_b)x_b^p\\ s.t. &\quad x_b^p \le x_b^e \ \forall b \in B,\\ &\quad x_b^p \in \{0,1\} \ \forall b \in B, \end{array}} \end{aligned}$$
(2.4)

and \({\tilde{g}} := \{{\tilde{g}}_b\}_{b \in B}\) represents the vector of random variables realizations of ore grades for every block \(b \in B\). Function \(\rho _{\alpha }: L^1 \rightarrow \mathbb {R}\) will be defined as a deviation or risk measure, which deals with the uncertainty of the ore grade in the second stage. The parameter \(\alpha \in {\mathbb {R}}\), whose range depends on the risk measure \(\rho _{\alpha }\) under consideration, represents the risk aversion level of the decision maker. We note that if \({\tilde{g}}_b\) has a deterministic value, then we recover formulation (2.2).

If \(\rho _\alpha (\cdot )\) is equal to the expected value \(\mathbb {E}[\cdot ]\), then we fall into the risk-neutral case, which was extensively studied in Marcotte and Caron (2013). One of their main results is that the objective value of \(\text {UP}_U\) is greater than or equal to the objective value of the deterministic problem when using the expected grade \(\mathbb {E}[{\tilde{g}}_b]\) of each block. This result, which can be seen as an application of Jensen’s inequality due to the convexity of function \(Q(x^e, {\tilde{g}}, r)\), shows that even under risk neutrality, the stochastic UP problem can provide a different—and more profitable—pit (in expected value) than its deterministic counterpart.

2.2 Nested pits in open pit mine planning

As mentioned before, nested pits are a well-known concept that serve as the basis for several software programs used by the industry. Given a sequence of revenue factors \(\beta _1< \beta _2< \ldots < \beta _M\), with \(\beta _i\ge 0\), we solve the deterministic problem \(\text {UP}(g,\beta _i\cdot r)\) for each \(i=1\ldots M\). The solution of these problems satisfies the nestedness property, where the resulting pit is obtained with revenue factor \(\beta _i\) is completely contained in the pit obtained by using revenue factor \(\beta _{i+1}\) (see Lerchs and Grossman 1965; Caccetta and Hill 2003; Marcos et al. 2015 for a formal proof of this property). Figure 2 shows the pits generated by solving a small instance of problem (2.2) with \(\beta _1< \beta _2< \beta _3 < \beta _4\).

Fig. 2
figure 2

Example of nested pits

Note that for \(\beta =0\), the resulting problem \(\text {UP}(g,\beta \cdot r)\) only considers the extraction and processing costs, so the optimal solution is to not extract any block. For \(\beta =1\), we recover the classical deterministic UP problem. Note that two consecutive revenue factors \(\beta _i<\beta _{i+1}\) do not necessarily result in different pits, implying that the number of different nested pits can be smaller than M.

Nested pits are not explicitly defined as a way to deal with uncertainty. They are used as input in a scheduling problem that defines the extraction and processing times for each block, incorporating other constraints such as maximum tonnage capacities per period, discounted cash flow over different periods, and other operational constraints. The idea behind nested pits approaches is that smaller pits (obtained using smaller revenue factors) are more profitable because they have a positive value, even if the ore grades (or mineral prices) are multiplied by these small factors. Finally, the interplay between revenue and discount factors can be interpreted as a way of risk aversion from the decision maker’s perspective. Parts of the mine that remain financially attractive even after being penalized by a revenue factor should be extracted first to avoid further devaluation by the discount factor.

In the last decade, several papers have questioned the use of nested pits for mine planing and production scheduling, replacing it by direct block scheduling of the problem (see Letelier et al. 2020 for a more detailed description of the state-of-the-art developments in this area). Nevertheless, the idea of nested pits is still broadly used due to its simplicity and practicality of having a sequence of pits as guidance for the successive steps of the mine planning process.

A key contribution of our work is to propose an adequate risk measure \(\rho _{\alpha }\) for mining problems and to be able to solve problems 2.32.4 efficiently. Given that we incorporate ore uncertainty, we want to derive conditions under which we can replicate the well-known idea of nested pits but consider different risk levels by varying the parameter \(\alpha \) to obtain these pits instead of using revenue factors.

3 Risk measurement

3.1 Defining risk under ore grade uncertainty

For the remainder of this paper, we focus on ore grade uncertainty. Nevertheless, most of the results presented can be easily extended to other sources of uncertainty. To select the appropriate risk measure for the UP problem, we need to define which we believe are the desirable characteristics of the pits generated by a given risk measure \(\rho _{\alpha }\) when we vary the parameter \(\alpha \).

As previously discussed, our intention is to study risk measures that provide nested pits with different risk aversion levels; that is, more conservative pits must be contained in riskier pits. The main justification for this property comes from the practical interpretation of the nested pits. Nested pits are widely used as a basis for the design of phases. In deterministic or risk-neutral settings, more profitable pits must be extracted earlier due to the discount factor of future cash flows of the mine. In a risk-averse context, a set of nested pits with an ordered level of risks can be used in the same way, providing a guideline on how to sequentially exploit the mine, from safer to riskier pits. If considering different risk levels provides different unrelated or unnested pits, then the intuition on how to exploit the deposit is unclear.

Our proposed definition of the risk nestedness property is as follows:

Definition 1

Assume that the level of risk aversion of the decision maker rises when the value of \(\alpha \in \mathbb {R}\) increases and that the ore grades are independently distributed. A risk measure \(\rho _\alpha \) is risk nested for the UP problem if for \(\alpha _1 > \alpha _2\) the set of extracted blocks in the optimal solution of problem (2.3) obtained by using \(\rho _{\alpha _1}\) is contained in the set of extracted blocks obtained using \(\rho _{\alpha _2}\).

Please note that ore grade independence between blocks is not a realistic assumption in a mining context. In contrast, geostatistical models for ore grades consider dependence between blocks that are physically close. However, we assume this independence to provide formal mathematical proofs of results that as we will see in the computational experiments, are also observed for realistic dependent distributions of the ore grades.

A second desirable property is related to the consistency of extracting decisions. Decisions over a set of blocks should not depend on the uncertainty of independent elements within the mine. More precisely, assuming independence of the ore grade distributions, if a pit \(U'\) is riskier than other pit \(U''\) and we add to both pits a completely unrelated block \(b\notin U'\bigcup U''\), then it is desirable that pit \(U'\bigcup \{b\}\) is still riskier than pit \(U''\bigcup \{b\}\). This property is known as the additive consistency property of a risk measure and can be defined as follows:

Definition 2

Let XY and Z be random variables, where Z is independent of both X and Y, and \(\rho _\alpha (\cdot )\) be a risk measure where \(\alpha \) is the risk level of the decision maker. If \(\rho _\alpha (X) < \rho _\alpha (Y)\) and \(\rho \) is additive consistent, then

$$\begin{aligned} \rho _\alpha (X+Z) < \rho _\alpha (Y+Z). \end{aligned}$$

We suggest these two properties as desired for a given risk measure. In the next subsection, we discuss about different risk measures in the literature, and we study whether they satisfy these desired properties.

3.2 Study of some risk measures for the stochastic UP problem

In the last 20 years, risk measures have been extensively studied, and their theoretical properties are well established. Motivated by applications in energy, finance, reliability, among others, several risk measures have been used in real-world problems, each of which has specific properties. One of the most popular in the optimization literature is the conditional value-at-risk (CVaR). Despite its popularity in finance and energy problems, we present somewhat surprising results that show that the CVaR does not satisfy risk nestedness or additive consistency.

In Cominetti and Torrico (2016), the authors proved that the only risk measure that is additive consistent is the entropic risk measure. Motivated by that result, we study this risk measure and discuss its suitability for the risk-averse UP problem.

3.2.1 The CVaR risk measure

The CVaR is a commonly used risk measure due to its tractability. Following (Rockafellar and Uryasev 1997), the CVaR can be defined as

$$\begin{aligned} \text {CVaR}_{\alpha }(X) = \min _{\eta \in {\mathbb {R}}} \left\{ \eta + \frac{1}{1-\alpha } \mathbb {E}\left[ (-X-\eta )^+\right] \right\} , \end{aligned}$$
(3.5)

where X is a random variable, \((\cdot )^+\) is the positive part function, and \(\alpha \in [0,1[\) is the risk level. A value of \(\alpha =0\) corresponds to the risk-neutral case (the expected value), and as \(\alpha \) approaches one, the risk measure protects against the worst-case realization of X.

CVaR has been previously used as a risk measure for the UP problem in Amankwah et al. (2013) and Lagos et al. (2011). However, the next proposition shows that it may not be the ideal risk measure for stochastic UP problems under our requirements.

Proposition 1

The CVaR risk measure is not risk-nested for the UP problem.

Proof

We will prove that the CVaR is not risk-nested with a counterexample: let B bet the set of blocks containing only three blocks \(b_X\), \(b_Y\) and \(b_Z\), and assume that to extract block \(b_X\) or \(b_Y\), we must extract block \(b_Z\) (i.e. \(b_Z\) is a precedence of \(b_X\) and \(b_Y\)).

We remark that if a random variable W is normally distributed with mean \(\mu \) and variance \(\sigma ^2\) (\(W \sim N(\mu , \sigma ^2)\)), then

$$\begin{aligned} \text {CVaR}_\alpha (W) = \mu + \frac{\sigma ^2}{(1-\alpha )\sqrt{2\pi }}e^{-z_{1-\alpha }^2/2}, \end{aligned}$$
(3.6)

where \(z_{1-\alpha }\) is the \((1-\alpha )\)-quantile of the standard normal distribution. Let the ore grade of the three blocks be independent and normally distributed, resulting in a final profit distribution given by the random variables \({X}\sim N(-2,18)\), \({Y}\sim N(1,0)\) and \({Z}\sim N(-4,10)\).

Table 1 Comparison of CVaR values for different pits

Table 1 shows the results of computing the CVaR using formula (3.6) for different values of \(\alpha \) and for different combinations of valid pits. It can be seen that with \(\alpha =0.1\), the optimal solution is to extract \(b_Y\) and \(b_Z\), but reducing the risk aversion level to \(\alpha = 0.05\), it is optimal to extract blocks \(b_X\) and \(b_Z\), which are not contained in the previous solution. This example shows that the CVaR is not risk nested. \(\square \)

3.2.2 The entropic risk measure

Let X be a random variable. The entropic risk measure of X at level \(\alpha \) is defined as

$$\begin{aligned} \rho ^{Ent}_\alpha (X) := {\left\{ \begin{array}{ll} \frac{1}{\alpha } \log {\mathbb {E}[e^{\alpha X}]} \quad \text {if } \alpha \ne 0,\\ \mathbb {E}[X] \quad \text {if } \alpha = 0. \end{array}\right. } \end{aligned}$$

If \(\alpha >0\), then the decision maker is risk-averse, \(\alpha < 0\) means risk-seeker, and \(\alpha = 0\) corresponds to the risk-neutral (expected value) case. We will focus on nonnegative values of \(\alpha \).

The entropic risk measure (\(\rho ^{Ent}\)) satisfies additive consistency (Cominetti and Torrico 2016), and it is a convex function, making it suitable for optimization problems. However, it is not a coherent risk measure (Artzner et al. 1999; Cominetti and Torrico 2016).

We claim that \(\rho ^{Ent}\) is a viable risk measure for the stochastic UP problem since it also satisfies the risk nestedness property:

Proposition 2

The entropic risk measure is risk-nested for the UP problem.

To prove this result, we use two technical lemmas from the entropic risk measure:

Lemma 1

Let \(\alpha _1,\alpha _2 \in {\mathbb {R}}\), with \(\alpha _1 > \alpha _2\). For any random variable X, we have that \(\rho ^{Ent}_{\alpha _1}(X) \ge \rho ^{Ent}_{\alpha _2}(X)\).

Lemma 2

The entropic risk measure is translation invariant; that is, if X is a random variable and \(a \ge 0\), then

$$\begin{aligned} \rho ^{Ent}_\alpha (X+a) = \rho ^{Ent}_\alpha (X)+a. \end{aligned}$$

Lemma 3

Let XY be two independent random variables, and \(0< \alpha < \infty \). We have

$$\begin{aligned} \rho ^{Ent}_\alpha (X+Y) = \rho ^{Ent}_\alpha (X) + \rho ^{Ent}_\alpha (Y). \end{aligned}$$

The proofs of Lemma 1 and Lemma 2 are presented in "Appendix A", and for the proof of Lemma 3, we refer the reader to Cominetti and Torrico (2016).

Proof of Proposition 2

Let \({\tilde{g}}_b\) be the random variable of the ore grade for each block \(b \in B\), assume they are independently distributed, and let \(\alpha _1, \alpha _2 \in \mathbb {R}\) be two risk levels where \(0< \alpha _2< \alpha _1 < \infty \). Let \(U_1\) and \(U_2\) be the resulting pits obtained by solving (2.3) using the entropic risk measure \(\rho ^{Ent}_{\alpha _1}\) and \(\rho ^{Ent}_{\alpha _2}\), respectively. We assume that \(U_1\) and \(U_2\) are minimal solutions for these problems.

By Lemma 2, we know that the objective function of problem (2.3) can also be written as \(\rho ^{Ent}_\alpha (\sum _{b\in B} c_b^e+x_b^e + Q(x^e,{\tilde{g}},r))\). To simplify the notation, we denote by \(\rho ^{Ent}_\alpha (U)\) the value of this objective function when evaluated on the indicator vector of U (that is, \(x^e_b=1\) if and only if \(b\in U\)).

By contradiction, assume that \(U_1 \nsubseteq U_2\). Since extracting an empty pit is a feasible solution for the UP problem and \(\rho ^{Ent}_\alpha (\emptyset )=0\) for any \(\alpha > 0\), we know that \(\rho ^{Ent}_{\alpha _1}(U_1) < 0\). Furthermore, since \(U_1\cup U_2\) is also a valid pit (but not optimal) for the UP problem with risk measure \(\rho ^{Ent}_{\alpha _1}\), Lemma 3 and the minimality of \(U_1\) implies that \(\rho ^{Ent}_{\alpha _1}(U_1 \backslash U_2) < 0\).

By Lemma 1, since \(\alpha _2 < \alpha _1\), then \(\rho ^{Ent}_{\alpha _2}(U_1 \backslash U_2) \le \rho ^{Ent}_{\alpha _1}(U_1 \backslash U_2) < 0\). However, since \(U_2 \cup (U_1 \backslash U_2)\) is also a feasible solution for the UP problem, by Lemma 3, we have that

$$\begin{aligned} \rho ^{Ent}_{\alpha _2}(U_2 \cup (U_1 \backslash U_2)) = \rho ^{Ent}_{\alpha _2}(U_2) + \rho ^{Ent}_{\alpha _2}(U_1 \backslash U_2) < \rho ^{Ent}_{\alpha _2}(U_2), \end{aligned}$$

which is a contradiction of the optimality of \(U_2\), completing the proof. \(\square \)

These results show that the entropic risk measure is a good candidate to be utilized for solving the stochastic UP problem. In the next section, we will exemplify our approach using different configurations of uncertainty in the ore grade of the blocks of a mine, and compare the pits obtained to the classical nested pits approach.

4 Computational results

We present two sets of results. First, we consider a small-sized mine that was designed to illustrate the properties and methodologies that we propose in this paper. The second mine is a real-life instance with a large number of blocks to test the applicability of our methodology in a more realistic environment.

We provide a comparison of the results using different risk measures and discuss the practical implications and managerial insights of the pits obtained with each method. We implemented all models in Python 3.7.2 using Gurobi 8.1 as a solver for the resulting optimization problems. All computations were performed on an Intel(R) i7 CPU 7700HQ @ 2.80 GHz and 12 GB of memory using Windows 10.

4.1 Small mine

4.1.1 Mine description and parameters of the problem

We will study a mine with 66 blocks in a two-dimensional configuration, whose topology is described in Fig. 3. The precedences are as follows: if a block is to be extracted, then the blocks on the top, right and left of the block on top must be extracted, emulating a 45-degree slope angle (see block b in Fig. 3). We assume that only a set B of blocks (indicated as blue blocks in the figure) can be extracted.

Fig. 3
figure 3

Topology of the small mine: blue blocks are part of set B, yellow arrows represent the precedences of block b, and white blocks are not being considered for extraction

In this simple example, to avoid approximations of the distribution and simplify the computations, we assume that the ore grade \({\tilde{g}}\) follows a multivariate normal distribution \(N(\mu , \varSigma )\), where \(\mu \) is the vector of mean values and \(\varSigma \) is the variance-covariance matrix. In this case, it can be proven that problem (2.3) for the entropic risk measure can be formulated as a mixed integer program with a convex quadratic objective function (see "Appendix A.3"), which can be solved by most modern optimization solvers.

The resulting equivalent formulation of the UP problem with \(\rho ^{Ent}\) is given by

$$\begin{aligned} \min _{x^e,x^p \in X^{EP}} \sum _{b \in B} c_b^e x_b^e + \sum _{b \in B} (c_b^p - r_b \bar{g}_b) x_b^{p}+ \frac{1}{2} \alpha \sum _{b\in B}\sum _{b' \in B} r_b x^p_b r_{b'}x^p_{b'} \varSigma _{bb'}, \end{aligned}$$
(4.7)

where \(X^{EP}\) is the same feasible set of solutions for problem (2.3).

The grade distribution \({\tilde{g}}\) of the blocks is illustrated in Fig. 4: all blue blocks are waste, with constant grade equal to zero, and the yellow block has a constant grade of 3. The grades of the red blocks follow a multivariate normal distribution, with means \((\mu _1,\mu _2,\mu _3,\mu _4)= (4, 4, 5, 6)\), from top to bottom, a variance \(\varSigma _{i,i}=\mu _i^2\) and covariances \(\varSigma _{i,j}=0.2 \cdot \mu _i\cdot \mu _j\) for all pairs (ij). The remaining parameters of the problem are \(c^e_b = 0.5\), \(c^p_b = 0\) and \(r_b = 1\) for all \(b \in B\).

4.1.2 Results for the stochastic UP using \(\rho ^{Ent}\)

Fig. 4
figure 4

Grade structure of a small mine: blue blocks are waste, yellow blocks have no variance, and red blocks have positive covariance

The different pits resulting from solving the UP problem with \(\rho ^{Ent}\) for (1) \(\alpha =0\), (2) \(\alpha =0.02\), (3) \(\alpha =0.04\), (4) \(\alpha =0.1\) and (5) \(\alpha = 0.2\) are shown in Fig. 5a. We recall that \(\alpha =0\) represents the risk neutral case, that is, maximizing the expected profit of the pit. It can be seen that in the risk-neutral case the whole mine is extracted. If we increase the risk-aversion level \(\alpha \), then the optimal UP becomes smaller and contained in the previous UP. Under high levels of risk aversion (\(\alpha =0.2\)) the pit only includes the four top-leftmost red blocks because this pit has a positive profit and zero variance, so there is no risk.

Fig. 5
figure 5

Optimal stochastic UP for \(\rho ^{Ent}_\alpha \) with different values of \(\alpha \) (in different colors) and covariance matrices

It is interesting to note the evolution of the pits when we multiply the parameter \(\varSigma \) by some positive constant: as variability grows, the optimal pit configuration changes for the same values of \(\alpha \). Figure 5b, c show the evolution of the pits when we multiply the covariance matrix \(\varSigma \) by 2 and 10, respectively. Note that if some numbers within the blocks are missing, it is because the pit is the same as the solution obtained by using a higher value of \(\alpha \), i.e., in Fig. 5b the optimal solutions using \(\alpha =0.1\) and \(\alpha =0.2\) are the same. In Fig. 5c we see that for any value of \(\alpha \ge 0.02\), all pits obtained extract the yellow block from Fig. 4 and avoid the red blocks.

4.1.3 Results for the classical nested pit approach

As a comparison, the same analysis can be performed using classical nested pit methodology, namely, solving problem (2.2) for different revenue factors \(\beta \). Note that this methodology assumes a deterministic grade for each block, so we use the mean values \(\mu \) as the ore grade of the random red blocks. Figure 6 shows some nested pits obtained for different ranges of revenue factors \(\beta \) as follows:

  1. 1.

    \(\beta \in [0.91\overline{6}, 1]\),

  2. 2.

    \(\beta \in [0.9, 0.91\overline{6}[\),

  3. 3.

    \(\beta \in [0.875, 0.9[\),

  4. 4.

    \(\beta \in [0.\overline{6}, 0.875[\),

  5. 5.

    No blocks are extracted: \(\beta < 0.\overline{6}\).

Fig. 6
figure 6

Classical nested pits for the small mine

As expected, it can be seen that using revenue factor \(\beta =1.0\), we obtain the same pit as in the risk-neutral when the red block is random. Interestingly, the resulting nested pits for other values of \(\beta \) do not consider extracting the less riskier pit on the upper left corner (those with the number 5 in Fig. 5). In contrast, for smaller revenue factors (\(\beta < 0.\overline{6}\)), the optimal solution is to not extract any blocks.

Note that this methodology does not consider the variance of the grades, so we will obtain the same optimal pits even in scenarios of high uncertainty. The classical nested pit approach is not designed for dealing with geological uncertainty, and given its invariance with respect to higher variability levels, it should not be used for this purpose. As we show in Fig. 6, under grade uncertainty, the risk-free profitable pit could not be detected by the nested pit approach, a pit that \(\rho ^{Ent}\) was able to identify.

We can clearly see how variability becomes an important factor with the \(\rho ^{Ent}\) approach and how it is controlled by the choice of \(\alpha \): conservative values aim for smaller but less variable expected profits. We also remark that similar to revenue factors for the classical nested pits, considering a very high level of risk aversion results in not extracting any block (or only a set of blocks that provides a positive revenue for any realization of its grade values). On the other extreme, a zero risk level is equivalent to the risk-neutral case, resulting in the classical ultimate pit problem. The entropic risk measure can be of practical use in the mine planning community especially due to its risk nestedness property, providing a potential alternative way for defining phases when the ore grade uncertainty is relevant for the problem. By inspecting the resulting nested pits, decision makers with different levels of risk aversion can define the mine contour that is more suitable to their risk tolerance. Since the risk-averse approach takes into account uncertainty in the ore grades, pits obtained for values of \(\beta \) close to zero, as shown in Fig. 2 might be left unexplored.

4.2 Andina case

The Andina is a copper mine that belongs to Codelco, Chile’s state-owned company that controls approximately 20% of global reserves of copper. It is located in Rio Blanco (approximately 80 km NE of Santiago) and is still active after 82 years. We will work with a sector of the mine represented by 26,400 blocks from the Sur-Sur open pit. The ore grade distribution is approximated by scenarios generated from a set of drilling holes via conditional simulations using the turning bands algorithm (Emery and Lantuéjoul 2006).

Since this mine does not follow a multivariate normal distribution of ore grades we cannot use the closed-form quadratic formulation (4.7). We attempted to solve the problem directly using nonlinear optimization solvers (MINOS) to check if we could handle the real mine with its 26,400 blocks. The solver was not able to close an optimality gap of 99.9% after 5 h running, even after reblocking the mine into less than 100 blocks. In summary, it is hopeless to try to attack the problem directly.

The presence of integer variables and the nonlinearity of function \(\rho ^{Ent}\) turns the risk-averse UP problem into an extremely challenging problem to solve for larger mines. However, given the convexity of the exponential function, we can use a piecewise linear approximation of the objective function. To accomplish this, we modify the objective function of (2.3) due to the translation invariance of \(\rho ^{Ent}\) (please refer to Lemma 2) and the monotonicity of the \(\log \) function. All necessary proofs and final model that approximate the value of problem (2.3) using \(\rho ^{Ent}\) can be found in "Appendix C".

4.2.1 Results

Given the computational cost of solving this problem with a large number of scenarios, we use 20 ore grade scenarios to solve the problem. Using more scenarios would lead to larger solving times, which requires specialized optimization techniques to solve the resulting model, and goes beyond the scope of this paper. Nevertheless, using this number of scenarios is not uncommon in the stochastic mine planning literature (Leite and Dimitrakopoulos 2007; Consuegra and Dimitrakopoulos 2010).

Our piecewise linear approximation was designed to be a uniform grid between integers \(-20\) and 20 with steps of 0.01 (4000 steps in total) because Gurobi treats values below \(10^{-8}\) as 0 and \(\exp (-20) \approx 2.06\cdot 10^{-9}\). We also computed a scaling factor and changed the NumericFocus parameter of Gurobi to its maximum value. With those changes, we gained more processing power in the numerical calculations and avoided numerical issues given the small numbers that will be used by our model. We use the mean absolute deviation (MAD) to check the performance of our approximation.

We denote by NP the resulting pits obtained by the nested pit methodology. To correctly approximate the true value of each solution, we evaluate the total profit of the resulting ultimate pit on a different out-of-sample set of 100 scenarios. This gives us a set of 100 values for each solution. We also compare these values with the lower bound (LB) obtained by solving problem (2.2) by fixing \(g_b\) to each scenario. LB will have a different pit and a negative profit per scenario, serving as a benchmark to the best possible result. Additionally, we present the coefficient of variation (CV), which is the ratio between the standard deviation and the average of the profits.

Table 2 shows the results for NP while changing the value of \(\beta \). If \(\beta \ge 0.7\), then the solutions obtained have less than 100 blocks extracted, and for \(\beta \ge 0.9\), no blocks are extracted, and we omit the results. A similar behavior was observed in the small mine example in Sect. 4.1: conservative values of \(\beta \) will avoid extraction altogether. However, as risk averseness grows, the CV becomes larger (6.9% on average).

Table 2 Table of descriptive statistics for NP of the out-of-sample results

Table 3 shows the results for \(\rho ^{Ent}\). We observe a larger number of blocks extracted in the pit of the neutral case, and using larger values of \(\alpha \), the pits obtained have fewer blocks. As risk averness grows, the value of the CV decreases (5.6% on average). The MAD for different values of \(\alpha \) are of order \(O(10^{-5})\) in the worst case, showing that our linear approximation is effective. The values for \(\alpha > 12.5\) are not present in the table since the model incurs numerical errors: the solution obtained was an empty pit.

Table 3 Table of descriptive statistics for \(\rho ^{Ent}\) of the out-of-sample results

Figure 7 shows a boxplot of the negative cost (therefore profits) distribution obtained by \(\rho ^{Ent}\), NP and LB. From a practical viewpoint we can see a higher dispersion of profits among different values of \(\beta \) and that NP with \(\beta \) equal to zero has a significant profit gap with respect to the lower bound (as noticed in Marcotte and Caron (2013)). The situation is different for the entropic: there is less dispersion among the different choices of the risk aversion level \(\alpha \), and a higher percentage of the maximum profit can be captured by \(\alpha \) values closer to zero.

Fig. 7
figure 7

Boxplots of the total cost distribution per model for the Andina case

Figure 8 shows a cross-sectional view of the pit, where we can observe the differences in the optimal solutions of \(\rho ^{Ent}\) and NP. We used \(\alpha =10\) and \(\beta =0.03\) to obtain pits with a similar number of blocks extracted (14,810 in \(\rho ^{Ent}\) versus 14,983 in NP). We can see that both models aim for different sections of the mine. In fact, NP extracts the block at the bottom of the central zone, which has a positive profit but with a high uncertainty. In contrast, \(\rho ^{Ent}\) avoids this zone, and prefers to go deeper in the eastern zone of the mine (left side of the figure).

Fig. 8
figure 8

Graphical comparison of the different resulting pits on a cross sectional view of the mine

The importance of these results is the effective application of an approximation method to prove \(\rho ^{Ent}\) as a viable tool for risk aversion in real life instances of mining optimization problems. In fact, these techniques are compatible with the current methodologies to solve large-scale instances of open-pit mine planning problems (Letelier et al. 2020). The error obtained in the approximation to the real value of the exponential function is of order \(O(10^{-5})\), which was within our acceptance tolerance, and Gurobi was able to solve the approximate formulation in an average time of under 5 minutes for all the values of \(\alpha \).

5 Conclusions

In this paper, we formulate and discuss how to solve a risk-averse version of the UP problem under geological uncertainty. By considering the block model, we incorporate ore-grade uncertainty and discuss the most appropriate risk-measure properties for the problem, especially if we are interested in replicating the use of nested pits for mine planning. We propose two properties we believe a risk measure should satisfy in this context: risk nestedness and additive consistency. Risk nestedness is a desirable property since modern mine scheduling algorithms use the idea of nested pits as an input for generating a sequence of block extractions over time. Additive consistency is also an intuitive property in the sense that it preserves preferences in the presence of independent waste blocks. An interesting consequence of our work is that one of the most popular risk measures, the conditional value-at-risk, fails to satisfy both properties even in the case of blocks with independent distribution, calling into question its use in mining.

In a small simulated mine, we contrast the results obtained using the revenue factor to obtain nested pits and the proposed entropic risk measure methodology. The former behaves as expected: we obtain different pits that vary from mining everything up to avoid working on the mine altogether in the most conservative cases. The proposed approach shows a different behavior: we obtain smoother pit transitions as we vary the level of risk aversion, with a focus on avoiding variance within the results. By generating higher variability scenarios in the entropic risk approach, we obtain smaller pits than low variability cases for the same value of \(\alpha \).

We apply our method to a real-life mine (Andina in Chile) with more than 25,000 blocks, which is challenging to solve directly. The entropic risk measure adds nonlinearity to the problem, which in most cases would make the problem intractable. We apply a linear approximation to efficiently solve the risk averse problem with the entropic risk measure, and our results show that the approximation errors were within tolerable margins, validating the linear approximation scheme. The pits generated by the entropic risk measure are completely different from those using revenue factors, and different choices of the risk parameter \(\alpha \) allow us to obtain better profits when compared to revenue factors.

Future work includes extensions to our numerical algorithm to be able to cope with larger mines with millions of blocks and the study of dynamic risk-averse models for mine scheduling problems.