As unconventional resources and piedmont high and steep structures gradually turn to focused fields of oil and gas exploration & development in China, the evaluation of reservoir anisotropy has become increasingly important^[1,2]. Anisotropy refers to the phenomenon that some physical parameter values of rock are directional, such as acoustic velocity, conductivity, and permeability etc., that is, the results measured along different directions are different. The directional arrangement of particles inside the rock and the external temperature and pressure it bears are responsible for anisotropy. The scale of rock anisotropy ranges from several microns to tens of kilometers. According to the classification system of crystal symmetry and the symmetry that wave physics can realize in underground media, the basic symmetry of real media is divided into 10 categories in geophysics, including monoclinic symmetry, triclinic symmetry, orthogonal symmetry and so on. The difference between them is mainly manifested in the difference of elastic coefficient matrix, among which the hex-agonal anisotropic medium is also called transverse isotropy (TI for short) medium. It is the most common medium in the earth, and also the most widely used medium model in seismic exploration. TI medium with vertical symmetry axis is called the VTI medium^[3,4,5].

This study takes the tight sands of Chang 7 Member of Triassic Yanchang Formation (hereinafter referred to as Chang7) in Ordos Basin as an example to study the anisotropy. The deep-water turbidite sand of Chang7 is mainly composed of fine siltstone, which is in close contact with the source rock, forming "sandwich biscuit" structure of frequent interbeds^[6,7]. In addition, there is a large amount of argillaceous tearing debris in the local part of the slumping, which further aggravates the anisotropy of this kind of tight reservoir. Previous studies have shown that there are a large number of weak structural planes, such as bedding and schistosity composed of clay minerals or argillaceous fine silt. Particles at these planes have lower cohesive force and contact degree between each other, which results in obvious differences in the mechanical, electrical and seepage characteristics in the directions perpendicular and parallel to the schistosity, namely, anisotropy^[8,9].

Due to the common existence of anisotropy in the formation, big errors would appear in the process of geophysical inversion and interpretation of various parameters of the formation if a simple isotropic model is used. These errors have drawn attention of the academic community, which can be traced back to the 1950s. The cause was that the existence of thin interbeds in offshore seismic exploration led to over 5% of error in time-depth conversion. By introducing the concept of anisotropy, we can successfully explain the phenomenon that the vertical and horizontal seismic velocities are inconsistent. Since the 1990s, more and more researchers delved deep into this issue. Liu Yunhe et al. systematically summarized the research history of anisotropy in electromagnetic exploration, and clearly pointed out that the anisotropy was observed in the 1950s, but quantitative extraction and inversion of the degree of anisotropy in electromagnetic measurement data did not become a research hotspot until the beginning of this century^[3]. Through the experimental test and simulation analysis, Tang Xinwei et al. found out the anisotropic mechanical characteristics and deformation and fracture rules of slate, granite and shale^[9,10,11,12]. Zhang Bing, Huang Xinrui, Ding Gongbo, Yin Xingyao, et al. focused on the rock physical modeling considering anisotropy in the geophysical inversion process^{[13,14,15,16]}. Liu Zhonghua, Xiao Jiaqi, Shen Jinsong, Xu Song, et al. investigated how to accurately calculate horizontal principal stress, three-dimensional conductivity and S-wave anisotropy coefficient in anisotropic formation^{[17,18,19,20]}.

All the above researches mainly focused on anisotropic characteristics of rock from the perspective of fracture experiment or single acoustic/resistivity logging, which included mainly array lateral or 3D scanning resistivity datum. The forward modeling theory of 3D scanning resistivity logging is relatively perfect, but the inversion of electrical anisotropy has ambiguous solutions and slow progress. Although array lateral resistivity logging can provide electrical anisotropy information, previous studies were limited to 2D or 2.5D inversion only considering single factors such as drilling fluid invasion or formation dip angle, while the actual logs is subject to the comprehensive influence of multiple factors, so the inversion results only considering one factor are not satisfactory.

In view of the current status and existing problems of quantitatively evaluating the anisotropy from logging information, this study analyzed the factors affecting the response of array lateral resistivity logging based on numerical simulation. Since the key factor affecting the inversion speed and accuracy is the setting of the initial value, the multi-parameters fast inversion method of array lateral resistivity logging with an accurate estimation of the initial value is proposed. Moreover, a new method of extracting anisotropic information of resistivity according to the characteristics of micro-resistivity scanning imaging logging has been proposed when there is no array lateral logging data available. At last, the results of the above methods were compared to find out the consistency of the acoustic and electrical anisotropies based on the origin of anisotropy in tight sands reservoirs. The study provides a new technical idea for the evaluation of the anisotropy of tight sands reservoirs and important parameters for seismic inversion and reservoir stimulation.

