In [1]:
import numpy
from matplotlib import pyplot
%matplotlib inline
from VortexPanel import Panel,solve_gamma_kutta,plot_flow,make_jukowski,make_circle
alpha = numpy.pi/16
N = 64
foil = make_jukowski(N)
solve_gamma_kutta(foil,alpha)
plot_flow(foil,alpha)
From the previous notebook we know the function march doesn't need the details of the geometry, but it does need:
The viscosity is obvious, but we'll need to get the other variables from the potential flow solution.
Next, let's get $s$. Note that a body will form two boundary layers, one on each side. We need to identify the starting point of these two flow regions.
This makes it straightforward to split the body into the two boundary layer sections:
In [2]:
# split panels into two sections based on the flow velocity
def split_panels(panels):
# positive velocity defines `top` BL
top = [p for p in panels if p.gamma<=0]
# negative defines the `bottom`
bottom = [p for p in panels if p.gamma>=0]
# reverse array so panel[0] is stagnation
bottom = bottom[::-1]
return top,bottom
In [3]:
foil_top,foil_bottom = split_panels(foil)
Note that we changed the direction of the bottom array so that it runs from the stagnation point to the trailing edge, in accordance with the flow direction.
Lets plot them to make sure we got it right:
In [4]:
# plot panels with labels
def plot_segment(panels):
pyplot.figure(figsize=(10,2))
pyplot.axis([-1.2,1.2,-.3,.3])
for i,p_i in enumerate(panels):
p_i.plot()
if i%10 == 0:
pyplot.scatter(p_i.xc,p_i.yc)
pyplot.text(p_i.xc,p_i.yc+0.05,
'panel ['+'%i'%i+']',fontsize=12)
In [5]:
plot_segment(foil_top)
In [6]:
plot_segment(foil_bottom)
In [7]:
# Pohlhausen Boundary Layer class
class Pohlhausen:
def __init__(self,panels,nu):
self.u_e = [abs(p.gamma) for p in panels] # tangential velocity
self.s = numpy.empty_like(self.u_e) # initialize distance array
self.s[0] = panels[0].S
for i in range(len(self.s)-1): # fill distance array
self.s[i+1] = self.s[i]+panels[i].S+panels[i+1].S
ds = numpy.gradient(self.s)
self.du_e = numpy.gradient(self.u_e,ds) # compute velocity gradient
self.nu = nu # kinematic viscosity
self.xc = [p.xc for p in panels] # x and ...
self.yc = [p.yc for p in panels] # y locations
def march(self):
# march down the boundary layer until separation
from BoundaryLayer import march
self.delta,self.lam,self.iSep = march(self.s,self.u_e,self.du_e,self.nu)
# interpolate values at the separation point
def sep_interp(y): return numpy.interp( # interpolate function
12,-self.lam[self.iSep:self.iSep+2],y[self.iSep:self.iSep+2])
self.s_sep = sep_interp(self.s)
self.u_e_sep = sep_interp(self.u_e)
self.x_sep = sep_interp(self.xc)
self.y_sep = sep_interp(self.yc)
self.delta_sep = sep_interp(self.delta)
A few implementation notes:
numpy.gradient function is used to get $u_e'$. Pohlhausen.march calls march from the last notebook and then interpolates linearly to get values at the separation point.
In [8]:
circle = make_circle(N) # set-up circle
solve_gamma_kutta(circle) # solve flow
top,bottom = split_panels(circle) # split panels
nu = 1e-5 # set viscosity
top = Pohlhausen(top,nu) # get BL inputs
u_e = 2.*numpy.sin(top.s) # analytic u_e
du_e = 2.*numpy.cos(top.s) # analytic du_e
# compare the boundary layer inputs
pyplot.xlabel(r"$s$",fontsize=16)
pyplot.plot(top.s,top.u_e, lw=2, label=r'Panel $u_e$')
pyplot.plot(top.s,u_e, lw=2, label=r'Analytic $u_e$')
pyplot.plot(top.s,top.du_e, lw=2, label=r"Panel $u_e'$")
pyplot.plot(top.s,du_e, lw=2, label=r"Analytic $u_e'$")
pyplot.legend(loc='lower left')
Out[8]:
Those look very good. Now lets march and look at $\delta$ and the separation point.
In [9]:
top.march() # solve the boundary layer flow
i = top.iSep+2 # last point to plot
# plot the boundary layer thicknes and separation point
pyplot.ylabel(r'$\delta$', fontsize=16)
pyplot.xlabel(r'$s$', fontsize=16)
pyplot.plot(top.s[:i],top.delta[:i],lw=2)
pyplot.scatter(top.s_sep,top.delta_sep, s=100, c='r')
pyplot.text(top.s_sep-0.6,top.delta_sep,
' separation \n s='+'%.2f' % top.s_sep,fontsize=12)
Out[9]:
Same answer as the previous notebook. Good.
Now that we know the code is working, lets write a function to set-up, solve, and plot the separation points for the boundary layer flow.
In [10]:
def solve_plot_boundary_layers(panels,alpha=0,nu=1e-5):
# split the panels
top_panels,bottom_panels = split_panels(panels)
# Set up and solve the top boundary layer
top = Pohlhausen(top_panels,nu)
top.march()
# Set up and solve the bottom boundary layer
bottom = Pohlhausen(bottom_panels,nu)
bottom.march()
# plot flow with separation points
plot_flow(panels,alpha)
pyplot.scatter(top.x_sep, top.y_sep, s=100, c='r')
pyplot.scatter(bottom.x_sep, bottom.y_sep, s=100, c='g')
return top,bottom
In [11]:
top,bottom = solve_plot_boundary_layers(circle)
The red and green dots mark the separation point for the top and bottom boundary layer, respectively.
Separation occurs soon after the flow begins to decelerate. Physically, the boundary layer loses energy to friction as it travels over the front of the body (remember how large $C_F$ was?) and can not cope with the adverse pressure gradient on the back of the body.
In [12]:
def predict_jukowski_separation(t_c,alpha=0,N=128):
# set dx to gets the correct t/c
foil = make_jukowski(N,dx=t_c-0.019)
# find and print t/c
x0 = foil[N/2].xc
c = foil[0].xc-x0
t = 2.*numpy.max([p.yc for p in foil])
print "t/c = "+"%.3f"%(t/c)
# solve potential flow and boundary layer evolution
solve_gamma_kutta(foil,alpha)
top,bottom = solve_plot_boundary_layers(foil,alpha)
# print message
print ("Separation at x/c = "+"%.3f"%
((top.x_sep-x0)/c)+" from the leading edge")
In [13]:
predict_jukowski_separation(0.2,alpha)
We can make sure the behavoir above is correct by validating against the analytic solution for simple geometries. Here is a summary figure from Chapter 3 of Hoerner's Fluid-Dynamic Drag
There are two Jukowski examples: $t/c=0.15$ which separates at $x/c\approx0.49$ from the leading edge, and $t/c=0.17$, which separates at $x/c\approx0.39$.
In [14]:
predict_jukowski_separation(t_c=0.15)
The $t/c=0.15$ case matches very well with Hoerner's picture.
In [15]:
predict_jukowski_separation(t_c=0.17)
In [16]:
def predict_ellipse_separation(t_c,N=128,alpha=0):
ellipse = make_circle(N,t_c)
print "t/c = "+"%.3f"%(t_c)
# solve potential flow and boundary layer evolution
solve_gamma_kutta(ellipse,alpha)
top,bottom = solve_plot_boundary_layers(ellipse,alpha)
# print message
print ("Separation at x/c = "+"%.3f"%
((top.x_sep+1)/2.)+" from the leading edge")
In [17]:
predict_ellipse_separation(t_c=0.5)
In [18]:
predict_ellipse_separation(t_c=0.25)
In [19]:
predict_ellipse_separation(t_c=0.125)
So I get the feeling Hoerner has a typo... that's the first one I've found.
Therefore, the drag coefficient is
$$C_D = \frac{-F_x}{\frac 12 \rho U^2_\infty A} = \frac1w\oint_{\cal S} c_p s_y ds$$where $A$ is the 2D projected area of the body (the width) and $s_y = -n_x$.
Using the vortex panel method we can determine the potential flow solution for $c_p$, but what does $c_p$ look like in a real flow with separation?
I've sketched the results for the flow around a circular cylinder above. The measured pressure at the front of a body in a viscous fluid is fairly well predicted by potential flow.
However, the pressure coefficient completely deviates from the potential flow prediction near the point of separation. Indeed it remains essentially constant in the separated flow region.
In [20]:
# you code here
Ignore the line below - it just loads the style sheet.
In [21]:
from IPython.core.display import HTML
def css_styling():
styles = open('../styles/custom.css', 'r').read()
return HTML(styles)
css_styling()
Out[21]: