Method for forming history of natural gas accumulation by using carbon isotopes by pyrolysis experiment

Abstract

The present invention provides a method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment. The method includes: obtaining activation energy distribution and a frequency factor of light carbon methane; carrying out carbon isotope kinetics simulation of natural gas in a study area by using a spreadsheet function of Excel to obtain activation energy, a mass fraction and a frequency factor of heavy carbon methane; establishing a burial history and a thermal history of the study area based on geological data; and combining the activation energy distribution and frequency factor of the heavy carbon methane with the burial history and thermal history of the study area, and establishing an instantaneous curve, a cumulative curve and a stage cumulative curve of natural gas under geological conditions on a geologic time scale.

Claims

1. A method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment, comprising the following steps: S1, performing a pyrolysis experiment on a selected kerogen sample to obtain activation energy, a mass fraction and a frequency factor of light carbon methane; S2, carrying out carbon isotope kinetics simulation of natural gas in a study area based on the activation energy, the mass fraction and the frequency factor of the light carbon methane, to obtain activation energy, a mass fraction and a frequency factor of heavy carbon methane; S3, establishing a burial history and a thermal history of the study area based on geological data; and S4, combining the activation energy, the mass fraction and the frequency factor of the heavy carbon methane with the burial history and the thermal history of the study area to obtain an instantaneous curve, a cumulative curve and a stage cumulative curve of the natural gas under geological conditions to form a history of natural gas accumulation; wherein step S2 comprises the following steps: S201: expanding the activation energy and the mass fraction of the light carbon methane, and setting expansion values as initial values of the activation energy and the mass fraction of heavy carbon methane; S202: calculating a frequency of the heavy carbon methane according to the frequency factor of the light carbon methane, wherein the frequency factor of the heavy carbon methane is expressed as follows: A.sup.13C=1.02×A.sup.12C wherein, A.sup.13C represents the frequency factor of the heavy carbon methane, and A.sup.12 represents the frequency factor of the light carbon methane; S203, calculating conversion rates and reaction rates of total reaction light carbon methane and heavy carbon methane precursors by the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and mass fraction of the heavy carbon methane and the frequency factor thereof; S204, calculating a carbon isotope cumulative amount according to the conversion rates and the reaction rates of the total reaction light carbon methane and heavy carbon methane precursors as well as the pyrolysis experiment data, wherein the carbon isotope cumulative amount is expressed as follows: δ.sup.13C (t)=[R.sub.0F(t)/F*(t)/R.sub.std−1]1000, where R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std, wherein δ.sup.13C(t) represents a carbon isotope cumulative amount at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; δ.sup.13C.sub.0 represents an initial carbon isotope amount, R.sub.std represents a standard ratio of Pee Dee Belemnite (PDB); S205, setting an activation energy difference ΔEa between the heavy carbon methane and the light carbon methane as a variable value; S206, obtaining a target value according to a carbon isotope measured amount from the pyrolysis experiment and the carbon isotope cumulative amount, wherein the target value X is expressed as follows: X=Σ.sub.t|δ.sup.13C*(t)−δ.sup.13C(t)| wherein, δ.sup.13C*(t) represents a carbon isotope measured amount at time t, and δ.sup.13C(t) represents a carbon isotope cumulative amount at time t; S207, establishing a constraint for obtaining the activation energy and mass fraction of the heavy carbon methane; and S208, fitting the target value by using a nonlinear interior point method according to the constraint; determining whether the target value is close to a preset limit value 0 during the fitting process; if yes, completing the fitting to obtain the activation energy, mass fraction and frequency factor of the heavy carbon methane; otherwise, returning to step S205.

2. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 1, wherein step S1 comprises the following steps: S101, performing the pyrolysis experiment on the selected kerogen sample to obtain a pyrolysis product of the selected kerogen sample and a yield thereof; and S102, carrying out hydrocarbon generation kinetics analysis based on the pyrolysis product and the yield thereof, to obtain the activation energy, the mass fraction and the frequency factor of light carbon methane.

3. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 1, wherein step S203 comprises the following steps: S2031, calculating reaction rate constants of all reaction light carbon methane and heavy carbon methane separately according to the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and the mass fraction of the heavy carbon methane as well as the frequency factor thereof; S2032, calculating unreacted amounts of all reaction light carbon methane and heavy carbon methane precursors separately according to the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and the mass fraction of the heavy carbon methane as well as the frequency factor thereof; S2033, calculating conversion rates of all the reaction light carbon methane and heavy carbon methane precursors separately according to the unreacted amounts of the light carbon methane and heavy carbon methane precursors; calculating reaction rates of all the reaction light carbon methane and heavy carbon methane precursors separately according to the conversion rates of the light carbon methane and heavy carbon methane precursors; S2034, calculating conversion rates and reaction rates of the total reaction light carbon methane and heavy carbon methane precursors separately according to the conversion rates and reaction rates of all the reaction light carbon methane and heavy carbon methane precursors.

4. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 3, wherein in step S2031, the reaction rate constant of each reaction light carbon methane or heavy carbon methane is expressed as follows:
k.sub.i=A exp(−Ea.sub.i/RT) wherein, Ea.sub.i represents the activation energy of the light carbon methane or the heavy carbon methane; A represents the frequency factor of the light carbon methane or the heavy carbon methane; R represents an ideal gas constant; T represents a temperature; k.sub.i represents a reaction rate constant of i-th reaction light carbon methane or heavy carbon methane; i=1, . . . , n, and n represents a total number of reactants.

5. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 3, wherein in step S2032, the unreacted amount of each reaction light carbon methane or heavy carbon methane precursor is expressed as follows: $w_i\left(t\right)=\exp \left[-\underset{t}{.Math.}\Delta U_i\left(t\right)\right]\Delta U_i\left(t\right)=\left[U_i\left(t\right)-U_i\left(t-1\right)\right]/\mathrm{Hr}{U}_i\left(t\right)=T\left(t\right)A\exp \left(-E{a}_i/RT\left(t\right)\right)\times \left\{1-\frac{{\left[E{a}_i/RT\left(t\right)\right]}^2+a_1\left[E{a}_i/RT\left(t\right)\right]+a_2}{{\left[E{a}_i/RT\left(t\right)\right]}^2+b_1\left[E{a}_i/RT\left(t\right)\right]+b_2}\right\}$ wherein, w.sub.i(t) represents an unreacted amount of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; ΔU.sub.i(t) represents an average cumulative amount of a light carbon methane or heavy carbon methane reaction product from time t to time t−1; Hr represents a constant heating rate; U.sub.i(t) represents a cumulative amount of the light carbon methane or heavy carbon methane reaction product at time t; R represents an ideal gas constant; T represents a temperature; k.sub.i represents a reaction rate constant of the i-th reaction light carbon methane or heavy carbon methane; Ea.sub.i represents the activation energy of the light carbon methane or heavy carbon methane; a.sub.1, a.sub.2, b.sub.1 and b.sub.2 represent constants for calculating the cumulative amount of the light carbon methane or heavy carbon methane reaction product at time t, a.sub.1=2.334733, a.sub.2=0.250621, b.sub.1=3.330657 and b.sub.2=1.681534.

6. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 3, wherein in step S2033, the conversion rate of each reaction light carbon methane or heavy carbon methane precursor is expressed as follows:
F.sub.i(t)=f.sub.o.sub.i[1−w.sub.i(t)] the reaction rate of each reaction light carbon methane or heavy carbon methane precursor is expressed as follows:
r.sub.i(t)=k.sub.i(t)[f.sub.o.sub.i−F.sub.i(t)] wherein, F.sub.i(t) represents a conversion rate of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; f.sub.oi represents a mass fraction of each reaction light carbon methane or heavy carbon methane; w.sub.i(t) represents an unreacted amount of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; r.sub.i(t) represents a reaction rate of the i-th reaction light carbon methane or heavy carbon methane at time t; k.sub.i(t) represents a reaction rate constant of the i-th reaction light carbon methane or heavy carbon methane at time t.

7. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 3, wherein in step S2034, the conversion rate of the total reaction light carbon methane or heavy carbon methane precursor is expressed as follows:
F*(t)=ΣF.sub.i(t)
F(t)=ΣF.sub.i(t) the reaction rate of the total reaction light carbon methane or heavy carbon methane precursor is expressed as follows:
r*(t)=Σr.sub.i(t)
r(t)=Σr.sub.i(t) wherein, F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors at time t; r(t) and r*(t) represent the reaction rates of the total reaction heavy carbon methane and light carbon methane precursors at time t; F.sub.i(t) represents the conversion rate of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; r.sub.i(t) represents the reaction rate of the i-th reaction light carbon methane or heavy carbon methane at time t.

8. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 1, wherein in step S207, the constraint is: a sum of the mass fractions f.sub.oi of all the reaction heavy carbon methane is set to 1, f.sub.oi≥0, and a range of the activation energy difference ΔEa between the reaction heavy carbon methane and the light carbon methane is set to 80-250 J/mol.

9. The method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment according to claim 1, wherein in step S4, the established cumulative curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C(t)=[R.sub.0F(t)/F*(t)/R.sub.std−1]1000
R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std wherein, δ.sup.13C.sub.cum(t) represents a carbon isotope cumulative amount of natural gas under geological conditions at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; δ.sup.13C.sub.0 represents an initial carbon isotope amount; R.sub.std represents a standard ratio of PDB; the established instantaneous curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C.sub.inst(t)={R.sub.0Σr.sub.i(t)/Σr.sub.i*(t)×[1−F(t)]/[1−F*(t)]/R.sub.std−1}1000
R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std wherein, δ.sup.13C.sub.inst(t) represents a carbon isotope instantaneous amount of natural gas under geological conditions at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; r.sub.i(t) and r.sub.i*(t) represent reaction rates of i-th reaction heavy carbon methane and light carbon methane precursors at time t; δ.sup.13C.sub.0 represents an initial carbon isotope amount; R.sub.std represents a standard ratio of PDB; the established stage cumulative curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C.sub.step(t)=δ.sup.13C.sub.inst(t−1)×{[F(t)−F(t−1)]/[F*(t)−F*(t−1)]}+{[F(t)−F(t−1)]/[F*(t)−F*(t−1)]−1}1000 wherein, δ.sup.13C.sub.step(t) represents a carbon isotope stage cumulative amount of natural gas under geological conditions at time t; F(t) represents a conversion rate of the total reaction heavy carbon methane precursor at time t; F*(t) represents a conversion rate of the total reaction light carbon methane precursor.

Description

BRIEF DESCRIPTION OF DRAWINGS

(1) FIG. 1 is a flowchart of a method according to the present invention.

(2) FIG. 2 shows an air source comparison between an instantaneous curve, a cumulative curve and a stage cumulative curve of natural gas and a carbon isotope measured amount as well as recovery of an accumulation history according to an example.

DETAILED DESCRIPTION

(3) The specific implementations of the present invention are described below to facilitate those skilled in the art to understand the present invention, but it should be clear that the present invention is not limited to the scope of the specific implementations. Various obvious changes made by those of ordinary skill in the art within the spirit and scope of the present invention defined by the appended claims should fall within the protection scope of the present invention.

(4) Limited by the simulation conditions and the complex simulation process, the application of the existing carbon isotope kinetics software based on a series of kinetic models is affected. The present invention conducts the carbon isotope kinetics simulation of methane by using a spreadsheet function of Excel. The method of the present invention is simple, and applies to carbon isotope kinetics simulation of natural gas in a wide range of areas, making the history of natural gas accumulation more authentic and reliable.

(5) As shown in FIG. 1, the present invention provides a method for forming a history of natural gas accumulation by using carbon isotopes by a pyrolysis experiment, including the following steps:

(6) S1, select a kerogen sample for a pyrolysis experiment to obtain activation energy, a mass fraction and a frequency factor of light carbon methane, including:

(7) S101, select a kerogen sample for a pyrolysis experiment to obtain a pyrolysis product of the kerogen sample and a yield thereof.

