**3. Calculation procedure**

424 Thermodynamics – Interaction Studies – Solids, Liquids and Gases

(15)

*ρ V A cons t* tan (16)

*<sup>A</sup>* (17)

(18)

0 0 0 *P ρ T P ρ T* 

The taking logarithm and then differentiating of relation (16), and also using of the relations

( ) *<sup>A</sup> dA F T dT*

1 1 ( ) ( ) ( ) 2 () *FT CT A P a T H T*

The integration of equation (17) between the critical state (*A\*, T\**) and the supersonic state (*A,* 

To find parameters *ρ* and *A,* the integrals of functions *Fρ(T)* and *FA(T)* should be found. As the analytical procedure is impossible, our interest is directed towards the numerical

The critical mass flow rate (Moran, 2007, Zebbiche & Youbi, 2005a, 2005b) can be written in

\* 00 . 0 0 \*

As the mass flow rate through the throat is constant, we can calculate it at the throat. In this section, we have *ρ=ρ\*, a=a\**, *M=1*, *θ=0* and *A=A\*.* Therefore, the relation (20) is reduced to:

0 0 0 0

*<sup>m</sup><sup>ρ</sup> a A <sup>ρ</sup> <sup>a</sup> <sup>ρ</sup> <sup>a</sup>*

00 0 ( ) ( )

The parameters *T, P, ρ* and *A* for the perfect gas are connected explicitly with the Mach number, which is the basic variable for that model. For our model, the basic variable is the temperature because of the implicit equation (10) connecting *M* and *T*, where the reverse

*a T γ T a γ T T* 

*m adA <sup>M</sup> Aa a A* 

 

*A*

*\**

The determination of the velocity sound ratio is done by the relation (5). Thus,

 

 

*\* Exp <sup>A</sup>*

calculation. *All parameters M, ρ and A depend on the temperature.* 

 \**<sup>T</sup> T*

 

cos

*\* \**

 

1 2 1 2 

*/ /*

*A dTTF*

*<sup>A</sup>* (19)

(20)

(21)

(22)

The mass conservation equation is written as (Anderson, 1988 & Moran, 2007)

(9) and (12), one can receive the following equation:

2

*T*) gives the cross-section areas ratio:

analytical expression does not exist.

non-dimensional form:

Where

In the first case, one presents the table of variation of CP and γ versus the temperature for air (Peterson & Hill, 1965, Zebbiche 2010a, 2010b). The values are presented in the table 1.


Table 1. Variation of *CP(T)* and *γ(T)* versus the temperature for air.

For a perfect gas, the *γ* and *CP* values are equal to *γ*=1.402 and *CP=*1001.28932 J/(kgK) (Oosthuisen & Carscallen, 1997, Moran, 2007 & Zuker & Bilbarz, 2002).. The interpolation of the *CP* values according to the temperature is presented by relation (23) in the form of Horner scheme to minimize the mathematical operations number (Zebbiche, 2010a, 2010b):

$$C\_p(T) = a\_1 + T(a\_2 + T(a\_3 + T(a\_4 + T(a\_5 + T(a\_6 + T(a\_7 + T(a\_8 + T(a\_9 + T(a\_{10}))))))))) \tag{23}$$

The interpolation (*ai i=1, 2, …, 10*) of constants are illustrated in table 2.


Table 2. Coefficients of the polynomial *CP(T)*.

Effect of Stagnation Temperature on Supersonic

0.00

derivative at temperature *TS.*

Where

Fig. 2. Variation of the function *FA(T)* in the interval [*TS,T\**] versus *T0*

condensation function has the following form (Zebbiche & Youbi, 2005a):

1 1

1 *<sup>i</sup>*

*N*

*i i*

Obtained *si* values, enable to find the value of *Ti* in nodes *i:*

0.01

0.02

0.03

0.04

For *T T* <sup>0</sup> , <sup>0</sup> ( ) *HT C T T <sup>P</sup>* For *T T* <sup>0</sup> ,we have two cases:

if : ( ) relation (24) *T T HT*

if : ( ) ( ) ( ) *T T HT C T T HT <sup>P</sup>*

following form:

Flow Parameters with Application for Air in Nozzles 427

Taking into account the correction made to the function *CP(T)*, the function *H(T)* has the

The determination of the ratios (14) and (19) require the numerical integration of *Fρ(T)* and *FA(T)* in the intervals *[T, T0]* and *[T, T\*]* respectively. We carried out preliminary calculation of these functions (Figs. 1, 2) to see their variations and to choice the integration method.

> 0 100 200 300 400 500 *T (K)*

Due to high gradient at the left extremity of the interval, the integration with a constant step requires a very small step. The tracing of the functions is selected for *T0*=500 K (low temperature) and *MS*=6.00 (extreme supersonic) for a good representation in these ends. In this case, we obtain *T\**=418.34 K and *TS*=61.07 K*.* the two functions presents a very large

A Condensation of nodes is then necessary in the vicinity of *TS* for the two functions. The goal of this condensation is to calculate the value of integral with a high precision in a reduced time by minimizing the nodes number. The Simpson's integration method (Démidovitch & Maron, 1987 & Zebbiche & Youbi, 2006) was chosen. The chosen

 *bz s b z b* 

<sup>1</sup> <sup>1</sup>

*<sup>i</sup> z iN*

2

tanh (1 ) (1 ) 1 tanh( )

2

(26)

*b* 

*i*

(27)

A relationship (23) gives undulated dependence for temperature approximately low than T 240 K . So for this field, the table value (Peterson & Hill, 1965), was taken

$$
\overline{\mathbf{C}}\_p = \mathbf{C}\_p \left( \overline{T} \right) = 1001.15868 \text{ J } / \text{ (kg K)}.
$$

Thus:

for *T T* , we have ( ) *CT C P P*

for *T T* , relation (23) is used.

The selected interpolation gives an error less than *ε=10-3* between the table and interpolated values.

Once the interpolation is made, we determine the function *H(T)* of the relation (8), by integrating the function *CP(T)* in the interval [*T, T0*]*.* Then, *H(T) is a function with a parameter T0 and it is defined when T≤T0.* 

Substituting the relation (23) in (8) and writing the integration results in the form of Horner scheme, the following expression for enthalpy is obtained

$$\begin{aligned} H(T) &= H\_0 \left[ \mathbf{c}\_1 + T(\mathbf{c}\_2 + T(\mathbf{c}\_3 + T(\mathbf{c}\_4 + T(\mathbf{c}\_5 + Tc\_6 + T(\mathbf{c}\_7 + T(\mathbf{c}\_8 + T(\mathbf{c}\_9 + T(\mathbf{c}\_{10})))))))) \right] \right] \end{aligned} \tag{24}$$

Where

$$\begin{aligned} H\_0 &= T\_0 \{ \mathbf{c}\_1 + T\_0 \{ \mathbf{c}\_2 + T\_0 \{ \mathbf{c}\_3 + T\_0 \{ \mathbf{c}\_4 + T\_0 \{ \mathbf{c}\_5 + T\_0 \{ \mathbf{c}\_6 + T\_0 \{ \mathbf{c}\_9 \} \} \} \} \\ T\_0 &\{ \mathbf{c}\_7 + T\_0 \{ \mathbf{c}\_8 + T\_0 \{ \mathbf{c}\_9 + T\_0 \{ \mathbf{c}\_{10} \} \} \} \} \} \end{aligned} \tag{25}$$

and

$$\mathcal{C}\_{i} = \frac{a\_{i}}{i} \quad \text{ (\$i = 1\$, \$2\$, \$3\$, \$..., \$10\$)}$$

Fig. 1. Variation of function *Fρ(T)* in the interval [*TS,T0*] versus*T0*.

A relationship (23) gives undulated dependence for temperature approximately low

*C CT . ) P p* 1001 15868 J / (kg K

The selected interpolation gives an error less than *ε=10-3* between the table and interpolated

Once the interpolation is made, we determine the function *H(T)* of the relation (8), by integrating the function *CP(T)* in the interval [*T, T0*]*.* Then, *H(T) is a function with a parameter* 

Substituting the relation (23) in (8) and writing the integration results in the form of Horner

[ ( ( ( ( ( ( ( ( )))))))))]

0 01 02 03 04 05 06

*H Tc Tc Tc Tc Tc Tc*

((((((

( 1, 2, 3, ..., 10) *<sup>i</sup>*

0 100 200 300 400 500 *T (K)*

(24)

(25)

1 2 3 4 5 6 7 8 9 10

*c T c T c T c T c Tc T c T c T c T c*

( ( ( ( ))))))))))

*<sup>a</sup> c i <sup>i</sup>*

0 7 0 8 0 9 0 10

*Tc Tc Tc Tc*

*i*

than T 240 K . So for this field, the table value (Peterson & Hill, 1965), was taken

Thus:

values.

Where

and

for *T T* , we have ( ) *CT C P P* for *T T* , relation (23) is used.

0

0.00

Fig. 1. Variation of function *Fρ(T)* in the interval [*TS,T0*] versus*T0*.

0.01

0.02

0.03

0.04

()

*HT H -*

scheme, the following expression for enthalpy is obtained

*T0 and it is defined when T≤T0.* 

Taking into account the correction made to the function *CP(T)*, the function *H(T)* has the following form:

For *T T* <sup>0</sup> , <sup>0</sup> ( ) *HT C T T <sup>P</sup>* For *T T* <sup>0</sup> ,we have two cases: if : ( ) relation (24) *T T HT* if : ( ) ( ) ( ) *T T HT C T T HT <sup>P</sup>*

The determination of the ratios (14) and (19) require the numerical integration of *Fρ(T)* and *FA(T)* in the intervals *[T, T0]* and *[T, T\*]* respectively. We carried out preliminary calculation of these functions (Figs. 1, 2) to see their variations and to choice the integration method.

Fig. 2. Variation of the function *FA(T)* in the interval [*TS,T\**] versus *T0*

Due to high gradient at the left extremity of the interval, the integration with a constant step requires a very small step. The tracing of the functions is selected for *T0*=500 K (low temperature) and *MS*=6.00 (extreme supersonic) for a good representation in these ends. In this case, we obtain *T\**=418.34 K and *TS*=61.07 K*.* the two functions presents a very large derivative at temperature *TS.*

A Condensation of nodes is then necessary in the vicinity of *TS* for the two functions. The goal of this condensation is to calculate the value of integral with a high precision in a reduced time by minimizing the nodes number. The Simpson's integration method (Démidovitch & Maron, 1987 & Zebbiche & Youbi, 2006) was chosen. The chosen condensation function has the following form (Zebbiche & Youbi, 2005a):

$$s\_i = b\_1 \ z\_i + (1 - b\_1) \left[ 1 - \frac{\tanh\left[ \begin{array}{c} b\_2 \cdot (1 - z\_i) \end{array} \right]}{\tanh(b\_2)} \right] \tag{26}$$

Where

$$z\_i = \frac{i-1}{N-1} \qquad 1 \le i \le N \tag{27}$$

Obtained *si* values, enable to find the value of *Ti* in nodes *i:*

Effect of Stagnation Temperature on Supersonic

**3.3 Supersonic nozzle conception**

possible mass of structure.

**3.4 Error of perfect gas model**

the following relationship:

prefer the first one.

(Peterson & Hill, 1965 & Zebbiche, Youbi, 2005b)

*P=PS*, *a=aS* and *A=AS.*

*N, b1* and *b2*.

Where

Flow Parameters with Application for Air in Nozzles 429

The ratios of pressures, speed of sound and the sections corresponding to *M=MS* can be calculated respectively by using the relations (15), (22) and (19) by replacing *T=TS*, *ρ=ρS,* 

The integration results of the ratios *ρ\* / ρ0, ρS/ρ0* and *AS/A\** primarily depend on the values of

For supersonic nozzle application, it is necessary to determine the thrust coefficient. For nozzles giving a uniform and parallel flow at the exit section, the thrust coefficient is

> *F <sup>F</sup> <sup>C</sup>*

The introduction of relations (21), (22) into (32) gives as the following relation:

0

*F E*

0

for the aerodynamic applications, the error should be lower than *5%*.

*y*

0 \*

\* \*

*<sup>E</sup>*

*a a C TM*

The design of the nozzle is made on the basis of its application. For rockets and missiles applications, the design is made to obtain nozzles having largest possible exit Mach number, which gives largest thrust coefficient, and smallest possible length, which give smallest

For the application of blowers, we make the design on the basis to obtain the smallest possible temperature at the exit section, to not to destroy the measuring instruments, and to save the ambient conditions. Another condition requested is to have possible largest ray of the exit section for the site of instruments. Between the two possibilities of construction, we

The mathematical perfect gas model is developed on the basis to regarding the specific heat *CP* and ratio *γ* as constants, which gives acceptable results for low temperature. According to this study, we can notice a difference on the given results between the perfect gas model and developed here model.The error given by the *PG* model compared to our *HT* model can be calculated for each parameter. Then, for each value *(T0, M)*, the *ε* error can be evaluated by

0

( ) ( )1 <sup>100</sup>

The letter *y* in the expression (35) can represent all above-mentioned parameters. As a rule

*y T , M <sup>ε</sup> T , M y T , M*

*PG*

*HT*

0

( )

(35)

0 \*0

*a a* 

*P A* (32)

*F mV mM a E EE* (33)

(34)

$$T\_i = \mathbf{s}\_i \left( T\_D - T\_G \right) + T\_G \tag{28}$$

The temperature *TD* is equal to *T0* for *Fρ(T)*, and equal to *T\** for *FA(T).* The temperature *TG* is equal to *T\** for the critical parameter, and equal to *TS* for the supersonic parameter. Taking a value b*<sup>1</sup>* near zero (b*1*=0.1*, for example*) and b*2*=2.0, it can condense the nodes towards left edge *TS* of the interval, see figure 3.


Fig. 3. Presentation of the condensation of nodes

### **3.1 Critical parameters**

The stagnation state is given by *M*=0. Then, the critical parameters correspond to *M*=1.00, for example at the throat of a supersonic nozzle, summarize by:

When *M*=1.00 we have *T=T\**. These conditions in the relation (10), we obtain:

$$2\,\,{\Im}\,{H}(T\_{\ast}) - a^{2}(T\_{\ast}) = 0\tag{29}$$

The resolution of equation (29) is made by the use of the dichotomy algorithm (Démidovitch & Maron, 1987 & Zebbiche & Youbi, 2006), with *T\*<T0*. It can choose the interval *[T1,T2]*  containing *T\** by *T1=0 K* and *T2=T0*. The value *T\** can be given with a precision *ε* if the interval of subdivision number *K* is satisfied by the following condition:

$$K = 1.4426 \,\mathrm{Log} \left( \frac{T\_0}{\varepsilon} \right) + 1 \tag{30}$$

If *ε*=10-8 is taken, the number *K* cannot exceed 39. Consequently, the temperature ratio *T\*/T0* can be calculated.

Taking *T=T\** and *ρ=ρ\** in the relation (14) and integrating the function *Fρ(T)* by using the Simpson's formula with condensation of nodes towards the left end, the critical density ratio is obtained.

The critical ratios of the pressures and the sound velocity can be calculated by using the relations (15) and (22) respectively, by replacing *T=T\*, ρ=ρ\*, P=P\** and *a=a\*,* 

### **3.2 Parameters for a supersonic Mach number**

For a given supersonic cross-section, the parameters *ρ=ρS, P=PS, A=AS,* and *T=TS* can be determined according to the Mach number *M=MS*. Replacing *T=TS* and *M=MS* in relation (10) gives

$$2\, H \text{(T}\_{\text{S}}\text{)}-M\_{\text{S}}^{2}\text{ }a^{2} \text{(T}\_{\text{S}}\text{)}=\text{0}\tag{31}$$

The determination of *TS* of equation (31) is done always by the dichotomy algorithm, excepting *TS<T\*.* We can take the interval *[T1,T2]* containing *TS*, by (*T1*=0 *K,* and *T2=T\*.* 

Replacing *T=TS* and *ρ=ρ<sup>S</sup>* in relation (14) and integrating the function *Fρ(T)* by using the Simpson's method with condensation of nodes towards the left end, the density ratio can be obtained.

The ratios of pressures, speed of sound and the sections corresponding to *M=MS* can be calculated respectively by using the relations (15), (22) and (19) by replacing *T=TS*, *ρ=ρS, P=PS*, *a=aS* and *A=AS.*

The integration results of the ratios *ρ\* / ρ0, ρS/ρ0* and *AS/A\** primarily depend on the values of *N, b1* and *b2*.
