Simulation of a Mining Value Chain with a Synthetic Ore Body Model : Iron Ore Example

Reconciliation of geological, mining and mineral processing information is a costly and time demanding procedure with high uncertainty due to incomplete information, especially during the early stages of a project, i.e., pre-feasibility, feasibility studies. Lack of information at those project stages can be overcome by applying synthetic data for investigating different scenarios. Generation of the synthetic data requires some minimum sparse knowledge already available from other parts of the mining value chain, i.e., geology, mining, mineral processing. The aim of the paper is to describe how to establish and construct a synthetic testing environment, or “synthetic ore body model” for data integration by using a synthetic deposit, mine production, constrained by a mine plan, and a simulated beneficiation process. The approach uses quantitative mineralogical data and liberation information for process simulation. The results of geological and process data integration are compared with the real case data of an apatite iron ore. The discussed approach allows for studying the implications in downstream processes caused by changes in upstream parts of the mining value chain. It also opens the possibility of optimising sampling campaigns by investigating different synthetic drilling scenarios including changes to the spacing between synthetic drill holes, composite length, drill hole orientation and assayed parameters. A synthetic deposit model can be a suitable tool for testing different scenarios for implementation of geometallurgical programs and also an educational tool for universities and companies.


Introduction
Geometallurgy aims to create a predictive model for mine-to-metal production chain through the combination of geological, processing and economic models.The highest returns on application of geometallurgy are expected to come from complex deposits with high variability in properties affecting its processing [1][2][3].The geometallurgical predictive model is to be established during development of a geometallurgical program [1].Establishment of the model requires involvement of various professionals from different disciplines and the availability of harmonised modelling and simulation tools.The high-fidelity orebody model [4] with a high resolution and utility of the information and possible moves from indicators to predictors is a spatial basis for the geometallurgical program.In a high-fidelity geometallurgical program, the critical methodology selection, such as sampling design, components included in the block model and their assaying methods, and geostatistical solutions should be based on scientific principles.The benefits and limitations of program design vary from case to case and thus no single deterministic approach can be developed.In real cases comprehensive study of different alternatives is difficult.The use of a synthetic orebody model combined with process Traditionally, the modelling of an ore body has been restricted to the geological domain (the definition of physical regions with homogeneous properties based on lithology, mineral grade and style of mineralisation [19,20], densities, and elemental grades).The processing properties have been almost totally neglected.However, a simulated synthetic ore body as defined in this paper must satisfy three conditions: Include processing properties; show realistic variability within synthesised data and impose spatial cohesion between data.Additional parameters, such as constraints by a mining method, processing performance, and economic response will produce more realistic output.
The selection of a modelling approach (Table 1) to simulate geological and propagate mineral processing properties in a three-dimensional physical space involves a trade-off between geological realism and conditioning capabilities.The number of data points, or number of samples, used in modelling is crucial when selecting the modelling method.The methods which are efficient with small number of data points can be used in geometallurgy, contrary to traditional resource estimate methods that work with large dataset.Prior knowledge about the studied phenomena and whether the modelling method is parametric or non-parametric may also have impact on number of data points [21].Another parameter is a smoothing, which can be defined as whether the model smooths predictions at sampling locations or not [22].Estimation methods such as kriging have a significant smoothing effect, while simulation techniques allow to produce values that are true to the fluctuations of the phenomenon [23].Simulation techniques provide higher realism than traditional geostatistics [24,25].A process-mimicking approach was selected for modelling in this study due to its relative simplicity, and that it does not require large number of data points has low smoothing and high realism.The algorithm was implemented in a MATLAB code and currently is not publicly available.Process (e.g., process based, process-mimicking) [24,31] Small Low High

Synthetic Ore Body Model
The synthetic ore body model includes two modules: The synthetic deposit module and the synthetic sampling module (Figure 1).Each module is comprised of a spatial component, which considers the location of each point in a physical space, and a database, which carries non-spatial quantitative and qualitative information for each point.Spatial and non-spatial information is connected by a unique ID number assigned to each point of the model in physical space and a corresponding record in the database.Both databases (deposit module and sampling module) have the same meta data and carry the same type of information (elemental composition, mineralogy, recoveries, throughput, mining cost, value etc.) about each block of the synthetic deposit and each segment of a synthetic drill core.

Synthetic Deposit Module
The Synthetic Deposit Model is described with three-dimensional voxel model and database.The Synthetic Deposit Model focus on outlining the deposit's borders and assigning geological, mining, processing and economic properties to each voxel of the voxel model.In practice, the uncertainty of the ore response in the process may also be impacted by extraction sequence, processing available, operators' team working at the plant, and batch with which the extracted block was blended during the process.However, those factors would have a lower impact on the process performance in iron ore mining than in, for example, copper or gold ores.

Spatial Data
The spatial part of the deposit module for the synthetic ore body is a three-dimensional voxel model, where each voxel corresponds to a minable block of the deposit.The size of the voxels can be set to any constant value.The location of each voxel is described by the coordinates x, y, z.The x and y coordinates are planar, and z represents height.The coordinate system follows the left-hand rule.Therefore, the centre of the coordinates lies in the lower left corner of the voxel model and values of the coordinate z increase upwards.

Database
A database of voxel properties describes geological, mineralogical, mining, mineral processing and economic properties of each voxel (Figure 2).Each voxel gets one entry in a database with complete information about its properties.The geological data in the database are derived from the three-dimensional voxel model.
The voxel model represents the extent of the mineral deposit, including both ore body and country rock.The geological domains, mineralogical and elemental properties are the key geological features chosen for modelling synthetic deposits.Generated mineralogical and chemical information is inputted to the mineral processing model implemented in HSC simulation software [38] for each block.Alternatively, other software may be used: CEET, FLEET [39]; Cycad Process [40]; ECS/CEMulator [41]; IDEAS Integrated [42]; JKTech's JKSimBlast, JKSimMet and JKSimFloat [43]; Microsoft Excel based simulator [44]; MinOOcad [45], MODSIM [46]; SimSci DYNSIM [47].Process parameters estimated for each block are constrained with a mining method, since sequencing of the Structure of the synthetic ore body model.

Synthetic Deposit Module
The Synthetic Deposit Model is described with three-dimensional voxel model and database.The Synthetic Deposit Model focus on outlining the deposit's borders and assigning geological, mining, processing and economic properties to each voxel of the voxel model.In practice, the uncertainty of the ore response in the process may also be impacted by extraction sequence, processing available, operators' team working at the plant, and batch with which the extracted block was blended during the process.However, those factors would have a lower impact on the process performance in iron ore mining than in, for example, copper or gold ores.

Spatial Data
The spatial part of the deposit module for the synthetic ore body is a three-dimensional voxel model, where each voxel corresponds to a minable block of the deposit.The size of the voxels can be set to any constant value.The location of each voxel is described by the coordinates x, y, z.The x and y coordinates are planar, and z represents height.The coordinate system follows the left-hand rule.Therefore, the centre of the coordinates lies in the lower left corner of the voxel model and values of the coordinate z increase upwards.

