Flow Across a Cylinder

In Chapter 13 we studied flow over a cylinder.  We showed simulations and experimental data that verified boundary layer separation and the formation of a wake behind the cylinder.  The images that we discussed were steady-state images like Figure1.  What we did not discuss in detail is that this flow is unstable as the Reynolds number increases and if we look at what happens at Reynolds numbers in excess of about 100, we see a situation like that shown in Figure 2.  We call the situation, which looks like a series of swirls downstream from the cylinder, a Kármán vortex sheet.  In Figure 2, the vortex sheet was generated in the atmosphere as winds passed over the Juan Fernandez Islands off the coast of Chile.

                   

                                                                        Figure 1                                                                                                                           Figure 2

Notice that the vortices are strongest as they first leave the vicinity of the cylinder and dissipate the further away they travel.  The vortices are shed from opposite sides of the obstructing object so that there is a symmetry line that can be drawn down the axis of the sheet.  Since only one vortex is shed at a time, the periodic vortex generation about the axis line changes the pressure distribution and can cause an object to vibrate if it is not held rigidly.

When considering something simple, like a long circular cylinder, the frequency of vortex shedding is approximated by the empirical relationship

 

                 (1)

where:  ƒ = vortex shedding frequency.


This formula will generally hold true for Reynolds numbers in the laminar the range 250 < Re < 2 × 10 5 . The dimensionless parameter, St = ƒd/Vmax is called the Strouhal number.  It is named after a Czech physicist who first investigated the singing of telegraph wires in the wind in 1878.  Figure 3 below is taken from a Comsol simulation related to this effect.

                                                                                                                               Figure 3

a)  Run the module and plot the velocity field and pressure field at various times.  Notice how the vortices develop and then dissipate as they travel further

      from the cylinder.

b)  Plot the total and viscous forces about the cylinder as a function of time.  How do the viscous and form drag behave?

c)  Run the module for a series of Reynolds numbers.  From plots of total force on half the cylinder as a function of time, compare the frequency of vortex

      shedding with equation (1).  The Reynolds number for the flow can be accessed from the module as Red0.

d)  The module can also be used to look at the flight of a curveball, or at least what rotating the cylinder can do to the vortex sheet.  The rotation rate is

      controlled by specifying the rotation direction (–) clockwise and (+) represents counterclockwise.  The strength of the rotational velocity is represented

      by the frequency of rotation, w.  Run the module for clockwise versus counterclockwise rotation and discuss what happens to the vortex sheet.  How

      does the rotation affect the form and viscous drag?

e)  Run the module in the clockwise direction and vary the rotation rate?  Comment on what you observe?  Can you solve the problem for any rotation rate? 

      Investigate the same phenomena when the rotation is counterclockwise.

 

Table of Contents

Model Properties

 

PropertyValue
Model name Cylinder
Author J.L. Plawsky
Company Rensselaer Polytechnic Institute
Department Chemical and Biological Engineering
Reference 
URL 
Saved dateJun 27, 2010 4:01:15 PM
Creation date 
COMSOL versionCOMSOL 3.5.0.608

 

 

Application modes and modules used in this model:

 

Constants

 

NameExpressionDescription
rho01density
Umax0.2free stream velocity
ds0.1cylinder diameter
T273.15reference temperature
P0101325reference pressure
w 0.5 cylinder rotation rate
circ pi*w*ds

cyclinder rotation speed  (–) clockwise

(+) counterclockwise

 

 

Geometry

Number of geometries: 1

Geom1

Point mode

The geometry is formed very simply from a large rectangle that represents the fluid and a small circle, cut out of the rectangle that rrepresents the cylinder.  Notice that the cylinder is not located along the centerline, y = 0.2.  This asymmetry is necessary to see the shedding of the vortices.  Had the cylinder been centered, we would only be able to view the steady-state, symmetric solution as we do in the stationary cylinder example.

 

Curve x1 x2 x3 y1 y2 y3 weight1 weight2 weight3
1 0 0   0 0.41   1 1  
2 0 2.2   0 0   1 1  
3 0 2.2   0.41 0.41   1 1  
4 2.2 2.2   0 0.41   1 1  
5 0.15 0.15 0.2 0.2 0.15 0.15 1 0.7071 1
6 0.15 0.15 0.2 0.2 0.25 0.25 1 0.7071 1
7 0.2 0.25 0.25 0.15 0.15 0.2 1 0.7071 1
8 0.2 0.25 0.25 0.25 0.25 0.2 1 0.7071 1

 

 

