RESEARCH ARTICLE | OCTOBER 27 2022

Dynamic mode decomposition analysis and fluid-mechanical

aspects of viscoelastic fluid flows past a cylinder in laminar
vortex shedding regime
F. Hamid ; C. Sasmal  ; R. P. Chhabra

Physics of Fluids 34, 103114 (2022)

https://doi.org/10.1063/5.0122103

CrossMark

 
View Export
Online Citation

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

Articles You May Be Interested In

Near-body vorticity dynamics of a square cylinder subjected to an inline pulsatile free stream flow
Physics of Fluids (September 2016)

Tensile properties of polypropylene (PP)/ virgin acrylonitrile butadiene rubber (NBRv)/Sago composites
under thermal aging
AIP Conference Proceedings (July 2021)
 Physics of Fluids ARTICLE scitation.org/journal/phf

Dynamic mode decomposition analysis

and fluid-mechanical aspects of viscoelastic
fluid flows past a cylinder in laminar vortex
shedding regime
Cite as: Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103
Submitted: 22 August 2022 . Accepted: 21 September 2022 .
Published Online: 27 October 2022

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

F. Hamid, C. Sasmal,a) and R. P. Chhabra

AFFILIATIONS
Soft Matter Engineering and Microfluidics Lab, Department of Chemical Engineering, Indian Institute of Technology Ropar,
Ropar, Punjab 140001, India

a)
Author to whom correspondence should be addressed: csasmal@iitrpr.ac.in

ABSTRACT
This study presents an extensive numerical investigation to understand the effect of fluid viscoelasticity on the flow dynamics past a
stationary cylinder in the laminar vortex shedding regime. The governing equations, namely, mass, momentum, and Oldroyd-B viscoelastic
constitutive equations, have been solved at a fixed value of the Reynolds number of 100 and over a range of values of the Weissenberg num-
ber as 0  Wi  2 and polymer viscosity ratio as 0:5  b  0:85. In particular, for the first time, this study presents a detailed analysis of
how the fluid viscoelasticity influences the coherent flow structures in this benchmark problem using the dynamic mode decomposition
(DMD) technique, which is considered to be one of the widely used reduced order modeling techniques in the domain of fluid mechanics.
We show that this technique can successfully identify the low-rank fluid structures in terms of the spatiotemporal modes from the time-
resolved vorticity field snapshots and capture the essential flow features by very few modes. Furthermore, we observe a significant difference
in the amplitude and frequency associated with these modes for Newtonian and viscoelastic fluids otherwise under the same conditions.
This, in turn, explains the differences seen in the flow dynamics between the two types of fluids in an unambiguous way, such as why the
fluid viscoelasticity suppresses the vortex shedding phenomenon and decreases the energy associated with the velocity fluctuations in visco-
elastic fluids than that in Newtonian fluids. However, before performing the DMD analysis, we also present a detailed discussion on the vari-
ous fluid-mechanical aspects of this flow system, such as streamline patterns, vorticity fields, drag and lift forces acting on the cylinder, etc.
This will ultimately set a reference platform for delineating the importance of the DMD analysis to get further insight into flow physics.
Published under an exclusive license by AIP Publishing. https://doi.org/10.1063/5.0122103

I. INTRODUCTION shedding phenomenon, etc., is now available for this geometry. Such
The flow past a bluff body is an area of intense research from the information is of utmost importance not only in the context of engi-
past several decades in the domain of fluid mechanics.1 It facilitates a neering applications but also in basic scientific exploration. In particu-
better and thorough understanding of various complex flow physics lar, the knowledge regarding the wake formation and the subsequent
phenomena occurring across multiple time and length scales in diverse vortex shedding phenomenon, which result in the formation of the
practical applications.2 Among several bluff bodies considered in these von Karman vortex street above a critical Reynolds number (ratio of
studies, the circular cylinder is probably the most used one owing to the inertia to that of the viscous forces) has significant practical impli-
its ubiquity in nature, simple configuration, and the absence of any cations in various industrial settings such as in the design of aircraft,
geometric singularities.3 Extensive studies comprising theoretical, offshore rigs, cooling towers, distillation columns, heat exchangers,
numerical, and experimental analyses have been conducted for this and so on.5,6
geometry in the literature.4 Suffice it to mention here that a significant Many factors can influence the flow physics past a cylinder, such
amount of information, such as forces acting on the cylinder, transi- as cylinder size, flow conditions, externally applied forces (magnetic or
tion to different flow regimes, wake formation mechanism, vortex electrical forces), etc. Additionally, the rheological behavior of the

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-1

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

working fluid, i.e., whether it is simple Newtonian or complex non- various flow variables are now becoming available both from the
Newtonian, can also significantly alter the flow phenomena. To date, numerical simulations and experimental measurements.24 As a result,
the number of studies considering the non-Newtonian behavior of big data has become a reality in the domain of fluid mechanics.25
fluid is far less than that available for Newtonian fluids. However, a However, much of these data are redundant as the fluid flow, rather
wide range of fluids (such as polymer solutions and melts, foams, than random, consists of many patterns and/or coherent flow structures
emulsions, paints, cosmetics, biological fluids like blood, saliva, etc.) representing its important characteristics.26 These are vital in calculating
frequently encountered in various industrial and natural processes are the quantities of engineering importance such as heat and mass transfer
non-Newtonian in nature.7,8 Among various non-Newtonian charac- rate.27 Currently, proper orthogonal decomposition (POD) and
teristics, viscoelasticity (a manifestation of the combined viscous and dynamic mode decomposition (DMD) are the most widely used mathe-
elastic properties) to various extents is exhibited by most complex matical techniques that aid in distilling these important spatial features
fluids.9 Adding a small amount of solid polymer into a Newtonian sol- of a flow field in terms of the so-called “modes.” However, DMD cap-
vent like water induces a significant extent of viscoelastic behavior in tures both the temporal and spatial coherent structures of a flow field,
the resulting polymer solution.10 This property of such complex fluids whereas POD captures only the spatial coherent structures. As a result,
has attracted much attention due to its practical implication in the DMD is often preferred over POD.28,29 Since its inception in the field of
drag reduction phenomenon.11–13 fluid mechanics, DMD has been successfully utilized to analyze various

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

