# Finite Element Analysis of the Dynamics of Power-Law Fluid around an Obstacle in a Channel

Authors: Rashid Mahmood; Mohammed M. M. Jaradat; Rifaqat Ali; Imran Siddique; Afraz Hussain Majeed; Asif Mehmood
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Control of fluid forces is an emerging area of research with numerous engineering applications. The uneven wake behind an obstacle causes undesirable structural oscillations, which can lead to fatigue or structural failure. Controlling the wake phenomena could directly benefit a wide range of engineering applications, including skyscrapers, naval risers, bridges, columns, and a few sections of airplanes. This study is concerned with the time dependent simulations in a channel in presence of an obstacle aiming to compute fluid forces. The underlying mathematical model is based on nonstationary Navier–Stokes equations coupled with the constitutive relations of power law fluids. Because the representative equations are complex, an effective computing strategy based on the finite element approach is used. To achieve higher accuracy, a hybrid computational grid at a very fine level is used. Theℙ2−ℙ1 elements based on the shape functions of the second and first-order polynomials were used to approximate the solution. The discrete nonlinear system arising from this discretization is linearized by Newton’s method and then solved through a direct linear solver PARADISO. The code validation study is also performed for Newtonian fluids as a special case, and then the study is extended to compute drag and lift forces for other cases of viscosity as described by the power law index. When looking at the phase plot, it can be seen that for the Newtonian case n = 1, there is only one closed orbit after the steady state is reached, whereas for n=0.5, there are multiple periodic orbits. Moreover, the effects of shear rate on the drag-lift phase plot are also discussed.


---

## Body

## 1. Introduction

The significance of flow around bluff bodies has been realized since the last couple of centuries, and hence it is considered as benchmark problem in fluid dynamics. Hence, researchers investigated different shapes of bluff bodies submerged in different fluids for their analysis. It has been proved both experimentally and analytically that hydrodynamic forces like drag and lift around obstacles play an important role in designing structures in many engineering applications. It is also worthy to mention that the interaction of fluid flow with such bluff bodies produces wake, and the formation of wake around these complex geometries is germane to understand due to its practical utilization. Some investigations are available in recent literature that realize the importance of such configurations in theoretical as well as practical aspects. Hussain et al. [1] considered a circular cylinder to investigate the incompressible flow behavior. The nonlinear fluid for specific values of n is used by Chhabra et al. [2] to investigate the blockage ratios. In [3], shear thinning and shear thickening behavior are observed for the range 100≤Re≤500. Drag and lift for a circular cylinder placed in a channel are studied in [4]. Lid-driven flow behavior is observed in [5]. A comparison of drag and lift for circular and square cylinders can be seen in [6]. Williamson [7] summarized some of the work about this benchmark. The results proposed by Schäfer et al. [8] consider Newtonian fluid at various Reynolds numbers for hydrodynamic forces. They justified the credibility of their work by constructing a comparison with experimental data and also using different computational approaches for attaining the solution. Fluid past a rotating and stationary cylinder is considered in [9, 10] to perform the numerical simulations. They computed instantaneous streamlines along with drag and lift. At low Reynolds number 0≤Re≤50, Tritton [11] performed analysis on the interaction of fluid flow with a circular obstacle and measured hydrodynamic forces. To demarcate the flow regime, Coelho and Pinho [12] used shear thinning fluids for elucidating the vortex shedding. Gupta et al. [13] investigated the steady flow features in Newtonian fluid for Reynolds number 5≤Re≤40. For sake of brevity, some more investigations across circular cylindrical obstacles are referred to in [14–18]. In [14], Bharti et al. presented a numerical study of heat transfer using the finite volume method. Some aspects of flow behavior around a circular cylinder are presented in [15]. The 2D and 3D computational analysis of fluid flow around a cylinder is discussed in [16]. Rajani et al. [17] used different values of Reynolds number to observe the impact of flow behavior around a circular cylinder. Very slow flow around a cylinder is numerically investigated in [18].A lot of numerical work has been accomplished to explore drag and lift forces on obstacles in the Newtonian flow field, but for non-Newtonian cases, the investigations are still at an embryonic stage to assess the relationship between drag and lift, especially for highly nonlinear viscosity. For nonlinear viscosity fluids, the flow around cylinders of different shapes is an interesting topic of research due to the investigation of the wake, recirculation zone length, and drag and lift characteristics. To analyse the laminar flow of non-Newtonian fluid around a circular cylinder, Chhabra et al. [2] used three values of Reynolds (1, 20, and 40). They used a finite difference scheme to observe the flow pattern. Mossaz et al. [19] used the most generic Herschel–Bulkley fluid model to evaluate the flow pattern around a circular cylinder. In [20], Nejat et al. used power law fluid for their numerical investigations. With the help of the Lattice Boltzman method, they computed drag and lift coefficients. In the low Reynolds number range 5≤Re≤40, the relationship between n and drag coefficient is investigated in [21]. Whitney and Rodin [22] proposed some results for shear thinning fluids. These results are compatible with some previous studies of flow around bluff bodies. The power law fluid past a heated cylinder is used in [23] to observe heat transfer characteristics. In order to eliminate the impact inlet and outlet on the results, some researchers [24–26] used very large computational domains to reduce such effects. The results show that for a circular cylinder, drag is inversely proportional to the computational grid. Abegunrin et al. [27] investigated a non-Newtonian fluid simulation on the upper horizontal surface of a paraboloid of revolution. It is deduced that when the flow is described as Newtonian fluid flow, the maximum velocity of the flow can be determined. Similarly, Patnana et al. [10] considered unsteady power-law (PL) fluid past a circular cylinder and calculated some results using the finite volume method.In the literature, there is only a small contribution towards non-Newtonian fluids moving around obstacles and the computation of fluid forces. There is a need to see the impact of shear-thinning and thickening behavior on drag and lift forces. The aim of this article is to elucidate the characteristics of unsteady nonlinear fluid flow in a channel over a cylinder and the estimation of fluid forces. The present work is organized as follows: a limited survey highlighting the importance of flow configuration is presented in Section1. In Section 2, the adopted research methodology is discussed with details of all solvers. The simulation results and analysis through graphs and phase-plots are performed in Section 3. The present attempt is concluded in Section 4.


