Category UNSTEADY AERODYNAMICS, AEROACOUSTICS AND AEROELASTICITY OF TURBOMACHINES

The modeled and measured compressor is a one-stage radial compressor. The compressor is driven directly by a 175 kW high-speed electric motor that is supplied with variable frequency AC-current by a 200 kW inverter. No gear­box is needed because of the single shaft construction, and the use of an in­verter allows variable speed control from zero RPM to full speed. Active mag­netic bearings are used in the rotor. The compressor has a mass flow rate of 1.574 kg/s, pressure ratio 2.01, and rotational speed 21640 RPM at the design operation point, and it produces a 50-kPa vacuum at the inlet.

The impeller has 7 full blades and 7 splitter blades with 30 ° backsweep at the trailing edge measured from the radial direction. The diameter of the
impeller is 309 mm. The diffuser is vaneless and the overhung volute has a circular cross-section. The area of the cross-section of the volute is designed for a constant static pressure along the circumference at the design operation point. More detailed information about the test facilities and the simulated compressor are presented in Ref. [6].

Teemu Turunen-Saaresti

Lappeenranta University of Technology Laboratory of Fluid Dynamics Teemu. Turunen@lut. fi

Jaakko Larjola

Lappeenranta University of Technology Laboratory ofFluid Dynamics Jaakko. Larjola@lut. fi

Abstract An unsteady experimental and numerical investigation was conducted to a cen­trifugal compressor with a vaneless diffuser. The unsteady static pressure was measured at the diffuser inlet and outlet at three different circumferential po­sitions. Also a time-accurate numerical simulation was carried out. The whole compressor was modeled. Turbulence was modeled with the low Reynolds num­ber k — e turbulence model. The measurements were made in operation points near the surge, near the choke and at the design operation point of the compres­sor. The numerical simulations were conducted at the design operation point and at the operation point near the choke. The gathered numerical and experimental unsteady pressures were compared with each other.

Keywords: Unsteady pressure, centrifugal compressor, vaneless diffuser

1. Introduction

The ft>w field inside a radial compressor is highly three-dimensional and unsteady. The wakes from the blades of the impeller and the non-axisymmetry of the volute causes a circumferentially varying ft>w field in the diffuser. The ft>w field in the vaneless diffuser is examined by measuring unsteady pres­sure at the inlet and outlet of the diffuser and analyzing the whole compressor numerically.

493

K. C. Hall et al. (eds.),

Unsteady Aerodynamics, Aeroacoustics and Aeroelasticity of Turbomachines, 493-503. © 2006 Springer. Printed in the Netherlands.

The simultaneous solution of the three-dimensional unsteady Navier-Stokes equations in the impeller and volute requires a large amount of a computa­tional resources. Previously Fatsis et al. [1] have carried out three dimensional unsteady fbw calculation of the impeller using Euler solver. Hillewaert and Van den Braembussche [2] have carried out the same simulation with the in­clusion of a steady volute calculation. Koumoutsos et al. [3] have conducted viscous three-dimensional unsteady flow calculation to the centrifugal com­pressor stage with vaned diffuser but the volute was not modeled.

In this paper, the time-accurate simulation was conducted to the whole com­pressor and the unsteady static pressure was measured at the diffuser inlet and outlet at three different circumferential positions. In the time-accurate simula­tion, the computational mesh describing the rotational part of the compressor was rotated. The connection between the stationary and the rotating part of the mesh was handled by using the sliding mesh technique. The grid lines between the impeller blocks and the stator block were discontinuous, so a mass con­serving interpolation was made at every time step [4]. At first, a time-averaged solution was calculated for the starting point of the time-accurate simulation. Then, using dual time-stepping, in which a pseudo-time integration was per­formed inside every physical time step, the time-accurate calculation was per­formed. This way the time history could be taken into account. Turbulence was modeled with the Chien’s k — є model [5]. A parallel computation with six pro­cessors was used to achieve a tolerable computational time. The convergence was speeded up with the multigrid method.

The unsteady static pressure measurements were made with Kulite XTC – 190 pressure transmitter and the data was collected with a Tektronix oscillo­scope connected to the computer. The unsteady pressure was measured within the time of ten rotations of the impeller. The measuring frequency was 1 MHz. Three sets of data were collected at each measuring point. The exact rotation speed of the rotor was measured at the same time.

Figure 1 shows the pressure contours of the steady state cascade ft>w (Ap = 333Pa). This ft>w field corresponds to the case of the back pressure of 75,000 Pa. The inlet Mach number is 1.2 and the exit Mach number is 0.605. The inlet ft>w angle is 57.6 deg., and the exit ft>w angle is 52.35 deg. At the leading edge we have oblique shocks, and in the cascade passage a normal shock. Interaction between a sonic boom and cascaded blades is studied based on this steady condition. In order to save computational region as well as time we chose the incident N wave of a short period of 8 ms. In reality the sonic boom induced by SST may have a longer period of 200 ms, we thought the interaction between the shock and cascaded blades should be essential, but the period of two shocks may not be so important in the phenomena. This tendency has been suggested by the semi-actuator disk analysis.

The maximum static pressure amplitude of the N wave given at the inlet boundary is selected to be 325 Pa. This value corresponds to the increase of relative total pressure of 1013 Pa, 1% of the atmospheric pressure. In this case the maximum amplitude of the absolute total pressure is 1.19%, and the static pressure is 0.78% of the steady value, respectively.

Hereafter we investigate the pressure variation history along the grid points indicated in Fig.2. The points from 1 to 15 locate in the upstream field, those from 15 to 25 locate in the cascade passage, and points from 25 and greater locate in the downstream field.

Figures 3 (a) to (e) show the results of pressure variation history for the inci­dent N-wave of period of 8ms. From the pressure history at the position 1, we can clearly recognize the incident N-wave followed by a little bit smaller N – wave refected back by the cascade at the leading edge plane. At positions 4, 6 and 8 the refected wave complicates the pressure signal significantly. Figures 3 (c) and 3 (d) show the pressure signals in front of and behind the passage shock, respectively. We see that high frequency components appear on the pressure signals behind the passage shock. Those components are probably due to var­ious types of resonance in the cascade passage. Further we go downstream the pressure signals seem to recover the N-shape as the transmitted wave.

Figure 4 shows the history of the non-steady aerodynamic force acting on blades during the passage of the N-wave. We see that at the passage of the front shock of the N-wave, the aerodynamic force once drops sharply then increases to have a peak. During the passage of the expansion portion of the N-wave, the

aerodynamic force increases gradually. At the passage of the end shock of the N-wave, the sharp variation of the aerodynamic force repeats again. The sharp variation can be explained as follows: First, the incident shock wave hits upon the suction side of blades reducing the aerodynamic force. Then, the shock wave refected at the suction side propagates downstream and hits upon the pressure side of blades increasing the aerodynamic force.

Figures 5 (a) to (d) show the blade surface pressure distribution along chord as measured by the perturbation from the steady state condition. We focus on the short period of the passage of the front shock, i. e., from 2.4ms to 3.0ms. The non-steady aerodynamic force experiences a sharp variation during this period. From 2.4ms to 2.6ms the suction side surface pressure increases be­cause of the incident shock wave. The pressure side surface pressure does not change during this time period. Then at 2.8ms the high pressure wave is diffracted around the trailing edge thus increasing the pressure side surface pressure near the trailing edge. In addition to this effect, at 3.0ms and after the pressure wave refected at the suction side of a blade arrives at the pressure

side of the adjacent blade increasing the whole level of the pressure side sur­face pressure. These behaviors of surface pressure clearly explain the variation of the non-steady aerodynamic force acting on blades.

Figure 5. Surface Pressure Distribution

The maximum amplitude of the non-steady aerodynamic force reaches 1.51% of the steady state aerodynamic force for the incident N-wave having the am­plitude of 1% of the steady state static pressure. Thus we can conclude that the amplitude of the non-steady aerodynamic force is of the same order of magnitude as the pressure amplitude of the incident N-wave.

2. Conclusions

Interaction between a sonic boom and a transonic cascade of blades is stud­ied numerically. The main results obtained can be summarized as follows.

For an incident N wave, the refbcted N wave and the transmitted N wave are clearly observed.

After passing through the passage shock the increase in the N wave am­plitude which has been suggested by the semi-actuator disk analysis is not clearly observed.

By the incidence of a sonic boom the aerodynamic force once drops sharply then rises to attain the peak pressure. The mechanism is the process such that first the incident shock wave hits upon the suction side of a blade then be refected to hit upon the pressure side of the adjacent blade.

