SciELO - Scientific Electronic Library Online

 
vol.89 issue220Optimization of urban land-cover classification workflow based on geographic-object analysis using very-high-resolution imageryPerformance control of root piles using a digital odometer author indexsubject indexarticles search
Home Pagealphabetic serial listing  

Services on Demand

Journal

Article

Indicators

Related links

  • On index processCited by Google
  • Have no similar articlesSimilars in SciELO
  • On index processSimilars in Google

Share


DYNA

Print version ISSN 0012-7353On-line version ISSN 2346-2183

Dyna rev.fac.nac.minas vol.89 no.220 Medellín Jan./Mar. 2022  Epub Sep 12, 2022

https://doi.org/10.15446/dyna.v89n220.97093 

Articles

Numerical model for predicting and evaluating sand production in weakly consolidated reservoirs

Modelo numérico para predecir y evaluar la producción de arena en yacimientos poco consolidados

Edson Felipe Araujo-Guerreroa 
http://orcid.org/0000-0002-8382-4878

Cristhian Bernardo Morales-Monsalveb 
http://orcid.org/0000-0003-0091-3240

Guillermo Arturo Alzate-Espinosaa 
http://orcid.org/0000-0001-6265-274X

Alejandra Arbelaez-Londoñoa 
http://orcid.org/0000-0003-0570-5125

a Universidad Nacional de Colombia, Medellin, Colombia, efaraujog@unal.edu.co, gaalzate@unal.edu.co, aarbelal@unal.edu.co

b Departamento de Ingeniería Civil y Medio Ambiente, Pontificia Universidade Católica do Rio de Janeiro, Brasil, cbmoralesm@puc-rio.br


Abstract

Sand production is a common phenomenon in oil and gas reservoirs, which occurs when reservoir fluids exert a sufficient drag force on reservoir rocks to erode the matrix. Numerical models for sand production have been used to understand the sanding mechanisms and forecast sand-production potential of formations to design well completion, optimize production, and prevent setbacks in future operations. This paper presents a mathematical model for defining the conditions of sanding onset as well as to predict and quantify the sand rate. We also introduce fluid-flow coupling and a geomechanical and sand-production model. By using the proposed model and a set of experimental data, sanding-related variables are analyzed, and a matching process for the simulated results and forecast analysis are performed. The results show that elastoplastic constitutive models are indispensable, and a clear relationship exists between the sanding and plastic strains.

Keywords: geomechanics; sanding; numerical simulation; elastoplasticity; mathematical modeling

Resumen

La producción de arena es un fenómeno común en los yacimientos de petróleo y gas. Este ocurre cuando los fluidos del yacimiento aplican una fuerza de arrastre suficiente sobre la roca yacimiento para erosionar la matriz. Los modelos numéricos para la producción de arena se utilizan para comprender los mecanismos de arenamiento y pronosticar el potencial de producción de arena de las formaciones con el fin de diseñar la terminación del pozo, optimizar la producción y evitar contratiempos en las operaciones futuras. Este artículo presenta un modelo matemático para definir las condiciones para el inicio del arenamiento, así como la predicción y cuantificación de la tasa de arena, acoplando un modelo de flujo de fluido, un modelo geomecánico y de producción de arena. Utilizando el modelo presentado y un conjunto de datos experimentales, se analizan las variables relacionadas con el arenamiento y se realiza un proceso de comparación de los resultados simulados y un ejercicio de pronóstico. Los resultados muestran que los modelos constitutivos elastoplásticos son indispensables y existe una clara relación entre la producción de arena y las deformaciones plásticas.

Palabras clave: geomecánica; arenamiento; simulación numérica; elastoplasticidad; modelamiento matemático

1. Introduction

Hydrocarbon reserves in unconsolidated or weakly consolidated formations represent 70% of world reservoir reserves [1]. During the production of these reservoirs, sanding frequently occurs, which increases the operational risks in wells and on surface facilities as well as the costs for remediation and cleaning operations. To prevent and control sanding, mechanical filters are used, which are costly and reduce productivity [2].

Sanding occurs when the drag forces exerted on a matrix by flowing fluids (during production) are greater than the resistive forces of the formation grains. Different factors influence these forces, such as rock and fluid properties, stress state around the well, and type of completion used [2]. In a model, the inclusion of all mechanisms and factors that influence sanding can result in a very complex numerical model. Thus, developing a model with lesser constraints and assumptions is necessary to analyze the most relevant mechanisms and to forecast the sand-production potential of formations in order to design a process for optimum completion of a well, optimize production, and avoid setbacks in future operations [2]. Therefore, to study sanding and its mechanisms, we need to use models that reproduce the actual behavior of the fluids and rocks [3].

The present study develops a mathematical and computational model. This model attempts to estimate the onset of sand production and quantify the rate of sand production for nonconsolidated or weakly consolidated formations. This study recognizes the importance of using an elastoplastic constitutive model to simulate the mechanical behavior of rocks. In addition, it demonstrates how the proposed model helps to perform better analysis of the sand-production phenomenon under laboratory conditions. This process is important for adjusting the parameters associated with elastoplastic behavior, including adjustment of the parameters related to porosity and permeability evolution.

Despite the current developments in sand-production modeling, the relationship between sand production and porosity-permeability of a material remains unclear. Furthermore, the mechanical behavior of sanding occurrence is not fully understood. Consequently, this study introduces a numerical model that includes current developments related to the effect of sand production on porous volume, porosity, permeability, and fluid-rate calculations.

2. Literature review

During oil and gas production, the reservoir formations are subjected to high seepage forces and sand is produced, along with the produced fluids. The literature provides two main types of models for modeling the sand-production phenomena, namely, models based on the continuum and discontinuum approaches. In addition, the solutions based on continuous mechanics consist of analytical and numerical solutions [2]. The current study focuses on the numerical models based on the continuum approach.

The numerical models use the following equations for the solutions: solid-mass balance, fluid-mass balance, and force equilibrium. These models result in a system of equations that represents the sanding phenomenology.

Early models for sanding quantification highlighted the importance of the erosional phenomena. One of the very first models stated that the volumetric solid production is a function of the erosional rate, and this erosional rate depends on the fluid velocity and solid concentration in the flow phase (eq. 1) [4].

where s is the sand-mass rate, p s is the solid density, 𝜙 is the porosity, X s is the concentration of solids in the moving fluid, X s cr is the critical solid concentration, and q i is the fluid-discharge rate.

An additional model defines the minimum fluid-velocity expression for the onset of solid production, which is defined as a function of the residual resistance and petrophysical properties [5]. The additional model includes elastoplasticity as a fundamental behavior of sanding mean, whereas the solid concentration is eliminated in eq. (1). The elastoplastic phenomena are included in the λ variable, which assumes two components: a minimum plasticization level to trigger solid production and the evolution of parameter λ proportional to the plastic strain until it reaches a maximum value [6]. Different researchers employ the minimum plasticization level for the sanding onset using experimental or numerical methods [7-9].

Several studies include sand production in their calculations and confirm that the sand-production level is related to the plasticized volume of a formation [10]. Additionally, three conditions that trigger the production of sand are proposed [11]. First, sanding occurs only at the cavity faces. Second, the element must have zero cohesion and be under tension condition. Third, the model removes the sanded element from the grid. Another assumption for the sanding models is that an element is eroded when it does not satisfy the equilibrium equations [12,13].

Another sand-production criterion is the requirement for a critical level of plastic deformation and fluid velocity to transport the detached grains [14]. The effect of the three-phase flow and capillary forces in the critical fluid velocity on sand production has also been investigated [15]. Another model introduces a pressure-gradient criterion at the pore scale for onset prediction [16]. However, this criterion is derived by mechanical analysis at the pore scale; thus, the geomechanical effect of stresses around the wellbore is not considered. In general, although, all previous models are based on force-equilibrium equations, fluid-mass balance, and solid-mass balance, they differ in the method in which sand is produced, effects of the produced sand, and involved constitutive equations.

The model presented in this research allows to asset the sand production onset and quantify the produced sand due to mechanical effects, as well as to study the effect of sand production on porosity-permeability relationship and its impact of hydrocarbons production rate, which leads to a general view of the causes and consequences of sand production.

3. Model specifications

3.1 Physical model

According to the problem's geometry, the simulation-domain geometry is a hollow cylinder with internal radius rw, external radius r ext , and thickness h (Fig. 1). This physical model is discretized in the r, 𝜃, and z directions depending on the number of required elements. In this domain, fluid flow and rock deformation simultaneously occur. The modeled formation is characterized as having high porosity, permeability, and deformability. This physical model incorporates boundary conditions similar to real conditions, such as stress at the inner radius, absence of flow at the top and bottom boundaries, and constant regional stresses at the borders of the reservoir.

Source: own elaboration.

Figure 1 Physical model and discretization.  

3.2 Differential and numerical models

The differential model consists of fluid-flow, coupled geomechanical deformation, and sand-production models.

3.2.1 Fluid-flow model

The fluid-flow model consists of the integration of the fluid- and solid-mass conservation equations with the monophasic diffusivity equation (Darcy's law) for a slightly compressible fluid. These equations, including the effects of porosity changes due to strain and solid production, result in the final equation for the fluid-flow model (eq. 2). The derivation of the fluid-flow model is shown in Annex A.

where p f is the fluid density, k is the permeability tensor, 𝜇 is the viscosity, p is the pressure, 𝜙 is the porosity, c f is the fluid compressibility, c bc is the bulk compressibility, c s is the solid compressibility, α is the Biot constant, t is the time, e v e is the elastic volumetric strain, e v p is the plastic volumetric strain, s is the mass rate of sand production, V b is the bulk volume, and p s is the solid density.

3.2.2 Geomechanical-deformation model

During the injection and production processes at the wellbore or the loading and unloading processes in the laboratory tests, the reservoir deforms. To develop a coupled model, three components must be integrated, namely, the stress-equilibrium equation (eq. 3), strain and displacement compatibility conditions (eq. 4), and constitutive model that integrates the stresses, strains, and pressure (eq. 5) [17].

where 𝝈 is the total stress tensor, Є is the total strain tensor, u is the displacement vector,G is the shear modulus, λ is the Lame coefficient, and I is the Kronecker delta. The geomechanical-deformation model for a porous medium with elastoplastic behavior (eq. 6) is obtained from the introduction of eq. (4)-(5) into eq. (3), which yields a model that considers the pressure and plastic strains.

The plastic strains (dЄp) are calculated during material yielding. Then, the displacements are recalculated so that the new calculated stresses satisfy the condition f s i ' N ) = 0. Among the different failure criteria, the Mohr-Coulomb criterion is the most conservative in its predictions, for this reason, Mohr-Coulomb criterion used in this model (eq. 7).

where σ1 and σ' 3 are the maximum and minimum principal effective stresses, respectively, φ is the internal friction angle, and c is the cohesion of the material.

3.2.3 Hardening and softening parameters

The presented model uses the shear-hardening parameter (e ps ) for the accumulation of the shear plastic strains to quantify the level of mechanical damage due to shearing.

where ė ps is the rate of e ps , J 2 is a tensor invariant, and Ė p is the rate of the plastic-strain tensor.

3.2.4 Sand-production model

This section defines the sanding onset conditions and the amount of produced sand. Sand onset occurs at a specific level of plastic strain, which is known as the critical level of plastic failure (ecr ps). ecr ps = 0 implies that sanding onset occurs immediately after "plasticity" starts, which is a common assumption of the analytical models.

This model assumes that sand is produced only from the cavity faces exposed to fluid flow because sand production is an erosive phenomenon and the eroded grains are large; thus, sand cannot come from inside the reservoir. An innovative element in this model is the relationship between the mass of produced sand and level of plastic failure e ps . This relationship suggests that the amount of sand available for production is proportional to the level of plastic failure associated developed in the porous medium. Finally, the model establishes a minimum flow velocity for the detached grains to be transported by the fluid, i.e., sand is only produced when the fluid velocity exceeds the critical value to transport it ( v f > v f min ). Considering the abovementioned conditions, the rate of sand production is defined as:

Where μ s is an experimental parameter. A f is the area of the sand-producing zone perpendicular to the flow and helps define a specific sand-production value (sand production per unit area), which can be used to extrapolate the results to cases such as wellbores and perforation tunnels.

In eq. (9) a high flow velocity does not imply that solids are being produced. In other words, a solid matrix must be subjected to a certain plastic state and a flow velocity for the solids to be produced. Under sand-production conditions, an increase in the level of plasticity (e ps ) represents an increase in sand production. In eq. (10), the new porosity (𝜙 n+1 ) is a function of the mean stress, pressure, strains, and sand-production changes.

3.3 Initial and boundary conditions

As initial conditions, the pressure at each calculation node is equal to the initial pressure, besides, the initial displacements and plastic strains are equal to zero (eq. 11) (condition of an intact reservoir). The initial-stress state distribution includes axial and horizontal stresses with an isotropic stress condition at the horizontal plane (eq. 12).

where is the initial displacement vector. In terms of the boundary conditions, the fluid-flow model denotes the pressure at the inner (p W f) and outer (p e ). Meanwhile, the solution uses no-flow boundary for the top and bottom boundaries. The boundary condition for the geomechanical model at the inner radius is σ rr = p wf , whereas that at the outer radius is Δσ rr = 0. Meanwhile, the top-boundary vertical displacement is set equal to zero.

3.4 Computational model

The proposed model is built in-house using the FORTRAN programming language. In addition, it is solved using the finite difference, fully coupled, and iterative method. The solution of the coupled system includes two convergence levels. In the first level, the solution of the pressure and displacement is obtained. When the first convergence level is finished, the stresses at time t n+1 are calculated, and the plasticity module evaluates the plastic deformations. In this manner, the second convergence level considers the plastic strains in the solution (see Fig. 2).

Source: own elaboration.

Figure 2 Flowchart of the coupled-solution system.  

4. Sand-production-phenomenon modeling

The fundamental step is the matching of the results of the proposed model with real data, which can be achieved using either the laboratory or field data. In this study, the results of an experimental sand-production test are used [18]. In this test, the mechanical properties for the mechanical behavior of a set of synthetic porous samples are characterized using an axial-deformation test. A second sample, similar to the first one except for the shape, is subjected to a hollow-cylinder test in which a pressure gradient through the sample and a progressive and controlled stress state are induced, which results in a fluid flow with progressive transport and production of sand grains in the fluid stream. The test results are used to determine and adjust the parameters of the sand-production model.

The model proposed in this paper is initially used to reproduce the mechanical behavior, pressure, and sand production recorded during the aforementioned sand-production test. The mechanical response of the material is adjusted according to the properties that control the elastoplastic behavior of the sample during the test. This study proposes to match the value of the initial parameters of cohesion and the internal friction of the material based on the axial and radial stresses recorded during the test. Subsequently, the sand production is adjusted by defining the parameters that control the proposed sand-production model, this to control the onset and amount of sand production.

5. Analysis of sanding tests

This section presents the analysis of the sand-production test performed to validate the proposed model. The sand-production test was performed on a hollow-cylinder sample, which was made from a mixture of water, sand, and cement [18]. The sample had external and internal diameters of 125 and 25.4 mm, respectively. Table 1 lists the mechanical and petrophysical properties.

Table 1 Sample mechanical and petrophysical properties-base case. 

Source: Nouri et al., 2006.

During this test, the axial and radial stresses were independently controlled, and the sample was subjected to a pressure gradient to induce a high fluid velocity. The fluid flowed from the external radius to the internal radius. The pressure in the internal diameter was constant and equal to the atmospheric pressure, whereas that in the external diameter was varied. The produced fluid and sand were quantified over the duration of the entire test. Fig. 3 shows the results recorded during the test, which included the pressure drop between the external and internal radii, applied axial and radial stresses, amount of produced sand, and axial deformation of the sample.

Source: Adapted from Nouri et al., 2006.

Figure 3 General results of the sand-production test.  

As part of this work, Fig. 3 shows the tangential stress at the inner diameter of the sample. In this case, the tangential stress was calculated at the internal radius of the sample by assuming that the material elastically behaved. This condition considered that σ r inn = 101 kPa, where σ r inn is the radial stress in the internal radius during the test.

To study the sand-production behavior during the test, the sand-production data were divided into four stages (E1, E2, E3, and E4). Each stage had initial and stabilized values of the sand-mass product, as shown in Fig. 4.

Source: Nouri et al., 2006.

Figure 4 Stages of sand production according to the behavior of the accumulated sand production.  

On the basis of our calculation of the tangential stress, the stresses in the inner diameter satisfied the condition that σ v > σ θ int > σ r int , where σ v is the vertical stress, σ θ int is the calculated internal tangential stress, and o r int is the internal radial stress. Using these stresses, we evaluated the Mohr-Coulomb failure function (f s ), i.e., eq. (7) in the inner diameter of the sample (using the data listed in Table 1, namely, pore pressure of 0.1 MPa and Biot coefficient of 1.0). The calculated Mohr-Coulomb failure function (f s ) is shown in Fig. 5. On the basis of the plasticity theory, the failure function cannot have positive values (only negative or zero). However, in this case, this function contained positive values because we evaluated it using the calculated stresses for a totally elastic-material behavior (not plastic). In this manner, the positive values of the f s function indicated that the material exhibited plastic behavior.

Source: own elaboration.

Figure 5 Calculus of the failure function in the internal face of the sample.  

In Fig. 5 the plastic behavior started very early in the test, and the sample then accumulated plastic strains and extended the plastic radius, causing in intense material disaggregation in the inner diameter and sand production when the fluid flowed transporting the disaggregated solids.

Fig. 6 shows the relationship between the accumulated sand production with the calculated failure function and the measured pressure drop in the four stages of sand production. Fig. 6 shows that regardless of the behavior of the pressure drop across the sample, a strong relationship existed between the amount of produced sand and the calculated failure function; however, no clear relationship existed between the pressure drop and accumulated sanding response. This condition indicated that the increase in the produced sand was associated with the increase in mechanical loading. This assertion indicated that for real-field applications with sanding events, efforts should be focused on determining the loading state of the formation by knowing the real stress path of the different operations, including those involving temperature changes.

Source: own elaboration.

Figure 6 Correlation of produced sand versus shear function.  

In contrast to the failure-function calculations shown in Fig. 5, the proposed numerical model calculated the plastic strains such that f s = 0 and f s = 0, which was the reason for the use of plastic strains to correlate the production of sand (instead of function f s , as shown in Fig. 6 ). In addition, this condition also explained the definition of a plastic-strain indicator for the sanding-rate criteria. As previously stated, we argued that the sand production in this test was a phenomenon that largely depended on the level of failure of the material near the internal cavity, whereas the moving fluid was responsible for transporting the released grains.

6. Results

6.1 Sand-production modeling

The experimental sand-production test is simulated using the proposed model. The simulated physical model is discretized using 10, 10, and 40 elements in the axial, tangential, and radial directions, respectively. The shape of the modeling elements is a section of a hollow circumference with a height of Δz. The discretization in radial direction is fine enough to get small variations in the calculated variables. In addition, a logarithmic radial discretization is used to obtain finer elements near the inner radius where larger changes in pressure and displacements occur. Table 2 lists the initial conditions for modeling, and the parameters described in Table 1 represent the properties of the base case for modeling. Table 3 lists the fluid properties for simulation.

Table 2 Initial condition and discretization parameters. 

Initial pressure (MPa) 0.1
Radial, tangential and axial stresses (MPa) 0.1
Blocks in radial, tangential and axial directions 40, 10, 10
Block thickness in axial direction [Δz] (m) 0.025

Source: own elaboration.

Table 3 Fluid properties for simulation. 

Property Value
Fluid density [pf] (g/cm3) 1.0
Fluid viscosity [μ] (cP) 1.0
Fluid compressibility [c f ] (kPa-1) 5.80 x 10-7

Source: own elaboration.

Fig. 7 shows the applied boundary conditions. At the top of the model, a displacement boundary condition is applied, thus the axial stress is a simulation result. At the internal radius, the radial stress and pressure are equal to atmospheric (0.1 MPa); at the external radius, the radial stress and pressure are equal to the experimental test data. Non-flow condition is used at the top and bottom faces.

Source: own elaboration.

Figure 7 Physical model and boundary conditions.  

The base case uses the basic parameters for the sand-test modeling as well as the initial and boundary conditions (Fig. 7 and Tables 2 and 3). The axial stress is used to compare and adjust the model results with the experimental results.

The model results for the base and elastic cases are shown in Fig. 8. For the elastic case, the plasticity module is deactivated so that the plastic strains are zero. Comparison of the elastic-case results with real data reveals two scenarios or regions within the test. The first region occurs when the axial strains are lower than 1%. Here the model reproduces the experimental data, and the mean relative error is around 5.2%. The second scenario occurs for axial strains greater than 1% where the simulated axial stresses overestimate the reported experimental results. This result reflects the limitation of the elastic model in representing the stress-strain behavior in highly complex processes where elastic behavior is not the dominant effect and a marked behavior of plasticity exists on the sample.

Source: own elaboration.

Figure 8 Modeling of the base and elastic cases.  

As presented, the reported base-case results demonstrate that the simulated rock appears to have a greater resistance than a real rock, which implies that the real rock loses its initial-strength properties during the plastic-deformation process. This result indicates that constant calibration of the constitutive-model parameters for field operations is critical in predicting the related effects of different operations and procedures over the zones near a wellbore. Therefore, to perform adequate modeling of a real rock, the weakening of the material must be taken into account.

To achieve the best fit between the simulated and recorded test results, a systematic process matching is performed using cohesion c, internal friction angle 𝜑, and dilatation angle 𝛹. This process also facilitates sensitivity analysis of the model, which yields the conclusion that the results are strongly dependent on the internal friction angle, very slightly dependent on the cohesion, and independent of the dilatation. This result emphasizes the importance of accurate determination of the friction angle and model calibration based on the factional angle. Using this methodology for matching parameters, we conclude that the sample also suffers weakening of its strength properties during the test. The best fit is obtained at = 413.7 kPa, 𝜑 = 27.5° and 𝛹 = -7°, for a relative error of 4% (Fig. 9).

Source: own elaboration.

Figure 9 Effect of friction angle on the simulated axial stress.  

6.2 Sand-production adjustment

Regarding the adjustment of the produced amount of sand in the test, the starting point is the adjustment of the stress-strain parameters. Following the same manner in which the correlation between the failure function (on the inner face of the synthetic sample) and sand production is obtained, a correlation between the shear-hardening parameter (e ps ) and experimental data of sand production is obtained by plotting the data of the accumulated sand production against parameter e ps of the test (Fig. 10). A linear correlation is generated for the data in the test, and a remarkable and decisive increase in sand production is verified. From this correlation, the base parameters that control the sand-production function presented in eq. (9) are obtained, namely, μ,A f ,e ps cr , and μ f min .

Source: own elaboration.

Fig. 10 Correlation between sand production and εps . 

Parameter μ s refers to the slope of the correlation, and parameter e ps cr is the value of the critical shear-hardening parameter that determines the sanding onset. Af is the internal area of the sample, and μ f min is the minimum flow velocity to transport the released sand grains. Parameter μf min is calculated using Darcy's law at the inner face of the sample. Table 4 lists the values of the mentioned parameters.

Table 4 Parameters of the sand-production module. 

Property Value
μs (kg/s) 0.348
A f (m2) 0.020
eps cr(adim) 0.056
μfmin (m/s) 0.00227

Source: own elaboration.

The values of the base parameters associated with the sand-production function are the starting set of the calculation and adjustment of the calculated data of the sand production, which are compared with the data recorded during the laboratory test. In the adjustment process, a sensitivity analysis of parameters μ s and μ f min is performed, which is demonstrated to significantly affect the different characteristics of the curve. Fig. 11 shows the results of the adjustment between the calculated sand production and actual production. The best adjustment is achieved with an error of 3% in the accumulated sand production at the following set of parameters: μ s = 0.45 kg/s, e ps cr = 0.056, μ f min = 0.00225 m/s, and A f = 0.020 m2.

Source: own elaboration.

Figure 11 Modeling cases of sand production versus a real case.  

6.3 Sand-production prediction

To verify the model behavior and based on the parameters and adjustment obtained for the axial stress and sand production against the actual sanding test, a sand-production forecast exercise is carried out. This forecast starts at the end of the real sand-production test (axial strain of 4%, see Fig. 11) and extends up to an axial strain of 5% using the same displacement velocity. Additionally, for comparison purposes, four simulation cases are defined (Table 5) where each case has a different increase rate in the external radial stress (based on the radial stress at an axial strain of 4%- 4654 kPa). In other words, at the end of the forecasting (axial strain of 5%) Cases 1-4 increase their external radial stress by 30%, 18%, 15%, and 8%, respectively. In all cases, the flow through the sample is modeled using the same pressure-drop program, which satisfies the criterion of critical flow velocity for sanding. This pressure-drop program models three stages of radial flow with interbedded stages under a no-flow condition. Fig. 12 shows the pressure-drop and radial stresses programs in each simulation case.

Table 5 Stress definitions for forecasting the cases. 

Case Increase of external radial stress at (ε α = 5%)
1 30%
2 18%
3 15%
4 8%

Source: own elaboration.

Source: own elaboration.

Figure 12 Radial stress in the sand-production forecast. 

Fig. 13 shows the forecast of sand production for each proposed case. As expected, the level of sand production is a function of the applied radial stress because of the greater concentration of stresses in the inner radius when the external radial stress is increased. In this manner, we can infer that the level of plasticization in Case 1 is more intense than those in the other cases, and this case generates greater sand production. In addition, the sand production is triggered by the fluid flow (sanding does not occur without a fluid flow), whereas the stress-load process (application of vertical displacement and radial stress) generates a constant release of sand. On the other hand, the change from a closed system to a system with a fluid flow (axial strain of 4.25%) generates production of additional sand, which is disaggregated due to the different radial stresses used in each simulation case. This result is reflected in the higher initial slope at the moment when sanding production steps up to an axial strain of 4.25%.

Source: own elaboration.

Figure 13 Sand production for each forecasting case. 

7. Conclusions

According to the laboratory results, using a poroelastic approach, sand production is related to resolving the changes in the yield function using the stress paths associated to the infield processes.

Despite its limitations, the developed model has predicted and quantified the sand production in a laboratory-scale sanding test and requires the adjustment of parameters such as the mechanical properties of a porous medium. Using the proposed model, comparison of the simulated and actual results achieves a relative error of almost 4%, which shows that the elastoplastic constitutive relationships with softening or hardening are essential for modeling weakly consolidated or nonconsolidated formations. In this analysis, the internal friction angle exerts the largest effect on the plastic behavior of the material, followed by cohesion and finally by the dilatation angle, which reflects the great significance of the friction angle in the sand-production predictions.

The elastic constitutive relationships are strongly limited when the behavior of sanding prone formations is modeled. However, the elastoplastic constitutive relationships are more adequate as long as good parameter determination is performed. In addition, modeling of the softening or hardening phenomenon is essential to achieve more coherence and representativeness of the results. This work demonstrates that the strength weakening that occurs in the triaxial tests also occurs in the hollow-cylinder test. This result should be considered in near-wellbore regions because the strength parameters decrease compared with their initial values during the production time and operating procedures.

A sand-production study is carried out, which considers different load-stress and fluid-flow schemes for a real-life experimental test. We find that sand production is directly affected by the level of applied radial stress. Further, this study finds that opening and closing the flow causes an increase in the level of produced sand, similar to the sanding cases found in cyclical injecting or producing wells where the production rate is stopped several times.

Finally, the application of this model in a field case enables quantification of the sanding level and understanding of the main mechanisms from a phenomenological perspective. Moreover, this model can help in accurately predicting possible events and avoiding them especially for unconsolidated or weakly consolidated formations.

This model fully couples fluid flow, geomechanical and sanding models, which allows not only the prediction of sanding onset but also sand quantification. The interaction enabled by this model between the different variables and phenomena is essential for correct management of the production rate and pressure in a wellbore if sanding is completely undesired or for proper design of a completion tool to control sand production.

Acknowledgments

The authors would like to thank to the Hydrocarbon National Agency of Colombia ANH and the Colombian science and technology department COLCIENCIAS (currently MINCIENCIAS) for their financial support through the contract number 015-2016.

References

[1] Nouri, A., Vaziri, H.H., Belhaj, H.A. and Islam, M.R., Comprehensive transient modeling of sand production in horizontal wellbores. SPE Journal, 12(04), pp. 468-474, 2007. DOI: 10.2118/84500-PA [ Links ]

[2] Rahmati, H., Jafarpour, M., Azadbakht, S., Nouri, A., Vaziri, H., Chan, D. and Xiao, Y., Review of sand production prediction models. Journal of Petroleum Engineering, 2013, art. 864981, 2013. DOI: 10.1155/2013/864981 [ Links ]

[3] Morita, N., Whitfill, D.L., Massie, I. and Knudsen, T.W., Realistic sand-production prediction: numerical approach. SPE Production Engineering, 4(01), pp. 15-24, 1989. DOI: 10.2118/16989-PA [ Links ]

[4] Vardoulakis, I., Stavropoulou, M. and Papanastasiou, P., Hydro-mechanical aspects of the sand production problem. Transport in Porous Media, 22(2), pp. 225-244, 1996. DOI: 10.1007/BF01143517 [ Links ]

[5] Detournay, C., Numerical modeling of the slit mode of cavity evolution associated with sand production. SPE Journal, 14(4), pp. 797-804, 2009. DOI: 10.2118/116168-PA [ Links ]

[6] Papamichos, E., Vardoulakis, I., Tronvoll, J. and Skjasrstein, A., Volumetric sand production model and experiment. International Journal for Numerical and Analytical Methods in Geomechanics, 25(8), pp. 789-808, 2001. DOI: 10.1002/nag.154 [ Links ]

[7] Wang, J., Walters, D., Wan, R.G. and Settari, A., Prediction of volumetric sand production and wellbore stability analysis of a well at different completion schemes. Alaska Rocks 2005, The 40th US Symposium on Rock Mechanics (USRMS), American Rock Mechanics Association, Anchorage, Alaska, USA, 2005. [ Links ]

[8] Gui, F., Khaksar, A., Zee, W.V.D. and Cadogan, P., Improving the sanding evaluation accuracy by integrating core tests, field observations and numerical simulation. SPE Asia Pacific Oil and Gas Conference and Exhibition, Society of Petroleum Engineers, Perth, Australia, 2016. DOI: 10.2118/182499-ms [ Links ]

[9] Santana, C. and Likrama, F., Workflow on incorporating thick-walled cylinder test results in finite element models of near wellbore for sanding prediction studies. 50th US Rock Mechanics/Geomechanics Symposium. American Rock Mechanics Association, Houston, TX., 2016. [ Links ]

[10] Wu, B., Choi, S.K., Denke, R., Barton, T., Viswanathan, C., Lim, S. and Jadid, M.B., A new and practical model for amount and rate of sand production estimation. Offshore Technology Conference Asia. Offshore Technology Conference, Kuala Lumpur, Malaysia, 2016. DOI: 10.4043/26508-ms [ Links ]

[11] Nouri, A., Kuru, E. and Vaziri, H., Elastoplastic modelling of sand production using fracture energy regularization method. Journal of Canadian Petroleum Technology, 48(04), pp. 64-71, 2009. DOI: 10.2118/09-04-64 [ Links ]

[12] Kim, S.H., Sharma, M.M. and Fitzpatrick, H.J., A predictive model for sand production in poorly consolidated sands. International Petroleum Technology Conference, Bangkok, Thailand, 2011. DOI: 10.2523/iptc-15087-ms [ Links ]

[13] Zhou, Z.Y., Yu, A.B. and Choi, S.K., Numerical simulation of the liquid-induced erosion in a weakly bonded sand assembly. Powder Technology, 211(2), pp. 237-249, 2011. DOI: 10.1016/j.powtec.2011.04.029 [ Links ]

[14] Wang, H., Cardiff, P. and Sharma, M.M., A 3-D poro-elasto-plastic model for sand production around open-hole and cased and perforated wellbores. 50th U.S. Rock Mechanics/Geomechanics Symposium, American Rock Mechanics Association, Houston, TX., 2016. [ Links ]

[15] Wang, H. and Sharma, M.M., A fully 3-D, multi-phase, poro-elasto-plastic model for sand production. SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, Dubai, UAE, 2016. [ Links ]

[16] Fan, Z. and Yang, D., Quantification of sand production by use of a pressure-gradient-based sand failure criterion. SPE Canada Heavy Oil Technical Conference, Calgary, Alberta, Canada, February, 2017. DOI: 10.2118/185009-ms. [ Links ]

[17] Chen, W.F. and Baladi, G.Y., Soil plasticity: theory and implementation (Vol. 38). Elsevier, 1985. DOI: 10.1016/c2009-0-09594-4 [ Links ]

[18] Nouri, A., Vaziri, H.H., Belhaj, H.A. and Islam, M.R., Sand-production prediction: a new set of criteria for modeling based on large-scale transient experiments and numerical investigation. In: SPE Annual Technical Conference and Exhibition. SPE J, 11, pp. 227-237, 2006. DOI: 10.2118/90273-PA [ Links ]

[19] Ghabezloo, S., Sulem, J. and Saint-Marc, J., Evaluation of a permeability-porosity relationship in a low-permeability creeping material using a single transient test. International Journal of Rock Mechanics and Mining Sciences, 46(4), pp. 761-768, 2009. DOI: 10.1016/j.ijrmms.2008.10.003 [ Links ]

How to cite: Araujo-Guerrero, E.F., Morales-Monsalve, C.B., Alzate-Espinosa, G.A. and Arbelaez-Londoño, A. Numerical model for predicting and evaluating sand production in weakly consolidated reservoirs.. DYNA, 89(220), pp. 54-63, January - March, 2022.

E.F. Araujo-Guerrero, received the Bs. Eng in Petroleum Engineering in 2012, the MS degree in Petroleum Engineering in 2015, both Universidad Nacional de Colombia - Medellín Campus (Unalmed), and is currently a PhD candidate in Energy Systems at the same university. Since 2012, he has worked as a researcher, at the Research group in Applied Geomechanics (GIGA) at National University, in several projects in applied geomechanics, such as hydraulic fracturing, sand production, heavy oil recovery, thermal recovery, mechanical earth model development, and in-situ stress state determination. His research interests include: geomechanics, reservoir simulation, constitutive modeling, sand production, and hydraulic fracturing, non-conventional reservoirs. ORCID: 0000-0002-8382-4878

C.B. Morales-Monsalve, is a Petroleum Engineer from Universidad Nacional de Colombia - Medellín Campus (Unalmed), with a master's degree in Petroleum Engineering from the same university and D.Sc. candidate at the Pontifical Catholic University of Rio de Janeiro (PUC-Rio). He has worked as a researcher in several projects in applied geomechanics, such as hydraulic fracturing, sand production, heavy oil recovery, wellbore stability and estimation/prediction of in-situ stresses. He currently works in the Petroleum Technology and Engineering Group (GTEP) where he has participated in the development of software for the evaluation of wellbore stability, in wellbore stability studies in Brazilian basins and in the evaluation of the mechanical behavior of shale and evaporite formations. ORCID: 0000-0003-0091-3240

G.A. Alzate-Espinosa, Bs. in Petroleum Engineering in 1989 at Universidad Nacional de Colombia - Medellín Campus (Unalmed). MSc in Petroleum Engineering at New Mexico Institute of Mining and Technology, New Mexico, USA in 2001. From 1988 to 1993, he works for oil and gas operating companies in Colombia and Venezuela, and since 1994 was enrolled as a professor at Universidad Nacional de Colombia, Medellin. Currently, he is associate professor in the Process and Energy Department, Facultad de Minas, Universidad Nacional de Colombia, Medellin. His research interests include: reservoir geomechanics, reservoir modelling and engineering, formation evaluation by openhole logging, sanding management and control, modelling of hydraulic fracturing, in-situ stress-state evaluation. ORCID: 0000-0001-6265-274X

A. Arbeláez-Londoño, received the Bs. Eng in Petroleum Engineering in 1998 from Universidad Nacional de Colombia - Medellín Campus (Unalmed); the MS degree in Petroleum Engineering in 2001 from the New Mexico Institute of Mining and Technology, Socorro, NM, US, and currently she is PhD candidate in Energy Systems at the Universidad Nacional de Colombia - Medellín Campus (Unalmed). From 2006 to 2017, she was an Assistant Professor in in Petroleum Engineering program at the Process and Energy Department, Facultad de Minas, Universidad Nacional de Colombia; also, from 2010 to date, she has worked as researcher at the Research group in Applied Geomechanics, GIGA for the Universidad Nacional de Colombia. Her research interests include: geomechanics, well stability, reservoir simulation, sand production, and hydraulic fracturing. ORCID: 0000-0003-0570-5125

Annex A-Fluid-flow model derivation

The fluid-mass and solid-mass balances are expressed by eq. (A1)-(A2).

To define the solid and porous volumes due to sand production, we assume that all the produced sand represents an increase in the porous volume (eq. A3).

Then, applying the material derivative and assuming a constant solid density, eq. (A2) results in:

With Darcy's law includes the solid velocities (eq. A5).

Applying the material derivative over eq. (A1) and including eq. (A5) results in:

After deriving and adding eq. (A4), eq. (A6) becomes

The porous-volume is defined by considering the elastoplastic effects [19]. Further, the term m s includes the effect in porous volume due to sanding.

The fluid is modeled as slightly compressible (eq. A9).

Finally, the flow model (eq. A10) is obtained by adding eq. (A8)-(A9) into eq. (A7), which accounts for the elastic, plastic, and sanding effects on the flow behavior.

where P f is the fluid density, k is the permeability tensor, μ is the viscosity, p is the pressure, 𝜙 is the porosity, Cf is the fluid compressibility, c bc is the bulk compressibility, c s is the solid compressibility, is the Biot constant, is the time, is the elastic volumetric strain, e v v is the plastic volumetric strain, s is the sand-production rate, V b is the bulk volume, V p is the porous volume, V s is the solid volume, p s is the solid density, v s is the solid velocity, and v f is the fluid velocity.

Received: July 04, 2021; Revised: November 26, 2021; Accepted: December 17, 2021

Creative Commons License This is an open-access article distributed under the terms of the Creative Commons Attribution License