Gas Turbine Simulation Taking into Account Dynamics of Gas Capacities

Abstract

The chapter considers one of the main dynamic factors of the turbine engine—the dynamics of gas capacities. Typically, the most influencing capacities in the turbine engine are combustion chamber, afterburner, mixing chamber, secondary duct of turbofan, and jet nozzle. Simulation of high-frequency transients in turbine engines needs taking into account this factor. For the needs of automatic control and parametric diagnostics, the equations of capacities must be combined with the equations of rotor dynamics and, sometimes, with the equations of a measurement system and actuators. The model complexity consists in two features. The first feature is in how many segments are used to simulate the capacity. The second feature is in what of three basic laws are taken into account at the gas motion description: the mass conservation law, the energy conservation law, and the momentum conservation law. This chapter includes the analysis of models of different complexity followed by the conclusions about their applicability. In the last part of the chapter, the real case of the engine dynamics analysis is considered when the designer does not need the simulation of the capacities’ dynamics in time, but needs estimating of the capacities’ ability to oscillate and in their natural oscillation frequencies.

1. Introduction

Engine development is known to include numerous stages and, among them, control systems and engine health management systems development. The development of these systems, however, includes conducting much experimental work, which is not a good choice, keeping in mind the cost of tests and time expenses. An alternative choice to be made is the involvement of the mathematical modeling into a development process [1]. One of the topics to be discussed in this chapter deals with the problem of considerable pneumatic volumes of a gas turbine: main and afterburning combustion chambers, bypass, exhaust nozzle, transition ducts, etc. Stationary gas turbines have extra volumes that must be considered when building up the mathematical model: intake with an air purification system and stack with a noise suppression device. In the models of a gas turbine-driven natural gas pumping compressor, the volumes of interest are input and output manifolds as well as a main pipeline. In all abovementioned cases, the designers must allow for the dynamics of pneumatic volumes when studying the transient behavior of an isolated engine or the engine as a component of a power plant.

The design process engages numerous kinds of models depending on the stage the development is at and the problems to be resolved by the particular model. But the “mother” of each and every model is a nonlinear component level model (thermodynamic model), which describes gas path variables using thermodynamic relations and performance maps of all main engine components, such as compressors, turbines, combustors, and input and output devices. When looking carefully at the transients described by this kind of models, one can conclude that the main factors affecting the engine transients are the inertia of rotors, the thermal inertia of engine parts, and the inertia of pneumatic volumes. Usually, designers are good to go with the model that eliminates the last two factors and considers only the rotor dynamics. Nevertheless, for some cases, it is good to have a model that is able to carefully simulate the processes, taking place in the pneumatic volumes, and the phenomenon of thermal inertia. Many researchers have paid their attention to the problem of pneumatic volumes within a total thermodynamic model. Next you will find a brief overview of their findings.

Almost five decades ago, Fawke and Saravanamuttoo proposed a method to simulate gas turbine inter-component volumes within the thermodynamic model [2, 3]. The proposed method found its niche in the field, but it has an increasing degree of differential equations and is too complicated for real-life calculations. The equations that describe the volumes in this method are also known for their low robustness. To cope with these challenges, the designer that is up to use this must take many assumptions to cope with these challenges. In this case, the equations can be simplified, which makes the solution more robust. But, unfortunately, the above studies do not provide any recommendations about a proper algorithm on how to compile the list of requirements for a particular engine. It is not even known whether the dynamics of all inter-component volumes should be simulated.

A few years later, a pretty complete analysis of the general problem of pneumatic volume simulation was made by Glikman in [4]. He described many methods to simulate the volume effect, but, unfortunately, did not pay enough attention to their comparison and highlighting the use cases. Moreover, the book does not cover the specific features of gas turbine engine simulation.

For the past two decades, many scientists have turned their sight to the problem of pneumatic volumes and its effect on engine transients [1, 5, 6, 7, 8, 9, 10, 11]. Most of the papers consider the pneumatic volume consisting of a single region with the performance described by the set of differential equations of mass and energy conservation. Conservation laws are added to the equation sets that describe the operation of components. The resulting set of the conservation laws and the equations of components’ operation compile the final set of equations, known as thermodynamic model. However, some of the works mentioned above still consider an isothermal process in the volume. Thus, only the pressure alternation is simulated keeping the temperature constant.

Wrapping up the analysis made above, the thermodynamic models including the models of pneumatic volumes take the assumptions listed below:

• In many cases the process in the volume is considered to be isothermal.

• The momentum conservation is omitted.

• The volume model can be called an “all-in-one-volume.”

One of the shortcomings of the methods overviewed before consists in omitting the momentum conservation law [12]; however, it plays a considerable role in an overall accuracy of the simulation.

The importance of the momentum conservation was also confirmed by Shi et al. in the study [13]. The authors deal with three engine models called “no volume effect,” “traditional simplified volume effect,” and “compressibility volume effect.” However, this study focuses only on the time delay of the transient because of the volume effect. The other parameters of the transient were not discussed in the study.

The gas turbine models allowing for the momentum conservation in pneumatic volumes are also presented in papers [14, 15, 16]. However, the authors do not indicate the range of tasks when it is essential to take into consideration the volume effect.

A new software PROOSIS for simulation in the area of propulsion allows more precise pneumatic volume description that includes the momentum conservation law. Henke et al. in their paper [17] introduce the PROOSIS capabilities in simulating gas turbine transients. One of their conclusions is that the time of transients caused by the volume effect is generally determined by the mass conservation. This is disputable, and it will be shown in the present chapter that indeed the transients are longer.

Wrapping up the overview of the existing methods to simulate the dynamics of the gas path with pneumatic volumes, one can draw the following conclusions:

• These volumes cause the delay of the transients.

• Most of the methods consider the phenomena of mass, energy, and momentum conservation, as well as hydraulic volume resistance; however, no method covers these phenomena simultaneously. Moreover, the used combinations of some phenomena are not compared.

• When the pneumatic volume algorithm is employed for control system design, it must be able to satisfy some specific needs related to the volume, such as oscillation analysis, time response evaluation, and analysis of the natural frequencies. The above tasks should be performed with minimum computation time expenses and, if possible, without simulation of the whole engine. However, none of the overviewed methods can solve these tasks.

• The authors of the overviewed studies did not focus on making faster the engine dynamic model with the algorithms describing pneumatic volumes.

The present chapter aims to overcome the mentioned bottlenecks and propose the best model of pneumatic volumes to be used as a component of either an engine dynamic model or an autonomous usage. Section 1 determines the set of differential equations to generally describe a pneumatic volume. This section introduces seven alternative volume models that adopt different simplifications. Section 2 provides numerical simulation of the volume by each model and, using the simulation results, the models’ comparison. The linearization algorithms are described in Section 3, and the linearization accuracy is studied in Section 4.

2. Mathematical models of the pneumatic volumes

2.1 Basic equations

To make a fast computational algorithm, the present approach assumes that a pneumatic volume has a constant transversal section. The proposed thermodynamic model deals with the gas path variables averaged over the radius and the circumference, i.e., it considers the averaged variables of a one-dimensional gas flow. Because most of the methods cited in the introduction use the total pressure and total temperature to characterize the flow (see [1, 2, 3, 6, 7, 8, 9, 10, 11, 12, 14, 15, 16, 17]), the same variables are employed in the present study. In the present approach, the difference between total and static parameters is ignored as the subject of the study is low-speed flow ( $\mathrm{M}<0.4$ ). For the same reason, static density is determined by the total pressure and temperature in the ideal gas state equation $\mathrm{p}=\mathrm{\rho RT}$ .

Let us introduce the volume of interest by Figure 1 and formulate the conservation laws for this volume. The mass conservation can be presented by

Figure 1.

Gas capacity design scheme.

$$\frac{\mathrm{dm}}{\mathrm{dt}}={\mathrm{W}}_{\mathrm{\text{in}}}-{\mathrm{W}}_{\mathrm{\text{out}}}$$ E1

Let us then use the ideal gas equation and express the mass as $\mathrm{m}=\mathrm{\rho V}=\frac{\mathrm{P}}{\mathrm{RT}}\mathrm{LA}$ , whence the mass conservation law transforms to

$$\frac{\mathrm{dm}}{\mathrm{dt}}=\frac{\mathrm{LA}}{\mathrm{R}}\left(\frac{1}{\mathrm{T}}\frac{\mathrm{dP}}{\mathrm{dt}}-\frac{\mathrm{P}}{{\mathrm{T}}^2}\frac{\mathrm{dT}}{\mathrm{dt}}\right)$$ E2

The internal energy conservation law for an adiabatic flow is written as

$$\frac{\mathrm{dU}}{\mathrm{dt}}={\mathrm{h}}_{\mathrm{\text{in}}}{\mathrm{W}}_{\mathrm{\text{in}}}-{\mathrm{h}}_{\mathrm{\text{out}}}{\mathrm{W}}_{\mathrm{\text{out}}}.$$ E3

After expressing the internal energy by the pressure, one obtains the energy conservation law that reflects the differential equation for the pressure:

$$\mathrm{U}={\mathrm{mc}}_{\mathrm{v}}\mathrm{T}=\frac{{\mathrm{LAc}}_{\mathrm{v}}}{\mathrm{R}}\mathrm{P},$$ E4

$$\frac{\mathrm{dU}}{\mathrm{dt}}=\frac{{\mathrm{LAc}}_{\mathrm{v}}}{\mathrm{R}}\frac{\mathrm{dP}}{\mathrm{dt}},$$ E5

$$\frac{\mathrm{dP}}{\mathrm{dt}}=\frac{\mathrm{\gamma R}}{\mathrm{LA}}\left({\mathrm{T}}_{\mathrm{\text{in}}}{\mathrm{W}}_{\mathrm{\text{in}}}-{\mathrm{TW}}_{\mathrm{\text{out}}}\right)$$ E6

Let us now get the differential equation for the temperature using Eqs. (6) and (2):

$$\frac{\mathrm{dT}}{\mathrm{dt}}=\frac{{\mathrm{RT}}^2}{\mathrm{PLA}}\left[\left(\mathrm{\gamma }\frac{{\mathrm{T}}_{\mathrm{\text{in}}}}{\mathrm{T}}-1\right){\mathrm{W}}_{\mathrm{\text{in}}}-\left(\mathrm{\gamma }-1\right){\mathrm{W}}_{\mathrm{\text{out}}}\right]$$ E7

The next equation to form presents the momentum conservation law expressed through the Darcy-Weisbach equation (it relates the pressure loss due to friction to the average velocity of the flow for an incompressible fluid):

$$\mathrm{\Delta }\mathrm{P}=\mathrm{\xi }\frac{\mathrm{RT}}{2{\mathrm{PA}}^2}{\mathrm{W}}^2.$$ E8

It helps to describe the relations between the variables of two elementary volumes situated to the left and to the right from the calculation point:

$$\frac{{\mathrm{dW}}_{\mathrm{\text{in}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}-\mathrm{\xi }\frac{{\mathrm{RT}}_{\mathrm{\text{in}}}}{2{\mathrm{P}}_{\mathrm{\text{in}}}{\mathrm{A}}^2}{\mathrm{W}}_{\mathrm{\text{in}}}^2\right);$$ E9

$$\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}-\mathrm{\xi }\frac{\mathrm{RT}}{2{\mathrm{PA}}^2}{\mathrm{W}}_{\mathrm{\text{out}}}^2\right),$$ E10

