Investigating the Influences of Pore-Scale Characteristics on Tight Oil Migration by a Two-Phase Pore Network Model

State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum, Beijing 102249, China College of Geosciences, China University of Petroleum, Beijing 102249, China The Key Laboratory of Gas Hydrate, Ministry of Natural Resources, Qingdao Institute of Marine Geology, Qingdao 266071, China Laboratory for Marine Mineral Resources, Pilot National Laboratory for Marine Science and Technology, Qingdao 266071, China Institute of Petroleum Engineering, Heriot-Watt University, Edinburgh EH14 4AS, UK


Introduction
Following the shale gas breakthrough in the United States, tight oil has become the next focus of unconventional petroleum exploration worldwide [1,2]. Tight oil from the Bakken Formation has accounted for more than 10% of the United States daily oil production in 2017 [3]. Unlike conventional reservoirs, tight oil reservoirs are usually characterized by low porosity, low permeability, complicated pore-throat structures, and strong heterogeneity ( [4][5][6]). The hydrocarbon migration is a process in which hydrocarbon expelled from low-permeability source rocks finds tortuous porous paths through carrier beds into the traps where hydrocarbon can accumulate to form reservoirs in the geological time-frame [7]. Of many factors that control this process and therefore the resultant reservoirs in terms of their sizes, fluid saturations, and fluid properties within them, which are of crucial importance for hydrocarbon exploration and exploitation [8,9], the tightness of pore space within carrier and reservoir rocks is of pivotal importance. Compared with conventional reservoirs, oil migration behaves distinctively in tight reservoirs in the following aspects [5,10]: (1) there is no clear demarcation between the primary migration and secondary migration process, (2) the migration distance is typically short which induces hydrocarbon inner source accumulation or near-source accumulation, and (3) the impact of capillary pressure is significantly important during the migration process while the effect of buoyancy is limited.

