Numerical domain and grid

The numerical domain and the grid used for this study are shown in Figure 1. The numerical domain is based on an isolated low-speed axial compressor (LSAC) rotor. The rotor blade of the LSAC is obtained through an inverse-design method using MISES (Youngren and Drela, 2008). The inverse design method is an iterative method used to obtain a desired blade shape from a prescribed blade loading. The aim of the inverse-design method is to obtain a tip-critical compressor rotor that is suitable for passive casing treatment study. Hence, the prescribed blade loading is obtained numerically from an axial compressor that is known to exhibit a spike-stall inception pattern. The 3D rotor blade is generated by radially stacking seven 2D profiles at as many spanwise positions. The hub and casing are parallel to one another resulting in a constant annulus height across the blade. The 3D blade design has a nearly constant exit flow angle (absolute) distribution in the spanwise direction. This agrees with the flow angle distribution of the tip critical compressor of which the inverse design is attempted. Aerodynamic design parameters of the LSAC are presented in Table 1.

Figure 1.

Isometric view of the numerical domain and grid of the LSAC rotor.

The single passage numerical domain and unstructured grid are generated automatically using ANSYS Turbogrid 17.1. The grids are generated by setting the target passage grid node count and the boundary layer grid node count near the walls. A grid independence study was conducted by varying the total number of elements, but with a fixed y+ value (≈2). The change in the converged mass flow rate value and the blade surface static pressure distribution at 98% span was found to be negligible for grid resolutions that are finer. The resulting blade domain has a total grid count of 3.2 million. This is made up of 143 elements in the spanwise direction (of which 30 are in the tip gap region), 72 elements in the circumferential direction and 310 elements in the chordwise direction (of which 147 are within the blade passage). Boundary layer meshing is used in the near wall regions over the blade surfaces (typically 30 layers normal to the blade) and the end wall surfaces.

Numerical setup and boundary conditions

The 3D steady RANS equations are solved using ANSYS CFX 17.1 code. The code uses a cell-centred finite volume approach to evaluate the numerical fluxes. Flux values at the cell-boundary are obtained using a high-resolution upwind scheme (ANSYS CFX Theory Guide, 2016). The high-resolution upwind scheme is a blended scheme that addresses the shortcomings of the first and second order upwind schemes. The flow is assumed to be fully turbulent and the standard k−ϵ model is selected for turbulence closure. A scalable wall function is required for resolving the laminar sublayer region (ANSYS CFX Theory Guide, 2016). The lessons learned from a previous study, Mustaffa and Kanjirakkad (2019, 2020), of a high-speed rotor with a similar tip loading suggest that the effect of wall bounded viscous phenomena such as separation plays a lesser role in the build-up of the near tip blockage. The above compressor was simulated both at 100% and at 60% (where no passage shock was present) speed. This makes the use of the k−ϵ turbulence model with wall functions justifiable in the current study.

At the inlet, the total pressure and total temperature profiles are prescribed. The steady calculation uses an inlet turbulence intensity of 5% along with an eddy viscosity ratio of 10. The total pressure profile used at the inlet boundary corresponds to that of a turbulent boundary layer with a displacement thickness of approximately 5% of the annulus height. A radial equilibrium exit static pressure condition is imposed at the outlet. Flow within the whole domain is solved in the rotating frame of reference, and all stationary walls are set to be counter-rotating. All walls are treated as smooth and adiabatic. The convergence rate is accelerated using an algebraic multi-grid (AMG) approach built into the solver.

Performance characteristics

The performance characteristics of the LSAC for the smooth casing case are shown in Figure 2. The abscissa of Figure 2 is the flow coefficient, ϕ, which is defined as the average axial velocity at the inlet, Vz,1, over the midspan blade speed, Umid. The total-to-static pressure rise coefficient, ψts,and isentropic efficiency, η, are calculated using Equations 1 and 2, respectively. The fluid density, ρ, in Equations 1 and 2 is calculated at the inlet. P and P0 are the static and total pressures, respectively and T0 is the total temperature. The static pressure in Equation 1 is its area averaged value at exit whereas the total pressure and total temperatures used in these equations are mass averaged at the appropriate stations.

Figure 2.

Performance characteristics of the LSAC. (a) Total-static pressure rise coefficient and (b) Isen-tropic efficiency.