where $\mathrm{\xi }=\mathrm{\varsigma }\frac{\mathrm{L}}{2\mathrm{D}}$ is the Darcy friction factor, $\mathrm{\varsigma }$ is the specific frictional resistance, $\mathrm{D}=\frac{4\mathrm{A}}{\mathrm{\Pi }}$ is the hydraulic diameter, and $\mathrm{\Pi }$ is the wetted perimeter of the cross section.

Eqs. (6), (7), (9) and (10) constitute a closed set of first-order nonlinear differential equations. The set describes four parameters P, T, W_in, and W_out that characterize the volume dynamics.

As proven before, the processes in the gas path pneumatic volumes can be simulated with different precision.

The simplest method only regards the mass conservation (see [18]). According to this method, the time derivative of temperature in Eq. (2) is considered negligibly small. Modification of this approach is used in the paper [19] for compressor dynamics simulation. The method described by Shevyakov [20] suggests the energy conservation to be considered only when deducing the set of equations. The problem of calculation accuracy is omitted. The method used by Dobryansky [21] already considers the mass and energy conservation. One of last publications on volume dynamic modeling [22] is based on the same suppositions. However, the relation between the internal energy and the temperature employs the heat capacity at constant pressure instead of the constant volume heat capacity. The method considered by Jaw and Mattingly [1] is also based on mass and energy conservation, but Eq. (2) has no time derivative of temperature. The method described by Gurevich [23] already takes into consideration the difference between the static and total parameters but still neglects the momentum conservation. Such diversity of the methods for simulating the pneumatic volumes of gas turbines implies the necessity to perform their comparative study.

The present chapter introduces and compares three groups of pneumatic volume models.

The first group unites all isothermal models. The assumption about the minuteness of the second item on the right side of Eq. (2) is equivalent to an assumption about the isothermal process in the volume. Keeping the temperature constant requires heat exchange with the ambiance, and hence the volume process cannot be adiabatic in this case. Indeed, the volumes in real engines are not absolutely adiabatic because the heat exchange is always present between the working substance and the construction elements surrounding the cavity. However, the characteristic time of the heat exchange is several orders greater than that of the mass and energy accumulation in the volumes. This fact proves the use of the adiabatic models. Although the isothermal models are not the best option for the volume effect simulation, the present chapter uses them for comparing the errors of different models.

The second group includes the models based on the mass and energy conservation. The volume process in these models is adiabatic.

The third group consists of adiabatic models considering all conservation laws (mass, energy, and momentum).

2.2 Model 1.1: Isothermal volume without hydraulic resistance

Given the assumptions that are taken for this model ξ = 0, $\frac{\mathrm{dT}}{\mathrm{dt}}=0$ , it follows that $\mathrm{T}={\mathrm{T}}_{\mathrm{\text{in}}}$ -const.

Let us differentiate Eq. (6) and then substitute the flow rate time derivatives that correspond to the case ξ = 0. As a result, we have

$${\mathrm{\tau }}_1^2\frac{{\mathrm{d}}^2\mathrm{P}}{{\mathrm{dt}}^2}+\mathrm{P}=\frac{1}{2}\left({\mathrm{P}}_{\mathrm{\text{in}}}-{\mathrm{P}}_{\mathrm{\text{out}}}\right),$$ E11

where ${\mathrm{\tau }}_1=\frac{1}{2}{\mathrm{\tau }}_0$ , ${\mathrm{\tau }}_0=\frac{\mathrm{L}}{\mathrm{a}}$ is a time required for the disturbance to pass through the volume, and $\mathrm{a}=\sqrt{\mathrm{\gamma RT}}$ is the sonic velocity in the cavity.

Thus, the lossless isothermal volume is modeled by a single second-order linear differential equation, whose solution depends on input disturbances and the time constant ${\mathrm{\tau }}_1$ .

2.3 Model 1.2: Isothermal volume with hydraulic resistance (momentum conservation is omitted)

Having applied the condition $\frac{{\mathrm{dW}}_{\mathrm{\text{in}}}}{\mathrm{dt}}=\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=0$ for Eqs. (9) and (10), we get the following:

$${\mathrm{W}}_{\mathrm{\text{in}}}=\mathrm{A}\sqrt{\frac{2{\mathrm{P}}_{\mathrm{\text{in}}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}\right)}{\mathrm{\xi }{\mathrm{RT}}_{\mathrm{\text{in}}}}};{\mathrm{W}}_{\mathrm{\text{out}}}=\mathrm{A}\sqrt{\frac{2\mathrm{P}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{\xi RT}}}.$$ E12

Then, the Eq. (6) is transformed to

$$\frac{\mathrm{dP}}{\mathrm{dt}}=\frac{\mathrm{\gamma }}{\mathrm{L}}\sqrt{\frac{2\mathrm{RT}}{\mathrm{\xi }}}\left[\sqrt{{\mathrm{P}}_{\text{in}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}\right)}-\sqrt{\mathrm{P}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right)}\right].$$ E13

Thus, in the case of isothermal volume with hydraulic losses, the volume is described by a first-order nonlinear differential equation.

2.4 Model 1.3: Isothermal volume with hydraulic resistance (momentum conservation is taken into account)

The volume is modeled by the following set of equations consisting of Eq. (6) modified for the constant temperature condition and Eqs. (9) and (10):

$$\frac{\mathrm{dP}}{\mathrm{dt}}=\frac{\mathrm{\gamma RT}}{\mathrm{LA}}\left({\mathrm{W}}_{\mathrm{\text{in}}}-{\mathrm{W}}_{\mathrm{\text{out}}}\right);$$ E14

$$\frac{{\mathrm{dW}}_{\mathrm{\text{in}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}-\mathrm{\xi }\frac{{\mathrm{RT}}_{\mathrm{\text{in}}}}{2{\mathrm{P}}_{\mathrm{\text{in}}}{\mathrm{A}}^2}{\mathrm{W}}_{\mathrm{\text{in}}}^2\right);$$ E15

$$\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}-\mathrm{\xi }\frac{\mathrm{RT}}{2{\mathrm{PA}}^2}{\mathrm{W}}_{\mathrm{\text{out}}}^2\right).$$ E16

2.5 Model 2.1 based on mass and energy accumulation in the volume without hydraulic resistance

As momentum loss is neglected, we can state that $\frac{{\mathrm{dW}}_{\mathrm{\text{in}}}}{\mathrm{dt}}=\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=0$ . Since the pressure loss in Eqs. (9) and (10) is equal to zero, we have $\mathrm{P}={\mathrm{P}}_{\mathrm{\text{in}}}={\mathrm{P}}_{\mathrm{\text{out}}}$ . But then it follows from Eq. (1) that $\frac{\mathrm{dm}}{\mathrm{dt}}=0$ , and we can transform Eq. (7) to

$$\frac{\mathrm{dT}}{\mathrm{dt}}=\frac{\mathrm{\gamma }{\mathrm{RWT}}^2}{\mathrm{PLA}}\left(\frac{{\mathrm{T}}_{\mathrm{\text{in}}}}{\mathrm{T}}-1\right).$$ E17

2.6 Model 2.2 based on mass and energy accumulation in the volume with hydraulic resistance

By substituting Eq. (12) into Eq. (6) and Eq. (7), we arrive to the following equations for the model under analysis:

$$\frac{\mathrm{dP}}{\mathrm{dt}}=\frac{\mathrm{\gamma T}}{\mathrm{L}}\sqrt{\frac{2\mathrm{R}}{\mathrm{\xi }}}\left[\sqrt{\frac{{\mathrm{P}}_{\mathrm{\text{in}}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}\right)}{{\mathrm{T}}_{\mathrm{\text{in}}}}}-\sqrt{\frac{\mathrm{P}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{T}}}\right];$$ E18

$$\frac{\mathrm{dT}}{\mathrm{dt}}=\frac{{\mathrm{T}}^2}{\mathrm{P}}\frac{1}{\mathrm{L}}\sqrt{\frac{2\mathrm{R}}{\mathrm{\xi }}}\left[\begin{array}{l}\left(\mathrm{\gamma }\frac{{\mathrm{T}}_{\mathrm{\text{in}}}}{\mathrm{T}}-1\right)\sqrt{\frac{{\mathrm{P}}_{\mathrm{\text{in}}}\left({\mathrm{P}}_{\mathrm{inl}}-\mathrm{P}\right)}{{\mathrm{T}}_{\mathrm{\text{in}}}}}-\\ -\left(\mathrm{\gamma }-1\right)\sqrt{\frac{\mathrm{P}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{T}}}\end{array}\right].$$ E19

2.7 Model 3.1 based on mass and energy accumulation and momentum conservation in the volume without hydraulic resistance

The set of equations, constituting this model, consists of Eqs. (6) and (7), and also Eqs. (9) and (10) changed for the lossless conditions (ξ = 0):

