A Novel Large-Scale Stochastic Pushback Design Merged with a Minimum Cut Algorithm for Open Pit Mine Production Scheduling

: Traditional optimization of open pit mine design is a crucial component of mining en-deavors and is inﬂuenced by many variables. The critical factor in optimization is the geological uncertainty, which relates to the ore grade. To deal with uncertainties related to the block economic values of mining blocks and the general problem of mine design optimization, under unknown conditions, the best ultimate pit limits and pushback designs are produced by a minimum cut algo-rithm. The push–relabel minimal cut algorithm provides a framework for computationally efﬁcient representation and processing of the economic values of mining blocks under multiple scenarios. A sequential Gaussian simulation-based smoothing spline technique was created. To produce push-backs, an efﬁcient parameterized minimum cut algorithm is suggested. An analysis of Indian iron ore mining was performed. The developed mine scheduling algorithm was compared with the conventional algorithm, and the results show that when uncertainty is considered, the cumulative metal production is higher and there is an additional increase of about 5% in net present value. The results of this work help the mining industry to plan mines in such a way that can generate maximum proﬁt from the deposits.


Introduction
Similar to other sectors, the mining sector uses production scheduling to determine the order in which mining blocks will be removed from the earth during various production periods during the life of the mine's existence to maximize the net present value (NPV).This is a well-known linear programming with mixed integers problem that is difficult to solve when the integer variables are large.In the process referred to as conventional, first, the optimal ultimate pit, i.e., final pit contour, is defined using the Lerchs-Grossman or L-G algorithm [1,2] so that the problem size can be reduced for mine production scheduling.To optimize the mine's annual profit while taking into account constraints such as pit slopes, the ultimate pit limit can be considered as the collection of blocks that must be removed.Selecting pit limits that maximize the gap between the overall mine value of extracted materials over the life of the mine and the total extraction cost is another way to think about maximizing the pit limits.The L-G technique applies graph theory to the problem and proves that the optimal solution to the ultimate pit problem equates to maximum graph closure.Although the purpose of the pit design decision is to reduce the size of the problem, in large mines, the reduced problem is occasionally unsolvable.
In practical scenarios, the large pit is therefore discretized into a series of small problems known as pushbacks so that mine production scheduling can be applied to those pushbacks.The generation of pushbacks is considered a long-term plan where certain constraints such as mining capacity constraints and plant capacity constraints should be respected to maximize profit from pushbacks.Pushbacks, a series of small pits, can be obtained by parameterization of the L-G algorithm [2,3].Picard [4] demonstrated that adopting a minimal cut technique to a related graph can solve the maximum closure problem of open pit mine planning.The final pit limit and pushbacks can be created more efficiently, as shown by Hochbaum and Chen [5], who used a minimum cut/maximum flow algorithm based on Goldberg's [6] push-relabel maximum flow/minimum cut approach.The results of Hochbaum's [7] detailed investigation of the L-G algorithm's theoretical complexity revealed that it is significantly more difficult than the push-relabel technique.However, these studies did not address how resource constraints are respected.In addition, the work is developed within the context of conventional optimization that does not consider uncertainty associated with different parameters; optimization is achieved assuming that all inputs are the actual values.This is a planning practice in the context of mining problems where block values are calculated from a limited sample of data and the price of metal is assumed to be constant over the mine's life.Due to uncertainty in orebody attributes, such as grades and metal in the ground, as well as market uncertainties, such as exchange rates, it has become clear that this assumption has major limitations and may lead to unrealistic assessments and decisions [8,9].Significant attempts have been undertaken recently to include uncertainty in mine design.To incorporate ore grade and metal uncertainty into mine planning, a set of equally likely orebody models is used in combination with stochastic optimization for the mine production scheduling problem.These efforts include the development of stochastic integer programming (SIP) models [10][11][12], intending to maximize NPV and minimize estimated deviations from production goals, as well as streamlining the use of a stochastic stockpile.Presently, the key limitation of the SIP approach is that it is computationally very demanding-an issue that is currently being addressed through research in parallelization.This is related to the work of Godoy and Dimitrakopoulos [13] and Leite and Dimitrakopoulos [14], in which scheduling is stochastically optimized using a simulated annealing algorithm.Although the optimization part of this approach is computationally very efficient, the approach requires a labor-intensive preparation of inputs, as discussed in detail by Albor and Dimitrakopoulos [12].
The mine design formulation under the stochastic scenario presented so far only takes into account geological uncertainty, assuming that the price of metal remains constant during the mine's life.As a result, mine designs based on a fixed pricing value may fail to meet the desired goal.The typical response to this situation is really to create a new modified long-term plan and invest more funds, resulting in inefficient capital allocation, low investor returns, and a lower net present value.As a result, it is determined that including price uncertainty into production planning is required and vital, resulting in optimized production scheduling throughout the mine's life.
To address the practical issues related to incorporating the metal price uncertainty, Meagher et al. [15] proposed a new minimal cut algorithm.Goldberg and Tarjan [16] proposed an algorithm for optimizing final pit limits that combines multiple simulated orebodies and metal price uncertainties, which generate sets of possible economic values for mining blocks.The advantage of this algorithm is that it only needs to generate a single directed graph for the minimum cut algorithm, irrespective of the number of the simulated orebodies or the sets of possible economic values of the mining blocks involved.It was shown in their paper that in a two-dimensional example, the approach can jointly handle grade/metal and market uncertainties, while other uncertainties relevant to the Systems 2022, 10, 159 3 of 14 mining block's economic value can also be integrated.To incorporate multiple sources of uncertainties and pushback generation, Asad and Dimitrakopoulos [17] proposed a Lagrangian relaxation-based approach with a minimum cut algorithm.In their approach, they incorporated price uncertainty by price simulation using the mean reverting price model [18].Since the mean reverting method follows the Markov properties, the major drawback of the approach is that for calculating the conditional probability of future price value, the algorithm only considers the present value, neglecting the past value.The price value, like other time series methods, is dependent on a previous set of price values [19].
The push-relabel minimal cut algorithm is used in this paper to generate the ultimate pit and pushback while taking into account geological uncertainty.The pushbacks are generated by parameterization of the arc capacities of a minimum cut graph.A case study of an Indian iron ore deposit demonstrates the method's complexities.