The maximum amplitude of the unsteady aerodynamic force reaches 1.51% of the steady state aerodynamic force for the incident N wave having the amplitude of 1% of the steady state static pressure.

References

Kaji, S., “Interaction of Sonic Boom and Engine Fan Rotor”, Proc. 9th ISUAAAT, Sep. 2000, pp.362-374.

Paynter, G. C., Clark, L. T., and Cole, G. L., “Modeling the Response from a Cascade to an Upstream Acoustic Disturbance”, AIAA J., Vol.38, No.8, Aug. 2000.

Freund, D. and Sajben, M., “Experimental Investigation of Outfbw Boundary Conditions Used in Unsteady Inlet Flow Computations”, AIAA paper 97-0610.

Dorney, D. J., “Unsteady Acoustic Wave Propagation in a Transonic Compressor Cas­cade”, AIAA paper 96-0246.

“Study on Advanced Aircraft Technology Development” Report No. 1204, SJAC, 2001.

A sonic boom is produced by a supersonic airplane because of its volume and aerodynamic forces. Its waveform is very much complicated near the air­plane, however, in the distant field the waves are naturally weakened and coa­lesced to be a so-called N-wave, i. e., the combination of a bow shock and an end shock. The detailed waveform and strength will be dependent, therefore, on the difference of fight altitudes of the SST and the subsonic transport in
question. Numerical simulation of the sonic boom of Concorde at various alti­tudes [5] indicates that at 40,000 ft the maximum amplitude of the sonic boom reaches as high as 100 Pa. This value corresponds to about 0.5strength should be increased as the speed and size of airplane increases. For example, as we go from Mach 2 to Mach 3 the sonic boom strength increases by 15 next genera­tion SST of 300 ton class, it increases by about 20 Therefore, the strength of a sonic boom encountered at altitudes of subsonic transport should be 1

In the present analysis a numerical approach is taken based on the two­dimensional Euler equation.

In the actual analysis, the basic equations are transformed into the generalized coordinate system £ = £(x, y),n = n(x, y).

To capture shock waves the upwind type TVD scheme by Harten-Yee is adopted. For time integration, the LU-ADI scheme is used. The Newton iter­ation method is applied to the non-steady calculation in order to ensure time accuracy. For the boundary conditions in case of steady calculation, the total pressure, the total temperature and the circumferential velocity are specified at the inlet, and the static pressure is fixed at the exit.

For the non-steady calculation, the nonrefecting boundary condition of the Thompson type is applied. For the incident N-wave, the physical conditions are specified at the inlet boundary. The finite amplitude, one-dimensional isen – tropic wave can be described by

where u is the velocity, a the velocity of sound, p the pressure, and p the density, and the subscript indicates the steady state condition before the wave arrives. The plus sign indicates the progressing wave and the minus sign rep­resents the retreating wave.

From the above equation, once the pressure is given, the density and the induced velocity are obtained. During the period of the incident N-wave, the total pressure and the total temperature are not fixed but the static pressure, density and velocity are given in accordance with the above relation of the N-wave. After the entrance of the N-wave, the total pressure and the total temperature are fixed again at the inlet. In order to simulate a realistic sonic boom, we introduce a start up time(or rising time) of 0.1ms for the shock wave to gain the peak pressure.

The combination of grids of H-type in the upstream and downstream fields and those of O-type near the region of cascaded blades is adopted. The blades have the tip section profile used in the NASA Quiet Engine Fan Program B. At the design condition the inlet Mach number is 1.2 and the exit Mach number is 0.853.

Shojiro Kaji

Department of Aerospace Engineering Teikyo University

1-1 Toyosato-dai, Utsunomiya 320-8551 Japan

kaji@koala. mse. teikyo-u. ac. jp

Takahiro Suzuki

Department ofAeronautics and Astronautics University of Tokyo

7-3-1 Hongo, Bunkyo-ku, Tokyo 113-8656 Japan

Toshinori Watanabe

Department ofAeronautics and Astronautics University of Tokyo

7-3-1 Hongo, Bunkyo-ku, Tokyo 113-8656 Japan

watanabe@aero. t.u-tokyo. ac. jp

Abstract Interaction between a sonic boom and cascaded blades is studied numerically.

The incident sonic boom is modeled by a finite amplitude N wave with a short rising time, and the time dependent transonic cascade fbw with a passage shock inside is solved. Detailed mechanism on refection and transmission of a sonic boom by the cascade and unsteady aerodynamic force exerted on blades are elucidated. The peak amplitude of the unsteady blade force is the same order of magnitude as the pressure amplitude of the incident sonic boom.