$$\frac{{\mathrm{dW}}_{\mathrm{\text{in}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left({\mathrm{P}}_{\mathrm{\text{in}}}-\mathrm{P}\right);$$ E20

$$\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right).$$ E21

This model contains a contradiction that can be illustrated by the following example. The change of pressure in a volume inlet results in the pressure drop between the volume inlet and outlet. So the pressure drop in its turn makes the gas flow to become transient (see Eqs. (20) and (21)). Since the volume of interest is lossless, the pressure at the inlet and outlet will eventually become equal when the transient comes to the steady state. However, this will never happen as the inlet pressure has already changed, and the outlet pressure will remain immutable forever. Thus, the transient will not stop within this model.

Despite the above contradiction, the considered model is not expelled from the study, because it still can be used autonomously to estimate the dynamic process in the volume.

2.8 Model 3.2 based on mass and energy accumulation and momentum conservation in the volume with hydraulic resistance

As it has been mentioned above, this model is the most comprehensive. It is based on Eqs. (6), (7), (9) and (10).

3. Simulation

The language VisSim was used for programming the algorithm of volume dynamic simulation with the above models. To ensure the accuracy, the simulation was performed by different integration algorithms (Euler, Runge Kutta 2nd order, Runge Kutta 4th order, adaptive Runge Kutta 5th order) and with different integration steps. The simulation results were trusted only when all integration algorithms provided similar results and the integration step was small enough not to influence them.

As proposed in the paper [12], to simulate the transients, we have chosen the volume with standard geometrical characteristics L = 1 m and A = 1 m². This volume is placed between two infinite capacities, whose parameters are P_in = 300 kPa, P_out = 150 kPa, and T_in = T_out = 300 K. The first capacity was simulated as a single volume. Having experimented with different integration technics and diverse integration steps, we arrived to proper computation conditions at which the integration method and integration step do not influence simulation results. Under these conditions, the computations were conducted with the seven models described above. The results are plotted in Figure 2 (the disturbing factor is a pressure drop at the inlet ΔP = 10 kPa) and Figure 3 (the disturbing factor is a temperature drop at the inlet ΔТ = 50 K).

Figure 2.

Reaction of the volume parameters on the perturbation in the inlet pressure (single volume model).

Figure 3.

Reaction of the volume parameters on the perturbation in the inlet temperature (single volume model).

Using these figures, let us firstly analyze the dynamic performance of each model and then study the effect of volume split-off on simulated parameters.

3.1 Model dynamic performances

Since Model 3.2 is the most comprehensive, we will employ its performance as a pattern to compare the performances of the other models with it.

Plots of model 1.1 express the results in the constant-amplitude continuous oscillations. The model response has this form because Eq. (11) is similar to a conservative element [24]. Thus, this model cannot be used to properly simulate the volume effect within the gas path model. However, this model represents the frequency of the parameter oscillations well enough and can be used in a volume Eigen frequency analysis.

Models 1.2 and 2.2 give the similar results when simulating the pressure. The transient lags in this case. Its time constant is very small (about 0.002 s), much smaller than the total transient simulation time. Models 2.1 and 2.2 output a very similar response when simulating the temperature (however, as mentioned above, model 2.1 cannot simulate the pressure).

The pressure simulation using models 1.2., 1.3, 2.2, and 3.2 for the case of an inlet pressure perturbation gives the same values to the end of the transient. The pressure transient computed by model 3.1 ends with the different value, which is equal to the average between the inlet and outlet pressure. The difference appears because models 1.2., 1.3, 2.2, and 3.2 consider hydraulic losses in contrast to model 3.1.

The parameters simulated by models 1.3, 3.1, and 3.2 change according to a damped oscillation law. Obviously, this is because these models take into account the momentum conservation law.

As shown in the figures, model 1.3 has a bit higher frequency of oscillations than model 3.2. In general, the frequencies of different models are close to each other. Hence, all these models can be used when estimating the amplitude-frequency characteristic of the volume.

As regards the oscillation decay time, it is two times greater for model 1.3 and five times greater for model 3.1 than the model 3.2 time.

Let us now consider the effect of volume split-off on the simulated parameters.

3.2 Volume split-off effect simulation

The most advanced model 3.2 has been chosen to study this effect. Three cases were considered: the entire cavity not split, the cavity split in the axial direction into three equal volumes, and the cavity split into five volumes. The section area remained the same.

The presented above equations, corresponding to the model chosen, were applied to each one elementary volume. The output parameters (pressure, temperature, and flow rate) of one volume were the input parameters of the next volume. The specific frictional resistance $\mathrm{\xi }$ of each elementary volume was determined in the way that results in the total pressure loss equal to that of the non-split cavity.

The simulation results for the three cases are plotted in Figures 4 and 5. As seen in these figures, the transient plots corresponding to these cases are pretty similar, i.e., all of them obey the damped oscillation law. When the pressure disturbance is considered (Figure 4), the rate of the damping is greater for the three-volume and five-volume models than for the single volume model. For the temperature disturbance, the damping rates of all the models are equal. The fundamental frequencies for different volume numbers are very close as well. As to the amplitude, the three-volume and five-volume models have approximately equal amplitudes that are about 20% greater than that of the single volume model.

Figure 4.

Reaction of the volume parameters on the perturbation in the inlet pressure: (1, )—Single volume, (3, )—Three volumes, and (5, )—Five volumes.

Figure 5.

Reaction of the volume parameters on the perturbation in the inlet temperature: (1, )—Single volume, (3, )—Three volumes, and (5, )—Five volumes.

4. Linear analysis of volume effect

The differential equations of the models described in the previous section can be generally presented by

$$\frac{\mathrm{dy}}{\mathrm{dt}}=\mathrm{f}\left({\mathrm{x}}_1,\dots ,{\mathrm{x}}_{\mathrm{n}},\mathrm{y}\right).$$ E22

The nonlinear form of the equations does not allow the direct use of the universal methods [24, 25] that have been specially developed for the dynamic analysis of the set of linear differential equations. Hence, let us linearize the equations describing the volume effect and transform them introducing small variations of arguments from their steady-state magnitudes denoted by the subscript “0.” We arrive to

$$\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{y}\right)}{\mathrm{dt}}=\frac{\mathrm{df}}{{\mathrm{dx}}_1}\mathrm{\Delta }{\mathrm{x}}_1+\dots +\frac{\mathrm{df}}{{\mathrm{dx}}_{\mathrm{n}}}\mathrm{\Delta }{\mathrm{x}}_{\mathrm{n}}+\frac{\mathrm{df}}{\mathrm{dy}}\mathrm{\Delta }\mathrm{y},$$ E23

where $\mathrm{\Delta }{\mathrm{x}}_{\mathrm{i}}={\mathrm{x}}_{\mathrm{i}}-{\mathrm{x}}_{\mathrm{i}0}$ and $\mathrm{\Delta }\mathrm{y}=\mathrm{y}-{\mathrm{y}}_0$ . The time derivatives were evaluated here in the point ( ${\mathrm{x}}_{10}$ ,…, ${\mathrm{x}}_{\mathrm{n}0}$ , ${\mathrm{y}}_0$ ).

Let us transform absolute deviations to relative deviations $\mathrm{\delta y}=\frac{\mathrm{\Delta }\mathrm{y}}{{\mathrm{y}}_0}$ , $\mathrm{\delta }{\mathrm{x}}_1=\frac{\mathrm{\Delta }{\mathrm{x}}_1}{{\mathrm{x}}_{10}}$ ,…, $\mathrm{\delta }{\mathrm{x}}_{\mathrm{n}}=\frac{\mathrm{\Delta }{\mathrm{x}}_{\mathrm{n}}}{{\mathrm{x}}_{\mathrm{n}0}}$ . Eq. (23) is then changed to

$${\mathrm{y}}_0\frac{\mathrm{d}\left(\mathrm{\delta y}\right)}{\mathrm{dt}}=\frac{\mathrm{\partial }\mathrm{f}}{\mathrm{\partial }{\mathrm{x}}_1}{\mathrm{x}}_{10}\mathrm{\delta }{\mathrm{x}}_1+\dots +\frac{\mathrm{\partial }\mathrm{f}}{\mathrm{\partial }{\mathrm{x}}_{\mathrm{n}}}{\mathrm{x}}_{\mathrm{n}0}\mathrm{\delta }{\mathrm{x}}_{\mathrm{n}}+\frac{\mathrm{\partial }\mathrm{f}}{\mathrm{\partial }\mathrm{y}}{\mathrm{y}}_0\mathrm{\delta y}.$$ E24

Using the linearization principle described above, linear differential equations have been formed for all the models under consideration, and their analytical solutions for volume pressure and temperature were derived (see Appendix). These solutions are determined by the totality of the physical laws and assumptions that are used in each of the considered models.

Let us now analyze basic properties of these solutions. Specifically, in the next section we will determine the order of equations, their parameters, and type of transients that they describe.

5. Analysis of the transients in the volume based on the linearized equations of each model

5.1 Model 1.1

Solution (51) of this model corresponds to the undamped harmonic oscillations with the angular frequency $\mathrm{\omega }=\frac{1}{{\mathrm{\tau }}_1}$ . When L = 1 m and a₀ = 500 m·s⁻¹, then ${\mathrm{\tau }}_1$ =0.001 s.

5.2 Model 1.2

Eq. (55) corresponds to an aperiodic system, whose dynamics is described by the time constant ${\mathrm{\tau }}_{\mathrm{p}}$ .

For the rough estimation of the time constant and the gain coefficient, we can neglect the Darcy friction factor and assume that $\frac{{\mathrm{P}}_{\text{in}0}}{{\mathrm{P}}_0}\approx \frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}\approx 1$ . As a result, we arrive to

$${\mathrm{\tau }}_{\mathrm{p}}=\frac{1}{2}{\mathrm{\xi \tau }}_0{\mathrm{M}}_0;{\mathrm{K}}_{\mathrm{\text{in}}}^{\mathrm{P}}=\frac{1}{2};{\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{P}}=\frac{1}{2}.$$ E25

Provided that L = 1, ξ = 0.02, and a = 500 ms⁻¹, the time constant is ${\mathrm{\tau }}_{\mathrm{p}}$ = 0.00006 s.

It is obvious that the model simulates the volume like an almost inertia-free object.

5.3 Model 1.3

Let us analyze Eq. (61) that presents this model. The aperiodicity conditions can be given by an inequality:

$${\left(\frac{{\mathrm{\tau }}_{\mathrm{P}}}{1+0.25\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\right)}^2>4\frac{{\mathrm{\tau }}_0^2}{2+0.5\mathrm{\xi \gamma }{\mathrm{M}}_0^2}.$$ E26

The Darcy friction factor and squared Mach number are minuscule. Hence, this condition can be simplified to

$${\mathrm{\tau }}_{\mathrm{p}}^2>2{\mathrm{\tau }}_0^2\text{or}{\mathrm{\xi }}^2{\mathrm{\tau }}_0{\mathrm{M}}_0^2>8.$$ E27

As we see, this condition is not fulfilled. Thus, the dynamic processes in the volume have an oscillatory nature. The coefficients in the right side of Eq. (61) are positive. Hence, the system is robust, i.e., the oscillations relax. The time constant ${\mathrm{\tau }}_2$ determines the intensity of the relaxation. It is inverse to the real root α of the characteristic equation corresponding to differential Eq. (61):

$${\mathrm{\tau }}_2=\frac{1}{\mathrm{\alpha }}=-\frac{2{\mathrm{\tau }}_0^2}{{\mathrm{\tau }}_{\mathrm{p}}}\frac{1+0.25\mathrm{\xi \gamma }{\mathrm{M}}_0^2}{2+0.5\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\approx -\frac{{\mathrm{\tau }}_0^2}{{\mathrm{\tau }}_{\mathrm{P}}}=-\frac{2{\mathrm{\tau }}_0}{\mathrm{\xi }{\mathrm{M}}_0}=-\frac{2\mathrm{\tau }}{\mathrm{\xi }},$$ E28

where $\mathrm{\tau }=\frac{\mathrm{L}}{{\mathrm{c}}_0}$ is the time needed by the flow to cross the volume.

The frequency is equal to an absolute value of the imaginary root:

$$\mathrm{\omega }=\frac{\sqrt{\frac{{\mathrm{\tau }}_{\mathrm{p}}^2}{{\left(1+0.25\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)}^2}-\frac{4{\mathrm{\tau }}_0^2}{2+0.5\mathrm{\xi \gamma }{\mathrm{M}}_0^2}}}{\frac{2{\mathrm{\tau }}_0^2}{2+0.5\mathrm{\xi \gamma }{\mathrm{M}}_0^2}}\approx \frac{\sqrt{2}}{{\mathrm{\tau }}_0}=\frac{\sqrt{2}{\mathrm{a}}_0}{\mathrm{L}}.$$ E29

The evaluated time constant is much greater than that from model 1.2. The frequency is $\sqrt{2}$ times lower than the frequency estimated by model 1.1. However, this difference is acceptable for rough estimation.

5.4 Model 2.1

Eq. (64) corresponds to an aperiodic system, which characteristic time of the transient is given by the time constant ${\mathrm{\tau }}_{\mathrm{T}}$ .

If L = 1 m, с₀ = 100 m·s⁻¹, and the ratio of specific heats is 1.4, then ${\mathrm{\tau }}_{\mathrm{T}}$ = 0.007 s. This constant is small, but it is considerably greater than the time constant obtained for model 1.2.

5.5 Model 2.2

The left sides of Eqs. (70) and (71) that represent this model have the same order and similar coefficients. This proves the dynamics of pressure and temperature in the volume to be equal.

The aperiodicity condition for these equations can be formulated as

$${\left(\frac{\mathrm{\tau }}{\mathrm{\gamma }}+\frac{\mathrm{\gamma }+1}{2\mathrm{\gamma }}{\mathrm{\tau }}_{\mathrm{P}}\right)}^2>4\frac{{\mathrm{\tau \tau }}_{\mathrm{P}}}{\mathrm{\gamma }}.$$ E30

As ${\mathrm{\tau }}_{\mathrm{P}}=\frac{1}{2}{\mathrm{\xi \tau }}_0{\mathrm{M}}_0$ , the condition ${\mathrm{\tau }}_{\mathrm{P}}<<\mathrm{\tau }$ is fulfilled, and the aperiodicity condition is transformed to $\mathrm{\tau }<4{\mathrm{\gamma \tau }}_{\mathrm{P}}$ . It is obvious that this condition is fulfilled. Hence, the transients have an aperiodic form. The dynamics of the volume is determined by its time constants:

$${\mathrm{\tau }}_1\approx \frac{\mathrm{\tau }}{\mathrm{\gamma }}\text{and}{\mathrm{\tau }}_2\approx {\mathrm{\tau }}_{\mathrm{P}}.$$ E31

5.6 Model 3.1

Let us analyze Eqs. (79) and (80) derived for model 3.1. For doing so, we must form a characteristic equation, which is common for both equations:

$${\mathrm{s}}^3+{\mathrm{as}}^2+\mathrm{bs}+\mathrm{c}=0,$$ E32

where $\mathrm{a}=\frac{\mathrm{\gamma }}{\mathrm{\tau }}$ , $\mathrm{b}=\frac{4}{{\mathrm{\tau }}_0^2}$ , and $\mathrm{c}=\frac{4}{{{\mathrm{\tau \tau }}_0}^2}$ .

Let us use the method proposed by Gerolamo Cardano [26]. For this we first check whether the volume dynamics is oscillatory. The condition of the oscillations is $\mathrm{Q}>0$ , where

$$\mathrm{Q}={\left(\frac{\mathrm{p}}{3}\right)}^3+{\left(\frac{\mathrm{q}}{2}\right)}^2;\mathrm{p}=-\frac{{\mathrm{a}}^2}{3}+\mathrm{b}=-\frac{{\mathrm{\gamma }}^2}{3{\mathrm{\tau }}^2}+\frac{4}{{\mathrm{\tau \tau }}_0{\mathrm{M}}_0};\mathrm{q}=2{\left(\frac{\mathrm{a}}{3}\right)}^3-\frac{\mathrm{ab}}{3}+\mathrm{c}=2{\left(\frac{\mathrm{\gamma }}{3\mathrm{\tau }}\right)}^3-\frac{4\mathrm{\gamma }}{3{\mathrm{\tau \tau }}_0^2}+\frac{4}{{\mathrm{\tau \tau }}_0^2}.$$

Since (1) ${\mathrm{\tau }}_0<\mathrm{\tau }$ , (2) $\frac{{\mathrm{\gamma }}^2}{3{\mathrm{\tau }}^2}<<\frac{4}{{\mathrm{\tau }}_0^2}$ and $\frac{4\mathrm{\gamma }}{3{\mathrm{\tau \tau }}_0^2}<\frac{4}{{\mathrm{\tau \tau }}_0^2}$ , (3) $\mathrm{p}>0$ , and (4) $\mathrm{q}>0$ , the condition $\mathrm{Q}>0$ is fulfilled. This obviously means that the transient has an oscillatory character.

The characteristic Eq. (32) has a single real root s₁ and two complex conjugate roots s₂ and s₃:

$${\mathrm{s}}_1=\mathrm{A}+\mathrm{B}-\frac{\mathrm{a}}{3}=\mathrm{A}+\mathrm{B}-\frac{\mathrm{\gamma }}{3\mathrm{\tau }};$$ E33

$${\mathrm{s}}_{2,3}=\mathrm{\alpha }\pm \mathrm{i\omega }=-\frac{\mathrm{A}+\mathrm{B}}{2}-\frac{\mathrm{\gamma }}{3\mathrm{\tau }}\pm \mathrm{i}\frac{\mathrm{A}-\mathrm{B}}{2}\sqrt{3},$$ E34

where $\mathrm{A}=\sqrt[3]{-\frac{\mathrm{q}}{2}+\sqrt{\mathrm{Q}}}$ ; $\mathrm{B}=\sqrt[3]{-\frac{\mathrm{q}}{2}-\sqrt{\mathrm{Q}}}$ .

Analyzing Eqs. (33) and (34), we can see that they have some infinitesimal summands that can be neglected. In this way, we get the simplified equations:

$$\mathrm{p}\approx \mathrm{b}=\frac{4}{{\mathrm{\tau }}_0^2};\mathrm{q}\approx -\frac{\mathrm{ab}}{3}\mathrm{b}+\mathrm{c}=\frac{4{\mathrm{M}}_0}{{\mathrm{\tau }}_0^3}\left(1-\frac{\mathrm{\gamma }}{3}\right);\mathrm{Q}\approx \frac{1}{{\mathrm{\tau }}_0^6}\left[\mathrm{c}+{\mathrm{x}}^2\right]\approx {\left(\frac{4}{3{\mathrm{\tau }}_0^2}\right)}^3;$$ E35

$$\mathrm{A}\approx -\frac{1}{{\mathrm{\tau }}_0}\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}};\mathrm{B}\approx -\frac{1}{{\mathrm{\tau }}_0}\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}},$$ E36

where $\mathrm{x}=2{\mathrm{M}}_0\left(1-\frac{\mathrm{\gamma }}{3}\right)$ ; $\mathrm{c}={\left(\frac{4}{3}\right)}^3$ ;

$${\mathrm{s}}_1\approx -\frac{1}{{\mathrm{\tau }}_0}\left(\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\frac{1}{3}\mathrm{\gamma }{\mathrm{M}}_0\right);$$ E37

$$\mathrm{\alpha }\approx -\frac{1}{{\mathrm{\tau }}_0}\left[\frac{1}{3}\mathrm{\gamma }{\mathrm{M}}_0-\frac{1}{2}\left(\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}}\right)\right];$$ E38

$$\mathrm{\omega }\approx \frac{\sqrt{3}}{2{\mathrm{\tau }}_0}\left(\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}}-\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}}\right).$$ E39

Taking into account the Mach number being less than 0.3, these formulas for the characteristic equation roots can be further simplified. Let us change the equations to linear relations of the following form:

$$\mathrm{y}=\mathrm{y}\left({\mathrm{M}}_0=0\right)+\frac{\mathrm{dy}}{\mathrm{dM}}\left(\mathrm{M}=0\right){\mathrm{M}}_0.$$ E40

The required derivatives are.

$$\frac{{\mathrm{ds}}_1}{\mathrm{dM}}\left({\mathrm{M}}_0=0\right)=-\frac{1}{{\mathrm{\tau }}_0};\frac{\mathrm{d\alpha }}{\mathrm{dM}}\left({\mathrm{M}}_0=0\right)=-\frac{\mathrm{\gamma }-1}{2{\mathrm{\tau }}_0};\frac{\mathrm{d\omega }}{\mathrm{dM}}\left({\mathrm{M}}_0=0\right)=0,$$ E41

that is why,

$${\mathrm{s}}_1\approx -\frac{{\mathrm{M}}_0}{{\mathrm{\tau }}_0}=-\frac{1}{\mathrm{\tau }};\mathrm{\alpha }\approx -\frac{\mathrm{\gamma }-1}{2\mathrm{\tau }};\mathrm{\omega }\approx \frac{1}{{\mathrm{\tau }}_0}.$$ E42

5.7 Model 3.2

The differential Eqs. (88) and (89) of this model have the common characteristic Eq. (32), where

$$\mathrm{a}=\frac{\mathrm{\gamma }+4\mathrm{\xi }}{\mathrm{\tau }},\mathrm{b}=\frac{2}{{\mathrm{\tau }}_0^2}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right],\text{and}\mathrm{c}=\frac{2\mathrm{C}}{{{\mathrm{\tau \tau }}_0}^2}.$$ E43

The following parameters were evaluated by the Cardano’s method:

$$\mathrm{p}=-\frac{{\mathrm{a}}^2}{3}+\mathrm{b}\approx \frac{2}{{\mathrm{\tau }}_0^2}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right],\mathrm{q}=2{\left(\frac{\mathrm{a}}{3}\right)}^3-\frac{\mathrm{ab}}{3}+\mathrm{c}\approx \frac{4{\mathrm{M}}_0}{3{\mathrm{\tau }}_0^3}\left(3-\mathrm{\gamma }-4\mathrm{\xi }\right),$$ E44

$$\mathrm{Q}={\left(\frac{\mathrm{p}}{3}\right)}^3+{\left(\frac{\mathrm{q}}{2}\right)}^2\approx {\left\{\frac{2}{3{\mathrm{\tau }}_0^2}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right]\right\}}^3;$$ E45

$$\mathrm{A}\approx -\frac{1}{{\mathrm{\tau }}_0}\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}};\mathrm{B}\approx -\frac{1}{{\mathrm{\tau }}_0}\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}},$$ E46

where $\mathrm{x}=\frac{2}{3}{\mathrm{M}}_0\left(3-\mathrm{\gamma }-4\mathrm{\xi }\right)$ ; $\mathrm{c}={\left\{\frac{2}{3{\mathrm{\tau }}_0^2}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right]\right\}}^3$ .

The expressions for the parameters of the characteristic equation roots are similar to that of Eqs. (37) and (38) (for model 3.1):

$${\mathrm{s}}_1\approx -\frac{1}{{\mathrm{\tau }}_0}\left(\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\frac{1}{3}\left(\mathrm{\gamma }+4\mathrm{\xi }\right){\mathrm{M}}_0\right);$$ E47

$$\mathrm{\alpha }\approx -\frac{1}{{\mathrm{\tau }}_0}\left[\frac{1}{3}\left(\mathrm{\gamma }+4\mathrm{\xi }\right){\mathrm{M}}_0-\frac{1}{2}\left(\sqrt[3]{\mathrm{x}-\sqrt{\mathrm{c}+{\mathrm{x}}^2}}+\sqrt[3]{\mathrm{x}+\sqrt{\mathrm{c}+{\mathrm{x}}^2}}\right)\right].$$ E48

Eq. (39) for $\mathrm{\omega }$ remains the same.

In the same way as it was done for Eq. (78), the linearization of these expressions over the Mach number allows their simplification:

$${\mathrm{s}}_1\approx -\frac{1}{\mathrm{\tau }};\mathrm{\alpha }\approx -\frac{\mathrm{\gamma }-1+4\mathrm{\xi }}{2\mathrm{\tau }};\mathrm{\omega }\approx \frac{1}{{\mathrm{\tau }}_0}.$$ E49

