Discrete fracture modeling by integrating image logs, seismic attributes, and production data: a case study from Ilam and Sarvak Formations, Danan Oilfield, southwest of Iran

Understanding the fracture patterns of hydrocarbon reservoirs is vital in the Zagros area of southwest of Iran as they are strongly affected by the collision of the Arabian and Iranian plates. It is essential to evaluate both primary and secondary (fracture) porosity and permeability to understand the fluid dynamics of the reservoirs. In this study, we adopted an integrated workflow to assess the influence of various fracture sets on the heterogeneous carbonate reservoir rocks of the Cenomanian–Santonian Bangestan group, including Ilam and upper Sarvak Formations. For this purpose, a combination of field data was used including seismic data, core data, open-hole well-logs, petrophysical interpretations, and reservoir dynamic data. FMI interpretation revealed that a substantial amount of secondary porosity exists in the Ilam and Sarvak Formations. The upper interval of Sarvak 1-2 (3491 m to 3510 m), Sarvak 1-3 (3530 m to 3550 m), and the base of Sarvak 2-1 are the most fractured intervals in the formation. The dominant stress regime in the study area is a combination of compressional and strike-slip system featuring reverse faults with a NW–SE orientation. From the depositional setting point of view, mid-ramp and inner-ramp show a higher concentration of fractures compared to open marine environment. Fracture permeability was modeled iteratively to establish a realistic match with production log data. The results indicate that secondary permeability has a significant influence on the productivity of wells in the study area.


A f
Cross-sectional area of the fracture network E(n) Possibility density function that defines the fractures amount F ij Fracture tensor F kk Is the trace of the fracture tensor F ij K Permeability k ij Permeability tensor N Fractures amount n i n j Is the elements of a unit normal to the fracture n n ik n jk Are the elements of a unit normal to the fracture k R m Resistivity of the drilling mud R xo Resistivity of the invaded zone T k Fracture transmissivity V Mesh volume V t Total volume V tf Is fracture total volume W Fracture aperture

