Научная статья на тему 'Application of the method of insignificant divergencies to diagnose the technical aircraft gas turbine engine state under the transient-state conditions of its operation'

Application of the method of insignificant divergencies to diagnose the technical aircraft gas turbine engine state under the transient-state conditions of its operation Текст научной статьи по специальности «Электротехника, электронная техника, информационные технологии»

CC BY
11
2
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
transient-state conditions / diagnostics / analytical models / aircraft gas turbine engines / classification of conditions / переходные режимы / диагностирование / аналитические модели / авиационные газотурбинные двигатели / классификация состояний

Аннотация научной статьи по электротехнике, электронной технике, информационным технологиям, автор научной работы — Oleg F. Mashoshin, Ilya G. Kharmats

The article deals with issues related to the use of parametric information of the transient-state gas turbine engines (GTE) operation conditions for diagnosing their technical condition during the operation. A review of general approaches to computational algorithms for the recognition and classification of the condition applicable to aircraft GTE has been carried out. The significance of analytical models in modern algorithms for assessing the technical GTE condition is emphasized. The construction of a linearized mathematical model for the transient-state condition of the generalized-scheme aircraft GTE operation has been considered. It represents a system of equations analytically combining the relative parameter divergences measured during the engine operation with the relative divergences of unmeasured thermogasdynamic parameters and geometric gas-air flow duct parameters allowing for the technical condition of gas-air channel elements to be classified. A method for constructing mathematical and diagnostic engine models, using the transient response data, has been formulated. The capability of employing a method of insignificant divergences, used to build linear (linearized) mathematical and diagnostic GTE models for the steady-state conditions of its operation, has been demonstrated as well. It is shown that, despite the structural similarity of linear models of the steady and transient-state processes, diagnostics by means of the stated above processes is based on completely different principles – under the steady-state condition, the classification of a technical condition is determined by the variation in the value of the group of controlled responses, and under the transient-state condition, this operation is based on correlating the change in the transient-state behavior. To ensure the versatility of employing proposed methods regarding various GTE designs installed on modern civil aircraft, a generalized-design aircraft GTE model – a three-shaft bypass turbojet engine with mixing flows in a common jet nozzle, has been considered.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Применение метода малых отклонений для диагностирования технического состояния авиационного газотурбинного двигателя на переходных режимах его работы

В статье рассмотрены вопросы, связанные с использованием параметрической информации переходных режимов работы газотурбинных двигателей (ГТД) для диагностирования их технического состояния в процессе эксплуатации. Проведен обзор общих подходов к вычислительным алгоритмам распознавания и классификации состояний применительно к авиационным ГТД. Показано место аналитических моделей в современных алгоритмах оценки технического состояния авиационных ГТД. Рассмотрено построение линеаризованной математической модели переходного режима работы авиационного ГТД обобщенной схемы – системы уравнений, аналитически связывающих относительные отклонения параметров, измеряемых в процессе работы двигателя, с относительными отклонениями неизмеряемых термогазодинамических параметров и геометрических параметров газовоздушного тракта, позволяющих классифицировать техническое состояние элементов проточной части газотурбинного двигателя. Сформулирован метод построения математической и диагностической моделей двигателя с использованием характеристик переходного процесса, а также показана возможность применения метода малых отклонений, используемого для построения линейных (линеаризованных) математических и диагностических моделей ГТД для стационарных режимов его работы. Показано, что, несмотря на структурное сходство линейных моделей установившегося и переходного процессов, диагностирование с их помощью базируется на совершенно разных принципах – на установившемся режиме классификация технического состояния определяется по изменению величины группы контролируемых откликов, а на переходном режиме эта операция основывается на сопоставлении изменения характера протекания переходного процесса. Для обеспечения универсальности применения предложенных методов к различным схемам ГТД, устанавливаемых на современных самолетах гражданской авиации, рассмотрена модель обобщенной схемы авиационного газотурбинного двигателя – трехвального двухконтурного турбореактивного двигателя со смешением потоков в общем реактивном сопле.

Текст научной работы на тему «Application of the method of insignificant divergencies to diagnose the technical aircraft gas turbine engine state under the transient-state conditions of its operation»

Vol. 26, No. 05, 2023

Civil Aviation High Technologies

МАШИНО СТРОЕНИЕ 2.5.12 - Аэродинамика и процессы теплообмена летательных аппаратов; 2.5.13 - Проектирование, конструкция и производство летательных аппаратов;

2.5.14 - Прочность и тепловые режимы летательных аппаратов;

2.5.15 - Тепловые электроракетные двигатели и энергоустановки

летательных аппаратов; 2.5.16 - Динамика, баллистика, управление движением летательных аппаратов

УДК 681.3.069:681.3.015

DOI: 10.26467/2079-0619-2023-26-5-81-95

Application of the method of insignificant divergencies to diagnose the technical aircraft gas turbine engine state under the transient-state conditions of its operation

O.F. Mashoshin1, I.G. Kharmats1

Moscow State Technical University of Civil Aviation, Moscow, Russia

Abstract: The article deals with issues related to the use of parametric information of the transient-state gas turbine engines (GTE) operation conditions for diagnosing their technical condition during the operation. A review of general approaches to computational algorithms for the recognition and classification of the condition applicable to aircraft GTE has been carried out. The significance of analytical models in modern algorithms for assessing the technical GTE condition is emphasized. The construction of a linearized mathematical model for the transient-state condition of the generalized-scheme aircraft GTE operation has been considered. It represents a system of equations analytically combining the relative parameter divergences measured during the engine operation with the relative divergences of unmeasured thermogasdynamic parameters and geometric gas-air flow duct parameters allowing for the technical condition of gas-air channel elements to be classified. A method for constructing mathematical and diagnostic engine models, using the transient response data, has been formulated. The capability of employing a method of insignificant divergences, used to build linear (linearized) mathematical and diagnostic GTE models for the steady-state conditions of its operation, has been demonstrated as well. It is shown that, despite the structural similarity of linear models of the steady and transient-state processes, diagnostics by means of the stated above processes is based on completely different principles - under the steady-state condition, the classification of a technical condition is determined by the variation in the value of the group of controlled responses, and under the transient-state condition, this operation is based on correlating the change in the transient-state behavior. To ensure the versatility of employing proposed methods regarding various GTE designs installed on modern civil aircraft, a generalized-design aircraft GTE model - a three-shaft bypass turbojet engine with mixing flows in a common jet nozzle, has been considered.