It follows from Eq. (49) that the hydraulic resistance of the volume mostly affects the oscillation damping rate but does not influence the frequency and the aperiodic component of the transient.

5.8 Simulation results

To verify that the linearization did not introduce big error and the obtained results can be trusted, we have compared them with the original nonlinear models. Figures 6 and 7 illustrate this comparison by plotting pressure for the transients caused by inlet pressure and inlet temperature perturbations. As seen in the figures, the linear model correctly simulates the main features and parameters of the transient: oscillatory nature, time of the transient, oscillation frequency, and the magnitude of the first overshoot. It is worth to mention that just these performances are the subject of the dynamics analysis for the development of ACSs. The comparison results for all the models confirm that the dynamic behavior of the linear models agree with the behavior of the nonlinear models. This allows recommending the obtained linear models and corresponding analytical solutions for practical usage.

Figure 6.

Reaction of the volume pressure on the perturbation in the inlet pressure (single volume cavity, nonlinear model, and linear model).

Figure 7.

Reaction of the volume pressure on the perturbation in the inlet temperature (single volume cavity, nonlinear model, and linear model).

The dynamic parameters obtained as a result of the linearization are presented in Table 1 for all the models. The numerical values correspond here to the input conditions of the example: L = 1 m, c = 100 m/s, a = 500 m/s, and ξ = 0.42.

Model	Transient	Time constants		Eigen frequency ω
		Formula	Value, s	Formula	Value, 1/s
1.1	Oscillatory	∞	∞	$\frac{2}{{\mathrm{\tau }}_0}$	1000
1.2	Aperiodic	${\mathrm{\tau }}_{\mathrm{P}}\approx \frac{1}{2}{\mathrm{\xi \tau }}_0{\mathrm{M}}_0$	0.000084	—
1.3	Oscillatory	$\frac{2\mathrm{\tau }}{\mathrm{\xi }}$	0.0476	$\frac{\sqrt{2}}{{\mathrm{\tau }}_0}$	707
2.1	Aperiodic	$\frac{\mathrm{\tau }}{\mathrm{\gamma }}$	0.00714	—	—
2.2	Aperiodic	$\frac{\mathrm{\tau }}{\mathrm{\gamma }}$ , ${\mathrm{\tau }}_{\mathrm{P}}$	0.00714, 0.000084	—	—
3.1	Oscillatory	$\mathrm{\tau }$ , $\frac{2\mathrm{\tau }}{\mathrm{\gamma }-1}$	0.01, 0.05	$\frac{1}{{\mathrm{\tau }}_0}$	500
3.2	Oscillatory	$\mathrm{\tau }$ , $\frac{2\mathrm{\tau }}{\mathrm{\gamma }-1+4\mathrm{\xi }}$	0.01, 0.00962	$\frac{1}{{\mathrm{\tau }}_0}$	500

Table 1.

Dynamic parameters of the volume.

6. Conclusions

The following conclusions can be drawn on the results of the carried-out research:

1. The use of the momentum conservation law makes a tangible contribution in the transient state simulation. Thus, it cannot be omitted, when simulating the engine transients by the engine model with the volume model integrated.

2. The volume effect can be accurately simulated by the single volume model. The simulation of big connected volumes (e.g., annular manifolds of gas pumping units or station or trunk pipelines) requires deeper understanding and further researching to prove the model applicability.

3. The isothermal models are not recommended to be integrated into the gas path models because they do not correspond to the operating conditions in the engines.

4. The time of transients evaluated by the conventional volume models 2.1 and 2.2 is significantly lower against the models that among other consider momentum conservation. Hence, when this time delay effect is the subject of simulation, it is reasonable to use model 3.2.

5. In some instances, it may become important to study the volume effect on the frequency responses of the engine. In this case, we once more recommend model 3.2, which consists of Eqs. (6), (7), (9) and (10). The momentum transformation in the cavity causes oscillations of the parameters. The frequency of oscillations depends on the velocity of sound and the volume length only.

6. The time of the transients depends on the aperiodic component duration and the oscillation decay time. The aperiodic component duration in its turn depends on the time during which the gas crosses the volume. The decay time may be greater than the time of the aperiodic process.

7. The hydraulic resistance mostly affects the oscillation decay time. High hydraulic resistance reduces it. If the hydraulic losses are negligible, then the oscillation decay time is about five times longer than the aperiodic process duration.

8. The analytical method to solve the equations of volume dynamics makes it possible to determine the main dynamic properties of the volume and to get simple equations for determining the dynamic parameters on the basis of known geometrical characteristics and gas properties.

9. The obtained analytical equations and solutions (Eqs. (88)–(93) are the most accurate) can be implemented when developing the combined algorithm of engine dynamic simulation with the volume effect integrated. Integration of differential equations of the rotor dynamics will be performed iteratively using the required integration step, and the volume effect will be computed analytically. The application of this method will allow significant reduction of the operational time.

Acknowledgments

This work has been carried out with the support of the Ministry of Education and Science of Ukraine (Research Project No. D203-3/2019-П).

A.1 Linearized equations’ derivation

A.1.1 Model 1.1

Differential Eq. (11) is already linear, and thus its structure is conserved despite switching to the relative deviations:

$${\mathrm{\tau }}_1^2\frac{{\mathrm{d}}^2\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^2}+\mathrm{\delta P}=\frac{1}{2}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}-\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right).$$ E50

The system, which behaves like this, is oscillatory. The solution of this equation in the case of inlet or outlet pressure step is changed by $\mathrm{\delta P}={\mathrm{A}}_{\mathrm{P}}$ must be found as

$$\mathrm{\delta P}\left(\mathrm{t}\right)=\frac{1}{2}{\mathrm{A}}_{\mathrm{p}}\left(1-cos\frac{\mathrm{t}}{{\mathrm{\tau }}_1}\right).$$ E51

A.1.2 Model 1.2

Let us start from linearizing Eq. (13):

$$\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{P}\right)}{\mathrm{dt}}=\frac{\mathrm{\gamma }}{\mathrm{L}}\sqrt{\frac{2\mathrm{RT}}{\mathrm{\xi }}}\left[\frac{\left(2{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0\right)\mathrm{\Delta }{\mathrm{P}}_{\mathrm{inl}}-{\mathrm{P}}_{\mathrm{inl}0}\mathrm{\Delta }\mathrm{P}}{2\sqrt{{\mathrm{P}}_{\mathrm{inl}0}\left({\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0\right)}}\right.-\left.\frac{\left(2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}\right)\mathrm{\Delta }\mathrm{P}-{\mathrm{P}}_0\mathrm{\Delta }{\mathrm{P}}_{\text{out}0}}{2\sqrt{{\mathrm{P}}_0\left({\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}\right)}}\right].$$ E52

One must consider then

$${\mathrm{W}}_0=\mathrm{A}\sqrt{\frac{2{\mathrm{P}}_{\mathrm{inl}0}\left({\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0\right)}{\mathrm{\xi RT}}}=\mathrm{A}\sqrt{\frac{2{\mathrm{P}}_0\left({\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}\right)}{\mathrm{\xi RT}}},$$ E53

whence

$${\mathrm{P}}_0=\sqrt{{\mathrm{P}}_{\mathrm{inl}0}^2+{\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_{\text{out}0}}{2}\right)}^2}-\frac{{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_{\text{out}0}}{2}.$$ E54

Then, when switching to the relative deviations, we get

$${\mathrm{\tau }}_{\mathrm{P}}\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}={\mathrm{K}}_{\mathrm{inl}}^{\mathrm{P}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+{\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{P}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}},$$ E55

where ${\mathrm{\tau }}_{\mathrm{P}}=\frac{{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\mathrm{\text{out}}0}}{\mathrm{\xi \tau }}_0{\mathrm{M}}_0$ is a time constant; ${\mathrm{M}}_0=\frac{{\mathrm{c}}_0}{\mathrm{a}}$ is a Mach number; and ${\mathrm{K}}_{\mathrm{inl}}^{\mathrm{P}}=\frac{2{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}}\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}$ and ${\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{P}}=\frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}}$ are gain coefficients.

The solution of this equation in the case of inlet pressure perturbation $\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}={\mathrm{A}}_{\mathrm{P}}$ can be presented by

$$\mathrm{\delta P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{P}}{\mathrm{K}}_{\mathrm{inl}}^{\mathrm{P}}\left(1-{\mathrm{e}}^{-\frac{\mathrm{t}}{{\mathrm{\tau }}_{\mathrm{p}}}}\right).$$ E56

A.1.3 Model 1.3

Let us linearize Eqs. (6), (9) and (10):

$$\frac{\mathrm{LA}}{\mathrm{\gamma RT}}\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{P}\right)}{\mathrm{dt}}=\mathrm{\Delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\Delta }{\mathrm{W}}_{\mathrm{\text{out}}};$$ E57

$$\frac{\mathrm{d}\left(\mathrm{\Delta }{\mathrm{W}}_{\mathrm{inl}}\right)}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left[\mathrm{\Delta }{\mathrm{P}}_{\mathrm{inl}}-\mathrm{\Delta }\mathrm{P}-\frac{\mathrm{\xi RT}}{2{\mathrm{A}}^2}\left(\frac{2{\mathrm{W}}_0}{{\mathrm{P}}_{\mathrm{inl}0}}\mathrm{\Delta }{\mathrm{W}}_{\mathrm{inl}}-\frac{{\mathrm{W}}_0^2}{{\mathrm{P}}_{\mathrm{inl}0}^2}\mathrm{\Delta }{\mathrm{P}}_{\mathrm{inl}}\right)\right];$$ E58

$$\frac{\mathrm{d}\left(\mathrm{\Delta }{\mathrm{W}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}=\frac{2\mathrm{A}}{\mathrm{L}}\left[\mathrm{\Delta }\mathrm{P}-\mathrm{\Delta }{\mathrm{P}}_{\mathrm{\text{out}}}-\frac{\mathrm{\xi RT}}{2{\mathrm{A}}^2}\left(\frac{2{\mathrm{W}}_0}{{\mathrm{P}}_0}\mathrm{\Delta }{\mathrm{W}}_{\mathrm{\text{out}}}-\frac{{\mathrm{W}}_0^2}{{\mathrm{P}}_0^2}\mathrm{\Delta }\mathrm{P}\right)\right].$$ E59

Next, we differentiate the equation for the pressure and substitute the derivatives of airflow. Then, considering $\frac{{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}}\mathrm{\Delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\Delta }{\mathrm{W}}_{\mathrm{\text{out}}}\approx \frac{\mathrm{LA}}{\mathrm{\gamma RT}}\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{P}\right)}{\mathrm{dt}}$ , we get

