The production curve of tight gas started in 1964 (Figure 4). The fitting curve of tight gas production maintains a trend of increasing first and then decreasing, and the increasing and decreasing periods are almost symmetrical about the peak value. This is consistent with the prediction curves of the Hubbert and Gauss models, and therefore, the study of tight gas production trend can refer to these two models.

As can be seen from the curve, one cycle of mining work has been completed from 1964 to 2020. However, in 2020, the potential of tight gas resources in the Sichuan Basin was re-evaluated, and the central Sichuan area was selected as the core production area for further exploration. At this time, only the production in 2021 and the historical production of the previous cycle are used as a reference. Therefore, the Ward model is used to predict the production peak and peak year of tight gas, and then the models are established to predict the change trend of tight gas production.

Before carrying out tight gas production predestining work, it is necessary to first estimate the numerical range of the ultimate recoverable reserves URR. Through geological exploration, it is found that the current resource of tight gas in the Sichuan Basin is 6.9 × 10 12 m 3 . By analyzing the exploration and development laws of many gas reservoirs in the Sichuan Basin, the range of the proven rate and recovery rate of tight gas is selected. According to the current tight gas exploration and mining technical conditions in the Sichuan Basin, the proven rate of tight gas is in the range of 10%–20%, That is, the range of cumulative proven reserves at the end of the life cycle is ( 0.69 − 1.38 ) × 10 12 m 3 . At present, the recovery factor of tight gas in the Sichuan Basin is about 20%, so the estimated URR range is ( 1380 − 2760 ) × 10 8 m 3 . After determining the range of URR, using it as the boundary condition, according to the principle of Ward model, the variation trend of tight gas production is predicted.

Select the proven rates of 10%, 12.5%, 15%, 17.5%, and 20%, respectively, and calculate the URR values corresponding to each proven rate (Table 1). The growth trend of tight gas production under these five scenarios with different proven rates is studied.

Substitute different URR values into Eq. 1 respectively, initial production and historical cumulative production data are known. The growth rate remains at 10% initially, and then gradually decreases (Höök et al., 2012), without considering the extreme cases of the proven rate, when URR= ( 1725 − 2415 ) × 10 8 m 3 , the production prediction curve and fitting analysis curve are shown in Figure 5.

In Figure 5A, the stable production period is roughly between 2038 and 2046 (between the black dotted lines), with peak production occurring in 2042 or 2043. In Figure 5B, for the prediction results under different URR, the best fit with historical production is URR= ( 2070 − 2415 ) × 10 8 m 3 , the higher the degree of fit, the higher the accuracy of the prediction. The final peak production result should choose the prediction value with a high degree of fit. Therefore, the predicted peak production is ( 70 − 91 ) × 10 8 m 3 / a , the predicted results can be used as reference values for the predictions of the Hubbert and Gauss models. Next, based on the above work, the production trend of tight gas is predicted under five different URRs, and different stages are analyzed.