Introduction
Faults and fractures directly affect reservoir characteristics and crucially control fluid production rates in carbonate reservoirs (Hanafy et al. 2017;Davarpanah et al. 2018;Ghalenavi et al. 2020). Natural fractures are mechanical failures in rock formations, usually resulting from geological stresses (Olson and Pollard 1989;Meng and Pater 2011;Bisdom et al. 2017;Liao et al. 2020). The term "fracture" is virtually used for all generic rock discontinuities, including a wide variety of geological failures, like cleats, vugs, slickensides, joints, or faults (National Research Council 1996). In this regard, such naturally occurring fractures in hydrocarbon reservoirs significantly impacts their fluid dynamics by increasing their permeability and production volumes through wellbores (Liu et al. 2020;Rashid et al. 2020). Fractured reservoirs are common worldwide, but their complex heterogeneous nature poses significant development and production challenges for the oil and gas industry (Cacas et al. 2001;Al-Ahmadi and Wattenbarger 2011;Salimi and Ghalambor 2011;Abushaikha et al. 2015;Smeraglia et al. 2021). Fracture systems usually happen in swarms and are generally irregular, mostly disconnected. Hydrocarbon traps in the pore matrix of reservoir rock, but open fractures have a great share in its primary flow to the well (Esmaeili et al. 2022;Hosseinzadeh et al. 2020). The fractures can play both positive and negative roles in fluid flow, either as a flow path or a barrier (Nelson 2001;Aydin 2000;Shalaby and Islam 2017). In some zones, they can provide preferential flow channels for water rather than oil or gas, thereby adversely impacting hydrocarbon production. Fracture permeability is mostly much greater than the rock matrix's permeability, particularly in carbonates, that is because the ability of fluid to flow in naturally fractured reservoirs is inclined to be governed by the fracture network (Guerriero et al. 2013;Carey et al. 2020;Lv et al. 2021;Ceccato et al. 2021). It is the fundamental characteristics of fractures and expansive fracture networks that control the flow of hydrocarbons in naturally fractured reservoirs (Al-Rubaye et al. 2021;Tiab and Donaldson 2015;Voorn et al. 2015;Guerriero et al. 2013;Afşar et al. 2014;Zambrano et al. 2021). Consequently, it is important that these characteristics are well defined, and their impacts on fluid dynamics thoroughly understood as early as possible during field development for optimal reservoir flow performance and resource recovery to be achieved.
The general characteristics of a naturally fractured reservoir are as follows (Bratton et al. 2006;Zhang and Sheng 2020;Jiang et al. 2020).
• Matrix permeability is normally less than fracture permeability and is usually heterogeneous throughout the reservoir, • The rate of field production is highly reliant on a wellconnected, open natural fracture network, • Natural fractures substantially impact fluid flow rate and recovered resources. • The main aims objectives of this study are to: • Estimate fracture intensity using borehole image data • Determine in situ stresses orientations using image log data • Match fracture properties with production data using an optimizer • Relate fractured zones to their sedimentary environments For fractured reservoir modeling, it is essential to recognize the kinds of fractures and their interaction with each other and the host rock (Maulianda et al. 2019). Fractures may be caused by stress conditions. A stress regime may be established before, during and after the creation of the reservoir. The paleo-stress is the stress condition associated with the fracture initiation time (khodaei et al. 2021). However, the impact of present stresses on reservoirs is most important since they govern the degree of opening of fractures with specific orientations and thereby affect the permeability (khadivi et al. 2022). Formations can be fractured due to some reasons such as stress regimes, according to which the four main groups of fractures are tectonic fractures, regional fractures, stylolites, and surface-related fractures (Gudmundsson 1987;Ju and Sun 2016;Laubach et al. 2004;Nelson 2001;Azim 2016;Peacock et al. 2018).
Tectonic fractures are caused by local-to regionalscale folding and faulting (Volland and Kruhl et al. 2004;Guo et al. 2016). Regional fractures usually extend over very large distances (Laubach et al. 2004;Lu et al. 2015). Stylolite-related fractures occur in carbonate and carbonatecemented sandstone reservoirs due to the solution (Marfil et al. 2005;Padmanabhan et al. 2015;Toussaint et al. 2018). Surface-related fractures appear where a formation's surface bears weathering at weak areas by definition and are not present at depth (Zhang and Chai 2020).
Assessing and modeling fracture distribution requires inputs from field observations, core analysis, well-log interpretation, fluid production trends, fluid losses during drilling, and wellbore images (Tingay et al. 2009, Lacombe et al. 2011Nwabia and Leung 2020;Liu et al. 2021;Khoshbakht et al. 2012;Rajabi et al. 2010, Ezati et al. 2018Yarmohammadi et al. 2020;Ashena et al. 2021). Valuable information about fracture systems is also available from seismic data (Al-Dossary and Marfurt 2006;Chopra and Marfurt 2007), including amplitude variations and details of wave scattering and elasticity (Malehmir et al. 2009;Tsvankin 2012) and seismic attribute analysis (Karimpouli et al. 2013;Kosari et al. 2015Kosari et al. , 2017Hosseinzadeh et al. 2019). Such various sources of information are best studied collectively; addressing them individually introduces unnecessary bias and uncertainty into the interpretations. The data sources complement each other, but individual data sources have certain limitations and resolution constraints. Seismic data information is substantially influenced by the quality and resolution of the recorded data. Seismic data generally cannot resolve individual fractures but simply gives an indication of the presence and orientation of fracture sets (Pérez et al. 1999;Protasov et al. 2016;Pan et al. 2021). On the other hand, the basic suite of well-logs typically recorded post-drilling is not able to resolve fractures, identify whether fractures present are open, or distinguish fracture porosity from matrix porosity (Mohebbi et al. 2007;Bahrami et al. 2012;Sun et al. 2014;Mazaheri et al. 2015;Aljuboori et al. 2022). Whole wellbore cores, if recovery is high, tend to provide good material for fracture identification and measurements (Khoshbakht et al. 2012). However, sidewall core plugs are more difficult to interpret and extrapolate from as they may or may not represent the zone sampled across its full vertical extent (Fornero et al. 2019;Stark et al. 2021). Moreover, fracture orientation information can usually not be resolved accurately from wellbore cores, and the collection of whole wellbore cores is costly and time-consuming during drilling (Zazoun 2013). Wellbore-derived data sources (i.e., cores and basic well-logs) are unable to provide, on their own, reliable interpretations on a three-dimensional spatial field scale. Borehole image logs are useful for determining some fracture properties, such as accurate orientations of the fractures, information on whether fractures are open or closed in the subsurface, as well as information on fracture apertures (Ozkaya and Mattner 2003;Aghli et al. 2016;Lai et al. 2017). However, these images do have some limitations. They are typically not able to resolve all the fractures visible in cores. Also, they do not provide details of fracture length and can only give indications of fracture height with substantial associated uncertainty (Olson et al. 2001;Özkaya 2003;Wilson et al. 2015;Sarkheil et al. 2021).
The stress regime of Zagros region is active due to the collision of Arabian plate and the Eurasia plate (Bordenave and Hegre 2005;Alavi 2007;Ghanadian et al. 2017;Mohajjel and Fergusson 2014). As a consequence, many of the oil fields in this basin have experienced faulting and fracturing over geological time. Wennberg et al. (2006) established that field studies can expressively enhance the recognition of the fracture propagation and its effect on permeability in hydrocarbon reservoirs. Specifically, outcrop studies can provide sufficient information about the relationship between geological deposition and fracture characteristics. This has been studied in a comprehensive geological assessment of the Asmari Formation on the northeast limb of the Khaviz Anticline in the Zagros in southwest Iran (Wennberg et al. 2006). The results showed that mud-supported deposits have higher fracture intensity than grain-supported textures. Pireh et al. (2014) analyzed naturally occurred fractures and the impact of deformation intensity on fracture density in the Garau Formation within two anticlines of the Zagros fold and thrust belt, Iran. They conducted an extensive study of outcrops of the Garau Formation in the Kabir Kuh area, identifying dominant fracture orientations of NW-SE and NE-SW. Nasseri et al. (2015) proposed an approach for improving the recognition and interpretation of faults and fractures exactly. By applying some fault detection attributes and pattern recognition techniques, they were able to improve faults and fracture identification. Their method applies the ant colony (AC) optimization algorithm to a seismic dataset from Dezful Embayment in the southwest region of Iran. However, this method cannot be applied to fracture detection for low-quality seismic data throughout the studied field. Yousefi et al. (2019) used image logs to analyze the fold and fault-related fracture network of the Asmari Formation in the Rag Sefid Anticline, southwest Iran. Although able to determine fracture orientation, this method was unable to establish fracture aperture, fracture porosity, and fracture permeability in the studied field. Aghli et al. (2020) studied reservoir heterogeneity and fracture parameter determination using electrical image logs and petrophysical data in carbonates of the Asmari Formation in the Zagros Basin, SW Iran. They evaluated FMI log data to determine fracture parameters, including fracture aperture, porosity, and permeability. However, they did not generate a 3D model to display these parameters across the field. Fracture modeling in this complex structural basin reveals how faults and fractures impact reservoir properties and production rates (Zhao et al. 2014;Xiaoxia et al. 2021).
The approach introduced in this study involves using multiple assessment scales, ranging from millimeters to kilometers, by combining drilling and core data, well-log petrophysical data (basic well-log suite as well as well-log imaging), and seismic attributes. This integrated approach provides a more representative and precise distribution model for fractures throughout the reservoir rocks (Ilam and Upper Sarvak) of the Bangestan group and the interaction between fractures and the matrix of these carbonate reservoirs. The results provide improved knowledge and understanding in relation to the studied field and its regional setting. Moreover, useful information is provided regarding the tectonic evolution of the field, the azimuth of horizontal stresses, and the propagation of the fractures with respect to the field's stress regime. This makes it possible to calculate fracture intensity using image log data and to quantify fracture reservoir properties, including porosity and permeability.