Operating point 1, in Figure 2, corresponds to a hub exit static pressure of 101,350 Pa. The characteristic is obtained by gradually increasing the hub static pressure by a value between 50–100 Pa. Operating point 4 and 8 represent the conditions at near design and near stall respectively. The last stable operating point (Operating point 10) is determined iteratively using a “bisection method” type approach. Here, the size of the hub static pressure increment is updated until any further increase of 5 Pa fails to produce a converged solution. The convergence criteria set for all simulations are as following:

The computational methodology as used here has been extensively validated in Mustaffa and Kanjirakkad (2019, 2020).

Tip region flow

As shown in Figure 2a, the stability limit occurs before the ψts characteristic reaches its peak (and hence the curve still has a negative slope). A spike-stall behaviour can therefore be inferred based on the findings in Camp and Day (1998). The spike-stall is found to initiate near the tip region when the critical incidence is exceeded. Figure 3 shows the normalised axial velocity contour for an axial plane located at 0.04cax,t for conditions near stall (operating point 8). The low axial velocity or blockage region (generally seen in Figure 3 as regions where Vz/V¯z,1<0.8) is predominantly located near the rotor casing. The blockage as seen in the hub region is caused by a separation on the suction side near leading edge and is largely confined radially and circumferentially to the LE suction side hub corner. The blockage region extends much further in the circumferential direction near the casing compared to that near the hub. The findings here further suggest that the LSAC is indeed a tip critical compressor that is suitable for studying casing treatment.

Figure 3.

Normalised axial velocity contour plotted for an axial plane located at 0.04 c_(ax,t) at operating point 8.

The build-up of blockage near the casing as the compressor approaches stall can be investigated by plotting the normalised axial velocity contour inside the tip gap for operating point 4 (near design) and 8 (near stall) as shown in Figure 4. In comparison to operating point 4, the blockage region due to the TLF at operating point 8 occupies a relatively larger portion of tip gap region within the blade passage. This is expected since the TLF is pressure-driven such that at a lower ϕ, the blade loading increases (close to the blade tip LE) with increasing ψts. This results in a stronger TLV and TLF with an increased pitch-wise velocity component. It is interesting to note that there are two regions of increased axial velocity deficit within the passage; one in the forward part of the blade in line with the leading edge tip leakage vortex and the second region close to the blade trailing edge. Although not shown here, these regions are consistent with the local blade (tip) loading distribution.

Figure 4.

Normalised axial velocity contour inside the tip gap (50% of tip gap height, τ) for (a) operating point 4 and (b) operating point 8.

In addition, at operating point 8, it can be seen that the interface of the TLF and the incoming flow is more aligned to the circumferential direction compared to operating point 4. This is similar to what has been described in Hoying et al. (1999) where, at near stall conditions, it is found that the upstream movement of the TLV eventually leads to a “spill forward” effect. The blockage near the blade LE of the suction side (SS) is likely to be caused by a stronger TLV and the increased TLF. The TLV trajectory can be found by extracting the normalised static pressure profiles along several constant pitch lines as shown in Figure 5.

Figure 5.

TLV trajectory (red line) extracted from static pressure contour at mid tip gap at operating point 4.

For each constant pitch line, the axial coordinate of the minimum static pressure is found through interpolation. The trajectory of the TLV can be found through linear fitting. The slope of the TLV trajectories is represented as TLV angle, ξ. The angle ξ is measured from the axial direction and is positive if measured in the direction of rotation. This procedure is repeated for all operating points and the TLV angles are plotted as shown in Figure 6. Here the variation of ξ is plotted against the pressure rise coefficient, ψts, with annotations for the respective operating points as found in Figure 2a). As the backpressure increases from operating point 4 to operating point 10, the change in ξ is about 4°. The increase in ξ and a stronger TLF contributes to the upstream movement of the TLV-main flow interface.

Figure 6.

TLV angles with respect to the axial direction for all operating points.

Tip leakage momentum

As mentioned earlier, the increased loading near the blade tip LE at conditions close to stall, contributes to the increase of the TLF pitch-wise velocity component which, in turn, appears as blockage due to reduction in the axial velocity component. The strength of the TLF that also affects the strength of the TLV can be quantified by calculating the TLF momentum inside the tip gap region as shown in Equation 3.