However, to date, very few studies are available (in comparison flow problems, such as the flow past solid boundaries, wakes, jet flows,
with that available for Newtonian fluids) in the literature for the flow of mixing layers, cavitating flow, etc.24,30,31 In his pioneering work,
viscoelastic fluids past a cylinder in the finite Reynolds number regime, Schmid32 demonstrated how the DMD technique could be utilized in
particularly in the periodic vortex shedding regime. Among a few, analyzing the flow field data of flame and water jets captured by a high-
Gadd12 performed an experimental investigation to delineate the effect speed camera. He explained the significance of different DMD modes
of the fluid viscoelasticity on the formation of the vortex street in flows and how they capture the underlying flow structures in detail.
past a cylinder at a fixed Reynolds number of 240. His study showed Thereafter, this technique has been widely used to get better insights
suppression in the vortex shedding phenomenon in viscoelastic fluids, into different flow systems, and one of them is the flow past a cylinder.
which again depended on the type of polymer additives used to induce In particular, the vortex shedding phenomenon past a cylinder has
the fluid viscoelasticity. This observation was further verified by a num- been widely used as a model phenomenon to test the ability and effi-
ber of other independent experimental studies in the literature.14–17 On ciency of the DMD technique and its different variants.33–36 For
the other hand, Oliveira18 numerically studied the laminar vortex shed- instance, Zhang et al.27 conducted a study to extract the coherent struc-
ding phenomenon past a circular cylinder using the FENE-CR (Finite tures for the flow past a single and two side-by-side cylinders using both
Extensible Non-linear Elastic-Chilcott and Rallison) viscoelastic fluid DMD and POD techniques. In the case of a single cylinder, it was
model at a fixed Reynolds number of 100. Similar to the experiments, observed that only one temporal DMD mode and two POD modes
he also observed attenuation of the vortex shedding frequency. (which are again contaminated with undesirable structures) are required
Furthermore, he also found that the maximum extensibility of the poly- to effectively capture the Karman vortex street. The same trend was
mer molecules significantly influences the vortex shedding phenomenon observed for the case of two side-by-side cylinders. Therefore, it suggests
along with the dimensionless Weissenberg number. Similar findings that the DMD technique efficiently captures the temporal and spatial
were also obtained in another numerical study by Sahin and Owens.19 flow structures by a lesser number of modes in comparison with the
They also conducted a linear stability analysis and found that the critical POD technique, for which the modes are again contaminated with
value of the Reynolds number for the onset of the vortex shedding phe- undesirable structures. Zhang et al.37 used the lattice Boltzmann method
nomenon increases with the maximum extensibility of the polymer (LBM) in combination with the DMD technique to perform the stability
molecules. In a very recent study, Minaeian et al.20 proposed correla- analysis of flow past a cylinder at a fixed Reynolds number of 50. Their
tions for these critical Reynolds numbers for the PTT (Phan–Thien– analysis was performed on two data sets: one obtained from the tran-
Tanner) viscoelastic fluid model. The three-dimensional numerical sient region (where the flow transits from steady to unsteady) and the
simulations of viscoelastic fluid flow past a cylinder were carried out by second from the periodic vortex shedding regime. The DMD analysis in
Richter et al.21 using the FENE-P (Finitely Extensible Non-linear the latter case generates steady modes, whereas, for the transient period,
Elastic-Peterlin) viscoelastic fluid model at two values of the Reynolds modes are associated with a growth rate, i.e., they are unstable as the
number, namely, 100 and 300. They also observed a stabilizing effect of flow is evolving. In a more recent study, Ping et al.26 performed the
the fluid viscoelasticity on the flow transition, particularly at high values DMD analysis on the flow past an oscillating cylinder submerged in a
of the polymer extensibility parameter. The corresponding study in Newtonian fluid. In addition to extracting the coherent modes in this
the steady flow regime at low Reynolds numbers using the Oldroyd-B flow system, another main aim of this study was to determine the struc-
and FENE-P viscoelastic fluid models was carried out by Sasmal et al.22 tures that have an impact on the hysteresis of the lift and drag forces act-
They found that the separation of boundary layers in flows past a cyl- ing on the cylinder, which was observed experimentally by Bishop and
inder was suppressed entirely in a constant viscosity Oldroyd-B visco- Hassan.38 In this work, it was successfully demonstrated that the DMD
elastic fluid, whereas it was accelerated in FENE-P viscoelastic fluid. technique could separate such essential features of complex flow phe-
Norouzi et al.23 performed the study using the Giesekus viscoelastic nomena responsible for the hysteresis effect, which are otherwise almost
fluid model and found an increase in the vortex shedding frequency impossible to assess by any other analysis.
with the mobility factor associated with this model. Therefore, based on the aforementioned discussion, two things
On the other hand, with the rapid advances in computational are clear. First, the fluid viscoelasticity can significantly influence the
power and measurement techniques, excessive high-fidelity data of vortex shedding phenomenon in the case of a cylinder, and second,

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-2

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

the DMD technique has the potential to capture the coherent struc- presence of polymer molecules (assumed to be linear elastic dumb-
tures of the flow field. Although experimental and numerical studies bells that can be stretched infinitely) in a viscous Newtonian sol-
on the flow of viscoelastic fluids past a cylinder were carried out in the vent. Therefore, the total deviatoric stress tensor in a viscoelastic
laminar vortex shedding regime; however, there is no systematic and fluid is a linear combination of the viscous term ss (the solvent
detailed analysis of the coherent flow structures present for this flow contribution) and the viscoelastic extra-stress term sp (the polymer
system, particularly using the DMD technique. The application of the molecules contribution), i.e.,
DMD technique is so far limited to the flows of simple Newtonian flu- s ¼ ss þ sp : (3)
ids. Therefore, the main aim of the present study is to provide a deep
understanding of the viscoelastic fluid flow phenomena past a cylinder The solvent contribution to the stress is given by
in the laminar vortex shedding regime by analyzing its coherent flow ss ¼ 2gs D; (4)
structures using the DMD technique. We believe this analysis will
where D is the rate of deformation tensor given by 12 ½ðruÞ þ ðruÞT 
unravel more underlying physics and provide more concrete explana-
and gs is the viscosity of the Newtonian solvent. The viscoelastic extra
tions for the phenomena in viscoelastic fluid flow past a cylinder, such
stress tensor due to the polymer molecules is calculated using the
as the suppression of vortex shedding. However, before performing
Oldroyd-B model as follows:
that, we first present and discuss the results on the fluid-mechanical
gp

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

aspects such as streamlines, vorticity, drag and lift forces acting on the sp ¼ ðC  dÞ; (5)
cylinder, etc. This will set a better platform for the importance of the k
DMD analysis to get further insights into the present flow physics. r
C þ k C ¼ d: (6)
The rest of the present paper is organized as follows. We first describe
the problem and write all the governing equations along with the asso-
In Eqs. (5) and (6), C is the non-dimensional conformation tensor, gp
ciated boundary conditions in Sec. II. The details of the numerical
is the zero-shear viscosity of the polymer, k is the polymer relaxation
solution methodology and the dynamic mode decomposition technique r
are presented in Sec. III. A detailed discussion of the present results, time, d is the identity tensor, C is the Oldroyd upper-convected deriva-
both on fluid-mechanical aspects and DMD analysis, is provided in tive of C given as follows:
Sec. IV. Finally, Sec. V summarizes the findings of the present study. r @C
C¼ þ u  rC  C  ru  ð ruT Þ  C: (7)
II. PROBLEM STATEMENT AND GOVERNING @t
EQUATIONS These governing and constitutive equations have been non-
The problem considered herein is the study of the effect of fluid dimensionalized using the following scaling parameters: d for all
viscoelasticity on the coherent flow structures in flows past an infi- lengths, Ud1 for time, U1 for velocity, q U1
2
for pressure, g0 dU1 for stress
nitely long circular cylinder of diameter d in the laminar vortex shed- tensor. Here g0 is the zero-shear rate viscosity of the polymer solution
ding regime, as shown schematically in Fig. 1(a). The incompressible given by the sum of the solvent viscosity ðgs Þ and the contribution of
viscoelastic fluid approaches the cylinder surface with a uniform veloc- the polymer ðgp Þ, i.e., g0 ¼ gs þ gp . The non-dimensional forms of
ity of U1 . Furthermore, the flow is assumed to be two-dimensional the mass, momentum, and Oldroyd-B viscoelastic constitutive equa-
and unconfined. The latter situation is realized here by placing the cyl- tions are as follows:
inder in a fictitious rectangular domain of fluid, as schematically
shown in Fig. 1(b). Due to the two-dimensional nature of the flow and r  u ¼ 0; (8)
the infinitely long size of the cylinder, it is safe to assume that uz ¼ 0 @u 
b ð1  bÞ 1  
and @ðÞ þ u  r u ¼ r p þ r  ss þ r  sp ; (9)
@z ¼ 0. The dimension of the rectangular domain in the stream- @t  Re Wi Re
wise direction is Lu þ Ld , where Lu ¼ 20d and Ld ¼ 60d are the
1b
upstream and downstream lengths of the domain, respectively. A total sp ¼ ðC  dÞ; (10)
height of Hd ¼ 40d is used in the spanwise direction with the cylinder Wi
r
placed in the middle. A typical grid structure used in this study is pre- C þ Wi C  ¼ d: (11)
sented in Fig. 1(c). The present flow is governed by the continuity and
momentum equations written in their dimensional forms as below: Here, Re ¼ q Ug1 d is the Reynolds number defined as the ratio of the
0
Continuity equation: inertial to that of the viscous forces, Wi ¼ k Ud 1 is the Weissenberg
r  u ¼ 0: (1) number defined as the ratio of the elastic to that of the viscous forces,
gs
and b ¼ g þg is the viscosity ratio. In the limit b ! 0, the fluid
s p
Momentum equation: behaves like a polymer melt, whereas it behaves like a Newtonian fluid
  in the limit of b ! 1. Finally, the following boundary conditions are
@u
q þ u  ru ¼ rp þ r  s; (2) prescribed to complete the problem formulation.
@t
Inlet: The fluid enters the computational domain with a uniform
where q is the fluid density fluid, u is the velocity vector, p is the stream-wise velocity in the x-direction, i.e., ux ¼ U1 and the span-
pressure, and s is the deviatoric stress tensor. The Oldroyd-B con- wise velocity component is set to zero, i.e., uy ¼ 0. The pressure gradi-
stitutive equation is used to model the fluid viscoelasticity.39 ent and the extra stresses due to polymeric contribution are set to zero
According to this model, fluid viscoelasticity arises due to the at this boundary.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-3

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 1. (a) Schematic of the present problem, (b) computational domain, and (c) a typical grid used in the present study.