Tectonic framework
The Zagros fold and thrust belt (ZFTB) is an NW-SE extended orogenic belt covering about 2000 km from Turkey through southwest Iran (Fig. 1). It is restricted to the north by the Main Zagros Fault, defined as the suture zone of the Neo-Tethys Ocean (Braud and Ricou 1971). The Zagros belt is created due to the collision between the continental Arabian and Eurasia plates (Takin 1972;Ghasemi and Talbot 2006). The first compressive movements occurred during the Late Cretaceous because of the obduction of ophiolites on the northeastern border of the Arabian continent. These activities got more extended following the continent-continent collision in Miocene time (Sherkati and Letouzey 2004). According to its surface structure, the ZFTB can be separated into the High Zagros Belt (HZB) and the Zagros Simply Folded Belt (ZSFB) (Alavi 2007). The Iranian ZSFB consists of different tectonostratigraphic areas based on areal facies properties: the Fars province or eastern Zagros; the Izeh zone; the Dezful Embayment or central Zagros (Fig. 1); and the Lurestan province or western Zagros. To the north and northeast, the Dezful Embayment is bounded by the Mountain Front Fault (MFF), a major structural and topographic boundary in the Zagros belt. In the north, the Dezful Embayment is separated from the Lurestan Salient by the Balarud Fault Zone, a transverse zone of the MFF with a sinistral strike-slip component. The Dezful Embayment is restricted by the Zagros Foredeep Fault or Zagros Frontal Fault (ZFF) in the south, which is located just south of the studied field (Sepehr et al. 2006;Sepehr and Cosgrove 2004) (Fig. 1).
The Zagros basin's sedimentary depositions comprise formations from Cambrian to Recent. The geological evidence suggests that this area was part of a passive continental margin, which subsequently experienced rifting in the Permian-Triassic, and collision in the Late Tertiary (Lalami et al. 2020).

Location and stratigraphy of the field
The Danan field is structurally located within the Dezful Embayment. The Dezful Embayment shows a foreland basin where sinking at the uplifting Mountain Front Fault' downhill has caused the deposition of thick Upper Miocene-Pliocene Aghajari and Bakhtiyari Formations (Farahzadi et al. 2019). In this region, the evaporates of the Miocene Gachsaran Formation play an excellent caprock role ‫ه‬in the Asmari reservoir. Generally, the subsurface evidence suggests low amplitude and large wavelength anticlines buried beneath the Gachsaran Formation.
The Danan anticline field has one thrust fault with a back thrust penetrating the evaporitic Gachsaran Formation and two main and five minor faults with the NW-SE trend transecting the reservoir zones. The main reservoir rocks of the field are shallow-marine carbonates of the middle Cretaceous Ilam and Sarvak Formations at a depth of approximately 3100 m. Also, the Oligo-Miocene Asmari Formation, located at approximately 2200 m depth, is another reservoir rock of the field (Bordenave and Hegre 2010).
The Ilam Formation, with a thickness of 100 m, consists mainly of off-white, yellowish brown, light gray, brownishgray argillaceous limestone with some gray to dark gray shale interbeds in its lower section. The Sarvak Formation, with a thickness of 750 m, consists mainly of off-white to light-cream, gray, brownish-gray limestone. It displays some oil staining, pyritic nodules, and fossiliferous bands. Thin layers of light-greenish gray to greenish gray and soft shale occur in the upper Sarvak Formation (modified after James and Wynd 1965) (Fig. 2).