Key words: transient-state conditions, diagnostics, analytical models, aircraft gas turbine engines, classification of conditions.

For citation: Mashoshin, O.F., Kharmats, I.G. (2023). Application of the method of insignificant divergences to diagnose the technical aircraft gas turbine engine state under the transient-state conditions of its operation. Civil Aviation High Technologies, vol. 26, no. 5, pp. 81-95. DOI: 10.26467/2079-0619-2023-26-5-81-95

Civil Aviation High Technologies

Vol. 26, No. 05, 2023

Применение метода малых отклонений для диагностирования технического состояния авиационного газотурбинного двигателя

на переходных режимах его работы

О.Ф. Машошин1, И.Г. Хармац1

1 Московский государственный технический университет гражданской авиации,

г. Москва, Россия

Аннотация: В статье рассмотрены вопросы, связанные с использованием параметрической информации переходных режимов работы газотурбинных двигателей (ГТД) для диагностирования их технического состояния в процессе эксплуатации. Проведен обзор общих подходов к вычислительным алгоритмам распознавания и классификации состояний применительно к авиационным ГТД. Показано место аналитических моделей в современных алгоритмах оценки технического состояния авиационных ГТД. Рассмотрено построение линеаризованной математической модели переходного режима работы авиационного ГТД обобщенной схемы - системы уравнений, аналитически связывающих относительные отклонения параметров, измеряемых в процессе работы двигателя, с относительными отклонениями неизмеряемых термогазодинамических параметров и геометрических параметров газовоздушного тракта, позволяющих классифицировать техническое состояние элементов проточной части газотурбинного двигателя. Сформулирован метод построения математической и диагностической моделей двигателя с использованием характеристик переходного процесса, а также показана возможность применения метода малых отклонений, используемого для построения линейных (линеаризованных) математических и диагностических моделей ГТД для стационарных режимов его работы. Показано, что, несмотря на структурное сходство линейных моделей установившегося и переходного процессов, диагностирование с их помощью базируется на совершенно разных принципах - на установившемся режиме классификация технического состояния определяется по изменению величины группы контролируемых откликов, а на переходном режиме эта операция основывается на сопоставлении изменения характера протекания переходного процесса. Для обеспечения универсальности применения предложенных методов к различным схемам ГТД, устанавливаемых на современных самолетах гражданской авиации, рассмотрена модель обобщенной схемы авиационного газотурбинного двигателя - трехвального двухконтурного турбореактивного двигателя со смешением потоков в общем реактивном сопле.

Ключевые слова: переходные режимы, диагностирование, аналитические модели, авиационные газотурбинные двигатели, классификация состояний.

Для цитирования: Машошин О.Ф., Хармац И.Г. Применение метода малых отклонений для диагностирования технического состояния авиационного газотурбинного двигателя на переходных режимах его работы // Научный Вестник МГТУ ГА. 2023. Т. 26, № 5. С. 81-95. DOI: 10.26467/2079-0619-2023-26-5-81-95

Introduction

With the increase in the operating time, the engine performance varies not merely under steady but also under the transient-state conditions of its operation. As the practice shows, a dynamic-response factor possesses high sensitivity to the varying of basic condition (factors) parameters of an object under study. Since up-to-date facilities to measure parameters allow for monitored parameters to be recorded multiple times within a short span of time and for the dynamics of their variation to be assessed, this makes it possible to utilize them to diagnose and predict the technical aircraft GTE condition.

For example, pressure, temperature, pressure-to-temperature ratio, flow velocity, oil and fuel consumption, flow areas of air-gas channel sections, thrust as well as the rotor speed are referred to the thermogasdynamic GTE parameters.

Time series of monitored parameters are used in various methods of the object (image) condition identification based on the mathematical object behavior description: analytical methods, statistical methods (including Bayesian hierarchical modeling [1, 2]), modeling methods based on the similarity (SBM) [3], machine learning methods (ML, Deep ML) [4-8], etc. Diagnosing the technical GTE condition, using the methods of mathematical modeling, lies in the solution to the problem of recognizing the EGT image

(condition) (classification), i.e., selecting a single variant of condition in the hierarchy (tree) of pre-determined possible GTE conditions. In many cases, a solution to the classification problem is supplemented by the calculation of the selected condition probability [9] and the limited parameter set, characterizing an ultimate limit state (for example, residual RUL life or synthetic divergence parameters from an object state classified as "normal").

In terms of calculations, the classification problem can be rationally solved based on:

• one-stage calculation - a result (the most probable object state) is defined during a single calculation cycle;

• two-stage calculation - the binary classification (the normal/abnormal state) is carried out in the first phase. In the second phase, the search for the most probable condition in the tree of abnormal states is executed provided that the result of the first calculation phase equals "the abnormal state";

• iterative calculation - the selection of a subset for variants of a solution, derived in the branch of the classification tree from the result of (i - 1)th calculation phase is the result of ith calculation phase. Iterations are completed while reaching the tree "leaf" or null decisions (for example, statistically imperceptible values for the probabilities of decisive subsets).

In modern GTE diagnosis algorithms, two-stage or iterative calculations are generally applied which allow for various mathematical models of an object, including analytical models, to be used at different stages of calculation.

The analytical GTE models are built on the functional dependencies of thermogasdynamic parameters derived from theoretical explorations [10, 11] and practical research of the thermodynamics and heat transfer processes. In the context of diagnostics terms, monitored parameters are attributed to signs, but parameters, classifying the GTE condition, are attributed to factors. Thus, thrust R, fuel consumption Gt, turbine inlet temperature T*r or turbine exhaust temperature T*T, parameters of the working substance along a duct as well as some parameters of fuel and oil systems can be referred to the engine condition signs. Burnouts of turbine air-

foils, flame tube of combustions, deformation of the flow section elements, etc. can serve as the examples of possible conditions.

The solution (classification of a state) is carried out by critical divergences of thermo-gasdynamic parameters. For example, the variation of turbine exhaust temperature is compared with a reference model. The reference model is built according to technical engine data. Temperature is monitored at a takeoff mode to which the reference turbine exhaust temperature complies with. In some cases, temperature ТТ as well as parameters TH and РН are used to calculate engine thrust, and it is compared with the thrust which should be generated under the conditions.

