1. Introduction

               

In nature, technological and transportation structures there are spatial branched systems characterized by certain parameters that vary in space as well as in time.

Typical examples of such systems in nature are circulatory systems of living organisms. A circulatory system has a branching, tree-like structure: firstly, large arteries branch to capillary vessels, and then the capillary vessels are embodied in the large venous vessels. The movement of blood through such a system is characterized by various parameters - pressure, velocity, flow rate, diameter and shape of the cross-section of a blood vessel, etc.

In industry, the most famous examples of such objects type are pipeline systems: there are trunk pipelines for transporting products for hundreds and thousands of kilometers; field pipelines for gathering product from production sites into the trunk pipelines; process pipelines for transporting substances within the plant; heat, water and gas networks ensuring population supply. All these pipeline systems usually have a branched system. The flow in the pipelines is characterized by pressure, temperature, velocity, composition, content of gas-vapor fraction, completeness of cross section isolation, and so on.

Finally, transport systems and roads are another example of branching systems. Their main characteristics are traffic intensity and vehicle speed.

All the examples listed above have one common characteristic: their length far exceeds their transverse dimension. For example, for trunk pipelines, the typical diameter of flow (pipe) is about 1 m, and the length is about 1000 km.

Today description of systems with such a ratio of spatial dimensions is most effectively carried out within the framework of one-dimensional approach by solving one-dimensional equations, when the parameters describing the state of the system vary only in the longitudinal direction. In the transverse direction, all the parameters are assumed to be unchanged. For example, when something flows in a pipeline, all flow parameters in any cross-sectional area are assumed to be constant, although in terms of multidimensional hydrodynamics the velocity on the axis of the pipeline differs from the velocity in the side-wall layer. Such an "equalization" by the over cross section speed is achieved by averaging.

Thus, when solving one-dimensional equations, one-dimensional distributions of proeng.files of various parameters are obtained depending on the distance. Such distributions are usually visualized in the form of one-dimensional graphs that are easily perceived and interpreted in practical implementations.

However, if the parameters in the branching system are described by one-dimensional distributions, then a problem arises. It is difficult to perceive numerous dependencies of parameters on the spatial variable not in one spatial interval, but in the whole set of basic and branching elements.

In the paper, this problem is studied by visualizing flow parameters in a branching pipeline system through which a slightly compressible liquid is transported under isothermal conditions.

 

2. The approach to visualization of data that describes spatio-temporal distribution of parameters in branching systems

               

To visualize flow parameters in branching systems, one can use the dependencies of various parameters on a single argument - distance. This can be done both in the "snapshot" mode for one moment of time and in the "video" mode for a certain time interval.

The situation when there is a single spatial interval with the parameters to be visualized on is the simplest situation. In this case, there is no problem either with the construction or with the perception of the visualized data: functions of one variable are the common data representation.

However, if there is a branched network, for example, such as shown in Fig. 1 (there are 23 linear sections in it) then the problem arises on how to distribute the parameters in the most convenient way.

The simplest and the most obvious way is to establish a correspondence between each linear section and the graph to illustrate the changes of the corresponding parameter along this linear section. With this approach, a set of graphical functions is obtained. These functions can be represented either in separate figures or in one figure. With a small number of branches, you can more or less easily see the distribution of parameters throughout the system. However, with more branching, this seems to be increasingly difficult. An example of such a representation form for the system shown in Fig. 1 is given in Fig. 2a where the dependence of a function f (x) corresponding to a certain instant of time is represented in each of the linear elements forming the branching system. The distance x in this graph is measured from the branch point (for element 1 - from the beginning of the system). In Fig. 2a, color matching was used to correlate various dependencies with the corresponding element of the branching system. It is hard to not acknowledge the fact that the perception of the data visualized in Fig. 2a is extremely difficult.

Another possible way of visualization is to use the following method: to display parameter values directly on the branch scheme. This can be done in various ways:

 - using a color scale on a plane network of a branched system;

- using the color scale on the spatial configuration of the branched system;

- using 3D graphics on the spatial configuration of a branched system.

The last version of the visualization is shown in Fig. 2b. This way of visualization, of course, simplifies the perception in comparison to to the version shown in Fig. 2a.