Analyzed data
The availability of fracture-related data from diverse data sources improves the reliability of the fracture model, which can reduce uncertainty in its interpretation. The data analyzed in this study are as follows.
• Conventional basic suite of well-logs. Conventional well-logs were available from 9 drilled vertical wells throughout the field, and their evaluation results for porosity, permeability, and lithology, which were used in creating both structural and property models for the field. • Borehole image (BHI) log. This log was used for fracture interpretation from well DA-07. This well is located  James and Wynd 1965) in the southern part of the Danan field (Fig. 3). It was drilled in 2009 with a small deviation (maximum deviation of 1 degree) through the Ilam and Sarvak Formations. Formation micro imager logs (FMI) provide an image from the borehole wall based on its resistivity property in water-based mud drilling with a resolution in the order of millimeters. Introduced by Schlumberger, they are now used broadly to distinguish geological features of the boreholes such as unconformities, type of fractures, reservoir properties, bedding plane dimensions, and stress conditions. Substantial natural fracture analysis in wellbores can be achieved with data collected by bore hole image logs. • Seismic data. Geophysical data, including 2D seismic amplitudes, were compiled. The quality of the available 2D seismic profiles (12 PSTM 2D lines with a total length of 234 km acquired in 2004) was deemed adequate for structural interpretation of the field. Figure 4 displays the seismic base map with the studied well location and available 2D seismic lines superimposed on surface topography. • Production data. A production log was run from the production services platform (PSP) soon after completion in July 2009 to evaluate the production profile of the open-hole sections. A summary of zonal oil production is shown in Table 1.
The data listed above are collectively critical for developing a reliable characterization of the fracture network in the studied field.
In order to analyze the structural evolution of the Danan anticline, different geological elements, including regional tectonics, have been taken into consideration. Also, to evaluate the structural growth of the field, a "flattening" technique has been applied to seismic lines with several horizons from deep to shallow being flattened. This approach makes it possible to follow the tectonic evolution in space and time. Flattened seismic sections at the horizons of interest are presented together with their interpreted structural aspects. These sections allow tracing, step by step, the structural growth from tectonic quiescence to ongoing activities.

Spatial organization of fractures
Diffuse fractures and fracture swarms are the two main groups of naturally occurring fractures.
Diffuse fractures are spread-out in the rock formation that tend to be concentrated in deposits with low ductility. These more brittle deposits may comprise one or more fracture sets with diverse spatial appearances. Diffuse fracturing can be identified by (1) high fracture intensities concentrated in significant formation beds and (2) considerable diffusion Fig. 3 The seismic base map for the Danan field area including studied well location and available 2D seismic lines superimposed on surface topography in fracture orientations stemming from the existence of some fracture sets in a similar bed.
Fracture swarms, on the other hand, are zones characterized by high fracture intensity and are typically attributed to the presence of faults and tend to be layer bound. Faults of all types tend to be associated with high fracture intensity in their vicinity. In fact, fracture swarms can be recognized as narrow zones with much higher fracture intensity than those present in the surrounding rocks. Fractures can be recognized as fracture swarms if (1) their emergence is not linked to the presence of geological deposits, (2) the dispersion of fracture orientation is small, and (3) the main fracture orientation follows the strike direction of neighboring faults.  Figure 6 shows the workflow applied in this research for discrete fracture network modeling.

Open fractures aperture
Fracture aperture can be computed for all conductive/ open fractures penetrated by wellbores using BHI log data. Its measurements are typically shown in centimeters in a  composite plot. Fracture aperture can be measured by Eq. 1. It shows that the fracture aperture is dependent on the resistivity of the drilling mud (R m ) and the formation resistivity of the invaded zone (R xo ) (Luthi and Souhaite 1990).
where W represents fracture aperture, A is access current associated with the FMI electrode; Rm is the mud resistivity; Rxo is the invaded zone resistivity; c and b are constants. Values for constants b and c were determined as 0.6 and 3.5, respectively, by calibrating with well-test data.

In situ stress assessment
Continental crust rarely remains at a hydrostatic stress state for extended periods. That hydrostatic stress condition refers to points in the subsurface where stresses are equal in all directions (σ1 = σ2 = σ3). However, stress situations are unusual in the subsurface due to the various local structural and/or regional tectonic movements that regularly occur and modify the prevailing stress regimes over time. Some of the more typical uneven stress states are a consequence of plate tectonic dynamics. These lead to the construction of regional stress patterns typically constrained within the boundaries of specific tectonic plates. Conversely, in some cases, regional tectonic stress regimes can be almost entirely overprinted because of localized movements and the stress conditions . The direction and magnitude of these local stress regimes can vary substantially and unexpectedly over small distances in any area (Zoback et al. 1985).
In this study, the breakouts' orientation was calculated according to compressively stressed zones observed in the FMI log in the well-07, which provides full-bore images.

Fracture permeability
Fracture permeability can be determined by considering collectively the fracture intensity, fracture interconnectedness, and fracture transmissivity data (Lee et al. 2011;Yuan et al. 2015). Oda (1985) presented a method for approximating the fracture permeability solution. The direction of each fracture in a grid cell is stated as a unit normal vector n in Oda's method. Integrating the fractures over all of the unit normal N, Oda gained a tensor N ij unfolding the mass moment of inertia of the normal fracture spread over a unit sphere (Eq. 2) (Golder 2012): where N is fractures amount in Ω , n i n j is the elements of a unit normal to the fracture n, E(n) is the possibility density Adopted workflow for studying the impact of fracture on production data on heterogeneous Ilam and Sarvak Formations, SW Iran function that defines the fractures amount whose unit vectors n are located within a small solid angle dΩ , and Ω is the entire solid angle equivalent to the surface of a unit sphere.
Discrete fractures for a sole mesh are biased by their surface and transmissivity to calculate an experimental fracture tensor (Eq. 3).
where F ij is the fracture tensor, V is the mesh volume, N is the whole quantity of fractures in the mesh, A k indicates the area of the fracture k , T k is the fracture transmissivity k , and n ik n jk are the elements of a unit normal to the fracture k.
By assuming that F ij indicates fracture flow as a vector along the fracture's unit normal, Oda's permeability tensor is achieved from F ij . However, F ij must be rotated into the planes of permeability if fractures are closed in an orientation corresponding to their unit normal (Eq. 4).
where k ij shows the permeability tensor, F kk is the trace of the fracture tensor F ij , and ij is Kronecker's delta.
The advantage of Oda's solution is that it may be determined without the involvement of flow simulations. However, it overestimates permeability and connectivity.

Fracture porosity
Fracture porosity is recognized as the fraction of overall fracture volume to unit volume. Fracture porosity is solely determined by fracture geometry and is computed based on fracture area for each mesh volume and fracture aperture (Eq. 5), where fracture area basically stems from fracture size and height.
where f is the fracture porosity, V tf is fracture total volume, and V t is total volume (m 3 ). Clearly, fracture porosity is subject to the fracture intensity, size, and aperture.

Faults detection using seismic data
To identify potential locations of intense fracturing, it is necessary to identify the locations and orientations of minor faults and fractures spatially across areas of interest. The approach applied in this study is to highlight structural discontinuities and distinguish faults with the aid of a suite of seismic attributes (Sigismondi and Soldo 2003;McDonnell (3) al. 2007;Marfurt and Alves 2015); specifically, amplitude, variance, dip, and semblance. Coherency is a seismic attribute that identifies variations in adjacent wave traces. Seismic signals are a consequence of wavelet responses to the geology they penetrate, thereby capturing abrupt structural variations in the subsurface. Therefore, coherency helps to detect faults associated with reflector displacements. Algorithms used to establish coherency include variance and semblance. Variance is a statistical measure of the spread of values across a distribution (Fisher 1922). Variance values distributions calculated for seismic signatures make it possible to delineate fault locations. Dip attributes are also useful for evaluating the curvature of dip-steered features such as faults (Basir et al. 2013).

Natural fracture characterization via image log interpretations
Substantial natural fracture analysis in wellbores can be achieved with data collected by borehole image logs. To extract quantitative information on conductive/open fracture characteristics such as aperture, the full-bore formation micro imager (FMI) log for well DA-07 was scaled with a shallow resistivity log and neutron log.

Classification of natural fractures
In the studied field, the image log-derived data, which have sufficient coverage throughout the Ilam and Sarvak Formations, were used to group fractures into fracture sets. Two types of natural fractures, including conductive and resistive fractures, based on their appearance on the wellbore images and continuity of their trace across borehole DA-07 were identified. Conductive fracture is identifiable as a sinusoidal shape with a dark color on FMI data, which indicates that the fracture is open and filled by mud in the borehole. However, the resistive fracture appears as a sinusoidal shape with bright color, which indicates the fracture is closed and is not conductive to fluid or mud in the wellbore. In addition, production services platform (PSP) data indicate clearly that fracturing is confined to the three uppermost reservoir zones.
In total, 282 open and conductive fractures were identified, which are further divided into three classes based on their nature and continuity within the wellbore. These three types are 1) major open and conductive, 2) medium open and conductive, and 3) minor open and conductive fractures (Fig. 7). Based on the interpreted image log, only two major open/conductive fractures existed, which seemingly have a larger aperture and clearly extend through the wellbore. Also, 43 medium open/conductive fractures by less apparent aperture traces across the wellbore are distinguished. Moreover, 237 minor open/conductive fractures are observed over the logged interval. Such fractures display limited or no aperture traces across the wellbore compared to the rest of the natural conductive/open fractures identified.
Only three clear indications of resistive fractures are observed over the interpreted interval. These are classified as resistive (closed) fractures (Fig. 8). (3346-3427 m; Fig. 9) The natural fractures identified in Ilam Formation consist of one major, 17 medium, and 100 minor conductive/open fractures. Their dominant dip azimuth direction is N45W with a corresponding strike of N45E-S45W. They have a mean dip inclination of 72 degrees. The relationship between the strike of these fractures and bedding plane orientations (with a strike of N65E-S65W) over the same zone shows that the identified fractures are slanted to the bedding plane (Figs. 6 and 10). Fig. 10) In total, one major, 26 medium, and 137 minor conductive/open frac-tures are observed within the Sarvak interval. They show a dominant strike of N25E-S25W and a 73-degree average dip inclination. It gives the impression that these fractures are transverse to oblique, corresponding to the predominant east-west strike direction of the bedding planes across this interval. Table 2 shows natural fracture orientations in every zone identified in the studied well. The statistics details of these fractures are provided in Fig. 11.