As previously mentioned, resistivity anisotropy refers to the directional measurement result. Strictly speaking, the characteristic of resistivity anisotropy can only be accurately characterized in the laboratory, which is finished by measuring plug samples drilled in different directions. Generally the resistivity result parallel to the bedding direction is defined as horizontal resistivity and recorded as R_h, the result in orthogonal direction is vertical resistivity and recorded as R_V. Then the degree of resistivity anisotropy is quantitatively characterized as:

Assuming that the formation with bedding is in horizontal (Fig. 1a), the resistivity along the bedding direction and perpendicular to the bedding direction must differ somewhat. How does the resistivity logging tool response to this common model? Let’s take the dual lateral logging (DLL) as example. As we know, the electric current mainly flows in the horizontal direction parallel to the bedding under the action of focusing electrode, so the measuring result can be approximately regarded as R_h. When an inclined well is drilled (Fig. 1b), there is a certain angle between the well axis or tool axis and the formation. Due to the objective existence of electrical anisotropy, the measurement result is affected to a certain extent by the vertical resistivity although forced focusing is used. The greater the θ, the more significant the effect of R_V is. Apparently the output of R_t and R_xo are synthetically determined by the conductivity of drilling fluid, flushed zone, undisturbed formation, surrounding rock (layer thickness) and degree of anisotropy. Although service companies provide plates that can eliminate borehole effects, for anisotropic study, even at given θ, these 4 parameters, λ, r_i, Rxo and Rt can’t be accurately inversed based on DLL alone.

There is a trend that more and more array sensors are introduced during the development of logging tools. In 1998, Schlumberger launched high-resolution array lateral logging tool HRLT, which has six measurement modes. Except mode 0 which mainly measures wellbore fluid, the instrument outputs five apparent resistivity logs (R_LA1, R_LA2, ... R_LA5) with different investigation depths and vertical resolution (Fig. 2), which make it possible to study the anisotropy. Many researchers, including some with Schlumberger, have carried out research on forward simulation and inversion process of HRLT logs, but so far there is no successful case of inversion algorithm covering these four parameters λ, r_i, R_xo and R_t synchronously. The key reasons include slow processing speed and ambiguous solutions of the inversion^[21,22].

As we know, successful inversion is based on forward modeling. 3-D finite element method has been used for HRLT logging response forward simulation in any inclined formation with electrical anisotropy, drilling fluid invasion, and wall rock influence successfully. There are many literatures in this field, which will not be covered in this paper^{[21, 23]}. Given a formation model with infinite thickness and R_h = 20.0 Ω·m, θ are 0° and 60° respectively, d = 20 cm, the influence of λ on HRLT outputs under different inclined conditions is investigated and results are shown in Fig. 2. Ordinarily sand shale sections usually have λ of no more than 2. In Fig. 2a, when λ increases from 1.0 to 2.0, the relative change of the array log doesn’t exceed 30%, even in the 60° deviated well (Fig. 2b), the relative change does not exceed 50%.

Under the same borehole environment, giving θ=0°, the infinite thickness formation has λ=2.0 that keeps same within the flushed zone and undisturbed formation, the formation resistivity decreases when invasion occurs. The simulation results are shown in Fig. 3: even if no invasion occurs in the electrical anisotropy formation (corresponding to r_i =0 in Fig. 3), these five apparent resistivity logs from HRLT have values higher than the R_th of the formation. When invasion occurs, the apparent resistivity decreases sharply with the increase of r_i. When r_i≥10 cm, these five logs turn to deflected in opposite direction. When r_i≥25 cm, the log value of R_LA5 and R_LA1 decreased by 25% and 50% respectively. When r_i≥50 cm, R_LA1 decreased by more than 70%. We all know that the actual invasion depth r_i mainly depends on the wellbore-formation pressure difference, soaking time and reservoir physical properties, but generally is more than 50 cm, even more than 120 cm in extreme cases^[22].

According to the above analysis, under the actual wellbore conditions, those major factors affecting the response of HRLT logs include invasion of drilling fluid, relative dip angle and electrical anisotropy coefficient, but their degrees of influence decrease in turn. To inversely model these four parameters, λ, r_i, R_xo and R_t from HRLT logs accurately, these factors must be considered comprehensively, but they should be treated differently and step by step.