Visualization with reference to the branching scheme enables the whole perception of the entire flow pattern, and the visual image of various sections with parameters in them is formed almost effortlessly and automatically. However, in this 3D approach it is more difficult to compare the parameters values at different points (compared to the 1D graphical dependence on one parameter). In addition, if we take into account significant variation in the length of the system in different directions, there will be an inevitable problem of perceiving different scales in three-dimensional visualization.

In this paper, we propose to solve this problem using graphs with a tree-like structure. As a research object, the branching pipeline system for the transportation of a slightly compressible liquid is selected for consideration. The main parameters characterizing the flow in such system are pressure, speed, flowrate, pipeline diameter, temperature and friction coefficient. The calculation of the pressure has a particular practical importance as it is the pressure excess of some values that lead to destruction of the pipeline.

               

               

2.1. Visualization using tree-like diagrams

               

Tree-like graphs represent a set of 1D dependencies of the corresponding parameter on the distance at a certain time moment. In this case, one graph shows dependencies for all linear sections in a single coordinate system, where the inlet of a pipeline system corresponds to the origin of the coordinate system, and the coordinate itself is the distance from a particular point to the origin of the coordinate system. Thus, for each branching system, a tree-like graph is constructed, completely corresponding to its topology. An example of this type of visualization is shown in Fig. 2c. This example also applies to the system shown in Fig. 1. Different colors are used to match some 1D dependence to a specific element of a branched system. However, it can be easily noticed that this practically does not affect the perception: in black and white, the visualization 2c would also remain sufficiently informative.

After comparison of these three visualization approaches presented in Fig. 2, we can conclude that the tree-like graph approach has the best visibility in combination with simplicity.

It should be noted that this form of representation has a certain advantage: in the branch points of those functions that are continuous, the graphs are branched like a tree. Therefore, the visualization of parameters occurs in the most natural form: the distribution of parameters in the system is presented in a form convenient for perception (dependence of the function on one parameter) with a reference to all branching points.

As a rule, certain integral characteristics of the flow are of interest when processes in branched systems are modelled. For example, when calculating flows in pipelines, it is important to not only know the current pressure values but also the maximum pressure across all the sections for all process time. It is also important to identify those sections where the permissible pressure values were exceeded.

So, the tree-like diagram can be supplemented by visualization of some integral characteristics:

- on the graph of a tree-like function corresponding to a certain time moment, can be plotted a chart of the maximum (or minimum) values of the corresponding parameter at a given point; for pipelines, this is usually the dependence on the coordinate of the maximum pressures reached in the corresponding time;

- on the graph of a tree-like function corresponding to a certain point in time сan be supplemented by a plot of the maximum permissible values of the corresponding value; for pipelines, such a value is usually set to the maximum allowable pressure values above which the operation of the pipeline can cause its destruction;

- on the graph of a tree-like function corresponding to a certain time moment, dangerous spatial intervals can be demonstrated; at these intervals the maximum permissible values have already been exceeded by this time moment; in case of consideration of the pipeline system, the places of exceeding the maximum allowable pressure at a given time moment can be highlighted, for example by color, directly on the line of the function.

In conclusion, it should be particularly noted that the presented approach to visualization using tree-like graphs is applicable not only to displaying at specified moments of time, but also in the form of animated films. In this case, in our opinion, the visibility of the proposed approach in the data representation increases due to the continuous perception of the entire space-time flow pattern.

In addition, animation visualization can be more visual in case of its implementation in real-time systems, when all the changes are reflected in real or advanced time, for example, in the control centers, from which the real pipeline systems are controlled. Animated visualization is indispensable in analyzing the appearances and development of emergencies.

 

2.2. The software environment used to visualize flow simulation results in the pipeline system

               

To visualize the flow parameters in the branched pipeline system, the software tool "TOXI+Water Hammer" [1] was used. The "TOXI+Water Hammer" software is developed for calculation of parameters in the pipeline systems of arbitrary configuration including all standard types of equipment necessary for pipeline operation: valves of various types, pumps, safety valves of various operating principles and parameters, check valves, branches, junction of pipes of various diameters, vessels. "ТОXI+Water Hammer" allows to describe stationary, transient and emergency processes during the pipeline transportation of stable slightly compressible liquids (oil, gasoline, water).

"TOXI+Water Hammer" does not require any special hardware - calculations with the necessary for practical purposes accuracy can be performed by standard personal computers.