Geom1

Space dimensions: 2D

Independent variables: x, y, z

Scalar Expressions

 

NameExpressionDescription
Red0ds*Umax/eta_nsReynolds number
uccw circ*sin(atan(y/x)) x-rotational velocity
vccw circ*cos(atan(y/x)) y-rotational velocity

 

 

Mesh

5.2.1. Mesh Statistics

The mesh was refined by solving as a stationary problem using the adaptive mesh routine.  Refinement was necessary to insure we capture the wake accurately.

 

Number of degrees of freedom61650
Number of mesh points6900
Number of elements13441
Triangular13441
Quadrilateral0
Number of boundary elements359
Number of vertex elements8
Minimum element quality0.515
Element area ratio0.003

 

Application Mode: Incompressible Navier-Stokes (ns)

Application mode type: Incompressible Navier-Stokes

Application mode name: ns

Scalar Variables

 

NameVariableValueDescription
visc_vel_factvisc_vel_fact_ns10Viscous velocity factor

 

 

Application Mode Properties

 

PropertyValue
Default element type Lagrange - P2 P1
Analysis typeTransient
Corner smoothingOff
FrameFrame (xy)
Weak constraintsOn
Constraint typeNon-ideal

 

 

Variables

Dependent variables: u, v, p, nxw, nyw

Shape functions: shlag(2,'lm1'), shlag(2,'lm2'), shlag(1,'lmp'), shlag(2,'u'), shlag(2,'v'), shlag(1,'p')

Interior boundaries not active

Boundary Settings

The boundary settings include a specified velocity at the inlet, an open boundary with specified pressure at the outlet, and either no-slip conditions on the sphere surface, or specified velocity that allows for clockwise or counterclockwise rotation.

 

Boundary12-34
TypeInletWallOutlet
x-velocity (u0) 000
Pressure (p0)00P0
Normal stress (f0)00(0)
Normal inflow velocity (U0in)Umax11
weakconstr000

 

 

Boundary 5 7 7 8
TypeWall Wall Wall Wall
x-velocity (uw) uccw uccw -uccw -uccw
y-velocity (vw) -vccw -vccw vccw vccw
weakconstr1 1 1 1

 

 

Subdomain Settings

The subdomain settings use property values for the fluid.  Temperature and pressure need to be specified to uniquely determine the properties.

 

Subdomain1
Shape functions (shape)shlag(2,'lm1') shlag(2,'lm2') shlag(1,'lmp') shlag(2,'u') shlag(2,'v') shlag(1,'p')
Integration order (gporder)4 4 2
Constraint order (cporder)2 2 1
Density (rho)rho(p[1/Pa],T[1/K])[kg/m^3] (Air)
Dynamic viscosity (eta)eta(T[1/K])[Pa*s] (Air)
Tuning parameter (delid)0.05
sdon0
cdon0

 

 

Subdomain initial value1
Pressure (p)P0

 

 

Materials/Coefficients Library

Air, 1 atm

 

ParameterValue
Heat capacity at constant pressure (C)Cp(T)
Dynamic viscosity (eta)eta(T)
Thermal conductivity (k)k(T)
Kinematic viscosity (nu0)nu0(T)
Density (rho)rho(p,T)

 

 

Functions

 

FunctionExpression
nu0(T)(-7.887E-12*T^2+4.427E-08*T+5.204E-06)/(1.013e5*28.8e-3/8.314/T)
Pr()eta/(k*Cp)
Cp(T)0.0769*T+1076.9
rho(p,T)p*28.8e-3/8.314/T
eta(T)-7.887E-12*T^2+4.427E-08*T+5.204E-06
k(T)10^(0.8616*log10(abs(T))-3.7142)

 

 

Solver Settings

Solve using a script: off

 

Analysis typeTransient
Auto select solverOn
SolverTime dependent
Solution formGeneral
SymmetricOff
Adaptive mesh refinementOff
Optimization/SensitivityOff
Plot while solving Off

 

 

 

Time Stepping

 

ParameterValue
Timesrange(0,0.3,15)
Relative tolerance0.001
Absolute tolerance0.00010
Times to store in outputSpecified times
Time steps taken by solverFree
Maximum BDF order5
Singular mass matrixMaybe
Consistent initialization of DAE systemsBackward Euler
Error estimation strategyInclude algebraic
Allow complex numbersOff