Database
A database of voxel properties describes geological, mineralogical, mining, mineral processing and economic properties of each voxel (Figure 2).Each voxel gets one entry in a database with complete information about its properties.The geological data in the database are derived from the three-dimensional voxel model.
The voxel model represents the extent of the mineral deposit, including both ore body and country rock.The geological domains, mineralogical and elemental properties are the key geological features chosen for modelling synthetic deposits.Generated mineralogical and chemical information is inputted to the mineral processing model implemented in HSC simulation software [38] for each block.Alternatively, other software may be used: CEET, FLEET [39]; Cycad Process [40]; ECS/CEMulator [41]; IDEAS Integrated [42]; JKTech's JKSimBlast, JKSimMet and JKSimFloat [43]; Microsoft Excel based simulator [44]; MinOOcad [45], MODSIM [46]; SimSci DYNSIM [47].Process parameters estimated for Minerals 2018, 8, 536 5 of 23 each block are constrained with a mining method, since sequencing of the blocks sent to the process depends on the production plan of the selected mining method.Process performance, e.g., concentrate quality, throughput, mining cost, dilution etc. provide inputs for the economical assessment of the designed mine.Geological, production and economic information about each block is stored in a database and linked to the block with a unique identification number.blocks sent to the process depends on the production plan of the selected mining method.Process performance, e.g., concentrate quality, throughput, mining cost, dilution etc. provide inputs for the economical assessment of the designed mine.Geological, production and economic information about each block is stored in a database and linked to the block with a unique identification number.

Domains: Physical space
The general term "domain" here refers to a volume or physical space with homogeneous properties (e.g., mineral distribution, rock density, grindability, commodity material recovery, or texture).In other words, voxels belonging to the same ore type, adjacent to each other and showing spatial continuity of any properties (e.g., geological, metallurgical properties) inside the ore body are referred to as domains [48,49].Ore types, alteration zones, mineral grade distributions, and weathering zones are called geological domains.Domains used to model process parameters are referred to as geometallurgical domains.
The complex shape of the domains is modelled (approximated) with simple geometrical bodies.Traditional geological modelling software, e.g., Maptek [50] and Dassault Systems [51] deal with 3D domains by using solids generated from wireframes or through implicit modelling.However, an ellipsoid (Figure 3 and Equation ( 1)) is chosen here, since it is a common shape used in geosciences, i.e., geostatistics.Other geometric shapes may be used in modelling domains as well (e.g., cubes, spheres and their sectors).It might be enough to use one ellipsoid to approximate each domain: Here (see also Figure 3) o = (x0, y0, z0) is the centre of the ellipsoid; a, b, c are its semi-axes.The spatial orientation can also be given by Euler rotation angles α, β, γ.

Domains: Physical space
The general term "domain" here refers to a volume or physical space with homogeneous properties (e.g., mineral distribution, rock density, grindability, commodity material recovery, or texture).In other words, voxels belonging to the same ore type, adjacent to each other and showing spatial continuity of any properties (e.g., geological, metallurgical properties) inside the ore body are referred to as domains [48,49].Ore types, alteration zones, mineral grade distributions, and weathering zones are called geological domains.Domains used to model process parameters are referred to as geometallurgical domains.
The complex shape of the domains is modelled (approximated) with simple geometrical bodies.Traditional geological modelling software, e.g., Maptek [50] and Dassault Systems [51] deal with 3D domains by using solids generated from wireframes or through implicit modelling.However, an ellipsoid (Figure 3 and Equation ( 1)) is chosen here, since it is a common shape used in geosciences, i.e., geostatistics.Other geometric shapes may be used in modelling domains as well (e.g., cubes, spheres and their sectors).It might be enough to use one ellipsoid to approximate each domain: Here (see also Figure 3) o = (x 0 , y 0 , z 0 ) is the centre of the ellipsoid; a, b, c are its semi-axes.The spatial orientation can also be given by Euler rotation angles α, β, γ. blocks sent to the process depends on the production plan of the selected mining method.Process performance, e.g., concentrate quality, throughput, mining cost, dilution etc. provide inputs for the economical assessment of the designed mine.Geological, production and economic information about each block is stored in a database and linked to the block with a unique identification number.

Domains: Physical space
The general term "domain" here refers to a volume or physical space with homogeneous properties (e.g., mineral distribution, rock density, grindability, commodity material recovery, or texture).In other words, voxels belonging to the same ore type, adjacent to each other and showing spatial continuity of any properties (e.g., geological, metallurgical properties) inside the ore body are referred to as domains [48,49].Ore types, alteration zones, mineral grade distributions, and weathering zones are called geological domains.Domains used to model process parameters are referred to as geometallurgical domains.
The complex shape of the domains is modelled (approximated) with simple geometrical bodies.Traditional geological modelling software, e.g., Maptek [50] and Dassault Systems [51] deal with 3D domains by using solids generated from wireframes or through implicit modelling.However, an ellipsoid (Figure 3 and Equation ( 1)) is chosen here, since it is a common shape used in geosciences, i.e., geostatistics.Other geometric shapes may be used in modelling domains as well (e.g., cubes, spheres and their sectors).It might be enough to use one ellipsoid to approximate each domain: Here (see also Figure 3) o = (x0, y0, z0) is the centre of the ellipsoid; a, b, c are its semi-axes.The spatial orientation can also be given by Euler rotation angles α, β, γ.However, cases that are more complex may require usage of combination of multiple ellipsoids.Interactions or conflicting overlaps between different domains are resolved by applying Boolean operations (Figure 4).Rock types, the same as other spatially continuous variables, are defined by modelling domains and encoding them with a unique identifier for the particular rock type.One voxel of the synthetic ore body model is occupied only by one rock type (ore type, ore class, weathering type etc.).

Mineralogy
Modal mineralogy is defined for the mineral grade model by considering mineralogical difference between geological domains and preserving relevant quantities, distribution, and composition of the minerals.Firstly, a complete list of minerals for the entire model is identified.Then an appropriate sub-list of minerals is defined for each geological domain.It is assumed that a mineral can occur in the orebody due to any one or more of several mineralization effects.Therefore, overlapping generations of the same mineral should be modelled as separate objects.Mineral grades from different populations of the same mineral can be summed up and given a new name (e.g., Magnetite_Sum) or kept separately in the database (e.g., Magnetite_1, Magnetite_2, Magnetite_3).
Mineral distribution for each voxel is modelled as a function of coordinates of the voxel in the block model (Figure 5).When extension of the mineral inside the block model is known, the physical space of the actual mineral distribution can be restricted by geological domains (ellipsoids) and scaled correspondingly.Mineral quantity distribution can also be complex since they represent several overlapping distributions with a drift.Initial hypothesis about mineral distribution is expressed as a generalised additive model (GAM) shown below by Equation (2).
= + + → = ( , , ), min = , max = = ( , , ), min = , max = = ( , , ), min = , max = where, is a mineral grade in a given voxel scaled to the interval of minimum and maximum values defined for the mineral; T is a trend function; S is a systematic error; N is a random error; T, S, N are functions of coordinates (x, y, z) of the voxel.Each ( , , ), ( , , ), ( , , ) can be chosen as follows, Equation (3), where each choice for the function defines a different hypothesis for mineral distribution and mineralisation.
where r is defined as the distance (e.g., Euclidian or Manhattan) in 1-, 2-or 3-dimensional space between the voxel of interest and the centre of the geological domain (ellipsoid): where n is any real number; x, y, z are coordinates of the voxels locations; and , , are cooridnates of the center of the geological domain (as shown in Figure 3).The physical location of each domain (ellipsoid) in relation to other domains is defined by aligning ellipsoid centres along an auxiliary design line or placed individually by the model user.The auxiliary design line is used as reference for the ellipsoids' locations and does not have any particular geological meaning.
Rock types, the same as other spatially continuous variables, are defined by modelling domains and encoding them with a unique identifier for the particular rock type.One voxel of the synthetic ore body model is occupied only by one rock type (ore type, ore class, weathering type etc.).

