Conjugate Solutions of Navier-Stokes Equation with deformed Pore Structure

Within the framework of block self-organizing of geological bodies with use of deformation theory the mathematical solution of a problem for effective final speed is proposed. The analytical and numerical integrated solutions of Navier-Stokes equation for deformable porous space were obtained. The decisions of multi-scaled regional problems «on a flow basis» were also presented: from lithology of rock space to a well and from a well to petro-physics. The evolutionary transformation of the linear solution of the equation on mass conservation up to the energetically stable non-linear solution of the equation on preserving the number of movements is also offered. Basing upon the analytical solution of Navier-Stokes equation and model of A.N. Kolmogorov we have obtained the energy model of turbulence pulsing controlled chaos, conjugated with risk stability of average well inflow and cluster structure of Earth de fluidization. Introduction Revolution with shale has brightly illustrated the possibilities of theoretical solutions for NavierStokes equation in the burning issue of the century – global geology correlated with global mathematics, geo-physics, geo-dynamics and global economy. Multi-scaled solutions for continuous medium of porous micro-structure, the equations based upon the conservation laws and phase trajectories of micro-world pulse moment are correlated by energy with conjugated solutions of Navier-Stokes equation and conditions of system equilibrium (porous matrix). Physical and mathematical solutions of inflow profiles with due consideration of geological risks in development of the reservoir with complex structure are well coordinated with geo-physical, seismic experiments, actual well production rates and cluster structure of the economics. Theoretical Substantiation Equations of Euler and Navier-Stokes describe the fluid movement in R n (n=2 or 3). These equations are solved using the unknown vector of velocity u(х ,t )=(u i(х , t ))1 ≤ i≤ п ∈ R n and pressure р(х,t)∈ R, determined by the position of х∈ R and time t ≥ 0. The equation of NavierStokes for conservation of number of movements has the following view ∂ ∂t ui + Σj=1 n uj ∂ui ∂xj = ν∆ui − ∂p ∂xi + fi (х, t), ( х ∈ R , t ≥ 0), (1) Bulletin of Mathematical Sciences and Applications Online: 2014-05-01 ISSN: 2278-9634, Vol. 8, pp 30-48 doi:10.18052/www.scipress.com/BMSA.8.3


Introduction
Revolution with shale has brightly illustrated the possibilities of theoretical solutions for Navier-Stokes equation in the burning issue of the century -global geology correlated with global mathematics, geo-physics, geo-dynamics and global economy. Multi-scaled solutions for continuous medium of porous micro-structure, the equations based upon the conservation laws and phase trajectories of micro-world pulse moment are correlated by energy with conjugated solutions of Navier-Stokes equation and conditions of system equilibrium (porous matrix). Physical and mathematical solutions of inflow profiles with due consideration of geological risks in development of the reservoir with complex structure are well coordinated with geo-physical, seismic experiments, actual well production rates and cluster structure of the economics.

Theoretical Substantiation
Equations of Euler and Navier-Stokes describe the fluid movement in R n (n=2 or 3). These equations are solved using the unknown vector of velocity u(х,t)=(u i (х,t)) 1 ≤ i ≤ п ∈ R n and pressure р(х,t)∈ R, determined by the position of х∈ R n and time t ≥ 0. The equation of Navier-Stokes for conservation of number of movements has the following view With initial conditions of u (x, 0) =u 0 (x), (х∈ R n ).
The final volume of liquid is placed into elastic shell that is in the state of equilibrium. Without any additional preconditions of elastic equilibrium the solution of Navier-Stokes equation becomes impossible. The equilibrium in saturated porous media having the changing density and elasticity may be described by a system of equations for classic continuous viscous-elastic media and generalized form of Hooke's law: σ q ij =µ q ·(ξ i , j +ξ j , i )+λ q δ ij· ξ i , i , σ q i j , j =ρ q ∂ 2 ξ i /∂t 2 (4) Where λ, µ, ρ -generalized parameters of Lame and density of porous media, q=1-Nnumber of a layer, i, j=1-3. At the boundaries of the intermediate layers of the porous media, the matrix and fluid: σ ij q -1 n j | x = h q =σ ij q n j | x = h q , ξ i q -1 | x = h q =ξ i q | x = h q (5) The external layer is immobile where nnormal to the boundary. Equations (1 to 6) define the equilibrium conditions for viscous and elastic deformations at the basis of fluid movement equations and porous frame with due consideration of inertial components. Normal and shear deformations ε are connected with relocation of ξ by geometrical correlations: ε ij =(∂ξ i /∂x j +∂ξ j /∂x i )/2. Excluding the relocations from (7) one may get the defoemation compatibility equations. The solution of elasticity equation system should satisfy the boundary conditions at the surface of the body (Fig. 1). Geometrical boundary conditions that define the character in fixing of a body at the external surface S u are overlapped by the re-location. Static boundary conditions of surface S σ define the character of body loading by external forces and are recorded as: t k =T k . Здесь T kprojections of loads upon axis, t i =σ i l i +τ ij l j +τ ik l k , where σ i ,τ ij -stress in boundary surface; l idirection cosines of the normal. Of all that's available, the actual status of the body within the limits of given loads and fixing conditions we may get this through defining the energy criterion of the equilibrium.