Specific capabilities are embedded in the diagnostic parameter "fuel consumption". The experience shows that the damage to the GTE airgas channel increases fuel consumption by 120...150 kg/h with the simultaneous variation of other thermodynamic parameters. Fuel consumption performance is well indicative of the technical condition of combustion chambers and turbine nozzles. However, to be based on the experience, the precise metering of consumption is hindered by errors of flow meters resulted from the necessity to take into consideration kerosine-based aviation fuel density at various temperatures.

Under certain conditions, GTE condition diagnostics can be also carried out by the pressure upstream the nozzles рт, but, in this case, measurement errors can be critical.

In order to ensure a reliable GTE condition classification by means of analytical models based on the results of gauged thermogasdynamic parameters, it is feasible to conduct the primary processing of time series. In addition to the standard processing (deleting bad-values, smoothing, recovery of missing values, etc.), parameter values should be brought to standard conditions and deterministic engine operation modes.

The results of research in the domain of GTE diagnostics by thermogasdynamic parameters allowed us to specify that the most sensitive and informational factor of the engine air-gas channel condition is adiabatic turbine efficiency цТ. Obviously, it is impractical to measure цТ. How-

ever, it can be expressed by rotor speed, extent of pressure ratio nK and turbine inlet temperature ТГ. This dependence will be empirical and specific referring to a particular engine.

A known dynamic analytical engine model is also proposed to diagnose. In this case, average integral values of monitored parameters over a definite span of time of load run-up and reset are used as a response.

Methods of the GTE performance calculation under the transient-state (not steady) conditions are considered in many theoretical explorations, applied research [12-15] and others. In this respect, investigations of the dynamic diagnostic GTE gas turbine model are given in [16].

Goal and objectives of research

The goal is the theoretical substantiation and investigation of capabilities for the use of ther-mogasdynamic parameters of the GTE operation recorded under some transient-state condition of its operation for the efficient technical diagnostics of major gas generator assemblies. Obtained models (besides compulsory requirements imposed towards similar models) must be simple, as far as feasible, to allow for a great array of monitored data to be efficiently processed for the purpose of obtaining practically significant results. The goals of research were as follows:

• an assessment of a principal capability to apply, for the transient-steady GTE condition, approaches and models used for the steady-state conditions;

• a construction of a linearized mathematical model of the transient GTE operation condition for the subsequent assessment of its technical condition using appropriate diagnostic models;

• an assessment of limitations and assumptions encompassing a practical field for the application of constructed models.

Methods and materials of research

While conducting the research, the methods of mathematical modeling of thermogasdynamic processes, linearization methods of second-kind

and greater model equations, methods of statistical processing of parameter arrays and the comparative assessment with modeling results.

As a starting point of research, a linear (linearized) mathematical model (MM), determining relationships among insignificant divergences of engine thermogasdynamic parameters, is under consideration. To be based on such a model, diagnostic GTE models are constructed, at the same time, a problem of determining variations of a series of parameters during the operation, which characterize its technical condition by the change of a specific group of its monitored parameters (factors), is solved. In order to calculate a diagnostic model of the given type, a method of parallel matrices, stated in [17], is employed. As a result, its solution is classified into the group of decisions of redundant subsystems of equations based on the same general MM. The derived diagnostic matrices are analyzed line by line on the limits of errors. Equations (lines) with high values of computed errors are excluded from the further analysis. The total calculation result represents an aggregate of parameter divergences according to which the technical engine condition is classified.

Insignificant divergences of argument parameters, involved in forming diagnostic models, are defined by comparing the results of their measurement at various time points using the same elements of measuring chains. Thus, the diagnosis accuracy by means of a diagnostic model is defined not only by an error of measurement but their frequency. The 1% error in the value of the relative divergence amounts barely a small portion of error percentage in the absolute parameter value.

The requirements concerning the scope of diagnosis outline the minimally required number of diagnostic sub-models (levels) under which an assigned degree of defect isolation will be reached. Notably, where appropriate, a list of structure-used sub-models of factors can be expanded not only by the calculation of an original linear MM, using additional equations, but also by means of a logical transfer to linear diagnostic models [18] not directly associated with the original linear MM.

It should be pointed out to the fact that the diagnostic models under consideration are poor-

Vol. 26, No. 05, 2023

Civil Aviation High Technologies

ly suited for the binary classification problems (normal/abnormal state) but can be applied in the two-stage or multi-stage diagnostic algorithms in phases of investigating individual tree branches of the condition classification.

Usage of the transient-state conditions models, obtained for the steady-state GTE operation conditions

While calculating the transient-state conditions, it is common practice to assume a quasi-steady state of the operation process, i.e., it is supposed that for the mathematical description of conditions, the most of relationships, used for the calculation of the steady-state (design) GTE operation condition, are applicable. Merely, the equations of power balance for the turbine and compressor located on the same shaft, in which, it is essential to take into consideration the timevariation of kinetic energy of spinning rotary masses, are subjected to significant variations.

One of the main engine performance figures is the power response which is determined by minimum time required for the transition from the idle mode for the maximum thrust mode, i.e., by the time of acceleration tPA3r from nMr to nB3^, as well as by the time of reverse transition from nB3^ to nMr - by the time of rotation speed reduction tCE.

Thus, let us assume that: • equations of continuity, consumption and its corollaries remain true for the transient response data;

balance of rotor power under the transient GTE operation condition is defined in conformity with the formula

nt - Nk.=de,

4m dt

(1)

where

E =

I a' 2

(2)

E is rotor kinetic energy; I is the moment of rotor inertia relatively the spinning axis; rn is an angle speed of rotor spinning.

Power can be defined according to the formula

N = LG = Ma.

(3)

Here L is operation; G is second consumption of working substance; M is torque. It follows from the equation (1):

d_

dt

f Ia2 1

v 2 у

= G

, - L - L

В T K

m

(4)

Here m = is the air-to-gas ratio. Suppo-

GK

sing that under the reset condition — LT 0 = LK

m0

it can be written in the form:

Gb

- T - T 1-G

t k i _ ^в

m )

lt lt0 i \lk lk0 )

m )

= gblk 0

—lt

m )

If to assume that under the short-term transient condition, the value m varies insignificantly, then

ô\ — LT | « 5LT . The latest expression will be as follows:

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

l m )

GB I — LT - LK I = GbLk 0 \SLT -5LK ]= Ia

