Confidence limits on kinetic models of primary cracking and implications for the modelling of hydrocarbon generation

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl

Confidence limits on kinetic models of primary cracking and implications for the modelling of hydrocarbon generation Seren B. Nielsen* Bergen University, Department of Chemistry, Allegaten 41, N-5007 Bergen, Norway

and Birger Dahl Norsk Hydro Research Centre, Lars Hillesgate 30, Bergen, Norway

Received 5 July 1990; revised 17 April 1991; accepted 21 April 1991 Kinetic models of primary cracking for use in hydrocarbon modelling are obtained through empirical calibration against laboratory measurements. This paper utilizes least-squares variance analysis to evaluate the uncertainty in the calibrated kinetic parameters and discusses some consequences for the prediction of hydrocarbon generation. It is substantiated that the sensitivity of the average kinetic parameters Ea* and logo(A)* to micropyrolysis errors can be parameterized by errors in pyrolysis Trnax. The absolute and relative errors in Trnax are modelled by normal (Gaussian) random variables and the variance analysis yields the joint normal probability density distribution of variability in --ca* and logo(A)*. It is found that Tma× errors reduce the ability to resolve independently the values of --ca* and logo(A)*: the estimation errors are highly positively correlated and repeated measurements on identical kerogens will produce values of Ea* and Ioge(A)* which tend to fall close to and along the line logo(A)* = EalRTp + constant, where Tp is an average pyrolysis temperature. This correlation is consistent with observations, and hence it is suggested that it is mainly an artifact which originates in the experimental errors of laboratory pyrolysis. When values of Ea* and Ioge(A)* are confined to the region of the confidence ellipsoid, the pyrolysis data are satisfied within their error bounds. However, the variability in predictions of hydrocarbon generation may be large depending on the particular geological temperature history. The uncertainties in the depth to the oil window may be several hundreds of metres, and the timing of peak hydrocarbon generation rate may vary considerably. The inclusion of basin geochemical data in the calibration drastically reduces the ambiguity in Ea* and logo(A)* and emphasizes the need to validate the laboratory calibrated kinetic models against basin observations. Keywords: primary cracking; kinetic modelling; hydrocarbon generation

Introduction It is common practice in hydrocarbon modelling to simulate the primary cracking of kerogen into liquid hydrocarbons by classical first-order reaction kinetics. The degradation of kerogen is described by a sum (the discrete model) or an integral (the continuous model) of independent parallel reactions with the temperature dependence of the respective rate constants governed by the Arrhenius law (Jungten and Klein, 1975; Tissot and Espitalie, 1975). This scheme allows the calculation of the hydrocarbon generation history of a source rock in response to an arbitrary geological temperature history and is widely used in computer programs for exploration purposes (Tissot and Welte, 1984; Tissot et al., 1987; Burnham et al., 1987; Ungerer 1990). The detailed chemical mechanisms of kerogen degradation are not considered in this approach. Instead, the model is empirical and the parameters of *Present address: Department of Earth Sciences, Geophysical Laboratory, The University of Aarhus, Finlandsgade 8, DK---8200 Aarhus N, Denmark

the independent parallel reactions merely provide a framework for calibration against observations, which usually derive from open or closed system pyrolysis experiments performed at elevated temperatures in the laboratory. In principle empirical models should not be used outside their calibration range (Larter, 1989). However, detailed studies have substantiated (Ungerer and Pelet, 1987; Ouigley and Mackenzie, 1988) that first-order kinetic models of primary cracking can simultaneously satisfy both laboratory and basin geochemical data, and that the values of the kimematic parameters in this instance are similar to those derived from laboratory experiments only. Furthermore, as the results of laboratory experiments involving heating rates differing by several orders of magnitude (Saxby et al., 1986; Freund and Kelemen, 1989) are consistent with the first-order reaction kinetic assumption, it seems justified to apply this principle in the translation of laboratory reaction rates to basin conditions. In the process of translation, the reaction rates are extrapolated over 10 to 12 orders of magnitude by the Arrhenuis rate law and the calibrated kinetic

0264-8172/91/040483-10 ©1991 Butterworth-Heinemann Ltd Marine and Petroleum Geology, 1991, Vol 8, November

