Examples of Vector VS Tensor Notation  

Quality  Vector Notation  Tensor Notation 
Scalar  a  a 
Vector  [latex]\bar{a}[/latex]  [latex]a_i[/latex] 
2nd order tensor  [latex]\overline{\overline{a}}[/latex]

[latex]a_{ij}[/latex] 
dot product (scalar)  [latex]\bar{a}\cdot\bar{b}[/latex]  [latex]a_{i} b_{i}[/latex] 
cross product (vector)  [latex]\bar{a}\times\bar{b}[/latex]  [latex]\varepsilon_{ijk}a_{j}b_{k}[/latex] 
del operator (vector)  [latex]\bigtriangledown[/latex]  [latex]\frac{\partial}{\partial x_i}[/latex] 
gradient (vector; tensor)  [latex]\bigtriangledown a;\bigtriangledown\bar{a}[/latex]  [latex]\frac{\partial{a}}{\partial {x_{i}}};\frac{\partial{a_{i}}}{\partial{x_{j}}}[/latex] 
divergence (scalar)  [latex]\bigtriangledown\cdot\bar{a}[/latex]  [latex]\frac{\partial{a_{i}}}{\partial{x_{i}}}[/latex] 
curl (vector)  [latex]\bigtriangledown\times\bar{a}[/latex]  [latex]\varepsilon_{ijk}\frac{\partial{a_k}}{\partial{x_j}}[/latex] 
Laplace operator (scalar)  [latex]\bigtriangledown^2[/latex]  [latex]\frac{\partial^2}{\partial{x_i}\partial{x_i}}[/latex] 
divergence of a 2nd order tensor  [latex]\bigtriangledown\cdot \overline{\overline{a}}[/latex]  [latex]\frac{\partial{a_y}}{\partial{x_i}}(vector)[/latex] 
dot product: vector & del  [latex]\bar{a}\cdot\bigtriangledown[/latex]  [latex]a_j\frac{\partial}{\partial{x_j}}(scalar)[/latex] 
dot product: vector & gradient  [latex]\bar{a}\cdot\bigtriangledown\bar{b}[/latex]  [latex]a_j\frac{\partial{b_i}}{\partial{x_j}}(vector)[/latex] 
Material Derivative  [latex]\frac{\partial}{\partial{t}}+\bar{u}\cdot\bigtriangledown or \frac{D}{Dt}[/latex]  
Vorticity: (pseudovector)  [latex]\omega_k=\bigtriangledown \text{x}\bar{u}=\left (\frac{\partial{u_i}} {\partial{u_j}}\frac{\partial{x_j}}{\partial{x_i}} \right )=\varepsilon_{ijk}\frac{\partial {u_k}}{\partial{x_j}}[/latex]  
Rate of strain tensor: (2nd order tensor)  [latex]\frac{1}{2}\left (\frac{\partial{u_i}} {\partial{u_j}}+\frac{\partial{u_j}}{\partial{u_i}} \right )[/latex] 
[latex]\frac{DP}{Dt}=\frac{\partial P}{\partial t}+u\frac{\partial P}{\partial x}+v\frac{\partial P}{\partial y}+w\frac{\partial P}{\partial z}[/latex]
(2.1)
Things to pay attention to in this equation are (i) the velocity components [latex]{u,v}[/latex] and [latex]{w}[/latex] are not vector quantities, but components of the velocity vector, V; (ii) the time variation is specified by the first term and is evaluated at a particular point in space (as are all of the other terms) but all other terms may in fact be time dependent; this equation is valid for scalars or tensors, and the tensor sense is determined by the quantity being evaluated (in Equation (2.1) it is a scalar since [latex]{P}[/latex] is a scalar). The material derivative of the velocity vector V, is:[latex]\frac{DV}{dt}=\frac{\partial V}{\partial t}+u\frac{\partial V}{\partial x}+v\frac{\partial V}{\partial y}+w\frac{\partial V}{\partial z}[/latex]
(2.2)
Since [latex]{V}[/latex] is a vector then [latex]{DV/Dt}[/latex] is also a vector, and each term on the right hand side is a vector and the components of each of these vectors is determined by inserting the component for the vector [latex]{V}[/latex]. In other words the x component is:[latex]\frac{Du}{Dt}=\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}+w\frac{\partial u}{\partial z}[/latex]
(2.3)
The above is for Cartesian coordinates. This can be written for any selected coordinate system. For example the result in cylindrical coordinates in [latex]{r,\theta}[/latex] is shown below.Material Derivative in Cylindrical Coordinates
We would like to find an expression for DV/DT in cylindrical coordinates that we can use to help interpret streamline coordinates. We will only examine a two dimensional situation, [latex]r, \theta[/latex] since [latex]z[/latex] is similar to Cartesian coordinates. Given the vector [latex]V = u_r\hat{r}+u_\theta\hat{\theta}[/latex] we write this as a material derivative: [latex]\frac{DV}{Dt}=\frac{\partial V}{\partial t}+u_r\frac{\partial V}{\partial r}+\frac{u_\theta}{r}\frac{\partial V}{\partial \theta}[/latex] When the angular position changes then there can be a change in both the unit vectors [latex]\hat{r}[/latex] and [latex]\hat{θ}[/latex] since their orientation changes. So we can write the following: [latex]\frac{\partial V}{\partial\theta}=\frac{\partial\left ( u_r\hat{r}+u_\theta\hat{\theta} \right )}{\partial\theta}=\hat{r}\frac{\partial u_r}{\partial\theta}+u_r\frac{\partial\hat{r}}{\partial\theta}+\hat{\theta}\frac{\partial u_\theta}{\partial\theta}+u_\theta\frac{\partial\hat{\theta}}{\partial\theta}[/latex] Now we must interpret the derivatives of the unit vectors with respect to θ. Refer to the figure below showing a blob of fluid moving [latex]\partial\theta[/latex] that results in a change of the unit vector [latex]\hat{r}[/latex], that we can write as: [latex]d\hat{r}=\hat{\left  r \right }d\theta\hat{\theta}=d\theta\hat{\theta}[/latex], which we can write as: [latex]\frac{\partial\hat{r}}{\partial\theta}=\hat{\theta}[/latex]. Similarly the change in the unit vector, [latex]d\hat{\theta}=rd\theta[/latex], as illustrated above on the right as well. This is negative because if the change in the unit vector is counterclockwise the change in the r unit vector is in the negative r direction. The result is [latex]\frac{\partial\hat{\theta}}{\partial\theta}=r[/latex]. All of this can now be combined and used in the material derivative above resulting in: [latex]\frac{DV}{Dt}=\frac{\partial V}{\partial t}+u_r\frac{\partial V}{\partial r}+\frac{u_\theta\partial V}{r\partial\theta} = \left ( \frac{\partial u_r}{\partial t}+u_r\frac{\partial u_r}{\partial r}+\frac{\partial_\theta\partial u_r}{r\partial\theta}\frac{u_\theta^2}{r} \right )\hat{r}+\left ( \frac{\partial u_\theta}{\partial t}+u_r\frac{\partial u_\theta}{\partial r}+\frac{u_\theta\partial u_\theta}{r \partial\theta}+\frac{u_r u_\theta}{r} \right )\hat{\theta}[/latex] This applies to Cylindrical Coordinates, but as we shall see helps to interpret streamline coordinates when there is curvature to the streamline. The net result is added terms that account for the change of the unit vectors due to the curving nature of the flow.[latex]\frac{Du_i}{Dt}=\frac{\partial u_i}{\partial t}+u_j\frac{\partial u_i}{\partial x_j}[/latex]
(2.4)
Note here that the vector, V is replaced with the tensor notation for a vector written as [latex]{u{}_{i}}[/latex]. Also the summation rule has been used in the second term on the right hand side, this implies that there are really three terms being represented by the one term shown. The first term on the right hand side is denoted as the "local acceleration" and physically represents the rate change of velocity at a fixed point in space. The second expression on the right hand side is denoted as the "convective acceleration" and physically represents how the velocity field changes over space. If a flow is steady the local acceleration is by definition zero. However a fluid mass as it moves through space where the velocity changes in space will experience acceleration (or deceleration). An example is steady flow through a nozzle, or diffuser, where the fluid velocity will increase, or decrease, along the flow direction. For these two cases a fixed mass of fluid will accelerate, or decelerate due to moving in space. The reader is reminded that the components of the acceleration depend on the component of velocity in the derivatives.The convective acceleration can be rewritten using a vector identity listed above. Namely, the convective acceleration in vector notation is [latex]\left(V\cdot \mathrm{\nabla }\right)V[/latex] and can be expressed as:[latex]\left(V\cdot \mathrm{\nabla }\right)V=\frac{1}{2}\mathrm{\nabla }\left(\mathrm{V}\mathrm{\cdot }\mathrm{V}\right)\mathrm{}\mathrm{V\times (}\mathrm{\nabla }\mathrm{\times V)}[/latex]. In tensor notation this becomes:[latex]u_j\frac{\partial u_i}{\partial x_j}=\frac{1}{2}\frac{\partial \left(u_ju_j\right)}{\partial x_j}{\varepsilon }_{ijk}\left(u_j{\varepsilon }_{klm}\frac{\partial u_m}{\partial x_l}\right)[/latex]
(2.5)
Consequently, the convective acceleration can be recast as the expression on the right hand side. Although this may not seem like it is of any benefit and just makes things more complicated we will show that it does indeed have advantages. The reader should pay attention to the subscripts in the above expression. Notice that each term is a vector with the free index of “[latex]{i}[/latex]”, all other indices are repeated and therefore summed. In the summation of vectors we basically need to sum components, this implies that they have the same index, in this case “[latex]{i}[/latex]”.[latex]a_s=\frac{\partial u_s}{\partial t}+u_s\frac{\partial u_s}{\partial s}=\frac{\partial u_s}{\partial t}+\frac{\partial {u_s}^2}{2}[/latex]
(2.6)
[latex]a_n=\frac{{u_s}^2}{{R}}[/latex] Here [latex]{R}[/latex] is the local value of the radius of curvature of the streamline. The value of [latex]{R}[/latex] goes to infinity for a straight streamline and there is no acceleration in the [latex]{n}[/latex] direction. Recall here that [latex]{u{}_{r}}[/latex] (or [latex]{u{}_{n}}[/latex]) is identically zero at the streamline and that the expression for [latex]a_n[/latex]only contains the last term associated with the changing unit vector [latex]{θ}[/latex], or in this case [latex]{s}[/latex]. Also we define in the streamline coordinates the direction [latex]{n}[/latex] to be inward towards the radius of curvature as shown in Fig. 2.1. This is opposite to the convention used in cylindrical coordinates where [latex]{r}[/latex] is radially outward, hence the sign is positive for the value of [latex]a_n[/latex]. That is, if there is a curving streamline the acceleration is inward towards the direction of the center of the radius of curvature. [caption id="attachment_251" align="aligncenter" width="548"] Fig. 2.1 Streamline coordinate illustration; $latex {s}$ is along the velocity vector direction and tangent to the streamline and $latex {n} $ is normal to $latex {s}$ while directed towards the center of radius, denoted as $latex {R}$.[/caption][latex]\frac{\partial u_j}{\partial x_i}=\frac{1}{2}\left(\frac{\partial u_j}{\partial x_i}\mathrm{\ }\mathrm{+}\frac{\partial u_i}{\partial x_j}\mathrm{\ }\right)+\frac{1}{2}\left(\frac{\partial u_j}{\partial x_i}\mathrm{\ }\mathrm{}\frac{\partial u_i}{\partial x_j}\mathrm{\ }\right)[/latex]
(2.7)
Notice here that the first term adds half of the transpose of the original tensor and the second term subtracts the same quantity completing this identity. The first term is symmetric that is by inverting the [latex]{i}[/latex] and [latex]{j}[/latex] subscripts there is no change in the value of the element. The second term is antisymmetric such that inverting the subscripts the value of the element has the opposite sign. The ½ multiplier is included in the definition of the symmetric and antisymmetric parts such that their sum results in the original tensor. In fluid mechanics applications the velocity gradient tensor is a measure of the velocity changes experienced in the fluid flow as one moves away (infinitesimally) from a given point. That is to say this tensor is a point function, having a value at any given point in the flow. We refer to that as a local function in space. It may be time dependent as well. So we can say that this tensor has a symmetric part and an antisymmetric part. As we will see when we discuss the viscous forces within a fluid the symmetric part is representative of the deformation rate, or strain rate, experienced by a fluid element. The antisymmetric part represents the rotation rate of a fluid element. If one examines the curl operator on [latex]{u{}_{j}}[/latex] given previously as [latex]{\varepsilon }_{kij}\frac{\partial u_j}{\partial x_i}[/latex], where the subscripts indicate that this is the [latex]{k}th[/latex] element, it can be shown by writing out the terms of the curl operation on [latex]{u{}_{j}}[/latex] that this is equivalent to twice the antisymmetric part of [latex]\frac{\partial u_j}{\partial x_i}[/latex] (or the part in parentheses of the second term in Eqn. (2.7)). Therefore the rotation part of the velocity gradient tensor is one half the vorticity of the flow. To be clear, the antisymmetric part of [latex]\frac{\partial u_j}{\partial x_i}[/latex] is a second order tensor, but from observation if i=j then the value of this is identically zero (all zero components along the diagonal of the tensor). The off diagonal terms are six in total but those with inverted indices are the negatives of the noninverted components ([latex]\frac{\partial u_j}{\partial x_i}=\frac{\partial u_i}{\partial x_j}[/latex]). Consequently there are only three independent values, so we identify this second order tensor as a pseudovector, having three components. It is a so called pseudovector since the sign is arbitrary. In other words do we set [latex]\frac{\partial u_1}{\partial x_2}[/latex] as a positive or negative quantity? The sign convention is typically selected so that a counter clockwise rotation is a positive value. This means that for [latex]\frac{\partial u_1}{\partial x_2}[/latex] being positive then as [latex]{u{}_{1}}[/latex] increases in the positive [latex]{x{}_{2}}[/latex] direction then the flow will tend to rotate clockwise and is a negative rotation direction. This is illustrated in Fig. 2.2. Note that if one considers a negative value of [latex]\frac{\partial u_1}{\partial x_2}[/latex] coupled with a positive value of [latex]\frac{\partial u_2}{\partial x_1}[/latex] such that$latex {\omega_3}=\left (\frac {\partial u_2}{\partial x_1} \frac {\partial u_1}{\partial x_2} \right)$ > $latex 0$
(2.8)
then rotation (or vorticity) about the 3 axis (out of the page in Fig. 2.2) is counterclockwise. Obviously, the signs of the two components of the velocity gradient that go into any of the vorticity components can be positive and/or negative. If in the above example for Eqn. 2.8 the magnitudes of [latex]\frac{\partial u_2}{\partial x_1}[/latex] and [latex]\frac{\partial u_1}{\partial x_2}[/latex] are equal and their signs are the same then the vorticity component [latex]{\omega }{_3}[/latex] will be zero. One can think of the net vorticity or rotation rate vector component to be the combined value from both velocity derivative components, which contribute to rotation about an axis. The same interpretation of the combined effects of the two velocity gradient elements can be made for the strain rate, or the symmetric part of the velocity gradient tensor. When the rotation rates of the vertical and horizontal axes are identical (both elements causing, say, counterclockwise rotation at the same rate) then the strain rate will be zero (the two elements in the symmetric part cancel) and the flow is in pure rotation, as a solid body. [caption id="attachment_496" align="aligncenter" width="464"] Fig. 2.2 Illustration of rotation rate as determined by the velocity gradient tensor components; the original square element has its sides experiencing rotation at a rate of $latex \frac{\partial u_2}{\partial x_1}$ for the horizontal axis and $latex \frac{\partial u_1}{\partial x_2}$ for the vertical axis.[/caption][latex]\frac{P}{\rho }+\frac{V^2}{2}+gz=constant[/latex]
(3.1)
The restrictions placed on the application of this equation are rather limiting, but still this form of the equation is very powerful and can be applied to a large number of applications. But since it is so restrictive care must be take in its application. The restrictions can be stated as:[latex]\frac{D\boldsymbol{V}}{dt}=\frac{1}{\rho }\frac{\partial P}{\partial \boldsymbol{x}}+\boldsymbol{g}[/latex]
(3.2)
Or in tensor notation as:$latex \frac{\partial u_i}{\partial t}+u_j\frac{\partial u_i}{\partial x_j}=\frac{1}{\rho }\frac{\partial P}{\partial x_i}+g_i$
(3.3)
The pressure term is representative of the net force caused by the compressive load of pressure along the direction of the vector component of interest, on a per mass basis. The body force term is again given as the “[latex]{i}[/latex]” vector component of the gravitational vector. The acceleration terms on the left can be recast using the vector identity for the convective acceleration given in Chapter 2.[latex] u_j\frac{\partial u_i}{\partial x_j}=\frac{1}{2}\frac{\partial \left(u_ju_j\right)}{\partial x_j}{\varepsilon }_{ijk}\left(u_j{\varepsilon }_{klm}\frac{\partial u_m}{\partial x_l}\right) [/latex]
(3.4)
The last term introduces the vorticity into the equation since the vorticity is defined as the curl of the velocity vector:[latex] {\omega }_k={\varepsilon }_{klm}\frac{\partial u_m}{\partial x_l} [/latex]
So the last term becomes the cross product of velocity and vorticity, or: [latex]{\varepsilon }_{ijk}\left(u_j{\varepsilon }_{klm}\frac{\partial u_m}{\partial x_l}\right)={\varepsilon }_{ijk}u_j{\omega }_k[/latex] Notice that the result of this operation is a vector in the “[latex]{i}[/latex]” direction, consistent with the other terms in the Euler equation. By introducing the vorticity into the convective acceleration term, this term can now can be considered to have two components. One is half of the gradient of the magnitude of the velocity vector squared (the first term), and the other is the cross product between the velocity vector and vorticity (a vector). This latter term has been identified as the "Lamb vector", after the applied mathematician, Horace Lamb (18491934). Euler's equation can now be rewritten as:[latex] \frac{\partial u_i}{\partial t}+\frac{1}{2}\frac{\partial \left(u_ju_j\right)}{\partial x_i}{\varepsilon }_{ijk}u_j{\omega }_k=\frac{1}{\rho }\frac{\partial P}{\partial x_i}+g_i [/latex]
(3.5)
Next we write the Lamb vector as a gradient of some unknown function:[latex] \frac{\partial \pi }{\partial x_i}={\varepsilon }_{ijk}u_j{\omega }_k [/latex]
(3.6)
In this expression [latex]{\pi}[/latex] must be a scalar quantity such that its gradient is a vector. This is done for convenience so that we can then integrate each term of our equation along any desired path, with elemental distance [latex]{ds}[/latex]. That is we take the projection of each term along [latex]{ds}[/latex] (take the dot product of each term with vector [latex]{ds})[/latex]. However, before we do this we make the following modification to the gravitational term such that we can integrate spatially along [latex]{ds}[/latex]. We can write the gravitational term as a "potential" as:[latex] g_i=g\frac{\partial h}{\partial x_i} [/latex]
(3.8)
where [latex]“{h}”[/latex] is a scalar and $latex g$ is the magnitude acceleration of gravity. The choice of the symbol [latex]“{h}”[/latex] is because it will have units of length and if gravity is vertically downward then [latex]“{h}”[/latex] is vertically upward and represents elevation above some arbitrarily chosen datum. Notice that if [latex]“{i}”[/latex] is horizontal then [latex]\frac{\partial h}{\partial x_i}[/latex] is the change in elevation over a differential change in the horizontal direction, consequently the value of [latex]\frac{\partial h}{\partial x_i}[/latex] is zero. If [latex]“{i}”[/latex] is in the vertical upward direction then [latex]\frac{\partial h}{\partial x_i}=1[/latex] and [latex]g_i=g[/latex]. If [latex]“{i}”[/latex] is at some angle [latex]{\theta}[/latex] to the horizontal, then [latex]\frac{\partial h}{\partial x_i}={\mathrm{sin}\theta}[/latex]. We are now able to write the integral along [latex]{ds}[/latex] of Euler's equation. Note that the dot product of each term results in a scalar that defines the magnitude of the change of the term along direction [latex]{ds}[/latex]: [latex] \int_s{\frac{\partial u_i}{\partial t}}ds_i+\int_s{\frac{1}{2}\frac{\partial \left(u_ju_j\right)}{\partial x_j}}ds_i+\int_s{\frac{\partial \pi }{\partial x_i}}ds_i+\int_s{\frac{1}{\rho }\frac{\partial P}{\partial x_i}ds_i+\int_s{g\frac{\partial h}{\partial x_i}ds_i=f(t)}} [/latex](3.9)
Since the integration is in space there may be a time dependence of each term that is not accounted for which is why the term on the right hand side appears  a general function of time can be added such that if one takes the spatial derivatives this function will vanish. Note that if there is no time dependence then [latex]{f(t)} = 0 [/latex], as well as the first term on the left. Bringing the time derivative outside of the integral in the first term, and performing the integration along the arbitrary line [latex]{ds}[/latex]:[latex] \frac{\partial}{\partial{t}}\int_s{u_i}ds_i+\Delta \left(\frac{\left(u_ju_j\right)}{2}\right)+\Delta \pi +\int_s{\frac{1}{\rho }\frac{\partial P}{\partial x_i}ds_i+\Delta \left(gh\right)=f(t)} [/latex]
(3.10)
where [latex]\Delta [/latex] represents the change of the variable, or expression, along the line [latex]{s}[/latex]. As shown we have two integrals that we can not evaluate at this point. The first represents the local acceleration and how it varies along the integration path. The second is the pressure term. Note however if the flow is incompressible, [latex]{\rho } [/latex] is constant, then we can say that: [latex]\int _S{\frac{1}{\rho }\frac{\partial P}{\partial x_i}ds_i=\frac{1}{\rho }\int_s{\frac{\partial P}{\partial x_i}ds_i}}=\Delta \frac{P}{\rho }[/latex] Equation (3.10) excludes frictional or viscous forces, but that is about its only limitation. As such it is a general form of the Bernoulli Equation. But considering incompressible and steady flow the result is:[latex] \Delta \left(\frac{\left(u_ju_j\right)}{2}\right)+\Delta \pi +\Delta \frac{P}{\rho }+\Delta \left(gh\right)=0 [/latex]
(3.11)
Consequently, the sum of these four terms which represent changes along any direction [latex]{s}[/latex] is zero, or[latex] \frac{\left(u_ju_j\right)}{2}+\pi +\frac{P}{\rho }+\left(gh\right)=constant [/latex]
(3.12)
To satisfy some curiosity, one would expect that by applying the four conditions listed above for the specialized form of Eqn. (3.10) will result in Eqn. (3.1). If density is constant the pressure integral term becomes [latex]{P/\rho}[/latex]. If the flow is steady then the time dependence terms are zero. If the integration is taken along a streamline then we can make the following argument. The velocity vector by definition is aligned with the streamline, therefore the direction of [latex]{ds}[/latex] and [latex]{u{}_{i}}[/latex] are identical. The vorticity may be in any arbitrary direction, but the cross product of [latex]{u{}_{i}}[/latex] and [latex]{\omega_j}[/latex] must be perpendicular to [latex]{ u{}_{i}}[/latex]. Consequently, since the integration is for the projection of each term along [latex]{ds}[/latex], and the cross product of $latex {u{}_{i}}$ and [latex]{\omega_j}[/latex] is normal to [latex]{ds}[/latex] (since it is normal to [latex]{u{}_{i}}[/latex]) then the net effect of this term is zero (has no component along [latex]{ds}[/latex] if [latex]{s }[/latex] is a streamline). So if integration is along a streamline we can delete the Lamb vector effect,[latex]\Delta \pi =0[/latex]. Combining all of these conditions we end up with Eqn (3.1) as we hoped. Other than along a streamline, another way in which [latex]\Delta \pi =0[/latex] is if the vorticity is zero along any chosen integration path. That is to say, the vorticity is zero throughout the flow field. Let's examine the consequences of this condition. If each component of vorticity is zero then we can write a set of conditions, one for each component of the vorticity,[latex] {\omega }_i=\left(\frac{\partial u_k}{\partial x_j}\frac{\partial u_j}{\partial x_k}\right)=0;\mathrm{where}\ i\neq j\neq k[/latex]
(3.13)
Now we define a "velocity potential" such that[latex]u_j=\frac{\partial \phi }{\partial x_j}[/latex] for all values of [latex]j[/latex]
(3.14)
This is a rather powerful condition  that a single scalar function,[latex]{\phi}[/latex], can be used to define the velocity field through its partial derivatives. If one takes the derivatives of [latex]{u{}_{k}}[/latex] and [latex]{u{}_{j}}[/latex] as shown in Eqn. (3.13) but rewrite this in terms of the velocity potential, the result is an identity for [latex]j\neq k[/latex]. That is to say, it is always = 0 since it is possible to reverse the order of differentiation. Therefore it is concluded that if one can replace the velocity vector with the derivatives of the scalar velocity potential defined in Eqn. (3.14) then the vorticity is zero everywhere. This can also be stated as: a scalar velocity potential exists and can be used to define the velocity field if the flow is vorticity free. Since vorticity can be defined as the degree of local rotation occurring in the flow we say the flow is "irrotational" if the vorticity is zero everywhere. Using the above definition of a velocity potential, that exists for irrotational flow, it can be said that irrotational flow results in a simplified form of the Bernoulli Equation since [latex]π[/latex] is zero when the vorticity is zero. Here we summarize the Bernoulli Equation and how it is modified for different conditions.[latex]\frac{\partial \phi }{\partial x_i}=u_i[/latex]
(4.1) This equation defines each component of the velocity in terms of the local spatial partial derivative in the direction of the velocity component. As stated in Chapter 2 this definition when inserted for velocity in the definition of the vorticity results in the identity that the vorticity is zero, hence irrotational flow.[latex] {\dot{m}}_{out,x_1}{\dot{m}}_{in,x_1}=\left[\rho u_1dx_2dx_3+\left(\frac{\partial }{\partial x_1}\left(\rho u_1\right)dx_1\right)dx_2dx_3\right]\left[\rho u_1dx_2dx_3\right]=\left(\frac{\partial }{\partial x_1}\left(\rho u_1\right)dx_1\right)dx_2dx_3 [/latex]
(4.2)
This is repeated for the [latex]{x{}_{2}}[/latex] and [latex]{x{}_{3}}[/latex] direction where changes of density times velocity are with respect to [latex]\partial x_2 \mathrm{and}\ \partial x_3[/latex] respectively, and [latex]{u{}_{2}}[/latex] and [latex]{u{}_{3}}[/latex] are used for the velocity, respectively, while using the area [latex]{dx{}_{1}x{}_{3}}[/latex] and [latex]{d{}_{1}xdx{}_{2}}[/latex], respectively. Summing all three of these net flow rates results in a scalar representation of the difference between the outflow and inflow where all possible flow paths in and out are included. Reversing the sign (to make it inflow minus outflow) this must equal to the change of mass within the volume element, [latex]{dx{}_{1}dx{}_{2}dx{}_{3}}[/latex].
This is expressed as: [latex]\frac{\partial (\rho )}{\partial t}dx_1dx_2dx_3=\frac{\partial }{\partial x_1}\left(\rho u_1\right)dx_1dx_2dx_3\frac{\partial }{\partial x_2}\left(\rho u_2\right)dx_1dx_2dx_3\frac{\partial }{\partial x_3}\left(\rho u_3\right)dx_1dx_2dx_3[/latex] or rearranging and dividing each term by [latex]dx_1dx_2dx_3[/latex]:[latex] \frac{\partial (\rho )}{\partial t}+\frac{\partial }{\partial x_1}\left(\rho u_1\right)+\frac{\partial }{\partial x_2}\left(\rho u_2\right)+\frac{\partial }{\partial x_3}\left(\rho u_3\right)=0 [/latex]
(4.3)
This is the continuity equation that must be satisfied to conserve mass. Notice that if the flow is steady the first term is zero. Also if the density is constant (incompressible) then the first term, or the partial derivative with respect to time, is zero, and density can be factored from each of the other terms and divided out of the equation. The result is:[latex] \frac{\partial u_1}{\partial x_1}+\frac{\partial u_2}{\partial x_2}+\frac{\partial u_3}{\partial x_3}=0=\frac{\partial u_i}{\partial x_i} [/latex]
(4.4)
where in the last term there is summation by tensor notation. This is the reduced form of continuity for incompressible flow. Notice that this form does not require the flow be steady (even though the unsteady derivative of density is no longer included). The velocity may in fact vary with time. If we now insert the definition of the velocity potential from Eqn. (4.1) for each of the velocity components in Eqn. (4.4) we end up with the following equation for [latex]{\phi}[/latex] for incompressible flow:[latex] \frac{{\partial }^2\phi }{\partial {x_1}^2}+\frac{{\partial }^2\phi }{\partial {x_2}^2}+\frac{{\partial }^2\phi }{\partial {x_3}^2}=0 [/latex]
(4.5)
This equation is the Laplace operation on the scalar velocity potential, [latex]{\phi}[/latex], and represents continuity (or conservation of mass) for an incompressible flow. Before moving on we write the continuity equation using the Material Derivative from Chapter 2. We combine the time derivative of density with the other three terms but notice that there is a difference in the three spatial derivative terms from those found in the Material Derivative, the velocity is included in the derivative in continuity. So if each of these terms is expanded, [latex]\frac{\partial \left({\rho u}_i\right)}{\partial x_i}=\rho \frac{\partial \left(u_i\right)}{\partial x_i}+u_i\frac{\partial \left(\rho \right)}{\partial x_i}[/latex] The we see that:[latex]\frac{D\rho }{Dt}+\rho \frac{\partial \left(u_i\right)}{\partial x_i}=0 [/latex]
And if the density is [latex]y[/latex] constant we obtain Eqn. (4.4) as expected.
[latex] u_1=\frac{\partial \psi }{\partial x_2}\ \mathrm{and}\ u_2=\frac{\partial \psi }{\partial x_1} [/latex]
(4.7)
Eqn. (4.7) represents the definition of the scalar [latex]\psi (x_1,x_2)[/latex] associated with each streamline. If we know expressions for the velocity components as a function of position then we can integrate Eqn. (4.7) to find the value of the streamfunction, [latex]\psi[/latex]. Similarly if we know the equation for the streamfunction then we can calculate the values of each velocity component through partial differentiation using Eqn. (4.7). Let's assume an incompressible flow so that the flow field follows the continuity equation given by Eqn. (4.4). Now insert for the derivatives of the velocities in terms of the derivatives of the streamfunction. The results is: [latex]\frac{{\partial }^2\psi }{\partial x_1\partial x_2}+\frac{{\partial }^2\psi }{\partial x_1\partial x_2}=0[/latex] This is an identity, in other words it is automatically true, so the existence of the streamfunction, by the given definition, automatically solves the continuity equation. Said another way, if the streamfunction exists by its definition of Eqn. (4.7) then the flow satisfies continuity for incompressible conditions. It is possible to take a given velocity field and construct a number of streamlines. At any given point there is a velocity vector and therefor a streamline that passes through it. The only time there can be two or more streamlines passing through a given point (intersecting at some random angle) is if the magnitude of the velocity is zero. Then both partial derivatives of Eqn. (4.7) are zero and the slope is not defined. A stagnation point is such an intersection of streamlines, as shown in Fig. (4.2) for flow over a cylinder. The streamfunction can be continuous up to the stagnation point and beyond, say following the cylinder surface, but it divides at the stagnation point one branch going up and another going down. Stagnation points don't have to be on surfaces they can be distributed within the flow field. [caption id="attachment_254" align="aligncenter" width="468"] Fig. 4.2 Streamlines illustrating a stagnation point streamline where flow separates away from this point, and the velocity at this point is zero.[/caption] Streamfunctions are valuable in that they can provide information on local flow rate conditions within a flow field. In general the flow rate (mass or volume) is determined by the velocity vector and an area through which the flow occurs. That is to say, the velocity vector only provides flow rate through an area if there is velocity vector component normal to the area. For a given area we define an outward normal unit vector, [latex]\hat{n}[/latex] as shown in Fig. (4.3). The mass flow rate through the area "A" with this outward normal is given as:[latex]\cdot{m}=\rho \boldsymbol{V}\bullet \hat{n}A[/latex]
(4.7)
The reader should check the units for this equation. Notice that [latex]\boldsymbol{V}\cdot \hat{n}[/latex] is the dot product between the velocity and outward normal that results in a scalar whose value represents the projection of the velocity vector in the [latex]\hat{n}[/latex] direction. To obtain the volume flow rate, [latex]\dot{Q}[/latex] this expression is divided by mass per volume, or the density:[latex]\cdot{Q}=\boldsymbol{V}\bullet \hat{n}A[/latex]
(4.8)
(4.9)
The reason there is a negative sign for [latex]u_2\Delta x_1[/latex]is that in determining the flow rate between points B and C we integrate along negative [latex]\Delta x_1[/latex] direction (or the change in [latex]\Delta x_1[/latex] is negative) Also we have used the definition of the streamfuction, Eqn. (4.7) to evaluate finite changes, on [latex]u_1=\bigtriangleup\psi/\bigtriangleup{x_2}[/latex] and [latex]u_2=\bigtriangleup\psi/\bigtriangleup{x_1}[/latex]. The interpretation then is that the flow rate [latex]\dot{Q}{'}[/latex] is equivalent to the change in streamfunction value between two points within the flow. [caption id="attachment_384" align="aligncenter" width="478"] Fig. 4.4 Flow between two different streamlines illustrating constant flow rate between the streamlines; control volume ABC is used to show relationships between flow rate and change in streamfunction values, Eqn. 4.9.[/caption] Interestingly, one can use streamline maps to qualitatively and quantitatively evaluate the velocity field. Imagine a wind tunnel test in a two dimensional flow over, say, a wing, as in Fig. (4.5). Smoke dye is injected at discrete points upstream separated by some vertical distance between each streamline. The lines of smoke travel downstream and over the wing. As the flow goes over the wing some of the streamlines diverge and some converge (the distance of separation between streamlines changes). Since we have shown that there is constant flow rate between streamlines when the distance between streamlines gets smaller the area of the flow decreases, so the velocity must increase. As streamlines diverge the velocity must decrease. The relationship between cross sectional area and velocity is linear as shown in Eqn. (4.8). Measuring the change of distance between adjacent streamlines provides a measure of the amount of increase or decrease of velocity. [caption id="attachment_287" align="aligncenter" width="309"] Fig. 4.5 Smoke flow visualization of flow over an inclined flat wing.[/caption] We now have a physical as well as mathematical interpretation for the streamfunction. Remember that this is a scalar field function representative of the local velocity. If we use the definition of streamfunction, Eqn. (4.7) and insert this into the definition of vorticity for a two dimensional flow in the [latex]x_1x_2[/latex] plane we obtain the following:[latex] {\omega }_3=\left(\frac{\partial u_2}{\partial x_1}\frac{\partial u_1}{\partial x_2}\right)=\left(\frac{\partial \psi_2 }{\partial {x_1}^2}\frac{\partial \psi_2 }{\partial {x_2}^2}\right)=\left(\frac{\partial \psi_2 }{\partial {x_1}^2}+\frac{\partial \psi_2 }{\partial {x_2}^2}\right) [/latex]
(4.10)
We see that the vorticity is equal to the negative of the Laplace of the streamfunction (Shown here in two dimensional flow, but is also the case in three dimensional flow). For [latex]{irrotational flow}[/latex] the vorticity is identically zero so: [latex]\left(\frac{\partial \psi_2 }{\partial {x_1}^2}+\frac{\partial \psi_2 }{\partial {x_2}^2}\right)=0[/latex] This shows that the Laplace of the streamfunction is zero for irrotational flow and follows the same results for the velocity potential for incompressible flow. So for irrotational, incompressible (ideal) flow the Laplace of both the velocity potential and streamfunction are equal to zero. This points to the ability to solve the Laplace equation for either of these quantities and from this solution determine the velocity field from the definitions of [latex]\psi[/latex] and [latex]\phi[/latex] in teams of the velocity. This will be the approach we take in the next Chapter. Since both [latex]\phi \mathrm{\ and}\ \psi [/latex] have the same functional form one might think that they are related to each other. We see this in comparing Eqns. (4.1) and (4.7), both are related to velocity derivatives. Notice that [latex] u_1=\frac{\partial \phi }{\partial x_1}=\frac{\partial \psi }{\partial x_2}\mathrm{\ and}\ u_2=\frac{\partial \phi }{\partial x_2}=\frac{\partial \psi }{\partial x_1}[/latex] The velocity is tangent to the constant streamfunction value, but the velocity is normal to the constant velocity potential value. Consequently lines of constant [latex]\phi [/latex] are normal to lines of constant [latex]\psi [/latex]. It is straightforward to show that the slope of constant potential lines is [latex]\frac{u_1}{u_2}[/latex] while the slope of constant streamfunction lines is [latex]\frac{u_2}{u_1}[/latex]. We will generate plots for specific flows in the next chapter to illustrate this, but shown in Fig. 4.6 is an illustration of the orthogonality of the velocity potential lines relative to the streamfunction lines. Velocity is always along the streamfunction line and normal to the potential lines. [caption id="attachment_257" align="aligncenter" width="300"] Fig. 4.6 Illustration of streamlines (lines of constant $latex \psi$) and velocity potential lines, $latex \phi$, indicating the orthogonal condition of the two lines.[/caption] In Table 4.1 are the two dimensional expressions in cylindrical coordinates for the various mathematical representations presented here. Note that [latex]{v{}_{r}}[/latex] is the radial velocity component and [latex]{v{}_{r}}[/latex] is the circumferential velocity component. The vorticity only has a z component, as shown, all others are identically zero for a two dimensional flow in [latex]{{r},\Theta}[/latex].[latex]\psi =Ux_2+C [/latex]
(5.1)
We can arbitrarily set [latex]\psi = 0[/latex] at [latex]{x{}_{2} = 0}[/latex] so that[latex]\psi =Ux_2[/latex]
(5.2)
If we proceed along the same line to solve for [latex]{\phi}[/latex] from its definition relative to the partial velocity derivative, Eqn. (4.1), the result is:[latex]\phi =Ux_1 [/latex]
(5.3)
where we have set [latex]{\phi}= 0[/latex] at [latex]{x{}_{1}} = 0[/latex], other wise there is an additive constant which is equal to the potential at [latex]{x{}_{1}} = 0[/latex]. For a condition of uniform flow, again with velocity of U, but at some angle [latex]{\alpha }[/latex] to the [latex]{x{}_{1}}[/latex] axis we have the following (the reader is encouraged to obtain this result by integrating the definition of the streamfunction): [latex]\psi =U\left(x_2cos\alpha x_1sin\alpha \right)+C [/latex](5.4)
The constant C is determined by a "datum" value of [latex]{\psi}[/latex] that passes through some point. For instance if [latex]\psi =0[/latex] at [latex](0,0)[/latex] then [latex]{C }= 0[/latex]. Going back to uniform flow (two dimensional) aligned in the [latex]{x_{1}}[/latex] direction we can find the volumetric flow rate per depth through some area say between [latex]{x{}_{2}} = 0[/latex] and [latex]{x{}_{2}}= 4[/latex] as: [latex]\dot{Q}{'}[/latex] = [latex] \Delta \psi =U(40)=4U[/latex][latex]v_r=\frac{\dot{m}{'}}{2\pi \rho }\frac{1}{r}=\frac{\dot{Q}{'}}{2\pi }\frac{1}{r} [/latex]
(5.5)
For a constant steady volumetric flow, [latex]\dot{Q}{'}[/latex], we see that the velocity decreases linearly as [latex]{ r}[/latex] increases. So this flow has straight streamlines all represented as radial lines from the center point. We see how the velocity increases rapidly in magnitude as [latex]{ r}[/latex] goes to zero, at which point it becomes infinite (or undefined). Based on the orthogonal condition between the velocity potential and streamfunction we see that lines of constant potential are circles (like the dashed line shown in Fig. 5.1.) We can combine the constants in Eqn. (5.5) [latex]{\mu }_s=\frac{\dot{Q}{'}}{2\pi }[/latex], where [latex]{{\mu}_{s}}[/latex] is the "source strength" with units of [latex]m{}^{2}/s[/latex] in SI units such that:[latex]v_{r,source}=\frac{{\mu }_s}{r} [/latex]
(5.6)
Eqn. (5.6) represents the flow from a source with strength [latex]{{\mu}_{s}}[/latex]. The source pumps fluid into the plane of flow with streamlines that are radially outward. The flow can be reversed such that the flow along the radial lines is inward. This implies flow is exiting the plane, and this is a "sink", with velocity inward, which has a negative [latex]{r }[/latex]component:[latex]v_{r,sink}=\frac{{\mu }_s}{r} [/latex]
(5.7)
The streamfunction is found by integrating the velocity based on the definition of the streamfunction using cylindrical coordinates: [latex]v_r=\frac{\partial \psi }{r\partial \theta }=\frac{{\mu }_r}{r}[/latex][latex]\psi ={\mu }_s\theta [/latex]
(5.8)
Velocity potential lines are found from the definition of the potential relative to the velocity: [latex]v_r=\frac{\partial \phi }{\partial r}=\frac{{\mu }_s}{r}[/latex] Integrating this results in:[latex]\phi ={\mu }_sln\ r[/latex]
(5.9)
[latex]\psi =Cln\ r[/latex]
(5.10)
Next we define the "circulation" [latex]\mathrm{\Gamma } [/latex], as the line integral of the velocity around a closed line, or a loop, as:[latex]\mathrm{\Gamma }\mathrm{=}\oint{v_{\theta }ds} [/latex]
(5.11)
This has units of [latex]m{}^{2}/s[/latex] in SI units. Inserting the expression for [latex] v_{\theta }[/latex] above and noting that [latex]ds=rd\theta[/latex] with the integration carried out between 0 and 2[latex]{\pi}[/latex], we get:[latex]\mathrm{\Gamma }{=}\oint{v_{\theta }ds}=2\pi C [/latex]
(5.12)
The result is that the circulation within the vortex is a constant that can be determined knowing [latex]{C}[/latex]. We can replace C with [latex]{rv}_{\theta }[/latex] and solve for the velocity as:[latex]v_{\theta }=\frac{\mathrm{\Gamma }}{2\pi r}=\frac{{\mu }_v}{r} [/latex]
(5.13)
where we define the "vortex strength" as:[latex]{\mu }_v=\frac{\mathrm{\Gamma }}{2\pi } [/latex]
(5.14)
[latex] \mathrm{\Gamma }\mathrm{=}\oint{v_{\theta }ds}=\int{\boldsymbol{\omega }dA} [/latex]
(5.15)
Consequently the circulation is the area integral of the vorticity in some selected area within the flow. Note that we are dealing with a two dimensional flow and dA is within the plane of the flow so the vorticity is aligned along a vector out of the plane. Now we have a bit of a dilemma. A free vortex has finite values of velocity [latex]v_{\theta }[/latex] that result is a certain circulation. This circulation is proportional to the vorticity within the flow. But we have said that the flow is irrotational, which means that the vorticity is zero. How can the circulation be nonzero while the vorticity is zero? The answer is that the vorticity (that drives the circular velocity) is concentrated at the center of the circle. Away from the center of the circle we know the vorticity is zero since the velocity we are using, [latex]v_\theta[/latex], was based on the vorticity being zero. If one calculates the circulation about a circular loop of radius [latex]r{}_{1}[/latex] (an area of [latex]\pi r{}_{1}{}^{2}[/latex]) and then repeats this calculation for a larger area of radius [latex]r{}_{2}[/latex], the result will be the same. This shows that in the area between [latex]r{}_{1}[/latex] and [latex]r{}_{2}[/latex] there is no added vorticity since the circulation remains the same. This can be done for any arbitrarily small radius [latex]r{}_{1}[/latex], showing that in the limit of r going to zero there is no vorticity within the flow except at r=0. A free vortex has a concentration of vorticity at r = 0. Also, the velocity must go to infinity, a velocity singularity. Interestingly this says that we can have isolated vorticity concentrations within an otherwise irrotational flow field.[latex]\psi ={\mu }_s{\theta }_1{\mu }_s{\theta }_2 [/latex]
(5.16)
where it is assumed the source and sink have equal but opposite strength and the angles are shown in Fig. 5.3. Similarly for the velocity potential:[latex]\phi ={\mu }_sln\ r_2{\mu }_sln\ r_1=\ {\mu }_sln\ \frac{r_2}{r_1}[/latex]
(5.17)
[latex]\psi ={\mu }_s\ arctan\ \left(\frac{2ax_2}{{x_1}^2a^2+{x_2}^2}\right) [/latex]
(5.18)
For example, the streamfunction at [latex]P = (3,4)[/latex] for {a} = 2 is [latex] \psi ={\mu }_s(37.3^o)[/latex]. At this point we have a source and sink separated by a distance of 2a. Extending this we can add a uniform flow in the positive [latex]{x{}_{1}}[/latex] direction to the combined source sink flow to obtain the following:[latex]\psi =\ Ux_2{\mu }_s\ arctan\ \left(\frac{2ax_2}{{x_1}^2a^2+{x_2}^2}\right)=Ur\ sin\ \theta +{\mu }_s\left({\theta }_1{\theta }_2\right) [/latex]
(5.19)
where the last term is the result written in cylindrical coordinates. If we set the streamfunction equal to some constant we can plot the associated streamline from Eqn. (5.19). In particular setting [latex]\psi =0 [/latex] an oval results as shown in Fig. 5.3 (b). This is known as the Rankine Oval. The characteristics of this oval can be adjusted by inserting different values for [latex]{\mu}_{s}[/latex] for a given [latex]{U}[/latex]. The characteristics of the oval are: Major axis [latex](x{}_{2}=0) =\ 2a{\left(1+\frac{2{\mu }_s}{a}\right)}^{1/2}[/latex]Minor axis [latex](x{}_{1}=0) =2a\left(cot\ \frac{h/2a}{2{\mu }_s/aU}\right)[/latex] (where [latex]{h }[/latex]= minor axis)
(5.20)
It is possible then to model the flow over an oval surface with an approach velocity of U by Eqn. (5.19). The velocity is determined by taking the derivatives of [latex]\psi[/latex] relative to [latex]x_1[/latex], and [latex]x_2[/latex]. The geometry of the oval can be adjusted by varying the strength of the source/sink as well as their locations, [latex]{a}[/latex].[latex]\psi ={\mu }_d\frac{x_2}{{x_1}^2+{x_2}^2} [/latex]
(5.21)
where [latex]{\mu}_{d}{}_{\ }= 2C[/latex] a constant that determines the "strength of the doublet". This equation can be rearranged: [latex]{x_1}^2+{\left(x_2+\frac{{\mu }_d}{2\psi }\right)}^2={\left(\frac{{\mu }_d}{2\psi }\right)}^2[/latex] Noting that for a constant value of [latex]\psi[/latex] the coefficient on the right inside the parenthesis, and the same term added to [latex]{x{}_{2}}[/latex] is a constant resulting in an equation of a circle. The center of the circle changes with changing values of [latex]\psi[/latex] as does the radius of the circle. The result is a series of streamlines for selected values of [latex]\psi[/latex], each of which is a circle centered along the [latex]{x{}_{2}}[/latex] axis as is shown in Fig. 5.4. Notice that the radius of the circle is [latex]\frac{{\mu }_d}{2\psi }[/latex]. [caption id="attachment_302" align="aligncenter" width="277"] Fig. 5.4 Streamlines associated with a doublet  a source and sink with a distance a approaching zero[/caption][latex]\psi =Ux_2\frac{{\mu }_dx_2}{{x_1}^2a^2+{x_2}^2} [/latex]
(5.22)
This equation can be recast in cylindrical coordinates where [latex]{r{}^{2} = x{}_{1}{}^{2}+x{}_{2}{}^{2}}[/latex], [latex]{x{}_{1} = r cos \theta}[/latex] and [latex]{x{}_{2} = r sin \theta}[/latex] and we define [latex]a\ =\ {\left(\frac{{\mu }_d}{U}\right)}^{1/2}[/latex], a constant for a given value of [latex]{U}[/latex], the result is:[latex]\psi =Ur\ \ sin\theta Ua^2\frac{sin\ \theta }{r}=Ur\ sin\theta \left(1\frac{a^2}{r^2}\right) [/latex]
(5.23)
It can be seen that in the limit of large values of [latex]{r}[/latex] the flow reverts back to uniform flow and the contribution for the doublet goes to zero. If we let [latex]\psi =0[/latex] then [latex]{r }= {a}[/latex], a constant. That is to say, the streamfunction of [latex]\psi =0[/latex] is a circle of radius a. So we have now constructed the flow field for uniform flow over a cylinder of radius [latex]{a}.[/latex] From this the velocity components in the [latex]{r,\theta}[/latex] coordinates are found to be: [latex]v_r=\frac{\partial \psi }{r\partial \theta }=U\ cos\theta \ \left(1\frac{a^2}{r^2}\right) [/latex][latex]v_{\theta }=\frac{\partial \psi }{\partial r}=U\sin\theta \ \left(1+\frac{a^2}{r^2}\right)[/latex]
(5.24)
and the velocity vector is [latex]{V{}^{2}}({r,\theta}) = {v_r}^2+{v_{\theta }}^2[/latex]. The streamline distribution is shown in Fig. 5.5. We are only interested in the flow on the outside of the circle, which represents a cylinder. It is important to note that at [latex]{r }= {a}[/latex], the velocity is not zero and depends on [latex]{\theta}[/latex]. So at different positions around the cylinder surface the velocity will change. The largest velocity occurs at [latex]{\theta} = 90{}^{o}[/latex] and [latex]270{}^{o}[/latex], which represents the top and bottom of the cylinder. This is where the streamlines converge the most indicative of high velocity flow. Also at [latex]{\theta} = 0{}^{o}[/latex] and [latex]180{}^{o}[/latex] the velocity is zero. These are stagnation points on the cylinder. Notice also that the streamlines are symmetric about the [latex]{x{}_{1}}[/latex] and [latex]{x{}_{2}}[/latex] axes. This has important implications on the forces that exist on the cylinder caused by the flow. [caption id="attachment_303" align="aligncenter" width="458"] Fig. 5.5 Streamline for flow over a cylinder formed from the superposition of a doublet with uniform flow; note stagnation on either side of the cylinder[/caption][latex]\psi =\ Ur\ sin\ \theta \ {\mu }_d\ \frac{sin\ \theta }{r}{\mu }_vln\ r\ +C[/latex]
(5.25)
where [latex]{C }[/latex] is a constant that can be determined by setting the value of [latex]\psi[/latex] at some point in the flow, as shown below. The coordinate system used here has its origin at the center of the circle generated by the uniform flow and doublet. The added vortex does not add a radial component of velocity, since its flow streamlines are all circles. The result is adding a [latex]{\theta}[/latex] component of velocity throughout the flow that depends on the radial location. In Eqn. (5.25) the direction of the added circumferential flow is clockwise (negative). To summarize the strengths of the elements we have [latex]{U}, {\mu }_d=Ua^2[/latex] and [latex]{\mu }_v=\ \frac{\mathrm{\Gamma }}{2\pi }[/latex]. We can think of these as adjustable parameters to the flow field.Rewriting Eqn. (5.25) by combining the first and second term as was done for flow over a cylinder and then setting [latex]{r = a}[/latex] , inserting the expression for [latex]{\mu }_d[/latex], as well as setting [latex]\psi =0[/latex], we have:
[latex]C=\mu_vln[/latex] [latex]a[/latex] This then sets [latex]\psi=0[/latex] on the circle with radius [latex]a[/latex]. Inserting this value of C in Eqn. (5.25) we obtain:[latex]\psi =\ U\ sin\ \theta \ \left(r\frac{a^2}{r}\right){\mu }_vln\ \frac{r\ }{a} [/latex]
(5.26)
This represents uniform flow over a rotating cylinder as shown in Fig. 5.7, where the streamline representing the cylinder has [latex]{\psi =0}[/latex]. Notice that since we have a clockwise (negative) circulation the cylinder is rotating clockwise. The velocity on the surface of the cylinder now has a contribution caused by the vortex in addition to the velocity found for a nonrotating cylinder. We see that the stagnation points have moved away from along the [latex]x_1[/latex] axis. Since we know that the velocity is zero at the stagnation points we can solve for their location (note that [latex]v_r=0[/latex] everywhere on the cylinder).[latex]v_{\theta }(r=a)=\frac{\partial \psi }{\partial r}=2Usin\ \theta \ +\frac{\mathrm{\Gamma }}{2\pi a} [/latex]
(5.27)
[latex]P(r=a)\ =P_{\infty }+\frac{1}{2}\rho U^2\frac{1}{2}\rho {v_{\theta }}^2\rho gh[/latex]
(5.28)
where [latex]{h }[/latex] is the local height on the circle above the centerline (datum) [latex]h=a\ sin\ \theta .[/latex] We insert Eqn (5.27) for the velocity on the surface, [latex]v_\theta[/latex], into (5.28) to obtain the expression for the local surface pressure. The integral of the pressure around the circle then determines its net force. First, we find the force component in the [latex]x_2[/latex] direction, this is denoted as the "lift force" per distance into the page, [latex]{F'}_L[/latex]: [latex]{F'}_L=\ \oint{P(r=a)\ \ sin\ \theta \ a\ d\theta }[/latex] Inserting the expression for [latex]{P }[/latex] and noting that integrals of odd powers of the sine functions around the entire circle are zero we are left with:[latex]{F'}_L=\ \rho U\mathrm{\Gamma }\ +\ \rho g\pi a^2[/latex]
(5.29)
The last term represents the net body force by the fluid on the volume of the circle (per distance into the page) or the buoyancy force. The lift force is usually defined without the buoyancy force included so we write:[latex]{F'}_L=\ \rho U\mathrm{\Gamma } [/latex]
(5.30)
A few things should be stated about this result. First, recall that the viscous force is not included so this is only due to the pressure distribution. Also, when the circulation is zero there is no lift force (there is still a buoyancy force however). We can conclude that the circulation provides the means to create asymmetric conditions around the circle so that a net pressure force occurs. This expression is often called the KuttaJoukowsky Law who showed this equation to hold for other shapes as well, and is often used in aerodynamics to determine the lift force on two dimensional airfoils based on the circulation associated with the flow around the wing. The sign convention is such that a counterclockwise rotation results in a downward force, and a clockwise rotation results in a upward force for flow along the positive [latex]{x{}_{1}}[/latex] axis. It is surprisingly accurate for real flows considering the restrictions on its application. This tends to indicate that viscous forces are small at best. It is only accurate for flows that have not separated from the object surface. We will discuss separation when we get into viscous flow affects. The KuttaJoukowsky Law (or theorem) as applied to an airfoil requires what is known as the Kutta condition. This is a condition on the flow at the trailing edge that says that the flow exits the airfoil on the top and bottom surfaces of the airfoil with equal velocity and pressure. This implies that the flow does not tend to wrap around the back end or trailing edge of the airfoil and is the boundary condition that allows the calculation of the lift using Eqn. (5.30). If interested see the web site: https://en.wikipedia.org/wiki/Kutta_condition[latex]F=\emptyset+i\psi[/latex]
(5.31)
The derivative of [latex]F[/latex] in terms of [latex]z[/latex] is defined as:[latex]W\left(z\right)=\frac{dF}{dz}[/latex]
(5.32)
Also: [latex]\frac{dF}{dz}=\frac{dF}{dz}\frac{\partial z}{\partial x}=\frac{\partial F}{\partial x}[/latex] Consequently, inserting the definition of [latex]F[/latex] we have[latex]W\left(z\right)=\frac{\partial\emptyset}{\partial x}+i\frac{\partial\psi}{\partial x}=u−iv[/latex]
(5.33)
This result shows that [latex]W(z)[/latex] represents the “complex velocity” of the flow and is determined by the real velocity components [latex]u,v[/latex]. It is often advantageous to use cylindrical coordinates, expressed as [latex]r,{\theta}[/latex] in two dimensions. The transformation from [latex]x,y[/latex] to [latex]r,{\theta}[/latex] is: [latex]u=u_rcos\theta−u_\theta sin\theta[/latex] [latex]v=u_rsin\theta+u_\theta cos\theta[/latex] The expression for the complex velocity is then (by direct substitution and using the identity of [latex]\cos{\theta}−isin\ \theta=e^{−i\theta}[/latex] :[latex]W=\left(u_r−iu_\theta\right)e^{−i\theta}[/latex]
(5.34)
The use of the complex variable representation in terms of [latex]F[/latex] can then be converted back into the physical space velocity components, [latex]u,v[/latex], through its derivative relative to the variable [latex]z[/latex]. The magnitude of the local velocity vector is found from [latex]W[/latex] by taking the square root of [latex]WW*[/latex], where [latex]W*[/latex] is the complex conjugate of [latex]W[/latex]. Once the velocity is determined then the pressure field is directly determined from the Bernoulli Equation. Some examples will help illustrate the use of complex variables to represent rather simple flows.[latex]F=Az^n=Ar^n(\cos{n\theta}+i\sin{n\theta)}[/latex]
(5.35)
Where “[latex]n[/latex]” represents a parameter to be specified, usually as a constant for a given flow, and “[latex]A[/latex]” is a constant for a specific flow, and shown below to be proportional to the velocity far from the sector, representative of the flow into the sector. An example of this class of flows is uniform flow where [latex]A=U[/latex] and [latex]n=1[/latex]. Note that uniform flow can be thought of as flow over a flat surface that is parallel to the uniform flow direction (there are no viscous forces so the no slip boundary condition does not hold). The flat surface can be thought of as a sector with an intersection of two lines with the lines being parallel. Or, the angle between the two lines is [latex]\pi[/latex]. The general expression above corresponds to the following expression for velocity potential and streamfunction: [latex]\emptyset=Ar^n\cos{(n\theta)}[/latex] [latex]\psi=Ar^n\sin(n\theta)[/latex] Using these two expressions we can easily locate lines of constant streamfunction can be found, and in particular when [latex]\psi=0[/latex]: when [latex]\theta=0[/latex] and [latex]\pi/n[/latex]. The flow between these two radial lines represents flow in a “sector”, as seen in the Fig. 5.8, below. [caption id="attachment_733" align="aligncenter" width="465"] Fig. 5.8 Twodimensional flow in a sector of arbitrary angle between two straight lines; the streamfunction is chosen to be zero along these two lines.[/caption] Using the definition of the complex velocity potential, [latex]W[/latex]: [latex]W=\frac{dF}{dz}=nAz^{n−1}=nAr^{n−1}e^{i(n−1)\theta}[/latex] Expanding the exponential reveals that the velocity components in [latex]r,\theta[/latex] coordinates are: [latex]u_r=nAr^{n−1}\cos(n\theta)[/latex] [latex]u_\theta=nAr^{n−1}\sin(n\theta)[/latex] By inserting different values for the parameter “[latex]n[/latex]” different flows can be simulated. Some examples are given below. [latex]n = 1:[/latex] [latex]u=A[/latex] and [latex]v = 0[/latex] (uniform flow) [latex]n = 2:[/latex] [latex]u=2Ar\cos{(\theta)}[/latex] [latex]v=2Ar\sin{(\theta)}[/latex] or in Cartesian coordinates: [latex]u=2Ax;\ v=2Ay[/latex] [latex]\psi=−Ar^2\sin{\left(2\theta\right)}=−Axy[/latex] or for a constant streamfunction (along a streamline) we have: [latex]xy = C[/latex] (where [latex]C[/latex] is a constant determined by the value of [latex]\psi[/latex] at a given location). This latter flow, for a range of values of [latex]C[/latex] yields flow into a right angled corner where the location ([latex]0,0[/latex]) is the position of the corner. Note that [latex]\psi=0[/latex] along the axes. The reader is encouraged to plot this function for different values of streamfunction to visualize the flow. The corresponding flow into a corner of arbitrary angle, [latex]\alpha[/latex], can be expressed as: [latex]F=Az^{\pi / \alpha}[/latex] Note that [latex]Az^n[/latex] is equivalent to [latex]Ar^n [cos(n\theta)+i \enspace sin(n\theta)][/latex] which can then be restated in terms of [latex]\pi / \alpha=n[/latex], so that for [latex]\alpha=\pi /3[/latex] we have [latex]n=3[/latex]. This then relates [latex]n[/latex] to the turning angle of the corner flow. For instance, if [latex]n=3/2[/latex] we have an angle of [latex]2\pi /3[/latex] which is greater than 90^{o}. We can also have flow over a flat plate parallel with the flow if [latex]n=1[/latex].[latex]F=Ae^{−i\beta}z^n=Ar^ne^{i(n\theta−\beta)}[/latex]
(5.36)
Noting that the real part of F is the velocity potential and the imaginary part is the streamfunction we have: [latex]\emptyset=Ar^n\cos(n\theta−\beta)[/latex] [latex]\psi=Ar^n\sin(n\theta−\beta)[/latex] This streamfunction expression shows us that when using the complex representation for the constant multiplying [latex]z^n[/latex] there is introduced a rotation of the streamlines through angle [latex]\beta[/latex] using the same coordinate system when compared to the streamfunction when the coefficient is real (not complex).(6.1)
where the integral is over the surface area of the flow (assuming two dimensional flow), S. This equation indicates that to solve for the velocity potential we must evaluate the integral on the flow boundaries (both the solid surface and infinitely far away). This can be shown to be a result of the application of the divergence theorem that says:$latex \int{n_ju_jdS}=\int{\frac{\partial u_i}{\partial x_i}d\forall } $
(6.2)
where $latex {u{}_{i}}$ is the velocity vector. This states that the rate of expansion of a volume of fluid, given by the divergence of the velocity vector, within a given volume equals the flux through the volume boundaries. If the fluid is incompressible then the right side is zero and the net flux through the entire boundary is zero (this is the general argument for the conservation of mass of steady flow). All this is nice and can be a powerful tool to find $latex \phi$ and therefore the velocity in the flow. But we really don't need to find the potential of the entire flow, what our goal is, is to set up a flow that we know satisfies the velocity boundary condition at the surface of the object (with no velocity component across the surface) and far from the surface where the velocity is known to be the freestream velocity. Once this is established the force can be found. That is to say we only need to evaluate the surface velocity and then the pressure on the surface. The general approach is to select a ``grid'' which is a series of ``panels'' that form the surface. Here we take the panels as straight flat surfaces arranged over the real surface. In the limit of very small panels the constructed surface will simulate the actual surface. On each panel we place a distribution of flow elements (like sources, sinks, vortices, etc) that when combined together will result in a flow field that will satisfy the surface boundary condition. There are lots of ways to identify which elements to use and how they may be distributed on the panels. Here we will use vortex elements, with one placed on each panel. The net flow is the result of superposition of the flow set up by each vortex on each panel element. So at each point in the field we add together the flow caused by all of the panel elements using the superposition rule. The panels that are far away from a given point will have less and less influence on the flow because the strength of the flow caused by a flow element decreases with distance from the element origin. For instance for a ``source'' the velocity decreases with increasing radial position because the flow is spreading out away from the source. But the influence never really goes to zero. [caption id="attachment_261" align="aligncenter" width="300"] Fig. 6.1 Illustration of a panel geometry; any three dimensional shape can be constructed; shown here is a surface of what could be a three dimensional object such as an entire airplane.[/caption] In placing a series of panels over the surface we first need to specify the size of each panel. We place a vortex of some strength some where on the panel (whose location is shown below) and we must identify points on the surface where we want to make sure that the velocity is zero across the surface. To be clear, individual points on each panel are used to evaluate the element (vortex) flow field it needs its own origin, or coordinate system, to write an equation for the flow generated by this element. We also only pick a point on the panel to check to make sure that the net sum of contributions from all elements results in zero flow across the surface. The fact that we only satisfy the condition at one point on each panel will be satisfactory if the panels are made to be reasonably small. These points are called ``collocation points'' on each panel. In the end each panel will have coordinates that define its location on the surface, coordinates for the element location on the surface and coordinates for the collocation point. [caption id="attachment_310" align="aligncenter" width="372"] Fig. 6.2 Geometry of panels over a two dimensional wing showing individual vortices, $latex \Gamma_j$ associated with each panel; $latex {n}$ is the outward normal for each panel; the wing is at an angle of attack to the freestream of $latex \alpha$.[/caption] The method of solution for the force on the object is the determination of the magnitudes or strengths of the elements on each panel. Once we have this distribution of strengths we can calculate the total lift on the surface that results from all of these elements. Recall that the lift experienced by a surface is really the component of the pressure at the surface integrated over the surface area in a direction normal to the approach velocity vector of the flow  it is not necessarily vertical, but normal to the freestream velocity. For illustration we are going to use flat plate panels with vortex elements, one per panel. To get an idea of how to define the vortex origins on each panel we can examine flow over a single flat plate at some angle of attack to the freestream velocity vector, $latex \alpha $. First we define the ``center of pressure'', $latex {x{}_{cp}}$. This is the location on the surface where the resultant lift force caused by the distributed load on the surface acts such that there is no net moment on the surface (see F.M. White, Chapt. 8, for some details on this). In general, if the moment on the surface about the leading edge is $latex {Mo}$, and $latex {L}$ is the lift force then using $latex {x}$ as the coordinate measured along the flat surface from the leading edge then: $latex {x{}_{cp\ }= Mo/L}$. Consider now flow over the same flat surface, the lift we have seen is $latex L= \rho {U} \Gamma$ note that we are going to define $latex \Gamma$ as positive clockwise  this is opposite to what we have done previously, but it helps get rid of some negative signs, this is not necessary, but convenient. To find out the value of $latex {x{}_{cp}}$ we use the theoretical evaluation of the lift force based on ``thin airfoil theory''. The interested reader can refer to F.M. White Chapt. 8 for details on this. First, assuming a two dimensional flow, we assume that there is some continuous distribution of vorticity over the surface (in contrast to discrete vortex elements), such that there is a local circulation per unit length of surface, $latex {\gamma (x)}$. The surface integral of this distribution results in the net total circulation: $latex \mathrm{\Gamma }=\int^c_0{\gamma dx}$ where ``c'' is the length of the flat surface, which we call the ``chord''. The Kutta condition is imposed where we want the flow to leave the surface flowing parallel with the back (trailing) edge. If this is a surface like an airfoil we want both the top and bottom flows to leave smoothly from the trailing edge. For this to happen then we don't expect there to be a pressure difference between the top and bottom of the surface at the trailing edge, which would cause the flow to deflect up or down. If there is no velocity or pressure difference across the flow at the trailing edge the local circulation at the trailing edge, $latex {(x=c)}$, must be zero. So now we have a boundary condition for the integration of $latex {\gamma (x)}$ to get the total circulation $latex \Gamma $, that is $latex {\gamma = 0}$ at $latex {x=c}$. The needed distribution has been figured out for a single flat surface at an angle of attack of $latex \alpha$, the angle between the oncoming flow and the surface. The function that works is: $latex \gamma =2 U \sin \alpha \ {(\frac{c}{x}1)}^{1/2}$. The integral of this expression for $latex {\gamma(x)}$ yields the value of $latex \Gamma $. Once the circulation is known then the lift, $latex {L,}$ can be found. Using this lift, and also using the above definition of $latex {x{}_{cp}}$, it turns out that the ``center of pressure occurs at $latex {x{}_{cp} = c/4}$. So by assuming the lift force occurs at $latex {x{}_{cp}}$ then there is no moment caused by the distributed pressure force. The lift force, $latex {L}$, as calculated as stated above, can be made nondimensional by divided by $latex {1/2 \rho U{}^{2}b}$, where b is the span of the surface (in and out of the page). The result is an expression for the lift coefficient:$latex C{}_{L} = 2 \pi \sin \alpha $
(6.3)
This result is the thin airfoil theory result for the lift on a surface at angle of attack, $latex \alpha $. So at this point, we take this result for a flat surface and place vortex elements at the center of pressure for each flat element, which will be at a point ¼ of the distance along the panel from the beginning edge of the panel. Although this is not absolutely required for the panel method it is a convenient choice. As panel size is reduced smaller and smaller the impact of this choice becomes small. The next step is to find the location where we want to impose the condition of zero velocity crossing the surface for each panel. For this we are using a coordinate system for the panel to be ($latex {x,z}$) where $latex {x}$ is, again, along the flat surface from the panel leading edge and $latex {z}$ is normal (upward) from the surface. To find these collocation points, we note that the the velocity set up by a vortex element at location $latex {r}$ from the origin is $latex v_{\theta }=\mathrm{\Gamma}/2\pi r$. Writing this in Cartesian coordinates where at the panel surface $latex v_{\theta }=\ w$ (where $latex {w}$ is the $latex {z}$ component of velocity) is:$latex w=\frac{\Gamma }{2\pi }\frac{(xx_o)}{{(zz_o)}^2+{(xx_o)}^2} $
(6.4)
where $latex (x_o,z_o)$ is taken as the center of pressure (where the vortex element is located). Also, since the collocation point is on the surface then both $latex {z}$ and $latex z{}_{o}$ are zero. To get the total normal component velocity at the panel surface we add together $latex {w}$ being the velocity normal to the surface induced by the vortex, with the component of the freestream flow normal to the surface which is at angle $latex {\alpha}$ to the surface, $latex U \sin \alpha$, and set this to zero:$latex w+U\ sin\ \alpha =0 $
(6.5)
We next insert Eqn. (6.4) for $latex {w}$, with $latex {z{}_{o} = 0}$ and $latex {x{}_{0} = c/4}$. And we define the location at which we evaluate the zero normal velocity component as the collocation point as: $latex {x= kc}$. Finally we solve for ``$latex {k}$'' to be $latex {k=3/4}$. This says that the $latex {x}$ location along any given panel where the boundary condition is to be satisfied is at $latex {x=3 \ c/4}$. Recall that the vortex element is located at $latex {x_{o} = c/4 }$. In summary, once the panels are set up on the surface, each has its own angle of attack, $latex \alpha_{\rho} $, depending on the local orientation of the surface. Each panel has a vortex element located at its own center of pressure and the normal velocity must be equal to zero at $latex {x=3 \ c/4} $ measured from the front of each panel. An equation can be written for the zero velocity condition at the collocation point of a panel by taking the sum of all of the contributions from all of the panels, each with a different vortex strength $latex \Gamma_{\rho}$, and setting this sum to zero. This will result in N equations (one for each panel) with N unknowns for $latex \Gamma_{\rho}$ for each panel. We can then solve for all $latex \Gamma_{\rho}$ using this set of equations. [caption id="attachment_564" align="aligncenter" width="432"] Fig. 6.3: A simple example of two panels along a flat surface of length c and angle of attack of $latex \alpha$.[/caption] Now we can set this up for a simple example to help pull all of this together, see Fig. 6.3. We use a thin flat surface representing a flat airfoil at an angle of attack of $latex \alpha $, with a uniform approaching flow of U. We use two panels for illustration. Note that the entire length of the two panels together is ``c'', NOT the length of each panel. Also, $latex {x} $ and $latex {z,}$ used previously for a single panel, are along and normal to a line from the start to the end of the surface, respectively. This provides a coordinate system for the ``system'' of panels. For this simple case the angle of attack for both panels are the same, but we allow individual panel vortex strengths, $latex {\Gamma}_{1}$ and $latex {\Gamma}_{2}$. In determining the vortex location for each panel, using the overall coordinate system we get (c/8,0) and (5c/8,0) for panel 1 and panel 2, respectively. The collocation points are (3c/8,0) and (7c/8,0), respectively. [caption id="attachment_566" align="aligncenter" width="646"] Fig. 6.4: Illustration of the coordinate system using the chord as the x axis; this is the coordinate system typically used to define $latex a_{ij}$.[/caption] Note that the each outward normal, $latex n_{i}$, points in the same direction for this case since both panels have the same orientation, but in general each panel could be at some angle $latex \beta_{i}$. Fig. 6.4 shows the general definition of the panel angle to the chord, the latter being along the $latex {x}$ axis. The outward normal to the panel is given as $latex n_{i}$. For the general case shown in Fig. 6.4 the outward normal is at an angle $latex \beta_{i}$ to the $latex {z}$ axis. The general equation for the boundary condition is Eqn. (6.5), which we write slightly differently as:$latex v_i\cdot n_i=U\cdot n_i$
(6.6)
Where [latex]v_i[/latex] is the induced velocity by the set of vorticies evaluated at each panel [latex]j[/latex] and [latex]n_i[/latex] is the panel outward normal. The left hand side is the dot product of the induced velocity with the outward normal. The right hand side is the component of the freestream velocity aligned with the outward normal. These two balance each other to yield zero velocity crossing the panel. The induced flow from all of the panel elements needs to be added together to get the value of $latex v_{i}$ in the above equation. The general form for the velocity, which has components identified as $latex {u,w}$,at any point in the flow with a vortex element located at ($latex x_o$, $latex z_o$) is obtained from above and written here as:$latex u=\frac{\mathrm{\Gamma }}{2\pi }\frac{(zz_o)}{{(xx_o)}^2+{(zz_o)}^2} $
(6.7)
$latex w=\frac{\mathrm{}\mathrm{\Gamma }}{2\pi }\frac{(xx_o)}{{(xx_o)}^2+{(zz_o)}^2} $
(6.8)
where we have transformed $latex v_{\theta }$ into the Cartesian components $latex {u,w}$. Writing this as a matrix we obtain a general expression for the induced velocity anywhere in the flow [latex]\left ( x, z \right )[/latex]:$latex \binom{u}{w}= \frac{\Gamma}{2 \pi r^2} \Big ( \frac{01}{10} \Big ) \binom{xx_o}{zz_o}$
(6.9)
Note here that $latex xx_o$ is the distance along the x axis from a vortex element to the point of interest where the velocity is being evaluated. Similarly for $latex {{z}z_o}$. Since the values for $latex \Gamma$ are not known, it is convenient to write this set of two equations for a ``unit value'' of $latex \Gamma $, that is $latex \Gamma=1$. Keep in mind that the flow caused by each vortex when evaluated on a surface results only in a component normal to the surface. We also are going to use values of $latex x_{o}$ and $latex {x}$ based on the location of the vortices and collocation points, respectively. Consequently, for our two panel example, where each panel is along the chord line, or the $latex {x}$ axis, we can write the velocity components given below, however here these velocities assume $latex {\Gamma}=1$ in the above equations. The notation used below is that the first subscript, $latex {i}$, represents the collocation point of interest at some panel, and the second subscript, $latex {j}$, identifies the vortex element at some panel that induces the velocity at the $latex {i}$ collocation point. That is to say these subscripts do NOT represent components of vectors in space in the conventional manner, and here [latex]u[/latex] and [latex]w[/latex] are the [latex]x[/latex] and [latex]z[/latex] component of the velocity at the surface.$latex (u_{11}, w_{11}) = \left (0, \frac{1}{\frac{2 \pi c}4} \right ) = (0, \frac{2}{\pi c})$
(6.10)
$latex (u_{21}, w_{21}) = \left (0, \frac{1}{\frac{2 \pi 3 c}4} \right ) = (0, \frac{2}{3 \pi c})$
(6.11)
$latex (u_{12}, w_{12}) = \left ( 0, \frac{2}{\pi{c}} \right )$
(6.12)
$latex (u_{22}, w_{22}) = \left ( 0,  \frac{2}{\pi{c}} \right )$
(6.13)
For instance, $latex u_{12}$ is the $latex {x }$ component of velocity at the collocation point of panel 1 caused by the vortex element on panel 2. Also, $latex w_{12}$ is the z component of velocity at collacation point 1 caused by the vortex element on panel 2. We can define, again for $latex \Gamma_{\rho} = {1}$, for each panel with vortex $latex {j}$, the following matrix $latex a_{ij}$:$latex a_{ij}={\boldsymbol{v}}_{ij}\cdot \boldsymbol{n} $
(6.14)
This represents the component of the induced velocity in the outward normal direction for a given panel [latex]\Gamma=1[/latex]. We need to be careful here with this notation, $latex {v}_{ij}$ in the equation above is the velocity vector caused by the vorticies with unit value circulation with the subscripts: $latex {i}$ represents the collocation point and $latex {j}$ the vortex location inducing that flow. Also, the dot product with $latex {n}$ (the outward normal for each panel), gives the projection of $latex {v}$ normal to the surface. For example for the above example [latex]n[/latex] is only in the [latex]z[/latex] direction.$latex a_{11} =\left ( 0,  \frac{2}{\pi c} \right ) \cdot (0.1)=\frac{2}{\pi c}$
(6.15)
so the elements $latex a _{ij}$ are quantities representing a measure of the velocity contribution normal to each panel from each vortex $latex {j}$ at each collocation point $latex {i}$ for a $latex \Gamma_{j}=1$ at each panel. Once we have identified the components $latex a_{ij}$, which, as shown above for the two panel example, are determined only by geometry, it is possible to use the impervious boundary condition to find each of the values for $latex \Gamma _{j}$. We have a general equation for our two element surface for panels $latex {i}=1,2$ that says the sum of the outward normal velocity at each panel with the freestream component normal to that panel must be zero. The outward normal velocity due to summation of the induced velocity over both panels is for some unknown distribution of circulations, $latex \Gamma_{j}$:$latex \sum^2_{j=1}{a_{ij\ }{\mathrm{\Gamma }}_j} $
(6.16)
So this expression is balanced by the contribution from the free stream velocity: $latex U\cdot n$ (this is the component of U normal to the surface). The final system of equations become:$latex \sum^2_{j=1}{a_{ij}{\mathrm{\Gamma }}_j=(U_x,U_z)\cdot (n_x,\ n_z)} $
(6.17)
Here the right hand side value of $latex {U}$ is represented as the vector ($latex {U}_{x},{U}_{z}$) indicating that $latex {U}$ is a vector with $latex {x}$ and $latex {z}$ components that depends on the angle of attack given by the chord line. The $latex {x,z}$ components of the outward normal will depend on the angle $latex {\beta}$ for each panel. The geometry of all of the panels fully determines $latex a_{ij}$, as shown above in the simple two panel example. Consequently, once the geometry is known the values of $latex a_{ij}$ and $latex (n_x,\ n_z)$ are all determined. Then knowing $latex {U}$ the set of Eqns. (6.17) provides the means to find all of the values of $latex {\Gamma_j} $, where the index represents each of the panel vortices. In the above two panel example there are two values of $latex {\Gamma}$. From the solution of $latex \Gamma_j $ we can find the lift on each panel, $latex \Delta L_j=\rho U \Gamma_j$ and then the total lift as the sum over all panels:$latex \sum^2_{j=1}{\Delta L_j}=C_L\left(1/2\rho U^2c\right)$
(6.18)
Recall that the lift is normal to the direction of the freestream velocity, so the lift of each panel, even though the panels may all be at different angles of attack, are all in the same direction. It is also possible at this point to calculate the velocity tangent to the surface at each panel. That is to say, instead of finding the component normal to the surface to satisfy the impermeable boundary condition find the component tangent to the panel, for each panel and combine this with the component of U tangent to the panel. The induced velocity tangent to the panel can be found using the dot product of the vortex induced velocity with the unit vector tangent to the panel. The tangent unit vector is given by its [latex]x[/latex] and [latex]z[/latex] components as [latex]cos\beta_i[/latex], [latex]sin\beta_i[/latex], respectively with [latex]\beta_i[/latex] being the angle of the panel relative to the chord line which is the [latex]x[/latex] axis, show in Fig 6.4. The freestream contribution is found from the value of U. This combined velocity yields the velocity vector of the flow at each panel and can be used in the Bernoulli equation to find the pressure at each panel. The Bernoulli equation is applied from a known upstream condition (known pressure and freestream velocity) to each panel in order to calculate the panel pressure.[latex]U\sin{\alpha}=\ \int_{0}^{c}{\frac{\gamma(x)}{2\pi(x_o−x)}dx}[/latex]
(6.19)
The above condition needs to be satisfied for the given distribution [latex]\gamma(x)[/latex] along the foil. However, we need a boundary condition at [latex]x = c[/latex] that specifies [latex]\gamma(c)[/latex]. This is done by introducing what is known as the Kutta Condition. This is based on the observed flow condition at the trailing edge of the foil, assuming the geometry is such that the flow from the top and bottom sides meet at a point (the trailing edge) and that flow exits the foil “smoothly” from the top and bottom surfaces such that there is no vorticity as the flow leaves the trailing edge. Physically this would occur if there is no rollup of the flow say from the bottom to the top surface. If the velocities on top and bottom at the trailing edge are identical, then the pressures are identical by the Bernoulli equation, and the flow has no tendency to swirl. So, this condition is expressed as [latex]\gamma(x=c) = 0[/latex]. It can then be shown that the integration in Eqn. (6.19) is found to satisfy the Kutta Condition for the following distribution of the vortex sheet:[latex]\gamma(x) = 2U \sin \alpha (c/x 1) [/latex]
(6.20)
(See the extra comments on the Kutta condition below) With this distribution of [latex]\gamma(x)[/latex] the total circulation, [latex]\Gamma[/latex], can be determined. Note that this value of [latex]\Gamma[/latex] depends on the angle of attack, [latex]\alpha[/latex], and the free stream velocity, [latex]U[/latex]. We can now use the previously given result that the lift force, [latex]L[/latex] , is dependent on the circulation that was shown for a rotating cylinder in a freestream, [latex]L= \rho U \Gamma[/latex]. This is further discussed below and is found to be valid for airfoils in irrotational flow. The integration of [latex]\gamma(x))[/latex] given in Eqn. (6.19) results in the following Lift expression for an airfoil: [latex]L=2 \pi \sin \alpha \rho U^2 c[/latex] Or by defining the nondimensional lift coefficient as the lift force divided by [latex]1/2 \rho U^2[/latex]:[latex]C_L=2 \pi \sin \alpha[/latex]
(6.21)
This surprisingly simple expression is a major finding of the thin air foil theory in the description of the lift force as determined by the foil angle of attack. This is based on a twodimensional steady flow (or infinitely wide foil) with very little camber with viscous effects ignored. Some additional features and effects will be added to this equation later. Figure 6.7, shows the lift coefficient versus angle of attack. At large angles of attack the results of eqn. (6.21) are not valid. This is because the angle is so large that flow over the leading edge cannot negotiate the needed turn to stay “attached” to the foil surface. Rather the flow “separates” from following the surface contour. This results in a recirculating flow region between the streamlines of flow and the body surface. This can be seen in Fig. 4.5 illustrating flow going over the leading edge of a thin airfoil at a fairly large angle of attack. The primary limitation on the application of thin airfoil theory is the assumption of inviscid flow (which is actually has been shown to be a fairly valid approximation since all of the friction induced vorticity is accounted for with the vortex sheet). Also, the foil must be thin (maximum thickness on the order of ten percent of the chord length). This is not necessarily that restrictive since even airfoils with thicknesses of over 10% of the chord follow thin airfoil theory fairly well. With some modifications to the basic lift result of Eqn. (6.19) thickness effects, relatively large camber and threedimensional effects can be accounted for, as given below. A generic plot of the lift coefficient, [latex]C_L[/latex], versus angle of attack, a, is given in Fig. 6.7, as well as results for a specific airfoil, the NACA 0012 compared with Eqn. (6.21). Very good agreement between data and theory are shown up to an angle of attack of about 16^{o}. [caption id="attachment_791" align="aligncenter" width="300"] Fig. 6.7. Illustration of a typical Lift coefficient versus angle of attack curve, left, and comparison with Eq1n. (6.19) and data for a NACA 0012 air foil.[/caption]$latex \frac{dP}{dx}=\frac{2\pi R}{\pi R^2}{\tau }_w=\frac{2{\tau }_w}{R}$
(7.1)
This indicates that the pressure decreases along the flow direction proportional to the magnitude of the wall shear stress caused by friction between the fluid and pipe wall. In this equation the values of the wall shear stress is taken as a positive number. External flows have a similar noslip condition at a surface, the difference being that as one moves away from the surface any other boundaries are ``infinitely'' far away in the sense that these other surfaces are not influencing the flow. This is shown in Fig. 7.1 (b) showing that the extent of the friction layer increases along the flow direction. However, there is still wall friction and fluidfluid friction up to some point far enough from the surface, which depends on the particular location along the surface. We will discuss this in detail when we examine boundary layer flows. Recall that in previous sections we neglected viscous effects for external flows and the results tend to be very accurate, this is because for the flows we looked at other forces that come into play (pressure and body forces) dominate the flow. However, even in these flows viscous forces can be important in that they lead to ``flow separation'' where streamlines no longer follow the general shape of the object. When this occurs the potential flow solutions tend to break down, even though pressure forces dominate the overall forces. This is because the velocity distribution is greatly altered by the flow separation. This topic of flow separation is beyond the scope of this book but we will discuss it to a limited extent in the next chapter. The reader is referred to ``Viscous Fluid Flow'' by F.M. White for a discussion of this and some of the many models used in these flows, to predict wall friction effects. The velocity profiles shown in Fig. 7.1 indicates a variation of velocity from the boundary condition at the surface of noslip to a maximum at the centerline for pipe flow or reaching the freestream velocity far from the surface in external flows. This implies a distribution of internal viscous stress caused by the local velocity distribution since frictional forces depend on the relative velocity variation, or velocity derivative, within the flow field. In order to determine the magnitude of these forces we will first define the viscous stress tensor. [caption id="attachment_267" align="aligncenter" width="373"] Fig. 7.1 Illustration of velocity profiles (a) steady fully developed pipe flow resulting in a balance of pressure forces and wall friction forces, (b) external flow over a flat surface indicating boundary layer profile.[/caption]$latex \frac{1}{l}\left ( \frac{\bigtriangleup{l}}{\bigtriangleup{t}} \right )=\frac{\bigtriangleup{l/l}}{\bigtriangleup t}=\frac{\partial u_1}{\partial x_1}$
(7.2)
$latex \frac{\partial u_i}{\partial x_i}={\dot{e}}_{ii}$
(7.3)
$latex \frac{\Delta l/l}{\Delta t}=\frac{\partial u_1}{\partial x_2}$
(7.4)
We can interpret this as a rotation, or angular strain rate, which in this case is about the $latex x{}_{3}$ axis. Similarly there can be a rotation about [latex]x_3[/latex] caused by rate, $latex \frac{\partial u_2}{\partial x_1}$. This can be repeated to obtain angular strain rates about each of the other axes as well. These results provide the basis for the definition of the strain rate tensor, $latex {\dot{e}}_{ij}$:$latex {\dot{e}}_{ij}=\frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)$
(7.5)
Notice that when i=j the strain rate tensor consists of components of the dilation rate, for instance $latex \frac{\partial u_1}{\partial x_1}$, etc. When $latex i\neq j$ we see the addition of the two terms causing rotation about any given axis. For instance if i=1 and j=2 we obtain rotation about the $latex x{}_{3}$ axis as shown above. But also, if i=2 and j=1 we also obtain rotation about the $latex x{}_{3}$ axis. In fact when added to calculate $latex {\dot{e}}_{ij}$ these resultant two cases are identical in value. This implies that the strain rate tensor is symmetric, $latex {\dot{e}}_{ij}={\dot{e}}_{ji}$. Keep in mind that when the velocity derivative component is negative we have rotation in the opposite direction. The sign convention is that counterclockwise rotation is positive. The inclusion of ½ is a matter of convenience as we will see later, but also can be thought of as providing an average of the two terms. The general strain rate tensor, $latex {\dot{e}}_{ij}$, consists of all possible combinations of velocity spatial derivatives, and contains 9 elements. The diagonal terms are ``linear deformation rates'' and the off diagonal terms are the ``shearing deformation rates''. Notice that there are really only six definitive terms, since this is a symmetric tensor. Recall that vorticity, is a measure of rotation rate as well. We need to distinguish this from the strain rate tensor. Consider the case of $latex \ {\omega }_3$:$latex {\omega }_3=\left(\frac{\partial u_2}{\partial x_1}\frac{\partial u_1}{\partial x_2}\right)$
(7.6)
First we notice the sign difference for the second term in the parenthesis when compared with the strain rate tensor. What is the consequence of this? If rotation caused by $latex \frac{\partial u_2}{\partial x_1}$ is counter clockwise because this derivative is positive, and the rotation caused by $latex \frac{\partial u_1}{\partial x_2}$ is negative when this derivative value is positive, then applying this to the above definition of $latex {\omega }_3$ results in a net positive rotation (counter clockwise) about the $latex x{}_{3}$ axis. Or said another way, the vorticity is the sum of the two contributions of positive rotation rate about a given axis. Notice that vorticity has only three terms needed to define it (hence it is a vector). The diagonal terms are identically zero and the off diagonal terms are just of opposite sign (an antisymmetric tensor). It is then possible to have zero vorticity (the two terms in the vorticity definition cancel) but still have a positive angular deformation. Consider the case of $latex \frac{\partial u_2}{\partial x_1}=\frac{\partial u_1}{\partial x_2}$ where the rotation rate is zero but the net angular deformation rate is not, and equal to twice the value of each component. This addition is shown schematically for a fluid element in Fig. 7.2 (c) where the net rotation about the $latex x{}_{3}$ axis would be zero if $latex \dot{\alpha }=\ \dot{\beta }$ , but we see deformation of the fluid element. We do one last manipulation. Suppose we identify the velocity gradient tensor as $latex \frac{\partial u_i}{\partial x_j}$. This has nine elements of all of the partial derivatives of the three components in the three coordinate directions. This second order tensor can be written as:$latex \frac{\partial u_i}{\partial x_j}=\frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)+\frac{1}{2}\left(\frac{\partial u_i}{\partial x_j}\frac{\partial u_j}{\partial x_i}\right)$
(7.7)
The first term on the right hand side is the deformation rate tensor and the second term is ½ of the vorticity, which (when including the ½) is identified as the rotation rate tensor. The latter can be thought of as the average rotation about any given axis. So the velocity gradient tensor consists of two parts, that due to the strain rate and that due to the vorticity (or rotation rate). At this point we will construct a fairly simple model for the relationship between viscous stress and deformation rate to be as follows: $latex {\tau }_{ij}\sim {\dot{e}}_{ij}$ This implies that the rate of deformation that occurs in a fluid element is directly proportional to the stress applied to that element. If we make the assumption that there is a proportionality factor equal to 2 $latex \mu$, then we can write this as:$latex {\tau }_{ij}=\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)$
(7.8)
It is obvious why the factor 2 was introduced. Notice that the relationship does not include any of the rotational components of the velocity gradient tensor, only the deformation rate components. The physical interpretation is that viscous stress is related to the degree of deformation occurring in a fluid element, not its rotation rate. Looking at the normal stress components, when $latex {i=j}$, we see that these terms result in linear deformation. As noted above these are elongation or contraction rates of a fluid element as it flows over time. Adding all three dimensions we get the net volume change of the fluid element over time. Consequently, according to the model of Eqn. (7.8) the viscous forces resulting in normal viscous stress components are responsible for these linear deformation rates. These stresses are the normal, or diagonal, components of $latex {\tau }_{ij}$. The shearing contributions of $latex {\tau }_{ij}$ are the off diagonal terms and are related through Eqn. (7.8) to the off diagonal terms of the strain rate tensor, $latex {\dot{e}}_{ij}$. These terms result in an angular deformation rate that distorts the original shape of the fluid element over time. Again, according to eqn. (7.8) this angular deformation rate, due to shearing viscous stresses, has a magnitude determined by the parameter, $latex \mu$. The parameter $latex \mu$ is a fluid property and can take on many forms depending on the characteristics of the fluid itself. This property is the fluid dynamic viscosity and may be a function of the state of the fluid (say its pressure or temperature), but it may also be a function of the stresses applied to the fluid. A Newtonian fluid is one that can have state varying viscosity (say like an oil that has decreasing viscosity with increasing temperature) but has a constant value of viscosity with changing stress within the fluid. These fluids have a linear relationship between viscous stress and deformation rate. A dilatant fluid is described as shear thickening. That is to say as a shearing stress is applied to the fluid the viscosity increases. This implies that increasing stresses will result in decreasing deformation rates, consequently not having a linear relationship. Mixtures of cornstarch and water can behave as a dilatant fluid. In contrast, a fluid with a decreasing viscosity with increasing applied stress is called pseudoplastic, or shear thinning, fluids. Examples of this type of behavior is ketchup, some paints, ink, and others. Fig. 7.3 illustrates the stressstrain rate relationship for various types of fluids where the slope is a measure of viscosity. For the remainder of this book we will deal with Newtonian fluids. Keep in mind, however, functional relationships can be used for nonNewtonian fluid viscosities, whereas we will use constant values in the problems we examine. [caption id="attachment_270" align="aligncenter" width="356"] Fig. 7.3 Illustration of the relationship between stress and strain rate for different fluids; shear thickening fluids have increasing viscosity (slope) with increasing stress (dilatant fluids) whereas shear thinning fluids have decreasing viscosity with increasing stress (pseudoplastics); Newtonian fluids have a constant viscosity with strain rate, or a linear relationship.[/caption]$latex {\sigma }_{ij}=P{\delta }_{ij}+{\tau }_{ij}$
(7.9)
We use the Kronecker delta, $latex {\delta }_{ij}$, which is equal to one when i=j and otherwise equal to zero. This conveniently indicates that P is a normal compressive stress component. The negative sign indicates that pressure is compressive. An additional force is introduced to account for forces that result in volume changes of a given fluid element. This is modeled as follows. As noted above the rate change of fluid volume per original volume is given by $latex \frac{\partial u_k}{\partial x_k}$, where there is summation over index k. Assuming that the force required to change the volume is a compressive normal force, and proportional to the volume change, we can write this force as: $latex \lambda \frac{\partial u_k}{\partial x_k}{\delta }_{ij} $ where $latex \lambda $ is a proportionality constant dependent on the fluid. This is added to the above stress expression to yield:$latex {\sigma }_{ij}=P{\delta }_{ij}+{\tau }_{ij}+\lambda \frac{\partial u_k}{\partial x_k}{\delta }_{ij}$
(7.10)
The net force per volume of fluid due to the stress can be written as:$latex \frac{\partial {\sigma }_{ij}}{\partial x_j}$
(7.11)
Recall that there is a summation over the index j indicating that this expression is in a vector in the i direction. Inserting Eqn. (7.8) for the stress term in Eqn. (7.10), and then insert this result into Eqn. (7.11) results in the net stress force on a fluid element. This is then added to the Euler equation keeping in mind that we have already included the pressure term with the stresses. The result is:$latex \rho \frac{Du_i}{Dt}=\frac{\partial {\sigma }_{ij}}{\partial x_j}+\rho g_i=\frac{\partial P}{\partial x_j}{\delta }_{ij}+\frac{\partial }{\partial x_j}\left(\lambda \frac{\partial u_k}{\partial x_k}\right){\delta }_{ij}+\rho g_i+\frac{\partial }{\partial x_j}\ \left(\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)\right)$
(7.12)
This equation has three velocity components and the pressure as unknowns. Two fluid parameters, $latex \mu$ and $latex \lambda$ are also included. The parameter $latex \lambda$ is known as the coefficient of bulk viscosity, and is representative of the force required to change the fluid volume. For an incompressible fluid, since $latex \frac{\partial u_k}{\partial x_k}=0,$ (with summation) this term disappears. For a compressible fluid Stokes hypothesis is that:$latex \lambda =\frac{2}{3}\mu$
(7.13)
This comes with some theoretical justification and has been shown over a wide range of conditions to be reasonably valid, but shown not to be valid for extreme pressures. Since in this book we restrict ourselves to incompressible flows, this is not an issue and we will neglect this term. Also, for incompressible flows with constant viscosity we can simplify the last (viscous) term in Eqn. (7.12) as:$latex \frac{\partial }{\partial x_j}\ \left(\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\right)\right)=\mu \frac{{\partial }^2u_i}{\partial {x_j}^2}+\mu \left(\frac{\partial }{\partial x_i}\left(\mu \frac{{\partial }u_j}{\partial {x_j}}\right)\right)=\mu \frac{{\partial }^2u_i}{\partial {x_j}^2}$
(7.14)
wherein the second term in the middle expression we have reversed the order of differentiation of $latex x{}_{i}$ and $latex x{}_{j}$ and then invoked the incompressible condition leaving only the first term, which is the Laplace of the velocity vector. The final form of the NavierStokes equation that we will use for incompressible constant viscosity (Newtonian fluid) flows is:$latex \rho \frac{Du_i}{Dt}=\rho \frac{\partial u_i}{\partial t}+\rho u_j\frac{\partial u_i}{\partial x_j}=\frac{\partial P}{\partial x_j}{\delta }_{ij}+\rho g_i+\mu \frac{{\partial }^2u_i}{\partial {x_j}^2}$
(7.15)
The left hand side has been expanded to show the local and convective acceleration terms and on the right hand side we have force contributions from pressure gradients, gravity and viscous forces (normal and shearing). This equation is coupled with the requirements for incompressible flow:$latex \frac{\partial u_k}{\partial x_k}=0$
(7.16)
Notice, since we have incompressible flow it is possible to rewrite the convective acceleration term:$latex \rho u_j\frac{\partial u_i}{\partial x_j}=\rho \frac{\partial \left(u_iu_j\right)}{\partial x_j} $
(7.17)
Sometimes this may be referred to as the ``conservative form'' of the NavierStokes equation since it takes advantage of the conservation of mass for incompressible flow. We can also rewrite the NavierStokes equation using the definition of vorticity. First the convective acceleration term invokes the vector identity mentioned, in an earlier chapter in discussion of the generalized Bernoulli equation, that is written as:$latex \rho u_j\frac{\partial u_i}{\partial x_j}=\ \rho \left(\frac{\partial \left(\frac{1}{2}u_ju_j\right)}{\partial x_i}+{\varepsilon }_{ijk}{\omega }_ju_k\right)$
(7.18)
The viscous terms can also use the following vector identity:$latex \mu \frac{{\partial }^2u_i}{\partial {x_j}^2}=\mu {\varepsilon }_{ijk}\frac{\partial {\omega }_k}{\partial x_j} $
(7.19)
This last term is the viscosity times the curl of the vorticity. Notice what happens when the flow is irrotational. The second term of Eqn. (7.18) is zero and the viscous stress term, Eqn. (7.19), is also zero. It is left for the reader to see that under these conditions one can arrive at the simplified Bernoulli equation when integrated between any two points within the flow field. If the flow is inviscid the last term of Eqn. (7.15) is zero, but the flow is not necessarily irrotational since the vorticity also enters into the equation in the convective acceleration term. So, an irrotational flow does imply an inviscid flow but an inviscid flow does not imply irrotational flow. Lastly, the viscous term vanishes if the curl of the vorticity is zero, it is not necessary that the flow be irrotational.$latex u_1=C\left(1\frac{{x_2}^2}{h^2}\right)$
(7.20)
where $latex C=\frac{h^2}{2\mu }\frac{dP}{dx_1}$ for pressure driven horizontal flow and $latex C=\frac{\rho gh^2}{2\mu }$ for gravity driven vertical downward flow. Given this velocity field the vorticity distribution can be determined and since the velocity is parabolic the vorticity distribution will be piecewise linear across the flow with zero value at the center, and maximum value at the surface, $latex {\omega }_3=\frac{2C}{h^2}x_2$. The same flow can be considered for flow in a round constant area horizontal pipe. This is often called HagenPoiseuille flow after two experimentalists who developed this relationship in the 1800s for steady pipe flow in the form of flow rate versus pressure drop. In this case the NavierStokes equation is written in cylindrical coordinates as shown in Table 7.1 and then reduced for steady, constant area, incompressible, constant property flow. Again, the flow has no acceleration for steady, constant area flow and the result is a balance of pressure and viscous forces in the [latex]z[/latex], or axial direction. The pressure gradient along z is a constant similar to the result for flow between parallel plates. The resulting simplified equation is integrated across the flow (r coordinate) using the noslip boundary condition at $latex r=R $ ($latex R$ is the pipe radius) and also noting that the velocity is a maximum at the centerline where the velocity derivative in $latex r $ is zero. The result is:$latex v_z=\frac{R^2}{4\mu }\frac{dP}{dz}\left(1\frac{r^2}{R^2}\right)$
(7.21)
For all of these internal flows the flow rate can be obtained by integrating the velocity distribution across the flow area. For instance, for a volumetric flow rate of $latex \dot{Q}$ and writing the pressure gradient $latex (dP/dx_1)$ as $latex \Delta P/L$ since it is constant (notice the sign change since [latex]\Delta P[/latex] is the upstream pressure minus the downstream pressure) and $latex L$ is the total length along $latex x{}_{1}$. We have for steady horizontal pressure driven flow between parallel plates:$latex \dot{Q}=\frac{2h^3\Delta P}{3\mu L}$l (per distance $latex x{}_{3})$
(7.22)
For steady flow in a pipe:$latex \dot{Q}=\frac{\pi R^4\Delta P}{8\mu L} $
(7.23)
It is straight forward to evaluate similar relationships for gravity driven flows by replacing the pressure gradient with $latex \rho g$. There are a number of other steady flows that can also be solved for analytically from the NavierStokes equation. For instance, steady axial flow in an annulus, either pressure or gravity driven; flow in the gap between two rotating cylinders that have different rotation rates; flow between parallel plates where one plate is moving relative to the other with or without a pressure gradient or gravity driving the flow, and a few others. All of these flows have no local or convective acceleration terms and the basic governing equation is a balance of forces. It is possible to solve some flows with acceleration as shown next.$latex \frac{\partial u_1}{\partial t}=\nu \frac{{\partial }^2u_1}{\partial x^2}$
(7.24)
where $latex \nu =\mu /\rho$ which is called the kinematic viscosity and has SI units of $latex (m{}^{2}$/s). This is obtained by dividing through by the fluid density. The associated boundary and initial conditions are: Noslip: $latex u_1=U_p\ at\ x{}_{2} = 0$ No flow far from the surface: $latex u_1=0\ at\ x_2\to \infty $ Initially no flow: $latex u_1=0\ at\ t=0\ for\ all\ x_2$ This is what is called a ``diffusion'' problem where velocity (or momentum per mass) of the fluid diffuses away from a source of momentum (the wall friction force). Fick's law of diffusion states that the diffusion rate, or flux of a quantity, is proportional to the local gradient of the quantity. As we will see this diffusion process is directly linked to the fluid viscosity. The above governing equation can be scaled in several different ways. First let's use dimensional analysis with the condition that the velocity at any point in space and time within the fluid near the plate is a function of the following list of variables$latex u=(U_p, x_2,\ t,\ \nu )$
(7.25)
where we have used u instead of $latex u{}_{1}$ since there is only one velocity component to consider. This list of variables should be obvious from the governing equation, (7.24) and the associated initial and boundary conditions. Now using dimensional analysis we arrive at the following three nondimensional Pi groups:$latex {\mathit{\Pi}}_1=\frac{u}{U_p}, {\mathit{\Pi}}_2=\frac{x_2}{\sqrt{\nu t}},{\mathit{\Pi}}_3=\frac{x_2}{U_p t}$
(7.26)
However, $latex {\mathit{\Pi}}_3$ provides no new information since all of the variables appear in the other two Pi groups. Hence we can conclude that$latex \Pi_1=f \left({\Pi_2}\right) \text{ or } \frac{u}{U_p} = f\frac{x_2}{\sqrt{\nu t}}$
(7.27)
This implies that the velocity scaling for this problem is [latex]U_p[/latex] and the length scale is $latex \sqrt{\nu t}$. Notice here that the length scale is a function of time. That is to say that the extent of the length scale increases over time, and is not a constant. We can think of this length scale as roughly the distance over which diffusion is expected to occur, which is time dependent. For convenience in the mathematics we will add an arbitrary ``2'' to the length scale and express it as:$latex \frac{u}{U_p}=f(\frac{x_2}{2\sqrt{\nu t}})$
(7.28)
We have two variables, $latex f=\frac{u}{U_p}$ and $latex \eta =\frac{x_2}{2\sqrt{\nu t}}$. And we seek a solution of $latex f(\eta )$. The next step is to substitute our nondimensional variables into the governing equation and initial and boundary conditions. This needs to be done carefully since one of our scales (length) is a function of time. The reader is encouraged to go through the following transformation using the chain rule: $latex \frac{\partial u}{\partial t}=\frac{\partial (fU_p)}{\partial \eta }\cdot \frac{\partial \eta }{\partial t}=U_p\frac{\partial f}{\partial \eta }(\frac{\eta }{2t}) $ $latex \frac{\partial u}{\partial x_2}=U_p\frac{\partial f}{\partial \eta }\cdot \frac{\partial \eta }{\partial x_2}=U_p\frac{\partial f}{\partial \eta }\frac{1}{2\sqrt{\nu t}} $ $latex \frac{{\partial }^2u}{\partial {x_2}^2}=\ \frac{\partial }{\partial x_2}\left(U_p\frac{\partial f}{\partial \eta }\cdot\frac{\partial \eta }{\partial x_2}\right)=U_p\frac{\partial f}{\partial \eta }\frac{{\partial }^2\eta }{\partial {x_2}^2}+U_p\frac{\partial \eta }{\partial x_2}\frac{{\partial }^2f}{\partial {\eta }^2}\frac{\partial \eta }{\partial x_2} $ $latex \frac{{\partial }^2u}{\partial {x_2}^2}=U_p{\left(\frac{1}{2\sqrt{\nu t}}\right)}^2\frac{{\partial }^2f}{\partial {\eta }^2} $ Inserting the first and the last of this set of equations in Eqn. (7.24) the result is:$latex \frac{{\partial }^2f}{\partial {\eta }^2}+2\eta \frac{\partial f}{\partial \eta }=0=f^{''}+2\eta f'$
(7.29)
where the last expression is used to write derivatives as primes $latex (f^{''}$ is the second derivative and $latex f'$ is the first derivative.) Notice here that the derivatives can actually be ordinary derivatives since as shown in (7.27) $latex f$ is only a function of $latex \eta$. Importantly here the second order partial differential equation in two variables is transformed into a second order ordinary differential equation, which in general has a much easier solution.The next step is to transform the boundary conditions as follows: $latex {\text {at}\ \ x}_2=0\ \ \eta =0\ \text{and}\ u=U_p\ \text{so}\ f\left(0\right)=1$ $latex \text {as}\ x_2\to \infty \ \ \ \eta \to \infty \ \text{and}\ u\to 0\ \text{so}\ \ f\left(0\right)=0 $ $latex \text{at}\ t=0\ \ \ \eta \to \infty \ \ \text{and}\ u\to 0\ \text{so}\ \ f\left(0\right)=0 $
(7.30)
The last two conditions are identical, so we end up with two boundary conditions for [latex]\eta[/latex] which is what we expect since we have a second order equation. We are left with the governing equation (7.29) and the boundary conditions (7.30). It turns out that there is an analytical solution. Letting: $latex g=f' $ then the governing equation is: $latex \frac{g'}{g}=2\eta $ this can be integrated to: $latex g=f'=C_1e^{{\eta }^2} $ which we integrate again to:$latex f=C_1\int^{\eta }_0{e^{{\eta }^2}d\eta }+C_2 $
(7.31)
The first term on the right hand side is a straightforward integral evaluation up to a variable value $latex \eta$. This can be replaced by the Gaussian Error Function, $latex erf(\eta)$ which is shown in Fig (7.6) as: $latex {erf \left(\eta \right)\ }=\frac{2}{\sqrt{\pi }}\int^{\eta }_0{e^{{\eta }^2}d\eta } $ so that: $latex f=\frac{\sqrt{\pi }}{2}C_1{erf \left(\eta \right)\ }+C_2 $ [caption id="attachment_273" align="aligncenter" width="300"] Fig. 7.6 Gaussian Error Function, erf(x), where x is the argument; the function is normalized to vary between 1 and +1 and passes through zero at $latex x = 0$.[/caption] We are left to evaluate the constants from the two boundary conditions for f.$latex \text{Since}\ erf (0) = 0 \text{ and } erf (\infty )=1$ then $latex C_1=\frac{2}{\sqrt{\pi }} \text{ and } C_2 = 1$ and the final form of the equation is:
$latex f=\frac{u}{U }=1erf(\eta) \text{ and } f' ={\frac{2}{\sqrt{\pi}}e^{^\eta}}^2$
(7.32)
The results of Eqn. (7.32) are plotted in Fig. 7.7. This solution is a ``similarity solution'' in that the results collapse, in this case, into a single variable, $latex \eta$, which depends on the two independent variables, $latex x{}_{2}$ and t. There are two parameters affecting the result of the local velocity, the boundary condition, [latex]U_p[/latex], and the fluid property, $latex \nu$. To evaluate a given set of conditions one does the following. (1) Identify the physical space and boundary conditions of interest, (2) transform variables such as values of space or time into associated values of the variable $latex \eta$. (3) Use Eqn. (7.32) to determine the value of [latex]f[/latex] versus [latex]\eta[/latex] and the desired $latex \eta$ value. (4) Transform f into the desired velocity $latex u$. [caption id="attachment_274" align="aligncenter" width="300"] Fig. 7.7 Nondimensional result of $latex f'$ versus $latex \eta$.[/caption] Besides determining the local velocity in the fluid these results are also useful in determining the transient wall shear stress, or force per unit area of the plate. Given in Eqn. (7.32) is the derivative, $latex f'$. This provides information of the velocity derivative, $latex \frac{\partial u}{\partial x_2}$. Since there is only one velocity component and it is only a function of $latex x{}_{2}$ we can write the shear stress as:$latex \tau =\mu \left(\frac{\partial u_2}{\partial x_1}+\frac{\partial u_1}{\partial x_2}\right)=\mu \frac{\partial u_1}{\partial x_2}=\mu \left(\frac{Uf'}{2\sqrt{\nu t}}\right) $
(7.33)
where $latex f'$ is the derivative $latex \frac{\partial f}{\partial \eta }$. The force on the surface of the plate with an area of LS, (where L is the distance along the plate in the direction of motion and S is the width, or span of the plate) due to its motion is:$latex F_p=\int^L_0{{\tau }_sdx_1S}$
(7.36)
where the shear stress, $latex {\tau }_s$, is evaluated at the surface where $latex x{}_{2}$ = 0. Evaluating (7.34) using (7.33) at $latex \eta =0$, with $latex f'(0)=\frac{2}{\sqrt{\pi } } $ results in:$latex {\tau }_s=\frac{\mu U_p}{\sqrt{\pi \nu t}}$ and the force $latex F_p=\frac{\mu U_p LS}{\sqrt{\pi \nu t}}$
(7.37)
Notice that the wall stress decreases with time as momentum diffuses away from the surface. The force therefore also decreases with time. At time t = 0 the equations have a problem, since we said the acceleration was initially infinite. We can define the diffusion distance way from the surface as the distance $latex x{}_{2}$ to where the velocity decays to 1% (small arbitrary value) of the plate velocity. So setting $latex f$ = 0.01 we see that $latex \eta$ = 1.85 and setting [latex]x_2=\delta[/latex] at this point we have:$latex \delta =3.7\sqrt{\nu t}$
(7.38)
Notice that the diffusion distance increases with time as well as kinematic viscosity. The greater the viscosity the further is the distance, $latex \delta$, at any given time. It is independent of the plate velocity. One last item on this seemingly simple problem that will help in the interpretation of boundary layer flows. This flow has vorticity, but only aligned with the $latex x{}_{3}$ direction, as can be seen by calculating the local vorticity distribution:$latex {\omega }_3=\left(\frac{\partial u}{\partial x_2}\right)=\frac{Uf'}{2\sqrt{\nu t}}$
(7.35)
The circulation is the area integral of the vorticity, so extending an integration from the surface at [latex]x_2=0[/latex] to [latex]\delta[/latex], for a distance L along $latex x{}_{1}$, then the total circulation per distance $latex x{}_{3}$ becomes:$latex \Gamma=\int ^\delta_0  \frac{U_p f'}{2 \sqrt{\nu t}} dx_2L= \int_{0}^{\infty }U_p f'd \eta L=U_p L $
(7.36)
Interestingly, the total circulation remains constant overtime, but the area grows since [latex]\delta[/latex] increases with time. The diffusion process merely spreads the vorticity over an increasing area, for increasing time and all of the vorticity is generated with the initial acceleration of the surface.]]>$latex {\tau }_s=\mu \ {\left(\frac{\partial u_2}{\partial x_2}\right)}_{x_2=0}$
(8.1)
where the subscript [latex]x_2=0[/latex] means that it is evaluated at the surface. Since the surface velocity derivative decreases along $latex {x}_{1}$ then so does the surface shear stress and the qualitative result is the trend of [latex]\tau_s[/latex] shown in Fig. 8.2. Notice that the surface shear stress does not reach a constant value, as in pipe flow, since the flow never becomes fully developed, or [latex]\frac{\partial{u_1}}{\partial{x_1}}\neq 0[/latex]. If the surface is perfectly flat then the freestream flow outside the boundary layer will be nearly parallel with the surface (in the $latex {x}_{1}$ direction) at velocity U, except for the small upward velocity component described earlier. For parallel straight flow the change of pressure along the flow, and normal to the flow, is zero as we saw in inviscid flow. If we anticipate that the thickness, $latex \delta$, of the boundary layer is relatively thin (we will show this to be the case later in this chapter when we develop the boundary layer equations) then the change of pressure due to hydrostatic effects between the outer edge of the boundary layer and the surface is small, even if the fluid is a liquid. However, even if hydrostatic pressure rise is not small across the boundary layer the change in pressure along the flow should in fact approach zero near the surface because the change in pressure along the flow is zero just outside of the boundary layer. This is explained by the following. The hydrostatic pressure change is the same between two vertical points evaluated at different locations along $latex {x}_{1}$. So the difference of pressure in the $latex {x}_{1}$ direction, $latex \frac{\partial P}{\partial x_1}$ will be identical for positions near the edge of the boundary layer and near the surface. So if $latex \frac{\partial P}{\partial x_1\ }=0$ just outside the boundary layer then it is also zero inside the boundary layer. The net effect of this is that for a flat surface the pressure gradient along the flow is essentially zero. We will show this to be true based on the governing equations later in this chapter. We could extend this further for a curved surface and deduce that the pressure gradient for flow over the curved surface is the same outside and inside the boundary layer if the boundary layer is thin. So a summary of the main attributes we will rely on for the analysis of boundary layers at this point are: (i) the boundary layer is thin, $latex \delta$ >> $latex L$, (ii), the wall normal velocity is small but finite, $latex {u}_{2}$ << $latex {u}_{1}$, (iii) for a (nearly) flat surface the pressure gradient is zero, $latex {\partial {P}}/{\partial {x}_{1}} = 0$, (iv) the outer edge of the boundary layer velocity is the freestream velocity, U, which is determined by inviscid flow conditions.$latex \rho \left(u_1\frac{\partial u_1}{\partial x_1}+u_2\frac{\partial u_1}{\partial x_2}\right)=\mu \left(\frac{{\partial }^2u_1}{\partial {x_1}^2}+\frac{{\partial }^2u_1}{\partial {x_2}^2}\right)$
(8.2)
we combine this with the differential conservation of mass for an incompressible flow:$latex \left(\frac{\partial u_1}{\partial x_1}+\frac{\partial u_2}{\partial x_2}\right)=0$
(8.3)
At this point we have twodimensional partial differential equations in two variables. It is worth noting here that the pressure gradient in the y direction, across the flow, is reduced to a very small number under the conditions that the $latex u{}_{2}$ component of velocity is very small compared with $latex u{}_{1}$, and that the boundary layer is thin so that there is essentially zero hydrostatic pressure variation across the flow. This can be shown by noting that all terms in the $latex {x}_{2}$ momentum equation containing $latex u{}_{2}$ become very small and the body force term is small since changes in $latex {x}_{2}$ are small, therefor the pressure gradient in $latex {x}_{2}$ must also be small. This process will become clearer in our discussions of scaling of the $latex {x}_{1}$ momentum equation below. Next we scale the variables in the governing equation. This means coming up with appropriate scaling for each variable (dependent and independent) that is on the order of the maximum value of the variable. For instance, the characteristic length scale for flow over a flat plate with length L would be just L. Similarly, the characteristic velocity scale for the $latex {x}_{1}$ velocity component would be the freestream velocity, U. Below is listed the characteristic scales for the various variables (note that the scales may not be constant values as is noted by the choice of scale for $latex {x}_{2}$):$latex x_1L$
(8.4)
$latex x_2\ \delta $ $latex u_1U_{\infty }$ $latex u_2unknown$ In general L can be any length of surface. We select $latex \delta $ to be the boundary layer thickness that coincides with this selected length in the analysis below. That is to say that the largest value of [latex]\delta[/latex] occurs at [latex]x_1=L[/latex]. But in general [latex]\delta[/latex] varies as we move in the $latex {x}_{1}$ direction. We do not know how large $latex u{}_{2}$ may be, but we can determine this by examining the conservation of mass, Eqn. (8.3). To do this we use the above scaling to form nondimensional variables, denoted by the superscript ``*'' by dividing each variable by its scale value:$latex {x_1}^*=\frac{x_1}{L}$
(8.5)
$latex {x_2}^*=\frac{x_2}{\delta }$ $latex {u_1}^*=\frac{u_1}{U}$ $latex {u_2}^*=\frac{u_2}{V_s}$ where $Latex V{}_{s}$ is an unknown scale for variable $latex u{}_{2}$ and is to be determined. Inserting these nondimensional variables into Equation (8.3), where $latex u=u^*U$, etc., we obtain:$latex \frac{U{}}{L}\left(\frac{\partial {u_1}^*}{\partial {x_1}^*}\right)+\frac{V_s}{\delta }\left(\frac{\partial {u_2}^*}{\partial {x_2}^*}\right)=0=\frac{\partial {u_1}^*}{\partial {x_1}^*}+\frac{V_sL}{\delta U}\left(\frac{\partial {u_2}^*}{\partial {x_2}^*}\right)$
(8.6)
where the last set of terms is a consequence of dividing the first set of terms through by $latex \frac{U}{L} $. The coefficient of the first term become one, the coefficient multiplying the second term becomes $latex \frac{V_sL}{\delta U}$. At this point we apply an order of magnitude analysis of the terms. Both terms have to balance each other in the equation, being of opposite sign. We use the nomenclature $latex \mathcal{O}(1)$ to indicate that the value of a term is ``order one''. By proper scaling of the variables we assume that, $latex \frac{\partial {u_1}^*}{\partial {x_1}^*}$, is $latex \mathcal{O}(1)$. That is to say, by proper scaling this derivative using scales that are approximately the maximum value of the variable then the resultant term will not be a very large or a very small number, but rather ``order one''. We make the same assumption for the derivative $latex \frac{\partial {u_2}^*}{\partial {x_2}^*}$. Consequently since neither of these derivatives are zero, as we have noted previously for the boundary layer, then for this equation to be valid we must have $latex \frac{V_sL}{\delta U}=\ \mathcal{O}(1)$ With this we can assign the value of the $latex u{}_{2}$ velocity scale to be:$latex V_s=U\frac{\delta }{L}$
(8.7)
By doing so then each term is $latex \mathcal{O}(1)$. Keep in mind that we are not trying to find exact values but at this point trying to determine approximately how large or small the scaling values should be. If we now take Equation (8.7) and insert this into Equation (8.6) we end up with:$latex \left(\frac{\partial u^*}{\partial x^*}\right)+\left(\frac{\partial v^*}{\partial y^*}\right)=0$
(8.8)
for the nondimensional representation of the conservation of mass with the scaling of the variables given by the relationships above. Using our new found scaling for $latex u{}_{2}$ we continue this same nondimensional process for the other terms in the NavierStokes Equation (8.2). That is, replacing all of the variables with the nondimensional variables, and then dividing the resulting equation by the coefficient of the first term, $latex \frac{{\rho U}^2}{L}$ resulting in:$latex {u_1}^*\frac{\partial {u_1}^*}{\partial {x_1}^*}+{u_2}^*\frac{\partial {u_1}^*}{\partial {x_2}^*}=\frac{\mu }{\rho UL}\left(\frac{{\partial }^2{u_1}^*}{\partial {{x_1}^*}^2}+\frac{L^2}{{\delta }^2}\frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2}\right)$
(8.9)
This nondimensional equation, that incorporates the scaling given above, leads to some profound conclusions. First note that the coefficient on the right hand side of the two terms in parentheses, which are the fluid friction terms, is the inverse of the Reynolds number that we define for boundary layer flows: $latex Re_L=\frac{{\rho U}_{\infty }L}{\mu }$. This Reynolds number is based on the freestream velocity and the total length of the surface. We have included a subscript L on the Reynolds number to inform us of the selected length scale. Next we have the coefficient $latex \frac{L^2}{{\delta }^2}$ of the last term in parentheses. Since we have assumed that the boundary layer is thin ($latex \delta \ll L$) then this coefficient will be large. Since each of the derivatives is assumed to be order one through proper scaling, then the entire second term in parentheses must be much larger than the first term, or:$latex \frac{{\partial }^2{u_1}^*}{\partial {{x_1}^*}^2}$ << $latex \frac{L^2}{{\delta }^2}\frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2}$
(8.10)
So the frictional forces in the boundary layer are dominated by the term containing $latex \frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2}$. This term is the net angular deformation caused by shearing stress. The first term is the linear deformation caused by the normal viscous stress. This result is consistent with the fact that the thin boundary layer will have very large derivatives of velocity in the $latex {x}_{2}$ direction when compared with the $latex {x}_{1}$ direction. So neglecting the normal viscous stress contribution the resulting nondimensional form of the boundary layer equation for steady flow over a flat surface is:$latex {u_1}^*\frac{\partial {u_1}^*}{\partial {x_1}^*}+{u_2}^*\frac{\partial {u_1}^*}{\partial {x_2}^*}=\frac{\mu }{\rho UL}\left(\frac{L^2}{{\delta }^2}\frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2}\right)$
(8.11)
There is one more conclusion that can be arrived at from this scaling. Noting that the shearing contribution derivative, $latex \frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2},$ is $latex \mathcal{O}(1)$ as are each of the nondimensional acceleration terms on the left hand side, then in order for this equation to be balance $latex \left(\frac{\mu }{\rho UL}\frac{L^2}{{\delta }^2}\right)$ should also be $latex \mathcal{O}\left(1\right)$:$latex \frac{\mu }{\rho UL}\left(\frac{L^2}{{\delta }^2}\right)\sim \ \mathcal{O}(1)$
or$latex Re{_L}\sim \frac{L^2}{{\delta }^2}$
(8.12)
This last condition implies that for a thin boundary layer to exist the Reynolds number should be large. So for large Reynolds numbers in boundary layer flows with a surface length of L the boundary layer thickness, $latex \delta$, is proportional to the square root of $latex Re{}_{L}$ which we can write as:$latex \frac{\delta }{L}={\frac{C}{\sqrt{Re}_{x_1}}}$
(8.13)
where C is some unknown constant inserted to make the equation quantitatively correct. Consequently, for very small values of [latex]L[/latex] the boundary layer conditions do not hold since $latex Re{}_{L}$ is small. The nature of the governing equation given by Eqn. (8.11) is that it is a parabolic partial differential equation. Since it only requires one boundary condition in $latex {x}_{1}$, which can be selected at $latex {x}_{1} = 0$ (the leading edge) where the velocity is given as [latex]U[/latex], then the solution is independent of any downstream boundary condition. This means it can be solved at any arbitrary position $latex {x}_{1}{}_{\ }$ along the surface. This also implies that the result (8.13) is valid by replacing L with location $latex {x}_{1}{}_{.\ \ }$. Doing this we have a Reynolds number based on any position $latex {x}_{1}$: $latex Re_{x_1}=\frac{\rho Ux_1}{\mu }$ and Equation (8.13) can be expressed as:$latex \frac{\delta }{x_1}={\frac{C}{\sqrt{Re}_{x_1}}}$
(8.14)
Equation (8.14) predicts that the boundary layer grows proportional to the distance $latex x{}_{1}{}^{1/2}$ (recall that there is a factor $latex x{}_{1}$ contained in the $latex {Re}_{x_1}$); also this will not be valid for very small values of [latex]x_1[/latex]. Using Eqn. (8.12) the final form of the boundary layer equation for flow over a flat surface is:$latex {u_1}^*\frac{\partial {u_1}^*}{\partial {x_1}^*}+{u_2}^*\frac{\partial {u_1}^*}{\partial {x_2}^*}=\frac{{\partial }^2{u_1}^*}{\partial {{x_2}^*}^2}$
(8.15)
From our discussion of the characteristics of the boundary layer we have the following boundary conditions: (i) $latex {x}_{2}= 0: u{}_{1}=0, u{}_{2} =0$ (noslip), (ii) $latex {x}_{2} = \delta: u{}_{1}$ = U (iii) $latex {x}_{1} = 0: u{}_{1}$ = U This equation seems to have eliminated all parameters (like viscosity, density and Re). However, recall that $latex {x}_{2}$ is nondimensionalized using $latex \delta$ and also [latex]u_2[/latex] is nondimensionalized using [latex]\delta[/latex]. In order to determine $latex \delta$ one needs to know the relationship between [latex]\delta[/latex] and Re. So even though we may be able to solve Eqn. (8.15) there is still work to do to interpret the results.$latex \rho \left(u_1\frac{\partial u_1}{\partial x_1}+u_2\frac{\partial u_1}{\partial x_2}\right)=\mu \left(\frac{{\partial }^2u_1}{\partial x_2}\right)$
(8.16)
with the conservation of mass as Equation (8.3). Following the general solution approach of Blasius we can introduce the streamfunction, and use it to replace the velocity components:$latex u=\frac{\partial \psi }{\partial x_2};\ \ \ \ v=\frac{\partial \psi}{\partial x_1}$
(8.17)
We have previously shown that this definition for two dimensional, incompressible flow satisfies conservation of mass, by inserting these expressions into Equation (8.3). Substituting these two expressions in terms of $latex \psi$ for the two velocity components in Equation (8.16) results in:$latex \frac{\partial \psi }{\partial x_2}\left(\frac{{\partial }^2\psi }{\partial x_1\partial x_2}\right)\frac{\partial \psi }{\partial x_1}\left(\frac{{\partial }^2\psi }{\partial x^2_2}\right)=\nu \left ( \frac{\partial^3\psi}{\partial x^3_2} \right ) $
(8.18)
Here we have also divided through by density to obtain the kinematic viscosity, $latex \nu = \mu / \rho$, as the coefficient of the viscous term which is the only parameter in the equation (a constant for a given fluid). Although the resulting equation seems formidable the use of the streamfunction has reduced the two velocities into a single dependent variable $latex \psi $. The streamfunction, with SI units of $latex m{}^{2}$/s, is nondimensionalized as:[latex]f\left(\eta \right)=\frac{\psi }{{\nu x_1\frac{1}{2}}}[/latex]
(8.19)
a nondimensional similarity variable,$latex \eta$, is defined to combine $latex x_{1}$ and $latex x$ into a varible:$latex \eta =\frac{x_2}{{\left(\frac{\nu x_1}{U}\right)}^{1/2}}=x_2{\left(\frac{U}{\nu x_1}\right)}^{1/2}$
(8.20)
The new variable f is assumed to be only a function of $latex \eta$. This will be verified by inserting these new variables into the original equation and see if the result is only [latex]\eta[/latex] dependent. Each of the terms in Eqn. (8.18) are transformed into the ($latex f,\eta$) variables using the chain rule similar to what was done in Chapter 7 to determine the governing equation for a suddenly accelerated flat surface. Similarly, the nondimensional stream function of Eqn. (8.18) can be written in terms of [latex]u_1[/latex] and [latex]u_2[/latex] and inserted into Eqn. (8.16). The result is the following, using the prime designation to represent the derivative with respect to $latex \eta$$latex u_1=\frac{\partial\psi}{\partial x_2}\frac{\partial \psi }{\partial \eta }\frac{\partial \eta }{\partial x_2}=\frac{\partial \psi }{\partial \eta }\left(\frac{U}{\nu x_1}\right)=f{'\left(\frac{U}{\nu x_1}\right)}$
(8.21)
$latex u_2=\frac{\partial \psi }{\partial x_1}=\left ( \frac{vU_\infty}{4x} \right )^{1/2}\left ( \eta f' f \right )$
(8.22)
$latex \frac{\partial u_1}{\partial x_1} =U f''\frac{\partial \eta}{\partial x_2} = \frac{uf'' \eta}{2x_1}$
(8.23)
$latex \frac{\partial u_1}{\partial x_2} Uf''\frac{\partial \eta }{\partial x_2}=Uf''{\left(\frac{U}{\nu x_1}\right)}^{1/2}$
(8.24)
$latex \frac{{\partial }^2u_1} {\partial {x_2}^2} = \frac{U^2f^{'''}}{\nu x_1} $
(8.25)
Substituting all of this into our boundary layer equation, (8.16), with a little bit of manipulation the result becomes:$latex f^{'''}+\frac{1}{2}ff^{''}=0$
(8.26)
This now is a single variable ordinary differential equation, which shows that $latex f = f(\eta)$ resulting in a similarity equation. To obtain a solution the boundary conditions need to be specified in terms of $latex \eta$. Specifically there are three: $latex x_2=0:u_1=u_2=0$ (no slip condition) so$latex f\left(0\right)=f'\left(0\right)=0$
$latex x_2\ge \delta :u_1=U$ so $latex f'\left(\eta \to \infty \right)=1.0$(8.27)
The no slip condition for $latex u_{2}$ results in $latex f\left(0\right)=0$ since $latex f'\left(0\right)=0$ for the no slip condition of $latex u_{1}$ so by Eqn. (8.22) $latex f\left(0\right)=0$. Also, we take the boundary condition at the edge of the boundary as $latex \delta \rightarrow \infty \text{ since } U$ remains constant for $latex x_2$ > $latex \delta $ . Physically this imposes an asymptotic limit of the solution for large $latex \delta$, or beyond the edge of the boudary layer. The solution to Equation (8.26) along with the boundary conditions of Equations (8.27) is rather straight forward using a RungeKutta technique. However, this was not the case when Blasius solved this equation numerical by hand back in the early 1900s at the University of Göttingen. The results of the numerical solution of Eqn (8.26) are expressed as $latex f'\left(\eta \right)=u_1/U_{\infty }$ and are shown in Figure (8.3). Obviously $latex u_{1}/U$ varies from zero at the surface ($latex \eta = 0$) and approaches 1.0 at $latex x_{2} = \delta$. The value of $latex \delta$ can be determined from this plot, however one must select a value of $latex f'\left(\delta \right)$ close to 1.0, and typically 0.99 is chosen. At this value of $latex f'\left(\delta \right) $ the value of $latex \eta \simeq 5.0$. Using the definition of $latex \eta$ given above and setting it equal to 5.0 when $latex x_{2} = \delta$ results in:$latex 5=\frac{\delta }{{\left(\frac{\nu x_1}{U}\right)}^{1/2}}$
Rearranging this yields an expression for the boundary layer versus $latex x_{1}$:$latex \frac{\delta }{x_1}=\frac{5}{{Re_{x_1}}^{1/2}}\$
(8.28)
This shows that $latex C = 5$ in Eqn. (8.13). Eqn. (8.28) allows for the evaluation of the boundary layer thickness, $latex \delta$, at any $latex x_{1}$ location for a given value of U and viscosity $latex \nu $. Using the solution in Fig. (8.3) also allows for the evaluation of the velocity components, $latex {u}_{1}$ and [latex]{u}_{2}[/latex] and any location $latex x_{1}$ and [latex]x_{2}[/latex] in the boundary layer. [caption id="attachment_279" align="aligncenter" width="300"] Fig 8.3 Blasius solution to steady laminar flow over a flat plate.[/caption] Other aspects of the flow can also be found from the solution. The total force on the surface caused by fluid flow friction can be found from the surface stress, $latex {\tau }_s=\mu \frac{\partial u_1}{\partial y}$, which requires evaluating the velocity derivative at the surface, $latex x_{2}$ = 0. The nondimensional form is used to evaluate $latex \frac{\partial u_1}{\partial y} \text{ at } x_2 =0$. Using Eqn. (8.24) we have:$latex {\left(\frac{\partial u_1}{\partial x_2}\right)}_{x_2=0}=f''(0)U^{3/2}\sqrt{\frac{\rho }{\mu x_1}} $
(8.29)
where $latex f''(0)$ is the second derivative of f evaluated at $latex \eta =0$. From the Blasius solution $latex f^{''}\left(0\right)=0.332$. After multiplying Equation (8.29) by the dynamic viscosity to obtain the surface stress and then dividing by $latex \frac{1}{2}\rho {U}^2$ which defines the surface skin friction coefficient, $latex c_f,$ the result is:$latex c_f=\frac{{\tau }_s}{\frac{1}{2}\rho U^2}=\ \frac{0.664}{{Re_{x_1}}^{1/2}}$
(8.30)
The above expression can be used to find the local surface stress at any [latex]x_1[/latex] location along the surface (recall $latex Re_{x_1}=\frac{\rho U_{\infty }x_1}{\mu }$). This then indicates that for a given freestream velocity of a given fluid the surface shear stress decreases as $latex {x_1}^{1/2}$. The final step to determine the total drag force, $latex F{}_{D}$, on the surface is to integrate [latex]\tau_s[/latex] over the entire surface using the results of Equation (8.30). [latex]F_d = \int_{0}^{L}\tau_s dx_1S[/latex] which is written as the drag coefficient, $latex {C}_{D}$:$latex C_D=\frac{F_D}{\frac{1}{2}\rho {U^2{(LS)}}^2}=\frac{1.328}{{Re_L}^{1/2}} $
(8.31)
$latex \frac{\partial P'}{\partial x_1}=\mu {\left(\frac{{\partial }^2u_1}{\partial {x_2}^2}\right)}_{x_2=0}={\left(\frac{\partial \tau }{\partial x_2}\right)}_{x_2=0} $
(8.32)
For the boundary layer, assuming it is thin, the pressure gradient can be evaluated outside the boundary layer as discussed previously. Also, in this expression we have absorbed the body force term into the pressure term: $latex P'=P+\rho gh$ where h is the height above some datum. This form of the pressure plus body force accounts for the hydrostatic pressure variation that may occur along the direction $latex x_1$ caused by changes in elevation given by changes in h. Eqn. (8.32) provides insight into the flow conditions very close to the surface for changing pressure gradients. For instance, if $latex \frac{\partial P'}{\partial x_1} < 0$ (a ``favorable pressure gradient'') then the pressure decreases in the flow (or $latex x_{1}$) direction. This implies from Eqn (8.32) that the surface stress decreases moving away from the surface as one might expect since it ultimately reaches zero at the edge of the boundary layer, for inviscid flow conditions. For conditions of very large favorable pressure gradients ($latex \frac{\partial P'}{\partial x_1} \ll 0$) the decrease of stress with $latex x_{2}$ is greater and the boundary layer becomes thinner, and the wall shear stress becomes larger. The distribution of stress is shown in qualitatively Fig. 8.4. [caption id="attachment_280" align="aligncenter" width="539"] Fig. 8.4 Stress distribution for favorable and adverse pressure gradients in boundary layer flows; flow over a cylinder illustrating favorable and adverse pressure gradients for flow over a cylinder.[/caption] If the pressure gradient has a positive sign, $latex \frac{\partial P'}{\partial x_{1}}$ > $latex 0$, (an "adverse pressure gradient'') the stress increases with increasing $latex x_{2}$ away from the surface again as shown in Eqn. (8.32). This implies that at some point the stress reaches a maximum, and then decays to zero at the edge of the boundary layer. See Fig. 8.4. As the adverse pressure gradient increases one can suspect that the gradient of stress near the surface also becomes larger, forcing the stress to approach zero or even go negative as shown in Fig. 8.4. Under these conditions, of negative stress, the implications are that the sign reversal causes the flow in a region near the wall to actually go in the negative, or upstream, $latex x_{1}$ direction. This is the definition of ``flow separation''. Separated flow no longer adheres to boundary layer flow as part of the flow near the surface is moving upstream and part of the flow further from the surface is moving downstream. This generally results in swirling flow near the surface. What can cause such a large adverse pressure gradient? One common situation is the existence of large surface curvature, such as flow over a cylinder of diameter D and length S into the page, shown in Fig. 8.5. In these cases the $latex {x}_{1}$ coordinate is taken as curving along the surface in the flow direction. On the upstream side of the cylinder the pressure gradient is negative as the flow accelerates with increasing velocity. The growth of $latex \delta$ is not as large as for the zero pressure gradient case, as noted above. As the flow continues over the top (or bottom) of the cylinder the pressure gradient becomes positive, increasing pressure in the flow direction. This is caused by the decrease of velocity for the flow beyond the top (bottom) of the cylinder. The inviscid flow over a cylinder has been studied previously (as a superposition of uniform flow and doublet). For this inviscid case the pressure can be found from the Bernoulli equation and the result can be expressed as the pressure coefficient as:$latex c_p=\frac{PP_{\infty }}{\frac{1}{2}\rho U^2}=1\frac{{v_s}^2}{U^2}=14sin^2\theta$
(8.33)
where $latex P_{\infty }$ is the pressure for upstream (or downstream) of the cylinder and $latex v_s$ is the surface velocity at any point around the cylinder (zero at the front and back stagnation points and maximum at the top and bottom). A plot of this inviscid pressure distribution around the cylinder is shown in Fig. 8.5. Notice that the pressure distribution is symmetric and the net force caused by pressure is zero. The maximum pressure is at the stagnation points, as expected, and equal to [latex]P_\infty + \frac{1}{2}pU^2[/latex] from Bernoulli equation. With the addition of viscous effects, for flow over a cylinder, as described above, the adverse pressure gradient on the rear portion of the cylinder results in flow separation when the freestream flow is sufficiently large. The flow then forms a ``wake'' across the back of the cylinder consisting of swirling, turbulent flow. The result is an asymmetric flow pattern around the cylinder and an asymmetric pressure distribution. In the wake the pressure is not as high as in the inviscid case. We see in the pressure coefficient plot that for viscous flows the pressure recovery on the downstream side is not as large as the inviscid case and the net affect is a downstream drag force. Turbulent flows tend to have a greater pressure recovery and a lower nondimensional drag coefficient ($latex C_D=\frac{F_D}{\frac{1}{2}\rho U^2DS}$) when compared with laminar flow. [caption id="attachment_281" align="aligncenter" width="876"] Fig. 8.5 Pressure distribution for flow over a cylinder and illustration of wake formation for laminar and turbulent flow.[/caption] ]]>$latex \dot{m}=\int^{\delta }_0{\rho u_1dx_2\ S}$
(9.1)
where $latex u{}_{1}$ is the $latex x{}_{1}$ component of the velocity and S is the span, or the distance into the page. The distribution of $latex u{}_{1}$ follows the two boundary conditions, $latex u{}_{1} = 0$ at the surface and $latex u{}_{1} = U$ at the edge of the boundary layer, $latex {\delta}$. [caption id="attachment_284" align="aligncenter" width="874"] Fig. 9.1 Illustration of the definition of the displacement thickness based on the total mass flow rate within the boundary layer; note that the mass flow rate changes along the flow in the $latex x{}_{1}$ direction consequently $latex {\delta}_{1}$ is a function of $latex x{}_{1}$.[/caption] The surface generated friction by the surface on the flow causes the velocity to decrease at any given $latex x{}_{2}$ position along the surface. This decrease across the entire boundary layer results in a decreasing flow rate within the boundary layer with increasing $latex x{}_{1}$. Therefore the amount of mass flow lost in the boundary layer due to surface friction can be thought of as local a mass flow deficit. This deficit can be calculated by comparing the actual velocity distribution to the ideal velocity without friction, which would be $latex {U}$ across the entire boundary layer. The local mass flow deficit at any $latex x{}_{2}$ location is then $latex d{\dot{m}}_{deficit}=\rho \left(Uu\right)dy\ S$ and the total mass flow deficit is obtained by integrating this across the flow as:$latex \int^{\delta }_0{d{\dot{m}}_{deficit}}={\dot{m}}_{deficit}=\int^{\delta }_0{\rho \left(Uu_1\right)dy\ S}$
(9.2)
The mass flow deficit can be equated to mass flow rate loss given as the product of density, $latex {U}$ and ($latex {\delta}_{1}$S), where $latex {\delta}_{1}$ is an unknown ``thickness''. So setting the integral term in Eqn. (9.2) to this product [latex]\left ( \rho U\delta_1 S \right )[/latex] and solving for $latex {\delta}_{1} $ results in:$latex {\delta }_1=\int^{\delta }_{\ 0}{\left(1\frac{u_1}{U}\right)dy}$
(9.3)
which assumes that the density is constant. Eqn. (9.3) is the definition of the displacement thickness and can be interpreted as follows. If the total mass flow rate lost in the boundary layer due to friction is divided by $latex \rho US$ the result is a length that would represent a distance normal to the surface through which all of the lost mass flow rate would pass if it traveled at velocity $latex U $(the velocity it would have if there were no friction.) Another way of looking at this distance is that it represents the vertical distance the surface would need to be moved upward (in the $latex x{}_{2}$ direction) to capture all of the lost mass flow rate if the velocity were uniformly at $latex {U}$, or frictionless. Either of these interpretations indicates that as $latex {\delta}_{1}$ grows larger along the flow there is an increase in the total mass fow deficit. Since surface friction slows more and more fluid, creating a larger boundary layer thickness in the x direction, then $latex {\delta}_{1}$ must increase in the downstream $latex x{}_{1}$ direction. At any given position along the surface the greater the local surface shear stress the greater the value of $latex {\delta}_{1}$. Consequently, a favorable pressure gradient which increases the local shear stress compared with zero pressure gradient, results in a larger relative value of displacement thickness. Momentum thickness, $latex {\delta}_{2}$ The momentum thickness, $latex {\delta}_{2}$, has somewhat of the same interpretation as $latex {\delta}_{1}$, but it has to do with the momentum rate loss by the frictional forces, rather than the mass flow rate loss. In this case we develop a measure for the momentum rate loss, or momentum deficit, generated by friction. Noting that the momentum rate lost crossing a given vertical plan across the boundary layer is given by:$latex \int^{\delta }_0{\rho u\left(Uu_1\right)dx_2\ S}$
(9.4)
where $latex \rho udx_2\ S$ is the mass flow rate within an elemental area $latex dx_2 S$ w and the expression $latex \left(Uu_1\right) $ is the momentum lost per mass flow rate. The product of these two terms integrated over the entire boundary layer $latex {\delta}$ is the momentum loss rate at any given $latex x{}_{1}$ position along the surface. Similar to what was done for the mass flow rate deficit, here we imagine that the momentum rate lost due to friction is equal to the mass flow rate through an area equal to $latex \left ( \delta_2S \right )$ with velocity $latex {U}$ or $latex \rho U{\delta }_2SU=\rho U^2{\delta }_2S=\int^{\delta }_0{\rho u\left(Uu_1\right)dx_2\ S}$ Dividing through by $latex \rho U^2S$ yields the definition of the displacement thickness:$latex {\delta }_2=\int^{\delta }_0{\frac{u}{U}\left(1\frac{u_1}{U_{\infty }}\right)dx_2} $
(9.5)
We can interpret this as the distance $latex {\delta}_{2}$ the surface would need to be moved in the $latex x{}_{2}$ direction into the flow if the velocity is uniformly at $latex {U}$ in order to reduce the momentum rate equivalent to the actual rate of momentum loss. We can also interpret this somewhat differently as follows. A flat surface is shown in Fig. 9.2 indicating the surface shear stress, $latex {\tau}_{s}$, acting on the fluid. A control volume with an elemental length along the surface, $latex x{}_{1}$, and extending upward to the edge of the boundary layer is shown in the figure. The sum of all of the forces acting on the fluid in the control volume consists of $latex {\tau}_{s}$, since there is zero pressure gradient along the flow, and we assume there is no frictional forces at the top edge of the boundary layer, by definition of the extent of the boundary layer. The momentum equation for this control volume becomes the sum of the forces equal to the change of momentum rate leaving compared to that entering:$latex {\tau }_sdx_1S={\left(\ \rho U^2{\delta }_2S\right)}_{x+dx_1}{\left(\ \rho U^2{\delta }_2S\right)}_{dx_1}=\left(\ \rho U^2d{\delta }_2S\right)$
(9.6)
where $latex d{\delta}_{2}$ is the change of the displacement thickness along d$latex x{}_{1}$, while noting that $latex \rho U^2S $ is constant along $latex x{}_{1}$. The negative sign is included to indicate that the stress on the fluid is in the negative $latex x{}_{1}$ direction at the surface. If we rearranging this equation and change the stress to be that on the surface caused by the fluid flow (change the sign of $latex {\tau }_s$) then:$latex {\tau }_s=\rho U^2\frac{d{\delta }_2}{dx_1} $
(9.7)
Note that this equation is valid for any particular velocity distribution that may exist in the boundary layer. Also, it indicates that the rate of change of the momentum thickness decreases along the flow direction since the surface stress decreases along the flow direction. If one can determine $latex \frac{d{\delta }_2}{dx_1}$ it is possible to find the wall surface stress, $latex {\tau }_s$.$latex u_1\frac{\partial u_1}{\partial x_1}+u_2\frac{\partial u_1}{\partial x_2}=\frac{1}{\rho }\frac{\partial P}{\partial x_1}+\nu \frac{{\partial }^2u_1}{\partial x_2} $
(9.8)
Note that $latex \nu =\frac{\mu }{\rho }$, and is the kinematic viscosity. Since the pressure gradient term can be determined from the freestream flow which is assumed to be inviscid the Euler equation can be used where we ignore gravitational effects such that:$latex \frac{1}{\rho }\frac{\partial P}{\partial x_1}=U\frac{dU}{dx_1} $
(9.9)
Note that body force terms could be included by replacing P with $latex P'=P+\rho gh$. Equation (9.8) is integrated across the boundary layer, in the $latex x{}_{2}$ direction, for a given value of $latex x{}_{1}$. The boundary conditions for this are noslip at $latex x{}_{2}$ = 0, and u = U at $latex x{}_{2} = \delta $. Before we complete this operation we examine the second term on the left hand side in Eqn. (9.8). First we use continuity for a two dimensional flow:$latex \frac{\partial u_2}{\partial x_2}=\frac{\partial u_1}{\partial x_1}$
(9.10)
which integrates to:$latex u_2=\int^{\delta }_0{\frac{\partial u_1}{\partial x_1}dx_2} $
(9.11)
and this expression is used in the second term of (9.8). Next each term in Eqn. (9.8) is integrated from [latex]y=0[/latex] to [latex]y=\delta[/latex]. For instance the second term becomes:$latex \int_{0}^{\delta}u_2\frac{\partial u_1}{\partial x_2}dx_2=\int^{\delta }_0{\frac{\partial u_1}{\partial x_2}\left(\int^{x_2}_0{\frac{\partial u_1}{\partial x_1}}dx_2\right)dx_2} $
(9.12)
This term is now integrated by parts: $latex \int^2_1{b\ da=ab{_1}^2\int^2_1{a\ db}}$ where we set $latex db=\frac{\partial u_1}{\partial x_1}dx_2, da=\frac{\partial u_1}{\partial x_2}dx_2=\partial u_1 $ and $latex \ a=u_1$ The result for this term using the boundary conditions for $latex u{}_{1}$ is:$latex U\int^{\delta }_0{\frac{\partial u_1}{\partial x_1}dx_2}+\int^{\delta }_0{u_1\frac{\partial u_1}{\partial x_1}}dx_2 $
(9.13)
We insert (9.13) into (9.8) while integrating all other terms and obtain:$latex \int^{\delta }_0{u_1\frac{\partial u_1}{\partial x_1}dx_2}U\int^{\delta }_0{\frac{\partial u_1}{\partial x_1}dx_2}+\int^{\delta }_0{u_1\frac{\partial u_1}{\partial x_1}}dx_2=\int^{\delta }_0{U\frac{dU}{dx_1}}+\int^{\delta }_0{\nu }\frac{{\partial }^2u_1}{\partial x_2}dx_2$
(9.14)
The last term is the viscous term and can be expressed in terms of the shear stress, $latex \tau =\mu \frac{\partial u_1}{\partial x_2}$ (the $latex x{}_{1}$ derivative of $latex u{}_{2}$ is taken to be small compared with the $latex x{}_{2}$ derivative of $latex u{}_{1}$, consistent with the boundary layer approximations.) So the last term is:$latex \int^{\delta }_0{\nu }\frac{{\partial }^2u_1}{\partial x_2}dx_2=\int^{\delta }_0{\frac{1}{\rho }}\frac{\partial \tau }{\partial x_2}dx_2=\frac{1}{\rho }\left[\tau \left(x_2=\delta \right)\tau \left(x_2=0\right)\right]=\frac{{\tau }_s}{\rho }$
(9.15)
We combine the first and third terms of Eqn. (9.14) as:$latex \int^{\delta }_0{2u_1\frac{\partial u_1}{\partial x_1}dx_2}=\int^{\delta }_0{\frac{\partial \left(u_1u_1\right)}{\partial x_1}dx_2}$
(9.16)
Then we rewrite the second term of Eqn. (9.14) by bringing U inside the integral since it is not a function of [latex]x_2[/latex], but may be a function of [latex]x_1[/latex]:$latex U\int^{\delta }_0{\frac{\partial u_1}{\partial x_1}dx_2=\int^{\delta }_0{\left(\frac{\partial \left(uU\right)}{\partial x_1}+u\frac{\partial U}{\partial x_1}\right)}}dx_2$
(9.17)
Getting closer to the final result we combine all these expressions in (9.14) to obtain:$latex \int^{\delta }_0{(\frac{\partial \left(u_1u_1\right)}{\partial x_1}\frac{\partial \left(u_1U\right)}{\partial x_1}+u_1\frac{dU}{dx_1}}U\frac{dU}{dx_1})dx_2 = \frac{\tau _s}{\rho }$
(9.18)
This is rearrange as:$latex \int^{\delta }_0{\frac{\partial \left(u_1(u_1U\right)}{\partial x_1}dx_2+\int^{\delta }_0{\frac{dU}{dx_1}}}\left(u_1U\right)dx_2=\frac{{\tau }_s}{\rho } $
(9.19)
Exchanging the order of integration and derivatives in the first term, changing sign of all the terms, yields:$latex \frac{\partial }{\partial x_1}\int^{\delta }_0{u_1(Uu_1)dx_2}+\frac{dU}{dx_1}\int^{\delta }_0{\left(Uu_1\right)dx_2=}\frac{{\tau }_s}{\rho }$
(9.20)
This now allows us to use the definitions of $latex {\delta}_{1}$ and $latex {\delta}_{2}$ to arrive at our final form as:[latex]\frac{d}{d x_1}\left(U^2{\delta }_2\right)+{\delta }_1U\frac{dU}{dx_1}=\frac{{\tau }_s}{\rho }[/latex]
(9.21)
We end up with a differential equation in terms of the variables $latex {\delta}_{1}$ and $latex {\delta}_{2}$ both of which are functions of $latex x{}_{1}$. Recall that U($latex x{}_{1}$) is determined by the pressure gradient through use of Eqn. (9.9), Euler's equation. When U = constant, as in flow over a flat plate this equation can be simplified to:$latex \frac{d{\delta }_2}{dx_1}=\frac{{\tau }_s}{\rho U^2} \text { or } 2\frac{d{\delta }_2}{dx_1}=\frac{{\tau }_s}{1/2\rho U^2}=c_f$
(9.22)
where the right hand side is the skin friction coefficient. The problem of determining the local surface stress distribution is now reduced to Eqn. (9.21) for flows over surfaces with a pressure gradient or (9.22) for a flat plate with no pressure gradient. Next, we will examine methods to arrive at solutions to this.$latex \eta =\frac{x_2}{\delta }$
(9.23)
being careful to note that this variable, $latex {\eta}$, is not the same as used in the Blasius definition. Next we assume that the nondimensional velocity, $latex f=\frac{u_1}{U}$, is some function of $latex {\eta}$ or:$latex f=f(\eta )$
(9.24)
One may say that it is necessary to precisely determine this function since the surface stress, $latex {\tau}_{s}$, depends on the velocity derivative at the surface. This certainly is the case for the Blasius approach. However, examining Eqns. (9.21) and (9.22) it can be seen that the surface stress is a function of the integral of the velocity distribution as expressed through the variables $latex {\delta}_{1}$ and $latex {\delta}_{2}$. Based on this we conclude that to get a solution we must obtain good values for the integrals expressed in these equations, rather than the surface derivative of the velocity. Taking this approach it is possible to assume a velocity distribution using the assumed nondimensional form in Eqn. (9.24) requiring it to satisfy the boundary conditions. One can use this distribution after transforming the integrals in (9.21) into the nondimensional form. We will use the example of flow over a flat plate, governing by Eqn. (9.22). This simplifies the details, but the same approach is used for flows with pressure gradients, but requires further input of a known pressure gradient. First looking at Eqn. (9.22) written in terms of the skin friction coefficient we have:$latex c_f=2\frac{d{\delta }_2}{dx_1}$
(9.25)
We can not evaluate $latex {\delta_2}$ since we do not know the value of $latex {\delta}$. We illustrate the process by assuming a polynomial velocity distribution based on Eqn. (9.24) variables:$latex f=a+b\eta +c{\eta }^2+d{\eta }^3+$ higher order terms
(9.26)
There are four parameters in this polynomial, $latex a,b,c \text { and } d$ that need to be determined based on boundary conditions on the velocity distribution. Notice that the higher the order of the polynomial the more conditions are needed to evaluate these coefficients. The boundary conditions written in terms of the variables, $latex f,{\delta}$ are: $latex f(0)\ =\ 0$ (no slip boundary condition at the surface) $latex f(1)\ =\ 1$ (velocity becomes the freestream value at $latex {\delta}$) $latex \frac{\partial f}{\partial \eta }=0$ (no shear at $latex {\delta}$) $latex \frac{{\partial }^2f}{\partial {\eta }^2}=0$ (no net viscous force at the surface) The last condition is a consequence of all terms in the NavierStokes equation becoming zero at the surface if the pressure gradient is zero, and [latex]u_1=u_2=0[/latex], and then transforming the viscous term into the nondimensional variables. These four conditions can now be used to determine the four coefficients in Eqn. (9.26). The results are: $latex a=0,\ \ b=3/2,\ \ c=0,\ \ d=1/2$ Inserting these into Eqn. (9.26) results in the velocity distribution:$latex f=\ \frac{3}{2}\eta \frac{1}{2}{\eta }^3$
(9.27)
This result is unique to the selection of a third order polynomial for the velocity distribution. Next, we rewrite the relationship of $latex c{}_{f}$ in terms of the nondimensional variables. First we define the momentum and displacement thicknesses in terms of $latex {\delta}$ as:$latex \frac{{\delta }_1}{\delta }=\int^1_0{\left(1\frac{u_1}{U}\right)d\eta }$
(9.28)
$latex \frac{{\delta }_2}{\delta }=\int^1_0{\frac{u_1}{U}\left(1\frac{u_1}{U}\right)d\eta }$
(9.29)
Here we have changed the integration from $latex dx_2$ to $latex d{\eta}$, where [latex]dx_2=\delta d\eta[/latex] and need to change the limits of the integral accordingly. Notice that once we have an expression for $latex f=\frac{u_1}{U}$ from our polynomial, Eqn. (9.26), then both $latex \frac{{\delta }_1}{\delta }$ and $latex \frac{{\delta }_2}{\delta }$ are just constants. The important point now is that we can rewrite Eqn. (9.25) as:$latex c_f=2\frac{d\delta_2}{dx_1}=2\frac{d}{dx_1}\left ( \left ( \frac{\delta_2}{\delta} \right )\delta \right )=2\left ( \frac{\delta_2}{\delta} \right )\frac{d\delta}{dx_1}$
(9.30)
So we have replaced having to know $latex {\delta}_2 = f(x_1)$ to needing to know $latex {\delta} = f(x_1)$. We find the latter relationship using our function [latex]f(\eta)[/latex] as: $latex {\tau }_s=\mu \frac{\partial u_1}{\partial x_2}=\mu \frac{\partial u_1}{\partial \eta }\frac{\partial \eta }{\partial x_2}=\frac{\mu U}{\delta }\frac{\partial f}{\partial \eta }=\frac{\mu U}{\delta }f'(0)$$latex c_f=\frac{{\tau }_s}{\frac{1}{2}\rho U^2}=2\frac{{\delta }_2}{\delta }\frac{d\delta }{dx_1}=\frac{\mu U}{\delta \frac{1}{2}\rho U^2}f'(0)$
(9.31)
This last relationship is a differential equation for $latex {\delta} = {f(x_1)}$ since $latex \frac{{\delta }_2}{\delta } $ and $latex f'(0)$ are both constants known from Eqn. (9.26) using the boundary conditions to obtain Eqn. (9.27). Notice that these two constants are dependent on the order of the polynomial selected. Now Eqn. (9.31) can be integrated to obtain:$latex \frac{\delta }{x_1}=\frac{{\left(2f'(0)\frac{\delta }{{\delta }_2}\ \right)}^{1/2}}{{Re_{x_1}}^{1/2}}=\frac{C}{{Re_{x_1}}^{1/2}} $
(9.32)
where $latex C={\left(2f'(0)\frac{\delta }{{\delta }_2}\ \right)}^{1/2}$ and is a constant dependent on the order of the polynomial selected. Knowing $latex {\delta = f(x_1)}$ we can find $latex \frac{d\delta }{dx_1}$ and then $latex c{}_{f}$ using Eqn. (9.31). The result is:$latex c_f=\frac{{\delta }_2}{\delta }\frac{C}{{Re_{x_1}}^{1/2}} $
(9.33)
$latex C_f=\frac{F_D}{\frac{1}{2}\rho U^2A_s}=\frac{1}{\frac{1}{2}\rho U^2A_s}\int^L_0{{\tau }_sdx_1w}=\frac{1}{L}\int^L_0{c_fdx_1=2c_f(x_1=L)}$
(9.35)
where the last designation, $latex 2c_f(x_1=L)$, means two times the value of the skin friction coefficient evaluated at $latex x{}_{1} = L$. The result assuming the third order polynomial for the nondimensional velocity profile given above is:$latex C_f=\frac{1.2928}{{Re_L}^{1/2}} $
(9.36)
The results when compared with experimental data for the appropriate Reynolds number conditions become very close to these predicted values. Consequently, this procedure is a fairly easy yet surprisingly accurate means to predict friction drag forces on flat surfaces. When pressure gradients come into play, such as surfaces with some degree of curvature it is more of a challenge to obtain a good velocity profile since for this case at the surface the net viscous forces are not zero but are balanced by the pressure gradient. The velocity profiles need to reflect this condition and there are a number of procedures available to address this. The interested reader is encouraged to see F.W. White's book, Viscous Fluid Flow for some examples.]]>(c)
[caption id="attachment_290" align="aligncenter" width="300"] Fig. 10.1. Illustrations of turbulence flow visualization; (a) smoke visualization of flow over an inclined wing indicating flow separation from the leading edge which transitions to turbulence; (b) smoke visualization from a turbulent jet injecting flow into a cross flow in a wind tunnel indicating large and small scale fluctuations; (c) particle image velocimetry (PIV) of flow in a porous media (through a random packing of spherical beads) showing instantaneous velocity vector distribution indicative of small scale swirling motion.[/caption] Shown in Fig. 10.2 is a typical time history signal of a measured velocity component say from a hot wire anemometer at a given location within the flow. Although this signal may seem random in nature, it has been learned that there is a lot of structure to these kinds of signals. One thing to notice is that there are time events happening over many different time scales. There are small scale (or short) events with relative small fluctuating values coupled with longer time varying, higher amplitude fluctuations. The short time events seem to be superimposed over the longer events. To help illustrate this variation of time scales Fig. 10.2 also illustrates the ``energy spectrum'' of the time series signal. Here the mean square value of the fluctuations are decomposed by using a Fourier transformation of the signal. The mean square fluctuations have units of velocity squared that can be interpreted as the kinetic energy per mass of fluid. This plot identifies the frequency dependency of the energy associated with the signal. We see that the larger amplitude fluctuations occur at the longer time scales (lower frequency) and the smaller amplitude fluctuations occur at the shorter time scales (higher frequency). Turbulent flows are characterized by a broad spectrum within the signal typically with large amplitude (or energy) fluctuations having low frequency and low amplitude high frequency components. Looking at Fig. 10.2 we see that the longer time scales have several orders of magnitude greater energy than the smallest scales. This distribution has implications concerning the mixing characteristics of turbulent flows such as in the ability to promote high reaction rates in combustion and other chemical reactions. In contrast the energy spectrum of a flow that has a single oscillating frequency would have an energy spectrum that would appear as a spike at the frequency of oscillation and would not be referred to as turbulent. [caption id="attachment_568" align="aligncenter" width="2250"] Fig. 10.2: (a) Typical velocity signal of a single component of velocity in turbulent flow measured as a given location, notice the large and small scale fluctuations that are superimposed on the signal. (b) energy spectrum on a loglog scale of the mean square velocity fluctuations verses frequency of fluctuations, notice the large range of energy and frequency that occurs.[/caption] The integration of the energy spectrum over all frequencies yields the total energy of the fluctuations in the signal. This is equivalent to the time average of the square of the fluctuations (the mean square value) $latex \overline{{u'}^2}$, where the fluctuating velocity is $latex u'=\left(u\overline{u}\right)$ with $latex \overline{u}$ being the time average value of the velocity and u is the instantaneous value of the velocity. This can be used to define the ``turbulent intensity'' as: $latex \frac{{\overline{\left({u'}^2\right)}}^{1/2}}{\overline{u}}$ = ``turbulent intensity'' Notice that $latex \overline{{u'}^2}$ is obtained for a single velocity component, say $latex u_{1}$, so this represents the turbulent intensity of one component. If we calculate the mean square value for each component the total energy, known as the turbulent kinetic energy (``tke'') is:$latex tke\ =\frac{1}{2}\ \left(\left(\overline{{{u'}_1}^2}\right)+\left(\overline{{{u'}_2}^2}\right)+\left(\overline{{{u'}_3}^2}\right)\right)$
(10.1)
The greater the tke typically the more intense is the turbulence. This measure of turbulence can be measured at each location within the flow to determine the spatial distribution of turbulence in a time averaged sense.$latex u=\left(u'+\overline{u}\right)$
(10.2)
where the primed quantity is the fluctuating value, which is time dependent, and $latex \overline{u}$ is the time average value. It should be noted that the time average of the fluctuating part is by definition zero. Each of the time dependent terms in the conservation of mass (or continuity equation) and NavierStokes equations can be replaced with the right hand side of Eqn. (10.2). Then each term is time averaged. For instance in the convective acceleration term:$latex \overline{u_j\frac{\partial u_i}{\partial x_j}}=\overline{\left({u_j}'+\overline{u_j}\right)\frac{\partial \left({u_i}'+\overline{u_i}\right)}{\partial x_j}}=\overline{{u_j}'\frac{\partial {u_i}'}{\partial x_j}}+\overline{u_j}\frac{\partial \overline{u_i}}{\partial x_j}$
(10.3)
Where the two terms containing products of the mean and fluctuating parts, when time averaged, are zero since the time average of the fluctuations in zero. The last term in Eqn. (10.3) is merely the convective acceleration using the time average velocities. Notice that the time averaging and the spatial derivatives can be interchanged. The first term on the very right hand side is the time average of the convective acceleration of the fluctuating velocities. This term is not necessarily zero for a term consisting of the product of two or more time varying signals. This is because the average of the products of time varying signals is not equal to the product of the averages of those signals. Let's examine this term a bit closer for incompressible flow we can write.$latex \overline{{u_j}'\frac{\partial {u_i}'}{\partial x_j}}=\frac{\partial \overline{\left({u_i}'{u_j}'\right)}}{\partial x_j}$
(10.4)
Where, for incompressible flow, we are able to bring the velocity $latex u_j$ inside the derivative with respect to $latex x_j$ from the continuity condition for incompressible flows. This term can be interpreted as the spatial derivative of the correlation function $latex \left(\overline{{u_i}'{u_j}'}\right)$ that is the time averaged product of the various velocity components. This is a second order tensor with six independent components since it is a symmetric tensor. This quantity, $latex \left(\overline{{u_i}'{u_j}'}\right)$ represents the contribution of turbulence to the mean momentum of the flow. That is to say, we have found a term that can change or adjust the mean momentum distribution in the flow caused by averaging the Reynolds decomposition of convective acceleration terms. We will see how this term affects the wall shear stress acting on the fluid for turbulent flow conditions. Based on our earlier introduction of the turbulent kinetic energy, tke, one might expect that increased levels of tke would be associated with an increased magnitude of $latex \left(\overline{{u_i}'{u_j}'}\right)$ and increase the effect on the mean momentum distribution. How this occurs and under what conditions it occurs is not a simply thing to determine. This is because the relationship of the energy content of the fluid fluctuations with the momentum transport properties of the turbulence is a set of complicated interactions. Never the less it is possible to observe trends of these interactions and to formulate models of how they work. We further manipulate the time averaged NavierStokes equations by bringing the new fluctuating term $latex \left(\overline{{u_i}'{u_j}'}\right)$ over to the right hand side of the equation and group it with the viscous stress term:$latex \rho \overline{u_j}\frac{\partial \overline{u_i}}{\partial x_j}=\frac{\partial \overline{P}}{\partial x_j}+{\rho g}_i+\left(\frac{\partial }{\partial x_j}\left(\mu \frac{\partial \overline{u_i}}{\partial x_j}\rho \left(\overline{{u_i}'{u_j}'}\right)\right)\right)$
(10.5)
Here we can see that the viscous term is combined with the momentum altering fluctuating term while all other terms are similar to the laminar flow equation except here we use the time average value of the variables. The last term has two contributions, the mean viscous stress term $latex \ \mu \frac{\partial \overline{u_i}}{\partial x_j}$, and the turbulence term. If we interprete the turbulence term as one that modifies the momentum distribution in the flow we can identify it as a ``stress'' and call it the ``turbulent stress'', $latex {\tau }_t $ and write it as:$latex {\tau }_t=\rho \left(\overline{{u_i}'{u_j}'}\right)$
(10.6)
In honor of Reynolds, who developed the Reynolds decomposition leading to this results this is often labeled the ``Reynolds stress''. We can write the total stress within the flow as:$latex {\tau }_{total}={\tau }_{viscous}+{\tau }_{turbulent}=\ \mu \frac{\partial \overline{u_i}}{\partial x_j}\rho \left(\overline{{u_i}'{u_j}'}\right) $
(10.7)
In a number of situations the viscous stress is actually a lot smaller than the turbulent stress and can be ignored except maybe in certain region of the flow, say near a surface or object where the fluctuating quantities are dampened.$latex \rho \overline{u_1}\frac{\partial \overline{u_1}}{\partial x_1}\ \ \ \rho \overline{u_2}\frac{\partial \overline{u_1}}{\partial x_2}=\frac{\partial \overline{P}}{\partial x_1}+\left(\frac{\partial }{\partial x_2}\left(\mu \frac{\partial \overline{u_1}}{\partial x_2}\rho \left(\overline{{u_1}'{u_2}'}\right)\right)\right)$
(10.8)
The effect on the mean momentum distribution, caused by turbulence, then is brought about by the cross correlation contained in the Reynolds stress. Referring to Fig. 10.3 we attempt to provide some rationale for how this correlation may impact the momentum distribution and consequently the wall shear stress. At an instant in time imagine a positive, upward value of $latex u{}_{2}$. This moves fluid with a lower value of $latex u{}_{1}$ into a region of higher $latex u{}_{1}$ on average. Therefore the instantaneous value of $latex {u_1}'{u_2}'$ will be a negative number. Similarly if there is a negative value of $latex u{}_{2}$ it will move high $latex u{}_{1}$ velocity fluid into a slower moving region again causing $latex {u_1}'{u_2}'$ to be a negative number. Assuming this to be the case, one may assume that on average the value of $latex \overline{{u_1}'{u_2}'\ }$ is a negative number which means the added turbulent stress is positive, see Eqn. (10.6). This added stress slows down higher speed fluid and speeds up lower speed fluid. The result is a much flatter velocity or momentum distribution across the flow. However, the noslip condition still holds true. This implies that on average the flatter velocity profile will have higher velocity gradients near the surface, such that the overall frictional force is larger. The net result is that a turbulent boundary layer is expected to have larger surface viscous stresses.
$latex {v_f}=\sqrt{\frac{{\tau }_s}{\rho }}$
(10.9)
Notice that it is not a specific velocity within the flow but has units of velocity and is proportional to the velocity gradient at the surface through the value of the wall shear stress, [latex]\tau_s[/latex]. Since it is reasonable to note that the surface stress depends on the freestream velocity, one can also anticipate that [latex]{v_f}=f(U)[/latex]. From this we form the nondimensional velocity components:$latex {u_1}^+=\frac{\overline{u_1}}{{v_f}}$
(10.10a)
$latex {u_2}^+=\frac{\overline{u_2}}{{v_f}}$
(10.10b)
As we will see in the next section, this scaling of velocity can be used to arrive at an accurate near wall velocity distribution in nondimensional variables. However, if we choose to try to determine the wall shear stress using an integral analysis similar to what was done in laminar flow, then we can use the following. First we assume a functional form:$latex \frac{\overline{u_1}}{U}={f}{(\eta)}$
(10.11)
In contrast to the polynomial distribution we used in laminar flow. We assume a power law instead of the polynomial (since it has been shown to match experimental data well) as:$latex \frac{\overline{u_1}}{U}=\ {C\eta }^{1/n}$
(10.12)
where n is a positive integer and C is a constant. However, if we note that $latex \frac{\overline{u_1}}{U}=1$ at $latex \eta$ =1 then C = 1. Using this expression in the definitions for $latex {\delta}_{1}$ and $latex {\delta}_{2}$ we get:$latex \frac{{\delta }_1}{\delta }=\ \frac{1}{n+1}$
(10.13)
$latex \frac{{\delta }_2}{\delta }=\ \frac{n}{\left(n+1\right)\left(n+2\right)} $
(10.14)
Experimental data shows that n is a weak function of Reynolds number, such that as Reynolds number increases the velocity profile becomes flatter and n increases. For moderately large Reynolds numbers n $latex \mathrm{\sim}$ 7 matches data very well. To go any further with this we need a relationship between $latex {\tau}_{s}$ and [latex]\delta[/latex], similar to what was done with laminar flow, based on the nondimensional velocity distribution. For this we use an alternative scaling with a power law assumption given by:$latex {u_1}^+=f_2(\frac{x_2{v_f}}{\nu })=\ C_2{\left(\frac{x_2{v_f}}{\nu }\right)\ }^{1/n} $
(10.15)
where $latex C{}_{2}$ is a constant and $latex \frac{\nu}{v_f}$ is a new length scale. At $latex x{}_{2} = \delta$, we know that $latex \frac{\overline{u_1}}{U}=1 $ so we can write:$latex \frac{U}{{v_f}}=\ C_2{\left(\frac{\delta {v_f}}{\nu }\right)\ }^{1/n} $
(10.16)
This provides the relationship between $latex {\tau}_{s}$ and [latex]\delta[/latex], since $latex {\tau}_{s}$ is contained in $latex {v_f}$, see Eqn. (10.9). After inserting the definition of $latex {v_f}$ and a little algebra we obtain:$latex \frac{{\tau }_s}{\rho U^2}={C_2}^{\frac{2n}{n+1}}{Re_{\delta }}^{\frac{2}{n+1}} $
(10.17)
where $latex Re_{\delta }=\frac{U\delta }{\nu }$ which represents a Reynolds number based on the boundary layer thickness. Now we use the momentum integral equation that was derived for laminar flow, which is also valid for turbulent flow that yields for a flat plate:$latex c_f=2\frac{{\delta }_2}{\delta }\frac{d\delta }{dx_1}=\frac{{\tau }_s}{\frac{1}{2}\rho U^2}$
(10.18)
Noting that $latex \ \frac{{\delta }_2}{\delta }$ is a constant for a given value of n, and eliminating [latex]\tau_s[/latex] between Eqn. (10.17) and (10.18) and solving for [latex]\delta[/latex] we obtain:$latex \frac{\delta }{x_1}=\ C_T{Re_{x_1}}^{\frac{2}{\left(n+3\right)}} $
(10.19)
where $latex C_T={\left[\frac{{{\left(\frac{n+3}{n+1}\right)C}_2}^{\frac{2n}{n+1}}}{{{\delta }_2}/{{\delta }}}\right]}^{\frac{n+1}{n+3}}$ that is fully determined by the constant $latex C{}_{2}$ and the value of n. Based on experimental data, much of it analyzed by Blasius for turbulent flow results in:$latex C_2=8.74 $
(10.20)
Now we have Eqn. (10.19) fully determined so that we can find the derivative $latex \frac{d\delta }{dx_1}$ and inserts this into the expression for the skin friction, Eqn. (10.18). For the case of n = 7 we have:$latex C_T=0.37$ $latex \frac{{\delta }_2}{x_1}=\ 0.036\ {Re_{x_1}}^{1/5}$ $latex c_f=\ 0.0577\ {Re_{x_1}}^{1/5} $
(10.21)
The determination of the total viscous drag force acting on the surface of length L and width S is obtained exactly as was done for laminar flow, by integration of the surface shear stress over the area. The results for the case of n = 7 is:$latex C_f==\ 0.072\ {Re_L}^{1/5} $
(10.22)
$latex F_D=C_f\frac{1}{2}\rho U^2LS $
(10.23)
$latex x_{transition}=\frac{\left(5x10^5\right)\nu }{U} $
(10.24)
Once the transition location is found then the integration over the entire surface becomes:$latex F_D=\int^{x_{transition}}_o{{\tau }_sdx_1\ w}\ +\int^L_{x_{transition}}{{\tau }_sdx_1 S}$
(10.25)
Each term is now divided by $latex \frac{1}{2}\rho {U_{\infty }}^2(LS)$ to convert this equation into a nondimensional equation for $latex C{}_{D}$. The integrand for each term will be the skin friction coefficient, $latex c{}_{f}$ and $latex dx{}_{1}$ becomes $latex d(x/L)$. There is a problem here though, the expression for ${c}_{f}$ for the turbulent part (the second integral) can not be solved directly because the skin friction expression found previously assumes the flow begins turbulent at the leading edge. The solution to this problem is to recast the above integral by adding and subtracting the contribution to the total force from $latex x = 0 \text{ to } x = x{}_{transition}$ based on turbulent flow. When this is added to the second term in Equation (10.25) the result is an expression for the force assuming turbulent flow from x = 0. The subtraction of this term in the first integral results in the difference between the laminar and turbulent generated force up until $latex x_{transition}$. Rearranging the terms results in$latex C_D=\ \int^L_0{c_{f,\ turbulent}d\left(\frac{x}{L}\right)}\int^{x_{transition}}_o{\left(c_{f,\ turbuent}c_{f,\ laminar}\right)d\left(\frac{x}{L}\right)}$
(10.26)
By using the transitional location given in Eqn. (10.24) the second integral in Eqn. (10.26) can be found directly using the expressions for $latex c{}_{f}$. Writing this using $latex Re{}_{L}$ (the Reynolds number evaluated at the end of the surface) results in:$latex C_D=\ \frac{0.072}{{Re_L}^{1/5}}\frac{1674}{Re_L}$
(10.27)
It should be noted that the above expression assumes n = 7. Also, if transition occurs at a slightly different Reynolds number the coefficient in the numerator of the second term in this equation will vary. Often, the numerator of the second term is 1750 to better fit data.$latex {\tau }_{total\ }=\mu \frac{\partial \overline{u_1}}{\partial x_2}\rho \left(\overline{{u_1}'{u_2}'}\right)$
(10.28)
We assume that this stress is a constant which can be written in terms of the friction velocity, $latex {v_f}$. Then we model this term based on an effective viscosity, $latex {\nu }_t$, that accounts for both the viscous and turbulent contributions which results in the total stress in terms of the friction velocity as:$latex \frac{{\tau }_{total\ }}{\rho }=\ {{v_f}}^2={\nu }_t\frac{\partial \overline{u_1}}{\partial x_2} $
(10.29)
The effective viscosity (written here as a kinematic viscosity) has units of velocity times length. Taking the friction velocity as the appropriate velocity scale and $latex x{}_{2}$ as the appropriate length scale (since the viscous effects can only extend away from the wall some amount $latex x{}_{2}$) then we model the effective viscosity as:$latex {\nu }_t=\kappa {v_f}x_2 $
(10.30)
where $latex \kappa $ is an unknown constant. Inserting Eqn. (10.30) into (10.29) we have:$latex {{v_f}}^2=\kappa {v_f}x_2\frac{\partial \overline{u_1}}{\partial x_2} $
(10.31x)
This is now integrated in the [latex]x_2[/latex] direction for a constant value of $latex {v_f}$ to obtain a velocity profile as:$latex \frac{\overline{u_1}}{{v_f}}=\frac{1}{\kappa }ln\ x_2+c$
(10.32)
This can easily be transformed into the following:$latex \frac{\overline{u_1}}{{v_f}}=\frac{1}{\kappa }ln\ \frac{x_2{v_f}}{\nu }+B$
(10.33)
where $latex \kappa $ and B are constants determined from experimental results and shown to be $latex \kappa =0.41$ and B = 5. Actually there are also theoretical arguments that can be made to show that $latex \kappa =0.41$. This constant stress layer result is only valid fairly close to the surface, but has been shown over a wide range of conditions (different fluids, different Reynolds numbers, etc.) to be very accurate. The nondimensional variable, $latex \frac{x_2{v_f}}{\nu }$, is often given the symbol $latex {x_2}^+\ or\ y^+$, and $latex \frac{\overline{u_1}}{{v_f}}=u^+$, resulting in a nondimensional velocity profile valid close to the wall that is expressed in terms of the wall shear stress contained in the variable [latex]u^t[/latex] and [latex]x_2^t[/latex].]]>This work is licensed by James Liburdy under a Creative Commons AttributionNonCommericalShare Alike 4.0 International License (CC BYNCSA)
You are free to: Share — copy and redistribute the material in any medium or format Adapt — remix, transform, and build upon the material
The licensor cannot revoke these freedoms as long as you follow the license terms.
Under the following terms: Attribution — You must give appropriate credit, provide a link to the license, and indicate if changes were made. You may do so in any reasonable manner, but not in any way that suggests the licensor endorses you or your use. NonCommercial — You may not use the material for commercial purposes. ShareAlike — If you remix, transform, or build upon the material, you must distribute your contributions under the same license as the original. No additional restrictions — You may not apply legal terms or technological measures that legally restrict others from doing anything the license permits.
]]>Version  Date  Change Made  Location in text 
1.0  9/16/2021  Publication  
1.01  12/09/2021  Updated Eq. 10.19  Introduction to turbulence effects 
[latex]F=\emptyset+i\psi[/latex]
(5.31)
The derivative of [latex]F[/latex] in terms of [latex]z[/latex] is defined as:[latex]W\left(z\right)=\frac{dF}{dz}[/latex]
(5.32)
Also: [latex]\frac{dF}{dz}=\frac{dF}{dz}\frac{\partial z}{\partial x}=\frac{\partial F}{\partial x}[/latex] Consequently, inserting the definition of [latex]F[/latex] we have[latex]W\left(z\right)=\frac{\partial\emptyset}{\partial x}+i\frac{\partial\psi}{\partial x}=u−iv[/latex]
(5.33)
This result shows that [latex]W(z)[/latex] represents the “complex velocity” of the flow and is determined by the real velocity components [latex]u,v[/latex]. It is often advantageous to use cylindrical coordinates, expressed as [latex]r,{\theta}[/latex] in two dimensions. The transformation from [latex]x,y[/latex] to [latex]r,{\theta}[/latex] is: [latex]u=u_rcos\theta−u_\theta sin\theta[/latex] [latex]v=u_rsin\theta+u_\theta cos\theta[/latex] The expression for the complex velocity is then (by direct substitution and using the identity of [latex]\cos{\theta}−isin\ \theta=e^{−i\theta}[/latex] :[latex]W=\left(u_r−iu_\theta\right)e^{−i\theta}[/latex]
(5.34)
The use of the complex variable representation in terms of [latex]F[/latex] can then be converted back into the physical space velocity components, [latex]u,v[/latex], through its derivative relative to the variable [latex]z[/latex]. The magnitude of the local velocity vector is found from [latex]W[/latex] by taking the square root of [latex]WW*[/latex], where [latex]W*[/latex] is the complex conjugate of [latex]W[/latex]. Once the velocity is determined then the pressure field is directly determined from the Bernoulli Equation. Some examples will help illustrate the use of complex variables to represent rather simple flows.[latex]F=Az^n=Ar^n(\cos{n\theta}+i\sin{n\theta)}[/latex]
(5.35)
Where “[latex]n[/latex]” represents a parameter to be specified, usually as a constant for a given flow, and “[latex]A[/latex]” is a constant for a specific flow, and shown below to be proportional to the velocity far from the sector, representative of the flow into the sector. An example of this class of flows is uniform flow where [latex]A=U[/latex] and [latex]n=1[/latex]. Note that uniform flow can be thought of as flow over a flat surface that is parallel to the uniform flow direction (there are no viscous forces so the no slip boundary condition does not hold). The flat surface can be thought of as a sector with an intersection of two lines with the lines being parallel. Or, the angle between the two lines is [latex]\pi[/latex]. The general expression above corresponds to the following expression for velocity potential and streamfunction: [latex]\emptyset=Ar^n\cos{(n\theta)}[/latex] [latex]\psi=Ar^n\sin(n\theta)[/latex] Using these two expressions we can easily locate lines of constant streamfunction can be found, and in particular when [latex]\psi=0[/latex]: when [latex]\theta=0[/latex] and [latex]\pi/n[/latex]. The flow between these two radial lines represents flow in a “sector”, as seen in the Fig. 5.8, below. [caption id="attachment_733" align="aligncenter" width="465"] Fig. 5.8 Twodimensional flow in a sector of arbitrary angle between two straight lines; the streamfunction is chosen to be zero along these two lines.[/caption] Using the definition of the complex velocity potential, [latex]W[/latex]: [latex]W=\frac{dF}{dz}=nAz^{n−1}=nAr^{n−1}e^{i(n−1)\theta}[/latex] Expanding the exponential reveals that the velocity components in [latex]r,\theta[/latex] coordinates are: [latex]u_r=nAr^{n−1}\cos(n\theta)[/latex] [latex]u_\theta=nAr^{n−1}\sin(n\theta)[/latex] By inserting different values for the parameter “[latex]n[/latex]” different flows can be simulated. Some examples are given below. [latex]n = 1:[/latex] [latex]u=A[/latex] and [latex]v = 0[/latex] (uniform flow) [latex]n = 2:[/latex] [latex]u=2Ar\cos{(\theta)}[/latex] [latex]v=2Ar\sin{(\theta)}[/latex] or in Cartesian coordinates: [latex]u=2Ax;\ v=2Ay[/latex] [latex]\psi=−Ar^2\sin{\left(2\theta\right)}=−Axy[/latex] or for a constant streamfunction (along a streamline) we have: [latex]xy = C[/latex] (where [latex]C[/latex] is a constant determined by the value of [latex]\psi[/latex] at a given location). This latter flow, for a range of values of [latex]C[/latex] yields flow into a right angled corner where the location ([latex]0,0[/latex]) is the position of the corner. Note that [latex]\psi=0[/latex] along the axes. The reader is encouraged to plot this function for different values of streamfunction to visualize the flow. The corresponding flow into a corner of arbitrary angle, [latex]\alpha[/latex], can be expressed as: [latex]F=Az^{\pi / \alpha}[/latex] Note that [latex]Az^n[/latex] is equivalent to [latex]Ar^n [cos(n\theta)+isin(n\theta)][/latex] which can then be restated in terms of [latex]\pi / \alpha=n[/latex], so that for [latex]\alpha=\pi /3[/latex] we have [latex]n=3[/latex]. This then relates [latex]n[/latex] to the turning angle of the corner flow. For instance, if [latex]n=3/2[/latex] we have an angle of [latex]2\pi /3[/latex] which is greater than 90^{o}. We can also have flow over a flat plate parallel with the flow if [latex]n=1[/latex].[latex]F=Ae^{−i\beta}z^n=Ar^ne^{i(n\theta−\beta)}[/latex]
(5.36)
Noting that the real part of F is the velocity potential and the imaginary part is the streamfunction we have: [latex]\emptyset=Ar^n\cos(n\theta−\beta)[/latex] [latex]\psi=Ar^n\sin(n\theta−\beta)[/latex] This streamfunction expression shows us that when using the complex representation for the constant multiplying [latex]z^n[/latex] there is introduced a rotation of the streamlines through angle [latex]\beta[/latex] using the same coordinate system when compared to the streamfunction when the coefficient is real (not complex).[latex]U\sin{\alpha}=\ \int_{0}^{c}{\frac{\gamma(x)}{2\pi(x_o−x)}dx}[/latex]
(6.19)
The above condition needs to be satisfied for the given distribution [latex]\gamma(x)[/latex] along the foil. However, we need a boundary condition at [latex]x = c[/latex] that specifies [latex]\gamma(c)[/latex]. This is done by introducing what is known as the Kutta Condition. This is based on the observed flow condition at the trailing edge of the foil, assuming the geometry is such that the flow from the top and bottom sides meet at a point (the trailing edge) and that flow exits the foil “smoothly” from the top and bottom surfaces such that there is no vorticity as the flow leaves the trailing edge. Physically this would occur if there is no rollup of the flow say from the bottom to the top surface. If the velocities on top and bottom at the trailing edge are identical, then the pressures are identical by the Bernoulli equation, and the flow has no tendency to swirl. So, this condition is expressed as [latex]\gamma(x=c) = 0[/latex]. It can then be shown that the integration in Eqn. (6.19) is found to satisfy the Kutta Condition for the following distribution of the vortex sheet:[latex]\gamma(x) = 2U \sin \alpha (c/x 1) [/latex]
(6.20)
(See the extra comments on the Kutta condition below) With this distribution of [latex]\gamma(x)[/latex] the total circulation, [latex]\Gamma[/latex], can be determined. Note that this value of [latex]\Gamma[/latex] depends on the angle of attack, [latex]\alpha[/latex], and the free stream velocity, [latex]U[/latex]. We can now use the previously given result that the lift force, [latex]L[/latex] , is dependent on the circulation that was shown for a rotating cylinder in a freestream, [latex]L= \rho U \Gamma[/latex]. This is further discussed below and is found to be valid for airfoils in irrotational flow. The integration of [latex]\gamma(x))[/latex] given in Eqn. (6.19) results in the following Lift expression for an airfoil: [latex]L=2 \pi \sin \alpha \rho U^2 c[/latex] Or by defining the nondimensional lift coefficient as the lift force divided by [latex]1/2 \rho U^2[/latex]:[latex]C_L=2 \pi \sin \alpha[/latex]
(6.21)
This surprisingly simple expression is a major finding of the thin air foil theory in the description of the lift force as determined by the foil angle of attack. This is based on a twodimensional steady flow (or infinitely wide foil) with very little camber with viscous effects ignored. Some additional features and effects will be added to this equation later. Figure 6.7, shows the lift coefficient versus angle of attack. At large angles of attack the results of eqn. (6.21) are not valid. This is because the angle is so large that flow over the leading edge cannot negotiate the needed turn to stay “attached” to the foil surface. Rather the flow “separates” from following the surface contour. This results in a recirculating flow region between the streamlines of flow and the body surface. This can be seen in Fig. 4.5 illustrating flow going over the leading edge of a thin airfoil at a fairly large angle of attack. The primary limitation on the application of thin airfoil theory is the assumption of inviscid flow (which is actually has been shown to be a fairly valid approximation since all of the friction induced vorticity is accounted for with the vortex sheet). Also, the foil must be thin (maximum thickness on the order of ten percent of the chord length). This is not necessarily that restrictive since even airfoils with thicknesses of over 10% of the chord follow thin airfoil theory fairly well. With some modifications to the basic lift result of Eqn. (6.19) thickness effects, relatively large camber and threedimensional effects can be accounted for, as given below. A generic plot of the lift coefficient, [latex]C_L[/latex], versus angle of attack, a, is given in Fig. 6.7, as well as results for a specific airfoil, the NACA 0012 compared with Eqn. (6.21). Very good agreement between data and theory are shown up to an angle of attack of about 16^{o}. [caption id="attachment_791" align="aligncenter" width="708"] Fig. 6.7. Illustration of a typical Lift coefficient versus angle of attack curve, left, and comparison with Eq1n. (6.19) and data for a NACA 0012 air foil.[/caption]