m

da dt

As

d- = d(ю-ш0) = ю0dt(5со) = ш0(ôco) =a0(bn) ,

Civil Aviation High Technologies

Vol. 26, No. 05, 2023

then (5) will be provided as follows:

gblk о I ш2

(bLj. -5LK) = (dш) = (dn) = S(5LT -5LK ).

(6)

Here S is the constant defined at the datum point of the transient-state condition.

S =

f gblk ^

v iz ш

(7)

Jo

The value S characterizes the ratio between the compressor power and kinetic energy of rotor spinning under the original steady-state condition.

Rotor acceleration at an assigned rotation frequency is defined by the value of fuel supply excess AGT. In real processes, fuel supply excess is conducted considering a series of limitations.

Under the steady-state condition, maximum turbine inlet temperature Trmax is reached at the maximum rotation speed nmax. During the acceleration, under n < nmax, a short-term temperature overheat Tr above the maximum value by 50...1200 is allowed. When there are no other limitations, fuel supply can be increased by 1.5.2 times at acceleration modes. Temperature Tr increase during acceleration due to the effect

of heat throttling causes the parameter

д(Лв )

increase which is identical to the ratio

GB

incease. The given circumstance results in to the fact that at the same rotation speed, a point of combined turbine-compressor operation shifts while accelerating towards the boundary of its steady operation. As a rule, under small and average values n, possible excesses of fuel supply are limited by the compressor stability, and under large values n are limited by temperture Tr. The rotation speed n and, in particular, thrust P vary insignificantly at the beginning of acceleration and increase abruptly at the end.

The rotor speed is reset by reducing the fuel supply. In this case, the turbine power becomes less than the power consumed by the compressor, and the rotor receives negative acceleration. The main factor limiting the reduction in fuel

supply when the rotation speed is reset is the limit of stable combustion chamber operation.

The specific features of the unsteady processes in a two-shaft GTE are associated with the fact that when the power balance on the tur-bocompressor shaft is disturbed, the high- and low-pressure rotor speed generally varies at different speeds, and the slip value C = n2 / n}, during the transient-state process, may differ substantially from its value both in the initial and final modes.

Thus, in a two-shaft turbojet engine (TJE), as the rotation speed decreases, the rotor slip increases, which is caused by variations of the angles of attack in the compressor spools, as well as by the nature of the redistribution of pressure differential between compressor spools. The same patterns of gas-dynamic relationships develop in these engines under the transient-state conditions.

However, during the acceleration process, the rotor slip increases, and during the reset process, it decreases compared to the steady-state conditions. This difference makes alterations to the nature of the position of the joint turbine-compressor operation line on the low-pressure compressor (LPC) characteristic.

During acceleration, the rotation speed n2 increases faster than n}. This reduces LPC throttling; the line of joint operating modes shifts to the right into the area of excessive air flow rates. When the rotation speed is reset, the increase in slip slows down and the high-pressure compressor (HPC) begins to throttle the LPC, which leads to a shift in the line of joint operating modes towards the boundary of stable operation.

The appearance of the curves, characterizing the joint turbine-compressor operation on the HPC characteristic, is similar to the appearance of these curves for a single-shaft TJE.

Two-shaft bypass turbojet engines retain the same features of the transient-state process as two-shaft TJE. The only difference is that the divergence of the line of operating modes on the fan characteristic from the line of the steady-

state conditions turns out to be less than that of a TJE; whereby, the greater the bypass ratio, the smaller the divergence. This issue is explained by the fact that the air flow through the fan is determined not only by the air flow through the HPC, but also by the flow through the external circuit, which insignificantly depends on the rotor sliding.

For the transient-state conditions caused by a slight variation of parameters compared to their values under the original steady-state condition, it is permissible to solve the problem under consideration in a linear formulation.

Assume that the power balance on the rotor shaft is disrupted as a result of some finite increment in fuel consumption SGT. This will cause an increase in n, Tr, ttK and other parameters, which will lead to a variation of LT and LK in accordance with the equations for the turbine-compressor operation in the initial mode [19] (hereinafter, for simplicity of calculations, formulas for a single-shaft engine are used):

SLT = 8Tr + K3STTt + 8r]T, (8)

5Lk = Kn5n + Kl5^k - 5-Hk, (9)

where K are the corresponding base coefficients. Taking into account (8) and (9), equation (6) will take the following form1:

(Sri)' = S(STr - K38XT - K„Sn - KlSttk). (10)

In the unsteady process under consideration, we also neglect the variation of gas mass in the engine air-gas channel volume (SGB = SGr = const). The compressor characteristics are assumed unchanged.

The continuity equations after simple transformations can be represented as follows:

(1 - Kw)8ttk - Km5n = 0.5Т, (11)

(1 - 0.5 K3K4 + Кб)5жт = Кбджк. (12)

Since power imbalance is caused by the variation of fuel consumption, we use an equation that relates SGT to other parameter variations. In this case, as already mentioned above, we neglect the combustion efficiency т]Г variation, the lag in heat generation with a sharp increase in fuel consumption, as well as the heat accumulation by the air-gas channel parts. According to the equation of energy conservation in the combustion chamber, let us write down:

SGt = SGb + К5дтг - (К - 1)STK. (13)

As equations that close the system, we use the expressions for air flow and temperature at the compressor outlet:

SGb = KmSn + KwSTTk, (14)

Tk = K2KnSn + KKSk. (15)

The system of equations (11)-(15) allows us to determine all the main engine parameters through Sn and SGT.

So, for example, it can be shown that:

Sttk = Kx SGt + Ky Sn, (16)

STr = Ka SGt - Kb Sn, (17)

Sttt = Кб Kz (Kx SGt + Ky Sn), (18)

where

Kx =

2(1 -К10)k5 -(Ks - 1)К2Kl + Кю'

Ky = ((К - ))К 2 К + Km(K -1) = KX Г(К5 -1) K2Kn +Km (2 К5 -1)) .

Y 2 (1 - Кю )К5 -(К5 - 1)К2Kl + Кю X ^ 5 У 2 " mV 5

1 We assume that for the transient-state condition дцт = 5цк = 0.

1

Civil Aviation High Technologies

Vol. 26, No. 05, 2023

Ka = — + -

( - 1)K2Kl + Ko

= 2Kx(1 - К ),