## 2. Research Methodology

### 2.1. Mathematical Formulation

The mathematical model for an incompressible, isothermal, and nonlinear viscosity fluid model is given by the following set of nonlinear equations, as shown in [30]:(1)∂ui∂xi=0,(2)∂ui∂t+∂∂xjuiuj+1ρ∂p∂xi=1ρ∂∂xjμγ˙∂ui∂xj+∂uj∂xi,i,j=1…3,we couple equations (1) and (2) with the well -known power-law model as(3)μγ˙=mγ˙n−1,where m represents the coefficient of fluid consistency; n is the power law exponent.


### 2.2. Problem Description and Numerical Approach

The computational domain is a channelΩ=0,2.2×0,0.41, as considered in [28–30]. The circular obstacle is placed at (0.2, 0.2) that is offset from the center to generate a nonzero lift. The top and bottom walls offer no slip conditions as u = v = 0. An inflow parabolic profile with maximum velocity as Umax=1.5 is exposed at the inlet of the channel to produce a periodic flow regime and a zero Neumann condition at the outlet. More elements are considered around the obstacle to accurately capture the hydrodynamic forces on the cylinder.The coarse grid for present simulations is shown in Figure1. To achieve higher accuracy, one element is mapped into four dotted elements of smaller size by joining the midpoints of the sides of the parent element through h-refinement. Following this strategy of refinement, Table 1 shows the details about the cardinality of elements and degrees of freedom (DOF) at various computational refinement levels.Figure 1
Computational grid at level-1.Table 1
Grid statistics at various levels.
Levels# Elements#DOF19627378215901200932426179184454032683567164761161179681603727194187374865288444873The incompressible generalized Navier–Stokes equations (1) and (2), as well as the rheological law (3), or power law, can be used to explain a variety of simulations. These equations help us understand the flow of materials in nature by describing particular transactions. We use finite element simulation for the numerical estimation of the governing equations because exact solutions to such problems are rare due to the model’s high nonlinearity. The conforming element pair ℙ2−ℙ1 is chosen for the velocity and pressure approximations in this direction. This element is a stable element pair [31–36] that meets the inf-sup condition. The inner linear subproblems are solved using a direct solver, and discrete nonlinear algebraic systems are solved using Newton’s method. The stopping criterion for the nonlinear iteration is set through a bound on the relative error as(4)ℵr+1−ℵrℵr+1<10−6,where ℵ represents one of the component of solution u,v,p.