Relationships between Pore-Scale Characteristics and Oil
Migration. The significance of pore-scale characteristics for the oil migration process has been assessed using experimental and geological methods in previous literature. Lai et al. [11] analyzed the pore-structure characteristics of tight sandstones in the Sichuan Basin, China, by using experimental data and basin-scale inferential methods. They concluded that the complex pore systems in tight sandstones are associated with the original depositional environments and subsequent diagenesis alteration; the coexistence of different types of pores results in the polymodal distribution of pore sizes. Relationships between the pore-structure characteristics and diagenesis alteration effects are concluded in their published review paper [12]. Similar work was reported by Xiong et al. [13]. They classified the pore system of tight sandstone into four types: intergranular pores, dissolution pores, clay-aggregated pores, and microfractures, the sizes of which are distributed in the different ranges varying from nano-to microscale. The flow potential during the oil migration process and the oil storage capacity in terms of oil migration and accumulation are contributed by the pores with different sizes. Cao et al. [14] further identified the key geological controlling factors of microscale oil distribution during oil migration and accumulation by geochemical approaches combined with mesoscale petrophysical methods, which includes porosity, pore throat radius, and hydrocarbon generation capacity. Zhang et al. [15] conducted and analyzed hydrocarbon charging experiments using natural tight sandstones samples. The experimental results show that the starting pressure of oil migration and terminal oil saturation values both show good logarithmic relationships with the sample permeability. However, the growth curves of oil saturation (S o ) with charging pressure are irregular and distinctive for different samples, such characteristics cannot be explained without sufficient pore-scale information of rock samples. To get a better understanding of the oil migration critical conditions and controlling factors at the pore scale, further investigations were completed in our previous works [16][17][18] by using a series of microscopic observation methods including CT (X-ray computed tomography), SEM (scanning electron microscopy), NMR (nuclear magnetic resonance), and casting thin section techniques. Several concluding remarks obtained from these works are as follows: (1) there exists a critical pore radius threshold for oil migration in tight cores, which ranges from 0.035 to 0.627 μm; (2) the oil-charging curves of tight reservoirs could be classified into four types, where the porestructure factors play crucial roles; and (3) the surface areas of effective pores/throats, pore quantity density, and coordination number are proved to be important factors controlling the oil migration process.
Despite some qualitative relationships and critical conditions obtained from these existing experiments, the influencing mechanisms of the pore-scale characteristics on the microscale oil migration behaviors are not clearly understood. Besides, due to the complex interactions of different factors and high costs of the forced fluid displacement experiments, the quantitative relationship between any single pore-scale parameter and oil migration cannot be achieved experimentally.
1.2. Pore Network Modeling of Physical Processes at the Pore Scale. PNM (pore network modeling) is referred to any analytical and/or numerical models of physiochemical processes prescribed on pore networks composed of pore elements that are connected through throat elements. Since pore networks can be constructed to resemble geometric and topologic characteristics of pore space of realistic porous sample, PNM has been proved to be able to capture the manifestation of pore-scale processes, through simulation, on sample scales and predict macroscopic transport properties in porous media [19,20]. Fatt [21] pioneered PNM and proposed a 2D regular network model for predicting capillary pressure and relative permeability of the primary drainage process in which the radii of pores were randomly assigned. Dullien et al. [22] extended that method to the 3D network model with more realistic pore-scale representation. The subsequent development of PNM has allowed the packing of grains with different sizes and shapes to be modeled, simulating a wide range of flow processes in different porous sedimentary rocks [23][24][25][26]. The readers are referred to the reviews [13,27] for a summary of the recent advances in PNM models applied to various processes. Recently, PNM has been applied to tight formations. Riepe et al. [28] presented a case history of combining CT imaging and PNM to establish a relationship between microscale structures and petrophysical properties of tight clastic rocks in Oman. That application shows the promise of an alternative approach for the evaluation of unconventional reservoirs. Ruspini et al. [29] introduced a multiscale pore network modeling workflow to compute the transport properties of clay-rich tight sandstone samples with wide pore size distributions. The representative network is composed of macropores, defined as completely resolved porosity, and micropores which are defined as unresolved porosity and treated as the equivalent continuum. Wang and Sheng [30,31] proposed a quasistatic PNM model to simulate the drainage process in shale and tight porous formations. The liquid non-Darcy flow mechanism was incorporated into the proposed PNM model using an empirical equation obtained from experiments [32]. Based on their simulations, Wang and Sheng concluded that the non-Darcy flow behaviors show significant influences on the absolute permeability predictions, but not on the relative permeability results.
PNM models have been developed for modeling the oil migration process in conventional formations, but not suitable for directly studying the oil migration process in tight formations. Compared with medium-or high-permeability media, the fluid flow in tight media shows obvious non-Darcy flow (NDF) behaviors [33,34], which is the flow velocity of fluid being not linearly dependent on the charging pressure gradient. Some researchers [35][36][37] proposed that these NDF behaviors were attributed to the existence of the boundary layer (BL) of water attached to the pore surface where the viscosity is much larger than the fluid in the pore center. The BL with different thicknesses in pores significantly reduces the effective seepage radius, which results in the NDF behaviors in tight media. In recent works, a parametric model of BL thickness of BL (see Equation (1)) was proposed by Cheng's group [14,38,39] based on the 2 Geofluids experiments where coefficients are fitted using the particle swarm optimization algorithm [40]. The calculation equation is given as where δ is the thickness of BL, μ w is the brine viscosity, and ∇p is the pressure difference between inlet and outlet per unit length. The three coefficients (h 1 , h 2 , and h 3 ) are determined to be 0.258, -0.261, and -0.419. The proposed method shows good predictions of BL thickness of numerous microtube experiment results conducted by Li et al. (2011) as shown in Figure 1. These experiments are conducted through simulating the pore-throat system with microfused silica capillary tubes, the radius of which ranges from 2.5 μm to 10 μm. The influence of BL thickness on the fluid NDF is analyzed by measuring the flow rate, pressure gradient, and viscosity experimentally.
A PNM for simulating the oil migration behaviors is developed in this work, considering NDF behavior due to BL of water. Equation (1) is employed to calculate BL thicknesses in the network elements to account for the NDF behaviors during tight oil migration. We are aware that the developed model employed the BL theory may be neither universally applicable nor accurate for all the cases of tight formations. The principal objectives of this work are as follows: (1) to extend the application of PNM to analyze the oil migration in tight formations; (2) to study the influences of pore-scale characteristics on oil migration using PNM method, which includes the radius of pores, average coordi-nation number, aspect ratio (the average diameter ratio of pores to throats), brine viscosity, and wettability conditions; and (3) to discuss the significances and limitations of the conducted PNM simulations.