Bulletin of Mathematical Sciences and Applications Vol. 8
Let's make the re-location so that they become continuous functions of coordinates and correspond to the fixing conditions of the considered body. There exist many systems of re-location making the battery of kinematically possible condition of the body. For each of these conditions the geometrical correlations enable to find the deformations and physical correlations define the stresses. In general these stresses will not satisfy the equations of equilibrium and static boundary conditions. The possible system of re-locations that satisfies equation (7) and condition t k =T k and is a genuine one. Within the frames of Cauchy problem at the external boundary of the viscous layer S σ we accept t k =P o , where P o -initial formation pressure, and at the external boundary S ucondition (6).
The method to calculate zonal non-uniform well inflow has been offered by M. Muskat (1937). The solution is considered in generalized form as: p=c o ln r+∑ r α (a α sin αθ+ b α cos αθ) +∑ r -α (c α sin αθ+ d α cos αθ).
The rate in view of static equilibrium conditions as per Masket is Q=2πk(P o ·s-P w )/[µ·ln(R o /r w )], where k is reservoir permeability, s -angle, P o , P w . R o , r w -pressures and radii accordingly at the well contour, µ -dynamic viscosity. However, this decision is not used in the analytical analysis of reservoir development, nor in geo-monitoring, nor in numerical modeling of porous media. That, actually, complicates interpretation of well geophysical studies, including seismic ones and deforms volumetric 3D multiphase structure of fluid movement and on self-organizing of geological structures in thermal, electro-magnetic and gravitational fields during the arrangement of deposits and geo-hydrodynamic processes of modern geology, exploration and field development.
First of all let's enumerate the hydro-dynamic problems and geological risks in resolving the equations of thermal conductivity type, used in commercial 3D simulators that are eliminated by generalized solutions of Navier-Stokes equations and viscous-elastic deformation of porous space.
 The uncertainty with shift migration and ambiguity with top and bottom sides, intra-well space during various cycles of oil and gas field development, seismic study, search and exploration.
 The final velocity of effect transfer, the value of drainage and reservoir oil productivity; dependence on porosity, permeability and reservoir heterogeneity.
 The integrity and combination of hydro-dynamic solutions with geological and geo-physical information, growth in geo-information of numerical models.

Increase in Geo-information Knowledge on Geological and Hydrodynamic Modeling of 3D
Filtration. In order to study the geo-dynamics of porous media, waves of various origin in layered heterogeneous sheared structures, seismic emission at the boundary of functional density, let's consider the energetically linked multi-layered model, described by the equations of conserving the number of movements and equilibrium of porous space. The stirred pulse movement of viscous non-compressible fluid is described by Navier-Stokes equation (1) for the mean average value of U and pulse velocity u: ∂u i /∂t+u j ·u i , j -‹u j u i , j ›+u j ·U i , j +U j ·u i , j =-1/ρ f ·∂p/∂x i +ν∆u i (9) and continuity of u i,i =0, νU'=u * 2 +‹u i ·u j ›, where u * =(τ/ρ f ) 1/2 -dynamic velocity. Disregarding quadratic terms of pulsing velocities the equation of movement for viscous noncompressible fluid (9) allows to have the solution as [1]: where k=(k 2 2 +k 3 2 ) 1/2wave number. In order to resolve this task we take Hilbertian phase space. The average values of type ‹u j ·u i,j ›=(u i,j ·u j * +u i,j * u j )/4, where * -complex conjugation, and it excludes the dependence of time. Substituting (10) into system (9) and neglecting quadratic terms of velocity pulsing, we get a system of ordinary differential equations (DE) for u and p:

Bulletin of Mathematical Sciences and Applications Vol. 8
The solution of the task dealing with the interaction of viscoelastic porous layer and viscous layer is determined with an accuracy of normalization factor с 1 , i.e. intensity of pulsing pressure at the external boundary of the marginal field. Substituting (25,26) into (23,24), we anticipate to get the following conditions at the viscoelastic surface: So, the solution of the joint task of inter-acting the viscoelastic surface with sheared layer of the constant stress is presented in a form of resolving a single task for oscillations of viscoelastic layer with modified boundary conditions of the type (27).
At the surface of the stiff layer the conditions of adhesions are presented in a form of: For the task σ 22 , σ 12 , ξ 1 , ξ 2 it's required, respectively, to define the unknown factors c i , i= (4-4N+3).
Resolving the joint task of viscous sub-layer oscillations and the oscillation of a layer with constant stress and their effect upon the deformed surface as a result of pulsing load actions and using the ideology of maximum energetic stability [2] principle of the conjugated solution we have performed the dispersion analysis and have disclosed the character of pulsing effective parameters that affect the average value of the velocity U (η). The resulted energy model for the controlled filtration chaos of pulses u is the conjugated risk stability for the average inflow in a cluster structure of Earth defluidisation organizing process [3].
We have calculated the profiles of pulsing energy E (η) = (u 2 +v 2 +w 2 )/2, diffusion -‹pv› as by the thickness of the sheared layer, the Reynolds stress τ for the average flow at varying the frequency ω and k. At the basis of Reynolds stress analysis that guides the stability of the average velocity we have done the following conclusion: at any value of phase velocity the Reynolds stress have the stab le for up to η~35. This means that the conjugation point with asymptomatic continuation of the profile is in the area of η~35.
The same calculations are performed for the sheared layer at the viscoelastic surface presented by the function of µ (t) =µ o +Σ n j=1 μ j e -t/τj jointly with the time relaxation spectrum τ and static viscosity μ o . Dynamic shear module for harmonic law of loadings: Dynamic volumetric module is equal to static, that's why The viscoelastic shear deformations reposition the conjugation point due to induced negative Reynolds stresses at the surface. In this case there is the reduction in parameter oe and the growth in the average velocity of the flow. (4), presented in a form of Helmholtz, in a cylindrical system of coordinates [4,5] has the type of Bessels' differential equation. Ther solution of the resulted system is written down in a form of longitudinal and transversal waves with complex wave number of k=α+iβ:

Natural Properties of Mixed Deformed Viscous-Elastic Layer. Equations
While satisfying the boundary conditionsm we get the transcendal system of equations of the 4 th order. This gives the characteristic equation to find the nominal frequencies of sheared layer In area of natural frequency ω к~1 the amplitude becomes several dozen times larger than thickness of a shift layer, forming transversal fluid inflows and attached lithology. Anisotropy of permeability becomes less than 1, increasing vertical displacement, water cross-flows, meniscus of inclined gas-oil and water-oil contacts, wave alluviums of carbon/micas clay particles at the borders of oil-containing rock. The profiles of transversal and attenuated normal stresses (Fig. 3) in the fixed layer with parameters: μ 0 2 /μ 0 1 =5, ρ 1 =ρ 2 , α̃ =0, 04 (1); 0, 5 (2); 1 (3); 2, 5 (4); 5 (5); 10 (6), where α=kh, show the structure of the layer and are used at inversion of 3D seismic. As is seen in Figure the greatest normal stress occurs at the surface, and inside the layer has the minimal value, i.e. the radial component of displacement naturally fades. For tangential stresses this dependence is of the other kind: the maximum values are at some depth of the cover or at the surface, i.e. there is the regeneration of a radial component of moving vector into a tangential one to turn fluid component round the axis, the normal stress is split and the velocity changes the phase. Fig. 4 presents the amplitude/frequency spectra obtained at the seismic stations located within the Voronezh crystalline massif at the ancient East-European platform. The studies [6] show that

Bulletin of Mathematical Sciences and Applications Vol. 8
within the frequency range 0.1 -0.3 Hz the micro-seismic noise is generated mainly by microseismic of decompresses migration channels incoming from the Atlantic Ocean. The seismic stations are located within the Voronezh crystalline massif for which the Atlantic Ocean is the nearest one. The analysis of average daily variations of micro-seismic noise within the frequency range of 1.0 -8.0 Hz have shown that the correlation of timely variations with micro-seismic noise may be allocated to the fact that geological structure and anthropogenic load produce the prevailing effect upon the arrangement of relatively high-frequency constituent.
With low frequencies we notice the distribution of a bended wave of channel migration with natural frequency ω o . The amplitude of the transverse wave at high frequencies asymptomatically attenuates. At frequencies more than 1 Hz the porous media generates numerous transversal waves with high geo-informatics.
In general it's possible to make a conclusion that micro-seismic noise and its variations have global character.