## 2.1. Mathematical Formulation

The mathematical model for an incompressible, isothermal, and nonlinear viscosity fluid model is given by the following set of nonlinear equations, as shown in [30]:(1)∂ui∂xi=0,(2)∂ui∂t+∂∂xjuiuj+1ρ∂p∂xi=1ρ∂∂xjμγ˙∂ui∂xj+∂uj∂xi,i,j=1…3,we couple equations (1) and (2) with the well -known power-law model as(3)μγ˙=mγ˙n−1,where m represents the coefficient of fluid consistency; n is the power law exponent.


## 2.2. Problem Description and Numerical Approach

The computational domain is a channelΩ=0,2.2×0,0.41, as considered in [28–30]. The circular obstacle is placed at (0.2, 0.2) that is offset from the center to generate a nonzero lift. The top and bottom walls offer no slip conditions as u = v = 0. An inflow parabolic profile with maximum velocity as Umax=1.5 is exposed at the inlet of the channel to produce a periodic flow regime and a zero Neumann condition at the outlet. More elements are considered around the obstacle to accurately capture the hydrodynamic forces on the cylinder.The coarse grid for present simulations is shown in Figure1. To achieve higher accuracy, one element is mapped into four dotted elements of smaller size by joining the midpoints of the sides of the parent element through h-refinement. Following this strategy of refinement, Table 1 shows the details about the cardinality of elements and degrees of freedom (DOF) at various computational refinement levels.Figure 1
Computational grid at level-1.Table 1
Grid statistics at various levels.
Levels# Elements#DOF19627378215901200932426179184454032683567164761161179681603727194187374865288444873The incompressible generalized Navier–Stokes equations (1) and (2), as well as the rheological law (3), or power law, can be used to explain a variety of simulations. These equations help us understand the flow of materials in nature by describing particular transactions. We use finite element simulation for the numerical estimation of the governing equations because exact solutions to such problems are rare due to the model’s high nonlinearity. The conforming element pair ℙ2−ℙ1 is chosen for the velocity and pressure approximations in this direction. This element is a stable element pair [31–36] that meets the inf-sup condition. The inner linear subproblems are solved using a direct solver, and discrete nonlinear algebraic systems are solved using Newton’s method. The stopping criterion for the nonlinear iteration is set through a bound on the relative error as(4)ℵr+1−ℵrℵr+1<10−6,where ℵ represents one of the component of solution u,v,p.


## 3. Results and Discussions