Keywords: sonic boom, cascaded blades, interaction, supersonic transport,

numerical analysis

483

K. C. Hall et al. (eds.),

Unsteady Aerodynamics, Aeroacoustics and Aeroelasticity of Turbomachines, 483-491. © 2006 Springer. Printed in the Netherlands.

1. Introduction

In this century, a new era of supersonic transport (SST) is in prospect where a great number of SST is introduced concurrently with subsonic transport. Then, there will be an increasing possibility for subsonic transport of encoun­tering a sonic boom of supersonic transport. In such an occasion shock waves due to a sonic boom of SST are swallowed in an engine inlet and interact with fan rotor blades and stator vanes. In another occasion a sonic boom may invade into the engine cowl from the rear side.

Another example of interaction of shock waves and cascaded blades can be found in an innovative concept of a pulse detonation engine (PDE) combustor combined with a turbine engine. In this case shock waves due to a detonation combustor hit upon turbine nozzle vanes and blades. The strength of shock waves for a PDE combustor reaches more than 200 dB in pressure amplitude. Again we have a strong motivation to investigate the interaction between shock waves and cascaded blades.

The purpose of the present paper is to investigate numerically the impact of a sonic boom produced by SST on the engine fan rotor of subsonic transport. The same problem was treated in the last symposium analytically by use of the semi-actuator disk model for cascaded blades [1]. The results obtained in the analysis showed that (1) the amplitude of the transmitted N wave across the cascade passage shock increases over the incident N wave, (2) the pressure fhctuation at the cascade passage shock has a discontinuous character, and (3) unsteady aerodynamic force exerted on blades by a sonic boom becomes strong in the direction of chord. The numerical approach taken in this paper enables to check the validity of the semi-actuator disk model and give more detailed interaction mechanism.

Similar analyses have been made from different view points. Paynter et al. [2] needed the rear end boundary condition for numerical analysis of an en­gine air intake. They solved the fbw field including the cascade to determine the refection coefficient of sound waves. Freund et al. [3] investigated exper­imentally the refected pressure fhcthation for the disturbance incident upon an engine compressor. Dorney [4] analyzed numerically the effect of sound waves incident upon a viscous transonic cascade.

3.2.1 Numerical-Experimental comparison on unsteady pressure

distribution. The amplitude (normalized by outlet value) and phase angle of the unsteady pressure distribution is plotted in figure 5 for a perturbation fre­quency of 100Hz. A pressure amplification of factor three can be observed downstream of the shock location for both experimental and RANS calcula-

tion results. It is interesting to note that this amplification is not observed in the Euler simulation and might thus originate viscous or turbulent effects, or even possibly the Shock Boundary Layer Interaction (SBLI).

The analysis of the phase angle distribution is facilitated by considering the behaviour of travelling pressure waves in duct fbws. Similarly to poten­tial interaction in turbomachines, outlet static pressure fhctuations propagate upstream at a relative velocity of |c-U|. As long as the propagating speed is unchanged, the slope of the phase angle also remains constant, which is the case in the outflow region. However, in the vicinity downstream of the shock, the phase angle stops decreasing and even increases, which would actually cor­respond to a downstream propagating pressure wave.

At higher perturbation frequency (Fp=500Hz on figure 6), the pressure am­plification for both experimental and 2D RANS simulation exhibits an atten­uation downstream of the shock whereas the phase angle distribution presents an important phase shift (about 160o) at the same location (x=95mm). Further­more the same "increasing phase angle” behaviour is still observed downstream of the shock. It is noteworthy that the Euler simulation differs both in the am­plitude and phase distribution and does not present the same characteristics.

Generally, a fairly good agreement was achieved between experiments and turbulent computations. However, although the phase angle distribution for experimental measurements and 2D RANS simulation are very similar, there seems to be an "offset" between the two distributions. This may be explained by the location of the experimental reference on the side wall. As the pressure perturbations are generated by the rotating ellipse and propagate upstream in the test section, there might be a possibility of a phase shift between the lower and side wall pressure measurements as the ft>w is non uniform in the diffusor downstream of the shock.