$$\frac{{\mathrm{L}}^2}{2\mathrm{\gamma RT}}\frac{{\mathrm{d}}^2\left(\mathrm{\Delta }\mathrm{P}\right)}{{\mathrm{dt}}^2}=\left(\frac{1}{2}+\frac{\mathrm{\xi }{\mathrm{RTW}}_0^2}{4{\mathrm{A}}^2{\mathrm{P}}_{\mathrm{inl}0}^2}\right)\mathrm{\Delta }{\mathrm{P}}_{\mathrm{inl}}+\frac{1}{2}\mathrm{\Delta }{\mathrm{P}}_{\mathrm{\text{out}}}-\left(1+\frac{\mathrm{\xi }{\mathrm{RTW}}_0^2}{4{\mathrm{A}}^2{\mathrm{p}}_0^2}\right)\mathrm{\Delta }\mathrm{P}-\frac{\mathrm{\xi }{\mathrm{RTW}}_0}{2{\mathrm{A}}^2{\mathrm{p}}_0}\frac{\mathrm{LA}}{2\mathrm{\gamma RT}}\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{P}\right)}{\mathrm{dt}}.$$ E60

In a relative deviations format

$$\frac{{\mathrm{\tau }}_0^2}{2+0.5\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\frac{{\mathrm{d}}^2\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^2}+\frac{{\mathrm{\tau }}_{\mathrm{P}}}{1+0.25\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}=\frac{2+\mathrm{\xi \gamma }{\mathrm{M}}_0^2}{4+\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\frac{2}{4+\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}.$$ E61

The transient process, which is initiated by the inlet pressure perturbation $\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}={\mathrm{A}}_{\mathrm{P}}$ , is expressed as

$$\mathrm{P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{P}}\left[1-{\mathrm{e}}^{-\frac{\mathrm{\xi }}{2\mathrm{\tau }}\mathrm{t}}\left(cos\mathrm{\omega t}+\frac{\mathrm{\xi }}{2\sqrt{2}}sin\mathrm{\omega t}\right)\right],$$ E62

where $\mathrm{\omega }\approx \frac{\sqrt{2}}{{\mathrm{\tau }}_0}$ .

A.1.4 Model 2.1

The linearization of Eq. (17) outputs the equation in the absolute deviations

$$\frac{\mathrm{PLA}}{\mathrm{\gamma RG}}\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{T}\right)}{\mathrm{dt}}=-{\mathrm{T}}_{\mathrm{inl}0}\mathrm{\Delta }\mathrm{T}+{\mathrm{T}}_0\mathrm{\Delta }{\mathrm{T}}_{\mathrm{inl}},$$ E63

which in the relative deviations has the following form:

$${\mathrm{\tau }}_{\mathrm{T}}\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}+\mathrm{\delta T}=\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}},$$ E64

where ${\mathrm{\tau }}_{\mathrm{T}}=\frac{\mathrm{PAL}}{\mathrm{\gamma }{\mathrm{RG}}_0{\mathrm{T}}_0}=\frac{\mathrm{\tau }}{\mathrm{\gamma }}$ is a time constant.

The transient process, which is initiated by the inlet temperature perturbation $\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}={\mathrm{A}}_{\mathrm{T}}$ , is described as

$$\mathrm{\delta T}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\left(1-{\mathrm{e}}^{-\frac{\mathrm{t}}{{\mathrm{\tau }}_{\mathrm{T}}}}\right).$$ E65

A.1.5 Model 2.2

Let us transform Eq. (18) and linearize it:

$$\frac{\mathrm{L}}{\mathrm{\gamma }}\sqrt{\frac{\mathrm{\xi }}{2\mathrm{R}}}\frac{\mathrm{dP}}{\mathrm{dt}}=\mathrm{T}\left[\sqrt{\frac{{\mathrm{P}}_{\mathrm{inl}}\left({\mathrm{P}}_{\mathrm{inl}}-\mathrm{P}\right)}{{\mathrm{T}}_{\mathrm{inl}}}}-\sqrt{\frac{\mathrm{P}\left(\mathrm{P}-{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{T}}}\right];$$ E66

$$\begin{array}{ll}\frac{\mathrm{L}}{\mathrm{\gamma }}\sqrt{\frac{\mathrm{\xi }}{2\mathrm{R}}}\frac{\mathrm{d}\left(\mathrm{\Delta }\mathrm{P}\right)}{\mathrm{dt}}=& \frac{{\mathrm{W}}_0}{2\mathrm{A}}\sqrt{\frac{\mathrm{\xi R}}{2}}\left(\mathrm{\Delta }{\mathrm{T}}_{\mathrm{inl}}-\mathrm{\Delta }\mathrm{T}\right)\\ & +\frac{\mathrm{A}}{2{\mathrm{W}}_0}\sqrt{\frac{2}{\mathrm{\xi R}}}\left[\left(2{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0\right)\mathrm{\Delta }{\mathrm{P}}_{\mathrm{inl}}+{\mathrm{P}}_0\mathrm{\Delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right.\left.-\left({\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}\right)\mathrm{\Delta }\mathrm{P}\right].\end{array}$$ E67

Next, we transform the coefficients and change the equation to the relative deviations:

$${\mathrm{\tau }}_{\mathrm{P}}^{\text{'}}\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}={\mathrm{K}}_{\mathrm{T}}^{\mathrm{P}}\left(\mathrm{\delta T}-\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)+{\mathrm{K}}_{\mathrm{inl}}^{\mathrm{P}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+{\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{P}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}},$$ E68

where ${\mathrm{\tau }}_{\mathrm{P}}^{\text{'}}=2{\mathrm{\tau }}_{\mathrm{P}}\frac{{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}}$ ; ${\mathrm{K}}_{\mathrm{T}}^{\mathrm{P}}=\frac{1}{2}\mathrm{\xi \gamma }{\mathrm{M}}_0^2\frac{{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}}$ .

In a similar manner we transform Eq. (19):

$$\frac{2\mathrm{\tau }}{\mathrm{\gamma }+1}\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}+\mathrm{\delta T}=\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}+{\mathrm{K}}_{\mathrm{inl}}^{\mathrm{T}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+{\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{T}}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}-{\mathrm{K}}_{\mathrm{p}}^{\mathrm{T}}\mathrm{\delta P},$$ E69

where ${\mathrm{K}}_{\mathrm{inl}}^{\mathrm{T}}=\frac{2\left(\mathrm{\gamma }-1\right)}{\mathrm{\gamma }\left(\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2}\frac{2{\mathrm{P}}_{\mathrm{inl}0}-{\mathrm{P}}_0}{{\mathrm{P}}_0}\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}$ ; ${\mathrm{K}}_{\mathrm{\text{out}}}^{\mathrm{T}}=\frac{2\left(\mathrm{\gamma }-1\right)}{\mathrm{\gamma }\left(\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2}\frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}$ ;

$${\mathrm{K}}_{\mathrm{P}}^{\mathrm{T}}=\frac{2\left(\mathrm{\gamma }-1\right)}{\mathrm{\gamma }\left(\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2}\frac{{\mathrm{P}}_{\mathrm{inl}0}+2{\mathrm{P}}_0-{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}.$$

Having combined (68) and (69), we will get the differential equations for pressure and temperature:

$$\begin{array}{ll}\frac{{\mathrm{\tau \tau }}_{\mathrm{p}}}{\mathrm{\gamma }}\frac{{\mathrm{d}}^2\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^2}+\left(\frac{\mathrm{\tau }}{\mathrm{\gamma }}+\frac{\mathrm{\gamma }+1}{2\mathrm{\gamma }}{\mathrm{\tau }}_{\mathrm{P}}\right)\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}=& \frac{\mathrm{\tau }}{2\mathrm{\gamma }}\left(\frac{\mathrm{d\delta }{\mathrm{P}}_{\mathrm{inl}}}{\mathrm{dt}}+\frac{\mathrm{d\delta }{\mathrm{P}}_{\mathrm{\text{out}}}}{\mathrm{dt}}\right)+\frac{1}{2}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)\\ & -\frac{\mathrm{\tau }}{4}\mathrm{\xi }{\mathrm{M}}_0^2\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{\mathrm{dt}}\end{array}$$ E70

$$\begin{array}{ll}\frac{{\mathrm{\tau \tau }}_{\mathrm{P}}}{\mathrm{\gamma }}\frac{{\mathrm{d}}^2\left(\mathrm{\delta T}\right)}{{\mathrm{dt}}^2}+\left(\frac{\mathrm{\tau }}{\mathrm{\gamma }}+\frac{\mathrm{\gamma }+1}{2\mathrm{\gamma }}{\mathrm{\tau }}_{\mathrm{P}}\right)\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}+\mathrm{\delta T}=& \frac{\mathrm{\gamma }+1}{2\mathrm{\gamma }}{\mathrm{\tau }}_{\mathrm{P}}\frac{\mathrm{d\delta }{\mathrm{T}}_{\mathrm{inl}}}{\mathrm{dt}}+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\\ & -\frac{\mathrm{\gamma }-1}{\mathrm{\gamma }}\frac{{\mathrm{\tau }}_{\mathrm{P}}}{\mathrm{\xi \gamma }{\mathrm{M}}_0^2}\left(\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}\right).\end{array}$$ E71

The transient that is initiated by the pressure perturbation ${\mathrm{A}}_{\mathrm{P}}$ is described as

$$\mathrm{\delta P}\left(\mathrm{t}\right)=\frac{{\mathrm{A}}_{\mathrm{P}}}{2{\mathrm{\tau }}_{\mathrm{P}}}\left({\mathrm{Ae}}^{\mathrm{at}}+{\mathrm{Be}}^{\mathrm{bt}}+{\mathrm{\tau }}_{\mathrm{P}}\right),$$ E72

where $\mathrm{A}=\frac{\mathrm{a}+\mathrm{d}}{\mathrm{a}\left(\mathrm{a}-\mathrm{b}\right)}$ ; $\mathrm{B}=\frac{\mathrm{b}+\mathrm{d}}{\mathrm{b}\left(\mathrm{b}-\mathrm{a}\right)}$ ; $\mathrm{d}=\frac{\mathrm{\gamma }}{\mathrm{\tau }}$ ;

$$\mathrm{a}=\frac{1}{{\mathrm{\tau \tau }}_{\mathrm{P}}}\left[-2\mathrm{\tau }-\left(\mathrm{\gamma }+1\right){\mathrm{\tau }}_{\mathrm{P}}-\sqrt{{\left(2\mathrm{\tau }+\left(\mathrm{\gamma }+1\right){\mathrm{\tau }}_{\mathrm{P}}\right)}^2-16\mathrm{k}{\mathrm{\tau \tau }}_{\mathrm{P}}}\right];$$

$$\mathrm{b}=\frac{1}{{\mathrm{\tau \tau }}_{\mathrm{P}}}\left[-2\mathrm{\tau }-\left(\mathrm{\gamma }+1\right){\mathrm{\tau }}_{\mathrm{P}}+\sqrt{{\left(2\mathrm{\tau }+\left(\mathrm{\gamma }+1\right){\mathrm{\tau }}_{\mathrm{P}}\right)}^2-16{\mathrm{\gamma \tau \tau }}_{\mathrm{P}}}\right];$$

$$\mathrm{\delta T}\left(\mathrm{t}\right)=\frac{{\mathrm{A}}_{\mathrm{P}}\left(\mathrm{\gamma }-1\right)}{\mathrm{\xi \gamma \tau }{\mathrm{M}}_0^2}\left({\mathrm{A}}_1{\mathrm{e}}^{\mathrm{at}}+{\mathrm{B}}_1{\mathrm{e}}^{\mathrm{bt}}\right),$$ E73

where ${\mathrm{A}}_1=\frac{\mathrm{a}}{\mathrm{a}\left(\mathrm{a}-\mathrm{b}\right)}$ ; ${\mathrm{B}}_1=\frac{\mathrm{b}}{\mathrm{b}\left(\mathrm{b}-\mathrm{a}\right)}$ .

The transient state, which is initiated by the temperature perturbation $\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}={\mathrm{A}}_{\mathrm{T}}$ , is described as

$$\mathrm{\delta P}\left(\mathrm{t}\right)=-\frac{{\mathrm{A}}_{\mathrm{T}}\mathrm{\xi \gamma }{\mathrm{M}}_0^2}{4{\mathrm{\tau }}_{\mathrm{P}}}\left({\mathrm{A}}_1{\mathrm{e}}^{\mathrm{at}}+{\mathrm{B}}_1{\mathrm{e}}^{\mathrm{bt}}\right);$$ E74

$$\mathrm{\delta T}\left(\mathrm{t}\right)=\frac{{\mathrm{A}}_{\mathrm{T}}\left(\mathrm{\gamma }+1\right)}{2\mathrm{\tau }}\left({\mathrm{A}}_2{\mathrm{e}}^{\mathrm{at}}+{\mathrm{B}}_2{\mathrm{e}}^{\mathrm{bt}}+\frac{2\mathrm{\tau }}{\mathrm{\gamma }+1}\right),$$ E75

where ${\mathrm{A}}_2=\frac{\mathrm{a}+{\mathrm{d}}_1}{\mathrm{a}\left(\mathrm{a}-\mathrm{b}\right)}$ ; ${\mathrm{B}}_1=\frac{\mathrm{b}+{\mathrm{d}}_1}{\mathrm{b}\left(\mathrm{b}-\mathrm{a}\right)}$ ; ${\mathrm{d}}_1=\frac{2\mathrm{\gamma }}{{\mathrm{\tau }}_{\mathrm{P}}\left(\mathrm{\gamma }+1\right)}$ .

A.1.6 Model 3.1

The linearized Eqs. (6), (7), (20), and (21) in the relative deviations format are

$$\frac{\mathrm{\tau }}{\mathrm{\gamma }}\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}=\mathrm{\delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\delta }{\mathrm{W}}_{\mathrm{\text{out}}}+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}-\mathrm{\delta T};$$ E76

$$\frac{\mathrm{\tau }}{\mathrm{\gamma }}\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}=-\mathrm{\delta T}+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}+\frac{\mathrm{\gamma }-1}{\mathrm{\gamma }}\left(\mathrm{\delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\delta }{\mathrm{W}}_{\mathrm{\text{out}}}\right);$$ E77