The periodic flow of incompressible and unsteady power law fluids is numerically investigated by implementing the finite element method. For the validation of the code, we first simulate the case forn = 1, Re = 100 which corresponds to the Newtonian case for which the reference data is available in [8]. The viscosity of time-independent nonlinear fluids is only affected by shear rate and temperature. Shear-thinning, shear thickening, and plasticization viscosity increase with increased shear rate. The visual impact of drag and lift with respect to time can be seen in Figures 2–4 for all three cases of n. The considered time interval 0<t≤10 seconds is enough for hydrodynamic forces to get stable. Overall, there is not much variation in the amplitude of drag and lift from 5 seconds onwards. Lift coefficient appeared to be symmetric CL=0 on the horizontal axis for shear thinning (n=0.5) and Newtonian case (n=1). In the shear thinning case (n=0.5), the drag oscillates between 2.3 and 4 (Figure 3), while the lift starts from zero and varies between −4 and +4. From Figure 4(a) and Figure 4(b) it is seen that both drag and lift coefficients converge to respective steady state constant values for n = 1.5 due to enhanced viscosity in shear-thickening case. The lift forces are most potent in the upward trend, as indicated by the negative CL. Even the obstacle is arranged as the lift coefficient is negative, resulting in a computational configuration like this.Figure 2
CDandCL values at n=1.0. (a) CD. (b) CL.
(a)(b)Figure 3
CDandCL values at n=0.5. (a) CD. (b) CL.
(a)(b)Figure 4
CDandCL values at n=1.5. (a) CD. (b) CL.
(a)(b)The velocity profiles at various time steps forn = 1, 0.5 and 1.5 are presented in Figures 5–7, respectively. The corresponding pressure plots are shown in Figures 8(a)–8(c). In Figure 5, using n = 1, we provided the time dependent fluid with a parabolic profile having velocity Umax=1.5. It can be seen that as the fluid strikes the circular cylinder, at t=0.1 it starts bifurcating and approximately at t=3 approach the outlet of the channel.Figure 5
Velocity profiles forn=1.0 at various values of time. (a) t=0.1. (b) t=0.5. (c) t=1. (d) t=2. (e) t=3. (f) t=4. (g) t=5. (h) t=6. (i) t=7. (j) t=8.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)Figure 6
Velocity profiles forn=0.5 at various values of time. (a) t=0.1. (b) t=0.5. (c) t=1. (d) t=2. (e) t=3. (f) t=4. (g) t=5. (h) t=6. (i) t=7. (j) t=8.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)Figure 7
Velocity profiles forn=1.5 at various values of time. (a) t=0.1. (b) t=0.5. (c) t=1. (d) t=2. (e) t=3. (f) t=4. (g) t=5. (h) t=6. (i) t=7. (j) t=8.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)Figure 8
Pressure distribution in a channel for variousn at t = 8. (a) n=0.5, (b) n=1.0, and (c) n=1.5.
(a)(b)(c)The outcome of the velocity profile for the shear thinning case is shared in Figure6. The flow is fully developed approximately for t≥3. Figure 7 reports the velocity variations for shear thickening fluid. Here the flow is fully developed and vortex formation does not occur due to high viscosity. At t = 8, the pressure distribution for each case of n is offered in Figure 8.The periodic flow atRe=100 has its impact on drag and lift patterns, as shown in Figure 9 in terms of phase-plots. It is observed that for the Newtonian case n = 1, there exists a single closed orbit after the steady state is achieved; however, there exist multiple periodic orbits for n = 0.5. The case of n = 1.5 is different where there is not any closed trajectory. Physically, the Reynolds number corresponds to a direct increase in the inertial forces capable of laminar flow for a fixed viscous fluid [37]; therefore, in future works, we would like to increase the Reynolds number to see its impact on CDandCL and on the Strouhal number elaborating the frequency of vortex shedding for shear thinning and thickening cases. [38, 39].Figure 9
Phase-portrait (drag versus lift) (a)n=0.5, (b) n=1.0, and (c) n=1.5.
(a)(b)(c)


## 4. Conclusions

We have considered the two-dimensional transient flow of an incompressible power law fluid over a circular cylinder inside a channel. A static circular cylinder of diameter 0.1 m is used as an obstacle. The two-dimensional flow equations solved using the finite element method (FEM). Power law indexn is considered for three particular values (n = 0.5, 1.0, and 1.5) at base Reynolds number Re = 100. For the time interval 0≤t≤10 seconds, the behavior of velocity and pressure is observed. Moreover, the effects of shear rate on the drag-lift phase plot are also discussed. The diagram of drag-lift phase provides the precise information about the variations of CDandCL. From the phase plots, it is concluded that for the Newtonian case n = 1, there is only one closed orbit after the steady state is reached, while for n=0.5, there are multiple periodic orbits. For the case of shear thickening n = 1.5, there is not any closed trajectory.



---

*Source: 3160951-2022-04-27.xml*