As we all know, the accuracy, speed and stability of the inversion algorithm largely depend on the selection of initial value. It is very important to give reasonable initial values of r_i, R_t and R_xo for the inversion process. In this study, a hierarchical flowsheet based on charts is proposed, which aims at accurately estimating the initial values of inversion parameters. The invasion is first consider, then the resistivity and λ. These four parameters in inclined anisotropic formation with invasion were modeled inversely (assuming that the flushed zone and the undisturbed formation possess same electrical anisotropy).

1.1.1. Determining the initial value of r_i inversion

Since the logs value of HRLT is highly sensitive to invasion degree, and the separation degree of 5 logs is directly related to r_i, we establish the charts shown in Fig. 4 by analyzing the forward simulation results at different invasion radii. In Fig. 4a, R_LA3/R_LA1 ratio of HRLT is taken as abscissa, R_LA5/R_LA1 ratio as ordinate, the red curve represents different r_i, the blue curve represents different R_t/R_xo ratios. For the layer in which we are interested, R_LA1, R_LA3 and R_LA5 values are picked respectively. According to their ratios and corresponding coordinate position, the initial value of r_i can be worked out between the two adjacent theoretical invasion radius values by interpolation method from Fig. 4a.

1.1.2. Determination of initial values of R_t and R_xo inversion

Similarly, we establish the chart as shown in Fig. 4b based on forward simulation analysis at different R_t/R_xo values under each r_i. In Fig. 4b, R_LA1 is choosen as the abscissa, R_LA5/R_LA1 is the ordinate, the blue curve in the longitudinal direction corresponds to different R_t/R_xo ratios, and the pink curve in the transverse direction corresponds to different R_xo. For the target interval, the inversion initial values of R_t and R_xo can be determined according to the positions of R_LA1 and R_LA5 curves and their ratio in Fig. 4b.

We establish a series of templates as shown in Fig. 4 when λ is from 1.0 to 3.0 with increment as 0.1. For the actual application, we assume that θ is known (determined by other logs), firstly we estimate a set of initial values r_i,R_xo and R_t and their corresponding λ value by the charts like Fig. 4, and then forward the theoretical curve with the 3D finite element method to match with the measured curve circularly., At last, we select the group with the smallest error as the output, and determine the electric anisotropy coefficient λ and other parameters, recorded as λ_HRLT. The error analysis shows that anisotropy and other parameters obtained by this method have relative errors of less than 10%, so it is applicative^[23].

Besides HRLT, three-component induction logging can also provide horizontal and vertical resistivity information ^[18], but it is rarely used in PetroChina, and domestic instruments haven’t been commercially used, so it is not discussed in this paper.

Micro-resistivity scanning imaging (shortened as micro-scanning-imaging) logging provides the micro-resistivity array data around the borehole. The change of image color in the vertical (along well axis) and horizontal direction (around the wellbore) represents the corresponding formation resistivity change. Therefore, the electrical anisotropy coefficient reflecting the resistivity difference in different directions can be extracted from the micro-scanning-imaging logging.

Fig. 5 is sketch of the data structure of the micro-scanning imaging logging. In Fig. 5a, there are N sample points at each depth point (take the FMI instrument of Schlumberger company as an example, N=192). M depth points are grouped into a group to form a 2-D matrix of M×N, corresponding to the scanning result of a thin-layer with thickness of H in the depth direction (Fig. 5b). To simplify the processing of 2-D matrix, the distribution frequency histogram of M×N elements of the matrix is counted, and its peak value is taken as the equivalent resistivity of thin layer H, which is recorded as R_t1. In the range of processing depth window length W, we assume that there are K matrices, corresponding to K thin layers, we can obtain R_t2, R_t3, ... R_tK. In order to study the electrical anisotropy, K thin layers within the window W can be regarded as series in the longitudinal direction, and their equivalent resistivity R_ev can be expressed as:

In the horizontal direction, K thin layers can be regarded as a parallel connection, and its equivalent resistivity R_eh can be expressed as:

When the equivalent resistivities in vertical and horizontal directions are obtained by equ. 2 and equ. 3, the quantitative evaluation of electrical anisotropy based on micro-scanning- imaging logging can be realized by using equ. 1, and the result is recorded as λ_FMI.

On the other hand, the azimuthal anisotropy of the formation resistivity can be analyzed according to the micro-scanning-image. For FMI tool, its pads are designed in pairs, each pair of pads has a 180° difference in spatial orientation, and the resistivity measured is along a certain orientation. According to Fig. 5, the peak value of the button resistivity distribution on a pair of pads within the processing window length W is taken as the equivalent resistivity of the position. If there are S pairs of pads (or all buttons are divided into S groups according to their directions), the resistivity R_tp1, R_tp2 … R_tpS of S directions can be extracted, the ratio of the maximum R_tpmax to the minimum R_tpmin represents the resistivity difference in different directions, hereinafter referred to as resistivity azimuthal anisotropy (recorded as λ_DANI), and the equation is:

Obviously, the ratio represents the difference of conductivity in different directions of the formation. The larger the value is, the more significant the difference is in different azimuth.

Cross-dipole sonic logging has been successfully used in the petroleum industry to evaluate the acoustic anisotropy for many years, and has been reported in many literatures^[24,25]. As we all know sonic wave propagates in porous rock through the vibration of particles, but crystal or mineral particles of sedimentary rock arrange differently in different directions. For the isotropic medium, the components of S-wave in all directions propagate at the same speed. But in anisotropic medium, S-wave splitting in which fast and slow S-waves are distributed along with different directions often occurs. Equ. 5 is commonly used to evaluate the anisotropy degree of fast and slow S-waves in acoustic anisotropic strata:

In addition to equ. 5, some other researchers introduce the difference of arrival time of fast and slow shear waves and the differences in energy and azimuth angle to evaluate the degree of anisotropy. Their principles are similar and needn’t repeat explanation here.

A large number of studies show that it is just because of the external stress field borne by the rock itself or the directional arrangement of particles in the rock that the natural rock has the inherent anisotropy. In 2002, Georgi D et al. reached the conclusion that the fluid seepage and electric current flow in rock had internal correlation through numerical simulation^[26]. They pointed out that on the micro-scale, it was exactly because of the pore structure difference in different directions caused by directional arrangement of particles, the rock had consistent anisotropy in permeability and resistivity. Teng Jiwen et al. pointed out that it is exactly because the shape of particles, preferred orientation of crystallization and fractures in the rock etc. lead to the anisotropy of the earth's internal medium and structure, the shear wave splitting occurs to form fast and slow shear waves^[27,28]. It can be deduced that the electrical anisotropy and acoustic anisotropy of rock are intrinsically related, and these two kinds of anisotropic coefficients as λ_FMI and λ_DANI representing the vertical-horizontal resistivity contrast and the resistivity contrast in different azimuths respectively are theoretically consistent with the anisotropic information of λ_HRLT and λ_SLOANI. Although different in the amplitude of numerical values, they should be comparable.

Fig. 6 shows an example of the electrical anisotropies from HRLT inversion and FMI image extraction. From Fig. 6, it is obvious that the invasion of the section measured in this case is relatively serious, and the 5 logs of HRLT(track 4) show regular amplitude deviation, and the amplitude difference relatively reduces in the corresponding mudstone sections less than 4250 m and 4261-4268 m depths. According to the inversion results in track 5, the section less than 4262 m deep is weaker in electric anisotropy, where the electric anisotropy coefficients λ_HRLT and λ_FMI are basically distributed in the range of background value. The degree of electrical anisotropy of 4262-4267 m section gets stronger, and it can be seen from the static image (track 3) that 1# and 3# pad show expansion in well diameter, which is the direct evidence of stress anisotropy leading to borehole collapse. Within the corresponding section, these two anisotropic coefficient values (track 5) also increase significantly, and reduce to the base value below 4267 m. Due to different information sources and calculation principles, these two anisotropic coefficient logs in this case differ in numerical distribution range and scale range, but they reflect the same trend of electric anisotropy, which shows that the electric anisotropy extracted from the FMI image is consistent with that from HRLT logging.

Fig. 7 shows the comparison of their consistency between λ_DANI and λ_SLOANI of two other cases. In these two examples, the electrical and acoustic anisotropies of the section are both weak generally, but still have some local fluctuations (Fig. 7). In Fig. 7a both λ_DANI and λ_SLOANI increase within the interval of 2345-2349 m in contrast with other section with background value. In Fig. 7b, these two logs overlap very well and show consistent trend of variation. This example demonstrates that the azimuthal anisotropy of resistivity extracted from FMI image by equ. 4 is consistent with the acoustic anisotropy from array acoustic logging and they are comparable.

Fig. 8 shows an example of a carbonate formation. For this case HRLT, FMI and cross-dipole acoustic logging are conducted. It can be seen from Fig. 8 that the trends of four anisotropic coefficients, as λ_HRLT, λ_FMI, λ_DANI and λ_SLOANI, have consistent change overall. Especially in the upper section of 4560-4580 m, λ_DANI and λ_SLOANI curves are very consistent. It means that the anisotropy caused by the internal structure of the rock makes the electrical and acoustic anisotropies consistent in strength. The anisotropy from FMI image can replace that from cross dipole acoustic logging for characterizing the reservoir anisotropy to some extent.