ð
Outlet: For all variables, the Neumann type boundary condition  
CL ¼ K FL ¼ K ðp d þ ss þ sp Þ  ns y dS; (13)
is used except for pressure which is assigned a zero value at this s
boundary.
Cylinder: At the cylinder surface, the standard no-slip boundary where FD and FL are the total drag and lift forces per unit length of the
condition, i.e., ux ¼ uy ¼ 0 is applied. The polymeric stresses are cylinder respectively, K is equal to q U22 d ; ns is the outward unit nor-
1
extrapolated linearly onto the surface, and the pressure gradient is set mal vector drawn on the cylinder surface, and S is the cylinder surface
to zero. area.
Top and bottom planes: At both these planes, symmetry bound- The Strouhal number (St) is another important dimensionless
ary condition is applied, which is given by @u
@y ¼ 0 and uy ¼ 0.
x number that is extensively used to analyze the flow physics in the
On solving these governing equations using the above- unsteady regime. It provides the frequency of the vortex shedding and
mentioned boundary conditions, the results are obtained in terms of is calculated by the expression Uf 1d , where f is the frequency whose value
velocity (u), pressure (p), and stress tensor ðsÞ. These are further post- can be obtained from the temporal evolution of the lift coefficient in
processed to calculate quantities like the hydrodynamic forces acting the periodic flow regime.
on the cylinder. It is a combination of the stream-wise drag and cross-
III. RESEARCH METHODOLOGY
stream lift components. In this study, these are discussed in terms of
the drag ðCD Þ and lift ðCL Þ coefficients, which are calculated using the A. Numerical methodology
following relations: While the governing equations described in Sec. II, namely, mass
ð and momentum equations, have been solved using the finite volume
  method (FVM) based open-source computational fluid dynamics
CD ¼ K FD ¼ K ðp d þ ss þ sp Þ  ns x dS; (12)
s (CFD) code OpenFOAM, the rheoFoam solver available in the

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-4

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

recently developed RheoTool40 has been used to solve the Oldroyd-B TABLE I. Details of different grids used in the grid independence study.
viscoelastic constitutive equation. Various numerical schemes available
in OpenFOAM have been used to discretize the different terms of the Grid Nxa Ny b DSc105 Nd
governing equations. For instance, all the advective terms have been
G1 355 180 5.14 74 100
discretized using the CUBISTA (Convergent and Universally Bounded
Interpolation Scheme for Treatment of Advection) scheme.41 The dif- G2 460 250 2.50 140 600
fusion term in the momentum equation was discretized using the G3 545 295 1.36 200 150
Gauss linear orthogonal interpolation scheme, which has an accuracy a
Number of grid elements in the x-direction.
of the order two. All time derivative terms have been discretized using b
Number of grid elements in the y-direction.
the Euler scheme. While the linear systems of the pressure and velocity c
Minimum cell area.
fields were solved using the Preconditioned Conjugate Solver (PCG) d
Total number of cells in the whole computational domain.
with DIC (Diagonal-based Incomplete Cholesky) preconditioner, the
stress fields were solved using the Preconditioned Bi-conjugate Gradient TABLE II. Values of time-averaged drag coefficient obtained with different grid
Solver (PBiCG) solver with DILU (Diagonal-based Incomplete LU) densities.
preconditioner. The pressure–velocity coupling was accomplished

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

using the SIMPLE algorithm. Here, the preconditioners are employed Grid CD;avg
to increase the efficiency and robustness of the solvers. Furthermore,
the simulations were stabilized using the log-conformation tensor G1 1.704
approach.42,43 The relative tolerance level for all the fields to be calcu- G2 1.695
lated was set to 1010 . Another important parameter that determines G3 1.690
the stability of the numerical solution is the Courant number, defined
as the ratio ofthe time step size ðD tÞ to the characteristic convective
Dx drag coefficient for different grids is presented in Table II. It can be
time scale u , i.e.,
seen that the results obtained with grids G2 (140 600 cells) and G3
uDt (200 150 cells) are almost indistinguishable from each other with a
Co ¼ :
Dx percentage difference of 0.29. Therefore, the grid G2 with 140 600 hex-
ahedral cells was selected in our study to carry out all the simulations.
To ensure the stability and accuracy of the solution (explicit numerical
Furthermore, a time step size of 0.001 Ud1 was selected based on a
schemes), the following condition must be satisfied:
detailed time independent study, likewise the domain and grid-
Co  1: independent studies.
In this study, we have maintained Co < 0.5 for all our simulations.
In addition to the interpolation schemes mentioned above, the B. Dynamic mode decomposition (DMD) analysis
size of the computational domain, grid density, and time step size also The data-driven dynamic mode decomposition (DMD) algo-
play a significant role in solving any CFD problem with the desired rithm44 extracts the coherent structures underlying any flow phenome-
precision and accuracy. Therefore, we have conducted extensive stud- non from the spatiotemporal data sequence of the flow fields obtained
ies to determine the best combination of grid fineness, domain, and either from numerical simulations or experiments. Essentially, it is an
time step sizes for efficiently utilizing computational resources without eigenvector/eigenvalue problem and involves the decomposition of a
any loss of precision. As already mentioned in Sec. II, the present com- linear model that maps two temporally consecutive datasets.32 The
putational domain is a rectangular-shaped domain with dimensions resulting eigenvectors, known as the Ritz vectors, physically represent
defined in terms of the upstream length Lu, downstream length Ld, the Eigen flow fields and are referred to as the DMD modes. The asso-
and height Hd. These dimensions have been chosen as Lu ¼ 20d, Ld ciated eigenvalues, known as the Ritz values, give a characteristic fre-
¼ 60d, and Hd ¼ 40d by performing the standard domain indepen- quency and a growth/decay rate.
dence study. These values obtained in this study are also in line with In the present study, time-resolved spatial vorticity data from
that used in earlier studies for simulating the same flow for Newtonian numerical simulations are acquired and converted into a matrix form
and viscoelastic fluids. After fixing the domain size, an extensive grid such that a column in the data matrix represents each state as follows:
independence study was performed to choose an optimum grid den-  
sity by running simulations on three different hexahedral grids, X1M ¼ xðt1 Þ xðt2 Þ    xðtM Þ 2 RNM ;
namely, G1, G2, and G3, with a different number of elements on the
cylinder surface and in the whole computational domain. A relatively where the subscript and superscript of X represent the initial and final
fine grid is used in the vicinity of the cylinder surface to capture the time steps, respectively. Using the Krylov technique, DMD assumes a
steep gradients of flow variables such as velocity or stress, as schemati- linear mapping A 2 RNN between the two consecutive snapshots,
cally shown in Fig. 1(c). The details of different grids used in the grid which is constant over the data sequence,32 i.e.,
independence study are shown in Table I. The simulations for differ- xjþ1 ¼ Axj : (14)
ent grid densities were run at Re ¼ 100 and at the highest value of the
Weissenberg number of Wi ¼ 2 and the lowest value of the polymer Using vector x1 and propagator A, the application of the Krylov
viscosity ratio of b ¼ 0:5 considered in this study. The time-averaged sequence results in the following:

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-5

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

X1M ¼ fx1 ; Ax1 ; A2 x1 ; …; AM1 x1 g: (15) /j ¼ Uvj : (24)

In this sequence, the snapshot vectors become linearly dependent on The eigenvalues kj provide the information about temporal dynamics
the preceding vectors, and the final snapshot can, therefore, be of the modes and are usually converted by logarithmic mapping into
expressed as a linear combination of the previous ones within the limit the following form:
of a large number of snapshots as follows:
uj ¼ logðkj Þ=Dt ¼ rj þ ixj : (25)
xM ¼ a1 x1 þ a2 x2 þ a3 x3 þ    þ aM1 xM1 þ r; (16)
Here, the real part rj represents the growth rate, and xj represents the
where r is the residual vector. Assembling the coefficients in the vector angular frequency of jth mode. The flow field can be, therefore, con-
form, we get structed as a linear combination of the modes as follows:
xM ¼ a X1M1 þ r; (17) X
r
xðtÞ  /j expðxj tÞbj ¼ U expðX tÞb; (26)
T
where a ¼ fa1 ; a2 ; a3 ; …; aM1 g. This system can be represented as j¼1