К5 2(1 -Kw)K5 -(K5 - 1)K2Kl + K1(

Kb = K- [((1 - K5 ) К2Kn + Kio ) KY - (K -1) К2Kn + Km ] = 2( - Kx (1 - Kio )),

1

Kz =

1 - 0.5K3 K4 + K6

Let us introduce an expression that establishes the relationship between the specified parameters (for example, using a fuel supply controller) in the form of a linear function:

(¿n)' + aSn = bSGT, (19)

where a and b are constant coefficients: a = M [Kb - Ky (KKKz - Kl) + K„], b = M [Ka + Kx (KKKz - Kl)].

The general solution of the linear differential equation (19) takes the form:

5n = e~at [Jb5GTeatdt + const . (20)

If the increment in fuel supply SGT (compared to fuel consumption in the original mode) is instantaneous and remains constant at the time of the entire transient-state process, then from (20) it follows that:

ôn = b SGT + const • e a

(21)

The integration constant is determined from the initial conditions. At t = 0 we have 5n = 0. Accordingly

b

const = —ôG

a

and

ôn = b(1 - e~at )ôGT. av '

(22)

Formula (19) determines the law of change in rotation speed in duration of the transient-state condition caused by the fuel consumption variation by the value SGT. By specifying different values of time t, you can determine the corresponding value Sn, and then, substituting it into the resulting expressions (11)-(15), find the value of the required parameters at any moment of the transient-state process. From (22) it follows that at t ^ œ:

ôn = b ôGT . a

(23)

This rotation speed value n characterizes the new final engine condition, which will be established as a result of the fuel consumption variation.

In the case when the initial disturbance is not instantaneous, but changes under the transient-state condition depending on any thermogasdy-namic parameter or time variation, coefficient b in equation (19) will be a variable value and the solution will take a more complex form.

Solutions can be obtained similarly if the imbalance of power on the shaft is caused by other factors, for example, a rapid change in the jet nozzle area 5FC or in the combustion mode in the afterburner.

To determine the behavior of a parameter variation in a two-shaft GTE during the transient-state condition caused by an imbalance of power on each compressor shaft, one should use the system of equations written for a two-shaft GTE and the general principles for calculating unsteady modes outlined above.

After a series of transformations similar to those carried out when composing a mathemati-

at

Vol. 26, No. 05, 2023

Civil Aviation High Technologies

Fig. 1. Schematic diagram of a three-shaft gas turbine engine

cal model of a single-shaft engine, instead of (10), we obtain a system of two ordinary linear inhomogeneous differential equations of the first kind of the form:

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

r

(n) + PM + 015n2 = ^5Gt ;

, (24) (5n2) + P28n1 + Q2 5n2 = N25GT ,

where P, Q, N are numerical coefficients.

An example of the construction and calculation of a linear mathematical model of the transient-state condition for a two-shaft PS-90A-type GTE is considered in detail in [19].

As a result of integrating system (24), the dependences of the turbocompressor rotation speed on the time and the change in fuel supply are determined:

§n1 = fi (t, §gt ); 5n2 = /2 (t, §gt ).

Using the linearized equations of the two-shaft GTE processes as obtained for its transient-state condition of its operation, from the values 5n1 and Sn2, if applicable, as before for a singleshaft one, to determine the laws of variation of all other necessary characteristics.

It should be noted that, within the meaning of the presented solution, the end effects retain the basic properties of linear equations, i.e., the proportional variation of the increments of all parameters at any moment of the transient-state

process depending on the magnitude of the disturbance (for example, from SGT). The joint influence of several disturbances is determined by the summation of partial parameter increments.

Modeling of the transient-steady processes of a generalized schematic diagram engine

We will consider the transient-state process that is caused by an insignificant parameter variation compared to their value in some original mode. To ensure the versatility of the approach to constructing a mathematical model, we will consider a generalized aircraft GTE schematic diagram. A mathematical model (MM) has been compiled for a three-shaft by-pass TJE with flow mixing in a common jet nozzle (fig. 1), since most engines in use in civil aviation are its particular variants and can be obtained by excluding some elements from it.

The construction of such a MM is practical on the basis of the universal linear GTE MM given in [20]. To do this, adjustments should be made to it related to the dynamic nature of the processes under consideration.

First of all, as shown earlier, it is necessary to take into account the power imbalance on the shafts of low, medium and high pressure turbo-compressors by introducing additive components S1(dn1)', S(dn)', and S2(dn2)'on the right side of the corresponding power balance equations.

Civil Aviation High Technologies

Vol. 26, No. 05, 2023

Here

( ) d (bn )

(bn )'= ^,

(bn2 ) = Ä dt

(25)

(26) (27)

are accelerations of the rotational motion of the low, medium, and high-pressure rotor, respectively, and the constant

W ^ I

П ) In

S =

S 2 =

30 n

1"1 0

GK10 LK10

30

n

In2

GK 0 LK 0 I n 2

2 20

(28)

(29)

(30)

characterize the relationship between power and kinetic energy of rotation of the corresponding rotor.

The expressions for Si, S and S2 comprise the inertia moments Ii, I and I2 of the low, medium, and high pressure turbocompressors.

Further neglect the working substance mass variation, combustion efficiency, lag in heat generation and its accumulation by engine parts (we assume that the engine is warmed up). Otherwise, the system of equations developed for the steady-state condition remains unchanged.

Some information about the original system of equations is given in [15, 20, 21].

This system of equations establishes a connection between parameters characterizing the engine (factors) condition, the condition (responses) characters and operating parameters.

The MM comprises 53 linear equations and 91 variables. As independent MM variables we take all the areas of the flow sections (FTB, FTC, Fth, Fci, Fc2, Fc), tightness coefficients (qi, q2,

qox, q£), efficiencies of the main components (1b, In, lk1, 1k, jjr, JJtb, Vtc, 17th, 1m1, 1m, 1M2), as well as total pressure conservation coefficients (obx, okc, &tb, ore, &th, 0c1, 0c2, 02, oCM, oC). Generally, the subscripts correspond to the sections shown in Figure 1.

Furthermore, there are independent operating parameters: rotor speeds n1, n, n2, fuel consumption GT (or time t), flight Mach number, as well as environmental parameters: PH and TH. Since there are m = 38 independent variables, the number of unknowns is equal to the number of equations (n = 53), and, accordingly, the above system is defined.

It should be noted that the stated universal linear MM comprises the most common list of thermogasdynamic parameters and standard formulas for bringing the parameters for the standard atmospheric conditions. If applicable to introduce additional parameters, control laws, formulas for given parameters, etc. for a specific type of gas turbine engine, additional equations may be included into the linear MM.

The main operations for calculating and analyzing diagnostic parallel matrices for the steady-state GTE operation conditions are described in [20, 22, 23].

Research results

Computational studies show that methods for calculating and analyzing diagnostic parallel matrices for the steady-state conditions are also suitable for calculating diagnostic models of the transient GTE operation conditions. At the same time, despite the external MM similarity of the steady and transient-state processes, diagnosing with their help is based on completely different principles. Under the steady-state condition, variations of the technical condition are determined by the variation of the value of the group of controlled responses, while under the transient-state condition, this operation is based on a comparison of changes in the transient behavior.

It should be noted that this circumstance allows the joint use of both models, which in this case, complement each other, thus, significantly expanding the diagnostic capabilities.

2

2

So, for example, after determining Зщ, 5?jK, 8т]ТВ, Srjm, 8т]Г and 5q under the steady-state condition, you can take their values as a first approximation, then substituting them into the MM of the transient-state process. This makes it possible to use it to calculate changes in a group of other factors, for example, 5<jkc, 5FTB and 5FTH. The first factor, along with 8т]Г, characterizes the combustion chamber condition, its operation efficiency, and records the appearance of burnout, while an increase in 5FTB or 5FTH may indicate burnout of one of the turbine nozzles. By substituting the obtained values of 8окс, 3FTB and 5FTH into the diagnostic model of the steady process, it is possible to determine Зщ, Stjk, 8т]ТВ, Sr]m, 8т]Г и 5q in the second approximation. Iterations should be repeated until the difference in the results of the ith and (i+1)st approximations satisfies the specified requirements for calculation accuracy [13].