The above four cases show that the electrical and acoustic anisotropies caused by the intrinsic factors of the tight reservoir are consistent in trends. On the premise that the quality of HRLT logs and FMI image are good, the methods pre-sented in this paper can provide effective electrical and acoustic anisotropies, and are comparable with the traditional anisotropic evaluation method based only on the cross dipole acoustic logging. Besides, the premise of equ. 2-4 is that all button electrodes in the processing window represent the resistivity information of the same sublayer. If the relative dip angle is large, and the data of the same sampling point corresponds to the resistivity information of different sublayers, then the influence of the relative dip angle needs to be considered, and the above method needs to be further improved.

Array lateral resistivity logging describes electrical anisotropy through combination of several detection modes. Numerical simulation shows that λ, θ, r_i are the major factors that determine its output, and the invasion has more significant effect in thick layers. By introducing the 3D forward modeling algorithm to build a new template to estimate the initial inversion value, and using the hierarchical inversion process, we can get accurate electric anisotropy coefficient and horizontal resistivity of the formation.

Micro-resistivity scanning image provides the information of resistivity distribution with azimuth. The equivalent resistivity in horizontal direction and vertical direction and their ratio can be used for characterizing the electric anisotropy of formation. The average resistivity of different directions can be extracted, and the comparison of the maximum and minimum average resistivity can represent the azimuthal anisotropy of formation resistivity.

The electric and acoustic anisotropies of the tight reservoir mainly depend on the internal microstructure of the tight reservoir and external temperature and pressure, but they are consistent. Case analysis shows that the anisotropic coefficients extracted from different logging information also reflect consistent degree of electrical and acoustic anisotropies. Therefore, array lateral resistivity, micro-scanning-imaging and array acoustic logging all can accurately evaluate the anisotropic strength of the tight reservoir. However, it should be pointed out that the main directions of sound wave propagation and current propagation in rock porous media are controlled by different mechanisms respectively, and the influence mechanism is more complex. This study can’t reflect the differences in the two rock physical properties, which needs to be explored with more experiments and research.

d—wellbore diameter, cm;

H—thickness of sublayers corresponding to M sampling points when processing micro-scanning imaging data, cm;

K—number of sublayers included in window W;

M—number of sampling points set when processing micro-scanning imaging data;

N—number of button electrodes of micro-scanning-imaging logging tool;

r_i—invasion radius of drilling fluid, cm;

R_eh—horizontal equivalent resistivity of formation corresponding to window length W, Ω•m;

R_ev—vertical equivalent resistivity of the formation corresponding to the window W, Ω·m;

R_h—resistivity value of electric anisotropic formation in the parallel bedding direction, Ω·m;

R_LA1, R_LA2, ..., R_LA5—apparent resistivity values of five measurement modes of array lateral logging, Ω·m;

R_t—true resistivity of electrically isotropic undisturbed formation, Ω·m;

R_th—true resistivity of undisturbed electrically anisotropic formation in parallel bedding direction, Ω·m;

R_ti—equivalent resistivity of the ith sublayer within window length W, Ω·m;

R_tpi—formation equivalent resistivity corresponding to the ith pair of polar plates in micro-scanning-imaging logging, Ω·m;

R_tpmax, R_tpmin—the maximum and minimum equivalent resistivity of different pairs of polar plates in micro-scanning-imaging logging, Ω·m;

R_V—the resistivity in the direction perpendicular to bedding of electric anisotropic formation, Ω·m;

R_xo—resistivity of flush zone of electric isotropic layer, Ω·m;

S—the number of polar plate pairs of micro-scanning imaging logging tool;

Δt_SF and Δt_SS—slowness of fast and slow S-waves, μs/m;

W—depth step selected each time when processing micro-scanning imaging data, m;

θ—angle between formation and instrument, °;

λ—resistivity anisotropy coefficient, shortened as electrical anisotropy coefficient, f;

λ_DANI—resistivity azimuthal anisotropy coefficient extracted from micro-scanning imaging logging, f;

λ_FMI—electrical anisotropy coefficient extracted from micro-scanning imaging, f;

λ_HRLT—electrical anisotropy coefficient obtained from inversion of array lateral resistivity logging, f;

λ_SLOANI—acoustic anisotropy coefficient calculated from time difference between fast and slow shear waves of cross dipole acoustic logging, f.