(8) S102, carry out hydrocarbon generation kinetics analysis by using KINETICS software based on the pyrolysis product and the yield thereof, to obtain activation energy, a mass fraction and a frequency factor of light carbon methane.

(9) S2, carry out carbon isotope kinetics simulation of natural gas in a study area by using a spreadsheet of Excel based on the activation energy, the mass fraction and the frequency factor of the light carbon methane, to obtain activation energy, a mass fraction and a frequency factor of heavy carbon methane.

(10) In a specific example, since heavy carbon methane occupies a very small proportion in methane while light carbon methane occupies a very large proportion, these parameters are assumed to be the frequency factor and activation energy distribution of light carbon methane, and the calculation process includes:

(11) S201: expand the activation energy and the mass fraction of the light carbon methane, and set expansion values as initial values of the activation energy and the mass fraction of heavy carbon methane.

(12) S202: calculate a frequency factor of the heavy carbon methane according to the frequency factor of the light carbon methane, where the frequency factor of the heavy carbon methane is expressed as follows:
A.sup.13C=1.02×A.sup.12C
where, A.sup.13C represents the frequency factor of the heavy carbon methane, and A.sup.12C represents the frequency factor of the light carbon methane.

(13) S203, calculate conversion rates and reaction rates of total reaction light carbon methane and heavy carbon methane precursors by the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and mass fraction of the heavy carbon methane and the frequency factor thereof, including:

(14) S2031, calculate reaction rate constants of all reaction light carbon methane and heavy carbon methane separately according to the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and the mass fraction of the heavy carbon methane as well as the frequency factor thereof.

(15) S2032, calculate unreacted amounts of all reaction light carbon methane and heavy carbon methane precursors separately according to the activation energy, the mass fraction and the frequency factor of the light carbon methane and the initial values of the activation energy and the mass fraction of the heavy carbon methane as well as the frequency factor thereof.

(16) S2033, calculate conversion rates of all the reaction light carbon methane and heavy carbon methane precursors separately according to the unreacted amounts of the light carbon methane and heavy carbon methane precursors; calculate reaction rates of all the reaction light carbon methane and heavy carbon methane precursors separately according to the conversion rates of the light carbon methane and heavy carbon methane precursors.

(17) S2034, calculate conversion rates and reaction rates of the total reaction light carbon methane and heavy carbon methane precursors separately according to the conversion rates and reaction rates of all the reaction light carbon methane and heavy carbon methane precursors.

(18) In a specific example, the entire parameter fitting process involves the following calculation formula:

(19) The process of kerogen pyrolysis to methane is represented by a series of parallel first-order reactions, and the rate constant of each parallel first-order reaction follows a semi-quantitative Arrhenius equation:
k.sub.i=A exp(−Ea.sub.i/RT)
where, t represents a temperature; R represents an ideal gas constant (J/mol/K); k.sub.i represents a reaction rate constant of i-th reaction light carbon methane or heavy carbon methane; i=1, . . . , n, and n represents a total number of reactants; A represents the frequency factor of the light carbon methane or the heavy carbon methane; Ea.sub.i represents the activation energy (J/mol) of the light carbon methane or heavy carbon methane.

(20) The experimental sample is type I kerogen, so it is further assumed that all parallel reaction methane precursors have the same initial carbon isotopic composition, and that each parallel reaction has the same frequency factor, and the activation energy is discretely distributed. The calculation is as follows:

(21) $w_i\left(t\right)=\exp \left[-{\int }_0^{t}k\left(t\right)\mathrm{dt}\right]$

(22) where, w.sub.i(t) represents an unreacted amount of the i-th reaction light carbon methane or heavy carbon methane precursor at time t.

(23) At a constant heating rate Hr, the above formula is expressed as:

(24) $w_i\left(t\right)=\exp \left[-\underset{t}{.Math.}\Delta U_i\left(t\right)\right]\Delta U_i\left(t\right)=\left[U_i\left(t\right)-U_i\left(t-1\right)\right]/\mathrm{Hr}{U}_i\left(t\right)=T\left(t\right)A\exp \left(-E{a}_i/RT\left(t\right)\right)\times \left\{1-\frac{{\left[E{a}_i/RT\left(t\right)\right]}^2+a_1\left[E{a}_i/RT\left(t\right)\right]+a_2}{{\left[E{a}_i/RT\left(t\right)\right]}^2+b_1\left[E{a}_i/RT\left(t\right)\right]+b_2}\right\}$
where, w.sub.i(t) represents an unreacted amount of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; ΔU.sub.i(t) represents a cumulative amount of a light carbon methane or heavy carbon methane reaction product from time t to time t−1; Hr represents a constant heating rate; U.sub.i(t) represents an average cumulative amount of the light carbon methane or heavy carbon methane reaction product at time t; R represents an ideal gas constant; T represents a temperature; k.sub.i represents a reaction rate constant of the i-th reaction light carbon methane or heavy carbon methane; Ea.sub.i represents the activation energy of the light carbon methane or heavy carbon methane; a.sub.1, a.sub.2, b.sub.1 and b.sub.2 represent constants for calculating the cumulative amount of the light carbon methane or heavy carbon methane reaction product at time t, a.sub.1=2.334733, a.sub.2=0.250621, b.sub.1=3.330657 and b.sub.2=1.681534; then, the following formulas are derived:
F.sub.i(t)=f.sub.oi[1−w.sub.i(t)]
r.sub.i(t)=k.sub.i(t)[f.sub.oi−F.sub.i(t)]
where, F.sub.i(t) represents a conversion rate of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; f.sub.oi represents a mass fraction of each reaction light carbon methane or heavy carbon methane; w.sub.i(t) represents an unreacted amount of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; r.sub.i(t) represents a reaction rate of the i-th reaction light carbon methane or heavy carbon methane at time t; k.sub.i(t) represents a reaction rate constant of the i-th reaction light carbon methane or heavy carbon methane at time t.

(25) Therefore, the sum of n parallel reactions is expressed as:
F*(t)=ΣF.sub.i(t)
F(t)=ΣF.sub.i(t)
r*(t)=Σr.sub.i(t)
r(t)=Σr.sub.i(t)
where, F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors at time t; r(t) and r*(t) represent the reaction rates of the total reaction heavy carbon methane and light carbon methane precursors at time t; F.sub.i(t) represents the conversion rate of the i-th reaction light carbon methane or heavy carbon methane precursor at time t; r.sub.i(t) represents the reaction rate of the i-th reaction light carbon methane or heavy carbon methane at time t.

(26) S204, calculate a carbon isotope cumulative amount according to the conversion rates and the reaction rates of the total reaction light carbon methane and heavy carbon methane precursors as well as the pyrolysis experiment data, where the carbon isotope cumulative amount is expressed as follows:
δ.sup.13C(t)=[R.sub.0F(t)/F*(t)/R.sub.std−1]1000
R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std
where, δ.sup.13C(t) represents a carbon isotope cumulative amount at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; δ.sup.13C.sub.0 represents an initial carbon isotope amount; R.sub.std represents a standard ratio of Pee Dee Belemnite (PDB).

(27) S205, set an activation energy difference ΔEa between the heavy carbon methane and the light carbon methane as a variable cell in a spreadsheet window of Excel.

(28) S206, obtain a target value in the spreadsheet function according to a carbon isotope measured amount from the pyrolysis experiment and the carbon isotope cumulative amount, where the target value X is expressed as follows:
X=Σ.sub.t|δ.sup.13C*(t)−δ.sup.13C(t)|
where, δ.sup.13C*(t) represents a carbon isotope measured amount at time t, and δ.sup.13C(t) represents a carbon isotope cumulative amount at time t.

(29) S207, establish a constraint for obtaining the activation energy and mass fraction of the heavy carbon methane, specifically: a sum of the mass fractions f.sub.oi of all the reaction heavy carbon methane is set to 1, f.sub.oi≥0, and a range of the activation energy difference ΔEa between the reaction heavy carbon methane and the light carbon methane is set to 80-250 J/mol.