Mineralogy
Modal mineralogy is defined for the mineral grade model by considering mineralogical difference between geological domains and preserving relevant quantities, distribution, and composition of the minerals.Firstly, a complete list of minerals for the entire model is identified.Then an appropriate sub-list of minerals is defined for each geological domain.It is assumed that a mineral can occur in the orebody due to any one or more of several mineralization effects.Therefore, overlapping generations of the same mineral should be modelled as separate objects.Mineral grades from different populations of the same mineral can be summed up and given a new name (e.g., Magnetite_Sum) or kept separately in the database (e.g., Magnetite_1, Magnetite_2, Magnetite_3).
Mineral distribution for each voxel is modelled as a function of coordinates of the voxel in the block model (Figure 5).When extension of the mineral inside the block model is known, the physical space of the actual mineral distribution can be restricted by geological domains (ellipsoids) and scaled correspondingly.Mineral quantity distribution can also be complex since they represent several overlapping distributions with a drift.Initial hypothesis about mineral distribution is expressed as a generalised additive model (GAM) shown below by Equation (2).
where, M is a mineral grade in a given voxel scaled to the interval of minimum and maximum values defined for the mineral; T is a trend function; S is a systematic error; N is a random error; T, S, N are functions of coordinates (x, y, z) of the voxel.Each f T (x, y, z), f S (x, y, z), f N (x, y, z) can be chosen as follows, Equation (3), where each choice for the function f defines a different hypothesis for mineral distribution and mineralisation.
where r is defined as the distance (e.g., Euclidian or Manhattan) in 1-, 2-or 3-dimensional space between the voxel of interest and the centre of the geological domain (ellipsoid): where n is any real number; x, y, z are coordinates of the voxels locations; and x o , y o , z o are cooridnates of the center of the geological domain (as shown in Figure 3).Describing minerals independently from each other might lead to a severe over or under estimation of mineral proportions for some voxels.This problem is avoided by using a simplified two-step approach: Commodity minerals are modelled independently from each other summing up to a value A%, the rest (gangue) minerals, are modelled as a fraction of available space remaining, defined as 100% − A% = B%.However, this poses another problem of constructing the elemental and mineral grades since it introduces spurious correlations among data.The simplification was made here in light of the deposit type under the study with high grades of commodity mineral, and thus high iron grades, relatively high geological homogeneity.
Mineral chemical composition can vary from deposit to deposit or even within the same deposit.Thus, the mineral composition is described with the elemental grades being a function of a mineral grain's spatial location in the voxel model.Modelling of the minerals distribution is done by considering the following aspects: i. Stationarity is insured by modelling the "stationarity ellipsoids".A stationarity ellipsoid is a geometry where the mineral distribution is the same for a given mineral throughout the portion of the voxel model enclosed by this stationarity ellipsoid.The algorithm of describing stationarity ellipsoids is identical to the one used for describing geological domains.Stationarity allow for describing anisotropy of the mineral distributions within the ellipsoids' dimensions and orientations.ii.Spatial conflicts between geological domains and stationarity ellipsoids are resolved with Boolean operations (Figure 4).iii.Minerals modelled with Equations ( 2)-( 4) may have compositions which are outside the normal range of values.Therefore, a certain minimax range can be imposed by rescaling (or truncating) the mineral distribution.iv.The total sum of mineral grades should be closed to 100%.Normalisation of values to a constant sum is also called a closure and forces negative correlations [52][53][54][55][56].According to Aitchison [57] log-ratios should be used when it is needed to maintain the constant sum constraint.The Describing minerals independently from each other might lead to a severe over or under estimation of mineral proportions for some voxels.This problem is avoided by using a simplified two-step approach: Commodity minerals are modelled independently from each other summing up to a value A%, the rest (gangue) minerals, are modelled as a fraction of available space remaining, defined as 100% − A% = B%.However, this poses another problem of constructing the elemental and mineral grades since it introduces spurious correlations among data.The simplification was made here in light of the deposit type under the study with high grades of commodity mineral, and thus high iron grades, relatively high geological homogeneity.
Mineral chemical composition can vary from deposit to deposit or even within the same deposit.Thus, the mineral composition is described with the elemental grades being a function of a mineral grain's spatial location in the voxel model.Modelling of the minerals distribution is done by considering the following aspects: i Stationarity is insured by modelling the "stationarity ellipsoids".A stationarity ellipsoid is a geometry where the mineral distribution is the same for a given mineral throughout the portion of the voxel model enclosed by this stationarity ellipsoid.The algorithm of describing stationarity ellipsoids is identical to the one used for describing geological domains.Stationarity allow for describing anisotropy of the mineral distributions within the ellipsoids' dimensions and orientations.ii Spatial conflicts between geological domains and stationarity ellipsoids are resolved with Boolean operations (Figure 4).iii Minerals modelled with Equations ( 2)-( 4) may have compositions which are outside the normal range of values.Therefore, a certain minimax range can be imposed by rescaling (or truncating) the mineral distribution.
iv The total sum of mineral grades should be closed to 100%.Normalisation of values to a constant sum is also called a closure and forces negative correlations [52][53][54][55][56].According to Aitchison [57] log-ratios should be used when it is needed to maintain the constant sum constraint.The significance of closure problem for environmental data was emphasised by Filzmoser et al. [58] and Reimann et al. [59], and for compositional geochemical data by Makvandi et al. [60].

Chemical composition
A voxel's chemical composition depends on the modal mineralogy and chemical composition of the minerals.Each mineral within a geological domain has a defined chemical composition.If the chemical composition of a mineral varies, then several mineral species may have to be modelled, e.g., amphibole 1, amphibole 2, amphibole 3. The vector of chemical composition (v) of a voxel is calculated from modal mineralogy and the chemical composition of minerals, Equation ( 5): where A is a matrix of chemical composition of minerals (sometimes also called mineral matrix), and x is a vector of mass proportions of minerals in a voxel.