If you do not expand the range of diagnosed parameters, you can increase the accuracy of their determination by increasing the number of independent calculation formulas.

To perform a comparative analysis of changing controlled parameters of the transient-state process, you can use their average integral values for a selected control period of time. From the resulting mathematical model, it is not difficult to obtain the corresponding calculation formulas.

A higher quality processing of the results of calculating diagnostic matrices for the transient-state operation conditions can be carried out taking into account the assumption of equal probabilities of damage to the GTE gas-air channel. Accordingly, when constructing adequate multisystem diagnostic models, there is a need to take into account the probabilities of occurrence of various gas-air channel malfunctions and their possible combinations in relation to a specific GTE type.

To take into account the operational GTE features when assessing its technical condition, the method of constructing multi-system diagnostic GTE models is employed, using elements of centralist hierarchical structures [20]. The essence of the method is to localize GTE defects

through sequential optimization of a set of factors used in the linear diagnostic models. The process of selecting factors (optimization) is reduced to constructing several levels of the linear diagnostic models, at each of which the following procedures are performed:

1) calculation of the linear diagnostic model of the ith level;

2) analysis of the reliability of the results of calculating factor divergences using the ith level model, analysis of the significance of the resulting deviations for each factor;

3) checking whether the boundary conditions have been achieved (if the result is positive, the optimization process is completed);

4) modification of the current set of factors

by:

• excluding factors with zero (insignificant) deviations from the model;

• determining the "leading" factors (factors with maximum divergences) and introducing additional factors into the model that characterize the same node as the "leading" factors;

• construction of a new set of factors, alternative to that calculated in the ith model;

5) construction of a linear diagnostic model of the (i + 1)th level.

The significance of the jth factor of the current model is determined both by its absolute value (absolute significance) and by the relatively obtained divergences of other factors (relative significance). When analyzing the resulting divergences, in most cases the relative importance of the factor should be considered priority.

Discussion of results

The completed studies provide a theoretical basis for the development and substantiation of practical methods and algorithms for diagnosing the technical condition of the GTE gas generator elements based on the values of thermodynamic parameters (responses) recorded directly during the GTE operation. The use of such methods will increase the validity of assessing the technical condition (if the methods are employed in con-

junction with other methods for diagnosing the condition of the flow part), as well as increase the efficiency of obtaining an assessment of the aircraft GTE condition under operational conditions.

Conclusion

Modern algorithms for assessing the technical GTE condition and early detection of faults use various methods and models of dependencies between the measured GTE parameters and the conditions classified by the algorithms. In order to achieve a reliable classification of the GTE condition by analytical models based on the results of measured thermogasdynamic parameters, it is necessary to carry out primary processing of time series. After standard processing, the values of the series are reduced to standard atmosphere conditions and deterministic engine operating modes.

Deterministic models for diagnosing aircraft GTE can be expressed through a system of equations for the engine condition, the solution of which makes it possible to assign the technical GTE condition under study to one or more classification elements in the state tree. Linear (linearized) diagnostic equations are a finite set of expressions constructed for the increment of air flow, turbine inlet temperature, specific flow and other thermogasdynamic parameters. The right side of these equations contains the parameter divergences, which are determined by comparing the current values with the reference values.

The number of diagnostic equations is determined by the classes of possible GTE conditions, as well as the nomenclature, frequency, and error of the measured parameters. Recently, for the GTE diagnosis, it has been proposed to use complex parameters that, in an analytical form, connect several parameters with each other and, thereby, most fully characterize the operating procedure occurring in the engine.

It should be noted that the universal linear MM specified in the article contains the most common list of thermogasdynamic parameters and standard formulas for reducing the parameters to the standard atmospheric conditions. If it is

necessary to introduce additional parameters, control laws, formulas for reducing the parameters, etc. for a specific type of aircraft GTE, additional equations may be included in the linear MM.

As a result of the research, a method for constructing mathematical and diagnostic models of an aircraft engine, using transient response data, was proposed, and the capability of employing the methods used to construct linear mathematical and diagnostic models for the steadystate operation condition was shown. As a result, the model makes it possible to extract and interpret diagnostic information from data series obtained during the process of loading and unloading, as well as to expand the amount of information obtained about the technical aircraft GTE condition.

References

1. Zaidan, M.A., Mills, A.R., Harrison, R.F., Fleming, P.J. (2016). Gas turbine engines prognostics using bayesian hierarchical models: A variational approach. Mechanical Systems and Signal Processing, vol. 70-71, pp. 120-140. DOI: 10.1016/j.ymssp.2015.09.014