Interaction Stability in Mean Characteristics of Mass Transfer and Pulsing Deformations.
To understand the thin equilibrium between a viscous layer and non-viscous internal area as well as to confirm a principle of maximal stability for mean characteristics of pulsing boundary layer, we should consider a task to simulate the pulsing current on a viscoelastic porous surface. Thus by an iterative way we take into account the non-linearity of a profile U (η) at the basis of the decision in quadrature.
The system of non-linear differential equations (9) to simulate the pulsing of a boundary layer at filtration is resolved. The profile of mean velocity U(η) is set on the basis of a principle of maximal stability for mean characteristics received after resolution of a system of linear differential equations, proceeding from approximation in viscous sub-layer U (η)= η (18-21).
The given system of non-linear differential equations is resolved numerically, by iterative approximation in final differences. To determine the function of viscosity θ(η) and pressure p(η) while resolving the equations (8-12) the value of pulsing v(η) is taken from the previous iterative step.
Conditions of the interface which have been written down in the form that doesn't contain the unknown parameters, look like: where τ= -ρ f ‹uv›, E= (u 2 +v 2 +w 2 )/2. First of them is a condition final normalizing. The requirement to have maximum curvature in a point of interface |E II (R)|, which provides the maximal stability of the considered profile and is a final precondition in defining the pulsing amplitudes and unknown parameters k, С, R.