X2M ¼ A X1M1 : (18) where U is a matrix of eigenvectors /j ; X ¼ diagðxj Þ, and b (a vector

of bj) gives the coefficient of each mode at t ¼ 0 (time corresponding

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

Here, the ordering is done in a fashion such that the corresponding to first snapshot) and is obtained from b ¼ U† x1 . This vector quanti-
columns in X1M1 and X2M are shifted by one time step. To approxi- fies the energy contribution made by each mode to the overall data set,
mate the unknown A matrix in Eq. (18), a companion matrix S 2
which can also be alternatively obtained from the norm of each mode
RðM1ÞðM1Þ is built based on Eq. (17) as follows:
as jj/j jj.
X2M ¼ A X1M1  X1M1 S; (19)
where S is given by IV. RESULTS AND DISCUSSION
2 3 A. Validation
0 0  0 a1
6 7 Extensive validation of the present numerical code was first
61 0  0 a2 7
6 7 carried out before investigating the present problem. Studies carried
6 .. .. .. 7
S¼6 7: (20) out by Bharti et al.46 and Soares et al.47 were chosen for the valida-
60 1 . . . 7
6. .. .. 7 tion of simple Newtonian fluids flow past a stationary cylinder. A
6 .. . . 0 aM2 7
4 5 comparison in terms of the frictional and pressure drag coefficients
0  0 1 aM1 is shown in Table III. The results are seen to be in good agreement
with each other, with a maximum error of less than 4%.
This matrix shifts the snapshots 1 through M  1 and approximates Furthermore, the length of the wake ðLw Þ behind the cylinder is also
Mth snapshot as a linear combination of preceding ones. To compute calculated, and the results are again found to be on good terms with
the eigenvectors and eigenvalues of the S matrix, the singular value the literature values with a maximum deviation of 3.6%, as shown
decomposition (SVD) method is preferred for its robustness. The pro- in Table IV. Additional validation for the flows of Oldroyd-B visco-
cedure for this method is shown as follows: elastic fluids in a lid-driven cavity has been carried out by compar-
ing the present results with that of Fattal and Kupferman.45 Both
X1M1 ¼ URW T ; (21)
the horizontal and vertical velocity profiles at different vertical and
where U is a unitary matrix and R is a diagonal matrix of singular val-
ues arranged in descending order of magnitude. The rank of the singu- TABLE III. Comparison of the present values of the frictional ðCDF Þ and pressure
lar matrix is usually very low, and hence eigenvectors and eigenvalues drag ðCDP Þ coefficients with Bharti et al.46 and Soares et al.47
of low dimensional representation of the S matrix are obtained.
Substituting Eq. (21) in Eq. (19) yields Present Bharti et al.46 Soares et al.47
~S ¼ U T X M WR1 ;
2 (22) Re CDFa CDPb CDF CDP CDP CDF
rr
where ~S ¼ U T AU 2 R is the projection of A onto low dimen- 5 3.7752 2.1470 3.8387 2.1345  
sional space of proper orthogonal decomposition (POD) modes of U. 10 2.4634 1.5709 2.5078 1.5723  
The eigenvectors vj and the corresponding eigenvalues kj of ~S matrix 15 1.9248 1.3460 1.9587 1.3497  
approximate some of the Ritz vectors and corresponding Ritz values of 20 1.6143 1.2200 1.6421 1.2244 1.6422 1.1900
A,32 which ultimately can be used to characterize the flow field. The
25 1.4066 1.1369 1.4300 1.1412  
DMD modes can now be computed from the eigen-decomposition of
~S as follows: 30 1.2555 1.0768 1.2754 1.0807  
35 1.1395 1.0310 1.1569 1.0341  
~Svj ¼ kj vj ; (23) 40 1.0466 0.9932 1.0633 0.9977 1.0400 0.9600
r
where vj 2 R are the eigenvectors for corresponding kj eigenvalues. a
Friction drag coefficient.
b
The dynamic modes /j are, therefore, given by Pressure drag coefficient.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-6

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

TABLE IV. Comparison of the present values of the wake length (Lw) with Bharti
et al.46 and Fornberg.48

Lw

Re Present Bharti et al.46 Fornberg48

5   

10 0.2400 0.2487 
15 0.5789 0.5821 
20 0.9150 0.9164 0.9100
25 1.2450 1.2486 
30 1.5586 1.5796 
35 1.8450 1.9081  FIG. 3. Visual description of the time instances for the streamlines acquisition (red
circles). While time instances (a) (t ¼ 0) and (c) (t ¼ T/2) correspond to the zero lift,
40 2.1450 2.2252 2.2400 (b) (t ¼ T/4) and (d) (t ¼ 3T/4) correspond to the minimum and maximum lifts,

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

respectively.

horizontal planes of the cavity exhibit nearly a perfect match, as can B. Analysis of fluid-mechanical aspects
be seen in Fig. 2.
After performing rigorous validation studies of the present In the present study, the flow is in the laminar vortex shedding
numerical solver, we now present and discuss our new results in regime; hence, the streamlines, vorticity, and polymer conformation
detail. First, the fluid-mechanical aspects, such as streamlines and fields have been generated at selected time instances of a shedding
vorticity fields, hydrodynamic forces (drag and lift) acting on the cycle of a total period of T, as shown in Fig. 3. These time instances
cylinder, polymer stretching, etc., are briefly presented and exam- correspond to the zero, maximum, and minimum lift forces acting on
ined. This initial discussion on the fluid-mechanical aspects will the cylinder. For two-dimensional flows, the velocity field is often ana-
shed light on how fluid viscoelasticity can influence the flow physics lyzed with the help of temporal and spatial variations of streamline
in viscoelastic fluids compared to Newtonian fluids. Finally, the patterns. Figure 4 represents the streamlines along with the velocity
DMD technique will be utilized to understand and examine the dif- magnitude plots both for Newtonian and viscoelastic fluids with a
ference between Newtonian and viscoelastic fluids in more detail by fixed Wi ¼ 2 and two values of the polymer viscosity ratio, namely,
analyzing the coherent flow structures obtained in these two fluids. b ¼ 0:85 and 0.5. For Newtonian fluids, at the time instance a (t ¼ 0),
In the present study, the results are obtained at a fixed value of the one can see that the flow is asymmetric with a small anti-clockwise
Reynolds number of Re ¼ 100 and over a range of the Weissenberg wake formed at the right bottom corner of the cylinder and a large
number as 0  Wi  2 and polymer viscosity ratio as clockwise wake formed almost one cylinder diameter away down-
0:1  b  0:5. These selected ranges of values of the governing stream of the cylinder. The velocity magnitude is seen to be always
parameters are in line with that used in earlier several stud- high at the top and bottom sides of the cylinder, whereas it is less (or
ies.18,19,21,22,49 Also, note that at these combinations of parameters, almost zero) just downstream of the cylinder where the wake is
the elastic instability is not expected to appear, which occurs when formed. At the next time instance b ðt ¼ T=4Þ of the cycle, the right
the effect of inertial forces is negligible in a flow field, i.e., in a creep- bottom corner vortex grows in size, whereas the large downstream vor-
ing flow.50,51 However, the shear flow instability will be present in tex disappears. Simultaneously, a small clockwise vortex is initiated at
the present system, which ultimately leads to the generation of the the top right corner of the cylinder. The anti-clockwise vortex detaches
von-Karman vortex street in both Newtonian and viscoelastic fluids. from the cylinder surface in the next time instance c ðt ¼ T=2Þ and

FIG. 2. Comparison of non-dimensional x-component (a) and y-component (b) velocities between the present results and Fatal and Kupferman45 at Wi ¼ 1.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-7

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 4. Streamline patterns for Newtonian and viscoelastic fluids with a fixed Weissenberg number of 2 and two values of the viscosity ratio, namely, 0.85 and 0.5, at different
time instances of a periodic cycle, (a) t ¼ 0, (b) t ¼ T/4, (c) t ¼ T/2, (d) t ¼ 3T/4, as schematically shown in Fig. 3.