2. Zaidan, M.A., Relan, R., Mills, A.R., Harrison, R.F. (2015). Prognostics of gas turbine engine: An integrated approach. Expert Systems with Applications, vol. 42, issue 22, pp. 8472-8483. DOI: 10.1016/j.eswa.2015.07.003

3. Marins, M.A., Ribeiro, F.M.L., Netto, S.L., Da Silva, E.A.B. (2018). Improved similarity-based modeling for the classification of rotating-machine failures. Journal of the Franklin Institute, vol. 355, issue 4, pp. 1913-1930. DOI: 10.1016/j.jfranklin.2017.07.038

4. Vaezipour, A., Mosavi, A., Seigerroth, U. (2013). Machine learning integrated optimization for decision-making. In: 26th European Conference on Operational Research, Rome. Available at: https://www.semanticscho lar.org/paper/Machine-learning-integrated-optimi zation-for-making-Vaezipour-Mosavi/c1ad5937 eecf961a3be64d18889d5ffeb888de50 (accessed: 27.12.2022).

5. Sina Tayarani-Bathaie, S., Khorasa-ni, K. (2015). Fault detection and isolation of gas turbine engines using a bank of neural networks.

Journal of Process Control, vol. 36, pp. 22-41. DOI: 10.1016/jjprocont.2015.08.007

6. Amozegar, M., Khorasani, K. (2016). An ensemble of dynamic neural network identifiers for fault detection and isolation of gas turbine engines. Neural Networks, vol. 76, pp. 106-121. DOI: 10.1016/j.neunet.2016.01.003

7. Zhernakov, S.V. (2006). Application of neural network technology for the diagnostics of the technical condition of aircraft engines. Intel-lektualnyye sistemy v proizvodstve, no. 2 (8), pp. 70-83. (in Russian)

8. Lu, F., Jiang, J., Huang, J., Qiu, X. (2017). Dual reduced kernel extreme learning machine for aero-engine fault diagnosis. Aerospace Science and Technology, vol. 71, pp. 742-750. DOI: 10.1016/j.ast.2017.10.024

9. Mashoshin, O.F. (2015). Evaluation of diagnostic information in solving the task of aircraft operation. Nauchnyy Vestnik MGTU GA, no. 219 (9), pp. 53-56. (in Russian)

10. Kotlyar, I.V., Gitelman, A.I. Yermol-chik, V.N. et al. (1973). Transient processes in gas turbine power plants, in Kotlyar I.V. (Ed.). Leningrad: Mashinostroyeniye, 256 p. (in Russian)

11. Kazandzhan, P.K., Alekseyev, L.P., Govorov, A.N., Konovalov, N.Ye., Necha-yev, Yu.N., Pavlenko, V.F., Fedorov, R.M. (1955). Theory of jet engines. Moscow: Voyen-izdat, 296 p. (in Russian)

12. Ivanova, E.V., Tsymbler, M.L. (2020). Overview of modern time series management systems. Bulletin of the South Ural State University. Series: Computational Mathematics and Software Engineering, vol. 9, no. 4, pp. 79-97. DOI: 10.14529/cmse200406 (in Russian)

13. Cherkez, A.Ya. (1975). Engineering calculations of gas turbine engines using the small deviation method. 3rd ed., pererab. i dop. Moscow: Mashinostroyeniye, 380 p. (in Russian)

14. Akhmedzyanov, D.A. (2006). Nonstable regimes of aviation GTE. Vestnik Ufimskogo Gosudarstvennogo Aviatsionnogo Tekhnicheskogo Universiteta, vol. 7, no. 1, pp. 36-46. (in Russian)

15. Dmitriyev, S.A. (1996). Diagnosing the flow part of a gas turbine engine in steady and transient modes of its operation: Abstract of

Technical Sc. Dissertation. Kyiv: KIIGA, 28 p. (in Russian)

16. Tsoutsanis, E., Meskin, N. (2017). Derivative-driven window-based regression method for gas turbine performance prognostics. Energy, vol. 128, pp. 302-311. DOI: 10.1016/j.energy. 2017.04.006

17. Urban, L.A. (1975). Parameter selection for multiple fault diagnostics of gas turbine engines. Journal of Engineering for Gas Turbines and Power, vol. 97, no. 2, pp. 225-230. DOI: 10.1115/1.3445969

18. Kharmats, I.G. (2004). On the issue of technical condition diagnosing of aircraft gas turbine engines using linear mathematical models. Nauchnyy Vestnik MGTU GA, no. 75, pp. 87-92. (in Russian)

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

19. Grin, V.T. et al. (Eds.). (1982). Mathematical modeling of non-stationary processes in power plants with a gas turbine engine and other aerodynamic devices: proceedings. Moscow: TsIAM, 128 p. (in Russian)

20. Mikhnenkov, L.V., Kharmats, I.G. (2017). Parametric diagnosis of jet engines during operation using mathematical modeling: Monography. Moscow: Izdatelskiy Dom Acade-mii Zhukovskogo, 112 p. (in Russian)

21. Runacres, T., Hong, G. (1995). Improving gas turbine engine condition assessment using a thermodinamic model. In: 2nd Pacific International Conference Aerospace Science and Technology: 6th Australian Aeronautic Conference, Melbourne, 20-23 March, pp. 41-46.

22. Akhmedzyanov, A.M., Dubrav-skiy, N.G., Tunakov, A.P. (1983). Diagnosis of the jet engines conditions by thermos gas dynamic parameters. Moscow: Mashinostroyeniye, 206 p. (in Russian)

23. Mashoshin, O.F. (2005). Prediction of the vibration state of aircraft engines from the standpoint of classification problems. Nauchnyy Vestnik MGTU GA, no. 85, pp. 39-45. (in Russian)

Список литературы

1. Zaidan M.A. Gas turbine engines prognostics using bayesian hierarchical models:

A variational approach / M.A. Zaidan, A.R. Mills, R.F. Harrison, P.J. Fleming // Mechanical Systems and Signal Processing. 2016. Vol. 70-71. Pp. 120-140. DOI: 10.1016/j.ymssp.2015.09.014

2. Zaidan M.A. Prognostics of gas turbine engine: An integrated approach / M.A. Zaidan, R. Relan, A.R. Mills, R.F. Harrison // Expert Systems with Applications. 2015. Vol. 42, iss. 22. Pp. 8472-8483. DOI: 10.1016/j.eswa. 2015.07.003