Volume 8
At each iterative stage the task is to pick up the unknown parameters of k and с for the pulsing movement, so that to make smooth interface in some point R with the profile of mean velocity down to U IV and by mean pulsing energy -down to E II with asymptotic dependences U=1/oe ln (η) +C oe , E=E 1 (1+B/η).was performed.
The parameters of movement in a buffer zone are under established by the formulas: Thus, the profile of mean velocity is determined in complete zone of viscoelastic properties' influence for the sheared layer of porous border. Distribution of pressure (by Reynolds) τ/ρ f =-‹uv› shows, that the non-linearity stabilizes the profile in relation to small perturbance both at a rigid and viscoelastic basis. The task then becomes correct. The calculation has shown such an important feature as the leap by 180° for the relative shear between pulsing pressure and tangential component of a pulsing velocity upon the transition through a layer of concurrence, where the local vector of velocity is equal to phase velocity of waves, i.e. a point of phase leap U (η к ) =C, that is determined by the initial profile U (η).
Being the first approximation to the subsequent non-linear iterative evolutionary specifications, the linear approximation of essential non-linear equations of Navier-Stokes with good quality and correctly reflects distribution of pulsing fields of velocities and their correlation moments. The point of interface being the border of viscous sub-layer and zone of generation for mean filtration channel, correctly defines even the very first linear iterative stage of the calculation. While evaluating the non-linearity U (η) the specific point η=R is not mixed and is equal to R = 32.
This confirms a physical solvency of linear approximation of a zone with constant stress η=0-50. The linear task arrangement for the interaction of viscoelastic boundary and pulsing boundary layer shows the integration way of seismic emission layer in the tasks to reduce the friction resistance, to save the process time in parallel calculations for the multi-sectional geological models.
Energy Balance of Pulsing Velocities at the Surface Deformed. Of late a great deal of attention was given to statistic studies of pulsing flow parameters in pliant surfaces. Specifically there are the data on the intensity of spatial and time correlative functions and spectra of pulsing velocities. Each of the named values characterizes the eventual field of velocities as one-sided. It's obvious to correlate and analyze the random information coming as a result of the test jobs. The basis for such summation may be presented as theoretical model that explains the effect of flow rate for the given radius of the pipe with the turbulent flow at a given pressure drop. This model is based upon the absorption of the flow pulsing energy by the surface of viscoelastic coating and the redistribution of characteristics in a pulsing boundary layer by thickness. The diffusional energy flow at the viscoelastic boundary is determined by the kinematic properties of the surface and the dissipative properties of the material deformed surface.
Basing upon the numerical solution of auto-modeling for the balance of pulsing energy and number of movements, recorded in the cylindrical system of coordinates for axial-symmetrical case, describing the multi-phase pulsing flow in a coated pipe and pulsing filtration in cavernous/channeltype porous space with low-permeable matrix: Equations (35,36) describe the interaction of generative process, transfer and dissipation of nonuniformal velocities and reflect such a simple fact that local balance of pulsing energy is made opf Bulletin of Mathematical Sciences and Applications Vol. 8 pulsing energy inflow from the average filtration movement, i.e. due to Reynolds stresses (the first element), diffusive flow of of pulsing energy (the second element) and dissipation of pulse energy into heat. The theoretical analysis of the set of equations (35,36), strictly speaking, is not possible as the system with any finite number of DE for the moments of hydro-dynamic fields in a pulsing flow is always not closed. To close the equations sooner or later we are to introduce some suppositions or to apply any test data.
Let's intropduce the factors of pulse exchange: kinematic factors of pulsing viscosity ε b and diffusion ε g through the ratios: We note that for the flow with transverser shift this ration (36) presents very simple transformations.
Pulsing flow is combination of numerous vortexes, inertia moments of by-passing that strike the flow and create the irregular quick-changing movement. The supposition of A.N. Kolmogorov includes the fact that the factor of exchange ε b , ε g and energy dissipation velocity ε t in this point of the flow is determined by the local average energy of the vortex and for this we may accept the value of an average pulse energy mass unit E, and local average size of vortex or proportional to it external scale of pulses L. The analysis of sizes result in standard correlations as ε g =A g LE 1/2 ; ε b =A b LE 1/2 . The velocity of dissipation as per the hypothesis of Rotta has the form: The values of A b =0,5; A g =0,01; c 2 =10 -5 ; c 1 =0,05/A b in kinematic factors or pulsing viscosity, diffusion and dissipation velocity (37) are chosen at the basis of collating the calibration calculations for the rigid smooth surface with experimental data of Conte-Bello and Clarke in channels and Laufer in a pipe.
In order to account for the additional effect with the change in kinematic factor of pulse diffusion against the turbulent flow in a rigid pipe that result from the previous calculations, we may record it as follows: where multiplier e -βη defines the attenuation of the viscoelastic layer effect as it extends from its surface. The presence of a dissipative bed under the sub-layer of pulsing flow violates its persistence when compared with the flow at a rigid surface and enables to change the balance of pulsing energy near the surface as well as symmetry of well inflow flitration profile. In order to approximize the shifting fuction L(η) we have used the cubic parabola with Van-Drist absorbing multiplier where η=1-R/R o , ae=0,4; ae o =0,14. Transferring to system of coordinates η=1-R/R o , after nondementionalization and basing upon the equation of movement numbers we get where Re * =u * R o /ν, u * 2 =τ/ρ f =1/ρ f │∂p/∂z│R o /2. The equation of balance for pulsing energy is recorded as follows Inside the internal area of the sheared layer all its members that are the parts of pulsing energy balance have significant meaning. That's why we should not disregard them as they characterize the diffusion of the pulsing energy. Following the data of Laufer for the internal area beside the member that originates the pulsing energy -̅̅̅̅̅̅ the positive role (in the pulsing energy) is brough by the transfer of the energy by viscosity ν∂ 2 E/∂y 2 , as well as by the pressure diffusion ( ⁄ ̅̅̅̅̅̅ ). The role of pressure diffusion in the zone of viscous sub-layer is at maximum. That's why while studying the pulse sheared layer at a deformed surface that influences primarily upon the internal area of the sheared layer, it is required to follow the changes with a memnber of pulsing energy balance equation that characterizes the diffusion of the pulse energy. Diffusive flow at the boundary of the viscoelastic layer increases the flow of the pulsing energy towards the external area of the sheared boundary layer. As a result of this the process is stabilized at lower levels of pulsing energy intensity (Fig. 5, a). Fig. 5, b contains the distribution of diffusive layer with several types of viscoelastic layers that are characterized by the diffusion factor ε g (0)=0 (1), ε g (0)=10 -3 (2), ε g (0)=10 -2 (3) and, respectively, by Reynolds stresses (1' -3') with the same parameters of the covering.
The growth in pulsing energy flow from the area of its generation towards the walls and in the direction of external flow with growing diffision factor is not the same in its meaning. The flow in the direction of the wall provokes the growth in pulsing energy dissipation, and in the direction from the surface it increases the pulsing energy intensity to the external area thus increasing the Reynolds shear stress. This growth in stresses stimulates the increase in growing the velocity of the boundary layer, the thickness of pulse loss thus resulting in a growth of friction resistance factor.
Consequently the analysis of equation (1) with the use of Kolmogorov's suppositions enables us to draw the conclusion on the mechanism of maintaining the dynamic equilibrium inside the The factor of resistance is Fig. 6. The average velocity in smooth pipes (1)(2)(3)(4)(5) and with damping area depending upon the number of Reynolds (а); resistivity factors with viscoelastic coating (1) and the ones calculated as per formulae of Blizius (2) and Nikuradze (3), as well as with roughness d/k=507-60 (4-7) (b) Fig. 6,а contains the calcuklation results for the profile of velocity in a rigi pipe and in a pipe with absorbing coating. Viscoelastic coating is characterised by kinematic factor of pulsing diffusion at the surface ε g (0)=10 -3 . Kinematic factor is chosen in such a way so as to arrange a diffusion flow at the boundary of viscoelastic layer that is characteritic for the polyurethane coating at given values of Reynolds number.
The results of numerical expertiment are verified by the experimental studies of kinetic flow conditions above the absorbing surfaces, where we get the change in spectrum density of pulsing energy and stresses by Reynolds in the internal part of the flow. Fig. 6,б shows the graph of resistance factors in a rigid smooth pipe that may be calculated as per the formula of Blazius C f =0,3164/Re 0,25 at Re <10 5 and as per formula of Nikuradze C f =0,0032+0,221/Re 0,237 at Re >10 5 . The effect of absorbing coating in a smooth pipe is proportional to Re number. The resistance factor in a pipe with absorbing coating at Re=1,5·10 6 , ε g (0)=2,5·10 -2 (curve 1) C f abs. =0,46C f. Bearing in mind the fact that in an instant sheared skin-layer of the well due to the change in velocity phase we lose around 50% of the pressure drop and the influence of absorbing reservoir surface (proppant) upon the friction coefficient will be at maximum. Due to the additive effect of these factors the production rate may be increased by 40-50% at a constant pressure drop.
Block-type Structure of Reservoir Arrangement, Smart Development of Hydrocarbons. Now the oil companies more and more face the problems with new well production rate reduction. In order to maintain and increase the levels of oil production it is necessary to locate the wells in challenging geological conditions (compacted low-permeable rock; reservoirs with high clay contents; decompressed unconsolidated zones). All these result in high uncertainty of calculations, both for initial production rates and for long-term well operation forecast. The same time we see the growth of hard-to-recover hydrocarbon reserves, like bitumen, etc. Recently a great deal of attention was given to shale gas and shale oil, their position in classical pools.
Thus the new technologies applied for the reservoirs, the development of means, technologies of control and production management as well as modeling the process of hard-to-recover field development had not found the adequate presentation in mass-scaled decisions inside the oil companies. It is related with the fact that very often there are no specialized physical/mathematical algorithms and techniques to be used for challenging cases of oil-saturated reservoir structures that finally allow improving the geophysical control and monitoring over the processes of their development.
The scientific support for the given problem is directly connected with the solution of important theoretical and applied tasks [7], such as:  Final transfer of influence with the due account of rock petro-physical parameters: oilsaturation, compressibility, porosity, etc. New fundamental decisions on productive wells and base production from oil-saturated and condensate reservoirs.
 Contact between fluid flow and porous media in space being one of basic problems of mathematical theory to make mass-scaled averaging of volumetric mass transfer and wave transfer of movement amount;  Interpretation of well test data with the use of conceptually new information on reservoir structure (dynamic skin factor, anisotropy, lamination, updating of knowledge on top and bottom location, etc.).