Considering the above observation on phase shift (increasing phase angle in the vicinity of the shock) it is believed that the unsteady pressure distribu­tion at this location results from the superposition of an upstream and down­stream travelling wave with respective varying amplitudes and phase angles. Assuming that only plane waves can propagate in the nozzle at 100Hz, a one dimensional acoustic decomposition (see equation 1) was performed on the 2D RANS numerical results using the isentropic flow velocity and sound speed on the bump surface.

respective amplitude 1.7 and 0.8. The relative phase angle distribution is more delicate to understand. For upstream propagating wave, the phase distribution tends to increase in the vicinity of the shock. This might be explained by the fact that the relative propagating velocity |c-U| tends toward zero in the vicin­ity of the shock. The wave length thus also tends toward zero locally, which is interpreted as an increase of the phase in the decomposition algorithm. Con­cerning the downstream propagating waves, it is not yet clear why the phase distribution is decreasing between x=80mm and x=210mm. It is assumed that the decomposition model is not fully accurate in the region and that other phe­nomena are not yet taken into account.

1.5.2 Unsteady pressure distribution at different channel heights.

Considering the fairly good agreement between experiments and 2D RANS simulations, a complementary analysis of the unsteady pressure distribution was conducted within the numerical domain at different positions in the chan­nel’s height. In addition to the upper and lower walls, three other locations were selected as illustrated in figure 2(c): the SBLI region (probe C), the mid­dle of the shock (probe B), and the region where the sonic line meets the shock (probe A).

At low perturbation frequency (100Hz), the highest pressure amplification occurs in the BL on the upper and lower walls. Although not located in the BL, a slight amplification occurs downstream of the shock at the location of probe A and C whereas the pressure amplification in the middle of the channel (probe B) is very similar to the one observed in Euler computation. Beside, the phase angle distribution exhibits the behaviour described previously con­cerning the increasing phase angle downstream of the shock exclusively for the location where a shock occurs. It is noteworthy that this phenomenon is most pronounced in the middle of the channel although no pressure amplification was noticed. It thus seems as the observed phenomenon on the phase angle is

rather due to the presence of the shock wave than to viscous or turbulent effects.

At higher perturbation frequency (500Hz), a pressure attenuation can be ob­served close to the bump surface whereas a slight amplification still occurs at other location in the channel. Surprisingly, the increasing phase angle be­haviour only occurs on the bump surface whereas the phase angle distribution refects a single upstream propagating wave behaviour at all other locations.

The particular amplitude and phase angle distribution is not yet clearly un­derstood and might be related to the lower amplitude of the motion of the shock, its inertia to oscillate, or the smaller wave length of the perturbation at higher frequencies.

Shock motion. The amplitude and phase angle of the shock mo­tion throughout the channel’s height are presented in figure 10 and 11 respec­tively for the perturbation frequencies 100Hz and 500Hz. A fairly good agree­ment is achieved between experiments and numerics regarding both the ampli­tude and phase distribution of the unsteady shock motion. For both frequen­cies, the amplitude of shock motion increases with the height of the channel. It is noteworthy that the same trend is observed also in the Euler simulation. The amplitude of motion of the shock is thus related to the mean fbw gradients rather than to the SBLI. Beside, the amplitude clearly decreases with the fre­quency. Assuming that the shock has a certain inertia to oscillate or to respond to a back pressure variation, the decreasing amplitude of motion of the shock at higher frequencies is then probably due to shorter perturbation wave length although the amplitude of perturbation remains the same.

On the other hand, the shock motion phase angle slightly differs both be­tween the two frequencies and the two numerical calculations. At low fre­quency, the phase angle seems to be constant throughout the channel height, which corresponds to a "rigid" motion between the foot and the top of the shock. At higher frequency however, there is an important phase lag (up to

90o for experiments) between the foot and the top of the shock. The "head" of the shock seems to oscillate with a certain time lag (delay) compared to the foot. This tendency can be observed both on experimental visualization and 2D RANS simulation, but not on Euler calculation where the phase angle is still constant throughout the channel height. This effect might be related to the higher propagation speed within the BL than in the free stream. As a result, the pressure perturbations reach the foot of the shock slightly in advance in the SBLI region. The same behaviour should theoretically also be seen at lower frequencies, but since the wave propagation speed is the same, the perturba­tions wave length is longer and thus the phase lag between perturbations in the BL and in the free stream is lower.

Figure 11. DFSD on shock motion over 2D bump for Fp =500Hz