3. Marins M.A. Improved similarity-based modeling for the classification of rotating-machine failures / M.A. Marins, F.M.L. Ribeiro, S.L. Netto, E.A.B. Da Silva // Journal of the Franklin Institute. 2018. Vol. 355, iss. 4. Pp. 1913-1930. DOI: 10.1016/j.jfranklin.2017. 07.038

4. Vaezipour A., Mosavi A., Seigerroth U.

Machine learning integrated optimization for decision-making [Электронный ресурс] // 26th European Conference on Operational Research. Rome, 2013. URL: https://www.semanticscho lar.org/paper/Machine-learning-integrated-optimi-zation-for-making-Vaezipour-Mosavi/c1ad5937 eecf961a3be64d18889d5ffeb888de50 (дата обращения: 27.12.2022).

5. Sina Tayarani-Bathaie S., Khora-sani K. Fault detection and isolation of gas turbine engines using a bank of neural networks // Journal of Process Control. 2015. Vol. 36. Pp. 22-41. DOI: 10.1016/j.jprocont.2015.08.007

6. Amozegar M., Khorasani K. An ensemble of dynamic neural network identifiers for fault detection and isolation of gas turbine engines // Neural Networks. 2016. Vol. 76. Pp. 106-121. DOI: 10.1016/j .neunet.2016.01.003

7. Жернаков С.В. Применение технологии нейронных сетей для диагностики технического состояния авиационных двигателей // Интеллектуальные системы в производстве. 2006. № 2 (8). С. 70-83.

8. Lu F. Dual reduced kernel extreme learning machine for aero-engine fault diagnosis / F. Lu, J. Jiang, J. Huang, X. Qiu // Aerospace Science and Technology. 2017. Vol. 71. Pp. 742-750. DOI: 10.1016/j.ast.2017.10.024

9. Машошин О.Ф. Оценка диагностической ценности информации при решении задач в области эксплуатации авиационной

техники // Научный Вестник МГТУ ГА. 2015. № 219 (9). С. 53-56.

10. Котляр И.В., Гительман А.И., Ер-мольчик В.Н. и др. Переходные процессы в газотурбинных установках / Под ред. И.В. Котляра. Л.: Машиностроение, 1973. 256 с.

11. Казанджан П.К. Теория реактивных двигателей / П.К. Казанджан, Л.П. Алексеев, АН. Говоров, НЕ. Коновалов, ЮН. Нечаев, В.Ф. Павленко, Р.М. Федоров. М.: Воениздат, 1955. 296 с.

12. Иванова Е.В., Цымблер М.Л. Обзор современных систем обработки временных рядов // Вестник ЮУрГУ. Серия: Вычислительная математика и информатика. 2020. Т. 9, № 4. С. 79-97. DOI: 10.14529/cmse200406

13. Черкез А.Я. Инженерные расчеты газотурбинных двигателей методом малых отклонений. 3-е изд., перераб. и доп. М.: Машиностроение, 1975. 380 с.

14. Ахмедзянов Д.А. Неустановившиеся режимы работы авиационных ГТД // Вестник УГАТУ. 2006. Т. 7, № 1. С. 36-46.

15. Дмитриев С.А. Диагностирование проточной части ГТД на установившихся и переходных режимах его работы: автореф. дисс. ... докт. тех. наук. Киев: КИИГА, 1996. 28 с.

16. Tsoutsanis E., Meskin N. Derivative-driven window-based regression method for gas turbine performance prognostics // Energy. 2017. Vol. 128. Pp. 302-311. DOI: 10.1016/j.ener gy.2017.04.006

17. Urban L.A. Parameter selection for multiple fault diagnostics of gas turbine engines // Journal of Engineering for Gas Turbines and Power. 1975. Vol. 97, no. 2. Pp. 225-230. DOI: 10.1115/1.3445969

18. Хармац И.Г. К вопросу о диагностировании технического состояния авиационных ГТД с использованием линейных математических моделей // Научный Вестник МГТУ ГА. 2004. № 75. С. 87-92.

19. Математическое моделирование нестационарных процессов в силовых установках с газотурбинным двигателем и в других аэродинамических устройствах: сборник статей / Под ред. В.Т. Гриня и др. М.: ЦИАМ, 1982. 128 с.

20. Михненков Л.В., Хармац И.Г. Параметрическое диагностирование ВРД в процессе эксплуатации с использованием математического моделирования: монография. М.: ИД Академии Жуковского, 2017. 112 с.

21. Runacres T., Hong G. Improving gas turbine engine condition assessment using a thermodinamic model // 2nd Pacific International Conference Aerospace Science and Technology: 6th Australian Aeronautic Conference. Melbourne, 20-23 March 1995. Pp. 41-46.

22. Ахмедзянов А.М., Дубравский Н.Г., Тунаков А.П. Диагностика состояния ВРД по термогазодинамическим параметрам. М.: Машиностроение, 1983. 206 с.

23. Машошин О.Ф. Прогнозирование вибросостояния авиадвигателей с позиций классификационных задач // Научный Вестник МГТУ ГА. 2005. № 85. С. 39-45.

Information about the authors

Oleg F. Mashoshin, Doctor of Technical Sciences, Professor, the Head of the Aircraft Engines Chair, Moscow State Technical University of Civil Aviation, o.mashoshin@mstuca.aero.

Ilya G. Kharmats, Candidate of Technical Sciences, Associate Professor, Associate Professor of the Technical Mechanics and Engineering Graphics Chair, Moscow State Technical University of Civil Aviation, kharmats@mail.ru.

Сведения об авторах

Машошин Олег Федорович, доктор технических наук, профессор, заведующий кафедрой двигателей летательных аппаратов МГТУ ГА, o.mashoshin@mstuca.aero.

Хармац Илья Григорьевич, кандидат технических наук, доцент, доцент кафедры технической механики и инженерной графики МГТУ ГА, kharmats@mail.ru.

Поступила в редакцию 12.02.2023 Received 12.02.2023

Одобрена после рецензирования 23.04.2023 Approved after reviewing 23.04.2023

Принята в печать 21.09.2023 Accepted for publication 21.09.2023

i Надоели баннеры? Вы всегда можете отключить рекламу.