Using the Hubbert and Gauss model, the peak production of tight gas is predicted and its variation trend is studied. The five URRs were respectively substituted into Eqs 5, 11 to establish the tight gas production-time relationship of Hubbert and Gauss model (Eqs 15, 16) and tight gas production prediction results based on the Hubbert and Gauss models (Figures 6A, 7A). Among them, Eq. 15 is the production-time relation of Hubbert model, and Eq. 16 is the production-time relation of Gauss model. The model parameters covered by these formulas are: annual peak production Q m , peak production time t m , peak slope b and standard deviation s . Eq. 15 corresponds to Figures 6A, and Eq. 16 corresponds to Figures 7A. Q = { 2 × 33.02 1 + cosh [ 0.0957 ( t − 2042 ) ] , U R R = 1380 × 10 8 m 3 2 × 49.29 1 + cosh [ 0.1143 ( t − 2042 ) ] , U R R = 1725 × 10 8 m 3 2 × 68.00 1 + cosh [ 0.1314 ( t − 2042 ) ] , U R R = 2070 × 10 8 m 3 2 × 86.03 1 + cosh [ 0.1425 ( t − 2042 ) ] , U R R = 2415 × 10 8 m 3 2 × 105.98 1 + cosh [ 0.1536 ( t − 2042 ) ] , U R R = 2760 × 10 8 m 3 , (15) Q = { 38.01 × e − ( t − 2042 ) 2 / ( 2 × 14.48 2 ) , U R R = 1380 × 10 8 m 3 54.43 × e − ( t − 2042 ) 2 / ( 2 × 12.64 2 ) , U R R = 1725 × 10 8 m 3 71.00 × e − ( t − 2042 ) 2 / ( 2 × 11.63 2 ) , U R R = 2070 × 10 8 m 3 88.00 × e − ( t − 2042 ) 2 / ( 2 × 10.95 2 ) , U R R = 2415 × 10 8 m 3 105.0 × e − ( t − 2042 ) 2 / ( 2 × 10.49 2 ) , U R R = 2760 × 10 8 m 3 . (16)

In order to deeply analyze the growth law of tight gas production with URR, it is necessary to study the relationship between the above model parameters and URR. Therefore, one hundred different URR values were uniformly sampled within the estimated range of U R R = ( 1380 − 2760 ) × 10 8 m 3 (Section 3.1). Substitute these URR values into Eqs 5–11 to calculate different production prediction results (Figures 6A, 7A).

It can be seen that the production prediction results in Figures 6A, 7A are highly similar. Therefore, both production prediction results can be analyzed simultaneously. The value of peak production is positively proportional to the value of URR, that is, the greater the URR is, and the greater the peak production is. The peak production always occurs in 2042, that is, t m ≡ 2042 , and the production growth curve is axisymmetric about t = 2042 .

According to the increase amplitude of tight gas production, the future production trend can be divided into two stages, that is, the production rising stage (2022–2038) and the production stable stage (2038–2046). In 2022–2038, production rises rapidly as the year increases. In 2038–2042, the production maintains steady growth and reaches the peak production Q m in 2042. In 2042–2046, production remains in a steady decline.

In order to compare the accuracy of production prediction results under the two different models, correlation analysis is required. This paper makes a fitting analysis diagram for the results in the early stage of prediction (Figures 6B, 7B). The fitting result represents the coincidence degree of predicted prediction and historical prediction. It can be intuitively observed from the figure that the better the fitting degree is, the closer the prediction is to the real situation, and that is, the more accurate the prediction result is.

As can be seen from the figure, for the Hubbert model, the best fit with historical production is URR= ( 2415 − 2760 ) × 10 8 m 3 , for the Gauss model, the best fit with historical production is URR= 2760 × 10 8 m 3 . And from the figure, the Hubbert model fits better than the Gauss model. After comparing the fitted graphs, this paper also calculates the correlation coefficient to determine the accuracy of the model from the perspective of data. The calculation formula is: r = n ∑ i = 1 n x i y i − ∑ i = 1 n x i ∑ i = 1 n y i n ∑ i = 1 n x i 2 − ( ∑ i = 1 n x i ) 2 n ∑ i = 1 n y i 2 − ( ∑ i = 1 n y i ) 2 . (17)

In the formula, r is the correlation coefficient, n is the number of data points, and x i and y i are the coordinate values of each point.

The correlation coefficient method can directly reflect the correlation degree between two groups of variables, which can be used as one of the judgment criteria for the accuracy of prediction results. The correlation coefficients of the production prediction results of the two models under different URRs are shown in Table 2. It can be seen from the table that the correlation coefficients of the two models are very close to 1, and the corresponding production prediction results are very accurate. However, the correlation coefficients of the predicted results of the Hubbert model under each URR are higher than those of the Gauss model. The higher the correlation coefficient is, the more accurate and representative the prediction result is. Therefore, the prediction data of the Hubbert model is selected as the prediction result of tight gas production. In addition, since the Hubbert model has the highest fitting degree under URR= ( 2415 − 2760 ) × 10 8 m 3 , the final peak range should also correspond to it.

The low degree of exploration and development of tight gas leads to many unknown factors. Therefore, the ultimate recoverable reserves URR is introduced into the study of production change trend as a main factor. The final range of peak production was calculated by defining the value range of URR, and Hubbert and Gauss model were established to realize the calculation. The accuracy of the results was judged by comparing the results of the correlation coefficients and referring to the degree of fitting of the predictions. Establishing a prediction model for the growth trend of tight gas production under different proven rates can provide a theoretical basis for tight gas exploration and development planning. The production results predicted by the Hubbert model are shown in Table 3.

Preliminary prediction results show that tight gas production in the Sichuan Basin will continue to maintain a rapid growth trend in the next 20 years. According to the fitting degree and correlation analysis between the predicted results and historical production, the production of tight gas will reach the peak range of ( 86 − 106 ) × 10 8 m 3 / a in 2042. At the end of the stable production period, the URR recovery was about 60%.

As described in Section 3.2, the production growth process for future time periods is divided into two stages. That is, the production rising stage (2022–2038) and the production stable stage (2038–2046). Therefore, it is necessary to calculate and analyze the production realization probability for different production growth stages.

In order to clearly describe the influence of URR on production prediction results, the variation trend of peak production and model parameters under different URRs was first analyzed and studied. As shown in Figure 8, with the increase of URR, the peak production Q m and the model parameter b both keep increasing trend on the whole. However, in a certain area where the URR does not change significantly, Q m and b show a negative correlation trend, that is, when Q m increases, b decreases. This law corresponds to Eq. 4. Then zoom in on the prediction results of the two production growth stages in Figures 6A, and change the abscissa of these figures from year to URR. The URR-production prediction results of each stage and year can be obtained (Figure 9).

It can be seen from Figure 9 that the change trend of annual production with URR is consistent with the change trend of peak production in Figure 8, regardless of whether it is a rising stage or a stable stage. This law corresponds to Eq. 4. Therefore, it can be considered that with the increase of URR, the annual production shows a step-like growth trend. At the same time, with the increase of URR, the interval between the production curves of different years is larger, which shows that the increase of URR also leads to the increase of annual production growth rate. This also explains that in Figure 6, the larger the URR, the greater the slope of the production prediction curve, and the greater the production increase.

To calculate the different production realization probabilities for each year of the two production growth stages, the Monte Carlo method described in Section 2.2.1 was applied. Take 2030 as an example to introduce the probability analysis process of tight gas production. Taking the Hubbert production calculation equation of Eq. 5 as the mathematical model of probability simulation, URR is the main independent variable affecting production. Since the value of URR is obtained by uniform sampling (Section 3.2). Therefore, uniformly distributed URR values are directly drawn at random multiple times. Set the number of URR extractions to 10,000 times. For each URR value extracted, calculate the corresponding Q m and t m values, and substitute them into Eq. 5 together with t = 2030 to obtain the production Q in 2030.

The production of a certain year is repeatedly calculated by cycle. After 10,000 cycles, the distribution probability of the target production can be achieved is obtained, and then the cumulative probability is calculated according to the distribution probability. That is, the probability at the minimum point of annual production is 1, the probability at the next point of production is equal to 1 minus the sum of the distribution probabilities of all previous production, and so on up to the maximum point of annual production. The distribution probability is the average probability of each production value point under the normal distribution, the sum is 1, and the cumulative probability is the realization probability of annual production.

Since the URR is uniformly distributed, the accuracy of the production probability statistics can be guaranteed. Figure 10 shows the production realization probability results for representative years in the 2 phases. As can be seen from Figure 10, for the distribution probability, its changing trend is similar to the normal distribution. That is, in the same year, the closer the production is to the extreme value, the lower the probability, and the closer to the middle value, the higher the probability. For cumulative probability, the lower the production value, the higher the realization probability.

Figure 11A is a graph of the production realization probability of each year during the stable production period. It can be seen from the figure that with the increase of the annual production prediction value, the realization probability gradually decreases. At the same time, the production corresponding to the realization probability of 10%–20% is the ideal production, and the ideal production in the stable production period under the prediction model is stable between ( 80 − 90 ) × 10 8 m 3 / a . Figure 11B shows the prediction of annual production under different probabilities (realization probability 10%–90%), where P90 indicates that the probability is 90%. It can be seen from the figure that the change trend of annual production under different probabilities is basically the same, and the production values of adjacent probabilities are also roughly the same. This indicates that the production change trend under the prediction model is relatively stable, and the prediction results are more accurate.

Table 4 represents the production values that can be achieved in several major years under different probability scenarios. The probability is the cumulative probability, so cumulative probability is also realization probability. In 2042, P 20 = 96.89 × 10 8 m 3 / a means that the probability of production reaching 96.89 × 10 8 m 3 / a is 20%. In this paper, the production of tight gas is calculated in the probability interval of 0%–100%, and the risk quantification of production in different probability intervals. Among them, the production corresponding to P80 is the guaranteed production, the production corresponding to P50 is the average production, and the production corresponding to P20 is the ideal production. The probability values calculated by Monte Carlo method represent the possibility of achieving tight gas production in different periods. This feasibility study has an important guiding role for the planning of tight gas production in the future.

In order to conduct risk quantification research on tight gas production, the risk matrix (Figure 3) in Section 2.2.2 needs to be introduced to evaluate the production risk level. According to the probability calculation method in Section 2.2, the production distribution probability and realization probability curve of each year in the two stages were obtained respectively. According to the distribution probability curve, the mean value μ and the standard deviation s of the annual production are obtained, so as to obtain the dispersion degree C of the annual production according to Eq. 14.

In the rising stage of production, 5 % < C ≤ 10 % . At this time, the annual production corresponding to P > 50 % is a risk level II, the annual production corresponding to 20 % ≤ P < 50 % is a risk level III, and the annual production corresponding to P < 20 % is a risk level IV.

In the stable production stage, 10 % < C ≤ 25 % . At this time, the annual production corresponding to P > 50 % is a risk level III, and the annual production corresponding to P ≤ 50 % is a risk level IV.

Since the dispersion degree is different in different stages, the risk level of production target in different stages and years can be obtained by combining the dispersion degree and realization probability. (Figure 12) (Guo et al., 2022).

It can be known from Figure 11 that the higher the target production is, the lower the corresponding realization probability value is, and the production -probability curve is also in the region with higher risk level. According to the production risk quantification results of different stages in Figure 12, the realization probability and risk level of different productions in each year can be directly obtained. Production risk level indicates how easy it is to achieve the production target, therefore, combined with the realization probability calculation and risk level evaluation, the risk quantification of production can be well studied, and the realization probability of tight gas target production can be calculated, which provides data support for feasibility analysis.

The quantitative research on tight gas production risk is based on the production prediction results. According to the uniformly distributed proved rate, URR of different degrees is calculated, and the probability of achieving the target production in different stages is studied to obtain the possibility of achieving the tight gas production in different years. Combined with the realization probability and the degree of dispersion, the risk grade analysis was used to evaluate the annual production in different stages of production growth. The target risk of tight gas production is comprehensively analyzed.

At present, the exploitation of tight gas in the Sichuan Basin is still in its early stage, and the production prediction and risk quantification of tight gas are also in the shallow stage. In this paper, URR was introduced as the influencing factor of production change trend by using the proved rate as the medium, without more comprehensive quantitative research considering other factors. Risk assessment is based on the results of production prediction and studies the degree of difficulty to achieve the target production. Therefore, the production prediction and risk assessment of tight gas need to carry out more research work, so as to provide a better reliable basis for tight gas production planning.