Literature Survey
Metals and geological uncertainties are addressed based on the number of equal probable scenario mining maps of grade, metal content, and geology using a spatial stochastic modeling approach.These techniques explicitly take into account the volumetric differences between the existing data and the assumed mining blocks and the spatial correlation among measurements from the drilling process [20][21][22].Dimitrakopoulos and Ramazan [23] introduced the concept of ore body risk discounting, which postpones the risk for later periods by delaying the extraction period of the high-risk blocks.The main drawback of the above assumptions are they previously assigned risk probabilities for each block, discarding the grade over simulations and ignoring the uncertainty that has been jointly evaluated as groups of blocks in mining period to period [23,24].
The first mathematical models that effectively illustrated the stochastic features of OPMPS were proposed by Ramazan and Dimitrakopoulos [25].They take into account a two-stage probabilistic model that employs multiple similar scenarios to represent uncertain geology.The decision-making process is divided into two stages: the first stage concerns the mining sequence, and the second step of each scenario addresses the production target deviation.The objective is to maximize the predicted net present value and reduce production target deviations, which decreases the chance of not fulfilling production goals.For long-term production planning in a stochastic framework, the simulated annealing concept was first presented by Godoy and Dimitrakopoulos [13] and further investigated by Leite and Dimitrakopoulos [14] and AlborConsuegra and Dimitrakopoulos [12].This methodology was expanded to incorporate multiple mines, stockpiles, and processing destinations [26].To maximize net present value while reducing the variation from production goals, Ramazan and Dimitrakopoulos [11] designed a stochastic integer programming model.Chatterjee and Dimitrakopoulos [27] presented a case study based on the stochastic formulation, under uncertainty.Due to this, incorporating geological uncertainty into the optimization process is more complicated but also more beneficial.These benefits were initially emphasized [11][12][13][14][28][29][30][31].