1.5.3 Relation between unsteady pressure and shock motion. The

phase lag between unsteady pressure distribution over the bump and the shock motion (closest from the bump) has been calculated and plotted in figure 12 as a function of the perturbation frequency. Experimental results show an increas­ing phase lag with the perturbation frequency, meaning that the shock response with a certain increasing time delay with the perturbation frequency. For nu­merical simulation however, the trend is not as clear and further calculations
should be performed in order to clearly define a tendency.

Additionally, the unsteady pressure phase angle jump through underneath the shock location was plotted as a function of the perturbation frequency in figure 13. A fairly good agreement was found between experimental pressure measurements and 2D RANS numerical simulation. This unsteady pressure phase shift, which seems to increase linearly with the perturbation frequency, is extremely important considering aeroelastic stability prediction. Indeed, the unsteady aerodynamic load on an airfoil is directly infhenced by the value of the phase shift and the overall stability of the airfoil might change from stable to unstable (and vice versa) for a certain value of this phase shift.

3.2.5 Unsteady separated region motion. The unsteady motion of the separated zone was evaluated as a function of the shock motion for the 2D RANS numerical simulation, and is presented in figure 14 for two perturba­tion frequencies (100Hz and 500Hz). The advantage to display the separation versus the shock motion is to be able to see both the amplitude of motion and the phase lag between the separation/reattachment and the shock oscillations. Clearly, the separation oscillates with the same amplitude and nearly in phase with the shock wave. This behaviour does not seem to change with the fre­quency and confirms the idea of a shock induced separation. On the other hand, the reattachment seems to oscillate with a much larger amplitude and a certain phase lag with the shock. Both the amplitude and phase lag seem to be related to the perturbation frequency. It is however not possible to state upon a clear tendency and further calculation as well as experimental measurement should be performed.

2. Summary

Unsteady pressure measurements and high speed Schlieren visualizations were conducted together with 2D RANS and Euler numerical simulations over a convergent divergent nozzle geometry in order to investigate the Shock Bound­ary Layer interaction. A fairly good agreement between experiments and nu­merical simulations was obtained both regarding the unsteady pressure distri­bution and shock motion.

Results showed that the unsteady pressure distributions, both on the bump and within the channel, result from the superposition of upstream and down­stream propagating waves. It is believed that outlet pressure perturbations propagate upstream within the the nozzle, interact in the high subsonic ft>w region according to the acoustic blockage theory, and are partly refected or absorbed by the oscillating shock, depending on the frequency of the pertur­bations. The amplitude of motion of the shock was found to be related to the mean flow gradients and the local wave length of the perturbations rather than to the shock boundary layer interaction. The phase angle between unsteady pressure distribution on the bump and the shock motion for experimental re­sults was found to increase with the perturbation frequency. However, no clear tendency could be defined for numerical results.

At last, but not least, the phase angle "jump" underneath the shock location was found to linearly increase with the perturbation frequency. The phase shift is critical regarding aeroelastic stability since it might have a significant im­pact on the phase angle of the overall aerodynamic force acting on the blade and shift the aerodynamic damping from stable to exciting.

1.5 Steady state results

The steady state shock wave in the 2D nozzle is presented in figure 3.1 for experimental visualization and viscous numerical simulations (at mid channel plane for the 3D RANS simulation). Although a fairly good agreement on the shock location and structure is achieved, the experimental shock position could only be matched by raising up the outlet static pressure value in the numerical simulations. For the same back pressure value, the simulations would position the shock more downstream in the diffusor. Although there might be a real probability that the k — и turbulent model underestimates the level of losses, it cannot, by itself, explain the differences in the shock location. A more prob­able explanation involves the thickening of the side wall BLs and the resulting change in the effective section area, which would act like a slight convergent and magnify the pressure gradient. As the back pressure is manually setup, the pressure right downstream of the shock is then higher and the shock moves upstream.

Using the continuity equation at the outlet (Qm = p2VX2 S2), the reduction of section area due to BL thickening can be estimated for numerical simula­tions by calculating the change of section necessary to obtain the experimental mass fbw under the same numerical outlet conditions. For the 2D RANS sim­ulation, in which no side wall BL is specified, the change of section area was estimated around 9.44cm2, equivalent to a BL with a displacement thickness of 3.9mm on each side wall. For the 3D RANS, which already features outlet side wall BLs, the change of section area was estimated around 7.63cm2 and is equivalent to an increase of the displacement thickness of 1.73mm on each side wall BL.

