Shale anisotropic rock physical modeling method, electronic equipment and medium
The proposed method addresses the complexity of shale rock physics by modeling the anisotropic properties of shale formations, improving the prediction of anisotropic parameters and reservoir characteristics through advanced elastic matrix calculations.
Patent Information
- Application Number
- CN202311507137.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2023-11-13
- Publication Date
- 2025-05-13
AI Technical Summary
Existing rock physics models for shale formations do not adequately consider the complex distribution of clay minerals, organic matter, and directional arrangement of particles, which complicates the analysis of anisotropic parameters and reservoir properties in shale oil and gas exploration.
A method for modeling the anisotropic rock physics of shale formations, involving the calculation of elastic matrices for clay and non-clay minerals, and incorporating these into a composite model using Backus and Hashin-Shtrikman theories to account for the directional arrangement and porosity, allowing for more accurate prediction of anisotropic parameters.
The method provides a more precise model for predicting anisotropic parameters and reservoir properties, enhancing the accuracy of seismic inversion and imaging in shale formations.
Smart Images

Figure CN119989603A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the field of rock physics technology, and more specifically, to a shale anisotropic rock physics modeling method, electronic equipment and medium. Background Art
[0002] Rock physics technology is gradually playing an important role in shale oil and gas exploration and development. The establishment of a rock physics model for marine Permian shale based on logging data and rock physics theory can provide methods and technologies for predicting key reservoir physical properties and anisotropic parameters, reservoir mechanical properties and oil-bearing evaluation, and provide a more accurate anisotropic velocity model for shale seismic forward inversion and imaging technology. Analysis of geophysical logging and core data shows that the microscopic characteristics of marine Permian shale are mainly non-clay minerals such as pyrite, carbonate rock and felsic distributed in a matrix composed of clay minerals. Among them, the clay matrix is a mixture mainly composed of illite and montmorillonite mineral particles. The clay mixture has a microscopically oriented structure, which makes the shale solid matrix show a horizontal bedding structure on a macro scale, and thus has VTI (Transverselyisotropy with a vertical symmetry axis) type of inherent anisotropy. The organic mixture is mainly composed of kerogen and hydrocarbon-saturated pores. Therefore, the spatial distribution of minerals such as clay, organic kerogen, pores and fractures in shale is complex. At this time, rock physics theories or modeling methods that can take these complex factors into account are needed. However, most of the current shale rock physics models do not take these situations into account or do not take them into account sufficiently.
[0003] Therefore, it is necessary to conduct rock physics modeling research on marine Permian shale based on the above characteristics and develop a shale anisotropic rock physics modeling method.
[0004] The information disclosed in the background technology section of the present invention is only intended to deepen the understanding of the general background technology of the present invention, and should not be regarded as acknowledging or suggesting in any form that the information constitutes the prior art already known to those skilled in the art. Summary of the invention
[0005] The present invention proposes a shale anisotropic rock physics modeling method, electronic equipment and medium, the purpose of which is to solve the difficult problems existing in the shale rock physics modeling research, such as the complex organic matter occurrence and directional arrangement of clay mineral particles in marine Permian shales, and provide a research basis for the analysis of key physical properties and anisotropic parameters of reservoirs, and the evaluation of reservoir mechanical properties and oil content.
[0006] In a first aspect, an embodiment of the present disclosure provides a shale anisotropic rock physics modeling method, comprising:
[0007] Calculation of the elastic anisotropy matrix C of illite and montmorillonite mineral particles cl0 , which is the VTI anisotropic elastic modulus C of the clay mixture cl ;
[0008] Calculation of elastic modulus C of non-clay minerals ncl ;
[0009] Calculation of VTI anisotropy C of a solid matrix composed of clay minerals and non-clay minerals 1 ;
[0010] Calculation of the elastic matrix C of porous kerogen using Krief formula k ;
[0011] Based on the elastic matrix C 1 and C k , the elastic matrix C of VTI matrix 2 is calculated using anisotropic equivalent constant theory 2 ;
[0012] Based on the elastic matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat .
[0013] Preferably, the elastic anisotropy matrix C cl0 for:
[0014]
[0015] in, C 12 =C 11 - <c 11 >+ <c 12 >, C 66 = <c 66 >, c 11 , c 12 , c 13 , c 33 , c 44 , c 66 are the elastic stiffness coefficients of the equivalent media of each component before Backus averaging; C 11 , C 12 , C 13 , C 33 , C 44 , C 66 They are five independent elastic parameters describing the VTI equivalent medium after Backus averaging; the symbol <> indicates the weighted average of elastic parameters of multiple components.
[0016] Preferably, the elastic modulus C of non-clay minerals is calculated ncl include:
[0017] The shale matrix is considered as a mixture of non-clay minerals and clay and minerals composed of pyrite, carbonate and felsic. The elastic modulus C of non-clay minerals is calculated using the Hashin-Shtrikman limit theory (HS). ncl .
[0018] Preferably, the elastic modulus C ncl for:
[0019] c 12 =c 11 -2c 44
[0020] in, c 44 =μ.
[0021] Preferably, the elastic matrix C of the kerogen with organic pores is present. k for:
[0022] c 12 =c 11 -2c 44
[0023] in, c 44 =μ dry .
[0024] Preferably, the elastic matrix C 2 for:
[0025]
[0026] Preferably, the inorganic aperture system stiffness matrix C f for:
[0027] C f =C m0 -C(f)
[0028] Among them, C m0 is the isotropic elastic tensor of the matrix,
[0029] Preferably, the elastic matrix C of the marine Permian shale sat for:
[0030] C sat =C 2 -C f .
[0031] In a second aspect, an embodiment of the present disclosure further provides an electronic device, the electronic device comprising:
[0032] A memory storing executable instructions;
[0033] A processor runs the executable instructions in the memory to implement the shale anisotropic rock physics modeling method.
[0034] In a third aspect, the embodiments of the present disclosure further provide a computer-readable storage medium, which stores a computer program, and when the computer program is executed by a processor, the shale anisotropic rock physics modeling method is implemented.
[0035] The methods and apparatus of the present invention have other features and advantages that will be apparent from, or will be described in detail in, the accompanying drawings and subsequent detailed descriptions incorporated herein, which together serve to explain the specific principles of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS
[0036] The above and other objects, features and advantages of the present invention will become more apparent through a more detailed description of exemplary embodiments of the present invention in conjunction with the accompanying drawings, wherein like reference numerals generally represent like components throughout the exemplary embodiments of the present invention.
[0037] Figure 1 A flow chart showing the steps of a shale anisotropic rock physics modeling method according to one embodiment of the present invention.
[0038] Figure 2a and Figure 2b The graphs respectively show the comparison between the model calculation result according to one embodiment of the present invention and the P-wave velocity curve and S-wave velocity curve of the actual logging data. DETAILED DESCRIPTION
[0039] The preferred embodiments of the present invention will be described in more detail below. Although the preferred embodiments of the present invention are described below, it should be understood that the present invention can be implemented in various forms and should not be limited to the embodiments set forth herein.
[0040] The present invention provides a shale anisotropic rock physics modeling method, comprising:
[0041] Calculation of the elastic anisotropy matrix C of illite and montmorillonite mineral particles cl0 , which is the VTI anisotropic elastic modulus C of the clay mixture cl ;
[0042] Calculation of elastic modulus C of non-clay mineralsncl ;
[0043] Calculation of VTI anisotropy C of a solid matrix composed of clay minerals and non-clay minerals 1 ;
[0044] Calculation of the elastic matrix C of porous kerogen using Krief formula k ;
[0045] Based on the elastic matrix C 1 and C k , the elastic matrix C of VTI matrix 2 is calculated using anisotropic equivalent constant theory 2 ;
[0046] Based on the elastic matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat .
[0047] In one example, the elastic anisotropy matrix C cl0 for:
[0048]
[0049] in, C 12 =C 11 - <c 11 >+ <c 12 >, C 66 = <c 66 >, c 11 , c 12 , c 13 , c 33 , c 44 , c 66 are the elastic stiffness coefficients of the equivalent media of each component before Backus averaging; C 11 , C 12 , C 13 , C 33 , C 44 , C 66 They are five independent elastic parameters describing the VTI equivalent medium after Backus averaging; the symbol <> indicates the weighted average of elastic parameters of multiple components.
[0050] In one example, the elastic modulus C of a non-clay mineral is calculated. ncl include:
[0051] The shale matrix is considered as a mixture of non-clay minerals and clay and minerals composed of pyrite, carbonate and felsic. The elastic modulus C of non-clay minerals is calculated using the Hashin-Shtrikman limit theory (HS). ncl .
[0052] In one example, the elastic modulus C ncl for:
[0053] c 12 =c 11 -2c 44
[0054] in, c 44 =μ.
[0055] In one example, the elastic matrix C of kerogen with organic pores is k for:
[0056] c 12 =c 11 -2c 44
[0057] in, c 44 =μ dry .
[0058] In one example, the elasticity matrix C 2 for:
[0059]
[0060] In one example, the stiffness matrix C of the inorganic aperture system is f for:
[0061] C f =C m0 -C(f)
[0062] Among them, C m0 is the isotropic elastic tensor of the matrix,
[0063] In one example, the elastic matrix C of the marine Permian shale sat for:
[0064] C sat =C 2 -C f .
[0065] Specifically, the elastic anisotropy C of illite and montmorillonite mineral particles is calculated by the anisotropic Backus theory. cl0 :
[0066] Backus (1962) proposed that in the long wavelength limit, multi-layer horizontal thin-layer media (single-layer thin-layer media are isotropic media) can be equivalent to transversely isotropic elastic media, and the stiffness coefficient matrix composed of five independent moduli, the elastic matrix C cl0 for:
[0067]
[0068] in, C 12 =C 11 - <c 11 >+ <c 12 >, C 66 = <c 66 >, c 11 , c 12 , c 13 , c 33 , c 44 , c 66 are the elastic stiffness coefficients of the equivalent media of each component before Backus averaging; C 11 , C 12 , C 13 , C 33 , C 44 , C 66 They are five independent elastic parameters describing the VTI equivalent medium after Backus averaging; the symbol <> indicates the weighted average of elastic parameters of multiple components.
[0069] Sevostianov et al. (2005) proposed anisotropic equivalent field theory. This method is different from the traditional self-consistent approximation (SCA) and differential effective medium (DEM) theory. Instead, it places the inhomogeneous medium inclusions in the equivalent stress field and calculates the elastic coefficient tensor of the equivalent medium by considering the constitutive relationship of the medium.
[0070] Where: tensor is the matrix elastic coefficient, is the elastic coefficient of the inclusion, c is the volume fraction of the inclusion, and P ijkl =∫ V G ik,lj (x-x')dx'| (ij)(kl), where: G(x) is the Green's function of the anisotropic medium.
[0071] The anisotropic equivalent field theory formula is used to calculate the effect of intergranular soft matter, assuming that the soft matter between clay particles is included, and substituting The VTI anisotropic elastic modulus C of the clay mixture was obtained cl , expressed as:
[0072]
[0073] The shale matrix is considered as a mixture of non-clay minerals and clay and minerals composed of pyrite, carbonate and felsic. The elastic modulus C of non-clay minerals is calculated using the Hashin-Shtrikman boundary theory (HS). ncl .
[0074] The Hashin-Shtrikman (HS) rock physics boundary theory can give the narrowest upper and lower limits of the elastic modulus of a mixture composed of two or more mineral phases. For mixtures containing more than two components, Berryman extended the theory to the following general form:
[0075] K HS+ =Λ(μ max )
[0076] K HS- =Λ(μ min )
[0077] μ HS+ =Γ(ζ(K max ,μ max ))
[0078] μ HS- =Γ(ζ(K min ,μ min ))
[0079] In the formula, The brackets <> indicate the weighted average of the volumes of the components in the mixture.
[0080] Calculate the elastic matrix C based on the bulk modulus K and shear modulus μ of non-clay ncl , the elasticity matrix C ncl for:
[0081] c 12 =c 11 -2c 44
[0082] in, c 44 =μ.
[0083] The improved anisotropic Backus theory is used to calculate the VTI anisotropy C of a solid matrix composed of clay minerals and non-clay minerals. 1 , expressed as:
[0084]
[0085] In the formula, C 11 , C 12 , C 13 , C 33 , C 44 , C 66 Respectively by C cl and C ncl The elastic stiffness coefficient is calculated.
[0086] Calculation of the elastic matrix C of porous kerogen using Krief formula k .
[0087] The Krief formula is as follows:
[0088]
[0089]
[0090] In the formula, K m and K dry are the bulk moduli of kerogen and kerogen rich in organic pores, μ m and dry are the shear moduli of kerogen and kerogen rich in organic pores, respectively. In addition, in is the porosity.
[0091] Elastic matrix C of kerogen with organic pores k for:
[0092] c 12 =c 11 -2c 44
[0093] in, c 44 =μ dry .
[0094] Based on the elastic matrix C 1 and C k , the VTI matrix ② elastic matrix C is calculated using the anisotropic equivalent constant theory 2 :
[0095]
[0096] Based on the elastic matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat :
[0097] Hudson (1980, 1981) analyzed the scattering theory of the average wave field of elastic solids containing thin coin-shaped ellipsoidal cracks or inclusions. By assuming that there are idealized coin-shaped cracks in the isotropic background medium and that the cracks are independent of each other and there is no fluid flow, he proposed the Hudson model to estimate the equivalent elastic modulus and attenuation of rocks. Schoenberg (1980) proposed the linear slip theory, ignoring the shape of the cracks, and equating the cracks to linear slip interfaces to conduct equivalent medium modeling research on fractured formations. Chapman (2003) proposed a multi-scale fracture theory to describe the jet flow mechanism of fractured porous media. In the model, the fracture-pore system is divided into spherical pores, microcracks, and directional cracks. Assume that the aspect ratio of microcracks and cracks is r, the radius of microcracks and spherical pores is the same as the rock particle size, both represented by a, and the radius of the cracks is much larger than that of spherical pores and microcracks, represented by a. f If , then the volume of spherical pores, microcracks and cracks can be expressed as:
[0098]
[0099]
[0100]
[0101] In the crack-pore system, microcracks are interconnected with each other, spherical pores are interconnected with each other, and each crack can be connected to multiple spherical pores and microcracks, but cracks are not connected with each other, and each microcrack and each spherical pore is connected to at most one crack. Therefore, the effective stiffness tensor expression in the Chapman model is:
[0102] C=C m0 -Φ p C p -ε c C c1 -ε f C c2
[0103] Where: C m0 is the isotropic elastic tensor of the matrix, C p , C c1 and C c2 are the disturbances caused by spherical holes, microcracks and cracks, respectively, and C p and C c1is a function of the Lamé parameter and the fluid, C c2 is a function of Lamé parameter, fluid, frequency and relaxation time; Φ p , ε c and ε f are functions of spherical pore porosity, microcrack density and fracture density respectively. The model is parameterized so that the above formula can intuitively describe the influence of each parameter on rock properties:
[0104] C(f)=C m0 (Λ,Γ)-Φ p C p (λ 0 ,μ 0 ,f)-ε c C c (λ 0 ,μ 0 ,f)-ε f C f (λ 0 ,μ 0 ,f)
[0105] Where: f is the frequency; Λ and Γ are Lame constants, which have no clear physical meaning and can be calculated from isotropic longitudinal and transverse wave velocities; the independent components of the elastic matrix C(f) are:
[0106]
[0107]
[0108]
[0109]
[0110]
[0111] In elasticity, the four subscripts of the stiffness tensor can be converted to two (Auld, 1990). In the above formula, C 1111 =C 11 ,C 2323 =C 44 ,C 3333 =C 33 ,C 1122 =C 12 ,C 1133 =C 13 , so the elasticity matrix C(f) can be expressed as:
[0112]
[0113] Stiffness matrix C of inorganic hole-slot system f It can be calculated by the following formula:
[0114] C f =C m0 -C(f)
[0115] Among them, C m0 is the isotropic elastic tensor of the matrix calculated by the HS limit theory assuming that the matrix is isotropic.
[0116] Marine Permian Shale Elastic Matrix C sat for:
[0117] C sat =C 2 -C f .
[0118] The present invention also provides an electronic device, which includes: a memory storing executable instructions; and a processor, which runs the executable instructions in the memory to implement the above-mentioned shale anisotropic rock physics modeling method.
[0119] The present invention also provides a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, the above-mentioned shale anisotropic rock physics modeling method is implemented.
[0120] To facilitate understanding of the solutions and effects of the embodiments of the present invention, three specific application examples are given below. Those skilled in the art should understand that the examples are only for facilitating understanding of the present invention, and any specific details thereof are not intended to limit the present invention in any way.
[0121] Example 1
[0122] Figure 1 A flow chart showing the steps of a shale anisotropic rock physics modeling method according to the present invention.
[0123] like Figure 1 As shown, the shale anisotropic rock physics modeling method includes: step 101, calculating the elastic anisotropy matrix C of illite and montmorillonite mineral particles cl0 , which is the VTI anisotropic elastic modulus C of the clay mixture cl ; Step 102, calculate the elastic modulus C of non-clay minerals ncl Step 103, calculate the VTI anisotropy C of the solid matrix composed of clay minerals and non-clay minerals 1 Step 104, using the Krief formula to calculate the elastic matrix C of the porous kerogen k ; Step 105, based on the elasticity matrix C 1 and C k , the elastic matrix C of VTI matrix 2 is calculated using anisotropic equivalent constant theory2 ; Step 106, based on the elasticity matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat .
[0124] The elastic anisotropy of illite and montmorillonite mineral particles was calculated by anisotropic Backus theory, and the influence of intergranular soft matter was calculated by anisotropic equivalent field theory. The anisotropic elastic modulus C of the VTI of the clay mixture was obtained. cl The shale matrix is considered as a mixture of non-clay minerals and clay and minerals composed of pyrite, carbonate and felsic, and the elastic modulus C of non-clay minerals is calculated using the Hashin-Shtrikman limit theory (HS). ncl ; The improved anisotropic Backus theory is used to calculate the VTI anisotropy C of a solid matrix composed of clay minerals and non-clay minerals. 1 ; The elastic matrix C of the saturated fluid in the organic pores in kerogen is calculated using the Krief formula k ; Based on the elasticity matrix C 1 and C k , the VTI matrix ② elastic matrix C is calculated using the anisotropic equivalent constant theory 2 ; Based on the elasticity matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat .
[0125] Figure 2a and Figure 2b The following graphs respectively show the comparison between the model calculation results and the P-wave and S-wave velocity curves of the actual logging data according to an embodiment of the present invention. The blue data points are the P-wave and S-wave velocities of the shale layer calculated based on the model, and the black curves are the P-wave and S-wave velocities of the actual logging data. By comparison, it is found that the model calculation results are in good agreement with the actual logging data curves.
[0126] Example 2
[0127] The present disclosure provides an electronic device, which includes: a memory storing executable instructions; a processor, which runs the executable instructions in the memory to implement the above-mentioned shale anisotropic rock physics modeling method.
[0128] An electronic device according to an embodiment of the present disclosure includes a memory and a processor.
[0129] The memory is used to store non-temporary computer-readable instructions. Specifically, the memory may include one or more computer program products, which may include various forms of computer-readable storage media, such as volatile memory and / or non-volatile memory. The volatile memory may, for example, include random access memory (RAM) and / or cache memory (cache), etc. The non-volatile memory may, for example, include read-only memory (ROM), hard disk, flash memory, etc.
[0130] The processor may be a central processing unit (CPU) or other forms of processing units having data processing capabilities and / or instruction execution capabilities, and may control other components in the electronic device to perform desired functions. In one embodiment of the present disclosure, the processor is used to run the computer-readable instructions stored in the memory.
[0131] Those skilled in the art should be able to understand that in order to solve the technical problem of how to obtain a good user experience, the present embodiment may also include well-known structures such as a communication bus and an interface, and these well-known structures should also be included in the protection scope of the present disclosure.
[0132] For detailed description of this embodiment, reference may be made to the corresponding descriptions in the aforementioned embodiments, which will not be repeated here.
[0133] Example 3
[0134] An embodiment of the present disclosure provides a computer-readable storage medium storing a computer program, which implements the shale anisotropic rock physics modeling method when executed by a processor.
[0135] According to the computer-readable storage medium of the embodiment of the present disclosure, non-transitory computer-readable instructions are stored thereon. When the non-transitory computer-readable instructions are executed by a processor, all or part of the steps of the above-mentioned methods of each embodiment of the present disclosure are executed.
[0136] The above-mentioned computer-readable storage media include, but are not limited to: optical storage media (e.g., CD-ROM and DVD), magneto-optical storage media (e.g., MO), magnetic storage media (e.g., magnetic tape or mobile hard disk), media with built-in rewritable non-volatile memory (e.g., memory card) and media with built-in ROM (e.g., ROM box).
[0137] Those skilled in the art should understand that the purpose of the above description of the embodiments of the present invention is only to exemplarily illustrate the beneficial effects of the embodiments of the present invention, and is not intended to limit the embodiments of the present invention to any given examples.
[0138] The embodiments of the present invention have been described above, and the above description is exemplary, not exhaustive, and is not limited to the disclosed embodiments. Many modifications and changes will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments.
Claims
1. A shale anisotropic rock physics modeling method, characterized in that: include: Calculation of the elastic anisotropy matrix C of illite and montmorillonite mineral particles cl0 , which is the VTI anisotropic elastic modulus C of the clay mixture cl ; Calculation of elastic modulus C of non-clay minerals ncl ; Calculation of VTI anisotropy C of a solid matrix composed of clay minerals and non-clay minerals 1 ; Calculation of the elastic matrix C of porous kerogen using Krief formula k ; Based on the elastic matrix C 1 and C k , the elastic matrix C of VTI matrix 2 is calculated using anisotropic equivalent constant theory 2 ; Based on the elastic matrix C 2 and the stiffness matrix C of the inorganic hole-slot system f , calculate the elastic matrix C of marine Permian shale sat .
2. The shale anisotropic rock physics modeling method according to claim 1, wherein: The elastic anisotropy matrix C cl0 for: in, C 12 =C 11 - <c 11 >+ <c 12 >, C 66 = <c 66 >, c 11 , c 12 , c 13 , c 33 , c 44 , c 66 are the elastic stiffness coefficients of the equivalent media of each component before Backus averaging; C 11 , C 12 , C 13 , C 33 , C 44 , C 66 They are five independent elastic parameters describing the VTI equivalent medium after Backus averaging; the symbol <> indicates the weighted average of elastic parameters of multiple components.
3. The shale anisotropic rock physics modeling method according to claim 1, wherein: Calculation of elastic modulus C of non-clay minerals ncl include: The shale matrix is considered as a mixture of non-clay minerals and clay and minerals composed of pyrite, carbonate and felsic. The elastic modulus C of non-clay minerals is calculated using the Hashin-Shtrikman boundary theory (HS). ncl .
4. The shale anisotropic rock physics modeling method according to claim 3, wherein: The elastic modulus C ncl for: c 12 =c 11 -2c 44 in, c 44 =μ.
5. The shale anisotropic rock physics modeling method according to claim 1, wherein: Elastic matrix C of kerogen with organic pores k for: c 12 =c 11 -2c 44 in, c 44 =μ dry .
6. The shale anisotropic rock physics modeling method according to claim 1, wherein: The elasticity matrix C 2 for:
7. The shale anisotropic rock physics modeling method according to claim 1, wherein: The inorganic pore system stiffness matrix C f for: C f =C m0 -C(f) Among them, C m0 is the isotropic elastic tensor of the matrix, 8. The shale anisotropic rock physics modeling method according to claim 1, wherein: The marine Permian shale elastic matrix C sat for: C sat =C 2 -C f 。 9. An electronic device, characterized in that: The electronic device comprises: A memory storing executable instructions; A processor, wherein the processor runs the executable instructions in the memory to implement the shale anisotropic rock physics modeling method according to any one of claims 1 to 8.
10. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the shale anisotropic rock physics modeling method described in any one of claims 1 to 8 is implemented.