To provide a representation of the results of fluid flow modeling in branched pipeline systems, a special subsystem is developed in the "TOXI+Water Hammer" software (Fig. 3). This subsystem allows to display various flow parameters (density, velocity, pressure, flow rate, etc.) at any time moment of simulation, or in dynamics, as animated movies with varying speed of view.

               

3. Mathematical model of the process

               

For a one-dimensional unsteady isothermal flow the flow of liquid in an elastic deformable pipe (in isothermal approximation) is described by a system of the following equations [2-5]:

continuity equation

(1)

momentum conservation equation

(2)

pressure-density relation (equation of state)

(3)

here r, p, u —density, pressure and velocity of the liquid averaged over the cross section; t —time; x — distance from the beginning of the pipeline; l(Re) — friction factor as a function of the Reynolds number Re=Du/ν; А — cross-sectional area of the pipeline, D — pipeline diameter, g — gravitational acceleration; b — the sine of the slope of the route which is determined by the elevation marks of the pipeline  h(х), ν - kinematic viscosity; ρ₀ — density of liquid at pressure p₀ and temperature of transportation (usually p₀=10⁵ Pa), с – perturbation propagation velocity.

In this model, the following factors are taken into account:

-  non-stationarity of the process;

- change in the cross-section of the pipeline during its deformation;