The steady state pressure distribution at mid channel (y=50mm) over the 2D bump surface is plotted in figure 3 for experimental results and viscous numer­ical simulations. Although the curves collapse fairly well regarding the shock location, they differ downstream of it. Indeed, experimental results present a smoother pressure recovery, which denotes a change of local curvatures (to­wards a more convex surface) usually due to a separated flow region. This phenomenon is even stronger closer to the wall (see pressure distribution at y=10mm) and denotes a large BL thickening or a separation of the fbw in the corners. Probably due to larger side wall BLs and the interaction with the shock, the pressure rise occurs more upstream in the region close to the side walls.

In figure 4 are presented the streamlines from experimental visualization and 3D RANS calculation. As mentioned previously, the side wall BLs start thickening right downstream of the shock and a 15 mm large and 60mm long separation appears in both corners. In comparison, the 3D RANS predicts a

Table 4. Separated region location

much lower separated region, both in the corner and at mid channel. Again, this can be an effect of mismatched inlet boundaries or an important underes­timation of the losses by the turbulent model. The size of the separated region, measured at mid channel, is presented in table 4. Whereas the 2D RANS calcu­lation presents a fairly good estimation of the separated region, it is noteworthy that the 3D RANS simulation actually gives a much worth prediction. A possi­ble reason might simply be the underestimation of the side wall BL thickening.

The shock motion was imposed by sinusoidal downstream static pressure plane fhctuations. The amplitude and frequency of the perturbations are sum­marized in table 2 and the corresponding reduced frequency are presented in table 3. The reduced frequency based on the BL thickness is considered small enough to justify a quasi-steady response of the turbulence that is compatible with the turbulent model used.

The data reduction consisted in performing an harmonic analysis on both the unsteady pressure distribution and the shock motion. Similarly to exper­iments, the amplitude of each harmonic was divided by the amplitude of the fundamental at the outlet, and the phase angle value at the outlet reference lo­cation was subtracted to all signals, for each harmonic respectively.

The computed configuration is the experimental 2D nozzle previously de­scribed. The numerical domain was however extended 70mm upstream and 164mm downstream of the bump in order to avoid numerical interaction with the boundaries. Steady state simulations were performed on both 2D and 3D structured H-grids in order to achieve a good understanding of the mean ft>w structures. The respectively meshes comprise 150×84 and 150x84x35 nodes with adapted grid density both around the shock location and in upper, lower and side wall BLs (containing respectively 33, 28 and 26 nodes). Unsteady simulations were however only performed on the 2D mesh due to computation time restriction. Both RANS and Euler unsteady computations were performed for comparison purposes.

1.4 Steady flow conditions

For RANS computations, the fhid is modelled as a viscous perfect gas. The specific heat ratio equals k=1.4 and the perfect gas constant is R=287 J/kg/K. The laminar dynamic viscosity and the thermal conductibility are assumed con­stant and respectively equal ^=1.81 10-5 kg/m/s and k=2.54 10-2 m. kg/K/s3. The inlet conditions in the free stream are such that the stagnation pressure, pmiet, and the stagnation temperature, Ttmlet, equal respectively 160kPa and 303K. A fully developed 7mm thick BL profile computed over a flat duct is specified as inlet condition. The outlet static pressure was adjusted in order to match the experimental shock configuration. The numerical operating condi­tions are summarized in table 2.

1.3 CFD tool

Simulations were performed using the computational model referenced as PROUST [12] and developed to simulate steady and unsteady, viscous and inviscid fbws. The fully three-dimensional unsteady, compressible, RANS equations are solved. The space discretization is based on a MUSCL finite vol­ume formulation. The convective flixes are evaluated using an upwind scheme based on Roe’s approximate Riemann solver, and the viscous terms are com­puted by a second order centered scheme. The turbulence closure problem is solved using Wilcox k-w two equations model and fully accounts for the ef­fect of the boundary layer (BL) separation which originates at the shock foot location. Compatibility relations are used to account for physical boundary conditions. One-dimensional numerical boundary conditions are implemented by retaining the equations associated to the incoming characteristics and fixing the wave velocity to zero to prohibit propagation directed into the computa-

tional domain. The resulting semi discrete scheme is integrated in time using an explicit five steps Runge-Kutta time marching algorithm.