The result of these works was the development of new knowledge on challenging geological media, on system organization [5] of multi-scale phase movements, development of concepts on auto wave phenomena, principles synergy and self-organizing of geospheres of the Earth. А.N. Dmitrievskiy has generated the concept «of vortex geo-dynamics» of the Earth and block self-organizing of hydrocarbons (HC). The system multi-scaled approach in increasing geoinformation knowledge during search, exploration and development of oil/gas fields, study of Palaeo-and Neo-tectonic evolution of lithosphere, effect from man-made processes upon biosphere and antroposphere is further development of V.I. Vernadskiy's ideas.
In spherical shells of the Earth there occur the wave processes reflecting energy carry-on and transformation. These are electromagnetic, gravitational and thermal fields, elastic seismic waves, influence of phase transitions and more. Between the shells there are reflecting borders that create energy barriers between them. During the arrangement of elements of auto-wave system the cells capable to keep their own energy are formed. Concentrating in these cells the energy in a resonant mode interacts with adjacent cells and in case with "energy over-emission" increases its amplitude. The change in the level of dynamic heterogeneity of the media results in spasmodic qualitative changes with zonal heterogeneity expressing in processes that change their macroscopic structure. The massive transitions in heterogeneous media are based upon the combination of the mathematical theory of averaging for "operators" with high-oscillating properties, asymptomatic methods and theory of filtration [7].
The planetary stresses of lithosphere form dynamic structures of stretching and compression as well as diagonal structures of shears. Thus they form the pools with various structural confinements 1 .
First, structures connected to erosive surfaces of the foundation and at their step-like immersing. Second, with lateral distribution of stresses at stretching, there appear finer structures for which the influence of sedimentological factor (con-sedimentalogical processes) is characteristic. Third, the stress of compression is distinctly seen arranging neotectonical structures.
One of the geological examples may be presented by block-type movement of Earth lithosphere. The situation with block-type movement is most obviously seen in formation of structural elements and hydrocarbon traps of the sedimentary cover and this is used while simulating the regions with ancient foundation (Fig. 7). Planetary stress in lithosphere forms dynamic structures with extension and compression as well as the diagonal structures of the shears. With this in view we have the formation of the pools with various structural allocations.
Inside the energy field of the Earth there is the process of polarization of porous space, of forming the layers with block-type self-organizing, that is obviously fixed by the change in thermal, concentration and other types of geophysical fields. Large gradient of velocity pulsing, their amplitudes and phase angles are severely changed at relatively short lengths, thus forming abnormally high porous pressure. The very important aspect is the one that within the frames of resonant pulse mechanism of interactions the role of convective pulse transfer in the force balance becomes significant.

