MSC 65R20
DOI: 10.14529/ m m p160405
ON APPLICATION OF THE STRUCTURE OF THE NONLINEAR EQUATIONS SYSTEM, DESCRIBING HYDRAULIC CIRCUITS OF POWER PLANTS, IN COMPUTATIONS
A.A. Levin, Energy Systems Institute SB EAS (Irkutsk, Russian Federation), [email protected],
V.F.Chistyakov, Institute for System Dynamics and Control Theory SB EAS (Irkutsk, Russian Federation), [email protected],
E.A. Tairov, Energy Systems Institute SB EAS (Irkutsk, Russian Federation), [email protected]
This paper presents the results of applying the methods of the hydraulic circuit theory in the interactive modeling of hydrogasdynamic and thermal processes that occur in the equipment of thermal power plants. The problem statement of flow distribution in the energy plants with different pressure drop laws in the branches of complex gas air and steam-water ducts is formulated. The research shows that the application of traditional methods of hydraulic circuit theory is challenging for such problems. Some aspects of solvability of the related systems of nonlinear equations are studied. The numerical methods for solving these systems as applied to the problems that require calculations in real time are considered. A computation scheme is proposed. The scheme makes it possible to reduce the initial statement of the problem to the classical scheme of the nodal pressure method. The method of decomposition of the hydraulic circuit configuration into interconnected circuits of smaller dimension is considered to reduce the computational effort. The results of tests that demonstrate high reliability of the method developed for the flow distribution calculation are presented.
Keywords: hydraulic circuits; mathematical models; ducts of power plants; systems of nonlinear equations; Newton method; fixed point iteration method.
Introduction
At the present time, a successful application of the decomposition principles in the modelling of complex power plants, as well as the success in implementing of the separate description of the slow heat-mass exchange and rapid hydrodynamic processes, facilitates development of universally applied techniques for constructing of mathematical models. Therefore, development of the generalized algorithm for the flow distribution calculation in ducts of power plants on the basis of the hydraulic circuits theory is of interest [1]. The hydraulic circuit is a planar graph of m nodes and n edges. The edges represent parts of the power plant through which various types of fluid (water, steam-water, gas, etc.) move. The ambiguity of mathematical descriptions for the elements of schemes for thermal power plants is mentioned in [2]. The authors argue that it is necessary to unify the mathematical description for the elements to develop a universally applied computational software. Thus, it is important to address the problem of developing a generalized algorithm for the flow distribution calculation in ducts of power plants, where the fluid moves along the pipelines of the various length and profile that have connections joints, control valves, and head sources (pumps, ventilators). The design project for the power generating unit
has hundreds of nodes and edges which results in high dimension systems of equations with respect to flows and pressures. Methods of the hydraulic circuits theory are designed to solve such systems. In the paper [3], we proposed a modified nodal pressure method which enabled us to take into consideration the impact of the fluid phase change on closing relations. The algorithm developed combined properties of the Newton methods and the fixed point iteration method and allowed the reduction of the problem to the common form of the nodal pressure method. Investigation of the duct modes different to the regular one was left behind the framework of the previous research. Therefore, the flow distribution calculation problem aims at solving the following sub-tasks:
1. Consideration of the pressure-drop law specificities in the duct sections with the dust-laden flow (in the pulverized-coal system).
2. Efficient consolidation of the calculation objects interconnected within the power generating unit (boiler, turbine, regeneration system).
1. The Hydraulic Circuit Model with Various Pressure-Drop Laws in Pipelines
Any hydraulic circuit satisfies the conditions
AX = Q, ATP = Y,
where A is an (m x n) incident matrix with the elements 0,1, — 1; T denotes matrix transposition, X = (x1, x2,..., xn)T is the pipeline flow vector; Q = (q1,q2, ■■■,Qm)T is the inflow in the nodes of the hydraulic circuit. The vector of pressures P has the desirable as well as the given values of pressures. The classic statement of the problem [1] assumes that the pressure is given at one of the nodes. This paper offers the opportunity to set pressure values in several nodes. Therefore, we consider the vector P as decoupled into P = (p1,p2, ■■■,Pi)T and P^ = iPi+i,Pi+2, ■■■,Pm)T with the desirable and the given
—T A
vector P. The vector Y = (y1,y2, ■■■,yn)T represents heads (pressure differences at the ends of the pipelines with a corresponding number).
Such vector equations reflect the first and the second Kirchhoff's laws: 1) at any junction, the sum of inflows is equal to the sum of outflows; 2) the sum of the pressure differences around any closed hydraulic circuit is zero. The first law entails that Yj=o Qj = 0.
In what follows, we need the matrix A, formed by throwing out the last m — l rows from A, and the vector Q, formed by throwing out the last m — l elements from the vector Q. Write down the new system
AX = Q, ATP + ATP*m = Y, (ATAT) = AT. (1)
The matrix A is full rank: rank A = l In the classic statement of the problem, the matrix A is obtained by deleting an arbitrary row from A.
The Kirchhoff's laws equations are complemented by the laws for the movement of various fluid along the pipelines. This paper assumes that the flows in the pipeline with the number v E [1, 2, • • • ,n] are connected with the pressure at the ends of the pipelines via
hv + fv (Pi (v),pj(v)) sv (Pi (v) , pj(v), xv )xv \xv \, (3)
where i(v),j(v) are nodes of the given pipeline, sv(Pi(v), Pj(v), xv) are the resistance functions (in many research works they are assumed to be constant: sv(pi(v),pj(v),xv) = sv = const). Note that the form of (3) provides solvability with respect to xv for any signs of hv + fv(pi(v),Pj(v)) and sv(pi(v),pj(v),xv) [1]. In applications, it is often assumed that
fv (Pi(v), Pj(v)) = Vv = Pi(v) - Pj(v). (4)
Recall that yv is a head in the given pipeline. We assume that the functions in (3) are continuous and have the following properties:
signfv (Pi(v),Pj(v)) = sign[Pi(v) - Pj(v)}, fv (Pi(v),Pj(v)) = 0 & Pi(v) - Pj(v) = 0; (5) there exist constants smin > 0, smax > 0, such that
Smin < sv (Pi(v),Pj(v),Xv) < smax ^Pi(v), Pj(v), Xv £ R, Vv £ [1, 2, ••• ,n}. (6)
In other words, the fluid moves from the node with a bigger pressure to the node with a smaller pressure, and the resistance functions are bounded for any values of the arguments.
Suppose that the elements of H are zero in the sections of form (4). Multiply the equations for the nonlinear sections by yv and divide by fv(^i(v),Pj(v))■ The closing relations will have the form
yv — sv (Pi (v), Pj(v), xv)xv\xv sv (Pi (v),'Pj(v) ,xv) f. / \sv (P
i(v), Pj(v) , xv ). (7)
Jv (Pi(v),Pj(v))
Note that sv(Pi(v),Pj(v),xv) > 0 when yv = 0 due to (5) and (6). To simplify further reasoning, suppose that
lim Sv(Pi(v),Pj(v),xv) = sv = const,
yv —
although analysis of the real-life models does not always justify it. Taking into account (1) and (7), we obtain
AX = Q, ATP + ATP*m = Y, Y + H = S(P, X)XX, (8)
where S(P,X) = diag{si(pi(i),pjW,Xi), §2(pi(2),Pj(2),x2), sn(pi (n) ,pj(n), xn)}) X = diag{|x1|, x!,---, \xn\}. We can omit the vector equations in (8), and as a result arrive at
AX = Q, ATP + H = S(P, X)XX, (9)
where H = P^ + H. System (9) has l + n equations with l + n variables, where n is the element of X and l is the element of P.
Hence, we obtained an analogue of the system that corresponds to the nodal pressure method in the hydraulic circuit theory [1]. Now investigate some properties of the system (9). In the monograph [1], the unique solvability of (9) is proved for the matrix S(P,X) with constant diagonal elements. In our situation, the following statement holds:
Lemma 1. Let (9) have no flows and no sources of head and pressures: Q = 0, H = 0. Then (9) has a unique zero solution.
X AX = 0 X
equation in (9). We obtain
(X, ATP — S(P, X)XX) =
n
= (P, AX) — (X,S(P,X)Xx) = — £ S(Pi(v),Pj(v),xv)x2v\xv\ =0. (10)
v=0
X=0
A ATP = 0
□
Theorem 1. Let system (9) have no inflows and no heads: Q = 0, H = 0, and the different pressure values be given only in the dangling nodes. Then (9) has at least one solution.
Proof. Choose an arbitrary vector S = (s1, s2, ■■■, sn) from the cube [smin,smax]n and consider (9) with the matrix S(P,X) = diag{s1 ,s2, ■■■,$„}. According to [1], this system has a unique solution, and the vector-functions (P(S),X(S)) are continuous with respect to S. Having searched through all the values from the cube [smin, smax]n, we derive bounded sets P, X, where the values of all solutions are located. The biggest and the smallest values of the pressures pmin, pmax are reached in the dangling nodes. Some tedious reasoning leads to the conclusion that the assumption of existence of a node with bigger or smaller values than pmin, pmax violates the first Kirchhoff's law in this node. Set K = [pmin, pmax]m, find S(P, X) and solve the resultant system. We conclude that all the values of the set P are located in the cube K. Correspondingly, we can point out a closed convex set that contains X
Brauer's theorem, the function corresponding to such a set has at least one fixed point in this set.
□
2. The Flow Distribution Models for the Power Plants Dynamic Calculations
The power plants models, which have been considered by the authors, have three ducts: water-steam duct, gas-steam duct, and fuel duct. Consider a specific hydraulic circuit for a water-steam duct:
Example 1.
In the node 14, the pump pumps water into the boiler and it is being heated along the duct from node to node 14 — 1 — 2 — 3 — 4 — 5 — 6 — 8 by hot gases and emission from the combustion unit. Water flows through the pipelines 1, 2, 3, 4, 7, 14; the water-steam mixture flows through 5 and 6; steam fills all the other pipelines. Steam goes from the boiler to the node 6 through the pipeline 21; the control valve is located in the pipeline 9. Steam goes from the node 7 to the turbine governing stage (pipeline 10). Afterwards, a portion of steam flows along the pipelines 12 and 11 to be heated again in the boilers 1 and 2. In the pipeline 11, steam is partly used to heat water in the high pressure reheater. Moving along the pipelines 12, 13, 15, 20 and being additionally
The principal design project: straight-through boiler-turbine
reheated, steam goes to the turbine stage 16 and then, after the necessary portion of steam gets into the low pressure reheater. it reaches the last turbine stage 17 and inflows into the capacitor (Iv). Steam turns into water in the capacitor, and water, reheated and purified in the regeneration system, flows into the node. The pipelines 2, 7, 14 have automatic regulators. The inflow q(t) is used to simulate the variation of the fluid density in the pipelines 5 and 6.
Figure shows the principal design project for the hydraulic circuit. To simulate the work of the real-life power equipment, the project must have hundreds of nodes and edges. The dynamics is defined by the fact that the data, in particular, the pressure in the boundary nodes, is a solution to some differential equations P^(t) = Z(P^(t), P,X, ...,t), ' = d/dt. During implementation, the differential equations are replaced by difference relations and
P, X
The step by step calculation of the flow distribution in dynamic models of power generating units implies splitting of the closing relations into four subsets, where the pressure-drop laws in the pipelines are defined by the formulas:
fv (Pi(v),Pj(v)) = Pi(v) - Pj(v) = sv \xv \xv ,v G Ii, (11)
fv (Pi(v),Pj(v)) = Pj(v) - Pj(v) = sv\xv\xv,v G Ij , (12)
fv (Pi(v),Pj(v)) = CoPl(v) + CiPi(v)Pj(v) + CjPj(v) = sv\xv\xv ,v G I3, (13)
fv (Pi(v),Pj(v)) = Pi(v) - Pj(v) = sv (Pi(v),Pj(v),xv ,...)\xv\xv ,v G I4, (14)
where c0 = —0, 09; ci = 1, 09; c2 = -1,sv(Pi(v),Pj(v),xv,...) = k(1 + x^(p'/p'')/p' is the resistance in the pipelines with boiling water; p',p'' stand for water and steam density at saturation; x is the flow quality; k, ^ are the coefficients reflecting the impact of the flow structure on the frictional head [4].
Here, (11) describes water flows in pipelines, (12) and (13) define steam movement through the turbine sections [5] and control valves [6]. Equation (14) describes the water-steam flow on heating surfaces and in pipelines with account of (4).
Modify (12) and (13) using (7). We obtain
fv(Pi(v),Pj(v)) = Pi(v) — Pj(v) =-Sv-\xv\xv,v e I2;
Pi(v) + Pj(v)
fv (Pi(v),Pj(v)) = Pi(v) — Pj(v) = --\xv \xv,v e I3;
Pi(v) — 0, 09Pj(v)
These equations were successfully employed in the nodal pressure method [3] when calculating the flow distribution in the ducts of the thermal power plant.
The flow distribution calculation for the pulverized-coal system possesses some features related to drastic changes of the fluid properties. In order to account for the nonhomogenity of the air/dust flow, we consider dust concentration [7]. Resistance of mills, downflow fuel drying, separators etc. has the form
s = so(1 + by),
where b is a coefficient, y is the dust concentration mixture defined in general as
= (1 — a A W)(1 — 6)KC y gi(1 + xKpr) + a A W, ( )
where a is a share of withdrawn humidity, W is fuel humidity, Kpr, Kc are some coefficients, x is a suction share at the section, g1 is the amount of drying agent per 1kg of humid fuel delivered to the inlet of the pulverized-coal system:
273 VQ AW) Yo
g1 = '
/ 273 Vg AW \ \273 + tg 1000Br — 0.804)
273 + tg 1000Sr 0.80 V 1 + Kpr
where Y0 is the mixture relative density [kg/m3], Vg is a gas flowrate. To identify the impact of the fuel flowrate dynamics on resistance, and therefore on the flow distribution in the pulverized-coal duct, consider (15) as applied to the pulverized fun. The fuel concentration required to calculate the pulverized fun performance is a function of the fuel flowrate Br, which can be readily seen from the relation [7]:
1000Br (100 — Wr )Kc
y =
(100 — WPl)Vj" '
Therefore, the resistance in the pulverized-coal section of the duct is, similarly to (14), a function of the fluid flow.
3. Numerical Methods and Application of Decomposition to Calculations of Complex Ducts
Since the elements of the matrix S(X, P) are not, in general, differentiate, we cannot apply the Newton method in its classic form. We suggest the following tactics. Consider a system
F (z) + tf(z )G(z) = 0, (16)
where F(z), G(z) are differentiate vector-functions, ^(z) is a non-differentiable matrix, z is a desirable vector from R^. Start an integrating process
zj+i — zj -
dF (zj) + ) dGj
dz + ) dz
[F (zj ) + Ф(^- )G(zj)], j — 0,1, 2,... . (17)
Method (17) is a combination of the Newton method and the fixed point iteration method. We have carried numerical experiments with model examples. Sometimes, the method preserved the quadratic convergence, but occasionally showed geometric convergence.
Theorem 2. Let in (9) the vector-functions F(z), G(z) and the matrix ^(z) be differentiable in Rv, and let (9) have a solution z* satisfying
det
provided that the vector-function
dF(z*) T . ,,dG(z*)
— 0, (18)
Ф^) — z -
dF (z>+*(z) dG(zГ-1
[F (z) + ¥(z)G(z)]
dz dz
possesses a contraction property in some neighborhood of 6 = {z : \\z* — z\\ < g}:
||$(zi) — $(z2)\\ < a\\zi — zv\\, a < 1, Nzx,z2 £ 6. Then, if z0 £ 6, method (17) converges geometrically.
Proof. Indeed, it follows from (18) that the matrix dF(z)/dz + ^(z)dG(z)/dz is invertible in some neighborhood of z*. The iterative process zj+1 = $(zj), j = 0,1, 2,..., converges geometrically [8] and coincides with the iterative process (17).
□
We applied this approach to systems of the form (9) and replaced the Jacobian by the matrix similar to the one in (17). The resultant iterative method had the form
(Xj+i\ = (Xj\ (S(Xj ,Pj )Xj AT\~ (S(Xj, Pj )Xj Xj + AT Pj — . U+J V Pj) V A oM AXj — Q ),j = 0,1'....
(19)
We continued the process until
max(\xvJ+i — Xvj\, v = 1,...,n; \pv,j+i — pv,j\, n = 1,l) > t,
where t is a given number. The initial approximation (X0, Po) was taken close enough to the solution, and then we could observe the almost quadratic convergence.
Consider another simplification for the matrix inversion in (19). When applying the method, at each step we solve a linear system
rjf)X T) * — P-pj%j-fpj - — (j + " (20)
Using the diagonal form of the block S(Xj, Pj)Xj, multiply the first block line in (20) by A[S(Xj, Pj)Xj]-1 and deduct it from the second block line. We obtain
fS(Xj ,Pj )X j AT_ \
\ 0 A[S(Xj,Pj)Xj]-1AT) -
j j j
The matrix A[S(Xj, Pj)Xj}-lA is called the Maxwell matrix and it is invertible due to
A Xj
the matrix from (20) is invertible, and therefore we can perform process (19). Moreover, if we apply this technique, the number of arithmetic operations required to solve (20) is proportional to l3 and not to (n + l)3.
This method was used to calculate the steam-water and gas-and-air ducts of the power generating unit with two straight-through boilers, two control valves, and a turbine with seven sections. The parameters of the hydraulic circuit were: n = 125, m = 85. The numerical results demonstrate stable convergence of the nodal pressure method.
nm
In this situation, we suggest decomposition of the original design project into graphs of smaller dimension with a few number of common nodes. As an example, consider a hydraulic circuit of two boilers and three turbines with 4 common nodes. The calculation algorithm has the form:
1. Set the initial approximation in the nodes common for the turbine and the boiler: Pi,Pi,P3-
2. Solve the following system for the turbine and the boiler with pi,pj,p3:
AiXi = Qi, AT Pi + Hi = Si(Pi,Xi)X iXi,
AjXj = Qj, ATPj + Hj = Sj(Pj, Xj)XjXj.
3. Write down the system using the first Kirchhoff's law:
(fi(Pi,Pj,P3) = xi,k (Pi,Pj,P3) - xi,t(Pi,Pj,P3) = 0 fj(Pi,Pj,P3) = xj,k (Pi,Pj,P3) - xj,t(Pi,Pj,P3) = 0 f3(Pi,Pj,P3) = x3,k (Pi,Pj,P3) - x3,t(Pi,Pj,P3) = 0
4. Set the increment and find the Jacobian
M = (ui,uj,u3),
where
ui = [F(Pi + AP,Pj,P3) - F(i,Pj,P3)}/ △ P Uj = [F(Pi,Pj + AP,P3) - F(i,Pj,P3)}/ △ P u3 = [F(Pi, ,Pj,P3 + △p) - F(i,Pj,P3)}/△ P
5. Solve the system
Mh = F(pi,'j,p3), h = (hi,hj, h3)T.
6. Find a new approximation
(Pi,Pj,P3) := (Pi + hi,'j + hj,'3 + h:i).
7. Verify the inequality \\F(pi,pj,p3)\\ < e, where e is a stopping criterion. If the inequality is not satisfied then move to Step 2.
The numerical solution of real-life problems showed that application of our algorithm reduces the number of arithmetic operations by a factor of 1.5 if we carry out calculations in the neighborhood of nominal conditions.
Acknowledgements. This work has been supported by the Russian Foundation for Basic Research, grant No. 15-01-03228-a.
References
1. Merenkov A.P. Hasilev V.Ya. Teoriya gidravlicheskikh tsepey [Theory of Hydraulic Circuits]. Moscow, Nauka, 1985. 278 p.
2. Shashkov O.K., Shashkov V.O. [Calculation of the Variable Mode for the Heat and Power Plants with Steam Turbines Using the Hydraulic Circuits Method]. Teploenergetika [Thermal Engineering], 2004, no. 4, pp. 67-71. (in Russian)
3. Levin A.A., Tairov A.E., Chistyakov V.F. [Load Flow Analysis in Power Supply Plants as in Hydraulic Circuits with Controlled Parameters]. Truboprovodnye sistemy energetiki: matematicheskoe modelirovanie i optimizaciya [Pipe System of Power Industry: Mathematical Modelling and Optimization], Novosibirsk, Nauka, 2010, pp. 115-124. (in Russian)
4. Peterson D.F., Saminskaya L.V., Khabenskiy V.B. [On the Subject of Hydraulic Calculation for the Straight Through Subscritical Boiler on the Electronic Digital Computer]. Trudy Tsentral'nogo kotloturbinnogo institute/, [Bulletin of the Central Boiler-and-Turbine Institute], 1969, issue 98, pp. 60-70. (in Russian)
5. Vul'man F.A., Khor'kov N.S. Teplovye raschety na ECVM teploehnergeticheskikh ustanovok [Thermal Calculation of Heat and Power Plants on the Electronic Digital Computer]. Moscow, Energiya, 1975. 200 p.
6. Ivanov V.A. Rezhimy moshchnykh paroturbinnykh ustanovok [Operation of Large Steam-Turbine Plants]. Moscow, Energiya, 1971. 280 p.
7. Raschet i proektirovanie pyleprigotovitel'nykh ustanovok kotel'nykh agregatov [Analysis and Engineering for Fuel Pulverizing Plants of Boiler Units]. Leningrad, Tsentral'nyiy kotloturbinnyiy institut, 1971. 310 p.
8. Krasnoselsky M.A., Vainikko G.M., Zabreiko P.P. et al. Priblizhennoe reshenie operatornykh uravneniy [Approximate Solution of Operator Equations]. Moscow, Nauka, 1669. 456 p.
Received May 18, 2016
УДК 519.711.3 DOI: 10.14529/mmpl60405
ОБ ИСПОЛЬЗОВАНИИ СТРУКТУРЫ СИСТЕМЫ НЕЛИНЕЙНЫХ УРАВНЕНИЙ, ОПИСЫВАЮЩИХ ГИДРАВЛИЧЕСКИЕ ЦЕПИ ЭНЕРГОУСТАНОВОК ПРИ ЧИСЛЕННЫХ РАСЧЕТАХ
A.A. Левин, В.Ф. Чистяков, Э.А. Таиров
В статье представлены результаты применения методов теории гидравлических цепей для интерактивного моделирования гидрогазодинамических и тепловых процессов в оборудовании тепловых электрических станций. Сформулирована постановка
задачи расчета потокораспределения в трактах энергоустановок с различными законами падения давления на ветвях в сложных газовоздушных и пароводяных трактах. Показано, что применение традиционных способов расчета потокораспределения методами теории гидравлических цепей затруднительно для задач. Исследованы некоторые аспекты разрешимости соответствующих систем нелинейных уравнений. Рассмотрены численные методы решения этих систем применительно к задачам, требующим выполнения расчетов в масштабе реального времени. Предложена расчетная схема, позволяющая свести исходную постановку задачи к классической схеме метода узловых давлений. Рассмотрен способ декомпозиции структуры гидравлической цепи на взаимосвязанные цепи меньшей размерности для уменьшения вычислительных затрат. Приведены результаты тестов, демонстрирующих высокую надежность разработанного способа расчета потокораспределения.
Ключевые слова: гидравлические цепи; математические модели; тракты, энергоустановок; системы нелинейных уравнений; метод Ньютона; метод прост,ой итерации.
Литература
1. Меренков, А.П. Теория гидравлических цепей / А.П. Меренков, В.Я. Хасилев. - М.: Наука, 1985. - 278 с.
2. Шашков, O.K. Расчет переменных режимов ТЭС с паротурбинными установками на основе метода теплогидравлических цепей / O.K. Шашков, В.О. Шашков // Теплоэнергетика. - 2004. - № 4. - С. 67-71.
3. Левин, A.A. Расчет потокораспределения в энергоустановках как гидравлических цепях с регулируемыми параметрами / A.A. Левин A.A., Э.А. Таиров, В.Ф. Чистяков // Трубопроводные системы энергетики: математическое моделирование и оптимизация. -Новосибирск: Наука, 2010. - С. 115-124.
4. Петерсон, Д.Ф. К разработке метода гидравлического расчета прямоточного котла до-критического давления на ЭЦВМ / Д.Ф. Петерсон, Л.В. Саминская, В.Б. Хабенский // Труды ЦКТИ. - 1969. - Вып. 98. - С. 60-70.
5. Пульман. Ф.А. Тепловые расчеты на ЭЦВМ теплоэнергетических установок / Ф.А. Вуль-ман, Н.С. Хорьков. — М.: Энергия, 1975. - 200 с.
6. Иванов, В.А. Режимы мощных паротурбинных установок / В.А. Иванов. - М.: Энергия, 1971. - 280 с.
7. Расчет и проектирование пылеприготовительных установок котельных агрегатов. - Л.: ЦКТИ, 1971. - 310 с.
8. Красносельский, М.А. Приближенное решение операторных уравнений / М.А. Красносельский, P.M. Вайникко, П.П. Забрейко, Я.Б. Рутицкий и др. - М.: Наука, Главная редакция физико-математической литературы, 1969. - 456 с.
Левин Анатолий Алексеевич, кандидат технических наук, ведущий научный сотрудник, Институт систем энергетики им. Л.А. Мелентьева СО РАН (г. Иркутск, Российская Федерация), [email protected].
Чистяков Виктор Филимонович, доктор физико-математических наук, главный научный сотрудник, Институт динамики систем и теории управления СО РАН (г. Иркутск, Российская Федерация), [email protected].
Таиров Эмир Асгадович, доктор технических наук, главный научный сотрудник, Институт систем энергетики им. Л.А. Мелентьева СО РАН (г. Иркутск, Российская Федерация), [email protected].
Поступила в редакцию 18 мая 2016 г.