Conductive/open natural fractures aperture
The appearance of the three different types of conductive fractures (major, medium, and minor fractures) on the FMI log is shown in Fig. 7. The detection of these three fracture types is based on the degree of opening size of the distinguished fractures, which are displayed as dark sinusoid features on the FMI log. Figure 11 (A) provides the varieties of aperture size for the conductive natural fractures identified on the BHI images for the Ilam and Sarvak reservoirs. Most of the identified natural conductive fractures are interpreted through Ilam 1, Sarvak 1-2, Sarvak 1-3, and Sarvak 2-1  (Table 1), these intervals have the main contribution in the production of the field.

In situ stress measurement
There are essentially two types of borehole failure, related to either tensile or shear failures, caused by the unbalanced stress systems present in the drilled wells. The drilled formations are able to withstand certain levels of compressive and shear stresses due to the replacement of drilled rock cuttings with the drilling mud. However, borehole fluids can only withstand compressive stresses. The radial stress experienced at the wellbore walls represents the drilling fluid's weight which is less than the prevailing pore pressure. This leads to the concentration of stress around the borehole walls, known as hoop stress. In low mud-weight conditions, the maximum hoop stresses tend to be substantially higher than the prevailing radial stresses (Fig. 12). As a result, the borehole walls are able to bear shear failures (borehole breakout), which are presented as elongated dark areas on the FMI images that are positioned 180 degrees from each other, and these failures are perpendicular to the maximum horizontal stress (σ Hmax ). Conversely, in the case of high mud-weight conditions, there is a growth in radial stress and a decline in hoop stress. In such circumstances, tension overcomes the strength of the wellbore walls and causes failure in rock formations adjacent to the wellbore. These failures are termed drilling-induced fractures and can be seen in the BHI logs (Fig. 13).
Drilling-induced fractures are consistently observed along the N25E-S25W direction of the wellbore that can be recognized as the maximum horizontal in situ stress direction (σ Hmax ) around well DA-07 (Fig. 13). This orientation of in situ stress around the study well matches the regional stress orientation of the Zagros fold belt.
Borehole breakouts occur because of the impacts of high stresses around the borehole. These features are distinguishable as wide vertical fractures around the wellbore on the image logs. Some borehole breakouts are observed on the FMI images (Fig. 14). The trend of these borehole breakouts, N45W-S45E, is almost perpendicular to drilling-induced fractures (Fig. 15). This also shows that the minimum horizontal in situ stress (σ hmin ) direction around well DA-07 is N45W-S45E. Comparing the fracture orientations and implied minimum and maximum horizontal stresses with the known Arabian Peninsula regional stress regimes (NE-SW for σ Hmax and NW-SE for σ hmin ) (Perotti et al. 2011) reveals that the two are similar (Fig. 15).

Identification of major and minor faults using seismic data
Flattening on Sarvak to Asmari horizons shows that apparently no active movements and seemingly no obvious tectonic activity occurred, except for periods of some minor movements probably related to basement activities (Fig. 16A). On the contrary, the flattening of post-Asmari to Gachsaran (Gs-7) Formations effectively shows that structural growth was triggered during the Miocene period (Fig. 16B). This would be coeval with regional tectonics of the Zagros (corresponding to the second Alpine phase). Additionally, thrust faults, which are a consequence of the compressional regime associated with basement-related faults as well as décollement levels, contributed to folding, and affected the structural geometry.
Available seismic data are used to determine the presence and orientation of major and minor faults/fractures within and surrounding the studied field. Various seismic attributes are available to identify potential faults. These attributes include cosine of phase, frequency, and phase shift, and they were analyzed in the seismic data available throughout the field. Based on these attributes, two sets of faults can be distinguished, some of them associated with major offsets, and their displacements are visually obvious, while others are associated with less visible minor displacements. The seismic dataset available enabled us to analyze a grid of Fig. 9 Dip angle, Dip azimuth, fracture intensity, and fracture cumulative logs in DA-07 for the Ilam Formation (CGR is the summation of potassium and thorium, and SGR is summation of thorium, potassium, and uranium); (WCS is lithology log of shale) nine strike lines and three dip lines throughout the field. The faults were identified by considering all these lines.
Fault interpretation reveals that the field structure is transected by a relatively long-wrench fault-oriented NW-SE (Fig. 17). This major fault is possibly deep-rooted (basement related), creating grabens with the same orientation into which Oligocene-Miocene layers have been downthrown (Fig. 17). However, this interpretation is somewhat uncertain due to the low resolution at the depth of the 2D seismic data available.
Based on Fig. 17, two main faults have influenced the formation of Danan anticline (green color) and, subsequently, have shaped the actual geometry of this field. One of these faults caused the Danan Oilfield to be separated from its neighborhood structures. Taking into account the orientation of this fault, it might be considered a branch of the well-known Balarud fault. The second reverse fault has been active in the Gachsaran Formation. Another fault, a priori with a different origin, has also had an influence; that is a wrench fault causing a minor vertical offset in the core of the structure. For two reasons, this fault is interpreted in the strike-slip category (light blue color). Firstly, this appears to be related to a compressional regime and is, therefore, not related to normal fault activity. Secondly, it has been sourced from the deeper part of the sedimentary column or may be attributed to a basement-related fault. The detailed characteristics of the faults (yellow color) influencing the geological interval (cap rock/Asmari horizons to Sarvak horizon) cannot be distinguished in detail because of the low quality of the seismic data. Figure 17 illustrates the status of the identified faults, including their NW-SE azimuth, and distinguishes major and minor faults.

