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.

- Model Properties
- Constants
- Geometry
- Geom1
- Materials/Coefficients Library
- Solver Settings
- Postprocessing
- Variables

Property | Value |

Model name | Cylinder |

Author | J.L. Plawsky |

Company | Rensselaer Polytechnic Institute |

Department | Chemical and Biological Engineering |

Reference | |

URL | |

Saved date | Jun 27, 2010 4:01:15 PM |

Creation date | |

COMSOL version | COMSOL 3.5.0.608 |

Application modes and modules used in this model:

- Geom1 (2D)
- Incompressible Navier-Stokes

Name | Expression | Description |

rho0 | 1 | density |

Umax | 0.2 | free stream velocity |

ds | 0.1 | cylinder diameter |

T | 273.15 | reference temperature |

P0 | 101325 | reference pressure |

w | 0.5 | cylinder rotation rate |

circ | pi*w*ds | cyclinder rotation speed (–) clockwise (+) counterclockwise |

Number of geometries: 1

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 |

Space dimensions: 2D

Independent variables: x, y, z

Name | Expression | Description |

Red0 | ds*Umax/eta_ns | Reynolds number |

uccw | circ*sin(atan(y/x)) | x-rotational velocity |

vccw | circ*cos(atan(y/x)) | y-rotational velocity |

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 freedom | 61650 |

Number of mesh points | 6900 |

Number of elements | 13441 |

Triangular | 13441 |

Quadrilateral | 0 |

Number of boundary elements | 359 |

Number of vertex elements | 8 |

Minimum element quality | 0.515 |

Element area ratio | 0.003 |

Application mode type: Incompressible Navier-Stokes

Application mode name: ns

Name | Variable | Value | Description |

visc_vel_fact | visc_vel_fact_ns | 10 | Viscous velocity factor |

Property | Value |

Default element type | Lagrange - P_{2} P_{1} |

Analysis type | Transient |

Corner smoothing | Off |

Frame | Frame (xy) |

Weak constraints | On |

Constraint type | Non-ideal |

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

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.

Boundary | 1 | 2-3 | 4 |

Type | Inlet | Wall | Outlet |

x-velocity (u0) | 0 | 0 | 0 |

Pressure (p0) | 0 | 0 | P0 |

Normal stress (f0) | 0 | 0 | (0) |

Normal inflow velocity (U0in) | Umax | 1 | 1 |

weakconstr | 0 | 0 | 0 |

Boundary |
5 |
7 |
7 |
8 |

Type | Wall | Wall | Wall | Wall |

x-velocity (uw) | uccw | uccw | -uccw | -uccw |

y-velocity (vw) | -vccw | -vccw | vccw | vccw |

weakconstr | 1 | 1 | 1 | 1 |

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

Subdomain | 1 |

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 |

sdon | 0 |

cdon | 0 |

Subdomain initial value | 1 |

Pressure (p) | P0 |

Parameter | Value |

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) |

Function | Expression |

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) |

Solve using a script: off

Analysis type | Transient |

Auto select solver | On |

Solver | Time dependent |

Solution form | General |

Symmetric | Off |

Adaptive mesh refinement | Off |

Optimization/Sensitivity | Off |

Plot while solving | Off |

Parameter | Value |

Times | range(0,0.3,15) |

Relative tolerance | 0.001 |

Absolute tolerance | 0.00010 |

Times to store in output | Specified times |

Time steps taken by solver | Free |

Maximum BDF order | 5 |

Singular mass matrix | Maybe |

Consistent initialization of DAE systems | Backward Euler |

Error estimation strategy | Include algebraic |

Allow complex numbers | Off |