Here, Mtip is the tip leakage flow momentum normalised by the blade tip momentum where Utip is the blade tip speed. V⊥ is the local velocity normal to the blade camber at the tip. The density, ρ can be omitted considering that at low-speed conditions, the flow is incompressible. The TLF momentum is calculated for operating point 4 and 8 as shown in Figure 7. The shaded region (between the tip LE and 20% of cax,t aft of tip LE) represents the area where the increase in Mtip due to the increased blade tip loading is significant. Within the shaded region it can be seen that Mtip at near stall conditions (operating point 8) is about 30% higher with respect to operating 4. The increased leakage momentum contributes to reduced axial velocity in the tip region. This coupled with the increase in back pressure at near stall conditions causes the upstream movement of the TLF-main flow interface. Interestingly, as the stall is approached, the tip leakage momentum is seen to decrease in the rear part of the blade (z~t>0.5). This is in agreement with the normalised axial velocity contours presented in Figure 4. The contours suggest that compared to the near design case, closer to stall, the low momentum region near the LE intensifies and moves upstream whereas that near the blade TE suction side becomes less intense. This is due to a reduction in the blade loading in the rear part of the blade. It is thought that the modification of the passage flow by the increased upstream blockage is responsible for this effect.

Figure 7.

Comparison of normalised tip leakage momentum at operating point 4 and 8.

Quantification of near casing blockage

The preceding discussion has shown that the upstream movement of the TLF-main flow interface can be linked to the increased blockage (axial velocity deficit) regions at conditions near stall. It is, therefore, thought to be a worthwhile exercise to quantify this blockage to further explore its link to the initiation of stall. The blockage quantification method is based on a mass flow overshoot method introduced by Sakuma et al. (2013) in order to identify and visualise “blocked” cells in an axial plane. Blocked cells are identified by sorting and summing the mass flow of each grid cell in that plane in a descending order. If cells with negative axial velocity (flow reversal) values exist, the sum of all positive axial velocity will overshoot the inlet mass flow value. The summation is stopped when the summed mass flow reaches the inlet mass flow value. Cells that are yet to be summed are considered as “blocked” and are assigned a binary blockage index (ψ) of 1. “Un-blocked” cells are assigned a ψ value of 0. The outcome of this method is shown in Figure 8 where the blue dots represents the location of the “blocked” cells (ψ=1) over an axial plane 0.04cax,t aft of the blade tip LE. Note that this is the same axial plane as in Figure 3 where the normalised axial velocity contours were presented.

Figure 8.

Location of the “blocked” cells at an axial plane location same as in Figure 3 at operating point 8 (looked from the front).

Based on the overshoot method described above, a non- dimensional blockage parameter, Ψm, is introduced to quantify the blockage in that plane. Ψm is the normalised absolute sum of the “blocked” cell mass flow rate. The definition of Ψm is shown in Equation 4.

Here, N is the total number of cells in the plane. ψ(i) and m(i) are the blockage index and the mass flow rate respectively that are associated with the ith cell on this plane.

Figure 9 shows the distribution of Ψm across the rotor blade for all the operating points that are shown in the performance map in Figure 2. The distribution is obtained by extracting data for about 100 axial planes across the blade. Here, the distribution of Ψm is only evaluated for the outer 20% span region. This is mainly because the blockage in the tip region is found to be more significant than at the hub region for a tip-critical compressor. The peak blockage located at approximately 0.15z~tnear design (operating point 4) grows rapidly as the compressor is throttled while gradually moving upstream. At the last numerically stable point (operating point 10) the peak blockage has grown almost eight times compared to the near design and has also moved upstream, albeit small amount, by approximately 3% of cax,t. The increased blockage is consistent with the TLF momentum increase at conditions near stall as shown in Figure 7. However, the relatively small upstream movement of the location of the peak blockage as the compressor approaches stall, seems to suggest that the peak is likely due to the build-up of blockage originating from the TLV. The TLV trajectory approximation exercise presented in Figure 6 also showed a moderate 4° change in TLV angle (ξ) from operating point 4 to 10. The rise in the peak value of the blockage at conditions near stall points to a stronger roll-up of the TLV that is a result of the increase in TLF momentum. This in turn causes an increase in the axial velocity deficit near the LE region. This increase of the axial velocity deficit caused by the TLV therefore results in the upstream movement of the TLF-main flow interface near stall.

Figure 9.

Distribution of Ψ_m across the blade for all operating points.

Near the TE (z~t>0.8), a small hump in the blockage distribution is visible for all operating points with flow coefficient above that of the NS case. This is consistent with the previous discussion on the appearance of a low momentum region near the TE with reference to Figures 4 and 7. As noted before, at conditions near stall, this low momentum region weakens, which in turn result in the flattening of the hump region as seen for operating points 7 through to 10.