Discrete fracture network (DFN) modeling
Where the porosity and permeability of fractures have a great share in the reservoir's characterization, the DFN modeling approach can be usefully applied to generate the static, geomechanics, and dynamic field reservoir models. This approach relies heavily on wellbore image log information and treats fractures as individual objects. Continuous fracture intensity modeling (CFM) is used as a basis for developing the DFN modeling, leading to the definition of discrete fracture planes (Figs. 18 and 19). The construction of reliable DFN models requires a representative number of BHI logs (i.e., from several wells) that can be used to provide a detailed definition and analysis of the in situ fracture patterns. However, an indicative DFN can be constructed based on general geological considerations but more limited well data. The indicative DFN provides useful provisional visualization of the likely fracture patterns impacting specific field structures (e.g., studied field with DFN based on image logs from well DA-07 as well as regional stress field and fault distributions).

Fracture intensity
One of the most vital parameters for DFN modeling is fracture intensity. It can be assessed in one, two, or three dimensions, according to Dershowitz (1984). Field measurements of fracture intensity in the formation, on the other hand, are normally determined in one dimension along a sampling wellbore. The fractures amount per unit length Fig. 13 Examples of drilling-induced fractures on the FMI images in the formation is, therefore, usually represented by onedimensional fracture intensity measurements (Dershowitz 1984;Pei et al. 2020). The use of these measures increases the accuracy of discrete fracture investigation and modeling greatly. Fracture intensity was calculated using image logs in this study. Based on the measured fracture intensity from image data, the concentration of the distributed fracture set for the Ilam-1 zone is very high, whereas that for the Sarvak-2 zone shows only a low concentration (Fig. 18).

Fracture orientation
Fracture orientation has also an important contribution to DFN models because it controls directional permeability and hence contributes to reservoir fluid flow. The inclination of the fracture surface from the horizontal represents the dip, and the vector normal to the fracture plane is referred to as the pole. Information regarding fracture orientation can be gleaned from various datasets. Wellbore images were used to assess the dip and azimuth of fractures intersecting the wellbore in this study (Figs. 7 and 8). Generally, the trend of all the natural fractures was identified to be mostly perpendicular to the faults (NE-SW) (Fig. 19).
In natural fracture prediction, tectonic modeling can provide some useful insight. Tectonic models can be configured to capture a series of one or more tectonic events, and fractures are the cumulative result of one or more successive tectonic events. To generate a realistic tectonic model, as a starting point, it is necessary to have detailed knowledge of the tectonic regime(s) and their corresponding relative stress magnitudes. Fractures are characterized spatially by their positions in the well (X, Y, Z, MD) and their orientations (dip angles and dip azimuths). Also, the geological fracture type needs to be specified, and that information comes from the well-log analysis already completed. The results of the tectonic modeling conducted to confirm that the dominant stresses for the field area belong to a normal regime in the reservoir section.
For DFN modeling, it was assumed that open and conductive fractures would be the dominant type of fractures, since  (Fig. 19). An orientation map can be used as input for DFN modeling, calculated based on the dip azimuth of the top reservoir surface.

Fracture aperture
The fracture aperture is the fracture's perpendicular width. The FMI log data were used to calculate fracture aperture in this study. The fracture aperture distribution follows a lognormal distribution, which is consistent with previous research (Keller 1996;Gale 1987;Olsen et al. 2007;Liu et al. 2021). The minimum and maximum values of fracture aperture in the study area are approximately 0.01 m and 1 m, respectively (Fig. 11).

Fracture size
Since it is impossible to measure the fracture size directly with any borehole tool, fracture size is perhaps the least known parameter. This is due to the fact that a wellbore typically intersects only part of each fracture. This makes the determination of fracture size from image logs quite challenging. Examining outcrop analogs can often provide a useful alternative source of fracture size information. In addition, dynamic data can offer a common sense about the fracture length (Li et al. 2013).
The relationship between fracture size and fracture aperture is well known. Since the exact practical relationship is unknown, a linear relationship is most likely (Ozkaya and Mattner 2003). A linear relationship assumption can be applied to determine fracture size in cases where a fracture aperture distribution is available. The fracture aperture is likely to be in the order of 1000 times smaller than the fracture length. In this study, that assumption suggests that the fracture length varies between 20 and 300 m based on Shows tectonic evolution of the field. A Flattened at Sarvak level. B Flattened at GS4 (fourth member of Gachsaran Formation, which play a seal role for Asmari reservoir) level measured fracture aperture from FMI data. Iterative facture property modeling has then been applied to establish the optimum value for fracture length.

Fracture porosity and permeability
Equations 2-5 are used to generate fracture porosity and permeability. According to Eq. 5, fracture porosity relies on fracture intensity, size, and aperture. The fracture intensity and aperture have been calculated from FMI data. The fracture size has been estimated by applying an iterative calculation to match production flow data over time. In this regard, the total volume of the fracture ( V tf ) is derived as ∑ � A f e � , which can be obtained using the aforementioned parameters derived from FMI data. In addition, the total volume ( V t ) value can be calculated based on the grid size assessed. Figure 20 shows the modeled fracture porosity throughout the Ilam and Sarvak Formations.
Fracture permeability is provided by Eq. 4, in which F ij is related to fracture geometry and is mostly the average of the void volume of the fractures weighted by the volume of the mesh, as expressed using Eq. 3. The Oda permeability generated represents the effective fracture permeability (Figs. 20,21). Based on the generated model, the Ilam 2 and Sarvak 1-2 zones are revealed as possessing high values of fracture porosity and permeability. These characteristics are in line with the fracture intensity and fracture aperture information obtained from image log data (Eq. 1) and correlated with production flow data.