$$\frac{1}{2}{\mathrm{\gamma \tau }}_0{\mathrm{M}}_0\left(\mathrm{\delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\delta }{\mathrm{W}}_{\mathrm{\text{out}}}\right)=\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}-2\mathrm{\delta P}.$$ E78

Let us transform Eqs. (76)–(78) to get the differential equations for the pressure and the temperature:

$$\begin{array}{c}\frac{1}{4}{\mathrm{\tau }}^2{\mathrm{\tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^3\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^3}+\frac{1}{4}{\mathrm{\gamma \tau \tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^2\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^2}+\mathrm{\tau }\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}=\\ =\frac{1}{4}{\mathrm{\gamma \tau \tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^2{\mathrm{T}}_{\mathrm{inl}}}{{\mathrm{dt}}^2}+\frac{1}{2}\mathrm{\tau }\left(\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}\right)+\frac{1}{2}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right).\end{array}$$ E79

$$\begin{array}{c}\frac{1}{4}{\mathrm{\tau }}^2{\mathrm{\tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^3\left(\mathrm{\delta T}\right)}{{\mathrm{dt}}^3}+\frac{1}{4}{\mathrm{\gamma \tau \tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^2\left(\mathrm{\delta T}\right)}{{\mathrm{dt}}^2}+\\ +\mathrm{\tau }\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}+\mathrm{\delta T}=\frac{1}{4}{\mathrm{\gamma \tau \tau }}_0{\mathrm{M}}_0\frac{{\mathrm{d}}^2{\mathrm{T}}_{\mathrm{inl}}}{{\mathrm{dt}}^2}+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}+\frac{\mathrm{\gamma }-1}{2\mathrm{\gamma }}\mathrm{\tau }\left(\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}\right).\end{array}$$ E80

The transient that is initiated by the temperature perturbation $\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}$ = ${\mathrm{A}}_{\mathrm{T}}$ is described as

$$\mathrm{\delta T}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{\mathrm{\gamma }}{\mathrm{\tau }}\left({\mathrm{A}}_3{\mathrm{e}}^{\mathrm{\alpha t}}sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_1\right)+{\mathrm{B}}_3{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+\mathrm{K}\right),$$ E81

where α and ω are expressed by (38) and (39); ${\mathrm{A}}_3=\frac{1}{\mathrm{\omega }}\sqrt{\frac{{\left({\mathrm{\alpha }}^2-{\mathrm{\omega }}^2+\frac{4}{{\mathrm{\gamma \tau }}_0^2}\right)}^2+{\left(2\mathrm{\alpha \omega }\right)}^2}{\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)\left[{\left(\mathrm{\alpha }-{\mathrm{s}}_1\right)}^2+{\mathrm{\omega }}^2\right]}}$ ; ${\mathrm{B}}_3=\frac{{\mathrm{s}}_1^2+\frac{4}{{\mathrm{\gamma \tau }}_0^2}}{{\mathrm{s}}_1\left[{\left({\mathrm{s}}_1-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2\right]}$ ; $\mathrm{K}=-\frac{4}{{\mathrm{\gamma \tau }}_0^2{\mathrm{s}}_1\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)}$ ; $\mathrm{\beta }=\mathrm{\text{arctg}}\frac{2\mathrm{\alpha \omega }}{{\mathrm{\alpha }}^2-{\mathrm{\omega }}^2+\frac{4}{{\mathrm{\gamma \tau }}_0^2}}-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }-{\mathrm{s}}_1}-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }}$ ;

$$\mathrm{\delta P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{\mathrm{\gamma }}{\mathrm{\tau }}\left({\mathrm{A}}_4{\mathrm{e}}^{\mathrm{\alpha t}}sin\left(\mathrm{\omega }\mathrm{t}+\mathrm{\beta }\right)+{\mathrm{B}}_4{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}\right),$$ E82

where ${\mathrm{A}}_4=\frac{1}{\mathrm{\omega }}\sqrt{\frac{{\left({\mathrm{\alpha }}^2-{\mathrm{\omega }}^2\right)}^2+{\left(2\mathrm{\alpha \omega }\right)}^2}{\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)\left[{\left(\mathrm{\alpha }-{\mathrm{s}}_1\right)}^2+{\mathrm{\omega }}^2\right]}}$ ; ${\mathrm{B}}_4=\frac{{\mathrm{s}}_1^2}{{\mathrm{s}}_1\left[{\left({\mathrm{s}}_1-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2\right]}$ .

The transient that is initiated by the pressure perturbation ${\mathrm{A}}_{\mathrm{P}}$ is described as

$$\mathrm{\delta T}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{2\left(\mathrm{\gamma }-1\right)}{{\mathrm{\gamma \tau }}_0^2}\left({\mathrm{A}}_5{\mathrm{e}}^{\mathrm{\alpha t}}sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_1\right)+{\mathrm{B}}_5{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}\right),$$ E83

where ${\mathrm{A}}_5=\frac{1}{\mathrm{\omega }}\sqrt{\frac{1}{{\left(\mathrm{\alpha }-{\mathrm{s}}_1\right)}^2+{\mathrm{\omega }}^2}}$ ; ${\mathrm{B}}_5=\frac{1}{{\left({\mathrm{s}}_1-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}$ ; ${\mathrm{\beta }}_1=-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }-{\mathrm{s}}_1}$ ;

$$\mathrm{\delta P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{2}{{\mathrm{\tau }}_0^2}\left({\mathrm{A}}_6{\mathrm{e}}^{\mathrm{\alpha t}}sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_2\right)+{\mathrm{B}}_6{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+{\mathrm{K}}_2\right),$$ E84

where ${\mathrm{A}}_6=\frac{1}{\mathrm{\omega }}\sqrt{\frac{{\left(\frac{{\mathrm{\tau }}_0^2}{2\mathrm{\tau }}+\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}{\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)\left[{\left(\mathrm{\alpha }-{\mathrm{s}}_1\right)}^2+{\mathrm{\omega }}^2\right]}}$ ; ${\mathrm{B}}_6=\frac{{\mathrm{s}}_1+\frac{{\mathrm{\tau }}_0^2}{2\mathrm{\tau }}}{{\mathrm{s}}_1\left[{\left({\mathrm{s}}_1-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2\right]}$ ; ${\mathrm{K}}_2=-\frac{{\mathrm{\tau }}_0^2}{2\mathrm{\tau }{\mathrm{s}}_1\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)}$ ; ${\mathrm{\beta }}_2=-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }+\frac{{\mathrm{\tau }}_0^2}{2\mathrm{\tau }}}-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }-{\mathrm{s}}_1}-\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\mathrm{\alpha }}$ .

A.1.7 Model 3.2

The model consists of Eqs. (6), (7), (9) and (10). Linearized Eqs. (6) and (7) are of the format presented in Eqs. (76) and (78). As a result of the linearization, we get the missing difference $\frac{{\mathrm{dW}}_{\mathrm{inl}}}{\mathrm{dt}}-\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}$ :

$$\begin{array}{c}\frac{{\mathrm{dW}}_{\mathrm{inl}}}{\mathrm{dt}}-\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=\frac{2}{{\mathrm{\gamma \tau }}_0{\mathrm{M}}_0}\{\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\frac{{\mathrm{P}}_0}{{\mathrm{P}}_{\mathrm{inl}0}}\right)\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}-\left(2+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)+\mathrm{\delta P}\\ -\mathrm{\xi \gamma }{\mathrm{M}}_0^2\left[2\left(\mathrm{\delta }{\mathrm{W}}_{\mathrm{inl}}-\mathrm{\delta }{\mathrm{W}}_{\mathrm{\text{out}}}\right)+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}-\mathrm{\delta T}\right]\}.\end{array}$$ E85

On the other hand, from Eq. (77) we get

$$\frac{{\mathrm{dW}}_{\mathrm{inl}}}{\mathrm{dt}}-\frac{{\mathrm{dW}}_{\mathrm{\text{out}}}}{\mathrm{dt}}=\frac{\mathrm{\gamma }}{\mathrm{\gamma }-1}\left(\frac{\mathrm{\tau }}{\mathrm{\gamma }}\frac{{\mathrm{d}}^2\left(\mathrm{\delta T}\right)}{{\mathrm{dt}}^2}+\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}-\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{\mathrm{dt}}\right).$$ E86

Let us determine $\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}$ from (76) and (78):

$$\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}=\frac{\mathrm{\gamma }}{\mathrm{\gamma }-1}\left[\frac{\mathrm{d}\left(\mathrm{\delta T}\right)}{\mathrm{dt}}+\frac{1}{\mathrm{\tau }}\left(\mathrm{\delta T}-\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)\right].$$ E87