the clockwise vortex formed at the top right corner subsequently grows top and bottom sides of the cylinder is always higher in viscoelastic
in size. Streamlines corresponding to the time instance d ðt ¼ 3T=4Þ fluids than in Newtonian fluids. This tendency gradually increases with
show the maximum size of the clockwise vortex, which detaches in the the Weissenberg number and with the decreasing polymer viscosity
repeated time instance of a (t ¼ T) where the periodic cycle starts again. ratio. This is because the stretching of polymer molecules at these two
The mechanisms of vortex formation, detachment, and subsequent sides of the cylinder (due to the presence of a high shearing flow field)
shedding downstream of the cylinder are well investigated by many creates a barrier for the fluid to move in the spanwise direction. This
earlier studies for Newtonian fluids with the help of both simulations results in an increased velocity magnitude at these two places and
and experiments. On the other hand, for viscoelastic fluids, the stream- increased vortex size downstream of the cylinder.
line features are similar to the Newtonian fluids at low values of the The corresponding vorticity fields are shown in Fig. 5 at the same
Weissenberg number and high values of the viscosity ratio (results are time instances of a periodic cycle, as depicted in Fig. 3. Regardless of
not shown here for the sake of brevity). However, as the Weissenberg the fluid type, the vorticity fields show the well-known von Karman
number increases and the viscosity ratio decreases, the vortex dynamics vortex street structure in all cases. This is expected as the results are
drastically changes. First and foremost, under otherwise identical con- obtained at Re ¼ 100, which lies in the laminar vortex shedding
ditions, the vortices in viscoelastic fluids are bigger in size than those regime. On close observation with the magnified view near the cylin-
seen for Newtonian fluids at the same instant of time. Furthermore, der surface, it can be seen that the vorticity is produced upstream of
they are now more stretched (as a result, the aspect ratio decreases) and the cylinder, which gradually diffuses into the boundary layer. Later, at
attached to the cylinder surface. Due to the increase in the vortex size, the separation point, this vorticity is convected into the associated
the distance of the saddle point (where two streamlines meet with each shear layers [positive shear layer (PS) shown in red color and negative
other) from the rear stagnation point of the cylinder also increases in shear layer (NS) shown in blue color], and consequently, vortex shed-
viscoelastic fluids. For instance, let us take the example of the time ding occurs due to the shear layer instability. It occurs both for
instance a. At this time, the anti-clockwise vortex formed at the right Newtonian and viscoelastic fluids irrespective of the values of the
bottom corner of the cylinder is seen to be smaller for Newtonian fluids Weissenberg number and viscosity ratio. However, the interaction
than that seen for viscoelastic fluids with Wi ¼ 2 and b ¼ 0:5. On the between the shear layers and base region vorticities, denoted by PBRV
other hand, the same trend is also observed for the clockwise vortex (Positive Base Region Vorticity) and NBRV (Negative Base Region
formed downstream of the cylinder. Moreover, this downstream vortex Vorticity), along with the growth and decay of the base region vortic-
is attached to the cylinder surface for viscoelastic fluids, whereas it is ities are different for the two fluids. These play a vital role in the vortex
already detached and away from the rear stagnation point of the cylin- shedding mechanism and ultimately result in the change of the shape
der in the case of Newtonian fluids. It suggests that the vortex detach- of the vortices being shed in the two fluids. For instance, at time
ment or the vortex shedding frequency should be suppressed in instance a, the NBRV at the right bottom corner appears to be at the
viscoelastic fluids compared to Newtonian fluids. This will be further maximum strength, whereas a minor PBRV appears at the top corner
examined in the subsequent discussion on the results of drag and lift for a Newtonian fluid. The NBRV then undergoes a decay in the next
coefficients and vorticity patterns. Our observations are in line with time instance b, whereas the PBRV grows and attains the maximum
that seen in earlier experiments.11,12,14,15 The velocity magnitude at the strength in the time instance c. The PBRV now undergoes a decay in

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-8

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 5. Vorticity fields for Newtonian and viscoelastic fluids with a fixed Weissenberg number of 2 and two values of the viscosity ratio, namely, 0.85 and 0.5, at different time
instances of a periodic cycle, (a) t ¼ 0, (b) t ¼ T/4, (c) t ¼ T/2, (d) t ¼ 3T/4, as schematically shown in Fig. 3. A magnified view of the vorticity fields near the cylinder surface is
shown in the inset figures.

the time instance d, and the NBRV (which was initiated in the last The stretching of polymer molecules in the present viscoelastic
time instance) is again seen to grow. In this way, the cycle repeats fluids has been quantified using the non-dimensional conformation
where the interactions of the PBRV with the NS layer and the NBRV tensor ðCÞ presented in Fig. 6. In all cases, strands of highly stretched
with the PS layer result in the formation of the von Karman vortex polymer molecules are formed upstream of the cylinder, which are
street. For viscoelastic fluids with high Weissenberg numbers and low then extended along both the top and bottom sides of the cylinder sur-
viscosity ratios, stretching of the shear layers can easily be seen as com- face and finally into the downstream wake. These strands become lon-
pared to that seen in Newtonian fluids; for example, see vorticity fields ger, wider, and extend more into the downstream wake as we increase
presented for Wi ¼ 2 and b ¼ 0:5 for all time instances. Furthermore, the Weissenberg number due to the increase in the polymer relaxation
the strength of the base region vortices notably decreases in viscoelastic time ðkÞ or decrease in the viscosity ratio due to the increase in the
fluids. This is again attributed to the stretched polymer molecules polymer contribution of the total viscosity. The formation of these
adhering to the cylinder surface. The reduced base region vorticity strands inhibits the vortex shedding as these polymer molecules sus-
indicates less interaction between the shear layers and the base region tain the stresses for a longer duration, which results in the elongation
vortices. This causes a delay in the onset of the shear layer instability of the vortices and delays the vortex shedding, as discussed above.
and hence delays the vortex shedding phenomenon in viscoelastic One can now expect that the difference in the vortex dynamics
fluids. seen between the Newtonian and viscoelastic fluids would also

FIG. 6. Stretching of polymer molecules for viscoelastic fluids at different Weissenberg numbers and viscosity ratios at different time instances of a periodic cycle, (a) t ¼ 0, (b)
t ¼ T/4, (c) t ¼ T/2, (d) t ¼ 3T/4, as schematically shown in Fig. 3.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-9

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

influence the hydrodynamic drag (FD) and lift (FL) forces acting on
the cylinder. These forces are presented in terms of the dimensionless
drag and lift coefficients represented by CD and CL, respectively, as
defined in Sec. II. For Newtonian fluids, these forces are made of two
components, namely, due to viscous stresses and due to pressure dif-
ference. For viscoelastic fluids, an additional contribution is present
due to the elastic stresses (because of the presence of polymer mole-
cules). To get an insight into how the fluid viscoelasticity influences
these forces acting on the cylinder, first, the temporal variation of the
dimensionless drag and lift coefficients are presented in Figs. 7 and 8,
respectively, for different Weissenberg numbers at a fixed value of
b ¼ 0:5. First, the trend in the temporal variation shown by both the
drag and lift coefficient is typical of the flow past a stationary cylinder
due to the vortex shedding phenomenon. The drag first decreases,
then increases, and finally oscillates within a range of values with a fre-
quency of 2fst (where fst is the natural vortex shedding frequency of a

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

particular fluid). The lift also oscillates with a frequency of fst after the
same initial transient period. Furthermore, the drag value increases,
whereas the lift value decreases with the Weissenberg number. This is
more clear in Fig. 9(a) wherein the time-averaged drag (CD;avg ) and FIG. 8. Temporal variation of CL at different values of Wi and at a fixed value of
the root mean square lift (CL;rms ) coefficients are plotted against the b ¼ 0:5.
Weissenberg number. It can be seen that the value of CD;avg increases
with Wi irrespective of the values of b. Note that here the results for
one value of b ¼ 0:5 is presented as the same trend is observed for
other values of it. This observation is in line with that seen in the previ-
ous studies.18,21 Another important observation from the drag coeffi-
cient’s temporal variation is that the amplitude of fluctuations
decreases with the Weissenberg number. This suppression of the fluc-
tuations in the temporal variation of the drag suggests the suppression
of the vortex shedding phenomenon. Similar behavior is also observed
in the case of the lift coefficient fluctuation, which is also evident in the
plot of CL;rms presented in Fig. 9(a). This behavior is due to the sup-
pression of the boundary layer separation and the presence of the elas-
tic extensional stresses due to the stretching of polymers.49,52,53

FIG. 9. Variation of (a) the time-averaged drag coefficient (CD;avg ) and RMS lift
FIG. 7. Temporal variation of CD at different values of Wi and at a fixed value of coefficient (CL;rms ) and (b) various components of the time-averaged drag coeffi-
b ¼ 0:5. cient with the Weissenberg number (Wi) at a fixed value of b ¼ 0:5.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-10

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