-  the convective motion of the medium (the second terms in the left-hand sides of (1) and (2);

- appearance and circulation of waves generated by stopping/starting pumps, closing valves (the first term on the right side of (2);

- the presence of friction on the walls of the tube (the second term on the right-hand side of (2);

-  effect of gravity on the flow during the passage of the pipeline route through a terrain with a complex relief (the third term on the right-hand side of (2).

Equations (1) - (3) are supplemented by initial and boundary conditions. The inlet and outlet pressures of the pipeline are set defined as the boundary conditions, these pressures correspond to the pressures of tanks or pumps located at the inlet and outlet of the pipeline. As initial data, the parameters of stationary flow are specified, these parameters  can be obtained analytically from the solution of the system (1) - (3).

When the valve is closed, the boundary condition is zero flow rate.

The system (1) - (3) describes the motion in one linear section. In case of a branched pipeline system, when individual linear sections are combined in a certain sequence, the system of equations (1) - (2) is recorded for each linear section. The boundary conditions at the junctions of two linear sections are set in such a way that the fluxes of mass and momentum from one section to the other one are the same.

To determine l(Re) we use the Colebroke-White relation [6], which links the friction coefficient l to the Reynolds number Re and the pipeline characteristics (diameter D and roughness k):

(4)

where k is the pipeline roughness

Thus, the system of equations (1) - (4) allows to consider fully real characteristics of the equipment operating on the pipeline: pipes, pumps, valves, etc.

To solve the system of equations (1) - (4), we used Godunov-type method previously described in [3, 7-8].

 Godunov-type method is applied to systems of equations written in a divergent form. For example, for a system of differential equations:

(5)

when its numerical solution is constructed its difference analogue is as follows:

(6)

here U – vector of unknown quantities representing a function of space x and time t,  Φ – flux of quantities U,  –  the discrete analogue of quantities U referred to the i-th discret cell of the space (the coordinate of its center is (i+1/2)∙∆x)  at the n-th time step (in n∙∆t time moment),  - discrete analogue of Φ, referred to the boundary of the i-th discret cell with the coordinate i∙∆x (right cell boundary),  и   - discrete steps in time and space.

Thus, using (6) and knowing value () on the n-th time step and the fluxes of this quantity through the boundaries of the corresponding cell ( and ), it is possible to find the values of U((n+1)∙∆t, (i+1/2)∙∆x)) () on each new n+1-th time step.

So, in Godunov's approach, the quantities  are found from the solution of Riemann`s problem, the problem of the discontinuity decay. In this problem at the initial moment of time, one part of the half-space is occupied by a medium with parameters , and the other one - with parameters . At subsequent time moments, two waves begin to propagate into two half-spaces: one wave into one half-space. And each of these two waves will convert the initial states  and  into the new states, and these new states are "stitched together" according to certain parameters (in gas dynamics these parameters are pressure and velocity). The problem of discontinuity delay can be solved analytically either accurately or approximately. Knowing the parameters in the «stitching» region one can calculate the flows  and , i.e. actually solve (6).

It should be noted that Godunov-type method today is one of the most widely used methods of numerical solution of systems of hyperbolic equations due to its high physical validity.

               

4. Problem statement

               

When the valves are closed and the pumps are stoped or in the process of the operating mode change, the product in the pipeline slows down or stops completely. In its turn, flow slowing-down upstream the valves or upstream the pump leads to the pressure increase, formation of  increased pressure area in the pipeline and the subsequent propagation of the increased pressure area along the pipeline route. Such propagation of compression waves in the pipelines is called water hammer [9].

A water hammer is a significant hazard to the integrity of pipelines since a pressure increase can lead to pipe rupture or formation of dangerous defects [10].

For example, a water hammer was one of the causes of the pipeline rupture accident at the Sheskharis tank farm (Novorossiysk) and oil release to the port water area [11].

The danger of a water hammer in pipeline systems is worsened by the fact that, as a rule, there are one-dimensional flows in pipelines. Since there is one-dimensional geometry, compression waves propagating through the pipeline are attenuated rather slowly what leads to formation of potentially dangerous zones over tens of kilometers.

Water hammer investigations have been carried out for dozens of years, beginning with the work of N.E.Zhukovsky [9] and ending with modern monographs, thesis and articles [1, 10, 12, 13, 14]. Typically, water hammer is considered in a single pipeline. However, real pipeline systems often have a more complex topology, including branching and looping (branching with a return to the main pipeline).

Branching structure of a pipeline can have a significant impact on the development of the water hammer flow.

The water hammer problem in a branching system is extremely important from a practical point of view. Indeed, the consequences of a water hammer are the most catastrophic in cases of quick valve shut-off. Such a quick (in 3-5 seconds) operation of the valve and complete cutoff of the flow often take place while shipping at sea terminals. It is important to cut off the flow quickly because the release can pollute the aquatic environment what leads to the huge environmental damage. Shipment of large volumes at high rate is often carried out through several pipelines in order to reduce the total shipment time. As a result, there may be an intense water hammer in a branched system.

Below we consider the problem of generation and development of a water hammer in a pipeline with one and several branches.

               

4.1 Single branch pipeline

               

Let us consider the following model configuration of the branching pipeline (see Figure 4). Distance from the inlet of the pipeline (point A) to the branch point (point O) is equal to 20 km, the pipe has a diameter of 1 m, at point O, the pipeline branches into two identical pipes 20 km long and 1/(2)^0,5 m in diameter each, that is, the cross-sectional area of the two branches is equal to the cross-sectional area of the main pipe. All the pipes have the same roughness of - 0.3 mm and are laid on flat terrain. A pressure at the inlet of the pipe is equal to 1.0 MPa. Pressures at the end of each branch pipes are equal to 0.3 MPa. Also, there are valves which shut off the flow completely within 0.1 sec

               

The scenario under which valve №2 cuts off the flow is considered.

In the calculation, it was assumed that the propagation celerity of the wave coincides with the speed of sound and is equal to 1300 m/s. Density of the transported product is equal to 840 kg/m³ at p₀=10⁵, and kinematic viscosity is equal to 2 10^-6 м²/с.

The problem was solved using a grid with a spatial step of 400 m.

               

4.2 Pipeline system with several branches

               

To solve the second problem the following configuration of a pipeline with three branches was considered (see Figure 5). The main pipeline consisted of 4 sections (№№ 1, 3, 5 and 7 in Fig. 5), each with a length of 10 km. The diameter of this main pipeline decreased from section to section:

- tube in section №1 had a diameter of 1 m;

- tube in section №3 had a diameter of 0,7071 m;

- tube in section №5 had a diameter of 0,5 m;

- tube in section №7 had a diameter of 0,3536 m;

In the system under consideration there were three branches each with a length of 5 km (see Figure 5):

- branch №2 departs from the main pipeline at the junction of sections №1 and №3; this branch has a diameter of 0.7071 m; pressure at the outlet of branch №2 is equal to 0.5 MPa;

- branch №4 departs from the main pipeline at the junction of sections №3 and №5; this branch has a diameter of 0.5 m; pressure at the outlet of branch №4 is equal to 0.3 MPa;

- branch №6 departs from the main pipeline at the junction of sections №5 and №7; this branch has a diameter of 0.3536 m; pressure at the outlet of branch №6 is equal to 0.1 MPa;

 

The pressure at the entrance of the pipeline system is equal to 1.0 MPa and the pressure at the end (at the outlet of section №7) - 0.1 MPa.

All pipes are laid on flat country.

All pipes have the same internal surface roughness - 0.03 mm.

Water is pumped through the pipeline at a temperature of 20°C. The initial density (at 0.1 MPa) was assumed equal to 1000 kg/m³ and kinematic viscosity - 8.9 10^-6 m²/s.

Propagation celerity of the wave equal to 1480 m/s was used in the calculation.

The following scenario of a water hammer was considered. With a steady flow in the pipeline a gate valve was activated at the end of the main pipeline (at the end of section №7). The time for complete cut off of the flow at the valve shut off was equal to 2 s.

The problem was solved by using a grid with a spatial step of 500 m.

               

5. Calculation Results

               

5.1 Water hammer calculations in the pipeline with one branch: operation of a single valve

               

Figure 6 shows the pressure proeng.files at different time moments. In this figure, the main pipe corresponds to the section from the inlet (point A) to valve №1, the branch in this figure corresponds to the section from point O to valve №2.

               

In Fig. 6a one can clearly see how the wave that was formed in the branch (see Fig. 6a, time moment - 10 s) enters the main pipeline (see Fig. 6a), the time moment is 20 s).

Figure 6b shows how after the stoppage of wave circulation in the system a slow increase in pressure occurs across all sections of the pipeline. Most significantly, the pressure rises from point O to valve №2, where the flow stops completely. In the other branch (from point O to valve №1) and in the main pipe (from the beginning to point O) the pressure increase is less significant (by the relative value).

Note that the maximum attainable pressures will also be affected by the change in the closing delay time of the valves on the pipeline branches. Only an accurate and effective visualization system can enable proper analysis of situations affected by various closing delay time of the pipeline branches valves.

               

5.2 Water hammer calculations in a pipeline with several branches

               

The results of the water hammer calculations in a pipeline system with several branches in the formulation described in 4.2 are presented below.

The distributions of pressure, velocity, density and mass flow rate for two stationary states of the system (before the valve is closed at the end of Section №7 and after it is closed and a new pumping regime is established) are shown in Fig. 7-10. In Fig. 7-10 it can be seen that after closing the valve the flow in section №7 is stopped and the pressure is equalized along its entire length. At the same time, the mass flow rate in the system drops slightly. This is well seen in Fig. 10 (a) and (b) where the flow rate of the whole system determined by the flow rate in section №1 falls from 1630 kg/s to 1610 kg/s. There is some redistribution of flows between the long ten-kilometer sections and the short five-kilometer sections. The flow rates in sections №3 and №5 decrease while the flow rates in sections №2, №4 and №6 increase (see Figure 10). In accordance with the change in flow rates, such parameters as velocity and pressure also change. The velocity in the short sections increases while in the long sections it decreases (see Figure 8). Gradients of pressure increase in short sections and decrease in long ones (see Figure 7).

From the presented data, it is also easy to see a slight change in the density of the flow (see Fig. 9).

Intermediate states in the system are shown in animations presented in Fig. 11 and 12 where pressures and velocities are visualized.

In these animations, the beginning of the valve shut off occurs in 290 s.

In Fig. 11 one can see that after 6 seconds after the closing of the valve, the wave from the end of the pipeline approaches the first branch from the end (at the end of section No. 5). This wave begins to unload in sections №5 and №6. The pressure on the valve at this time moment reaches 1.8 MPa. Approximately 9 seconds after the valve closes, the wave reaches the end of section №6 and begins to unload through this outlet. By this time, the pressure on the closed valve reaches almost 1.9 MPa and at the first branch from the end - 0.7 MPa.

The compressive wave that goes through section №5, reaches the second branch from the end and begins to unload in sections №3 and №4 after about 12 seconds after the valve shut off starts. By this moment, the pressure on the closed valve has reached 1.9 MPa (10 seconds after the valve shut off starts) and begins to decrease due to the mass outflow from section №7 (the velocity for Section №7 is less than zero, see Figure 12). The pressure at the first branch from the end is equal to 1.1 MPa. At the same time, due to the larger (compared to section №6) diameter of section №5 the pressure in section №5 increases substantially up to 1.17 MPa.

After 18 seconds after the valve shut off starts, the water hammer wave begins to unload into sections №1 and №.2. At the same time, pressure in Section №3 increases, and pressures in sections №5 and №7 decrease. A traveling wave is formed in the main pipeline. The local maximum pressure in this wave is equal to approximately 0.86 MPa. This maximum is located at 14 km from the beginning of the pipeline at 24 seconds after the valve shut off starts.

The minimum pressure on the valve is reached approximately in 31 seconds after the valve shut off starts. After that, the outflow of mass from section №7 is replaced by its inflow in this section from section №3 (the velocity in section №7 grows higher than zero, see Figure 12) and the pressure begins to increase in section №7 on the closed valve. Also, pressure is rising throughout the pipeline.

The pressure increase on the valve continues up to 51 s after the start of its closing. The pressure on the valve reaches 0.76 MPa. Then the pressure on the valve begins to decrease. With some time delay (about 3 s), the pressure begins to decrease in the whole pipeline system.

After reaching the next minimum pressure (0.12 MPa) on the valve (73 s after the start of its closing), another pressure increase begins.

Such oscillations in the flow are observed along the entire length of the pipeline. The greatest amplitude of oscillations takes place on a closed valve. In the main part (sections № 1, 3, 5, 7) the amplitude of the oscillations increases with the approach to the pipeline beginning. Pressure peaks in damped pulsations on the valve are achieved after 94 s, 133 s, 153 s, 181 s and 213 s (after the valve starts to operate). They are 0.72 MPa, 0.54 MPa, 0.46 MPa, 0.50 MPa and 0.49 MPa respectively. Pressure minimums are reached after 118 s, 146 s, 164 s, 200 s and 222 s (after the valve has started to operate) and they are 0.31 MPa, 0.44 MPa, 0.41 MPa, 0.42 MPa and 0.45 MPa.

Approximately 240 s after the closure of the valve, the pressure and velocity oscillations in the system are practically damped.

All these pressure oscillations at the end of Section №7 are shown in the animation in Fig. 13.

 

 

 

 

 

 

 

 

Fig. 11. Time variation of pressure in the branched system.

Fig. 12. Time variation of velocity in the branched system.

Fig. 13. Time variation of pressure at the closed valve latch (the end of section №7).

               

6.  Conclusions

               

An approach that allows visualization of flow parameters in branched systems is proposed in this paper. Graphs with a tree-like structure (one-dimensional dependencies of the corresponding parameters on the distance from the beginning of the system) are used.

Such visualization combines the simplicity of the traditional representation of one variable functions and the convenience of perceiving these functions for branched systems.

The proposed approach to visualization is illustrated by examples from the field of pipeline transport, for example, the solution of the water hammer problem in a branched system.

The results of calculations with visualization in the form of a tree-like structure allow:

- to observe the propagation of compression/rarefaction waves along a branching system including dynamic visualization ("video" mode);

- to observe the unloading of compression waves at the branching points into the pipelines connected at the branching point and at the inlet/outlet of the pipeline system;

- to perceive in the optimal way the interaction of circulating waves, including compression.

The last one is especially important since it allows to find out quickly the cause of the emergence of zones with high pressures appearance. The proposed visualization (especially in the dynamic mode of "video") allows to identify interacting waves in real time and to understand instantly which events and what combination of factors has led to the appearance of high-pressure areas.

Based on the approach to modeling and visualizing flows in the branched pipeline systems that was presented above "TOXI+Water hammer" software has been developed. It helps to simulate flows in pipeline systems using a visualization subsystem based on tree-like graphs.

Thus, the modeling of flows in pipeline systems shows that only when visualization of the pattern of wave propagation in pipes is available it is possible to perceive the flow pattern directly. Besides, this perception is the most convenient when tree-like graphs are used.

The work was supported by the Russian Science Foundation, grant No. 16-19-00188

               

References

               

1. Sverchkov A. M., Agapov A. A., Sumskoi S.I., Egorov A. F., Sof’in A.S. Modelirovanie avarijnyh utechek na magistral'nyh nefteprovodov [Modeling  of emergency releases at oil trunk pipeline] //Bezopasnost truda v promyshlennosti, № 9, 2014, pp. 50-53 [in Russian].
2. Lurie M.V. Matematicheskoe modelirovanie processov truboprovodnogo transporta nefti, nefteproduktov i gaza [Mathematical Modeling of Gas, Oil and Oil Product Transport].  Moscow: Gubkin University Press, 2003 [in Russian].
3. Sumskoi S. I., Sverchkov A. M., Lisanov M. V., Egorov A. F.   Modelling of non-equilibrium flow in the branched pipeline systems // Journal of Physics: Conference Series, Vol. 751, №1, 2016, 012022 (doi:10.1088/1742-6596/751/1/012022)  (http://iopscience.iop.org/article/10.1088/1742-6596/751/1/012022 )
4. Sumskoi S. I., Sverchkov A. M., Lisanov M. V., Egorov A. F.  Simulation of systems for shock wave/compression waves damping in technological plants // Journal of Physics: Conference Series, Vol. 751, №1, 2016, 012023 (doi:10.1088/1742-6596/751/1/012023) (http://iopscience.iop.org/article/10.1088/1742-6596/751/1/012023)
5. Sumskoi S. I., Sofin A. S., Lisanov M. V. Developing the model of non-stationary processes of motion and discharge of single- and two-phase medium at emergency releases from pipelines// Journal of Physics: Conference Series, Vol. 751, №1, 2016, 012025 (doi:10.1088/1742-6596/751/1/012025) (http://iopscience.iop.org/article/10.1088/1742-6596/751/1/012025)
6. Mirkin A.Z., Ucinsh V.V. Truboprovodnye sistemy [Pipeline systems]. M: Khimia, 1991. [in Russian]
7. Gubin S.A., Gubina T.V., Sumskoi S.I. Modelirovanie perehodnyh i avarijnyh processov v magistral'nyh truboprovodah s pomoshh'ju metoda S.K.Godunova [Modeling of Transient and Emergency Processes in Trunk Pipeline Systems using Godunov’s approach] // Bezopasnost truda v promyshlennosti, №10, 2013, pp.66-71. [in Russian]
8. Sumskoi S.I., Sverchkov A.M. Modeling of Non-equilibrium Processes in Oil Trunk Pipeline Using Godunov Type Method // Physics Procedia, Vol. 72, 2015, рр. 347–350.
9. Zhukovsky N. E. O gidravlicheskom udare v vodoprovodnyh trubah [About water hammer in water pipes]. Moscow:-Leningrad: State Publishing of Technical and Theoretical Literature, 1949. [in Russian]
10. Gusein-zade M.A. Osobennosti volnovogo techenija v trubah. Gidravlicheskij udar. [Features of wave propagation flows in pipeline. Water hammer]. Moscow: Oil and Gas, 1999. [in Russian]
11. Guschin V.V. Avarija na truboprovode № 9 neftebazy «Shesharis» AO «Chernomortransneft'» [The Accident at Pipeline #9 of  “Shisharis’ Oil Tank Farm”] // Bezopasnost truda v promyshlennosti, № 1,  1998, pp. 8–9. [in Russian]
12. Truboprovodnyj transport nefti [Oil Pipeline Transportation]. In 2 volumes // Ed/ S.M.Winestock: Nedra - Business Center, 2002-2004. [in Russian]
13. Arbuzov N. S. Obespechenie tehnologicheskoj bezopasnosti gidravlicheskoj sistemy morskih neftenalivnyh terminalov v processe naliva sudov u prichal'nyh sooruzhenij (teorija i praktika) [Providing of technological safety of the hydraulic system of marine oil terminals in the process of loading vessels at mooring facilities (theory and practice)], Moscow: Thesis for the degree of Doctor of Technical Sciences, 2014.  [in Russian]
14. Sumskoi S. I., Sverchkov A. M., Lisanov M. V., Egorov A. F.  Simulation of compression waves/shock waves propagation in the branched pipeline systems with multi-valve operations // Journal of Physics: Conference Series, Vol. 751, №1, 2016, 012024 (doi:10.1088/1742-6596/751/1/012024) (http://iopscience.iop.org/article/10.1088/1742-6596/751/1/012024)