Setting-up a Problem and Methods to Resolve It.
Multi-scale "abnormal" behavior of "matrixfluid" pair with allocation of the main filtration flows is the generalized concept of Pollard-Pierson's double porous space ∆p(t)=Σ 3 j=1 a j exp(-α j t). Here α 1 >>α 2 >>α 3 , correspondingly, for fractures, porous matrix and viscous sheared layer. It's presented as the further integrated conjugation of stress interaction in porous space and have the generalized model of visco-elasticity ∆p(t)=Σ n j=1 a j μ j e -t/τj , with flow, displayed by normal pressure gradients. Along with that the concept of "self-organizing geological bodies" gives us the understanding of various-in-space structures and time, repeatability of processes in geological time, modeling in predictable "technogenic time" effect upon the geological body. With a stationary fluid movement along the stream-line the equations of movement and its continuity are presented by the equations in partial derivative of the velocity and pressure in view of the micro-structure scale h and relaxation time τ [8]: ρ f ·(u∂u/∂l+τu∂/∂l(u∂u/∂l))=-(∂p/∂l+h 2 /6∆∂p/∂l), (54) ∂/∂x k (u lk )+h 2 /6∆∂/∂x k (u lk )=0 (55) Differential equations (54,55) describe the relaxation processes in velocity and pressure delaying. The structure of porous space and composite content of continuous body will be presented as the microstructure. In classic constant-active geological/hydro-dynamic models of oil-saturated reservoir we don't account for the time required to bypass the compacted micro-structure. The modification of DE in view of instant micro-structure is presented in papers by Е. Koser (1903), L.I. Sedov, V.N. Nikolaevskiy [9].
The dynamic non-equilibrium sheared layer with tangential stresses possesses the maximum level of energy accumulation and dissipation. The stress-deformed status of deep geosphere arranges the asymmetric tensors of deformation and matrix permeability. Paleo-and neo-tectonic fracturing possessing the anisotropic direction (depending upon dynamics of lithosphere in a specified geological age) strengthens the lateral and vertical hydro-dynamic connection between regional oil and gas bearing complexes. The presence of porous media radically changes the hydro-thermo-dynamics of multi-phase mixtures [10]. The continuous medium of isothermal well inflow having radius r w is presented in a form of viscous block with thickness H, porosity ф, average velocity (as per Darcy) U=Q/(2πRHф) and viscoelastic media of saturated porous space with structure h (Fig. 8).

3D Modeling of Block-type Organization with Capillary-gravitational Multi-phase Fluid and
Porous Space Equilibrium. The obtained energetically steady solution of mean profile of inflow allows to create multiphase simulator of a new generation to calculate the net-type models of developed geological objects with the volume of hundred millions cells by a method of splitting the calculations by physical processes using sides of cells to transmit wave pulses. Thus the energy regulating parameter is presented by the energy of pulse or velocity.
The equations of modern classical three-phase (black-oil) hydrodynamic simulators have the form of: Here, p α -pressure in phase α, S α -saturation, k -permeability, b αvolumetric factor, k rα (S α )relative phase permeability (RPP); q α , m 3 /(m 3 ·s)velocity for drain unit volume, D -depth, α=o,w. The ratio of hydrodynamic and capillary forces is expressed by dimensionless parameter of capillary number N c =μu/σ, where σ -surface tension. With non-stationary filtration we take the dependence of RPP from N c :