(30) S208, fit the target value by using a nonlinear interior point method of Excel according to the constraint; determine whether the target value is close to a preset limit value 0 during the fitting process; if yes, complete the fitting to obtain the activation energy, mass fraction and frequency factor of the heavy carbon methane; otherwise, return to step S205.

(31) As shown in FIG. 2, in a specific example, the carbon isotope kinetics simulation of methane is conducted by using the spreadsheet function of Excel based on the carbon isotope amount measured by the pyrolysis experiment, to obtain the activation energy distribution and frequency factor of the heavy carbon methane. The parameter fitting is as follows:

(32) The activation energy and mass fraction of the light carbon methane are expanded and assumed to be the initial values of the kinetic parameters of the heavy carbon methane. The above calculation formulas are entered into Excel, and the activation energy difference ΔEa between the heavy and light carbon methane is set as a variable cell. A function value of a sum of errors between the carbon isotope measured amount from the pyrolysis experiment and the calculated carbon isotope cumulative amount is set as a target value, and the target value is set to 0. Then the constraint of the calculation process is established. Finally, the activation energy difference ΔEa and the mass fraction are solved by using a nonlinear interior point method of Excel. The error is minimized through repeated adjustments, thereby obtaining the activation energy difference ΔEa. If the error of the fitted parameter is still large, the mass fraction is continuously adjusted by using the above fitting method until the fitting result is satisfactory and the error is minimized.

(33) S3, establish a burial history and a thermal history of the study area based on geological data.

(34) S4, combine the activation energy, the mass fraction and the frequency factor of the heavy carbon methane with the burial history and the thermal history of the study area to establish an instantaneous curve, a cumulative curve and a stage cumulative curve of the natural gas under geological conditions to form a history of natural gas accumulation, thereby completing the formation of the natural gas accumulation history by using carbon isotopes by the pyrolysis experiment.

(35) The established cumulative curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C(t)=[R.sub.0F(t)/F*(t)/R.sub.std−1]1000
R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std
where, δ.sup.13C(t) represents a carbon isotope cumulative amount of natural gas under geological conditions at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; δ.sup.13C.sub.0 represents an initial carbon isotope amount; R.sub.std represents a standard ratio of PDB.

(36) The established instantaneous curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C.sub.inst(t)={R.sub.0Σr.sub.i(t)/Σr.sub.i*(t)×[1−F(t)]/[1−F*(t)]/R.sub.std−1}1000
R.sub.0=(δ.sup.13C.sub.0/1000+1)R.sub.std
where, δ.sup.13C.sub.inst(t) represents a carbon isotope instantaneous amount of natural gas under geological conditions at time t; R.sub.0 represents an initial carbon isotope ratio of the methane precursors; F(t) and F*(t) represent the conversion rates of the total reaction heavy carbon methane and light carbon methane precursors, respectively; r.sub.i(t) and r.sub.i*(t) represent reaction rates of i-th reaction heavy carbon methane and light carbon methane precursors at time t; δ.sup.13C.sub.0 represents an initial carbon isotope amount; R.sub.std represents a standard ratio of PDB.

(37) The established stage cumulative curve of natural gas under geological conditions is expressed as follows:
δ.sup.13C.sub.step(t)=δ.sup.13C.sub.inst(t−1)×{[F(t)−F(t−1)]/[F*(t)−F*(t−1)]}+{[F(t)−F(t−1)]/[F*(t)−F*(t−1)]−1}1000
where, δ.sup.13C.sub.step(t) represents a carbon isotope stage cumulative amount of natural gas under geological conditions at time t; F(t) represents a conversion rate of the total reaction heavy carbon methane precursor at time t; F*(t) represents a conversion rate of the total reaction light carbon methane precursor.

(38) Through the above design, the present invention provides a new way of recovering the history of natural gas accumulation, breaking the limitation of static models established only based on thermal simulation experimental data or measured carbon isotope data. The present invention improves the accuracy of the natural gas accumulation history, and applies to carbon isotope kinetics simulation of natural gas in a wide range of areas, making the natural gas accumulation history more authentic and reliable.