The method of characteristics described in Section 13.2 is applicable only to supersonic flows; the characteristic lines are not defined in a practical fashion for steady, subsonic flow. Also, the particular finite-difference method outlined in Section 13.4 applies only to supersonic flows; if it were to be used in a locally subsonic region, the calculation would blow up. The reason for both of the above comments is that the method of characteristics and the steady flow, forward-marching finite-difference technique depend on the governing equations being mathematically “hyperbolic.” In contrast, the equations for steady subsonic flow are “elliptic.” (See Reference 21 for a description of these mathematical classifications.) The fact that the governing equations change their mathematical nature in going from locally supersonic to locally subsonic flow has historically caused theoretical aerodynamicists much grief. One problem in particular, namely, the mixed subsonic-supersonic flow over a supersonic blunt body as described in Section 9.5, was a major research area until a breakthrough was made in the late 1960s for its proper numerical solution. The purpose of this section is to describe a numerical finite-difference solution which readily allows the calculation of mixed subsonic-supersonic flows—the time-dependent method—and to show how it is used to solve supersonic blunt-body flows. Time-dependent techniques are very common in modem computational fluid dynamics, and as a student of aerodynamics, you should be familiar with their philosophy.
Consider a blunt body in a supersonic stream, as sketched in Figure 13.9a. The shape of the body is known and is given by b = b(y). For a given freestream Mach number Moo, we wish to calculate the shape and location of the detached shock wave, as well as the flow-field properties between the shock and the body. The physical aspects of this flow field were described in Section 9.5, which you should review before progressing further.
The flow around a blunt body in a supersonic stream is rotational. Why? Examine Figure 13.10, which illustrates several streamlines around the blunt body. The flow is inviscid and adiabatic. In the uniform freestream ahead of the shock wave, the entropy is the same for each streamline. However, in crossing the shock wave, each streamline traverses a different part of the wave, and hence experiences a different increase in entropy. That is, the streamline at point a in Figure 13.10 crosses a normal shock, and hence experiences a large increase in entropy, whereas the streamline at point b crosses a weaker, oblique shock, and therefore experiences a smaller increase in entropy, Sb < sa. The streamline at point c experiences an even weaker portion
Figure 13.9 Blunt-body flow field in both the
physical and computational planes.
of the shock, and hence sc < sh < sa. The net result is that in the flow between the shock and the body, the entropy along a given streamline is constant, whereas the entropy changes from one streamline to the next; that is, an entropy gradient exists normal to the streamlines. It can readily be shown (see chapter 6 of Reference 21) that an adiabatic flow with entropy gradients is rotational. Hence, the flow field over a supersonic blunt body is rotational.
In light of the above, we cannot use the velocity potential equation to analyze the blunt-body flow. Rather, the basic continuity, momentum, and energy equations must be employed in their fundamental form, given by Equations (7.40), (7.42a and b), and (7.44). With no body forces, these equations are
Notice the form of the above equations; the time derivatives are on the left, and all spatial derivatives are on the right. These equations are in the form necessary for a time-dependent finite-difference solution, as described below.
Return to Figure 13.9a. Recall that the body shape and freestream conditions are given, and we wish to calculate the shape and location of the shock wave as well as the flow field between the shock and body. We are interested in the steady flow over the blunt body; however, we use a time-dependent method to obtain the steady flow. The basic philosophy of this method is as follows. First, assume a shock-wave shape and location. Also, cover the flow field between the shock and body with a series of grid points, as sketched in Figure 13.9a. At each of these grid points, assume values of all the flow variables, p, u, v, etc. These assumed values are identified as initial conditions at time t = 0. With these assumed values, the spatial derivatives on the right sides of Equations (13.59) to (13.62) are known values (obtained from finite differences). Hence, Equations (13.59) to (13.62) allow the calculation of the time derivatives dp/dt, du/dt, etc. In turn, these time derivatives allow us to calculate the flow properties at each grid point at a later instant in time, say, At. The flow properties at time t = At are different from at t = 0. A repetition of this cycle gives the flow-field variables at all grid points at time t — 2At. As this cycle is repeated many hundreds of times, the flow-field properties at each grid point are calculated as a function of time. For example, the time variation of ut J is sketched in Figure
13.11. At each time step, the value of u, j is different; however, at large times the changes in Uij from one time step to another become small, and Uij approaches a steady-state value, as shown in Figure 13.11. It is this steady-state value that we want; the time-dependent approach is simply a means to that end. Moreover, the shock-wave shape and location will change with time; the new shock location and shape at each time step are calculated so as to satisfy the shock relations across the wave at each of the grid points immediately behind the wave. At large times, as the flow-field variables approach a steady state, the shock shape and location also approach a steady state. Because of the time-dependent motion of the shock wave, the wave shape is a function of both t and у as shown in Figure 13.9a, s = s(y, t).
Given this philosophy, let us examine a few details of the method. First, note that the finite-difference grid in Figure 13.9a is curved. We would like to apply our finite differences in a rectangular grid; hence, in Equations (13.59) to (13.62) the independent variables can be transformed as
£ =———- and г/ — у
where b = b(y) gives the abscissa of the body and і = s(y, t) gives the abscissa of the shock. The above transformation produces a rectangular grid in the computational plane, shown in Figure 13.9b, where the body corresponds to £ = 0 and the shock corresponds to $ — 1. All calculations are made in this transformed, computational plane.
The finite-difference calculations themselves can be carried out using MacCor – mack’s method (see Section 13.4) applied as follows. The flow-field variables can be advanced in time using a Taylor series in time; for example,
[13.63]
In Equation (13.63), we know the density at grid point (i, j) at time t; that is, we know pij(t). Then Equation (13.63) allows us to calculate the density at the same grid point at time t + At, that is, Pij(t + At), if we know a value of the average time derivative [(dp/df), j]ave. This time derivative is an average between times t and t + At and is obtained from a predictor-corrector process as follows.