Density
The density of a voxel point is based on the densities of each of the separate minerals as a weighted harmonic mean assuming zero porosity of the rock, Equation ( 6): where, ρ i is the density of the mineral i, wt i , % is the weight fraction of the mineral i.The non-zero porosity can be considered by multiplying by the corresponding coefficient a ∈ (0, 1].Other geological (e.g., alterations, mineral textures), mineralogical (e.g., modal mineralogy, mineral liberation, grain size distribution), geophysical (e.g., magnetic susceptibility, porosity, dielectric permittivity, electrical conductivity) and geotechnical (e.g., Poisson's ratio) properties might be considered in the future.

Mining
A mining plan defines the ore sequence prior to coming to the processing plant.The mining plan is determined by the mining method.To enable the mining module, the voxel model is transformed into the resource block model by giving dimensions to voxels.For the simplicity, and for preserving the maximum resolution of the model, the conversion is made by assigning one voxel to one block.However, for low variability models, the properties of one block can be retrieved by averaging the properties of several voxels.
Mining constraints for the synthetic ore body model can be introduced by applying a mining method, which may be either surface mining or underground mining.However, the algorithm for implementation of underground mining remains to be developed.The mining constraints can also be omitted (thereby assuming zero ore dilution, and zero ore losses).The time constraint of the mining production can be modelled by applying, for example, a series of pushbacks in an open pit mining (Figure 6).Each following pushback is represented by a cone of a larger size with preserved slope angle (ratio between cone's diameter and height).Blocks enclosed in mining area of the synthetic deposit are extracted at once.Ore can be sequenced based on metal grade, mineral grade, processing cost, commodity recovery in the separation process, or by any other desirable property.For the mine scheduling and open pit optimization, there are several methods available such as floating cone [61], Lerchs-Grossman algorithm [62], network flow approaches [63], Dagdelen-Johnson Lagrangian Parametrization [64].Newman et al. [65] and Meagher et al. [66] provide a review of those and some more methods for the open pit optimization and mine planning.In this study, the optimization step was done by approximating open pit by a cone enveloping the whole ore body.The pushbacks are simply the sequential cones allocated with equal intervals from the bottom to the top of the ore body.

Processing
The process model is implemented through process simulation in the HSC Sim 7.1 process simulator [38].The simulation is done by considering the liberation and distribution of the mineral particles, and the process simulator is capable of handling liberation information and multiphase particles.The process model treats each voxel separately.The voxel information is retrieved from the synthetic ore body model and the time aspect is controlled by the mining model.The information on the plant feed, gathered through geological information from each voxel, enables the utilisation of particle-based models [67][68][69].As the material enters the plant, it is converted to particles by applying the liberation distribution of the corresponding geometallurgical domain.The modal composition of a block and the geometallurgical domain may be different; the mass proportions of particles in the particle population are adjusted using Equation ( 7) [70]: where ̂ is the iteratively adjusted mass proportion of the mineral grades; is the mass proportion of the mineral grades before adjustment; ( ) is the mass proportion of mineral in a particle; L is a total number of minerals; and the correction factor is calculated for each mineral i before each iteration round as follows, Equation ( 8 For the mine scheduling and open pit optimization, there are several methods available such as floating cone [61], Lerchs-Grossman algorithm [62], network flow approaches [63], Dagdelen-Johnson Lagrangian Parametrization [64].Newman et al. [65] and Meagher et al. [66] provide a review of those and some more methods for the open pit optimization and mine planning.In this study, the optimization step was done by approximating open pit by a cone enveloping the whole ore body.The pushbacks are simply the sequential cones allocated with equal intervals from the bottom to the top of the ore body.

Processing
The process model is implemented through process simulation in the HSC Sim 7.1 process simulator [38].The simulation is done by considering the liberation and distribution of the mineral particles, and the process simulator is capable of handling liberation information and multiphase particles.The process model treats each voxel separately.The voxel information is retrieved from the synthetic ore body model and the time aspect is controlled by the mining model.The information on the plant feed, gathered through geological information from each voxel, enables the utilisation of particle-based models [67][68][69].As the material enters the plant, it is converted to particles by applying the liberation distribution of the corresponding geometallurgical domain.The modal composition of a block and the geometallurgical domain may be different; the mass proportions of particles in the particle population are adjusted using Equation (7) [70]: where pj is the iteratively adjusted mass proportion of the mineral grades; p j is the mass proportion of the mineral grades before adjustment; χ(i) j is the mass proportion of mineral in a particle; L is a total number of minerals; and the correction factor κ i is calculated for each mineral i before each iteration round as follows, Equation (8): where M(i) is a mineral grade in the sample; n is number of particles; and p refers to the mass proportion of particle in a size class.As such, the denominator is the mineral grade back-calculated from the liberation data.
The unit process models contain a description of the behaviour of each particle based on one or several properties: Density, size, mineral composition, and shape [71].For each processed voxel or block the process model returns product quantities and qualities (elemental grades, mineral grades, particle size distribution (PSD), mineral liberation information) with processing information (time spent on processing the block, processing costs and consumables quantities required for processing the block.).
The process model uses the geological information modelled in the synthetic ore body model from the previous steps (see Section 2.3.2.1 Geology) as an input (in other words, feed stream).The geological information of each voxel is fed to the beneficiation simulation separately (no blending).The output of the beneficiation simulation includes the composition of all the process streams; recoveries, mass pulls of the separation processes; throughputs and energy consumption calculated for each voxel.Beneficiation information produced at this step can be added to the database and be linked to the spatial part of the model by a unique ID.
For the beneficiation simulation done at an elemental level, an input would include elemental composition of the voxel, i.e., Fe, Si, Al etc.At a mineralogical level, such input would include modal mineralogy of the voxel, i.e., commodity mineral (e.g., magnetite) and gangue minerals (e.g., biotite, quartz etc.).Beneficiation simulation by size is possible if mineral distributions in the synthetic ore body model are simulated by size.Alternatively, a fragmentation model may be applied over the synthetic data, c.f. [48,72].
In order to design a beneficiation simulation model, process performance information is needed for the similar ore types.Such information can include lab scale tests (e.g., WLIMS, Davis tube, flotation, Bond work index test), pilot plant tests, or plant survey.

Economics
The final performance of every mining project is always estimated economically.Conclusions can be drawn by comparing costs related to the implementation of the project and operating mining production, and revenues obtained from the sold final material.Capital and operating costs of production can be estimated from cost models (e.g., InfoMine [73] and Sayadi et al. [74]) and revenues-from the historical commodity prices taken at any date (e.g., Kitco Metalc Inc. [75] and LME [76]).While revenues are directly connected to the recoverable commodity grade, operating cost cannot be estimated so easily.Firstly, mining operating cost is linked to the ease of extraction of each block and the cost of transportation to the processing plant.For each mining method there is a mineral recovery and dilution; and mining layout also has an influence (Table 2).Depending on the mining method (Table 2), the cost of separate block extraction can have a connection to the cost of extraction of neighbouring blocks.Hardness and mass of the volume unit would be another parameter to control.Process operating costs will depend on the energy used in comminution and the chemicals and energy costs for separation processes.
The mining cost solution implemented in the synthetic ore body model accounts for the depth of the mining and density of the extracted material.The processing cost per tonne is considered to be a constant value.Present value (PV), future value (FV) and net present value (NPV) of the material (concentrate or metal) produced by the process are calculated as shown by Equation (9).The results from the economic model can be used for comparison between different scenarios, or for estimating the impact of managerial decisions.
where j is a period [1,N]; d is a discounting rate, %; P is the metal price; s is the sales cost; Q r is the material recovered in mining and process, units per period; c is the milling cost; Q c is the processed amount of material, units per period; m is the mining cost; Q m is the mined amount of material, units per period.

Sampling Module
Sampling is the main source of information on the deposit's geology, mineralogy and processing properties.Synthetic sampling is implemented by simulating synthetic drill cores.Mathematically they can be represented as infinitesimally thin cylinders, or as lines.
The sampling module (Figure 7) enables the investigation of different drilling patterns, orientations, and sampling densities.The spatial part of the sampling module includes x, y, z coordinates of each sample.For sampling with drill holes, coordinates x, y, z correspond to the location of the section of the drill core.The coordinate system used in a sampling module is the same as for the Deposit module described in Section 2.3.1.Spatial data.Depending on the mining method (Table 2), the cost of separate block extraction can have a connection to the cost of extraction of neighbouring blocks.Hardness and mass of the volume unit would be another parameter to control.Process operating costs will depend on the energy used in comminution and the chemicals and energy costs for separation processes.
The mining cost solution implemented in the synthetic ore body model accounts for the depth of the mining and density of the extracted material.The processing cost per tonne is considered to be a constant value.Present value (PV), future value (FV) and net present value (NPV) of the material (concentrate or metal) produced by the process are calculated as shown by Equation (9).The results from the economic model can be used for comparison between different scenarios, or for estimating the impact of managerial decisions.
where j is a period [1,N]; d is a discounting rate, %; P is the metal price; s is the sales cost; is the material recovered in mining and process, units per period; c is the milling cost; is the processed amount of material, units per period; m is the mining cost; is the mined amount of material, units per period.

Sampling Module
Sampling is the main source of information on the deposit's geology, mineralogy and processing properties.Synthetic sampling is implemented by simulating synthetic drill cores.Mathematically they can be represented as infinitesimally thin cylinders, or as lines.
The sampling module (Figure 7) enables the investigation of different drilling patterns, orientations, and sampling densities.The spatial part of the sampling module includes x, y, z coordinates of each sample.For sampling with drill holes, coordinates x, y, z correspond to the location of the section of the drill core.The coordinate system used in a sampling module is the same as for the Deposit module described in Section 2.3.1.Spatial data.The difference between a drill core, and reverse circulation (RC) chips may be considered by selecting the sample sizes (half core, quarter core) and sample preparation methods.Sample size in the drill holes is defined by the sample or composite length.A smaller composite length allows better capturing of the variability of the studied parameter, e.g., magnetite grade, Bond Work Index, etc.A larger composite size would smooth the data.The advantage of larger sample size is the lower total number of samples needed to assay or test all the drill cores and thus lower cost.The synthetic sampling module does not have any sample size limitation and allows simulation of infinitely small samples.This will lead to a larger number of records in the database of the sampling module.
The difference between sampling methods is implemented through the errors ( ) added to the sampled values Equation (10).This allows for studying different sampling strategies.The number of vertical drill holes is limited by the number of voxels covering the horizontal cross-section of the The difference between a drill core, and reverse circulation (RC) chips may be considered by selecting the sample sizes (half core, quarter core) and sample preparation methods.Sample size in the drill holes is defined by the sample or composite length.A smaller composite length allows better capturing of the variability of the studied parameter, e.g., magnetite grade, Bond Work Index, etc.A larger composite size would smooth the data.The advantage of larger sample size is the lower total number of samples needed to assay or test all the drill cores and thus lower cost.The synthetic sampling module does not have any sample size limitation and allows simulation of infinitely small samples.This will lead to a larger number of records in the database of the sampling module.
The difference between sampling methods is implemented through the errors (ε) added to the sampled values Equation (10).This allows for studying different sampling strategies.The number of vertical drill holes is limited by the number of voxels covering the horizontal cross-section of the block model.A single sample can consist of one or several voxels.Non-vertical drill holes are described with dip and azimuth in addition to the total depth and collar coordinates.Synthetic drill cores are extracted as composites where composite length can be defined depending on the variability in the deposit and purpose of the study.
Synthetic element composition in a drill core sample is obtained by transferring values from the nearest voxel in a voxel model to the segments of the synthetics drill cores.Those chemical compositions do not account for the error of the assaying methods, i.e., XRF.Therefore, chemical compositions of the synthetic drill core samples are converted to chemical assay values by applying an error model.The error model is based on the precision and accuracy information for chemical analyses, e.g., for X-ray fluorescence (XRF), Equation (10): where G XRF is a component's grade analysed by synthetic XRF; G true is a synthetic value of the component grade of the sample; ε is the measurement error which can be described by a normal distribution with standard deviation σ and expected value 0. The standard deviation for each elemental assay (σ El ) is computed as a product of elemental grades (G true ) and relative standard deviation (RSD El ) as shown in Equation (11).

Malmberget Case Study
A synthetic iron ore body model is created based on results from previous extensive and careful characterisations of the structural and mineralogical information from Malmberget iron ore deposit in Northern Sweden.The deposit, upon which our model is based, consists of several tabular to stock-shaped ore bodies of massive to semi-massive magnetite and/or hematite deformed into a synformal shape [70,78,79].
The massive ore is defined by its high Fe content (55-60%) and low SiO 2 content.In the eastern part of the deposit the massive ore is surrounded by semi-massive mineralisation, characterised by a lower Fe grade (<55%) and higher SiO 2 content.The semi-massive zone can be several tens of meters thick, occurring as rims or as inclusions in the massive ore.Mineralogically, the ore is composed of magnetite and hematite as the main minerals and apatite and amphibole-pyroxene as typical gangue minerals.The semi-massive ore contains various proportions of silicates, i.e., feldspars (albite and orthoclase), amphibole, quartz, and biotite, which display a broad variation of more or less complex mineral-texture relationships [78].
The metallurgical test work was conducted on five different sample batches of >100 kg which reflects and represents the main ore types of Malmberget deposit.A textural classification was established for the massive and semi-massive ore to include both mineral-(mineral phases, mineral chemistry, and modal mineralogy) and textural information (grain size, shape, associating mineral).Due to the overall relatively high Fe-oxide content in the ore, the main gangue minerals for each ore type (normalised value) were used as a key feature to identify a variation between the textural types (each mineral selected, was related to the ore processing).The challenge was to create a geological model that offers quantitative information to be used in a process model.As mineral ore/textural type is usually descriptive, it is assumed that there is a close relationship between mineral liberation and mineral textural type.Therefore, classification of the ore handles the textural information from a processing point of view, in this case, incorporating the particle liberation distribution for each textural type and sample (by size) (Table 3).The feldspar rich textural type (Fsp) is typical for the lower grade magnetite ore (semi-massive), the amphibole rich textural type (Amph-(Ap-Bt) is widespread in the massive ore together with another common apatite-bearing textural type (Ap-(Amph) (Figure 8), details in Lund et al. [70].Abbreviations: Mgt-magnetite; Ab-albite; Act-actinolite; Ap-apatite; Bt-biotite; Qtz-quartz; Amph-amphiboles.The process model in Lund et al. [78] for each sample was made using a one-unit concentration model to quantify the mineral processing performance (grades for Fe, Si, P and Fe-recovery) in Table The process model in Lund et al. [78] for each sample was made using a one-unit concentration model to quantify the mineral processing performance (grades for Fe, Si, P and Fe-recovery) in Table 3.The process model uses textural type and modal mineralogy of the sample as inputs.The output of the model by Lund et al. [78] (Figure 9) was the forecast of elements' and minerals' recoveries to concentrate, and the concentrate's chemical and mineralogical composition [80].
The estimated error of the outputted elements was the lowest for iron and did not exceed 3% for gangue elements (Table 4).3. The process model uses textural type and modal mineralogy of the sample as inputs.The output of the model by Lund et al. [78] (Figure 9) was the forecast of elements' and minerals' recoveries to concentrate, and the concentrate's chemical and mineralogical composition [80].
The estimated error of the outputted elements was the lowest for iron and did not exceed 3% for gangue elements (Table 4).Table 4. Relative standard deviations (RSD) of the XRF assay analysis, % [80].

Geological Model
The synthetic ore body model was created in MATLAB.The deposit has three textural types, i.e., Amph-(Ap-Bt), Ap-(Amph) and Fsp (Figure 10) were spatially modelled based on characterisation described earlier for the Malmberget case study.Spatial distribution of the textural types was arranged according to the distribution suggested by Lund [80].The textural types were modelled with ellipsoids in several steps and were aligned along the centre line such that the centre of the ellipsoids were lying on the centre line.The centre line of the deposit was selected with a 15° dip that crosses the centre point of the voxel model (x = 25, y = 25, z = 25) mimicking the model proposed by [80].This model is layered (Figure 10) with textural types changing with depth.The same pattern was mimicked in the synthetic ore body model with 20 layers: 7 layers for Ap-(Amph), 7 layers for Amph-(Ap-Bt) and 6 layers for to Fsp-type.The shape of the ore body (comprised of the layered textural types) was chosen to be ellipsoidal, with the widest part in the middle (x = 25, y = 25, z = 25) and gradually narrowing with distance away from the centre.Where

Geological Model
The synthetic ore body model was created in MATLAB.The deposit has three textural types, i.e., Amph-(Ap-Bt), Ap-(Amph) and Fsp (Figure 10) were spatially modelled based on characterisation described earlier for the Malmberget case study.Spatial distribution of the textural types was arranged according to the distribution suggested by Lund [80].
Minerals 2018, 8, x FOR PEER REVIEW 14 of 23 3. The process model uses textural type and modal mineralogy of the sample as inputs.The output of the model by Lund et al. [78] (Figure 9) was the forecast of elements' and minerals' recoveries to concentrate, and the concentrate's chemical and mineralogical composition [80].
The estimated error of the outputted elements was the lowest for iron and did not exceed 3% for gangue elements (Table 4).Table 4. Relative standard deviations (RSD) of the XRF assay analysis, % [80].

Geological Model
The synthetic ore body model was created in MATLAB.The deposit has three textural types, i.e., Amph-(Ap-Bt), Ap-(Amph) and Fsp (Figure 10) were spatially modelled based on characterisation described earlier for the Malmberget case study.Spatial distribution of the textural types was arranged according to the distribution suggested by Lund [80].The textural types were modelled with ellipsoids in several steps and were aligned along the centre line such that the centre of the ellipsoids were lying on the centre line.The centre line of the deposit was selected with a 15° dip that crosses the centre point of the voxel model (x = 25, y = 25, z = 25) mimicking the model proposed by [80].This model is layered (Figure 10) with textural types changing with depth.The same pattern was mimicked in the synthetic ore body model with 20 layers: 7 layers for Ap-(Amph), 7 layers for Amph-(Ap-Bt) and 6 layers for to Fsp-type.The shape of the ore body (comprised of the layered textural types) was chosen to be ellipsoidal, with the widest part in the middle (x = 25, y = 25, z = 25) and gradually narrowing with distance away from the centre.Where  The textural types were modelled with ellipsoids in several steps and were aligned along the centre line such that the centre of the ellipsoids were lying on the centre line.The centre line of the deposit was selected with a 15 • dip that crosses the centre point of the voxel model (x = 25, y = 25, z = 25) mimicking the model proposed by [80].This model is layered (Figure 10) with textural types changing with depth.The same pattern was mimicked in the synthetic ore body model with 20 layers: Minerals 2018, 8, 536 15 of 23 7 layers for Ap-(Amph), 7 layers for Amph-(Ap-Bt) and 6 layers for to Fsp-type.The shape of the ore body (comprised of the layered textural types) was chosen to be ellipsoidal, with the widest part in the middle (x = 25, y = 25, z = 25) and gradually narrowing with distance away from the centre.Where textural types modelled in the synthetic ore body model overlap, the assignment of a voxel to a specific textural type was made by applying priorities through the use of Boolean operations (Figure 4).The highest priority in spatial modelling was given to the Ap-(Amph), then to Amph-(Ap-Bt) and then to Fsp-type.The lowest priority was given to the country rock.For instance, when ellipsoids representing Ap-(Amph) and Fsp-types were overlapping, all the voxels inside the overlapped space were assigned to the Ap-(Amph) textural type.In reality the location of the geometallurgical domains does not always follow the regular pattern, i.e., domains centres are not aligned along a line and are not necessarily oriented and dipping in the exactly same direction.Therefore, noise was added to the dimensions and locations of the individual geometallurgical domains: ±1 voxel to the location of the domain centres, ±1 voxel to the size of semi-axis (a, b, c) and ±10 • to the rotation of each geometallurgical domain in xz and yz planes.The resulting spatial distribution of textural types in the generated synthetic ore body is shown in Figure 10 and can be compared to the Figure 8 from Malmberget [80].The synthetic ore body mimics reasonably well the actual variation of textural type in Malmberget.

Mining Model
The synthetic ore body is modelled as a surface deposit, therefore open pit mining is used as the mining constraint, which is approximated with a cone model.Mining and time constraints are modelled as 13 sequential cones, representing pushbacks in a synthetic open pit mine.The overall slope angle of this open pit model is 40 • .Mined cones set the time variable in the later process simulation, when blocks are sent to the processing according to the mining plan.A production plan was generated from the mining model (Figure 11).The mining plan includes textural ore composition and the total number of minable blocks extracted by every pushback.Synthetic blocks can be used for process simulation only after extraction from the synthetic mine.The applied mining constraint allows for consideration of dilutions from the country rock and losses due to mining method limitations.However, in order to convert the block model into a fully time indexed sequence of blocks fed to the processing plant, some other aspects (i.e., blasting, loading, muck-piling, blending, stockpiling management) need to be considered.In Figure 11, an interesting pattern emerges, where mining of the synthetic ore body shows considerable variations in the ore fed to the concentrator.This is based only on the assumed geological model and textural information.

Process Model
The process model used for this case study (Figure 12) is based on [80] and is comprised of two main sections: dry processing and wet processing.Dry processing includes cobbing and size reduction in a cone crusher.Wet processing has three stages of WLIMS and two stages of grinding  The number of mined blocks of ore increases during the early production stages, before pushback 5 (Figure 11), due to the larger volume included in each sequential pushback.In addition, the amount of waste rock in those pushbacks is relatively low compared to the amount of ore.The number of blocks of ore extracted with each pushback starts to decrease after pushback 5 (Figure 11), despite the increasing volume of the following pushbacks.The number of the ore blocks extracted over time decreases because the slope angle of the pushbacks is constant, the orebody is dipping vertically, and the slope angle of mineralization is lower than slope angle of the pushbacks.Thus, the extracted volumes of the pushbacks also contain waste rock blocks.

Process Model
The process model used for this case study (Figure 12) is based on [80] and is comprised of two main sections: dry processing and wet processing.Dry processing includes cobbing and size reduction in a cone crusher.Wet processing has three stages of WLIMS and two stages of grinding and dewatering.Here, the final concentrate is the feed to the pellets plant and recycled water is returned to the head of the process.

Process Model
The process model used for this case study (Figure 12) is based on [80] and is comprised of two main sections: dry processing and wet processing.Dry processing includes cobbing and size reduction in a cone crusher.Wet processing has three stages of WLIMS and two stages of grinding and dewatering.Here, the final concentrate is the feed to the pellets plant and recycled water is returned to the head of the process.The simulation used to create the process performance forecasts was run for 12,231 blocks, using five iterative calculation rounds.The final outcome of the simulation model was an amendment to the stream file, and included information on the concentrate: modal mineralogy, chemical composition, recoveries of main commodity minerals and elements, processing time, concentrate tonnage per hour and feed tonnage per hour.These values were used as an input to the economic model.The simulated magnetite recoveries (Figure 13) show that there are three distinctive populations of the process responses influenced by textural differences (Amph-(Ap-Bt), Ap-(Amph) and Fsp).The dispersion of the simulated recoveries is due to the texture differences of the different textural types.The process performance metrics are listed in (Table 5).recoveries of main commodity minerals and elements, processing time, concentrate tonnage per hour and feed tonnage per hour.These values were used as an input to the economic model.The simulated magnetite recoveries (Figure 13) show that there are three distinctive populations of the process responses influenced by textural differences (Amph-(Ap-Bt), Ap-(Amph) and Fsp).The dispersion of the simulated recoveries is due to the texture differences of the different textural types.The process performance metrics are listed in (Table 5).Validation of the synthetic ore body model was done by visual comparison of the geological structures (textural types) proposed by [80], as shown in Figure 8, and the textural types generated from the synthetic ore body model, shown in Figure 10.Comparison between the predicted and measured commodity recoveries was done based on Figures 9 and 13.Additional possible validation could be done using geostatistical methods, e.g., comparing variograms.Variograms of the real case study and for the synthetic ore body could be compared for the entire ore body, for selected zones (e.g., oxidation zones, ore types), and for drill holes.However, this is not possible since the exact locations of the Malmberget samples' origin in the ore body are not known.

Synthetic Sampling
The sampling module of the synthetic ore body model can be used to optimise the drilling campaign, depending on the purpose of the drilling.Drilling for resource estimation can use the cost of additional drilling as an optimisation parameter.The benefits of conducting synthetic drilling campaigns are the possibility to generate any number of drill holes, with any desirable orientation and length, and added accuracy from more dense sampling campaigns would suggest the optimum Validation of the synthetic ore body model was done by visual comparison of the geological structures (textural types) proposed by [80], as shown in Figure 8, and the textural types generated from the synthetic ore body model, shown in Figure 10.Comparison between the predicted and measured commodity recoveries was done based on Figures 9 and 13.Additional possible validation could be done using geostatistical methods, e.g., comparing variograms.Variograms of the real case study and for the synthetic ore body could be compared for the entire ore body, for selected zones (e.g., oxidation zones, ore types), and for drill holes.However, this is not possible since the exact locations of the Malmberget samples' origin in the ore body are not known.

Synthetic Sampling
The sampling module of the synthetic ore body model can be used to optimise the drilling campaign, depending on the purpose of the drilling.Drilling for resource estimation can use the cost of additional drilling as an optimisation parameter.The benefits of conducting synthetic drilling campaigns are the possibility to generate any number of drill holes, with any desirable orientation and length, and added accuracy from more dense sampling campaigns would suggest the optimum drilling density.Drilling for variability studies would rather focus on sample size (composite length) and cost of samples characterisation.Increasing the number of samples, thus decreasing the composite length, would allow finding the optimum sample size for mineral grade estimate or process characterisation.Synthetic drilling might also help to resolve the issue of change of support, since variability studies of geological and process variability may require different sample sizes.In both cases, the cost of additional actions (drilling, characterisation) compared with potential gain is a border condition for optimisation.Error studies on the synthetic resource estimate based on synthetic drilling is a good way for project risk assessment and checking whether the project fulfils the accuracy requirements of its development stage, e.g., pre-feasibility, feasibility or production stage.
Sampling optimisation is not the aim of this study.Instead, the ease of generating synthetic drilling data is highlighted here.Multiple parameters can be selected and tested in synthetic drilling easily, fast and with no additional cost: Spacing between collars, regular or irregular patterns, orientation for individual drill cores (azimuth, dip), preferential composite length, assayed properties (e.g., commodity elements, gangue minerals, grain size, recovery etc.).Therefore, two synthetic sampling campaigns were produced to illustrate the flexibility of the synthetic drilling in terms of drill hole spacing and assayed parameters.In both campaigns, drill cores are dipping vertically and have regular pattern of collars locations.The first (Figure 14A) comprised a 25 drill holes pattern, where thicker composites represent higher iron grade.The second pattern (Figure 14B) comprised of 81 drill holes, with thicker composites representing higher magnetite grade and a colourbar showing the actinolite grade.The composite length for both campaigns was chosen as 1/3 of the voxel size (25/3 = 8.33 m).Higher iron/magnetite grades identify the extent of the ore body thus showing efficiency of the sampling campaign in outlining the deposit.The sampling database contains the same information as deposit database (Figure 1): Mineralogical, mining, process and economic information about each part of the drill core.All those parameters can be easily included and visualised by the sampling campaign.
Minerals 2018, 8, x FOR PEER REVIEW 18 of 23 drilling density.Drilling for variability studies would rather focus on sample size (composite length) and cost of samples characterisation.Increasing the number of samples, thus decreasing the composite length, would allow finding the optimum sample size for mineral grade estimate or process characterisation.Synthetic drilling might also help to resolve the issue of change of support, since variability studies of geological and process variability may require different sample sizes.In both cases, the cost of additional actions (drilling, characterisation) compared with potential gain is a border condition for optimisation.Error studies on the synthetic resource estimate based on synthetic drilling is a good way for project risk assessment and checking whether the project fulfils the accuracy requirements of its development stage, e.g., pre-feasibility, feasibility or production stage.Sampling optimisation is not the aim of this study.Instead, the ease of generating synthetic drilling data is highlighted here.Multiple parameters can be selected and tested in synthetic drilling easily, fast and with no additional cost: Spacing between collars, regular or irregular patterns, orientation for individual drill cores (azimuth, dip), preferential composite length, assayed properties (e.g., commodity elements, gangue minerals, grain size, recovery etc.).Therefore, two synthetic sampling campaigns were produced to illustrate the flexibility of the synthetic drilling in terms of drill hole spacing and assayed parameters.In both campaigns, drill cores are dipping vertically and have regular pattern of collars locations.The first (Figure 14A) comprised a 25 drill holes pattern, where thicker composites represent higher iron grade.The second pattern (Figure 14B) comprised of 81 drill holes, with thicker composites representing higher magnetite grade and a colourbar showing the actinolite grade.The composite length for both campaigns was chosen as 1/3 of the voxel size (25/3 = 8.33 m).Higher iron/magnetite grades identify the extent of the ore body thus showing efficiency of the sampling campaign in outlining the deposit.The sampling database contains the same information as deposit database (Figure 1): Mineralogical, mining, process and economic information about each part of the drill core.All those parameters can be easily included and visualised by the sampling campaign.

Implications and Limitations
The paper presents a synthetic ore body modelling algorithm which allows for data integration and potentially may be developed for supporting decision making in planning, sampling, testing, and beneficiation model building.The approach is not a substitute for the existing geostatistical approaches or process models but is a tool for testing the sensitivity of the project to constraints at

Implications and Limitations
The paper presents a synthetic ore body modelling algorithm which allows for data integration and potentially may be developed for supporting decision making in planning, sampling, testing, and beneficiation model building.The approach is not a substitute for the existing geostatistical approaches or process models but is a tool for testing the sensitivity of the project to constraints at different stages of the mining value chain.The impact of the variability in an ore deposit can be traced through the whole mining value chain, down to the concentrate and tailings quality.Potentially, the impact can be traced to the smelter and environmental parameters may also be considered.
The procedure presented here is aimed at geologists, geostatisticians, mining engineers, process engineers, mine planners and economists.The main benefit for geologists is the possibility to develop several hypotheses on ore classification and the drilling program.Selections of the grade estimation method such as kriging, machine learning, multivariate statistics or domaining is easier with synthetic data.Therefore, resource estimates made by geostatistics can be validated and supported.This may be done by designing several synthetic drilling scenarios with different drilling density and composite length, performing new grade estimates for each scenario, and comparing newly estimated grades with grades known from the synthetic ore body model.The mining engineers can test economic dependences of the project from mining methods and scheduling methods.Additionally, a guess can be made on possible ore losses and dilutions based on the block model resolution chosen for mine planning.Knowledge about rock properties can be used to model cracks and rock mechanical properties, which is important in blast design and rock support planning.The developed procedure is a step forward to improve blending strategies in mine planning and make an estimate of approximate grades and mineral associations to be sent to the process plant.The economic performance of the project may be estimated at different production stages, including pre-feasibility, feasibility studies, and production strategies.Finally, the sampling module may be a good proxy tool for drilling optimisation for the mining projects.
Uncertainty of the geological model (i.e., errors in sampling, measurements and geological modelling) was not addressed here.Applying simplified uncertainty or soft contacts (e.g., described by random error) may make the Figures 8 and 10 look more realistic; however, the model would lack physical realism.The elemental and mineral grades were constructed considering an apatite iron ore deposit type.Up-scaling of geometallurgical attributes was not considered since, the breakage model implemented in the simulation is specific to the ore types, and the magnetic separation process does not face an upscaling problem [81].The mining method accounted only for the simplified open pit mining.

Conclusions
All data used in modelling come from the real-life case study of an iron ore.A complete chemical, mineralogical, mining, and processing information data set was generated with respect to the spatial constraints.Applying an algorithm for spatial constraints, including textural type allocation and mining production, makes this method more realistic compared to other studies, where only mineralogical and process parameters are considered for data integration.The uniqueness of this study is that the geological model is generated on a mineralogical level.Thus, modal composition and textural information are described quantitatively for the entire deposit.The resolution of the model is defined by the model designer and is reflected in the number of records in the deposit or sampling database.Thus, highly detailed models can be developed and studied with no additional costs.
The proposed approach allows for investigating various production strategies considering the impact of geological variability, mining constraints caused by mining method, and time constraints caused by mining plan on the process variability.Error properties, traced through the mining value chain simulated with this approach, can highlight the need for additional sampling and suggest techniques to be used in mineralogical characterisation.The approach can be used as a valuable tool in risk assessment at different stages of the mining projects or as a platform for testing different scenarios for the future geometallurgical program, or for educational purposes.
Future work should include creation of a deposit library covering VMS and porphyry deposits.More parts of the mining value chain can be included, i.e., geophysics, process metallurgy, environmental assessment, economic simulations.Adding a breakage model in blasting will improve production models.

Figure 1 .
Figure 1.Structure of the synthetic ore body model.

Figure 2 .
Figure 2. Structure of the deposit module of the synthetic ore body model.

Figure 3 .
Figure 3. Ellipsoid used to model domains.

Figure 2 .
Figure 2. Structure of the deposit module of the synthetic ore body model.

Figure 2 .
Figure 2. Structure of the deposit module of the synthetic ore body model.

Figure 3 .
Figure 3. Ellipsoid used to model domains.

Figure 3 .
Figure 3. Ellipsoid used to model domains.

Figure 4 .
Figure 4. Examples of Boolean operations.The physical location of each domain (ellipsoid) in relation to other domains is defined by aligning ellipsoid centres along an auxiliary design line or placed individually by the model user.The auxiliary design line is used as reference for the ellipsoids' locations and does not have any particular geological meaning.Rock types, the same as other spatially continuous variables, are defined by modelling domains and encoding them with a unique identifier for the particular rock type.One voxel of the synthetic ore body model is occupied only by one rock type (ore type, ore class, weathering type etc.).

Figure 5 .
Figure 5. Applying geological domains to model rock types and property (e.g., mineral grade, elemental grade, density etc.) inside the rock type.

Figure 5 .
Figure 5. Applying geological domains to model rock types and property (e.g., mineral grade, elemental grade, density etc.) inside the rock type.

Figure 7 .
Figure 7. Structure of the sampling module and its connection to the synthetic ore body model.

Figure 7 .
Figure 7. Structure of the sampling module and its connection to the synthetic ore body model.

Figure 10 .
Figure 10.Spatial distribution of textural types in cross-sections generated from the synthetic ore body.A-in a plane xz, B-in a plane yz, C-in a plane xy.

Figure 10 .
Figure 10.Spatial distribution of textural types in cross-sections generated from the synthetic ore body.A-in a plane xz, B-in a plane yz, C-in a plane xy.

Figure 10 .
Figure 10.Spatial distribution of textural types in cross-sections generated from the synthetic ore body.A-in a plane xz, B-in a plane yz, C-in a plane xy.

23 Figure 11 .
Figure 11.Feed forecast generated from the mining model.

Figure 11 .
Figure 11.Feed forecast generated from the mining model.

Figure 11 .
Figure 11.Feed forecast generated from the mining model.

Figure 12 .
Figure 12.Process model design of the Malmberget simplified flow sheet in HSC SIM used to populate the synthetic ore body model with process parameters.

Figure 12 .
Figure 12.Process model design of the Malmberget simplified flow sheet in HSC SIM used to populate the synthetic ore body model with process parameters.

Table 1 .
Listing of spatial modelling approaches.

Table 3 .
The geometallurgical framework of how to simulate and forecast the metallurgical variation of ore deposit variables.

Table 3 .
The geometallurgical framework of how to simulate and forecast the metallurgical variation of ore deposit variables.

Table 5 .
Basic statistics of the modal mineralogy and elemental composition of the produced concentrate.

Table 5 .
Basic statistics of the modal mineralogy and elemental composition of the produced concentrate.