Discussion
Fractures grow in formations at different periods and due to different reasons. The type of fractures depends upon the situation of a specific rock and its surrounding environment. As the FMI log is unable to measure the length of fractures, and also the resolution of the 2D seismic data was not good enough to use for that purpose, we draw on published analysis to provide insight to fracture lengths in the region of the study area. Such available analysis is mostly based on outcrop studies (Pireh et al. 2015;Aghli et al. 2020;Wennberg et al. 2006Wennberg et al. , 2007Nasseri et al. 2015;Yousefi et al. 2019). Additionally, production data can offer good guidance about the fractures' length (Li et al. 2013). Fracture length is important in calculating fracture permeability. As permeability has a key role in a reservoir's fluid production potential, it is useful to estimate the optimum length of the fractures by iteratively adjusting the discrete fracture model to match the dynamic PSP data ( Table 1). The fracture length used to generate DFN models for the Ilam and Sarvak Formations is provided ( Table 2).
The results suggest that fracture-dependent reservoir quality in the upper interval of the Ilam Formation leads to more suitable hydrocarbon fluid production compared to the lower interval of the Ilam Formation.
In the studied area, the upper interval of Sarvak 1-2 (3491 m to 3510 m), Sarvak 1-3 (3530 m to 3550 m), and There are no fracture characteristic contrasts between the Ilam and Sarvak Formations. Also, FMI analysis of sedimentary features, including bedding, fracture types, and orientations, indicate that the Ilam and Sarvak Formations throughout the studied well have behaved as a single structural unit. Moreover, all faults detected on seismic data transect the Ilam and Sarvak Formations simultaneously. Consequently, it is apparent that these two adjacent formations have related sedimentary and tectonic histories and their fracture characteristic are the same.
Fracture permeability generally plays a greater role in fluid flow through a reservoir than matrix permeability. Fracture permeability can be measured using production data, like PSP. This approach calculates transmissibility, by applying the permeability and the aperture, which provides an indication of flow capacity. PSPs record cumulative oil production across a specified depth interval influenced by all the fractures (and matrix porosity) connected to the wellbore in that interval. PSP data are used to tune fracture parameters, in particular fracture lengths, in order to improve fracture property estimates, including fracture porosity and permeability (Table 3). Fracture intensity has a substantial impact on fracture permeability because it determines the interconnectivity between matrix porosity and fractures, thereby influencing fluid transmissibility (Golder 2012). In this study, the fracture intensity is calculated using image log data and propagating it by applying a deterministic model (Kriging method). DFN model outputs, including  Figure 20 reveals that the Ilam 1 zone displays higher average fracture permeability than Ilam 2 zone. Additionally, the Sarvak 1-2 zone displays higher fracture permeability than the Sarvak 1-3 zone, which according to production data (Table 1), these intervals have the main contribution to the production of the field.
The Ilam and Sarvak Formations of the Upper Cretaceous Bangestan group have been deposited on a homoclinal carbonate ramp platform. Schematic depositional models of   Fig. 22. These formations have strongly been affected and deformed by a Miocene tectonic events (Rajabi et al. 2010). The Ilam Formation is classified into Ilam 1 (upper zone) and Ilam 2 (lower zone), and Sarvak interval includes Sarvak 1-1, Sarvak 1-2, Sarvak 1-3, Sarvak 2-1, and the Sarvak 2-2 reservoir zones from top to the base of the formation. According to Fig. 22A, the most fractured interval of Ilam and Sarvak Formations deposited in the inner-ramp and middle-ramp depositional environment, and the outer-ramp of the both formation lacks of significant amount of fracture.

Conclusions
In this study, an integrated geological modeling of the Ilam and Sarvak Formations was conducted integrating seismic data, core data, conventional and advanced open-hole logs, petrophysical and dynamic data to provide valuable insight to reservoir properties and fluid flow characteristics. Considering the results of this study, followings are concluded.
• Interpreted FMI logs and secondary porosity indicators calculated from open-hole logs reveal that a substantial amount of secondary porosity exists in the reservoir zones of the Ilam and Sarvak Formations. • Stress regime investigation of both formations, based on maximum horizontal stress orientation and stress ratio, shows that the main stress regimes impacting the field are reverse faulting with a strike-slip component. • Integration of the modeled secondary permeability with interpreted production log data shows that the vertical fluid inflow profile for the wells correlates with secondary permeability. It is concluded that secondary permeability exerts a significant influence over the productivity of wells in this field. • Characterization of the Ilam and Sarvak reservoir properties indicates that the Ilam 1, Sarvak 1-1 to Sarvak 1-3, and the base of Sarvak 2-1 zones possess the most favorable properties and contribute most to fluid production.
Funding This research received no specific grant from any funding agency in the public, commercial, or not-for-profit sectors.

Conflict of interest
On behalf of all the co-authors, the corresponding author states that there is no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.