Stochastic Mine Planning under Uncertainties
Uncertainty can be measured, modeled, and quantified by geostatistical methods.This can be achieved by generating the multiple equal probability scenario of orebody realization using any geostatistical simulation method.By including uncertainty in the decision-making approach, risk can be reduced.This will allow the mine design and production planning team to generate higher profit value and a better risk management strategy.However, the mining industry is aware of the risks and uncertainties.Because stochastic models are more relevant and realistic in how they handle uncertainty and risk, they attract more attention.A common strategy for solving stochastic challenges is to evaluate the profit under various scenarios and aim to minimize production target deviations [14,[28][29][30][31].
Geological uncertainty is considered to be the main reason for the unrealized forecast of cash flow by mining companies as it has a direct impact on the supply of ore and metals.The question arises of how to define the availability of ore for processing.It is not an easy question to answer, because it not only depends on the spatial distribution of the ore but also depends on the extraction sequence over time.Different extraction sequences will produce different ore supplies from the same deposit.The definition of ore is variable over time because it is a function of economic parameters and is time-dependent.The definition of available ore supply has traditionally been evaluated in the context of assumptions of fixed technical and economic constraints in time and space.In the conventional mine design and production scheduling framework, the average type of ore grade model is used in conjunction with geotechnical, economic, and environmental constraints to define the extraction order that produces the greatest economic return.The use of risk-free scenarios has shown a significant difference between forecasting economic returns and real economic returns.Some researchers have studied the impact of geologic uncertainty on the economics of a project and used conditional simulations to conduct a risk analysis of the performance parameters of the mine plan [8,9,27,28].They agreed that the use of smooth images or representations of geological reality (i.e., average type ore grade mode) when applying nonlinear transfer functions such as mining plans leads to unpredictable underestimation or overestimation of the performance parameters under consideration.To minimize the misclassification of resources, conditional simulation techniques are often used to estimate geological information at non-sampled locations based on stochastic model estimation techniques.This is done by interpolating data from a small number of exploration samples [27].The stochastic framework takes multiple simulations or similarly possible grade tonnage curves by taking geological uncertainty into account, and while taking into account a single (or constant and known) grade tonnage curve, conventional techniques neglect geological uncertainty [26,29,30].The creation of innovative scheduling systems that include simulated geological uncertainties in mine planning processes is facilitated by the existence of uncertainty modeling methods.This concept will be extended to a more general overview of the mining business when the global optimization of the mining complex is integrated.For this optimization, important factors such as hybrid demand, different processing flows, transportation systems, and product sales are included in the optimization [26].A thorough review of long-term mine production planning with uncertainty is provided by Koushav et al. [30].The stochastic graph closure problem of the original pushback design is a relaxed problem where the resource constraints are not included.To generate pushbacks, resource constraints have to be imposed.There are several approaches described in the literature.One strategy involves relaxing the resource constraints, putting them into the target function as a Lagrangian relaxation problem, and solving iteratively until a feasible solution is found [17].Although the Lagrangian relaxation indeed provides a reasonable upper bound, finding a feasible solution is not assured.The convergence of the algorithm is not also guaranteed, and the algorithm is also sensitive to many parameters [32].An alternative approach is relaxing the constraints from integer to real, solving the linear programming problem (LP-relaxation), and rounding the solution to an integer one.The limitation of the LP-relaxation algorithm is that the constraints may be violated during the rounding operation; a feasible solution is not guaranteed.It takes some postprocessing steps to reach a realistic solution.The third and most popular approach in the mining application is parameterized algorithm where the economic block value is parameterized with a multiplier λ.The parameterization algorithm of the L-G algorithm is well documented in various works of literature [2,3,33].By scaling the economic values for all blocks using a multiplier parameter λ, a series of "nested" pits can be generated.The structure of the article is as follows: The stochastic pushback design is presented mathematically.An undiscounted graph closure formulation with geological scenarios is discussed.A discounted graph closure formulation using the minimum cut algorithm and solution approach for push-back designs are presented.A case study of an Indian iron deposit is presented, followed by conclusions.

Stochastic Pushback Design Formulation
A three-dimensional block model is used to represent a mineral deposit.Each block has a set of data that describes its location in the three-dimensional model, the amount of ore or waste present, the grade, and so on.The goal of pushback design is to define the open pit mine's contour while keeping resource limits in mind.In other words, the goal is to locate a group of blocks that can be recovered from the deposit while keeping all limits in mind to maximize profit.The problem can be defined as an ultimate pit limit problem if resource limitations are removed from the pushback design formulation.The mine pushback design challenge is presented mathematically as follows: subject to where Γ i is the set of blocks which have to be removed before removing block i, c i,s is the amount of revenue generated or loss incurred by removing block i from scenario s, a i,s k is the amount of kth resource generated from block i from scenario s, b k is maximum allowable kth resource constraints, and K is the total number of resource constraints.
Objective (1) is to maximize the total profit from the mine over all scenarios.Equation ( 2) is a set of constraints for pushback design that must be satisfied.Equation ( 3) is a set of slope constraints such that for extracting block i from a deposit, all blocks in a set Γ i have to be extracted.Equation (4) ensures that a block is either extracted or remains in the deposit.
The c i,s of the individual blocks is calculated as follows: where MC is mining cost, PC is processing cost, Ti is block i tonnage, gi is block i grade, REC is recovery, ps is price simulation scenario, and SC is selling cost.It is clear from the above formulation of both single scenarios and multiple scenarios that pushback design is an NP-hard problem.However, if the constraints of Equations ( 2) and ( 3) are eliminated from these problems, they are efficiently solvable by any graph closure algorithm.In this paper, a push-relabel minimum cut approach is used for solving the graph closure problem and then a parameterized approach is applied to impose resource constraints for pushback design.

Stochastic Graph Closure with Multiple Scenarios
A minimum graph cut issue can be used to model pushback design without resource limitations.Using Equation (5), the block economic value (c i,s ) is determined.Blocks with positive values (c i,s > 0) are ore blocks, while those with negative or zero values (c i,s ≤ 0) are waste blocks.The directed graph is considered to have nodes at each block.Two special nodes are present: a source (s) and a sink (t).All nodes identified as ore blocks are linked to the source node, and the economic values of those blocks are represented by the capacity of those arcs (c i,s ).However, all nodes in the model that are labeled as waste blocks are linked to the sink node, and the capacities are equal to the absolute value of the waste blocks' economic values (|c i,s |).
A group of underlying blocks must first be extracted to mine a specific block, i.The slope constraints stated in Equation ( 3) can be maintained to achieve this.While removing the target block, it is necessary to determine the underlying blocks that must be removed to meet slope constraints.As a result, arcs are drawn from block i to the overlying blocks that have capacities of infinity (∝).Because slope constraint arcs have an infinite value, they will never be included in the minimum cut.
After forming the directed graph using multiple models as shown in Figure 1, the minimum cut algorithm can be implemented to solve pushback design without resource constraints.If the minimum cut algorithm is implemented with the above-described directed graph with S number of scenarios, It is not guaranteed that the same block will fall on the same side of the minimum cut in all scenarios.However, the decision should be binary (1, 0); i.e., a block is either within the closure or outside the closure, regardless of the number of scenarios.In the graph, an additional constraint is required to ensure that a provided block falls on the identical side of the cut in all scenarios.This constraint can be implemented by placing bidirectional arcs with infinite capacity (∝) among the same block i from different scenarios.Since these bidirectional arcs have infinite capacity, they will never be in the minimum cut, which ensures that there will be no such situation when the same block from different scenarios will fall on different sides of the minimum cut.Therefore, the pit generated by the minimum cut algorithm with multiple scenarios of orebody models is a valid pit where all constraints are respected except resource constraints.These nodes can be combined into a single node because the identical block from the many scenarios will fall on the same side of the minimal cut.It is possible to combine the arcs from the source node to the merging nodes from various scenarios into a single arc, and the capacity of the resulting arc will be ∑ S s=1 c i,s if c i,s > 0. Similarly, a single arc can be formed from a merge to the sink node(t), and the capacity of the arc will be ∑ S s=1 c i,s if c i,s ≤ 0. This concept is presented in Figure 2, where the graph is developed by merging the three scenarios of orebody models from Figure 1.In this figure, each block is assigned two values.The top-left corner value is the sum of the economic values of a block over all scenarios if c i,s > 0. The bottom-right corner value is the sum of the absolute block economic value over all scenarios if c i,s ≤ 0. The value in the top-right corner represents the arc capacity from the source to that node, and the value in the bottom-left corner represents the arc capacity from that node to the sink.For example, the top-left block in Figure 2 has an arc capacity from the source node to the node of 3, which is a sum of 2 and 1-the block To generate pushbacks, an iterative parameterization algorithm is applied until resource constraints are respected.The authors have applied a parameterized approach and proposed an approach to update the λ value in iteration.
After applying the parameterization of the block economic values, the problem can be rewritten as where A parameterization factor has a bound within the λ multiplier (0, 1).When λ is equal to 1, the problem is solved using the objective function Equations ( 6) and ( 9) without the need for pushback constraints.The goal of the parameterization strategy is to choose an appropriate λ value that will enable the problem posed in Section 3 to have a large upper-bound solution.The minimum cut algorithm covered in Section 3.1.1can be used to solve the objective function of Equation ( 6) for a particular value of λ.
Pushbacks are significantly influenced by the parameter λ, which is a multiplier of the directed arc capacity.The random selection of λ might lead to pushbacks with wide gaps, which would be unfavorable for mining.Any minimum cut algorithm can be used to maximize the optimization problem for a specific λ value once the merged graph has been constructed.The achieved solution, however, will ignore pushbacks' resource limitations.An iterative technique was used to choose the λ value to satisfy two resource limitations.The process begins with λ = 0 and progresses through the update of the λ value by λ + ∇λ as far as resource constraints are satisfied for all scenarios S.

3.
Assign best production k s = 0 at λ = 0 for all resource limitations k and scenario s.
Thus, both for resolution and computational efficiency, selecting the step size ∇λ is crucial.If ∇λ has been chosen too small, the algorithm will take a long time to converge.On the other hand, if the ∇λ value has been chosen too high, a large gap may be observed between the two pushbacks.Therefore, instead of choosing a fixed ∇λ value, a ∇λ value is updated in each iteration depending on how far away from the convergence point the solution is.By choosing a monotonically decreasing function for variable ∇λ, the large step size can be selected when the algorithm is far from the convergence point, and the small step size can be selected when the algorithm is approaching the convergence point.Equation ( 10) is used for calculating the value of ∇λ.
where w i,s are weights associated with k th resource constraints for grade scenario s and ε is a small constant number.Since the values of b k are constant and the values of best production k s are increased after each iteration, Equation (10) ensures that ∇λ value will have a monotonically decreasing function.When the algorithm of the pseudo-code is close to violating the resource constraints, the ∇λ value will be incremented by a small magnitude.

Results and Discussion
Composited iron sample analyses were statistically examined.Table 1 displays the statistical analysis for the composited data.The iron grade is skewed towards the right, although with relatively less variance, as the table demonstrates.The frequency distribution of the dataset was also prepared using the histogram plot (Figure 3).It clearly shows the skewness property of the dataset.The grade was estimated using ordinary kriging [34].The variogram analysis was used to determine the data's spatial correlation [35].For testing the anisotropy, both directional and omnidirectional variograms were examined.The variogram analysis of the deposit was exceptionally complicated in a proposed mine study due to the presence of many different lithotypes.Seven different lithotypes are present in this mine.The modeling of the deposit without classifying lithology may lead to a large error in estimation.Some of the lithotypes show very less variance compared to the others as well as compared to the original data.A better variogram model may be obtained for some of the lithotypes by modeling each lithology separately.Both approaches were used to estimate the deposit.First, without considering lithology, an omnidirectional variogram for iron was constructed.A spherical variogram model was fitted for modeling spatial continuity.The range, sill, and nugget of the variable are 255.28, 31.26, and 38.08, respectively.After fitting the variograms, the search for anisotropy was carried out.The variograms were calculated in eight different directions with 22.5 • tolerance.The directional variogram demonstrated that the deposit has directional anisotropy with a major axis along the N-S direction.The major to semi-major and major to minor ratios are 3:1 and 5:1, respectively.In the second approach, an attempt was made to develop a separate variogram for each lithotype.Unfortunately, it was not possible to obtain an experimental variogram due to a smaller number of available data for individual lithology.Therefore, it was decided to perform resource evaluation with composited data without considering lithology in this study mine.
The block model of the deposit was estimated using the ordinary kriging technique.Because of the size of the selective mining unit (SMU), the estimated block size is 20 m × 20 m × 10 m [36].A total number of 14 benches is considered for this deposit.Figure 4 depicts the iron ore grade map of one of the fourteen benches.Following that, the estimated block model was used as input for mine production scheduling.
The proposed models were implemented in a MATLAB platform and solved with the CPLEX solver [37].Table 2 presents the specific economic parameters that were utilized in this study.The final result is presented in Tables 3 and 4 for the traditional model and Tables 5 and 6 for the proposed model.Figure 5 depicts the production schedule of the proposed model.A period-wise production schedule is shown in Tables 3 and 5. Ore production (million tons), metal production (kilotons), and NPV (USD 1 million) generated by the traditional and proposed model are calculated in Tables 4 and 6.It was found that when these two models were compared, they produced different amounts of ore, waste, and metal content over the length of the mine's life.However, the proposed model generates comparatively more NPV than the traditional model because of the different extraction sequence.Due to the production of more ore and metal in a shorter amount of time and the consequent generation of greater discounted cash flow as compared to the traditional model, the suggested model's NPV is determined to be higher.The traditional approach requires less calculation time than the suggested model; due to the increase in the number of constraints and decision factors, this is quite obvious.Additionally, the proposed model is more reliable for efficient production scheduling over the mine's life.The mining sector is looking for more advanced solutions to address the challenges of large-scale production scheduling in many unpredictably occurring scenarios.Mine planners typically employ optimization strategies that produce the highest undiscounted profits.Different companies of mining software are becoming aware of the drawbacks of conventional mine planning methods and the necessity of including uncertainty in the mine plan.The findings of this case study demonstrate how incorporating uncertain parameters enables mine planners to better understand the risk associated with not attaining targets and the efficacy of their cash flow analysis.Mine planners have excellent potential to maximize returns on their investment with a high degree of confidence and a shorter mine life by using these optimization techniques.

Conclusions and Future Scope
In this paper, an uncertainty-based ultimate pit limit and pushback design algorithm was developed using a graph cut algorithm to integrate geological-based models.The suggested approach demonstrates that it is comparatively easy to create a minimum cut network flow model that incorporates uncertainty, and it can deal with uncertainty in the block economic value of the mining blocks being sequenced.Different scenarios from various production periods were merged with the stochastic framework as opposed to traditional pit optimization, which only considered one constant grade value for the life of the mine.The smoothing spline model was used to estimate trends, and a sequential Gaussian simulation approach was used to simulate the residual values.The directed graph was constructed using a simulated grade, and pushbacks were generated by parameterizing the directed graph.The key benefit of the suggested technique is that it is computationally very quick, making it possible for large-scale deposits to integrate grade uncertainties in the optimization process.The proposed method was applied to an Indian iron ore deposit.Results show that a total number of nine nearly equal size pushbacks were generated.When compared with a traditional model, it was observed that the proposed stochastic optimization algorithm generates 5% more NPV.It is anticipated that the applicability and clarity of the solution structure given in this research would make it simpler for the mining industry to incorporate.However, there are still vital areas for future study that need to be investigated further.First, while grade-related uncertainty is taken into account in our models, the actual range of uncertainty is far wider and displays itself in a variety of ways, most significantly in the economical or operational constraints.Therefore, it is expected that including market uncertainty inside the suggested framework will be beneficial-given that doing so can be done without falling into "computational complexity".

Figure 1 .
Figure 1.Graph constructed from multiple scenarios of models.

Figure 2 .
Figure 2. Graph constructions after merging three block economic models generated using multiple scenarios.

Figure 4 .
Figure 4. Map showing the ore grade for one of the deposit's fourteen benches.

Figure 5 .
Figure 5. Production schedule of the proposed model year wise.

Table 2 .
Specific details on the values of certain parameters.

Table 3 .
Traditional model block extraction year-wise.

Table 4 .
Traditional model of calculated ore, waste, metal, and NPV per year.

Table 5 .
Proposed model block extraction year-wise.

Table 6 .
Proposed model of average calculated ore, waste, metal, and NPV per year.