483

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl parameters, which govern the temperature of each cracking process cannot be resolved dependence. In addition to considerations about the independently: if one parameter is in error then almost chemical validity of a translation over this range, it is certainly the other will also be in error. When the clear that the results in terms of the prediction of kerogen model determined in this manner is translated hydrocarbon generation are sensitive to the precise into sedimentary basin conditions, the induced values of the kinetic parameters. uncertainty in the depth to the top of the oil window The ojective of this paper is the quantification of this typically amounts to about _+400 m for a temperature sensitivity for the micropyrolysis (e.g. Rock-Eval) gradient of 30°C k m - i and a heating rate of 2°C Ma ~. calibration of primary cracking models. To achieve this It turns out that a similar situation holds for the it is assumed that the laboratory calibration of kinetic calibration procedure based on the interpretation of models of primary cracking and the subsequent micropyrolysis reaction rate profiles, which is the focus translation to basin conditions is possible in principle. of this paper. On this basis the problems are addressed of (1) how When the admissible class of kerogen models has well the kinetic parameters can be determined from been quantified by variance analysis, problem (2) is pyrolysis data which have small errors, and (2) how the dealt with by selecting some members of this class and errors in the kinetic parameters affect the results of the performing hydrocarbon modelling in the usual way. extrapolation in terms of the accuracy of predictions of The variability in the results as seen in, for example, hydrocarbon generation in basins. the depth to the oil window or the timing of peak In particular the latter point should be of great generation rates yields an indication of the sensitivity of interest to explorationists, who face the problem of the hydrocarbon generation of a particular geological assessing the hydrocarbon potential of sparsely drilled temperature history to the inevitable errors in the areas. However, since problem (2) arises because of (1) kinetic kerogen model. it is necessary to consider both. In the following section a simple methodology is Problem (1) is dealt with within the well established developed which deals with problem (1) for primary framework of variance analysis. Qualitatively stated, cracking calibration using the R o c k - E v a l $2 peak. this analysis technique quantifies the problems of determining how large a class of different models Sensitivity of kinetic parameters to pyrolysis errors satisfy a certain data set, which has observational errors. For example, assume that a kinetic model of How sensitive are primary cracking kinetic parameters primary cracking has been found, which satisfactorily derived from the R o c k - E v a l $2 peaks to small errors in describes the pyrolysis data of Figure 1A. Intuitively it the data? Before we attempt to answer this question, is clear that certain small perturbations to the kinetic the concept of 'sensitivity' must be given a more explicit model parameters could probably result in an equally formulation, and some basic observations regarding good description of the pyrolysis data within their error this calibration procedure need to be established. bounds, i.e. many solutions (i.e. sets of kinetic parameters) are more or less equally valid. The Calibration by the R o c k - Eval method and the variance analysis technique quantifies which importance of T,,,,,.~ perturbations to the kinetic parameters result in kerogen models that describe the pyrolysis data within As the continuous kerogen model reduces to the their error bounds. discrete model by sampling, it is sufficient in this discussion to deal with the discrete case. Nielsen and Berth (1991) used variance analysis to The calibration (e.g. Tissot et el., 1987; Burnham et study this problem for isothermal hydrous pyrolysis. el.. 1987) of a discrete kerogen model against They found that the allowed perturbations to the micropyrolysis reaction rate profiles is illustrated in average activation energy and the average Arrhenius Figure 1: reaction rate profiles are obtained by open factor of the processes of primary (kerogen ~ bitumen) system pyrolysis at different heating rates (Figure 1A). and secondary (bitumen ---, gas) cracking are strongly The relative weights of a number of parallel reactions positively correlated (9 > 0.99). The correlation means with fixed and equidistantly sampled activation that the activation energies and the Arrhenius factors energies are adjusted (Figure 1B), together with a shared pre-exponential factor, until a satisfactory fit of model to data has been obtained. The temperature of -. ..... Observed maximum generation rate (Tm~×) exhibits a .~, Modelled 0 3 °C/Min characteristic shift with heating rate, which is of special significance. Figure I is an artificial example, where the 7 ~ 0 °C/Min ~ 3 0 0 pyrolysis data are satisfied exactly by the model. This general fitting technique is a refinement of the ~z001 5 i I, w ' technique of heating rate variation analysis (HVA) C g ! m (Redhead, 1962). H V A allows the determination of the 4 ' i El004 activation energy and the pre-exponential factor of a single first-order reaction solely from Tm,x and the shift 200 240 280 of Tm~× with heating rate. Figure 2 illustrates how E~, E (kJ/mole) and logo(A) are determined from the slope and A intercept of a straight line through the points [(l/T, . . . . , T ~ ~ 800 900 1000 200 300 400 500 600 700 loge(edT,.... 2) for different associated values of heating Temperature (°C) rate c~ and T,..... . Also shown in Figure 2 is a point deFigure 1 Calibration of p r i m a r y cracking m o d e l s f r o m reaction fined by a geological heating rate of 2°C M e - 1 and rate profiles obtained at different heating rates. Artificial the associated geological T,..... of 134°C. In essence the example

484

Marine and Petroleum Geology, 1991, Vol 8, November

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl 250

I0

50

~--~_ 15 J

J

~

~

q)

°C/Min

.-~ 240O

4 °C:Min .3 *C/Min

~'~230

\ -4o~

~ 220 N

C/Ms

;z3 210-

-45- E =, R'Slope , .. 50- L°ge(A) = I n t e r e e p t + L n ( E / R ) 1,0

1.15

2.0

( 1/Tmax)*

-

2.'5

3.0

1000

Figure 2 Heating rate variation analysis allows the determination of the activation energy, Ea, and the logarithmic Arrhenius factor, Iog,(A), from the slope and intercept of a straight line through points [1/Tma×, Ioge(a/Tmax2)]

o "~ 2 0 0 .

A SAM

<~190-

[] IFP

180

kinetic calibration problem is one of extending the straight line through the three laboratory heating rate points into the regime of geological heating rates. The refinement of the calibration technique to include the entire production rate profile and multiple parallel reactions has not decreased the importance of Tmax and the shift of Tm~x with heating rate to the determination of kinetic parameters. This is apparent from the fact that application of H V A to Tmax values derived from general kerogens with multiple parallel reactions yields good estimates of the average kinetic parameters (Braun and Burnham, 1987). To illustrate this the following procedure was undertaken: Tmax values for the discrete kerogen models available to this study were determined by heating (computer stimulation at rates of 0.3, 4 and 50°C min -1. The calculated Tmax values were then used to determine the kinetic parameters E, and log~(A) of a single first-order reaction by the technique of Figure 2. Also calculated were the average kinetic parameters Ea* and logo(A)* of the measured discrete kerogen models according to the formulae E.* = ]~kWkEak and log~(A)* = Zkwkloge(Ak), where w, is the relative weight of the kth reaction. The results are displayed in Figure 3, which shows the ability of Ea to predict E,* (Figure 3A and the ability of loge(A) to predict log~(A)* (Figure 3B). It is apparent how E a and log~(A), although based solely on the information represented by Tmax and the shift of Tmax with heating rate, are very capable of explaining the variation present in E.* and log~(A)*. This accentuates the point that Tm~x and the shift of Trn~x with heating rate represent essential information regarding the determination of kerogen kinetic parameters.

Analysis of uncertainty As Tma x is important to the determination of kinetic parameters, it is clear that small errors in Tmax are also important. The question posed is what size of perturbations in the average kinetic parameters are allowed without affecting the Tmax values beyond their acceptable error limits. To arrive at a quantitative formulation we proceed as follows. Define the kerogen model Kerogen to be any discrete or continuous kinetic model of primary cracking, which has been determined by calibration

i

1

180

19o

i

2oo

i

i

21o

Activation

24o

Energy

i

A

i

2~o'24o

250

(ka/mole)

/

38

36

"~'34 v

~-32 [] o

A

[] IFP

28

26

SAM

B 26

'

218

'

310

'

3f2

'

314

'

316

'

38

Loge(A ) ( l / s ) Figure 3 The ability of (A) Ea and (B) Ioge(A) derived from heating rate variation analysis to predict values of E~* and Iog,(A)* obtained by averaging discrete kinetic kerogen models. SAM are data supplied by Alan Samoun and IFP are data supplied by Norsk Hydro

against R o c k - E v a l data or by other means. Let the quantities AE and Aloge(A) be perturbations of the kinetic parameters of Kerogen. The effect of AE is to shift all the activation energies of Kerogen by the amount AE, and the effect of Aloge(A) is to shift the value of logo(A) of Kerogen by the amount Alog~(A). Define further the operation TMAX[Kerogen, AE, Alog~(A), et] to be the process of calculating pyrolysis /'max of Kerogen at heating rate 0¢ with perturbations AE and Alog¢(A) to the kinetic parameters of Kerogen. Let T m a x l , Tmax2 , and Tmax3 be the T m a x values obtained by heating Kerogen on the computer at rates o:1, 0c2and 0~3,respectively, with no perturbations to the kinetic parameters of Kerogen [i.e. AE = Aloge(A) = 0]. The allowed perturbations AE and AIoge(A) of the kinetic parameters of Kerogen due to errors in Tmax

Marine and Petroleum Geology, 1991, Vol 8, November

485

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and Bo Dahl hence formally are constrained by the system of around the true absolute value of Tmax. The accuracy equations: errors at two different heating rates are not uncorrelated. This is the statistical formulation of the TMAX[Kerogen, AE, Alogc(A), ~l] = statement that when Tm,x is biased by, say, 5°C at a Tmax I -4- Ep I + t~a1 heating rate of 1°C min -~. then it is probably biased by TMAX[Kerogen, AE, Alog~(A), 0~2] = a comparable amount at comparable heating rates. Tmax2 + Ep2 -l- Ea2 The investigation of Equation (1) is performed by TMAX[Kerogen, AE, Alog~(A), 0~3] = least-squares variance analysis. In this approach (e.g. Tarantola, 1987) this requires that the perturbations in Tmax3 + Ep3 -F Ea3 the kinetic parameters and the errors in Tm,,x can be (1) understood as normal (Gaussian) random variables and where Ep, and e,i, i = 1, 2, and 3, represent the that Equation (1) is reasonably linear, which is true precision and accuracy errors of the respective Tin,× (Nielsen, 1991). values. Equation (1) is a formal statement of the These models of errors in Tm~,x appear in the requirement that the allowed perturbations in the covariance matrices of the variance analysis formalism. kinetic parameters of Kerogen are those which keep the In the examples given here, the accuracy errors in T,n~x Tm~x values within their acceptable error bounds. are represented by a covariance matrix with elements In the special case when Kerogen contains the given by efexp[-loge(OCi/%)2/y, where y is a factor parameters of one single first-order reaction, Equation which regulates the degree of correlation of the (1) is the formal statement of the H V A technique accuracy error between different heating rates c~i and illustrated in Figure 2, with consideration of the errors %, and e, is the standard deviation of the accuracy in Tm~x. Hence, Equation (1) can be understood as a error. This type of covariance function ensures that generalization of the classical H V A to kerogens with realizations of the absolute error in Tma×vary smoothly multiple parallel reactions. Note that Equation (1) with the heating rate, as would be expected to hold for signals no intention to determine the distribution of R o c k - E v a l pyrolysis. The precision errors in Tma× are activation energies of Kerogen, which is assumed to be represented by a covariance matrix with zero offknown. diagonal elements and a diagonal consisting of the square of the precision error, ep, in each Tm,x. The Errors in Tmax combined effect of accuracy and precision errors is Accuracy errors are rooted in the absolute temperature obtained by adding the two matrices. Further details of calibration of the pyrolysis equipment, whereas the analysis may be found in Nielsen and Barth (1991) precision errors are due to the finite experimental and Nielsen (1991). reproducibility. Precision errors have been considered in detail by Results of the variance analysis Espitalie (1986). He finds that operating conditions, The analysis results relevant to this study are the five mineral matrix effects and heavy oil accumulations are parameters which describe the joint behaviour of the probably responsible for the -4-5oc dispersion of Tmax normal random variables AE and Aloge(A), i.e. which which is frequently observed during the pyrolysis of describe the allowed variability in the average kinetic samples from the same depth interval and of the same parameters. The parameters are the standard deviation source rock quality. However, he states that under Oe of AE, the standard deviation oa of Alog~(A) and the normal operating conditions the reproducibility of Tmax coefficient of correlation between the two, [- The measured on the IFP reference rock is + l°C. In the remaining two parameters, the expected values ft c and present approach, it is assumed that the precision ft, of AE and Alog~(A), respectively, are zero in the errors in Tm,x may be modelled by uncorrelated zero present application of Equation (1) because Tin,×1, mean normal (Gaussian) variables. The validity of this Tm~x2 and Tm~x3 have been calculated from a given approach may be tested by performing many identical kerogen model. The standard deviations can be pyrolysis experiments and generating a histogram of considered as the amplitude of variability in the the observed Tma~ values. The assumption is good parameters. However, only approximately two times when a normal probability density function (p.d.f.) out of three can we expect the value of AE to lie in the approximately describes the shape of this histogram. interval from - % to + % . p is a number between - 1 Accuracy errors are of an elusive nature. For and +1 which describes the degree of relatedness example, if a calibration error is known to exist, the between values of AE and Aloge(A). When p = 0 the calibration is immediately modified to eliminate the two variables take values completely independently of error. However, following the recalibration it is each other. When P = +1 (or - 1 ) there is perfect probably not wise to state that calibration errors have coupling between the two parameters. Further now been completely eliminated. It is less obvious that interpretations are given in the following examples. the accuracy errors in Tm~x may be modelled by zero For the sake of illustration, the variance analysis is mean normal variables. For a particular py'rolysis performed for two different measurement procedures. equipment the accuracy error is systematic, i.e. ideally Laboratory A uses three different heating rates (0.3, 4 it does not change until the temperature calibration is and 50°C rain-l), spanning more than two orders of changed. However, when comparing different pyrolysis magnitude, and laboratory B uses four different equipment this error is also random. The requirement heating rates (5, 10, 20 and 50°C rain-t), spanning one is that individual pieces of equipment exhibit systematic order of magnitude. absolute temperature errors in Tmax which, if they To illustrate the effect of calibration errors, the could be revealed for many pieces of equipment analysis is performed for three classes of accuracy performing identical pyrolysis experiments, would generate a histogram shaped like a normal p.d.f. errors given by the values e~, = 0, 3 or 5°C and with

486

Marine and Petroleum Geology, 1991, Vol 8, November

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl 4O constant ~, -- 20. It is emphasized that each class of 66.7~ C o n f i d e n c e Level errors has the expected value of zero, representing the belief that each set of equipment is calibrated to the best of our knowledge. Nevertheless there will be 36 calibration errors, differing from one set of equipment to the other, but with an overall standard deviation of 0, 3 or 5°C. The value of ,/ = 20 ensures that the realization of the calibration error for a given set of equipment varies smoothly with heating rate. If this o=30-1 was not so (e.g. ,/ = 0), the calibration error would have a very severe effect on the possibility of measuring 28 kinetic parameters.

t

-.-t-/.'_

26

Joint behaviour of A E and Aloge(A) To illlustrate the type of variability in AE and Aloge(A) which always results from the variance analysis, the analysis is performed with the above specified classes of accuracy errors and for precision errors of ep = 0, 1, 2 and 4°C. The parameters of Kerogen are those of Figure lB. However, the results in terms of sensitivity depend mainly on the measurement conditions and not much on the particular kerogen investigated. The parameter values obtained are listed in Table 1. The values of/zc, /z,, oc, o, and p define the joint normal p.d.f, of AE and Alog~(A). We define the

A

24 160

tso

26o

220

2~o

2~o

280

Activation Energy (k J / m o l e ) 4O

66.7~ C o n f i d e n c ; lLiv21 38-

36"~34~32

-

Table 1 Results of variance analysis Ep

Oe

Oa

(°C)

(kJ mol 1)

(s-l)

Excluding accuracy e r r o r s scale: Laboratory A

Laboratory B

(E a = O ° C ) "

0~30 -

p

28-

dispersion on a relative 26-

0 1 2 4

0.0 3.62 7.23 14.47

0.0 0.636 1.271 2.542

-0.9987 0.9987 0.9987

0 1 2 4

0.0 7.13 14.25 28.50

0.0 1.220 2.440 4.878

-0.9998 0.9998 0.9998

Including accuracyerrors ~a = 3 ° C a n d y = 20): dispersion on an absolute sca~: Laboratory A

0 1 2 4

9.33 10.01 11.82 17.23

1.620 1.745 2.065 3.019

0.9969 0.9969 0.9973 0.9980

Laboratory B

0 1 2 4

7.07 13.12 18.02 30.55

1.182 2.228 3.071 5.222

0.9937 0.9977 0.9986 0.9993

Includingaccuracyerrors ~a =5~Candy = 2 ~ : d i s p e r s ~ n on an absolu~ sca~: Laboratory A

0 1 2 4

15.54 15.97 17.17 21.26

2.700 2.778 2.992 3.718

0.9969 0.9968 0.9969 0.9974

Laboratory B

0 1 2 4

11.78 19.69 23.25 33.90

1.968 3.334 3.951 5.784

0.9949 0.9973 0.9979 0.9989

The + one standard deviation kinetic parameters o f Kerogen are given by Low values: ~t E a , = E a , - oekJ mol 1

J High

A = exp[Ioge(A) - p%] s -~

values: ~ Eak = Ea, + oe kJ m o l - '

J

A = exp[Ioge(A) + p o a ] s -1

24 160

z~o

2~o

220

240

2~o

Activation Energy (kJ//mole)

B

280

Figure 4 Confidence ellipsoids for the laboratories A and B measurement procedures, en = 2°C, ea = 3°C and y = 20. The

examples are listed in Table 1

associated confidence ellipsoid to be the iso-probability contour bounding the area associated with probability 66.7%. It is hence analogous to the well known one-dimensional normal confidence interval. These confidence regions are shown in Figure 4 for the laboratories A and B cases of Table 1 with ep = 2°C, I~a = 3°C and y = 20. Also indicated in Figure 4 are the locations in the E - l o g e ( A ) plane of the observed kerogen model and the two extremes obtained by shifting the activation parameters of Kerogen to the pointed ends of the confidence ellipsoid. Consider, for instance, Figure 4B. The interpretation of the confidence ellipsoid is that, for a given pyrolysis experiment, the true average kinetic parameters of the kerogen will, with a probability of 66.7%, be located within the ellipsoid. Alternatively, if many identical sets of kerogen samples were measured by laboratory B on equipment with a precision error of 2°C and with the class of accuracy errors given by e, = 3°C and y = 20, then the scatter in the results in terms of average kinetic parameters would occupy a region in the E-loge(A) plane approximately given by the confidence ellipsoid, but with one-third of the points falling outside, mainly in the elongated direction. It is apparent that the confidence regions of Figure 4,4 and B are very elongate. This is the expression of the high degree of correlation between AE and

Marine and Petroleum Geology, 1991, Vol 8, November

487

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl correlation, this effect ought to be visible in ~10 experimental results. Perturbed =d 9 Observed It has long been realized that average values of E, 0 3 C/Mill and logo(A) of kinetic kerogen models plot close to and 8 ]/~ 4 "C/Min 35 , along a straight line (Wood, 1988). Ungerer (1990) ~,~ ] 66.7 t~ Confidence Level ? presented a cross-plot of peak activation energy and 2 ~ 6 associated loge(A) value for a number of kinetic 5 kerogen models derived from micropyrolysis reaction O rate profiles from different laboratories (his Figure 8). The points plot close to and along a straight line, ~ 3 indicating that the two variables are linearly correlated. At first sight this correhition could be the expression of c~ 1 some fundamental chemical kinetic conditions of ~ 0 sedimentary organic matter, which then would require 200 300 400 500 600 700 800 900 lO00 an explanation. Temperature (~C) However, from knowledge about how the interpretational procedure (probably the OPTIM computer program) reacts to small shifts in the reaction ~10[ .~ : Perturbed rate profiles, Ungerer (1990) further notes that a major Observed part of this correlation is probably due to ,,03 C/Min 8 micropyrolysis errors. The present variance analysis A 4 °C/Min 35 1 . . . . . . . . . 66 7~; Confidence Level 7 supports this assessment. Figure 6 shows points of E,* /,50 C,,'Min i J a m and logo(A)* for the discrete kerogen models available to this study, together with the confidence ellipsoids of Figure 4. Also shown is the correlation line of Ungerer (1990). It is apparent that a major part of the correlation could have been induced by pyrolysis errors rather than by real differences in the kinetic parameters :~

2i

~o~

~00

/Vii 1 :

2

,10 :

300

400

500

600

Temperature

700

800

900

1{}00

i

('C)

/J

(orr Line of Ltl~eter ([990)

;~I

/ /

Figure 5 Variability in laboratory production rate profiles associated with (A) Ep = 1°C and (B) ~p = 2°C. The broken lines correspond to perturbations 1, 2, 3 and 4 on the confidence ellipsoid. The solid lines represent the observed kerogen

',36 i 1

34

Alog~(A) as emphasized by the values of p close to unity (Table 1). Repeated identical pyrolysis experiments on the same equipment reveal the precision error. The accuracy error only appears when the results from different sets of equipment are compared. Figure 5 shows how the effect of precision errors given by ep -I°C (Figure 5A) and ep = 2°C (Figure 5B) would show up in repeated identical pyrolysis measurements. The average kinetic parameters of the kerogen have been shifted in the E-loge(A) plane according to the points numbered 1-4 on the confidence ellipsoids, and the production rate profiles according to the laboratory A heating rates have been calculated. It is apparent that the variability of the average kinetic parameters on the confidence ellipsoid induces only minor changes in the pyrolysis response of the kerogen. This supports the assertion of the present variance analysis model that the quantification of micropyrolysis errors by errors in Tm,~ is relevant and useful.

o

A SA. 28

/\// i

E6

I

./~"

[J I F P

/"

a4

,

1

Activation

}':f]{:!r'~y ( k J / m o l e )

40 :38

Corr Lln~ of /ngerez (1990}

•

{36

f'

- . :34 ~32 zA O~£30 -2

SAM

Marine and Petroleum

Geology,

IFP

/

E6

Experimental evidence of the correlation between Ea and log~(A)

488

/// /

2g

The main result of the previous section is that pyrolysis errors, as they show up in absolute and relative errors in Tm~,~, tend to induce a correlation between the average values of E~, and log~(A) of kinetic kerogen models. Considering the size of the induced

A!

i

160

160

1B0 200 AcLivation

7J:~{/ F2nergy

240

26'0 (k J/mole)

B/

280

Figure 6 Correlation between E a and IOge(A) seen by discrete kinetic models of primary cracking. Also shown are the confidence ellipsoids of Figure 4 and the correlation line of Ungerer (1990)

1 9 9 1 , V o l 8, N o v e m b e r

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl 1 i set of pyrolysis equipment with unchanged temperature Incl. calib, errors I Exel. calib, errors calibration. In this instance the scatter in results 66.7% Conf. L e v e l 66.7% Conf. Level depends solely on the ability of the equipment to 22 reproduce measurements. However, the results are biased on an absolute scale by an amount determined ~3 '-~'--. by the characteristics of the particular accuracy error• ""::-~ . . . . . . . The effect of increasing the precision error is not the same in the two instances because the induced uncertainty in the average kinetic parameters grows Q9 approximately propomonally t o Ea[1 -F (ep/l~a) ] . T h i s 5 expression does not change very much when ep is small compared with ea. k~

•

0

0J

7

71

o.o 'o.~ 'o.'4 'o.'~ 'o.'8 Transformation

1.o

Ratio

B

o.o 'o.'z 'o.'4 'o.'8 'o.'8 Transformation

1.o

Ratio

Figure 7 Variability of the t r a n s f o r m a t i o n ratio as a funct i on of depth f o r the l a b o r a t o r y A procedure. The solid line corresponds to o b s e r v e d kerogen. Broken lines correspond to: En = 0, 1, 2 and 4°C; ~ = 3°C and ¥ = 20 w h e n calibration errors are included; and ~, = 0°C w h e n calibration errors are excluded. The e x a m p l e s are given in Table 1

of kerogens. This characteristic condition should be considered when average kinetic parameters are used for the classification of kerogens; differences between kerogens in the direction perpendicular to the strike direction of the confidence ellipsoid are much more likely to be significant than those in the direction parallel to the strike.

Consequences for the prediction of hydrocarbon generation In this section we investigate some aspects of how the derived statistical uncertainty in kinetic models of primary cracking influences the ability to predict hydrocarbon generation•

Depth to the oil window Figure 7A and B show the transformation ratios as a function of depth for calibration by the procedure of laboratory A, including and excluding accuracy errors, respectively• The oil window is largely the depth interval where the transformation ratio changes from 0 to 1. The kerogen models are obtained by perturbation of the kinetic parameters of the kerogen of Figure 1B according to the ep = 0, 1, 2 and 4°C parameters of Table 1. The geological temperature history is heating at a rate of 2°C Ma -~, the temperature gradient is 35°C km-~ and the surface temperature is 0°C. The regions in Figure 7A between symmetrically placed broken lines indicate, for each value of ep, the variability in modelling results which would arise if the kinetic parameters of identical kerogens were measured by the procedure of laboratory A on different pyrolysis equipment characterized by the above-specified class of accuracy errors• The regions are the 66•7% confidence interval for the depth to the oil window. The regions in Figure 7B between symmetrically placed broken lines indicate, for each value of ep, the variability in modelling results which would arise if the kinetic parameters of identical kerogens were measured by the procedure of laboratory A on a single

21,5

.

Generation history The possible hydrocarbon generation patterns associated with more general geological temperature histories are usually investigated by Monte Carlo simulations• The input parameters of the modelling problem such as thermal conductivities, sediment compaction properties, surface temperatures, background heat flow, thicknesses removed during erosion episodes and the kinetic parameters are varied within 'reasonable limits', and a forward calculation is performed for each new set of parameters• The different results obtained for many parameter sets, together with the likelihood of each parameter set, allow a statistical assessment of the prospectivity of a region or a particular play. The statistical error bounds derived in this paper offer a straightforward definition of 'reasonable limits' which should be used for the kinetic parameters in Monte Carlo simulations• Observing these error bounds ensures that kerogen models still satisfy the important constraint given by pyrolysis Tmax- More arbitrary definitions of variability in the kinetic parameters could easily produce a larger variability in the outcome of Monte Carlo simulations than is strictly necessary• To further illustrate the character of the dispersion in modelling results, which is induced by the suggested variability of the kinetic parameters within the confidence ellipsoids, the transformation ratios and generation rate histories were calculated for the observed kerogen and for four alternative kerogens obtained by perturbing the kinetic parameters of the observed kerogen according to the parameters of the confidence ellipsoid. The perturbations are numbered 1-4. Perturbations 1 and 2 are located at the pointed ends of the confidence ellipsoid as indicated in Figure 8B. Perturbations 3 and 4 are located on the intersections with the confidence ellipsoid of a vertical line through the average activation energy• The perturbations hence sample extreme cases on the confidence ellipsoid• The example geological temperature history is characteristic of a possible lower Upper Jurassic source rock in an area of the North Sea Central Trough. During Kimmeridgian time large thicknesses of the Farsund Formation (time equivalent to the Kimmeridge Clay Formation) were deposited and caused the initial rapid temperature rise, which is visible in Figure 8A. During Cretaceous time the temperature rose slowly, initially mainly due to fading of the transient thermal regime induced by the rapid Upper Jurassic sedimentation, later mainly due to deeper burial beneath Upper Cretaceous sediments•

Marine and Petroleum Geology, 1991, Vol 8, November

489

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl 200 . . . . I 0 During the Palaeogene the rate of sedimentation was ........... P e r t u r b e d Kerogens ~ J relatively low until Late Palaeogene time, when the { - Measured Kerogen .~ d sedimentation rate started to increase gradually. At ] . . . . . Temperature History ~/ /,;~L-0.u o present the rate of sedimentation and the heating rate // are at their maximum since Upper Jurassic time, and the modelled present day source rock temperature is believed to be around 150°C. Figure 8A and B show the transformation ratio and 1 . /./ / 04 generation rate histories for ~p = I°C, respectively. The t ~ " /./ .'" measurement procedure is that of laboratory A with ~, = 0°C, i.e. the variability in modelling results induced by calibration errors is not considered. It is apparent how the transformation ratio histories of the perturbed 160 120 80 40 0 Age (Ma) kerogens always bracket (although not symmetrically) the transformation ratio history of the measured kerogen. This rule does not, however, apply to the 35 66.7Z Confidence . . . . . . . . . . . . . 7 004 generation rate histories, which generally show much greater complexity. It seems to be usual that perturbations 3 and 4 are always very close to the evolutions of the measured kerogen, whereas perturbations 1 and 2 may exhibit large deviations. This can be explained by consideration of Figure 2 and i ! ' observing that perturbations 3 and 4 involve a vertical shift of the line without changing the slope (the average .<,,,.o,.> ., /I o activation energy is constant), whereas perturbations 1 and 2 involve tilting of the line. From Figure 2 it is clear that, for given error bars on the Tm,~ values, the tilting perturbations influence the predictions in the regime of 160 120 80 40 0 geological heating rates more than the constant slope Age (Ma) shift perturbations. In Figure 9 a calibration error of ~, = 3°C and ,! = 20 Figure 9 Evolution of (A) transformation ratio and (B) normalized production rate for a measured kerogen and for the has been added to the Ep = l°C precision error, and the kerogen perturbations 1-4 indicated on the confidence ellipsoid. Ep = 1°C, ~a = 3°C and 7 = 20. Compare with Figure 5,4 ]

200

10

....

Perturbed

Kerogens

j

Measured Kerogen

i

~

Uº5O'" . . . . Temperature History-y

~"

o

/

4

0;

.

1/ j /

~

160

.

.

~

.

{30

40

Age (Ma) .

.

.

.

.

.

.

33

o~

o~

i:04

"~

A0o

r

r

120

.

/

0

.

o,o4

35]'66.7~ Confidence Level

:

2

"~"

0.03

0.02

i 2 9 ~ 200

- -~ ......... 220

g

~

4

1

/

'~

~ 0.01

g (ka/rnole)

~-~,/J

j

,"7/

120

80

40

b~

~ B

160

o

8

!

l:~

.

~

0.00 0

Age (Ma)

Figure 8 Evolution of (A) transformation ratio and (B) normalized production rate for a measured kerogen and for the kerogen perturbations 1-4 indicated on the confidence ellipsoid. ~p = 1°C and ~ = 0°C. Compare with Figure 5A 490

scatter in the modelling results hence indicates the variability on an absolute scale. Whereas from Figure 9,4 there is little doubt that at present the source rock has released most of its hydrocarbon potential, it is much more uncertain when this happened. Figure 9B shows that perturbation 1 yields relatively high generation rates during the time interval from about 80 to 50 Ma whereas perturbation 2 yields relatively high generation rates from about 30 Ma until the present. In this instance, the potential trap for the hydrocarbons existed at least since the Early Cretaceous. The difference between the generation histories hence suggests the possibility of two extremes: (1) early hydrocarbon generation with a long residence time of the hydrocarbons in the trap and, therefore, the risk of leakage, or (2) late hydrocarbon generation and vigorous generation at present possibly enhancing the likelihood of a plentiful reservoir, if the permeability still exists. Figure 10 shows an example in which basin geochemical information has been included in the determination of kerogen kinetic parameters. The laboratory A measurement parameters are ep = 2°C, el, = 5°C and y = 20, i.e. reflecting a relatively large, but possibly realistic, uncertainty in pyrolysis measurements. The basin geochemical information consists of the geological temperature (134 _+4°C) of maximum generation rate at a geological heating rate of 2°C M a - '. This kind of information can be obtained in rare instances when the transformation ratio profile of a homogeneous source rock has been mapped as a function of depth and when the geological temperature

M a r i n e a n d P e t r o l e u m G e o l o g y , 1991, V o l 8, N o v e m b e r

Confidence limits on kinetic models of primary cracking: S. B. Nielsen and B. Dahl 200 10 laboratory and basin geochemical data, and its - Perturbed Kerogens / ) numerical handling is simple and efficient. However, M e a s u r e d Kerogen ~ only rarely are useful basin geochemical data available . . . . . ... . . . Tempera-ture History -o.8 o o~ 150 ~ for calibration of the kinetic models due to the sparsity ~ . - ~ ~ / o8 o= and heterogeneity of the sample material obtained for 10o ~ different levels of thermal maturity in boreholes. 04 E It is therefore common to base the calibration solely , ~ on laboratory measurements. This paper has discussed the fundamental problem of parameter sensitivity in 50 ,,"3 -o.z ~ the calibration of primary cracking models by micropyrolysis reaction rate profiles, and some aspects 0 A 00 of how these problems influence the ability to predict 180 1~0 8'0 4'0 hydrocarbon generation. Age (Ma) The basic reasoning was as follows: micropyrolysis measurements have small errors. Consequently, the 0o4 kinetic parameters cannot be determined perfectly a5 88~zconfid.... L~v,~ from such measurements. Translation of the laboratory ~ ,~3a

t

iiiiii/~

f

.~3~ / ~

~

,z~ o

16o

z~o

do

2,

8'0 Age (Ma)

4'0

0.0a ~ ~ a~ i-0.0z~ o ~ \ 0.ol ~ o z B 0.00 o

Figure 10 Evolution of (A) transformation ratio and (B) normalized production rate for a measured kerogen and for the

kerogen perturbations 1-4 indicated on the confidence ellipsoid. Ep 2°C, ea = 5°C and y = 20. The kinetic parameters have been further constrained by inclusion of knowledge about the peak generation rate at a heating rate of 2°C Ma-;. Compare with FigureSB =

history is well constrained. It is apparent from Figure 10 that this kind of additional information greatly reduces the variability in modelling results and makes useful pyrolysis measurements which would otherwise yield useless kinetic parameters, Figures 8B, 9B and 10B further show examples of the fact that the measured kerogen model, as compared to the symmetrical perturbations 1 and 2, may define intermediate generation rates (140-60 Ma), as well as high (60-30 Ma) or low generation rates (20-10 Ma). This is a non-linear effect, which implies that there is no guarantee that hydrocarbon modelling based exclusively on the measured kerogen kinetic parameters yield the statistically most likely generation rate history, a quantity which would probably be useful in hydrocarbon prospecting. This would only be the outcome if the evolution of the measured kerogen was always bracketed symmetrically by the evolution of the perturbations. This emphasizes that the statistically most likely generation rate or transformation ratio histories, as well as their variances and other statistical properties, should be obtained by statistical simulations, Conclusions The first-order reaction kinetic model of primary cracking is a sensible basis for hydrocarbon modelling, It has been demonstrated to b e able to satisfy both

reaction rates to the thermal regime characteristic of sedimentary basins makes use of the kinetic parameters: they determine the temperature dependence of the reaction rates. As the translation occurs over a large range of heating rates and temperatures, the predictions of hydrocarbon modelling have the potential of becoming very uncertain. This paper has demonstrated that, with respect to the sensitivity of kinetic parameters, micropyrolysis errors can be quantified by errors in pyrolysis Tmax. This strategy works well because Tmax and the shift of Tmax with heating rate represent essential information regarding the determination of the average kinetic parameters of primary cracking models. It was further demonstrated that small random errors in Tmax induce a very characteristic type of error in the kinetic model parameters: the errors in the average values of loge(A) and Ea are highly positively correlated with a coefficient of correlation typically in excess of 0.99. This means that repeated measurements on identical kerogens will tend to produce results which fall close to and along a straight line in a diagram of average loge(A) and Ea. It is believed that this result explains a large part of the observed correlation between average log~(A) and Ea. The line is given by log~(A) = Ea/RT p + constant, where Tp is an average pyrolysis temperature. More specifically, when the errors in Tmax are normally distributed, the variability in average log¢(A) and Ea follows in a good approximation a joint normal probability density distribution. The parameters of the distribution were found by applying least-squares variance analysis to a generalization of the classical heating rate variation analysis technique. This analysis does not resolve the distribution of activation energies, which is assumed to be known. However, for a given distribution it imposes a constraint on the location of this distribution in the average log~(A) and Ea plane by means of the, say, 66.7% confidence ellipsoid. The established behaviour of the variability in average loge(A) and Ea defines a class of kerogen models which satisfy the constraints of pyrolysis Tmax and which should be useful in Monte Carlo simulations of hydrocarbon generation. Members of any other class of kerogen models are likely to violate the pyrolysis constraint and thereby possibly give rise to greater variability in modelling results than is strictly necessary. When identical kerogens are calibrated by one set of

M a r i n e a n d P e t r o l e u m G e o l o g y , 1991, V o l 8, N o v e m b e r

491

C o n f i d e n c e l i m i t s o n k i n e ti c m o d e l s o f p r i m a r y c r a c k i n g : S. B. N i e l s e n a n d B. D a h l

pyrolysis equipment with unchanged calibration, the variability in hydrocarbon modelling results depends, in this analysis, solely on the ability of the equipment to reproduce measurements. However, the results are biased by an amount determined by the characteristics of the particular calibration error. When identical kerogens are measured by different sets of equipment, the individual calibration errors contribute to the dispersion in the results. How severely depends on the characteristics of these errors. As hydrocarbon modelling aims to predict the absolute amounts of hydrocarbons generated, it is the latter case, referring to variability on an absolute scale, which is relevant. If the subject is to evaluate the prospectivity of a play relative to one which is known, the former relative scale may possibly be relevant. The induced variability in modelling results was investigated in both situations by selecting members of the determined classes of kerogen models and performing hydrocarbon modelling. The example geological temperature history illustrated the evolution of a possible lower Upper Jurassic source rock in a region of the North Sea Central Trough. The examples showed that on the scale of relative variability, reproducibility (precision) errors in Tm,x of about 1°C and less, together with more than two orders of magnitude difference in pyrolysis heating rates, yield a satisfactory resolution of the hydrocarbon generation history as seen through the evolution of the transformation ratio and the normalized generation rate. The further inclusion of a correlated calibration (accuracy) error of standard deviation 3°C considerably increases the variability in modelling results. A possible way of evading this problem is the inclusion of basin geochemical information in the kerogen calibration procedure. In this formulation this reduces to supplying the temperature of peak hydrocarbon generation rate and its uncertainty at one geological heating rate. Such additional information, which is rarely available, greatly reduces the ambiguity in kerogen models and thus also the dispersion in the hydrocarbon modelling results. For the example studied it was shown that the supplementary information furnished by 'peak generation rate at 134 +4°C at a heating rate of 2°C Ma -I (corresponding to a depth uncertainty of +0.114 km for a temperature gradient of 35°C k m - ' ) made the pyrolysis results obtained with a precision error of 2°C and a correlated accuracy error of standard deviation 5°C useful. This level of pyrolysis uncertainty, which may well be realistic, would otherwise yield kinetic models which in most instances are useless. In general, the variance analysis yields the class of primary cracking models consistent with associated values of heating rate and temperature of peak generation rate. If furnishes a straightforward means of evaluating the reliability of measured ,kinetic parameters, when the quality of the pyrolysis measurements in terms of reproducibility and calibration has been assessed. It further allows for combining pyrolysis data with knowledge about the occurrence of peak generation rate under geological conditions, i.e. to calibrate kerogen models to comply with known geological generation conditions. In any instance, the admissible class of kerogen models is

492

parameterized in a simple way which should be useful in general Monte Carlo simulations of hydrocarbon generation.

Acknowledgements This study was undertaken while S. B. Nielsen visited the Department of Chemistry, University of Bergen, Norway, under a Nordic Ministry Council research fellowship in oil geology. Thanks are due to Alan Samoun (Lab Instruments Inc.) and Norsk Hydro, Bergen, for permission to use these kinetic data. Tanja Barth (University of Bergen), Steve Larter (University of Newcastle) and Erik Thomsen (Danish Geological Survey) are thanked for stimulating discussions.

References Braun, R. L. and Burnham, A. K. (1987). Analysis of chemical reaction kinetics using a distribution of activation energies and simpler models Energy Fuels 1, 153 Burnham, A. K., Braun, R. L., Gregg, H. R. and Samoun, A. M. (1987) Comparison of methods for measuring kerogen pyrolysis rates and fitting kerogen kinetic parameters Energy Fuels 1,452-458 Espitalie, J. (1986) Use of Tm~x as a maturation index for different types of organic matter. Comparison with vitrinite reflectance. In: Thermal Modelling in Sedimentary Basins (Ed. J. Burrus), Editions Technip, Paris, pp. 475-496 Freund, H. and Kelemen, S. R. (1989) Low temperature pyrolysis of Green River kerogen Am. Assoc. Petrol GeoL Bull. 73, 1011-1017 Juntgen, H. and Klein, J. (1975) Entstehung von Erdgas aus kohligen Sedimenten Erdol Kohle 2, 65-73 Larter, S. (1989) Chemical models of vitrinite reflectance evolution. GeoL Runds. 78, 349-359 Nielsen, S. B. (1991) Characterization of the ambiguity in kinetic models of primary cracking Proceedings of the NPF Meeting on Basin Modelling, Stavanger, Norway, 13-15 March 1991, in press Nielsen, S. B. and Barth, T. (1991) An application of least-squares. Inverse analysis in kinetic interpretations of hydrous pyrolysis Math. GeoL 23, 565-582 Quigley, T. M. and Mackenzie, A. S. (1988) The temperature of oil and gas formation in the sub-surface Nature (London) 333, 549-552 Redhead, P. A. (1962) Thermal desorption of gases Vacuum 12, 203-211 Saxby, J. D., Bennett, A. J. R., Corcoran, J. F., Lambert, D. E. and Riley, K. W. (1986) Petroleum generation: Simulation over six years of hydrocarbon formation from torbanite and brown coal in subsiding basins Org. Geochem. 9, 69-81 Tarantola, A. (1987) Inverse Problem Theory - - Methods for Data Fitting and Model Parameter Estimation Elsevier, Amsterdam, 613pp Tissot, B. P., Pelet, R. and Ungerer, Ph. (1987)Thermal history of sedimentary basins, maturation indices, and kinetics of oil and gas generation Am. Assoc. Petrol GeoL Bull. 71, 1445-1466 Tissot, B. P. and Welte, D. H. (1984) Petroleum Formation and Occurrence, Springer-Verlag, Berlin, 699pp Tissot, B. P. and Espitalie, J. (1975) L'evolution de la materiere organique des sediments: applications d'une simulation mathematique Rev. Inst, Fr. Petr. 30, 743-777 Ungerer, P. (1990) State of the art of research in kinetic modelling of oil formation and expulsion. In: Proceedings of the 14th International Meeting on Organic Geochemistry, Paris, France, 18-22 September 1989 (Eds B. Durand and F. Behar), Org. Geochem. 16, 1-25 Ungerer, P. and Pelet, R. (1987) Extrapolation of kinetics of oil and gas formation from laboratory experiments to sedimentary basins Nature (London) 327, 52-54 Wood, D. A. (1988) Relationships between thermal maturity indices calculated using the Arrhenius equation and Lopatin method: implications for petroleum exploration Am. Assoc. Petrol GeoL Bull. 72, 115-134

M a r i n e a n d P e t r o l e u m G e o l o g y , 1991, V o l 8, N o v e m b e r