Having equalized the right sides of Eqs. (85) and (86), derived the obtained equation, and substituted the derivative (87) in it, we get a differential equation for the temperature in the volume:

$$\begin{array}{c}\frac{{\mathrm{\tau \tau }}_0^2}{2\mathrm{C}}\frac{{\mathrm{d}}^3\left(\mathrm{\delta }\mathrm{T}\right)}{{\mathrm{dt}}^3}+\frac{{\mathrm{\tau }}_0^2}{\mathrm{C}}\left(0.5\mathrm{\gamma }+2\mathrm{\xi }\right)\frac{{\mathrm{d}}^2\left(\mathrm{\delta }\mathrm{T}\right)}{{\mathrm{dt}}^2}+\frac{\mathrm{\tau }}{\mathrm{C}}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right]\frac{\mathrm{d}\left(\mathrm{\delta }\mathrm{T}\right)}{\mathrm{dt}}+\mathrm{\delta }\mathrm{T}\\ =\frac{{\mathrm{\gamma \tau }}_0^2}{2\mathrm{C}}\frac{{\mathrm{d}}^2\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{{\mathrm{dt}}^2}+\frac{\left(\mathrm{\gamma }+1\right){\mathrm{\xi \tau }}_0{\mathrm{M}}_0}{\mathrm{C}}\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\\ +\frac{\mathrm{\tau }\left(\mathrm{\gamma }-1\right)}{\mathrm{C}}\left[\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{{\mathrm{P}}_{\mathrm{\text{out}}0}}{{\mathrm{P}}_0}\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}\right]\end{array}$$ E88

where $\mathrm{C}=2+\mathrm{\xi \gamma }{\mathrm{M}}_0^2$ .

Let us use the Laplace transform and transfer functions ${\mathrm{W}}_{\mathrm{T}}^{\mathrm{p}}\left(\mathrm{s}\right)$ , ${\mathrm{W}}_{{\mathrm{p}}_{\mathrm{inl}}}^{\mathrm{T}}\left(\mathrm{s}\right)$ , ${\mathrm{W}}_{\mathrm{\text{out}}}^{\mathrm{T}}\left(\mathrm{s}\right)$ , and ${\mathrm{W}}_{{\mathrm{T}}_{\mathrm{inl}}}^{\mathrm{T}}\left(\mathrm{s}\right)$ to obtain the equation for the pressure. The transfer functions we will get from Eqs. (87) and (88). The final result is

$$\begin{array}{c}\frac{{\mathrm{\tau \tau }}_0^2}{2\mathrm{C}}\frac{{\mathrm{d}}^3\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^3}+\frac{{\mathrm{\tau }}_0^2}{\mathrm{C}}\left(0.5\mathrm{\gamma }+2\mathrm{\xi }\right)\frac{{\mathrm{d}}^2\left(\mathrm{\delta P}\right)}{{\mathrm{dt}}^2}+\frac{\mathrm{\tau }}{\mathrm{C}}\left[2+\left(2\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0^2\right]\frac{\mathrm{d}\left(\mathrm{\delta P}\right)}{\mathrm{dt}}+\mathrm{\delta P}=\\ =\frac{{\mathrm{\gamma \tau }}_0^2}{2\mathrm{C}}\frac{{\mathrm{d}}^2\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{{\mathrm{dt}}^2}+\frac{{\mathrm{\tau }}_0\mathrm{\xi \gamma }{\mathrm{M}}_0}{\mathrm{C}}\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{\mathrm{\tau }}{\mathrm{C}}\left[\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}\right)}{\mathrm{dt}}+\frac{{\mathrm{P}}_{\mathrm{\text{out}}0}}{{\mathrm{P}}_0}\frac{\mathrm{d}\left(\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right)}{\mathrm{dt}}\right]+\\ +\frac{1}{\mathrm{C}}\left[\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)\mathrm{\delta }{\mathrm{P}}_{\mathrm{inl}}+\frac{{\mathrm{P}}_{\text{out}0}}{{\mathrm{P}}_0}\mathrm{\delta }{\mathrm{P}}_{\mathrm{\text{out}}}\right].\end{array}$$ E89

The transient state, which is initiated by the temperature perturbation $\mathrm{\delta }{\mathrm{T}}_{\mathrm{inl}}$ , (the magnitude of the step is ${\mathrm{A}}_{\mathrm{T}}$ ) is described as

$$\mathrm{\delta T}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{\mathrm{\gamma }}{\mathrm{\tau }}\left(\begin{array}{l}{\mathrm{A}}_6\cdot 1\left(\mathrm{t}\right)-\left({\mathrm{A}}_6+{\mathrm{B}}_6\right){\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+\\ +\frac{{\mathrm{B}}_6\sqrt{{\left(\frac{{\mathrm{K}}_3}{{\mathrm{B}}_6}-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}}{\mathrm{\omega }}{\mathrm{e}}^{\mathrm{\alpha t}}\cdot sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_3\right)\end{array}\right),$$ E90

where α and ω are expressed by (48) and (39): ${\mathrm{A}}_6=-\frac{{\mathrm{b}}_1}{{\mathrm{s}}_1\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)}$ ; ${\mathrm{B}}_6=\frac{\left(2\mathrm{\alpha }{\mathrm{s}}_1-{\mathrm{s}}_1^2\right){\mathrm{A}}_6-{\mathrm{s}}_1-{\mathrm{a}}_1}{{\mathrm{\alpha }}^2+{\mathrm{\omega }}^2-{\mathrm{s}}_1\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right)}$ ; ${\mathrm{K}}_3=1+{\mathrm{s}}_1{\mathrm{A}}_6-\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right){\mathrm{B}}_6$ ; ${\mathrm{\beta }}_3=\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\frac{{\mathrm{K}}_3}{{\mathrm{B}}_6}-\mathrm{\alpha }}$ ; ${\mathrm{a}}_1=\frac{2\left(\mathrm{\gamma }+1\right)\mathrm{\xi }{\mathrm{M}}_0}{{\mathrm{\gamma \tau }}_0}$ ; ${\mathrm{b}}_1=\frac{2\mathrm{C}}{{\mathrm{\gamma \tau }}_0^2}$ ;

$$\mathrm{\delta P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{T}}\frac{\mathrm{\gamma }}{\mathrm{\tau }}\left({\mathrm{A}}_7{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+\frac{{\mathrm{B}}_7\sqrt{{\left(\frac{{\mathrm{K}}_4}{{\mathrm{B}}_7}-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}}{\mathrm{\omega }}{\mathrm{e}}^{\mathrm{\alpha t}}\cdot sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_4\right)\right),$$ E91

where ${\mathrm{A}}_7=-{\mathrm{B}}_7$ ; ${\mathrm{B}}_7=-\frac{{\mathrm{a}}_2+{\mathrm{s}}_1}{{\mathrm{\alpha }}^2+{\mathrm{\omega }}^2-{\mathrm{s}}_1\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right)}$ ; ${\mathrm{K}}_4=\frac{\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right){\mathrm{A}}_7-{\mathrm{a}}_2}{{\mathrm{s}}_1}$ ; ${\mathrm{\beta }}_4=\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\frac{{\mathrm{D}}_2}{{\mathrm{C}}_2}-\mathrm{\alpha }}$ ; ${\mathrm{a}}_2=\frac{2\mathrm{\xi \gamma }{\mathrm{M}}_0}{{\mathrm{\tau }}_0\mathrm{\gamma }}$ .

The transient that is caused by the pressure perturbation ${\mathrm{A}}_{\mathrm{P}}$ is described as

$$\mathrm{\delta T}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{P}}\frac{2\left(\mathrm{\gamma }-1\right)}{{\mathrm{\tau }}_0^2}\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)\left({\mathrm{A}}_8{\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+\frac{{\mathrm{B}}_8\sqrt{{\left(\frac{{\mathrm{K}}_5}{{\mathrm{B}}_8}-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}}{\mathrm{\omega }}{\mathrm{e}}^{\mathrm{\alpha t}}\cdot sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_5\right)\right),$$ E92

where ${\mathrm{A}}_8=-{\mathrm{B}}_8$ ; ${\mathrm{B}}_8=-\frac{1}{{\mathrm{\alpha }}^2+{\mathrm{\omega }}^2-{\mathrm{s}}_1\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right)}$ ; ${\mathrm{K}}_5=\frac{\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right){\mathrm{A}}_8-1}{{\mathrm{s}}_1}$ ; ${\mathrm{\beta }}_1=\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\frac{{\mathrm{K}}_5}{{\mathrm{B}}_8}-\mathrm{\alpha }}$ .

$$\mathrm{\delta P}\left(\mathrm{t}\right)={\mathrm{A}}_{\mathrm{P}}\frac{2}{{\mathrm{\tau }}_0^2}\left(\frac{{\mathrm{P}}_{\mathrm{inl}0}}{{\mathrm{P}}_0}+\mathrm{\xi \gamma }{\mathrm{M}}_0^2\right)\left({\mathrm{A}}_9\cdot 1\left(\mathrm{t}\right)-\left({\mathrm{A}}_9+{\mathrm{B}}_9\right){\mathrm{e}}^{{\mathrm{s}}_1\mathrm{t}}+\frac{{\mathrm{B}}_9\sqrt{{\left(\frac{{\mathrm{K}}_6}{{\mathrm{B}}_9}-\mathrm{\alpha }\right)}^2+{\mathrm{\omega }}^2}}{\mathrm{\omega }}{\mathrm{e}}^{\mathrm{\alpha t}}\cdot sin\left(\mathrm{\omega }\mathrm{t}+{\mathrm{\beta }}_6\right)\right),$$ E93

where ${\mathrm{A}}_9=-\frac{{\mathrm{b}}_3}{{\mathrm{s}}_1\left({\mathrm{\alpha }}^2+{\mathrm{\omega }}^2\right)}$ ; ${\mathrm{B}}_9=\frac{\left(2\mathrm{\alpha }{\mathrm{s}}_1-{\mathrm{s}}_1^2\right){\mathrm{A}}_9-1}{{\mathrm{\alpha }}^2+{\mathrm{\omega }}^2-{\mathrm{s}}_1\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right)}$ ; ${\mathrm{K}}_6={\mathrm{s}}_1{\mathrm{A}}_9-\left(2\mathrm{\alpha }-{\mathrm{s}}_1\right){\mathrm{B}}_9$ ; ${\mathrm{\beta }}_6=\mathrm{\text{arctg}}\frac{\mathrm{\omega }}{\frac{{\mathrm{K}}_6}{{\mathrm{B}}_9}-\mathrm{\alpha }}$ ; ${\mathrm{b}}_3=\frac{1}{\mathrm{\tau }}$ .

Nomenclature

a	velocity of sound
c v	specific heat capacity at constant volume condition
c	gas velocity
D	hydraulic diameter
A	constant section area
W	gas flow rate
h	enthalpy
ɣ	ratio of specific heats
L	length of the volume
M	Mach number
m	mass
P	pressure
R	gas constant
T	temperature
U	internal energy
V	volume
ρ	static density
ζ	specific frictional resistance
ξ	friction factor
ω	Eigen frequency
τ	time constant
Indexes
in	inlet
out	outlet
0	initial static value