Volume 8 k ro (N c ,S o )=a o N c [(S o -S ro (N c ))/(1-S ro (N c ))] εo(Nc) .
(57) The equations of mass conservation for phase α (e.g. for oil) in a shift layer accepts a kind of pulse conservation equation Here, Diffusive conductivity of displacement front λ ∑ = о λ 2 о |u о | is proportional to velocity, reverse proportional to surface tension, volumetric factor and has obvious zonal heterogeneity. Conductivity of the boundary for phase block α attenuates as per quadratic law with reduction in permeability k, RPP, and compressibility ∂∕∂x(1/B α ). The growing conductivity of the freely penetrated boundary leads to self-organizing channel with capillary jam in a stagnant area of the sheared non-permeable layer and to high well water-cut. Multi-scaled Non-standard Well Studies. As an example we have considered the wells that are operating the Famennian stage of Devonian age at Yuzhno-Orlovskoye field. Well № 001 has crossed the zone of rock resonant-wave decompression and in conditions of practically impermeable matrix (rock porosity of 3%) is operated with a rate of 1000 m 3 /day, producing waterless oil (see Fig. 9). Dynamics of RPP has shown, that the real oil stocks there make something like 15 MT, whereas the ones calculated using a volumetric method (as based upon the porous space of matrix) were around 2.5 MT.
Hydrodynamic Solutions in Scaled Filtration. Papers of the authors and the Report to SPE Conference (SPE-166893), EAGE [11] present the theoretical studies and physical/mathematical modeling to develop complex-structured fractured-cavernous reservoir using the 3D three-phase software for porous media like CMG STARS, ECLIPSE, FLORA [12], that take into account the dynamic changes of RPP. With this in mind we have proposed the procedures with horizontal wells to be used in developing these complex-structured reservoirs.
Within the instantaneous multi-scale structure of inflow in self-organized geo-systems there are "convective-diffusional counter-flow drains and "skin-effects". Low energy segments are not included into development thus arranging the oppressed filtration and low final inflow into a horizontal part of a well (Fig. 10). The growth in horizontal well-bore length with similar reservoirs does not result in growth of production rate, since only high-permeable zone of the reservoir is operated. The efficiency to develop the similar complex (cluster) reservoirs directly depends on diagnostics of inflow profile by geophysical and hydrodynamic methods, allocation of highpermeable channels, isolation and restriction of the cross-flow negative influence, effect of the increased pressure upon the segments low-rate diffusional drain (Fig. 11).
Experimental measuring (Fig. 11) verifies the substantiation of multi-scaled filtration model. Theoretical and lab dependencies of growing residual oil-saturation at dropping the filtration velocity are well-coordinated and supported by numerous studies [14,15].
While simulating the inflow profile for the block-type segmental porous space with complexorganized boundaries it is required to make conjugation of boundary conditions in pressure and in rates at the basis of relative phase permeability (RFP) that depends upon the velocity. The conjugation of dynamic phase velocities provokes the arrangement of energetically heterogeneous segments: choke-type, channel-type, radial-type and convective-diffusional types of flow, of elasticgravitational, concentrated counter-flow impregnation, capillary jamming and filtration depression.
The method of dynamic zoning for RFP was used during 3D hydro-dynamic simulation with an objective to study the self-organizing phenomena in high-permeability channels with breakthroughs and in stagnant non-drained areas in zones with super-low filtration velocities, to define the energy scaling depending the filtration stage and period of development, Fig. 6. The multi-scaled geological phenomena in search, exploration and development are diagnosed, monitored and controlled by a system of electro-magnetic, hydro-dynamic and geophysical sensors, momentum energy principles of thermal, electric and magnetic resonance occurrence. The adjusted for critical energy velocity of choke gas and water breakthrough, for diffusive filtration, dynamic counter-flow valves weaken or complete stop the displacing agent supply into water-encroached  Physical and mathematical grounding of porous space block structure and resonant wave effects of natural and technogenic character opens new possibilities in search, exploration and pool modeling, in evaluating the potential flow-rates, diagnostics and control of well fluid inflows, economic optimization of pool development in view of geological risks.