The suppression of the vortex shedding frequency with the TABLE V. Calculated shedding frequencies (s1 ) for viscoelastic fluids with different
increasing Weissenberg number has been quantified by calculating the values of the Weissenberg number and polymer viscosity ratio.
power spectral density (PSD) of the temporal variation of the lift coef-
ficient, as presented in Fig. 10. The shedding frequency shifts away b/Wi 0.5 1 1.5 2
from the Newtonian value of 0.167 s1 to lower values (e.g., 0.12 s1 0.85 0.166 0.159 0.153 0.147
at Wi ¼ 2 and b ¼ 0.5) as we gradually increase the Weissenberg num-
0.70 0.166 0.158 0.147 0.133
ber irrespective of the value of b. However, as the value of b gradually
decreases, the trend becomes more prominent; for instance, see the 0.50 0.160 0.147 0.133 0.120
results presented at b ¼ 0:5 in Fig. 10(c). Furthermore, the magnitude

of the PSD also shows a decreasing trend, which suggests that the vor-
tex strength also decreases. The shedding frequencies calculated for
different viscoelastic fluids are listed in Table V. Furthermore, the vari-
ation of the individual components of the drag forces acting on the
cylinder is shown in Fig. 9 as a function of the Weissenberg number.

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

The pressure drag (CDP) component increases with Wi, whereas the
other two components, namely, friction drag (CDF) and polymer drag
(CDS), decrease.
The temporal variations of both the drag and lift coefficients at
different values of b are presented in Figs. 11 and 12, respectively at a
fixed value of the Weissenberg number of 2. As the value of b
decreases from 1 (Newtonian) to 0.5, the drag (CD) coefficient
increases but the lift coefficient decreases. These are also reflected in
the variations of both the time-averaged drag coefficient and RMS lift
coefficient presented in Fig. 13(a). The temporal fluctuations of both
the drag and lift coefficients flatten out with the decreasing value of b.
This implies that the vortex shedding frequency is suppressed with the
decreasing value of b, i.e., with an increase in the polymer concentra-
tion in the solution. The variation of different components of the drag
coefficient is plotted at different values of b in Fig. 13(a). From this
plot, it is clear that the two components, namely, the pressure drag
and the polymeric stress drag, increase, whereas the friction drag com-
ponent decreases with the increase in the polymer concentration in
the liquid, i.e., with the decreasing value of b.

FIG. 10. Variation of the power spectral density (PSD) of the lift coefficient (CL) with
the Weissenberg number at different values of the polymer viscosity ratio, namely, FIG. 11. Temporal variation of CD at different values of b and at a fixed value of Wi
(a) b ¼ 0:85; (b) b ¼ 0:7; and (c) b ¼ 0:5. of 2.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-11

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

knowledge, this analysis has not been attempted before. This analysis
will also further facilitate understanding of the results on the fluid-
mechanical aspects, which have already been presented in Sec. IV B. A
total of M ¼ 150 snapshots of the z-component of the vorticity field
were generated at a regular non-dimensional time interval of 0.2 units
in the periodic shedding regime for each case to perform the DMD
analysis. This number of snapshots was utilized to ensure the conver-
gence of the results. The analysis was carried out using the procedure
already outlined in Sec. III B. To begin with, the Ritz values of the first
21 modes are shown in Fig. 14. The results for both Newtonian and
viscoelastic fluids with different Weissenberg numbers and viscosity
ratios are included in the same figure. It can be seen that for all the
cases, the Ritz values fall on the unit circle of jkj j ¼ 1. This indicates
that these modes neither decay nor grow with time, i.e., these are stable
modes. The mean mode corresponding to the Ritz value of 1 and of
f d
zero frequency, i.e., Stj ¼ Uj 1 ¼ 0 (where Stj is the non-dimensional

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

frequency and subscript j indicates the mode number) has been
pointed out in Fig. 14. The zero frequency of the mean mode implies
that it does not vary throughout the flow. The rest of the modes are
FIG. 12. Temporal variation of CL at different values of b and at a fixed value of Wi of 2. then numbered from first to jth on the circle in the anti-clockwise
direction. These modes are obtained in complex conjugate pairs,
C. Dynamic mode decomposition (DMD) analysis where the complex conjugates of the modes are encountered when we
In the present study, the dynamic mode decomposition (DMD) traverse the circle in the clockwise direction, Fig. 14. Furthermore, the
analysis was employed to analyze the difference in the coherent flow modes have been recreated using the real part of /j .
structures of Newtonian and viscoelastic fluids. To the best of our The mean modes both for Newtonian and viscoelastic fluids with
different values of Wi and b are visualized in Fig. 15. For a Newtonian
fluid, it can be seen that a concentrated structure is formed in the
vicinity of the front, top, and bottom sides of the cylinder. Two shear
layer type structures are clearly visible from this plot, which are
extended to some distance downstream of the cylinder. Apart from
this, two additional small structures are also seen on the rear side of
the cylinder; see the zoomed-in figure presented in Fig. 15. This is in
accordance with the already stated mechanism for the vorticity

FIG. 13. Variation of (a) the time-averaged drag (CD;avg ) and RMS lift coefficient
(CL;rms ) and (b) different components of the time-averaged drag coefficient with the FIG. 14. Ritz values kj both for Newtonian and viscoelastic fluids with different val-
viscosity ratio (b) at a fixed value of Wi ¼ 2. ues of Wi and b. A total of 21 Ritz values are shown here.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-12

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 15. Visualization of the mean DMD modes of the vorticity fields both for (a) Newtonian and (b)–(e) viscoelastic fluids with different values of Wi and b. The magnified
view (on the right side of each image) shows the contours in the vicinity of the cylinder.

generation in the flow of an incompressible fluid past a bluff body, see trend is also seen for other modes. This is because, in viscoelastic flu-
Sec. IV B. On the other hand, for viscoelastic fluids, this mean mode is ids, the presence of polymer molecules in the fluid develops normal
seen to be more extended downstream of the cylinder than that seen elastic stresses (which leads to high elongational viscosity in these flu-
for Newtonian fluids. This tendency is further accentuated as we ids), suppressing the flow fluctuations and hence reducing the mode
increase the value of the Weissenberg number and decrease the value energy. Based on the energy ranking, we have further plotted the first
of the viscosity ratio. Therefore, this result again confirms that a signif- mode of our DMD analysis in Fig. 17. This mode captures the von
icant stretching in the vortex street occurs in viscoelastic fluids. Karman vortex shedding pattern with alternate signed bubble-like
However, the strength of the minor structures in the rear side of the structures in the cylinder wake both for Newtonian and viscoelastic
cylinder reduces as the fluid viscoelasticity increases. For instance, see fluids. The frequency associated with this mode, i.e., St1, is the natural
Fig. 15(e) for the results at Wi ¼ 2 and b ¼ 0.5, wherein the strength
of these two small structures of the mean mode is almost negligible.
This signifies the suppression of the base region vorticity in viscoelastic
fluids, which was also observed during the discussion of the vorticity
patterns in Sec. IV B. Therefore, much of the flow physics observed
during the fluid-mechanical aspects is further verified from the analy-
sis of the mean DMD mode itself.
Next, the corresponding magnitudes (jj/j jj) and/or energy con-
tributions of the mean mode along with other few modes are pre-
sented in Fig. 16 as a function of the associated non-dimensional
frequencies (St). This plot shows that the mean mode captures the
maximum energy (more than 75%) of the flow field both for
Newtonian and viscoelastic fluids. This testifies to the statement often
made regarding the DMD analysis, i.e., most of the flow phenomena
can be explained based on low rank extracted coherent structures
from the high dimensional flow data.24,44 One can see that the energy
captured by the mean mode of Newtonian fluids is higher than that of
viscoelastic fluids. As the extent of the fluid viscoelasticity increases
(either by increasing the value of the Weissenberg number or decreas- FIG. 16. The variation of the magnitude (and/or the energy contribution) corre-
sponding to the mean mode and first few modes with the associated non-
ing the value of the viscosity ratio), this energy contribution further dimensional frequencies ðStj Þ both for Newtonian and viscoelastic fluids. The mean
decreases. The maximum energy difference between the Newtonian mode captures the maximum energy in all fluids, as can be seen from the chart at
and viscoelastic fluids ðWi ¼ 2; b ¼ 0:5Þ is nearly 30%. A similar St0 ¼ 0.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-13

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 17. Visualization of the DMD mode 1 of the vorticity fields both for (a) Newtonian and (b)–(e) viscoelastic fluids with different values of Wi and b.