A Quasistatic Two-Phase PNM with BL for
Modeling the Oil Migration in Tight Formation In this section, a quasistatic two-phase pore-network model with BL is developed for modeling oil migration in tight formation, by adapting a well-established quasistatic PNM model for the primary drainage process. This newly developed model is implemented to support a simulation workflow to be discussed in a later section.
2.1. The Model Development. The oil migration into initially water-filled water-wet pore networks occurs in forms of piston-like displacement [41]. When the driving force of oil phase exceeds the capillary entry pressure (P entry ) of the neighboring water-filled network elements, oil will displace water in pores/throats. We add a BL of immobile water of finite thickness, δ, to each network element whose radius and volume available to oil phase are reduced proportionally. So, as a result, we should expect higher capillary entry pressure at any oil pressure and higher oil relative permeability and lower water relative permeability at any water saturation. The capillary entry pressure, for a cylindrical throat, is calculated as follows: where σ is the interfacial tension between the oil and water, θ is the contact angle that quantifies the wettability property, and r is the pore radius. For an angular pore element, the capillary entry pressure can be adapted from that of Øren [42] as follows: The parameter G in Equations (3) and (4) is the shape factor for the network element, which can be calculated through Equation (6): where A is the cross-sectional area and S is the perimeter of the element.

Geofluids
To calculate permeability of each connected phase from the inlet to outlet faces of a PNM model, the single-phase flow between two neighboring pore elements is given by a Poiseuille-type law [43]: where P i and P j are the pressure values of pore i and pore j and g ij is the conductance between the two pores, which can be determined through the following Equation (8).
Here, g i , g j , and g t are the conductances of pores i and j and the throat connecting the pores.
In this work, the conductance of network elements is calculated following the functional form in the "3Rs" approach [44,45] but takes into consideration BL, and it is expressed as where A is the conductance constants, λ is the conductance exponent, L is the length of the element. μ in Equation (9) is the average viscosity of the fluids in the pores and calculated by where f o and f w are the fractions of oil and water in the calculation element. Based on the calculations of Equation (7), the pressure field of the network model can be calculated by applying the mass conservation at each pore element: When the pressure field for a phase in the network is obtained, we can calculate the relative permeability of each phase (oil or water phase) as follows: The new model is implemented on an earlier version of the open-source PNM code (https://github.com/ ahboujelben/numSCAL_basic) by extending its quasistatic two-phase model. For technical details in that code, the reader is referred to Regaieg et al. [46] and Boujelben et al. [47].
2.2. Simulation Procedure. Simulations are carried according to the procedure illustrated by Figure 2. For a given pore network that may be extracted from X-ray CT images of real rock samples or stochastically generated, (1) it is initialized 100% water saturation (S w ) at a reference pressure while wettability condition is applied to each network element; (2) explore the minimum and maximum radius (Min R and Max R) of elements, then calculate the corresponding maximum and minimum entry capillary pressure (Max Pc and Min Pc) with considering the BL effect. The calculations of Max Pc and Min Pc were realized by the iterative method; (3) set an incremental charging pressure (P n ) which increases from Min Pc to Max Pc. The incremental step is set by (Max Pc-Min Pc)/simulation steps; (4) explore for and determine the elements to invade. An invasion event occurs when the charging pressure overcomes the entry pressure of element. The entry pressure of each element is calculated by Equations (3)-(5). (5) Update phase saturation and flow rate with the consideration of BL effect. (6) Check the terminal condition whether the charging pressure is larger than Max Pc. If no, increase the charging pressure by step and do loop for steps (4) to (6). If yes, output the phase saturation and calculate relative permeability for each phase.

Model Validations
Model validation has been conducted by making comparisons between the predicted relative permeability and experimental results for a Berea sandstone sample and two tight sandstone samples (Samples S5 and S6). The relative permeability data measured on the Berea sample and its matched CT image data can be downloaded from the Imperial College database (the website address is http://www .imperial.ac.uk/earth-science/research/research-groups/perm/ research/pore-scale-modelling/). Samples S5 and S6 are taken from the Chang-7 members of the Yangchang Formation in the Ordos Basin, a typical tight oil basin in northwest China.

Geofluids
Oil charging experiments are carried out by using Samples S5 and S6 to simulate the physical process of oil migration in tight media. Experimental apparatus, experimental procedures, and the properties of used water/oil have been described in our previous work [15,48]. Both tight samples are scanned to obtain CT images using Zeiss Xradia-500 Versa with a resolution of 0.99 μm, at the Stata Key Laboratory of Petroleum Resources and Respecting, China. Raw CT images are preprocessed and reconstructed to 16-bit grey images using a multithreading software package "tomoRecon" [49]. The image filtering and phase segmentation are completed by using Avizo 9.7 software. An interactive thresholding module is used to segment the phases by selecting values of image intensity ranges for each phase [50]. The properties of the three experimental samples are shown in Table 1. It shows that the porosity, permeability, and pore connectivity of the tight sandstone samples (S5 and S6) are much smaller than those of the Berea sample. The pore size ranges of the tight samples are distributed in smaller and narrower scope compared to the Berea sample. Note that the lower limit of pore radius for tight samples is limited by the solution of CT images which means some pores with radius smaller than 1 μm may not be detected by CT scanning. More information concerning the pore-structure features in the tight formations of the Ordos Basin can be found in our previous work [17,18]. A modified maximal ball (MB) method developed by Dong and Blunt [51] was applied to extract pore networks from 3D CT images. The modified MB method explore and build an inscribed sphere at each void voxel of images first by a two-step searching algorithm, then remove those included in others. The rest spheres (MBs) are defined as pores and throats followed the clustering rule of family trees sorted by their size and rank. Geometrical properties are calculated for each pore and throat, including the size of pore and throat, volume, and shape factor. For more details of the pore network extraction, the readers are referred to Dong and Blunt's original publication [51]. In our work, pore networks of tight sandstone samples were extracted using an open-source software package (website for download: https://github.com/aliraeini/pnextract) that implements those mentioned techniques. The flow modeling of these three samples is conducted followed the steps in Sections 2.1 and 2.2. For the three samples, the predicted relative permeability curves from the newly developed model with consideration of BL effect are compared with those without BL effect, while both the model predictions are compared to experimental data (Figures 3 and 4). Figure 3 indicates that both the PNM models with and without BL consideration show good agreements with the experimental data, which in turn partly proves that Equation (1) for BL calculation is basically reasonable for the pore size distribution in natural samples. Figure 4 shows that the predictions of the BL model are in better agreement with the experimental data than that without BL effect in tight sandstone samples. Compared to the model without BL effect, the BL model shows lower permeability predictions for each phase of which locates in a narrower and higher S w range. The differences between NBL and BL models are mainly attributed to the decrease of pore radius and connectivity in BL model caused by the existence of BL. Specific impacts of the single influencing factor on relative permeability predictions are further discussed in Sections 4 and 5.

Sensitivity Analyses of Five Pore-Scale Factors
Stochastically synthetic pore networks were constructed for the sensitivity analysis, with respect to five featured porescale factors (i.e., the radius of pores, average coordination number, aspect ratio, brine viscosity, and wettability), on their impacts on the oil migration. The changing scopes of these factors are determined by considering the porestructure characteristics of tight formations and referring to the measured results in several published works [4,10,53]. The resultant relative permeability curves are employed for showing the relative fluid flow capability during the oil migration, while oil charging curves for explaining the oil content accumulated. Technically, the water in the BL is simplified as dead water without flowing but count for S w in the calculation of relative permeability and S o . Note that each pore network is generated as a regular lattice network first and throats were then removed randomly to match the  compared to the experimental data by Oak [52] (square). The oil relative permeability data are marked with red, while the water relative permeability data are marked with blue. Data modeling NBL/BL means that the calculations using the model without/with BL consideration, respectively.

Geofluids
prescribed average coordination number. The size of pore elements follows a uniform distribution in the prescribed size range. Pore networks and their parameters set for sensitivity analysis are given in Table 2.
Runs 1~3 are used to assess the impact of radius on relative permeability and oil charge where the coordination number remains the same but the upper limits of pore size ranges increase from 1 to 5 μm with a constant lower limit of 0.1 μm. The calculated absolute permeabilities of the models are 0.001 mD, 0.03 mD, and 0.9 mD. Figure 5 shows that with the decrease of upper limit pore radius range the maximums of the oil/water phase relative permeability decrease. Meanwhile, the endpoints of K rw at the S w axis increase from 3.5% to 20.8%. The predicted K ro with small pores at the low S w (less than 35%) is higher than that with large pores but turns to be lower at a high S w condition. The predicted K rw shows an opposite trend within the corresponding saturation range. The S w at the intersection points of K rw and K ro curves stayed constant for Runs 1~3. The oil charging curves in Figure 6 show that the calculated S o of models increases gradually with increasing charging pressure and terminates at the S o of 96.5%, 86.5%, and 79.2% for Runs 1~3, respectively. The required charging pressure for 0.1~1 μm model is much larger than that with 0.1~5 μm.   Figure 4: Predicted relative permeability of S5 and S6 tight sandstone (lines) compared to the experimental results (square). The oil relative permeability data are marked with red, while the water relative permeability data are marked with blue. Data modeling NBL/BL means that the calculations using the model without/with BL consideration, respectively. Geofluids the connectivity of the networks [54,55]. A higher coordination number means a better connection condition of the pores in the network. Runs 4~6 set with different coordination numbers (i.e., Z = 2:5, Z = 3, and Z = 3:5) are used to evaluate the impact of coordination number on oil migration behaviors. The calculated absolute permeabilities of the models are 0.01 mD, 0.03 mD, and 0.06 mD. Figure 7 shows that the maximums of oil/water phase relative permeability decrease with the decreasing Z values. Meanwhile, the endpoints of K rw curves at the S w axis are kept constant (S w = 13:5%) for different simulations. The calculated K ro in the poor connectivity model (Run 4) is slightly higher than that in the better connectivity models (Runs 5 and 6) at a low S w condition but turns to be lower at a high S w condition. The predicted K rw curve shows an opposite trend. The corresponding X-axis values at the intersection points of K rw and K ro curves stayed constant while the corresponding Y -axis values decrease with the decrease of Z values. The oil charging curves in Figure 8 show that the calculated S o of all three models increase with the charging pressure and terminates the same S o = 86:5%. However, the required charging pressure reaching the maximum S o decreases with the increasing Z of the networks.   Figure 9 shows that the influence of aspect ratio on relative permeability is relatively weak; the predicted K ro of Run 9 is slightly lower than others (Runs 7~8) within a limited water saturation range of 0~35%. The oil charging curves in Figure 10 show that the maximums of S o in Runs 7~9 are 77.8%, 86.5%, and 89.2%, respectively. The required charging pressure reaching the maximum S o for three models are equal.
4.3. The Impact of Brine Viscosity. Runs 10~12 with different viscosity (μ) settings are used to assess the impact of brine viscosity on oil migration behaviors. The calculated relative permeability results in Figure 11 show that the predicted K ro and K rw values both decrease with the increasing brine viscosity; under the condition of high water saturation, the impact of brine viscosity on relative permeability is more obvious which induces that the intersections point of K rw and K ro shift towards the bottom-right area of the coordinates. The oil charging curves in Figure 12 show that the S o of all three models increase with the charging pressure and terminates the maximum S o of 63.2%, 74.6%, and 86.5%. The required charging pressure reaching the maximum S o for three models are equal. Run7: aspect ration=1 Run8: aspect ration=2 Run9: aspect ration=4 K ro K rw rw Figure 9: Effect of aspect ratio on relative permeability. Run7: aspect ration=1 Run8: aspect ration=2 Run9: aspect ration=4 Figure 10: Effect of aspect ratio on oil charging curves.  8 Geofluids wettability condition before the hydrocarbon migration process is assumed to be water wet. In our simulations, contact angle (CA) is used for qualifying the wettability of rock surface [58,59]. Three models with different wettability settings (CA = 0°, Run 13; CA = 30°, Run 14; and CA = 60°, Run 15) are modeled to assess the influence of wettability on oil migration, where CA is assigned to be uniformly distributed in the pores and throats. Figure 13 shows that the calculated K rw becomes larger when CAs increase from 0°to 60°, while the wettability setting has little influence on the K ro predictions. Figure 14 shows that S o of the three networks grow and terminate at the same maximum S o value of 86.7% with the increase of charging pressure. The required charging pressure reaching the maximum S o decreases with the increasing CA under the limited water-wet condition.

The Applicability and Limitations of the Developed PNM
Model. Numerous models have been developed in the past decades to simulate the fluid flow processes in porous media. The forces set in these models mainly involve capillary pressure, viscous force, and gravity force. These network models simulating the two-phase flow system can be divided into two types: quasistatic model and dynamic model. The quasistatic model is suitable to simulate the fluid flow when the capillary pressure is dominating while the viscous force is negligible. When the nonnegligible viscous force is present; however, a dynamic model needs to be considered. When the capillary number is less than 10 4 , the simulations of the quasistatic model are believed to be accurate enough [19]. In conventional reservoir formations, oil migration can also be divided into two types: steady oil migration and episodic oil migration [60]. The oil flow under the steady oil migration mode is assumed to be governed by capillary pressure, and the flow in the episodic migration is governed by viscous pressure while the natural buoyancy plays a role in the fluid flow. However, for the low-permeability reservoir (e.g., tight reservoir and shale reservoir), whether both the steady oil migration and episodic oil migration can be identified to explain the flow phenomena in low-permeability media is still a question. A common consensus is that the effect of buoyancy on fluid flow can be neglected because of the extremely small size of pores in tight formations [61,62]. Due to the same reason, the effect of capillary pressure becomes significantly important, which means the fluid flow in tight media is more likely to follow the steady migration mode. This is the reason why the quasistatic method is employed in our present work. However, some other scholars [63] also proposed that present day reservoirs may not represent the conditions of the reservoir during the time of hydrocarbon expulsion and migration, indicating that the viscous force may have been the main controlling factor, and buoyancy may have been an important force forming the hydrocarbon distribution in tight oil/gas migrates and accumulates. Some other limitations of the PNM methods should be noticed here. (1) The fluid inside the networks is considered incompressible and immiscible. (2) The real pore-throat system of the rock is simplified as simple 3D geometric shapes following the pore-network generation rules. (3) The chemical interactions between the fluids and solid matrix are neglected. (4) The wettability conditions are set as constantly distributed, which may not exactly match the real condition in nature [64]. (5) The gravity forces and influence of microfracture in the pore-throat system are not considered in the developed model. For these limitations, the PNM model may not provide precise predictions for the actual geological conditions.
The permeability of the tested tight samples may be slightly larger than the permeability scope of tight sandstone by generally speaking. But the upper limit of tight media permeability is not standard. In some literature, tight oil is defined as oil resource reserved in very low permeability formation which cannot be produced at economic rates unless the formation is simulated by large hydraulic fracture Oil saturation (%) Figure 14: Effect of wettability conditions on oil charging curves. 9 Geofluids treatments or produced by the use of horizontal wellbores [65]. The main reason we select these two samples with larger permeability is that the forced displacement experiments are hard to implement when the permeability of the sample is smaller than 0.5 mD. For the same reason, a complete series of experiments using different samples cannot be completed. Furthermore, the BL theory is employed in this work to explain the non-Darcy flow phenomenon in tight media. This theory is merely one of the numerous non-Darcy theories. Some other theory such as threshold pressure gradient theory can be also used for explaining the non-Darcy flow [66]. It may be acceptable to use the BL theory under most conditions, although the range of application, prediction accuracy, and limiting conditions of this theory remains to be further investigated.

The Simulation Results and Their
Significances. Based on the reservoir exploration and production experiences, previous literature has concluded that the migration of tight oil behaves differently from the conventional reservoirs (as described in the Introduction). Some of these behaviors could be explained using the obtained simulation results in this work. For example, the initial water saturation (S wi , the measured water saturation before the well drilling) in tight reservoirs is reported to be much higher than that in the conventional reservoirs [63]. Seen from the simulation results in Sections 4.1 and 4.2, the S wi increases dramatically with the decrease of pore size and pore connectivity. The calculated S wi with the pore size of 0.1~1 μm (absolute permeability K = 0:01 mD) is 20.8% and the calculated S wi with the pore size of 0.1~5 μm (absolute permeability K = 0:9 mD) is only 3.5%. Note that the low S wi here does not mean most of water could be replaced by the injected oil during the oil migration process in the natural system, which is only an idealized calculation value under the ideal hypothesis conditions.
The influential mechanisms of several factors on relative permeability and oil saturation during the oil migration are analyzed as follows. The decrease of pore radius leads to the decreases of both oil and water relative permeability and the increase of BL thickness. This increase of BL thickness further results in a higher irreducible water saturation and a lower effective seepage radius. Thus, the two-phase seepage zone in the relative permeability curves decreases with the decreasing pore size; the maximum of oil saturation which corresponds to the irreducible water saturation in the oil saturation curves also shows a decreasing trend with the decrease of pore radius. Under the condition of the same pore radius distribution, the scope of two-phase seepage remains unchanged with the change of pore connectivity. However, the decrease of coordination number will reduce the seepage area of the media which leads to the decrease of relative permeability for oil and water. Besides, under the condition of poor connectivity, water phase is easy to stay in the corner of the networks. This phenomenon may cause relative permeability to move towards the high water saturation direction. The changes of connectivity of the pore systems may change the pore space volume and water/oil content by absolute values but show no influence on the irreducible water saturation. Therefore, the maximum of oil saturation would keep constant at a high injection pressure. The differences between the relative permeability curves of Berea sample and tight samples (Figures 3 and 4) are mainly attributed to these aforementioned two factors. The impact of brine viscosity on oil migration is relatively simple. The thickness of the BL would increase with the increase of brine viscosity which causes the decrease of effective seepage radius for oil and water phases as well as the increase of irreducible water saturation. Meanwhile, a higher brine viscosity would directly reduce the permeability of water phase, which causes the decreasing trend of water relative permeability in Figure 11 to be much more obvious than that of oil with the increasing brine viscosity. Overall, the simulated five factors all show influence on the relative permeability and S o , but with different sensitivities. Compared with aspect ratio and hydrophilia of the rock, the pore radius, pore connectivity, and brine viscosity play more important roles in controlling the oil migration behaviors.
The improved oil migration PNM provides a new solution for analyzing the effects of complex pore-scale characteristics on tight oil migration and accumulation. Combining the quantitative results obtained from models and traditional geological analysis methods may be more helpful for studying the mechanisms of tight oil reservoir formation and evolution.

Conclusions
To better understand and quantitatively analyze the relationship between the oil migration characteristics and pore structures in tight formations, a two-phase PNM with consideration of the BL effect is developed in this work. Several conclusions are drawn as follows: (1) This work proves that PNM method with consideration to the BL effect is applicable and effective for tight oil migration analyses, which shows better prediction accuracy than that without the BL effect (2) Five important factors including pore radius, coordination number, aspect ratio, brine viscosity, and wettability show influences on the tight oil migration, but sensitivities of these factors on oil migration are different (3) Based on the simulations, the influences of pore radius, coordination number, and brine viscosity on oil migration are more pronounced while the influences of aspect ratio and hydrophilia of the rock surface are relatively weak (4) With the decrease of pore radius, coordination number or the increase of brine viscosity in the oil-phase relative permeability is weakened, which induces less oil accumulation in the pore system

Data Availability
The experimental and simulation data used to support the findings of this study have been deposited in the Mendeley 10 Geofluids repository (the download address is https://data.mendeley .com/datasets/j3467gsr3k/1). The data is also available on request.

Conflicts of Interest
The author(s) declare(s) that they have no conflicts of interest.