vortex shedding frequency as shown in Fig. 18. For instance, the value than that in Newtonian fluids. This was already seen in the PSD analy-
of St1 ¼ 0.167 for Newtonian fluids. This value is the same as that cal- sis of the temporal variation of the lift coefficient presented in Fig. 10.
culated earlier using the PSD of the temporal variation of the lift coeffi- It should be mentioned here that these frequencies of the DMD modes
cient (CL) in Sec. IV B. Likewise, these values are also the same for are calculated from the global flow field data, whereas the PSD utilizes
viscoelastic fluids as that presented in Table V. This suggests that the the time variation of the local lift forces acting on the cylinder.
dominant frequency in the wake downstream of the cylinder is the Therefore, the DMD analysis can provide information on different
corresponding natural shedding frequency of the respective fluids. It flow phenomena (like the vortex shedding phenomenon in this case)
can be seen from Fig. 18 that the frequency of all DMD modes utilizing only the visual representation of the global flow field using
decreases as the fluid viscoelasticity increases. Therefore, it naturally noninvasive measurements. This shows one of the potential applica-
shows that the frequency of the vortex shedding phenomena (which is tions of the DMD analysis in the domain fluid mechanics.24
captured by the DMD mode 1) in viscoelastic fluids should be less Furthermore, the suppression of the vortex shedding phenomena
and the base region vorticity can also be visualized from these modes;
for instance, see the results presented in Fig. 17 for mode 1. As we
gradually move from Newtonian to viscoelastic fluids, the bubble-like
structures become substantially larger, resulting in less number of
structures at a given distance downstream of the cylinder. Additionally,
small structures attached to the rear surface of the cylinder also vanish
as the fluid viscoelasticity increases. Therefore, it can be seen that the
modes of Newtonian and viscoelastic fluids differ in both near bubble-
like structure and their intensity in the downstream wake. The changes
in the near structure are caused due to the vortex stretching phenome-
non, as already discussed above. The high intensity of the DMD struc-
tures in viscoelastic fluids far away downstream indicates that the
dissipation of the vortices is slower in viscoelastic fluids than that seen
in Newtonian fluids, for instance, see Figs. 17(a) and 17(e) for the
results of Newtonian and viscoelastic fluids with Wi ¼ 2 and b ¼ 0:5,
respectively. It implies that the vortex street in viscoelastic fluids
extends to a larger distance downstream than Newtonian fluids under
otherwise identical conditions. It was not that much of prominent in
FIG. 18. The normalized magnitudes corresponding to the most dominant modes
(excluding the mean mode) are plotted against the non-dimensional frequency of the vorticity analysis presented in Sec. IV B. However, it is now more
each mode. The red dashed line represents the natural vortex shedding frequency clear because of the DMD analysis. Figure 19 presents the energy
(fst ¼ 0:167 s1 ) of Newtonian fluid. magnitude and the total energy as a function of the number of modes.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-14

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

FIG. 19. Variation of the energy contribution (a) and cumulative energy contribution (b) with the number of modes.

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

FIG. 20. Visualization of the DMD mode 2
of the vorticity fields both for (a)
Newtonian and (b)–(e) viscoelastic fluids
with different values of Wi and b.

It is clearly visible that most of the energy of the flow field is captured
TABLE VI. Summary of the flow features extracted from different DMD modes in the
by the first few modes, as already pointed out earlier. In fact, the mean case of flow of viscoelastic fluids.
mode and the first mode capture more than 90% of the energy. The
higher modes are just the harmonics of the first mode and do not Mode Features identified
show additional new flow physics. For example, see mode 2 presented
in Fig. 20, which has an anti-symmetric bubble-like structure, and the Mean mode Stretching in the shear layer type flow structure
frequency associated with it is the first harmonic of the natural vortex Suppression of the base region structure con-
shedding frequency. Once again, the DMD structures formed near the centration
cylinder and far away downstream of the cylinder have the same trend Suppression of the fluctuations in the flow field
as that seen for mode 1 presented in Fig. 17. A summary of different Mode 1 Convected vortex street and associated frequen-
flow features extracted from different DMD modes is presented in cies
Table VI. Suppression of the Karman vortex shedding fre-
quency
V. CONCLUSIONS Suppression of the vortex dissipation far down-
This work has discussed the critical changes in the coherent stream of the cylinder
flow structures that occur in the flow past a circular cylinder in the Mode 2 and Anti-symmetric nature of the vortex street
laminar vortex shedding regime as a consequence of introducing a higher modes Higher harmonics and finer-scale structures
minute amount of solid polymers into a Newtonian solvent like water.

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-15

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

8
In performing so, extensive numerical simulations have been con- H. A. Barnes, A Handbook of Elementary Rheology (University of Wales,
ducted for Oldroyd-B viscoelastic fluids over a range of Weissenberg Institute of Non-Newtonian Fluid Mechanics, Aberystwyth, 2000), Vol. 1.
9
A. Y. Malkin and A. I. Isayev, Rheology: Concepts, Methods, and Applications
numbers as 0  Wi  2 and polymer viscosity ratios as 0:5  b 
(Elsevier, 2017).
0:85 at a fixed value of the Reynolds number of 100. The dynamic 10
R. B. Bird, R. C. Armstrong, and O. Hassager, Dynamics of Polymeric Liquids.
mode decomposition (DMD) technique has been used in this study Vol. 1: Fluid Mechanics (John Wiley and Sons Inc., New York, 1987).
11
for post-processing the results, which is considered one of the widely G. Gadd, “Reduction of turbulent friction in liquids by dissolved additives,”
used reduced order modeling (ROM) techniques to analyze the coher- Nature 212, 874–877 (1966).
12
G. Gadd, “Effects of long-chain molecule additives in water on vortex streets,”
ent flow structures in various flow dynamical problems. This study has
Nature 211, 169–170 (1966).
shown that very few DMD modes can effectively capture the gradual 13
E. Burger, L. Chorn, and T. K. Perkins, “Studies of drag reduction conducted
change in the von Karman vortex street pattern with the increasing over a broad range of pipeline conditions when flowing Prudhoe Bay crude
values of fluid viscoelasticity. In particular, the stretching of shear layer oil,” J. Rheol. 24, 603–626 (1980).
14
structures and suppression of the vortex shedding frequency in visco- O. Cadot and M. Lebey, “Shear instability inhibition in a cylinder wake by local
elastic fluids have been observed from the DMD analysis. The energy injection of a viscoelastic fluid,” Phys. Fluids 11, 494–496 (1999).
15
O. Cadot and S. Kumar, “Experimental characterization of viscoelastic effects
associated with the corresponding DMD modes decreases with the on two-and three-dimensional shear instabilities,” J. Fluid Mech. 416, 151–172
fluid viscoelasticity, which is attributed to the suppression of the veloc- (2000).

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

16
ity fluctuations in viscoelastic fluids compared to Newtonian fluids. J. Cressman, Q. Bailey, and W. Goldburg, “Modification of a vortex street by a
Furthermore, suppression of the base region vorticity occurs, and the polymer additive,” Phys. Fluids 13, 867–871 (2001).
17
corresponding DMD structures disappear as the fluid viscoelasticity V. Kalashnikov and A. Kudin, “Karman vortices in the flow of drag-reducing
polymer solutions,” Nature 225, 445–446 (1970).
gradually increases. A detailed discussion on the fluid-mechanical 18
P. J. Oliveira, “Method for time-dependent simulations of viscoelastic flows:
aspects, such as streamline patterns, vorticity fields, drag and lift forces Vortex shedding behind cylinder,” J. Non-Newtonian Fluid Mech. 101,
acting on the cylinder, etc., is also presented in this study. This facili- 113–137 (2001).
19
tates the importance of the DMD analysis to get further insight into M. Sahin and R. G. Owens, “On the effects of viscoelasticity on two-
flow physics. dimensional vortex dynamics in the cylinder wake,” J. Non-Newtonian Fluid
Mech. 123, 121–139 (2004).
20
AUTHOR DECLARATIONS A. Minaeian, M. Nili-AhmadAbadi, M. Norouzi, and K. C. Kim, “Effects of vis-
coelasticity on the onset of vortex shedding and forces applied on a cylinder in
Conflict of Interest
unsteady flow regime,” Phys. Fluids 34, 013106 (2022).
21
The authors have no conflicts to disclose. D. Richter, G. Iaccarino, and E. S. Shaqfeh, “Simulations of three-dimensional
viscoelastic flows past a circular cylinder at moderate Reynolds numbers,”
J. Fluid Mech. 651, 415–442 (2010).
Author Contributions 22
C. Sasmal, M. B. Khan, and R. Chhabra, “Combined influence of fluid visco-
Faheem Hamid: Data curation (lead); Formal analysis (lead); elasticity and inertia on forced convection heat transfer from a circular cylin-
Investigation (lead); Methodology (lead); Software (lead); Validation der,” J. Heat Transfer 142, 041801 (2020).
23
(lead); Writing – original draft (lead). Chandi Sasmal: Conceptualization M. Norouzi, S. Varedi, M. J. Maghrebi, and M. Shahmardan, “Numerical inves-
tigation of viscoelastic shedding flow behind a circular cylinder,” J. Non-
(lead); Investigation (equal); Supervision (lead); Writing – review & Newtonian Fluid Mech. 197, 31–40 (2013).
editing (lead). Raj Chhabra: Formal analysis (supporting); Investigation 24
J. N. Kutz, S. L. Brunton, B. W. Brunton, and J. L. Proctor, Dynamic Mode
(supporting); Writing – original draft (supporting); Writing – review & Decomposition: Data-Driven Modeling of Complex Systems (SIAM, 2016).
25
editing (supporting). S. L. Brunton, B. R. Noack, and P. Koumoutsakos, “Machine learning for fluid
mechanics,” Annu. Rev. Fluid Mech. 52, 477–508 (2020).
26
H. Ping, H. Zhu, K. Zhang, D. Zhou, Y. Bao, Y. Xu, and Z. Han, “Dynamic
DATA AVAILABILITY mode decomposition based analysis of flow past a transversely oscillating cylin-
der,” Phys. Fluids 33, 033604 (2021).
The data that support the findings of this study are available 27
Q. Zhang, Y. Liu, and S. Wang, “The identification of coherent structures using
within the article. proper orthogonal decomposition and dynamic mode decomposition,”
J. Fluids Struct. 49, 53–72 (2014).
28
Z. Wu, S. L. Brunton, and S. Revzen, “Challenges in dynamic mode decom-
REFERENCES
position,” J. R. Soc. Interface 18, 20210686 (2021).
1 29
M. M. Zdravkovich, Flow Around Circular Cylinders: Volume 1: Fundamentals H. Eivazi, H. Veisi, M. H. Naderi, and V. Esfahanian, “Deep neural networks
(Oxford University Press, 1997), Vol. 2. for nonlinear model order reduction of unsteady flows,” Phys. Fluids 32,
2
M. M. Zdravkovich, Flow Around Circular Cylinders: Volume 2: Applications 105104 (2020).
(Oxford University Press, 1997), Vol. 2. 30
Y. Liu, J. Long, Q. Wu, B. Huang, and G. Wang, “Data-driven modal decompo-
3
G. E. Karniadakis and G. S. Triantafyllou, “Frequency selection and asymptotic sition of transient cavitating flow,” Phys. Fluids 33, 113316 (2021).
states in laminar wakes,” J. Fluid Mech. 199, 441–469 (1989). 31
Y. Liu, B. Huang, H. Zhang, Q. Wu, and G. Wang, “Experimental investigation
4
C. H. Williamson, “Vortex dynamics in the cylinder wake,” Annu. Rev. Fluid into fluid–structure interaction of cavitating flow,” Phys. Fluids 33, 093307
Mech. 28, 477–539 (1996). (2021).
5 32
S. Krishnamoorthy, S. Price, and M. Paidoussis, “Cross-flow past an oscillating P. J. Schmid, “Application of the dynamic mode decomposition to experimen-
circular cylinder: Synchronization phenomena in the near wake,” J. Fluids tal data,” Exp. Fluids 50, 1123–1130 (2011).
33
Struct. 15, 955–980 (2001). K. K. Chen, J. H. Tu, and C. W. Rowley, “Variants of dynamic mode decompo-
6
J. Carberry, J. Sheridan, and D. Rockwell, “Controlled oscillations of a cylinder: sition: Boundary condition, Koopman, and Fourier analyses,” J. Nonlinear Sci.
Forces and wake modes,” J. Fluid Mech. 538, 31–69 (2005). 22, 887–915 (2012).
7 34
R. P. Chhabra and J. F. Richardson, Non-Newtonian Flow and Applied M. S. Hemati, M. O. Williams, and C. W. Rowley, “Dynamic mode decomposi-
Rheology: Engineering Applications (Butterworth-Heinemann, 2011). tion for large and streaming datasets,” Phys. Fluids 26, 111701 (2014).

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-16

Published under an exclusive license by AIP Publishing
 Physics of Fluids ARTICLE scitation.org/journal/phf

35 45
K. Taira, M. S. Hemati, S. L. Brunton, Y. Sun, K. Duraisamy, S. Bagheri, S. T. R. Fattal and R. Kupferman, “Time-dependent simulation of viscoelastic flows
Dawson, and C.-A. Yeh, “Modal analysis of fluid flows: Applications and out- at high Weissenberg number using the log-conformation representation,”
look,” AIAA J. 58, 998–1022 (2020). J. Non-Newtonian Fluid Mech. 126, 23–37 (2005).
36 46
M. H. Naderi, H. Eivazi, and V. Esfahanian, “New method for dynamic mode R. P. Bharti, R. Chhabra, and V. Eswaran, “Steady flow of power law fluids
decomposition of flows over moving structures based on machine learning across a circular cylinder,” Can. J. Chem. Eng. 84, 406–421 (2006).
47
(hybrid dynamic mode decomposition),” Phys. Fluids 31, 127102 (2019). A. Soares, J. Ferreira, and R. Chhabra, “Flow and forced convection heat trans-
37
W. Zhang, Y. Wang, and Y.-H. Qian, “Stability analysis for flow past a cylinder fer in crossflow of non-Newtonian fluids over a circular cylinder,” Ind. Eng.
via lattice Boltzmann method and dynamic mode decomposition,” Chin. Phys. Chem. Res. 44, 5815–5827 (2005).
48
B 24, 064701 (2015). B. Fornberg, “A numerical study of steady viscous flow past a circular cylin-
38
R. E. D. Bishop and A. Hassan, “The lift and drag forces on a circular cylinder der,” J. Fluid Mech. 98, 819–855 (1980).
49
oscillating in a flowing fluid,” Proc. R. Soc. London, Ser. A 277, 51–75 (1964). M. Alves, F. Pinho, and P. Oliveira, “The flow of viscoelastic fluids past a cylin-
39
J. G. Oldroyd, “On the formulation of rheological equations of state,” Proc. R. der: Finite-volume high-resolution methods,” J. Non-Newtonian Fluid Mech.
Soc. London, Ser. A 200, 523–541 (1950). 97, 207–232 (2001).
40 50
F. Pimenta and M. Alves, see https://github.com/fppimenta/rheoTool for P. Pakdel and G. H. McKinley, “Elastic instability and curved streamlines,”
“Rheotool” (last accessed April 22, 2022). Phys. Rev. Lett. 77, 2459 (1996).
41 51
M. Alves, P. Oliveira, and F. Pinho, “A convergent and universally bounded S. S. Datta, A. M. Ardekani, P. E. Arratia, A. N. Beris, I. Bischofberger, G. H.
interpolation scheme for the treatment of advection,” Int. J. Numer. Methods McKinley, J. G. Eggers, J. E. L opez-Aguilar, S. M. Fielding, A. Frishman et al.,
Fluids 41, 47–75 (2003). “Perspectives on viscoelastic flow instabilities and elastic turbulence,” Phys.

Downloaded from http://pubs.aip.org/aip/pof/article-pdf/doi/10.1063/5.0122103/16574089/103114_1_online.pdf

42
R. Fattal and R. Kupferman, “Constitutive laws for the matrix-logarithm of the Rev. Fluids 7, 080701 (2022).
52
conformation tensor,” J. Non-Newtonian Fluid Mech. 123, 281–285 (2004). M. B. Khan, C. Sasmal, and R. Chhabra, “Flow and heat transfer characteristics
43
F. Pimenta and M. Alves, “Stabilization of an open-source finite-volume solver for of a rotating cylinder in a FENE-P type viscoelastic fluid,” J. Non-Newtonian
viscoelastic fluid flows,” J. Non-Newtonian Fluid Mech. 239, 85–104 (2017). Fluid Mech. 282, 104333 (2020).
44 53
P. J. Schmid, “Dynamic mode decomposition of numerical and experimental H. H. Hu and D. D. Joseph, “Numerical simulation of viscoelastic flow past a
data,” J. Fluid Mech. 656, 5–28 (2010). cylinder,” J. Non-Newtonian Fluid Mech. 37, 347–377 (1990).

Phys. Fluids 34, 103114 (2022); doi: 10.1063/5.0122103 34, 103114-17

Published under an exclusive license by AIP Publishing