Numerical verification of sharp corner behavior for Giesekus and Phan-Thien–Tanner fluids

Evans, J. D.; Palhares Junior, I. L.; Oishi, C. M.; Ruano Neto, F.

doi:10.1063/5.0125940

We verify numerically the theoretical stress singularities for two viscoelastic models that occur at sharp corners. The models considered are the Giesekus and Phan-Thien–Tanner (PTT), both of which are shear thinning and are able to capture realistic polymer behaviors. The theoretical asymptotic behavior of these two models at sharp corners has previously been found to involve an integrable solvent and polymer elastic stress singularity, along with narrow elastic stress boundary layers at the walls of the corner. We demonstrate here the validity of these theoretical results through numerical simulation of the classical contraction flow and analyzing the $270 °$ corner. Numerical results are presented, verifying both the solvent and polymer stress singularities, as well as the dominant terms in the constitutive equations supporting the elastic boundary layer structures. For comparison at Weissenberg order one, we consider both the Cartesian stress formulation and the alternative natural stress formulation of the viscoelastic constitutive equations. Numerically, it is shown that the natural stress formulation gives increased accuracy and convergence behavior at the stress singularity and, moreover, encounters no upper Weissenberg number limitation in the global flow simulation for sufficiently large solvent viscosity fraction. The numerical simulations with the Cartesian stress formulation cannot reach such high Weissenberg numbers and run into convergence failure associated with the so-called high Weissenberg number problem.

I. INTRODUCTION

The construction of approximate solutions through asymptotic techniques is widely adopted in fluid mechanics, for example, van Dyke (1964), Lagerstrom (1988), Ockendon and Ockendon (1995), and Schlichting and Gersten (2000). It is particularly useful in identifying and quantifying singular behavior, such as boundary layers and stress singularities. In a mathematical context, asymptotic analysis is a powerful technique, often providing accurate approximate solutions where numerical schemes can encounter difficulty or even fail. Here, we consider such a situation that arises for viscoelastic fluids in a contraction domain that possesses sharp corners, where both stress singularities and boundary layers are encountered.

A variety of fluid flow problems include sharp corners at the boundaries of the flow domain. Examples of such cases include flows around forward- or backward-facing steps (Demuren and Wilson, 1994; Wang, 1994), cavity flows (both over and inside) (Gatski and Grosch, 1985), and flows in expanding or contracting channels (Fearn et al., 1990; Shapira et al., 1990; Battaglia et al., 1997; Alleborn et al., 1997; and Drikakis, 1997). These situations have received much attention for Newtonian fluids, where the sharp change in the slope of the boundary at the corner causes singularities in both the vorticity and stress. The identification of the singularity goes back to the work of Dean and Montagnon (1949) and Moffatt (1964). For Newtonian fluids, these corner singularities not only dominate the flow behavior near the corner, but also can extend to the entire domain, where reversed flow regions and flow separation can occur. Such sharp corner singularities can also create serious numerical difficulties in accurately resolving the flow field around the corner. Convergence of numerical schemes is often suppressed, unless specialized techniques, such as mesh refinement or singular elements, are used. Moreover, these difficulties may affect the accuracy of the overall numerical simulation in the whole computational domain. The approach of incorporating the stress singularity within an analytical expansion (Phillips, 1989) or numerical scheme (Holstein and Paddon, 1981; Floryan and Czechowski, 1995; Hawa et al., 2002; Egger et al., 2014; and Deliceoğlu et al., 2019) for Newtonian fluids has been particularly successful.

The flow of viscoelastic fluids in such flow geometries though has posed far more problems. So much so that the contraction flow problem became a benchmark problem for the testing of numerical schemes (Hassager, 1988) back in 1987. A comprehensive discussion and reference for early work can be found in Owens and Phillips (2002), while the more recent literature can be found in Niethammer et al. (2018), Alves et al. (2021), Afonso et al. (2011), Fernandes et al. (2017), and Davoodi et al. (2022). The source of the difficulties can be apportioned to (1) a lack of knowledge of the stress singularity at sharp corners, and (2) the more challenging mathematical nature of the viscoelastic constitutive equations.

First, in regard to the stress singularity, it has only been relatively recently, since the work of Hinch (1993) and subsequent authors (Renardy, 1995, 1997; Rallison and Hinch, 2004; and Evans, 2005a; 2005b; 2008a; 2008b; 2010a; 2010b), that the asymptotic behavior has been determined for the UCM model, Oldroyd-B, PTT, and Giesekus. These models all share the common feature of narrow stress boundary layers at the walls, upstream and downstream of the singularity, which is a feature that Newtonian fluids do not possess. The Newtonian asymptotic stress behavior is uniformly valid around a sharp corner, a far simpler structure than the three-region asymptotic structure of the viscoelastic models. Second, the viscoelastic constitutive equations when coupled with the momentum equation give a mixed elliptic–hyperbolic system of equations (Joseph, 1990; Owens and Phillips, 2002; and Gerritsma and Phillips, 2001; 2008). This hyperbolic aspect of the viscoelastic equations ensures that the effect of any degradation in the numerical solution does not remain local, unlike Newtonian fluids (Renardy, 1989). This accounts for the immense difficulty that early approaches encountered. The problems of the stress singularity and narrow stress boundary layers are exacerbated as the relaxation time (elastic memory) in the viscoelastic models increases. This is notoriously referred to as the high Weissenberg number problem (HWNP) (Joseph et al., 1985; Keunings, 1986; and Renardy, 2000a) and has been a major stumbling block in computational rheology for several decades (see, e.g., the reviews Keunings, 2000; Keunings and Euler, 2001; Walters and Webster, 2003; and Owens and Phillips, 2002). It refers to the inability of numerical schemes to resolve the stress singularities and boundary layers, causing computations to break down at relatively low Weissenberg numbers. The observed critical Weissenberg number at which convergence fails varies with both the numerical method, and geometry and rheological model. However, a significant step forward was made through the stabilization approach of Fattal and Kupferman (2004), who suggested that the HWNP was primarily a numerical issue of stiffness in capturing steep spatial gradients and their subsequent convection. The numerical consequence of inaccurately approximating transported stress profiles is the local multiplicative growth of the errors, resulting in a numerical blowup (Fattal and Kupferman, 2005). The approach to alleviate this numerical stiffness was to introduce the matrix logarithm of the stress (Fattal and Kupferman, 2004). Further discussion on this approach and subsequently other stabilization approaches can be found in Afonso et al. (2011; 2012) and Niethammer et al. (2018). However, this approach does not totally resolve the HWNP issue, and in order to make further progress, the stress singularity and boundary layers need to be adequately addressed by any numerical scheme.

In this current work, we consider two realistic viscoelastic models of Giesekus (1982) and Phan-Thien–Tanner (PTT) (1977), which capture the shear thinning behavior and non-zero normal stress differences of many common polymeric fluids (Renardy, 2000b). The stress singularity and boundary layer structures have been theoretically determined for these models in Renardy (1997), Hagen and Renardy (1997), and Evans (2010a; 2010b). We consider these two models together, since their constitutive equations differ by the presence of a term quadratic in the elastic stress [when the PTT model is considered in its “simplified” form of upper convected derivative and linear stress function (Poole et al., 2019)]. However, the difference in this term leads to both different stress singularities and boundary layer structures for the two models. We also remark that these viscoelastic models are currently used in complex computational modeling (e.g., Minaeian et al., 2022; Bayat et al., 2022; Jeyasountharan et al., 2022; and Maqbool et al., 2022) evidencing the need of further study in the context of comparisons between theoretical and numerical results.

The goal of the work is to verify whether the theoretical stress singularity and boundary layer structures of the Giesekus and PTT models can be captured numerically in full numerical simulations. This would then confirm that the theoretical results are borne out in practical simulations and are the relevant structures to address in future work. We contrast two formulations of the constitutive equations. One in which the stress tensor is expressed in a fixed basis, which we term the Cartesian stress formulation (CSF), and another in a basis aligned with the velocity flow field, termed the natural stress formulation (NSF). We demonstrate that the NSF more accurately captures both the stress singularity and boundary layer structures, and moreover is not subject to a limiting Weissenberg number (for sufficiently large values of the solvent viscosity fraction). In addition, we have included a study concerning the dominance analysis of the constitutive equations around the reentrant corner. This investigation leverages further interpretations and constructions of data-driven methods (Callaham et al., 2021) for viscoelastic models in complex flows.

We remark that because the stresses and deformation rates at a reentrant corner are infinite, the local Weissenberg number at the corner should also be viewed as infinite (Renardy, 2000a), regardless of the value of the global Weissenberg number. This consequently justifies performing simulations at order one Weissenberg numbers, where the main features of interest will still be manifest. The layout of the paper is as follows: in Secs. II and III, we summarize the models and asymptotic results, and then, in Sec. IV, we present the numerical simulations confirming the theoretical asymptotic behaviors.

II. GOVERNING EQUATIONS

The dimensional governing equations in the present work are the mass and momentum equations given, respectively, by

\nabla \cdot v = 0, ρ (\frac{\partial v}{\partial t} + v \cdot \nabla v) = - \nabla p + \nabla \cdot τ,

(2.1)

where v is the fluid velocity, p is the pressure, ρ is the density, and $τ$ is the extra-stress tensor. The extra-stress is written as the sum of Newtonian solvent and elastic (polymer) contributions in the form

τ = 2 η_{s} D + τ^{p},

(2.2)

where η_s is the solvent viscosity, $D = \frac{1}{2} (\nabla v + {(\nabla v)}^{T})$ is the rate-of-strain tensor, and $τ^{p}$ the elastic extra-stress. As constitutive equations for the elastic extra-stress, we adopt the PTT (the affine and linear stress function versions) and Giesekus models in the form

τ^{p} + λ_{p} (\overset{▽}{τ^{p}} + \frac{κ}{η_{p}} g (τ^{p})) = 2 η_{p} D,

(2.3)

where λ_p is the relaxation time, $\overset{▽}{τ^{p}}$ is the upper-convected derivative of the elastic extra-stress, η_p is the polymeric viscosity, and κ is a model dependent parameter [often termed ϵ for PTT and influences the elongational behavior (Phan-Thien and Tanner, 1977), while it is termed the mobility factor α for Giesekus] controlling the influence of the quadratic stress terms

g (τ^{p}) = {\begin{matrix} (tr τ^{p}) τ^{p}, & PTT, \\ {(τ^{p})}^{2}, & Giesekus . \end{matrix}

(2.4)

We non-dimensionalize using a characteristic length L and flow speed U as follows:

\begin{matrix} x = L x^{*}, u = U u^{*}, t = \frac{L}{U} t^{*}, p = \frac{η_{0} U}{L} p^{*}, τ^{p} = \frac{η_{0} U}{L} T^{p *}, \end{matrix}

(2.5)

with $η_{0} = η_{s} + η_{p}$ being the total viscosity, which give (upon dropping ^*s) the dimensionless governing equations

\nabla \cdot v = 0,

(2.6)

Re (\frac{\partial v}{\partial t} + v \cdot \nabla v) = - \nabla p + β \nabla^{2} v + \nabla \cdot T^{p},

(2.7)

T^{p} + Wi (\overset{▽}{T^{p}} + \frac{κ}{(1 - β)} g (T^{p})) = 2 (1 - β) D,

(2.8)

g (T^{p}) = {\begin{matrix} tr (T^{p}) T^{p}, & PTT, \\ {(T^{p})}^{2}, & Giesekus, \end{matrix}

(2.9)

with dimensionless parameters of the Reynolds number $Re = ρ U L / η_{0}$ ⁠, the Weissenberg number $Wi = λ_{p} U / L$ (the dimensionless relaxation time), and retardation parameter $β = η_{s} / η_{0} \in [0, 1]$ (a dimensionless retardation time or solvent viscosity fraction, so that $η_{p} / η_{0} = 1 - β$ ⁠). Explicitly, the upper-convected stress derivative is defined as

\overset{▽}{T^{p}} = \frac{\partial T^{p}}{\partial t} + (v \cdot \nabla) T^{p} - (\nabla v) T^{p} - T^{p} {(\nabla v)}^{T} .

(2.10)

The system of Eqs. (2.6)–(2.8) is considered here for planar flows, for which there are two notable representations of the constitutive equations (2.8), when expressed in component form (see, e.g., Evans and Oishi, 2017; Evans et al. 2019b; 2019a; 2020). The first and most common gives the Cartesian stress formulation (CSF), where the elastic extra-stress tensor is expressed with respect to the usual fixed Cartesian orthonormal basis, so that the Cartesian stress components for (2.8) satisfy

Wi [\frac{\partial T_{11}^{p}}{\partial t} + u \frac{\partial T_{11}^{p}}{\partial x} + v \frac{\partial T_{11}^{p}}{\partial y} - 2 (T_{11}^{p} \frac{\partial u}{\partial x} + T_{12}^{p} \frac{\partial u}{\partial y}) + \frac{κ}{(1 - β)} g_{11}] + T_{11}^{p} - 2 (1 - β) \frac{\partial u}{\partial x} = 0,

(2.11)

Wi [\frac{\partial T_{12}^{p}}{\partial t} + u \frac{\partial T_{12}^{p}}{\partial x} + v \frac{\partial T_{12}^{p}}{\partial y} - \frac{\partial v}{\partial x} T_{11}^{p} - T_{22}^{p} \frac{\partial u}{\partial y} + \frac{κ}{(1 - β)} g_{12}] + T_{12}^{p} - (1 - β) (\frac{\partial u}{\partial y} + \frac{\partial v}{\partial x}) = 0,

(2.12)

Wi [\frac{\partial T_{22}^{p}}{\partial t} + u \frac{\partial T_{22}^{p}}{\partial x} + v \frac{\partial T_{22}^{p}}{\partial y} - 2 (T_{22}^{p} \frac{\partial v}{\partial y} + T_{12}^{p} \frac{\partial v}{\partial x}) + \frac{κ}{(1 - β)} g_{22}] + T_{22}^{p} - 2 (1 - β) \frac{\partial v}{\partial y} = 0,

(2.13)

where

\begin{matrix} g_{11} = {\begin{matrix} (T_{11}^{p} + T_{22}^{p}) T_{11}^{p}, & PTT, \\ {(T_{11}^{p})}^{2} + {(T_{12}^{p})}^{2}, & Giesekus, \end{matrix} g_{12} = {\begin{matrix} (T_{11}^{p} + T_{22}^{p}) T_{12}^{p}, & PTT, \\ (T_{11}^{p} + T_{22}^{p}) T_{12}^{p}, & Giesekus, \end{matrix} \\ g_{22} = {\begin{matrix} (T_{11}^{p} + T_{22}^{p}) T_{22}^{p}, & PTT, \\ {(T_{12}^{p})}^{2} + {(T_{22}^{p})}^{2}, & Giesekus . \end{matrix} \end{matrix}

(2.14)

An alternative stress formulation introduces the conformation tensor A through the transformation

T^{p} = \frac{(1 - β)}{W i} (A - I),

(2.15)

so that (2.8) becomes

Wi \overset{▽}{A} + (A - I) + κ g (A - I) = 0.

(2.16)

The conformation tensor may then be decomposed using the velocity field and an orthogonal vector to express the elastic part of the extra-stress as

A = λ v v^{T} + μ (v w^{T} + w v^{T}) + ν w w^{T},

(2.17)

with

v = {(u, v)}^{T}, w = \frac{1}{| v |^{2}} {(- v, u)}^{T} .

This is referred to as the natural stress formulation (NSF) and results in the components of (2.17) in (2.16) satisfying

Wi [\frac{\partial λ}{\partial t} + (v \cdot \nabla) λ + \frac{2 λ}{| v |} \frac{\partial | v |}{\partial t} + \frac{2 μ}{| v |^{4}} (v \frac{\partial u}{\partial t} - u \frac{\partial v}{\partial t}) + 2 μ \nabla \cdot w] + (λ - \frac{1}{| v |^{2}}) + κ g_{λ} = 0,

(2.18)

Wi [\frac{\partial μ}{\partial t} + (v \cdot \nabla) μ + (λ - \frac{ν}{| v |^{4}}) (u \frac{\partial v}{\partial t} - v \frac{\partial u}{\partial t}) + ν \nabla \cdot w] + μ + κ g_{μ} = 0,

(2.19)

Wi [\frac{\partial ν}{\partial t} + (v \cdot \nabla) ν - \frac{2 ν}{| v |} \frac{\partial | v |}{\partial t} + 2 μ (u \frac{\partial v}{\partial t} - v \frac{\partial u}{\partial t})] + (ν - | v |^{2}) + κ g_{ν} = 0,

(2.20)

where

\begin{matrix} g_{λ} = {\begin{matrix} (λ | v |^{2} - 2 + \frac{ν}{| v |^{2}}) (λ - \frac{1}{| v |^{2}}), & PTT, \\ {(λ - \frac{1}{| v |^{2}})}^{2} | v |^{2} + \frac{μ^{2}}{| v |^{2}}, & Giesekus, \end{matrix} \\ g_{μ} = (λ | v |^{2} - 2 + \frac{ν}{| v |^{2}}) μ, \\ g_{ν} = {\begin{matrix} (λ | v |^{2} - 2 + \frac{ν}{| v |^{2}}) (ν - | v |^{2}), & PTT, \\ {(ν - | v |^{2})}^{2} \frac{1}{| v |^{2}} + μ^{2} | v |^{2}, & Giesekus \end{matrix} \end{matrix}

(2.21)

and

\nabla \cdot w = \frac{1}{| v |^{4}} [(v^{2} - u^{2}) (\frac{\partial v}{\partial x} + \frac{\partial u}{\partial y}) + 4 u v \frac{\partial u}{\partial x}] .

It is notable that the coupling of the terms in the upper convective derivative part (square brackets) of the CSF equations (2.11)–(2.13) takes place through the spatial velocity gradients. In the equivalent terms in the NSF equations (2.18)–(2.20), the coupling takes place through the temporal velocity gradients. This emphasizes that the two formulations are significantly different and that the upper convected derivative of the natural stress variables achieves significant decoupling for steady flows.

III. REVISITING THE ASYMPTOTIC RESULTS

In this section, we present a summary of the asymptotic results for the reentrant corner viscoelastic flows of the PTT and Giesekus models (Renardy, 1997; Evans, 2010a; 2010b). They are presented for the specific case of the $270 °$ corner angle, which arises in the contraction geometry shown in Fig. 1(a). The numerical simulations presented later in Sec. IV are obtained for the specific benchmark contraction case of 4:1.

FIG. 1.

View large Download slide

(a) Contraction geometry and (b) zoomed-in corner region.

The key feature of the asymptotic stress behaviors at a reentrant corner for both PTT and Giesekus models is that the Newtonian solvent stress dominates the elastic polymer stress. The velocity field is thus Newtonian in character and allows the classical works of Dean and Montagnon (1949) and Moffatt (1964) to be used to construct its behavior local to the reentrant corner. Consequently, away from the walls, but close to the corner, in a region we shall call the core region (Hinch, 1993), we can write the stream function and pressure, respectively, as

ψ = C_{0} r^{n + 1} f_{0} (θ), p = C_{0} \frac{β}{(n - 1)} r^{n - 1} ({(n + 1)}^{2} f_{0}^{'} (θ) + f_{0}^{″'} (θ)),

(3.1)

where n = 0.5445 and $(r, θ)$ are the polar coordinates centered at the reentrant corner and

f_{0} (θ) = \frac{A_{0}}{C_{0}} cos ((n + 1) θ) + \frac{B_{0}}{C_{0}} sin ((n + 1) θ) + cos ((n - 1) θ) + \frac{D_{0}}{C_{0}} sin ((n - 1) θ),

(3.2)

with

\frac{A_{0}}{C_{0}} = 0.2030 \dots, \frac{B_{0}}{C_{0}} = 0.5430 \dots, \frac{D_{0}}{C_{0}} = - 0.3738 \dots,

(3.3)

and C₀ is a constant, which is determined by the properties of the flow away from the reentrant corner. In these expressions, we have used the numerical value of the first eigenvalue that controls the dominant behavior in the separable solution (3.1). This flow field gives the order of magnitude estimates determined in Evans (2010a; 2010b) and summarized in Table I for the main flow variables of velocity, velocity gradient (and pressure), elastic (polymer) stress, and natural stress variables.

TABLE I.

Asymptotic behavior of flow properties as $r \to 0$ ⁠.

Model	v	$\nabla v$ (and p)	$T^{p}$	B.L. thickness	λ	μ	ν
PTT	$r^{0.5445}$	$r^{- 0.4555}$	$r^{- 0.3286}$	$r^{1.1518}$	$r^{- 1.4176}$	$r^{- 0.1268}$	$r^{1.1638}$
Giesekus	$r^{0.5445}$	$r^{- 0.4555}$	$r^{- 0.2796}$	$r^{1.2278}$	$r^{- 1.3686}$	$r^{0.0}$	$r^{1.3686}$

Model	v	$\nabla v$ (and p)	$T^{p}$	B.L. thickness	λ	μ	ν
PTT	$r^{0.5445}$	$r^{- 0.4555}$	$r^{- 0.3286}$	$r^{1.1518}$	$r^{- 1.4176}$	$r^{- 0.1268}$	$r^{1.1638}$
Giesekus	$r^{0.5445}$	$r^{- 0.4555}$	$r^{- 0.2796}$	$r^{1.2278}$	$r^{- 1.3686}$	$r^{0.0}$	$r^{1.3686}$

The interesting features of the viscoelastic models are that the results in Table I for the elastic stresses do not hold uniformly around the corner. Rather, close to the walls viscometric elastic stresses dominate that are appropriate to shearing behavior. Consequently, stress boundary layers arise that reconcile the two behaviors. These stress boundary layers are the same that occur in the high Weissenberg limit discussed in Hagen and Renardy (1997), which appear to manifest themselves at stress singularities even in Weissenberg order one flows. The asymptotic estimates of the boundary layer (B.L.) widths are given in Table I. Figure 2 illustrates these boundary layer widths for the PTT and Giesekus models, along with three selected streamlines that pass close to the singularity. The PTT model is noted to have a slightly wider boundary layer than the Giesekus model, since its radial exponent for the B.L. thickness is smaller than that for Giesekus. A similar behavior in boundary layer thickness is seen in the stick-slip problem (Evans et al., 2017).

FIG. 2.

View large Download slide

A plot of selected streamlines using (3.1) near the singularity, along with the boundary layer curves for all two models at both the upstream and downstream boundary layers.

The actual boundary layer equations for both models and formulations are now recorded below. These sets of equations assume a set of local Cartesian axes, orientated with x along the solid boundary and y normal into the fluid.

Boundary layer equations for PTT-CSF are as follows:

u \frac{\partial T_{11}^{p}}{\partial x} + v \frac{\partial T_{11}^{p}}{\partial y} - 2 (T_{11}^{p} \frac{\partial u}{\partial x} + T_{12}^{p} \frac{\partial u}{\partial y}) + \frac{κ}{(1 - β)} {(T_{11}^{p})}^{2} = 0,

(3.4)

u \frac{\partial T_{12}^{p}}{\partial x} + v \frac{\partial T_{12}^{p}}{\partial y} - \frac{\partial v}{\partial x} T_{11}^{p} - T_{22}^{p} \frac{\partial u}{\partial y} + \frac{κ}{(1 - β)} T_{12}^{p} T_{11}^{p} = \frac{(1 - β)}{Wi} \frac{\partial u}{\partial y},

(3.5)

u \frac{\partial T_{22}^{p}}{\partial x} + v \frac{\partial T_{22}^{p}}{\partial y} - 2 (T_{22}^{p} \frac{\partial v}{\partial y} + T_{12}^{p} \frac{\partial v}{\partial x}) + \frac{κ}{(1 - β)} T_{11}^{p} T_{22}^{p} = 2 \frac{(1 - β)}{Wi} \frac{\partial v}{\partial y} .

(3.6)

Boundary layer equations for PTT-NSF are as follows:

(v . \nabla) λ + 2 μ \frac{\partial (1 / u)}{\partial y} + \frac{κ}{Wi} λ^{2} u^{2} = 0,

(3.7)

(v . \nabla) μ + ν \frac{\partial (1 / u)}{\partial y} + \frac{κ}{Wi} λ μ u^{2} = 0,

(3.8)

(v . \nabla) ν + \frac{κ}{Wi} λ u^{2} (ν - u^{2}) = 0.

(3.9)

Boundary layer equations for Giesekus–CSF are as follows:

u \frac{\partial T_{11}^{p}}{\partial x} + v \frac{\partial T_{11}^{p}}{\partial y} - 2 (T_{11}^{p} \frac{\partial u}{\partial x} + T_{12}^{p} \frac{\partial u}{\partial y}) + \frac{κ}{(1 - β)} {(T_{11}^{p})}^{2} = 0,

(3.10)

u \frac{\partial T_{12}^{p}}{\partial x} + v \frac{\partial T_{12}^{p}}{\partial y} - \frac{\partial v}{\partial x} T_{11}^{p} - T_{22}^{p} \frac{\partial u}{\partial y} + \frac{κ}{(1 - β)} T_{12}^{p} T_{11}^{p} = \frac{(1 - β)}{Wi} \frac{\partial u}{\partial y},

(3.11)

u \frac{\partial T_{22}^{p}}{\partial x} + v \frac{\partial T_{22}^{p}}{\partial y} - 2 (T_{22}^{p} \frac{\partial v}{\partial y} + T_{12}^{p} \frac{\partial v}{\partial x}) + \frac{κ}{(1 - β)} ({(T_{12}^{p})}^{2} + \frac{{(1 - β)}^{2}}{{Wi}^{2}}) - \frac{(1 - β)}{{Wi}^{2}} = 2 \frac{(1 - β)}{Wi} \frac{\partial v}{\partial y} .

(3.12)

We remark that Eqs. (3.11) and (3.12) are more commonly written using ${\bar{T}}_{22}^{p} = T_{22}^{p} + \frac{(1 - β)}{W i}$ and then take the more compact form

u \frac{\partial T_{12}^{p}}{\partial x} + v \frac{\partial T_{12}^{p}}{\partial y} - \frac{\partial v}{\partial x} T_{11}^{p} - {\bar{T}}_{22}^{p} \frac{\partial u}{\partial y} + \frac{κ}{(1 - β)} T_{12}^{p} T_{11}^{p} = 0,

(3.13)

u \frac{\partial {\bar{T}}_{22}^{p}}{\partial x} + v \frac{\partial {\bar{T}}_{22}^{p}}{\partial y} - 2 ({\bar{T}}_{22}^{p} \frac{\partial v}{\partial y} + T_{12}^{p} \frac{\partial v}{\partial x}) + \frac{κ}{(1 - β)} {(T_{12}^{p})}^{2} = \frac{(1 - β)}{{Wi}^{2}} (1 - κ) .

(3.14)

Boundary layer equations for Giesekus–NSF are as follows:

(v . \nabla) λ + 2 μ \frac{\partial (1 / u)}{\partial y} + \frac{κ}{Wi} λ^{2} u^{2} = 0,

(3.15)

(v . \nabla) μ + ν \frac{\partial (1 / u)}{\partial y} + \frac{κ}{Wi} λ μ u^{2} = 0,

(3.16)

(v . \nabla) ν + \frac{κ}{Wi} (μ^{2} - \frac{(1 - κ)}{κ}) u^{2} = 0.

(3.17)

In both models, the elastic stresses of the CSF and NSF are linked in the boundary layers through the equations

T_{11} = \frac{(1 - β)}{Wi} λ u^{2}, T_{12} = \frac{(1 - β)}{Wi} (λ u v + μ), T_{22} = \frac{(1 - β)}{Wi} (λ v^{2} + 2 μ \frac{v}{u} + \frac{ν}{u^{2}} - 1) .

(3.18)

A convenient summary of the terms that dominate in both the core and boundary layer regions, for both models and both formulations, is given in Table II. These are the leading-order terms that hold in the respective regions and it should be noted that not all terms in the quadratic stress expressions and the rate-of-strain term are present, as recorded in the explicit boundary layer Eqs. (3.4)–(3.17). The subscript on the upper convected derivative refers to the appropriate component and explicitly given in Eqs. (3.4)–(3.6) or (3.10)–(3.12). Moreover, it is noted for the Giesekus model that the constant terms in (3.12) arise from the leading order behavior $T_{22}^{p} \sim - (1 - β) / Wi$ in the linear and quadratic stress terms. Purely for convenience, we absorb the additional constant term from the linear stress component into the g₂₂ quadratic stress term for the purposes of Table II. Such behavior of the T₂₂ component in the boundary layer is peculiar to the Giesekus model and does not occur for PTT. These results hold for small or order one Reynolds number, order one (or large) Weissenberg numbers and $0 < β < 1$ ⁠. We would expect these results to breakdown as the solvent viscosity vanishes. The case β = 0 thus needs separate consideration and for which currently no complete asymptotic results are known (Evans and Sibley, 2008; 2009).

TABLE II.

Boundary layer equations according to the stress formulation, viscoelastic model, and regions.

Formulation	Downstream and upstream BL	Core
CSF	$Wi [\overset{▽}{T^{p}_{11}} + \frac{κ}{1 - β} g_{11}] = 0$
	$Wi [\overset{▽}{T^{p}_{12}} + \frac{κ}{1 - β} g_{12}] = 2 (1 - β) D_{12}$	$\overset{▽}{T^{p}} = 0$
	$Wi [\overset{▽}{T^{p}_{22}} + \frac{κ}{1 - β} g_{22}] = 2 (1 - β) D_{22}$
$g_{11} \sim {(T_{22}^{p})}^{2}, g_{12} \sim T_{11}^{p} T_{12}^{p}, g_{22} \sim {\begin{matrix} T_{11}^{p} T_{22}^{p}, & PTT \\ {(T_{12}^{p})}^{2} + \frac{{(1 - β)}^{2}}{{Wi}^{2}} \frac{(κ - 1)}{κ} & Giesekus \end{matrix}$
NSF	$Wi [(v \cdot \nabla) λ + 2 μ \nabla \cdot w] + κ g_{λ} = 0$	$(v \cdot \nabla) λ = 0$
	$Wi [(v \cdot \nabla) μ + ν \nabla \cdot w] + κ g_{μ} = 0$	$(v \cdot \nabla) μ = 0$
	$Wi [(v \cdot \nabla) ν] + κ g_{ν} = {\begin{matrix} 0, & PTT \\ \| v \|^{2}, & Giesekus \end{matrix}$	$(v \cdot \nabla) ν = 0$
$g_{λ} \sim λ^{2} \| v \|^{2}, g_{μ} \sim λ_{μ} \| v \|^{2}, g_{ν} \sim {\begin{matrix} λ \| v \|^{2} (ν - \| v \|^{2}), & PTT \\ (μ^{2} + \frac{κ - 1}{κ}) \| v \|^{2}, & Giesekus \end{matrix}$

Formulation	Downstream and upstream BL	Core
CSF	$Wi [\overset{▽}{T^{p}_{11}} + \frac{κ}{1 - β} g_{11}] = 0$
	$Wi [\overset{▽}{T^{p}_{12}} + \frac{κ}{1 - β} g_{12}] = 2 (1 - β) D_{12}$	$\overset{▽}{T^{p}} = 0$
	$Wi [\overset{▽}{T^{p}_{22}} + \frac{κ}{1 - β} g_{22}] = 2 (1 - β) D_{22}$
$g_{11} \sim {(T_{22}^{p})}^{2}, g_{12} \sim T_{11}^{p} T_{12}^{p}, g_{22} \sim {\begin{matrix} T_{11}^{p} T_{22}^{p}, & PTT \\ {(T_{12}^{p})}^{2} + \frac{{(1 - β)}^{2}}{{Wi}^{2}} \frac{(κ - 1)}{κ} & Giesekus \end{matrix}$
NSF	$Wi [(v \cdot \nabla) λ + 2 μ \nabla \cdot w] + κ g_{λ} = 0$	$(v \cdot \nabla) λ = 0$
	$Wi [(v \cdot \nabla) μ + ν \nabla \cdot w] + κ g_{μ} = 0$	$(v \cdot \nabla) μ = 0$
	$Wi [(v \cdot \nabla) ν] + κ g_{ν} = {\begin{matrix} 0, & PTT \\ \| v \|^{2}, & Giesekus \end{matrix}$	$(v \cdot \nabla) ν = 0$
$g_{λ} \sim λ^{2} \| v \|^{2}, g_{μ} \sim λ_{μ} \| v \|^{2}, g_{ν} \sim {\begin{matrix} λ \| v \|^{2} (ν - \| v \|^{2}), & PTT \\ (μ^{2} + \frac{κ - 1}{κ}) \| v \|^{2}, & Giesekus \end{matrix}$

IV. UNSTEADY PLANAR 4:1 CONTRACTION FLOW

We now use the full numerical simulation of the mass and momentum conservation equations (2.6) and (2.7) in combination with the constitutive equations for CSF Eqs. (2.11)–(2.13) or NSF Eqs. (2.18)–(2.20). The numerical methodology for solving the set of equations is based on our previous work (Evans et al., 2019b; 2019a; 2020). In summary, a Cartesian non-uniform mesh is employed to discretize the domain of the contraction flow geometry, with a more refined mesh close to the corners, as described in Fig. 3. The equations are discretized using finite-difference approximations, while a projection scheme is employed in order to uncouple the velocity and pressure fields in Eqs. (2.6) and (2.7). Explicit time discretizations are used for solving the constitutive equations of CSF and NSF, while the convective terms presented in these equations are approximated by an upwind methodology.

FIG. 3.

View large Download slide

Contraction geometry and zoomed-in interest regions.

In the simulations, we have fixed some parameters as $κ = 0.25$ and $Re = 10^{- 2}$ ⁠. The other parameters β and $Wi$ vary, as described in this section. Two different mesh sizes were used: M1 and M2 with the minimum size and channel lengths given in Table III.

TABLE III.

Meshes details.

Mesh	$Δ_{x_{min}} = Δ_{y_{min}}$	Channel lengths
M1	$3.5 \times 10^{- 3}$	$40 L \times 100 L$
M2	$1 \times 10^{- 4}$	$40 L \times 40 L$

A. Numerical verification: Corner vortex analysis and Couette correction

Preliminary to investigating the asymptotic behaviors at the sharp corners, it is necessary to confirm first that we have complete attached flow at such corners. The asymptotic results presented in Sec. III are premised on the absence of a separating streamline at the reentrant corner, either due to the presence of a lip vortex or the extension of the salient corner vortex. This is confirmed for both models in Figs. 4 and 5, for a range of Weissenberg values with $β = 1 / 2$ and 1/9, and fixed $κ = 0.25$ ⁠. No presence of lip vortices was detected (in keeping with the numerical results for 4:1 contractions in Alves et al. (2003; 2004), and it is clear that the larger value $β = 1 / 2$ confines the salient corner vortex more than the traditional benchmark value of $β = 1 / 9$ ⁠. Further, it is noteworthy that our simulations have no upper Weissenberg limit in the NSF formulation for either β value. However, the restriction of the salient corner vortex and time required to reach steady-state makes the $β = 1 / 2$ solvent viscosity fraction case optimum to investigate the reentrant corner behavior and which we exclusively use in Secs. IV B and IV C.

FIG. 4.

View large Download slide

Comparison of the vortex size using $β = \frac{1}{2}$ (first column) and $β = \frac{1}{9}$ (second column) and mesh M₁ for the PTT model and CSF. (a) $Wi = 0.5$ ⁠; (b) $Wi = 0.5$ ⁠; (c) $Wi = 2$ ⁠; (d) $Wi = 2$ ⁠; (e) $Wi = 10$ ⁠; (f) $Wi = 10$ ⁠.

FIG. 5.

View large Download slide

Streamlines for the planar contraction flow using the PTT (first column) and Giesekus (second column) viscoelastic models and the NSF with mesh M₂ and $β = \frac{1}{2}$ ⁠. (a) $Wi = 1$ ⁠; (b) $Wi = 1$ ⁠; (c) $Wi = 10$ ⁠; (d) $Wi = 10$ ⁠; (e) $Wi = 100$ ⁠; (f) $Wi = 100$ ⁠; (g) $Wi = 1000$ ⁠; (h) $Wi = 1000$ ⁠; (i) $Wi = 10 000$ ⁠; (j) $Wi = 10 000$ ⁠.

Our first simulation results are listed in Table IV and compare against the benchmark data of Alves et al. (2003) for linear PTT (no benchmark data being available for the Giesekus model). Recorded are the size of the corner vortex X_R and the Couette correction coefficient C for the pressure, both of which are defined in Alves et al. (2003). For consistency, we use the same minimum mesh spacing as in the benchmark case (Alves et al., 2003), this relatively coarse mesh being denoted as M₁, with the same model parameter values $β = 1 / 9$ and $κ = 0.25$ ⁠. The CSF simulation results are in good agreement with the benchmark data for Weissenberg values up to 100, the simulations then ceasing to converge at higher Weissenberg numbers. This breakdown signifies the classical HWNP, which all numerical schemes using the CSF have historically encountered.

TABLE IV.

Benchmark results for the linear PTT with $κ = 0.25$ and $β = \frac{1}{9}$ and mesh M₁.

$Wi$	X_R (CSF)	X_R (NSF)	X_R (Alves et al., 2003)	C (CSF)	C (NSF)	C (Alves et al., 2003)
0	1.4922	1.4922	1.5002	0.3753	0.3753	0.3741
0.5	1.4922	⋯	1.5060	0.1696	⋯	0.1672
1	1.5309	⋯	1.5420	0.0980	⋯	0.0951
2	1.6112	⋯	1.6390	0.0284	⋯	0.0261
10	2.1799	⋯	2.1310	−0.1543	⋯	−0.1113
50	2.7240	1.4173	2.4930	−0.1653	−1.0352	−0.1444
100	2.7240	1.4921	2.570	−0.1168	−0.7405	−0.0917
1000	⋯	1.4543	⋯	⋯	−0.1847	⋯

$Wi$	X_R (CSF)	X_R (NSF)	X_R (Alves et al., 2003)	C (CSF)	C (NSF)	C (Alves et al., 2003)
0	1.4922	1.4922	1.5002	0.3753	0.3753	0.3741
0.5	1.4922	⋯	1.5060	0.1696	⋯	0.1672
1	1.5309	⋯	1.5420	0.0980	⋯	0.0951
2	1.6112	⋯	1.6390	0.0284	⋯	0.0261
10	2.1799	⋯	2.1310	−0.1543	⋯	−0.1113
50	2.7240	1.4173	2.4930	−0.1653	−1.0352	−0.1444
100	2.7240	1.4921	2.570	−0.1168	−0.7405	−0.0917
1000	⋯	1.4543	⋯	⋯	−0.1847	⋯

Focusing upon $β = 1 / 2$ ⁠, Table V gives the vortex size and Couette correction for both PTT and Giesekus. The CSF and NSF results are comparable for Weissenberg numbers for which both schemes converged. However, an upper Weissenberg limit was encountered in the CSF, which the NSF did not encounter for either model. In fact on the finest mesh M₂, relatively low Weissenberg numbers only were obtainable with the CSF, with the Giesekus model being more restrictive.

TABLE V.

Estimates of the vortex size X_R and the Couette correction C for both viscoelastic models and formulations using the mesh M₂ and $β = \frac{1}{2}$ ⁠.

	PTT				Giesekus
Wi	CSF		NSF		CSF		NSF
	X_R	C	X_R	C	X_R	C	X_R	C
1	1.5059	0.2306	1.5059	0.2282	1.5059	0.2668	1.5059	0.2646
10	1.5682	0.2305	1.5682	0.2285	⋯	⋯	1.6330	0.3275
100	⋯	⋯	1.5059	0.2662	⋯	⋯	1.5059	0.3493
1000	⋯	⋯	1.5059	0.3599	⋯	⋯	1.5059	0.3641
10 000	⋯	⋯	1.5059	0.3640	⋯	⋯	1.5059	0.3648

	PTT				Giesekus
Wi	CSF		NSF		CSF		NSF
	X_R	C	X_R	C	X_R	C	X_R	C
1	1.5059	0.2306	1.5059	0.2282	1.5059	0.2668	1.5059	0.2646
10	1.5682	0.2305	1.5682	0.2285	⋯	⋯	1.6330	0.3275
100	⋯	⋯	1.5059	0.2662	⋯	⋯	1.5059	0.3493
1000	⋯	⋯	1.5059	0.3599	⋯	⋯	1.5059	0.3641
10 000	⋯	⋯	1.5059	0.3640	⋯	⋯	1.5059	0.3648

Figure 6 records the behavior of the velocity components, pressure, and Cartesian extra-stress components with radial distance from the reentrant corner along the ray $θ = π / 2$ ⁠. In each plot, the two cases $β = 1 / 9$ and $β = 1 / 2$ are presented for comparison, and this is done for three selected order one Weissenberg numbers with $κ = 0.25$ kept fixed. The theoretical value of the radial exponent for these variables is stated in each plot. The velocity components are a reasonable fit over the range of Weissenberg numbers, as are the Cartesian stress components (although these appear deteriorate as Weissenberg increases). However, the pressure behaviors do not capture the expected theoretical asymptotics.

FIG. 6.

View large Download slide

Comparison of the slopes using $β = \frac{1}{2}$ and $β = \frac{1}{9}$ and mesh M₁ for PTT along the ray $θ = π / 2$ ⁠. (a) u and v velocity components with $Wi = 0.5$ ⁠; (b) u and v velocity components with $Wi = 2$ ⁠; (c) u and v velocity components with $Wi = 10$ ⁠; (d) pressure p with $Wi = 0.5$ ⁠; (e) pressure p with $Wi = 2$ ⁠; (f) pressure p with $Wi = 10$ ⁠; (g) elastic extra-stress components with $Wi = 0.5$ ⁠; (h) elastic extra-stress components with $Wi = 2$ ⁠; (i) elastic extra-stress components with $Wi = 10$ ⁠.

B. Numerical verification: Sharp corner asymptotics for $W i = O (1)$

In this section, we present results for the parameter case $W i = 1, β = 1 / 2, κ = 0.25$ for both models. Figures 7–9 give the limiting behaviors of the velocity, pressure, Cartesian and natural stress components at the reentrant corner. Shown are the radial behaviors along selected fixed angles, uniformly spaced by $π / 4$ around the $3 π / 2$ corner, together with the theoretical slopes of the variables from Table I. Figure 7 confirms the asymptotic Newtonian behavior of the velocity field for both CSF and NSF, with the NSF giving better convergence behavior. Less convincing is the pressure behavior, which we note is better for NSF, but still seems to require a far finer mesh for confirmation. Figure 8 gives the Cartesian stress components using the CSF, with the results indicating slower convergence and loss of resolution at the end points, than compared to the natural stress components in Fig. 9. Notable for the CSF in Fig. 8 is that the convergence rates to the limiting behaviors vary significantly with the angle of approach, with both models generally behaving similarly apart from the T₂₂ component along the downstream ray $θ = π / 4$ ⁠. It is also noteworthy that the convergence of the natural stress components is significantly less affected by the ray angle of approach and showing good convergence properties uniformly around the reentrant corner.

FIG. 7.

View large Download slide

Numerical behavior of the velocity components and pressure along different radial distances $θ = \frac{π}{4}, \frac{π}{2}, \frac{3 π}{4}$ ⁠, π, and $\frac{5 π}{4}$ with $Wi = 1$ ⁠. (a) PTT: u velocity component; (b) Giesekus: u velocity component; (c) PTT: v velocity component; (d) Giesekus: v velocity component; (e) PTT: pressure p; (f) Giesekus: pressure p.

FIG. 8.

View large Download slide

Numerical behavior of the stress components along different radial distances $θ = \frac{π}{4}, \frac{π}{2}, \frac{3 π}{4}$ ⁠, π, and $\frac{5 π}{4}$ with $Wi = 1$ ⁠. (a) PTT: $T_{11}^{p}$ elastic extra-stress component; (b) Giesekus: $T_{11}^{p}$ elastic extra-stress component; (c) PTT: $T_{12}^{p}$ elastic extra-stress component; (d) Giesekus: $T_{12}^{p}$ elastic extra-stress component; (e) PTT: $T_{22}^{p}$ elastic extra-stress component; (f) Giesekus: $T_{22}^{p}$ elastic extra-stress component.

FIG. 9.

View large Download slide

Numerical behavior of the natural stress components along different radial distances $θ = \frac{π}{4}, \frac{π}{2}, \frac{3 π}{4}$ ⁠, π, and $\frac{5 π}{4}$ with $Wi = 1$ ⁠. (a) PTT: λ component; (b) Giesekus: λ component; (c) PTT: μ component; (d) Giesekus: μ component; (e) PTT: ν component; (f) Giesekus: ν component.

Next, we examine the dominant terms in the constitutive equations close to the reentrant corner. To enable this, we record in Table VI the radial sizes of the component terms occurring in both the CSF and NSF constitutive equations. These may be deduced directly from Table I and help order the sizes of terms in the constitutive equations in the asymptotic limit $r \to 0$ as the radial distance from the corner vanishes. First, we discuss the CSF for both models and present in Fig. 10, for each component equation in (2.11)–(2.13) the term whose modulus is a maximum. Table VI records the asymptotic sizes (with respect to the radial distance from the corner) of the terms in the CSF equations, which should hold in the core region near the corner but away from the walls. The plots do confirm the dominance of the upper convected derivative of stress for each component. However, Table VI indicates that the quadratic stress terms followed by the rate-of-strain terms are not that much smaller, and their appearance in Fig. 10 is not surprising. To see the relative ordering of the terms, we plot in Fig. 11 the modulus of the individual groups of terms in the constitutive equations along a given streamline that passes close to the reentrant corner (and shown in the plots in Fig. 10). The orientation of the polar coordinates is such that θ = 0 is the downstream wall and $θ = 3 π / 2 \approx 4.7$ is the upstream wall. The intention is to determine whether the theoretical asymptotic balances summarized in Table II hold for the boundary layers at the walls and the core region away from the walls. Figure 11 gives the sizes of the groups of terms in (2.11)–(2.13), plotted by component for the CSF. Away from the walls, it is clear that the upper convective derivative dominates in the 11- and 12-component equations, but that the quadratic stress terms are comparable in the 22-component equation. Near the upstream wall (⁠ $θ \approx 4.7$ ⁠), the linear stress terms do appear smaller than the quadratic stress terms in the 11- and 12-component equations, and remain so along the rest of the streamline. The 22-component equation suggests near the upstream wall the linear stress term is comparable with the quadratic stress term, which is consistent with the balance for this component in the Giesekus equation [see Table II and Eq. (3.12)]. Figure 12 gives the behavior of the terms along the same streamline but now parameterized with the radial distance from the corner and starting from near the upstream wall. It confirms and reinforces the previous remarks.

TABLE VI.

Asymptotic sizes of the terms in the steady CSF and NSF equations.

CSF

\overset{▽}{T^{p}}

\frac{κ}{(1 - β)} g (T^{p})

T^{p} / Wi

2 (1 - β) D / Wi

PTT

r^{- 0.7841}

r^{- 0.6572}

r^{- 0.3286}

r^{- 0.4555}

Giesekus

r^{- 0.7351}

r^{- 0.5592}

r^{- 0.2796}

r^{- 0.4555}

NSF

v \cdot \nabla λ

2 μ \nabla \cdot w

\frac{κ}{Wi} g_{λ}

\frac{1}{Wi} (λ - \frac{1}{| v |^{2}})

PTT

r^{- 1.8743}

r^{- 1.6713}

r^{- 1.7462}

r^{- 1.4176}

Giesekus

r^{- 1.8241}

r^{- 1.5445}

r^{- 1.6482}

r^{- 1.3686}

NSF

v \cdot \nabla μ

ν \nabla \cdot w

\frac{κ}{Wi} g_{μ}

\frac{1}{Wi} μ

PTT

r^{- 0.5823}

r^{- 0.3807}

r^{- 0.4554}

r^{- 0.1268}

Giesekus

r^{- 0.4555}

r^{- 0.1759}

r^{- 0.2796}

r⁰

NSF

v \cdot \nabla ν

⋯

\frac{κ}{Wi} g_{ν}

\frac{1}{Wi} (ν - | v |^{2})

PTT

r^{0.7083}

r^{0.7604}

r^{1.089}

Giesekus

r^{0.9131}

r^{1.089}

r^{1.089}

FIG. 10.

View large Download slide

Dominant terms using the CSF with $Wi = 1$ ⁠. (a) PTT: 11 component; (b) Giesekus: 11 component; (c) PTT: 12 component; (d) Giesekus: 12 component; (e) PTT: 22 component; (f) Giesekus: 22 component.

FIG. 11.

View large Download slide

Profile of the dominant terms in the constitutive equations (CSF) along the streamline $3.6 \times 10^{- 5}$ with $Wi = 1$ ⁠. (a) PTT: 11 component; (b) Giesekus: 11 component; (c) PTT: 12 component; (d) Giesekus: 12 component; (e) PTT: 22 component; (f) Giesekus: 22 component.

FIG. 12.

View large Download slide

Profile of the dominant terms in the Cartesian stress equations along the streamline $3.6 \times 10^{- 5}$ parameterized with the radial distance r with $Wi = 1$ ⁠. (a) PTT: 11 component; (b) Giesekus: 11 component; (c) PTT: 12 component; (d) Giesekus: 12 component; (e) PTT: 22 component; (f) Giesekus: 22 component.

Figures 13–16 repeat the previous CSF analysis but for NSF. The color map in Fig. 13 shows the dominant terms in the component equations (2.18)–(2.20). The convective term of the NS variable λ does dominate for both models in a region close to the corner. However, the quadratic stress terms appear to intrude very close to the corner, and this is confirmed when we follow the sizes of terms along the specified streamline close to the corner in Figs. 14 for PTT and 15 for Giesekus. In these two figures, results using the full numerical simulation code are denoted by “Full.” For comparison, we solve the NSF equations along the theoretical streamline given in (3.1) and present the same groups of terms in the plots denoted “Simplified” as a consistency check. Table VI gives the asymptotic sizes of the respective terms, and although the convective terms of the natural stresses dominate, it is seen that the quadratic stress terms are only slightly smaller. This would suggest that a much finer mesh is needed before the convective NS terms clearly dominate. It is clear that the streamline plots for PTT are supporting that the first three terms in the λ and μ NS equations and the first two terms in the ν NS equation in Table VI are balancing, with notably the terms $(λ - 1 / | v |^{2}), μ$ and $(ν - | v |^{2})$ being consistently smaller. We note that this is also the expected balance in the wall boundary layers as recorded in Table II. A similar behavior is seen in Fig. 15 for the Giesekus model, with notably the $(ν - | v |^{2})$ term now no longer being small and consistent with the predicted balance in Table VI (again also being the boundary layer balance in Table II). Nevertheless, despite the unclear dominance of all the convective NS variables at the reentrant corner, their behavior is still sufficiently close to the required theoretical asymptotic behaviors as evidenced by the NS variables in Fig. 9. Noticeably in Fig. 9, the μ variable has more difficulty converging for the Giesekus model, which may be understandable given the lack of dominance of its convective derivative in Figs. 13(d) and confirmed by Figs. 15(c) and 16(d). We note though that the behavior along the theoretical streamline in Fig. 15(d) would suggest the convected μ term should dominate strongly, yet its real behavior in Fig. 15(c) appears subdued. In a similar way of Fig. 12, we have described the behavior of the terms along the same streamline in Fig. 16 in order to further clear the above discussion.

FIG. 13.

View large Download slide

Dominant terms using the NSF with $Wi = 1$ ⁠. (a) PTT: λ component; (b) Giesekus: λ component; (c) PTT: μ component; (d) Giesekus: μ component; (e) PTT: ν component; (f) Giesekus: ν component.

FIG. 14.

View large Download slide

Profile of the dominant terms in the constitutive equations (NSF) of the PTT fluid along the streamline $3.6 \times 10^{- 5}$ with $Wi = 1$ ⁠. (a) Full: λ component; (b) simplified: λ component; (c) full: μ component; (d) simplified: μ component; (e) full: ν component; (f) simplified: ν component.

FIG. 15.

View large Download slide

Profile of the dominant terms in the constitutive equations (NSF) of the Giesekus fluid along the streamline $3.6 \times 10^{- 5}$ with $Wi = 1$ ⁠. (a) Full: λ component; (b) simplified: λ component; (c) full: μ component; (d) simplified: μ component; (e) full: ν component; (f) simplified: ν component.

FIG. 16.

View large Download slide

Profile of the dominant terms in the NS equations along the streamline $3.6 \times 10^{- 5}$ parameterized with the radial distance r with $Wi = 1$ ⁠. (a) PTT: λ component; (b) Giesekus: λ component; (c) PTT: μ component; (d) Giesekus: μ component; (e) PTT: ν component; (f) Giesekus: ν component.

C. Numerical study at high Weissenberg numbers

In this section, we explore the effect of increasing the Weissenberg number for both models. The parameter values $β = 1 / 2, κ = 0.25$ are kept fixed, while the Weissenberg number is increased from unity to 10⁴. The NSF was used exclusively for the simulations, due to the CSF breaking down for Weissenberg around 10². Behavior of the velocity components and pressure are given in Fig. 17 and confirm the Newtonian velocity field over the range of Weissenberg numbers. However, the pressure behavior is not confirmed. Convergence for the natural stress variables is confirmed in Fig. 18 over the range of Weissenberg numbers, and it is evident that convergence to the limiting behavior slows as the Weissenberg number increases.

FIG. 17.

View large Download slide

Asymptotic behavior of the velocity components and pressure for different values of $Wi = 1, 10, 100, 1000, 10 000$ ⁠. (a) PTT: u velocity component; (b) Giesekus: u velocity component; (c) PTT: v velocity component; (d) Giesekus: v velocity component; (e) PTT: pressure p; (f) Giesekus: pressure p.

FIG. 18.

View large Download slide

Asymptotic behavior of the natural stress variables for different values of $Wi = 1, 10, 100, 1000, 10 000$ ⁠. (a) PTT: λ component; (b) Giesekus: λ component; (c) PTT: μ component; (d) Giesekus: μ component; (e) PTT: ν component; (f) Giesekus: ν component.

Figures 19 and 20 show the behavior of the groups of terms in the NSF equations for large Weissenberg numbers. The balances in Table II again hold along the chosen streamline $ψ = 3.6 \times 10^{- 5}$ passing close to the reentrant corner. It is clear that the sizes of the terms grow extremely large at the walls, where there is now theoretically the occurrence of high Weissenberg boundary layers (Hagen and Renardy, 1997) that are present along all the walls of the domain. The thicknesses of these boundary layers are $W i^{- 1 / 3}$ for PTT and slightly thinner $W i^{- 1 / 2}$ for Giesekus. It is noticeable for both models that the behavior along the streamline is very different for $W i = 10^{3}$ and $W i = 10^{4}$ from $W i = 10^{2}$ ⁠, where the terms involving $\nabla \cdot w$ now dominate the quadratic stress in the λ and μ component equations. This is a direct consequence of the Weissenberg number being high. The transition to this behavior occurs between 10² and 10³. Figure 21 illustrates that the dominance of the convected natural stress variables is more evident than the $Wi = 1$ case, where small regions near the corner in the Giesekus case [panels (d) and (f) in Fig. 21] have started to emerge.

FIG. 19.

View large Download slide

Profile of the dominant terms in the constitutive equations (NSF) along the streamline $3.6 \times 10^{- 5}$ for the PTT viscoelastic model and $Wi = 100, 1000, 10 000$ ⁠. (a) λ component and $Wi = 100$ ⁠; (b) λ component and $Wi = 1000$ ⁠; (c) λ component and $Wi = 10 000$ ⁠; (d) μ component and $Wi = 100$ ⁠; (e) μ component and $Wi = 1000$ ⁠; (f) μ component and $Wi = 10 000$ ⁠; (g) ν component and $Wi = 100$ ⁠; (h) ν component and $Wi = 1000$ ⁠; (i) ν component and $Wi = 10 000$ ⁠.

FIG. 20.

View large Download slide

Profile of the dominant terms in the constitutive equations (NSF) along the streamline $3.6 \times 10^{- 5}$ for the Giesekus viscoelastic model and $Wi = 100, 1000, 10 000$ ⁠. (a) λ component and $Wi = 100$ ⁠; (b) λ component and $Wi = 1000$ ⁠; (c) λ component and $Wi = 10 000$ ⁠; (d) μ component and $Wi = 100$ ⁠; (e) μ component and $Wi = 1000$ ⁠; (f) μ component and $Wi = 10 000$ ⁠; (g) ν component and $Wi = 100$ ⁠; (h) ν component and $Wi = 1000$ ⁠; (i) ν component and $Wi = 10 000$ ⁠.

FIG. 21.

View large Download slide

Dominant terms using the NSF with $Wi = 1000$ ⁠. (a) PTT: λ component; (b) Giesekus: λ component; (c) PTT: μ component; (d) Giesekus: μ component; (e) PTT: ν component; (f) Giesekus: ν component.

V. DISCUSSION

The focus of this work was to verify the theoretical asymptotics for the PTT and Giesekus models at reentrant corners. This was investigated here for the $270 °$ corner that arises in the benchmark 4:1 contraction problem. The two models both contain quadratic stress terms that give them quite different behavior as compared to the Oldroyd-B model. While all three models share a similar asymptotic structure at a reentrant corner, their stress and boundary layer behaviors are very different. The key behavior for the PTT and Giesekus models is the dominance of the solvent stress over the polymer stress in contrast to Oldroyd-B, with thicker wall boundary layers. These features make the simulation of such models more attainable than Oldroyd-B.

The numerical simulation results of Sec. IV do verify the predicted theoretical behavior for the velocity and stress fields, which has been demonstrated over a large range of Weissenberg numbers. Not only has the theoretical stress and velocity gradient singularities of Table I been confirmed, the boundary layer balances of Table II are also seen to hold. Care was taken to ensure that a separating streamline was not present at the reentrant corner that may occur due to a lip vortex or extension of the salient corner vortex. This ensures that the reentrant corner remains sharp with the full corner angle producing a consistent stress singularity. This was one of the motivating reasons for choosing $β = 1 / 2$ rather than the benchmark value of $β = 1 / 9$ ⁠. For the parameter values considered here of solvent viscosity fraction $β = 1 / 2$ and the model-specific quadratic stress parameter $κ = 0.25$ ⁠, neither PTT nor Giesekus appears to form a lip vortex at the reentrant corner. This is noteworthy, as such a structure is present in the Oldroyd-B model and grows with the Weissenberg number (Alves et al., 2003), which would make the corner singularity verification more difficult.

In the derivation of the theoretical asymptotic results in Evans (2010a; 2010b), an artificial small parameter ϵ based on the length scale was used as the asymptotic expansion parameter and order the sizes of terms. An estimate for its value can be deduced from Fig. 9 at around $10^{- 3}$ (as an upper bound) for $Wi = 1$ ⁠, which is the length scale on which the NSF asymptotics begin to appear. We remark that far smaller length scales are required for all components of the CSF to capture their asymptotics, as seen in Fig. 8, where convergence of some components at some angles has yet to attain their limiting behaviors. Further, the boundary layer thicknesses would be at most $3.5 \times 10^{- 4}$ for PTT and $2 \times 10^{- 4}$ for Giesekus using Table I, which are quite narrow to resolve, even for the finest mesh M2 used here. The color maps Figs. 10 and 13 seem to confirm this.

The simulations for PTT and Giesekus did not encounter an upper Weissenberg number value using the NSF equations. In contrast, the CSF equations broke down at Weissenberg around 10², as found also in Alves et al. (2003). It appears that the NSF is not only advantageous over the CSF near-stress singularities but also for global simulation of these models. Future work will examine the properties of the NSF and compare with the log conformation approach for the CSF.

One aspect of the numerical results that need improving is the pressure variable. While the velocity field behaviors are acceptable, accurate computations of the pressure variable seem more difficult to obtain as evidenced in Figs. 6, 7, and 17, and moreover suffer with increasing Weissenberg number. This may be due to the numerical treatment of the momentum equation as a diffusion equation. Certainly as the solvent fraction β reduces and/or the Weissenberg number increases, a different approach to the projection scheme for the velocity and pressure may be required. One feature we noted was that the pressure slopes improved dramatically for the NSF when truncated channel lengths (⁠ $8 L \times 8 L$ ⁠) were used, as can be seen in Fig. 22. According to this figure, it is possible to see a better pressure convergence of the numerical scheme as compared with CSF for both fluid models. Therefore, high-order projection schemes for viscoelastic fluid flows could be considered for future work.

FIG. 22.

View large Download slide

Slopes for the pressure in a reduced geometry of the 4:1 contraction flow (⁠ $8 L \times 8 L$ ⁠). Simulation parameters: $Re = 0.01, We = 1, β = 0.5, κ = 0.25$ ⁠. The coarse and refined meshes have as $Δ_{x_{min}} = Δ_{y_{min}}, 5 \times 10^{- 4}$ and $1 \times 10^{- 4}$ ⁠, respectively. (a) PTT; (b) Giesekus.

ACKNOWLEDGMENTS

The authors would like to thank the financial support given by SPRINT/FAPESP Grant No. 2018/22242–0, The Royal Society Newton International Exchanges Grant No. 2015/NI150225, and the Núcleo de Processamento de Alto Desempenho of the Universidade Federal do Rio Grande do Norte—NPAD/UFRN to allow us to access their computer facilities. This research was developed with computational resources from Centro de Ciências Matemáticas Aplicadas à Indústria (CeMEAI) supported by FAPESP (Fundação de Amparo à Pesquisa do Estado de São Paulo) Grant No. 2013/07375–0. C. M. Oishi acknowledges CNPq (Conselho Nacional de Desenvolvimento Científico e Tecnólogico), Grant No. 307459/2016–0 and FAPESP Grant No. 2021/13833–7. F. Ruano Neto acknowledges the financial support of FAPESP Grant No. 2021/05727–2.

AUTHOR DECLARATIONS

Conflict of Interest

The authors have no conflicts to disclose.

Author Contributions

Jonathan Evans: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Supervision (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Irineu Lopes Palhares Junior: Data curation (equal); Formal analysis (equal); Software (equal); Validation (equal); Visualization (equal); Writing – review & editing (equal). Cassio Oishi: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Supervision (equal); Writing – original draft (equal); Writing – review & editing (equal). Fabiano Ruano Neto: Data curation (equal); Investigation (equal); Software (equal); Visualization (equal); Writing – review & editing (equal).

DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

References

1.

Afonso

,

A. M.

,

Oliveira

,

P. J.

,

Pinho

,

F. T.

, and

Alves

,

M. A.

, “

Dynamics of high-Deborah-number entry flows: A numerical study

,”

J. Fluid Mech.

677

,

272

–

304

(

2011

).

https://doi.org/10.1017/jfm.2011.84

Google Scholar

Crossref

2.

Afonso

,

A. M.

,

Pinho

,

F. T.

, and

Alves

,

M. A.

, “

The kernel-conformation constitutive laws

,”

J. Non-Newtonian Fluid Mech.

167–168

,

30

–

37

(

2012

).

https://doi.org/10.1016/j.jnnfm.2011.09.008

Google Scholar

Crossref

3.

Alleborn

,

N.

,

Nandakumar

,

K.

,

Raszillier

,

H.

, and

Durst

,

F.

, “

Further contributions on the two-dimensional flow in a sudden expansion

,”

J. Fluid Mech.

330

,

169

–

188

(

1997

).

https://doi.org/10.1017/S0022112096003382

Google Scholar

Crossref

4.

Alves

,

M. A.

,

Oliveira

,

P. J.

, and

Pinho

,

F. T.

, “

Benchmark solutions for the flow of Oldroyd-B and PTT fluids in planar contractions

,”

J. Non-Newtonian Fluid Mech.

110

(

1

),

45

–

75

(

2003

).

https://doi.org/10.1016/S0377-0257(02)00191-X

Google Scholar

Crossref

5.

Alves

,

M. A.

,

Oliveira

,

P. J.

, and

Pinho

,

F. T.

, “

On the effect of contraction ratio in viscoelastic flow through abrupt contractions

,”

J. Non-Newtonian Fluid Mech.

122

(

1–3

),

117

–

130

(

2004

).

https://doi.org/10.1016/j.jnnfm.2004.01.022

Google Scholar

Crossref

6.

Alves

,

M. A.

,

Oliveira

,

P. J.

, and

Pinho

,

F. T.

, “

Numerical methods for viscoelastic fluid flows

,”

Annu. Rev. Fluid Mech.

53

,

509

–

541

(

2021

).

https://doi.org/10.1146/annurev-fluid-010719-060107

Google Scholar

Crossref

7.

Battaglia

,

F.

,

Tavener

,

S. T.

,

Kulkami

,

A. K.

, and

Merkle

,

C. L.

, “

Bifurcation of low Reynolds number flows in symmetric channels

,”

AIAA J.

35

(

1

),

99

–

105

(

1997

).

https://doi.org/10.2514/2.68

Google Scholar

Crossref

8.

Bayat

,

J.

,

Emdad

,

H.

, and

Abouali

,

O.

, “

A mechanical model of partially liquefied vitreous dynamics induced by saccadic eye movement within a realistic shape of vitreous cavity

,”

Phys. Fluids

34

(

2

),

021905

(

2022

).

https://doi.org/10.1063/5.0079194

Google Scholar

Crossref

9.

Callaham

,

J. L.

,

Koch

,

J. V.

,

Brunton

,

B. W.

,

Kutz

,

J. N.

, and

Brunton

,

S. L.

, “

Learning dominant physical processes with data-driven balance models

,”

Nat. Commun.

12

,

1016

(

2021

).

https://doi.org/10.1038/s41467-021-21331-z

Google Scholar

Crossref

PubMed

10.

Davoodi

,

M.

,

Zografos

,

K.

,

Oliveira

,

P. J.

, and

Poole

,

R. J.

, “

On the similarities between the simplified Phan-Thien–Tanner model and the finitely extensible nonlinear elastic dumbbell (Peterlin closure) model in simple and complex flows

,”

Phys. Fluids

34

(

3

),

033110

(

2022

).

https://doi.org/10.1063/5.0083717

Google Scholar

Crossref

11.

Dean

,

W. R.

, and

Montagnon

,

P. E.

, “

On the steady motion of viscous liquid in a corner

,”

Math. Proc. Cambridge Philos. Soc.

45

(

3

),

389

–

394

(

1949

).

https://doi.org/10.1017/S0305004100025019

Google Scholar

Crossref

12.

Deliceoğlu

,

A.

,

Çelik

,

E.

, and

Gürcan

,

F.

, “

Singular treatment of viscous flow near the corner by using matched eigenfunctions

,”

Proc. Inst. Mech. Eng., Part C

233

(

5

),

1660

–

1676

(

2019

).

https://doi.org/10.1177/0954406218772603

Google Scholar

Crossref

13.

Demuren

,

A. O.

, and

Wilson

,

R. V.

, “

Estimating uncertainty in computations of two-dimensional separated flows

,”

J. Fluids Eng.

116

(

2

),

216

–

220

(

1994

).

https://doi.org/10.1115/1.2910257

Google Scholar

Crossref

14.

Drikakis

,

D.

, “

Bifurcation phenomena in incompressible sudden expansion flows

,”

Phys. Fluids

9

,

76

–

87

(

1997

).

https://doi.org/10.1063/1.869174

Google Scholar

Crossref

15.

Egger

,

H.

,

Rüde

,

U.

, and

Wohlmuth

,

B.

, “

Energy-corrected finite element methods for corner singularities

,”

SIAM J. Numer. Anal.

52

(

1

),

171

–

193

(

2014

).

https://doi.org/10.1137/120871377

Google Scholar

Crossref

16.

Evans

,

J. D.

, “

Re-entrant corner flows of Oldroyd-B fluids

,”

Proc. R. Soc. A

461

(

2060

),

2573

–

2603

(

2005a

).

https://doi.org/10.1098/rspa.2004.1410

Google Scholar

Crossref

17.

Evans

,

J. D.

, “

Re-entrant corner flows of UCM fluids: The initial formation of lip vortices

,”

Proc. R. Soc. A

461

(

2062

),

3169

–

3181

(

2005b

).

https://doi.org/10.1098/rspa.2005.1491

Google Scholar

Crossref

18.

Evans

,

J. D.

, “

Re-entrant corner flows of UCM fluids: The Cartesian stress basis

,”

J. Non-Newtonian Fluid Mech.

150

(

2–3

),

116

–

138

(

2008a

).

https://doi.org/10.1016/j.jnnfm.2007.10.018

Google Scholar

Crossref

19.

Evans

,

J. D.

, “

Re-entrant corner flows of UCM fluids: The natural stress basis

,”

J. Non-Newtonian Fluid Mech.

150

(

2–3

),

139

–

153

(

2008b

).

https://doi.org/10.1016/j.jnnfm.2007.10.019

Google Scholar

Crossref

20.

Evans

,

J. D.

, “

Re-entrant corner behaviour of the PTT fluid with a solvent viscosity

,”

J. Non-Newtonian Fluid Mech.

165

(

9–10

),

527

–

537

(

2010a

).

https://doi.org/10.1016/j.jnnfm.2010.01.011

Google Scholar

Crossref

21.

Evans

,

J. D.

, “

Re-entrant corner behaviour of the Giesekus fluid with a solvent viscosity

,”

J. Non-Newtonian Fluid Mech.

165

(

9–10

),

538

–

543

(

2010b

).

https://doi.org/10.1016/j.jnnfm.2010.01.010

Google Scholar

Crossref

22.

Evans

,

J. D.

, and

Oishi

,

C. M.

, “

Transient computations using the natural stress formulation for solving sharp corner flows

,”

J. Non-Newtonian Fluid Mech.

249

,

48

–

52

(

2017

).

https://doi.org/10.1016/j.jnnfm.2017.08.012

Google Scholar

Crossref

23.

Evans

,

J. D.

, and

Sibley

,

D. N.

, “

Re-entrant corner flows of PTT fluids in the Cartesian stress basis

,”

J. Non-Newtonian Fluid Mech.

153

(

1

),

12

–

24

(

2008

).

https://doi.org/10.1016/j.jnnfm.2007.11.006

Google Scholar

Crossref

24.

Evans

,

J. D.

, and

Sibley

,

D. N.

, “

Re-entrant corner flow for PTT fluids in the natural stress basis

,”

J. Non-Newtonian Fluid Mech.

157

(

1–2

),

79

–

91

(

2009

).

https://doi.org/10.1016/j.jnnfm.2008.09.008

Google Scholar

Crossref

25.

Evans

,

J. D.

,

Palhares Junior

,

I. L.

, and

Oishi

,

C. M.

, “

Stresses of PTT, Giesekus, and Oldroyd-B fluids in a Newtonian velocity field near the stick-slip singularity

,”

Phys. Fluids

29

(

12

),

121604

(

2017

).

https://doi.org/10.1063/1.4993782

Google Scholar

Crossref

26.

Evans

,

J. D.

,

Cuminato

,

J. A.

,

Palhares Junior

,

I. L.

, and

Oishi

,

C. M.

, “

Numerical study of the stress singularity in stick-slip flow of the Phan–Thien Tanner and Giesekus fluids

,”

Phys. Fluids

31

(

9

),

093101

(

2019a

).

https://doi.org/10.1063/1.5100730

Google Scholar

Crossref

27.

Evans

,

J. D.

,

França

,

H. L.

, and

Oishi

,

C. M.

, “

Application of the natural stress formulation for solving unsteady viscoelastic contraction flows

,”

J. Comput. Phys.

388

,

462

–

489

(

2019b

).

https://doi.org/10.1016/j.jcp.2019.02.045

Google Scholar

Crossref

28.

Evans

,

J. D.

,

França

,

H. L.

,

Palhares Junior

,

I. L.

, and

Oishi

,

C. M.

, “

Testing viscoelastic numerical schemes using the Oldroyd-B fluid in Newtonian kinematics

,”

Appl. Math. Comput.

387

,

125106

(

2020

).

https://doi.org/10.1016/j.amc.2020.125106

Google Scholar

Crossref

29.

Fattal

,

R.

, and

Kupferman

,

R.

, “

Constitutive laws for the matrix-logarithm of the conformation tensor

,”

J. Non-Newtonian Fluid Mech.

123

,

281

–

285

(

2004

).

https://doi.org/10.1016/j.jnnfm.2004.08.008

Google Scholar

Crossref

30.

Fattal

,

R.

, and

Kupferman

,

R.

, “

Time-dependent simulation of viscoelastic flows at high Weissenberg number using the log-conformation representation

,”

J. Non-Newtonian Fluid Mech.

126

,

23

–

37

(

2005

).

https://doi.org/10.1016/j.jnnfm.2004.12.003

Google Scholar

Crossref

31.

Fearn

,

R. M.

,

Mullin

,

T.

, and

Cliffe

,

K. A.

, “

Nonlinear flow phenomena in a symmetric sudden expansion

,”

J. Fluid Mech.

211

,

595

–

608

(

1990

).

https://doi.org/10.1017/S0022112090001707

Google Scholar

Crossref

32.

Fernandes

,

C.

,

Araujo

,

M. S.

,

Ferrás

,

L. L.

, and

Nóbrega

,

J. M.

, “

Improved both sides diffusion (iBSD): A new and straightforward stabilization approach for viscoelastic fluid flows

,”

J. Non-Newtonian Fluid Mech.

249

,

63

–

78

(

2017

).

https://doi.org/10.1016/j.jnnfm.2017.09.008

Google Scholar

Crossref

33.

Floryan

,

J. M.

, and

Czechowski

,

L.

, “

On the numerical treatment of corner singularity in the vorticity field

,”

J. Comput. Phys.

118

(

2

),

222

–

228

(

1995

).

https://doi.org/10.1006/jcph.1995.1094

Google Scholar

Crossref

34.

Gatski

,

T. B.

, and

Grosch

,

C. E.

, “

Embedded cavity drag in steady laminar flow

,”

AIAA J.

23

(

7

),

1028

–

1037

(

1985

).

https://doi.org/10.2514/3.9034

Google Scholar

Crossref

35.

Gerritsma

,

M. I.

, and

Phillips

,

T. N.

, “

On the use of characteristic variables in viscoelastic flow problems

,”

IMA J. Appl. Math.

66

(

2

),

127

–

147

(

2001

).

https://doi.org/10.1093/imamat/66.2.127

Google Scholar

Crossref

36.

Gerritsma

,

M. I.

, and

Phillips

,

T. N.

, “

On the characteristics and compatibility equations for the UCM model fluid

,”

Z. Angew. Math. Mech.

88

(

7

),

523

–

539

(

2008

).

https://doi.org/10.1002/zamm.200700058

Google Scholar

Crossref

37.

Giesekus

,

H.

, “

A unified approach to a variety of constitutive models for polymer fluids based on the concept of configuration-dependent molecular mobility

,” Rheol. Acta

21

,

366

(

1982

).

Google Scholar

38.

Hagen

,

T.

, and

Renardy

,

M.

, “

Boundary layer analysis of the Phan–Thien–Tanner and Giesekus model in high Weissenberg number flow

,”

J. Non-Newtonian Fluid Mech.

73

(

1–2

),

181

–

189

(

1997

).

https://doi.org/10.1016/S0377-0257(97)00035-9

Google Scholar

Crossref

39.

Hassager

,

O.

, “

Working group on numerical techniques

,”

J. Non-Newtonian Fluid Mech.

29

,

337

–

346

(

1988

).

https://doi.org/10.1016/0377-0257(88)85060-2

Google Scholar

Crossref

40.

Hawa

,

T.

,

Rusak

,

Z.

, and

Ruban

,

A. I.

, “

Numerical-asymptotic expansion matching for computing a viscous flow around a sharp expansion corner

,”

Theor. Comput. Fluid Dyn.

15

,

265

(

2002

).

https://doi.org/10.1007/s001620200055

Google Scholar

Crossref

41.

Hinch

,

E. J.

, “

The flow of an Oldroyd fluid around a sharp corner

,”

J. Non-Newtonian Fluid Mech.

50

,

161

(

1993

).

https://doi.org/10.1016/0377-0257(93)80029-B

Google Scholar

Crossref

42.

Holstein

,

H.

, and

Paddon

,

D. J.

, “

A singular finite difference treatment of re-entrant corner flow. I. Newtonian fluids

,”

J. Non-Newtonian Fluid Mech.

8

,

81

(

1981

).

https://doi.org/10.1016/0377-0257(81)80007-9

Google Scholar

Crossref

43.

Jeyasountharan

,

A.

,

D'Avino

,

G.

, and

Del Giudice

,

F.

, “

Confinement effect on the viscoelastic particle ordering in microfluidic flows: Numerical simulations and experiments

,”

Phys. Fluids

34

,

042015

(

2022

).

https://doi.org/10.1063/5.0090997

Google Scholar

Crossref

44.

Joseph

,

D. D.

,

Fluid Dynamics of Viscoelastic Liquids

(

Springer

,

New York

,

1990

).

Google Scholar

Crossref

45.

Joseph

,

D. D.

,

Renardy

,

M.

, and

Saut

,

J. C.

, “

Hyperbolicity and change of type in the flow of viscoelastic fluids

,”

Arch. Ration. Mech. Anal.

87

(

3

),

213

–

251

(

1985

).

https://doi.org/10.1007/BF00250725

Google Scholar

Crossref

46.

Keunings

,

R.

, “

On the high Weissenberg number problem

,”

J. Non-Newtonian Fluid Mech.

20

,

209

–

226

(

1986

).

https://doi.org/10.1016/0377-0257(86)80022-2

Google Scholar

Crossref

47.

Keunings

,

R.

, “

A survey of computational rheology

,” in

13th International Congress on Rheology

,

2000

.

Google Scholar

48.

Keunings

,

R.

, and

Euler

,

B.

, “

Advances in the computer modeling of the flow of polymeric liquids

,”

Comput. Fluid Dyn.

9

,

449

–

458

(

2001

).

Google Scholar

49.

Lagerstrom

,

P. A.

,

Matched Asymptotic Expansions

(

Springer

,

New York

,

1988

).

Google Scholar

Crossref

50.

Maqbool

,

K.

,

Siddiqui

,

A. M.

,

Mehboob

,

H.

,

Maqbool

,

K.

,

Mehboob

,

H.

, and

Jamil

,

Q.

, “

Study non-Newtonian synovial fluid flow by a recursive approach

,”

Phys. Fluids

34

,

111903

(

2022

).

https://doi.org/10.1063/5.0121918

Google Scholar

Crossref

51.

Minaeian

,

A.

,

Nili-AhmadAbadi

,

M.

,

Norouzi

,

M.

, and

Kim

,

K. C.

, “

Effects of viscoelasticity on the onset of vortex shedding and forces applied on a cylinder in unsteady flow regime

,”

Phys. Fluids

34

(

1

),

013106

(

2022

).

https://doi.org/10.1063/5.0080956

Google Scholar

Crossref

52.

Moffatt

,

H. K.

, “

Viscous and resistive eddies near a sharp corner

,”

J. Fluid Mech.

18

(

1

),

1

(

1964

).

https://doi.org/10.1017/S0022112064000015

Google Scholar

Crossref

53.

Niethammer

,

M.

,

Marschall

,

H.

,

Kunkelmann

,

C.

, and

Bothe

,

D.

, “

A numerical stabilization framework viscoelastic fluid flow using finite volume method general unstructured meshes

,”

Int. J. Numer. Methods Fluids

86

,

131

(

2018

).

https://doi.org/10.1002/fld.4411

Google Scholar

Crossref

54.

Ockendon

,

H.

, and

Ockendon

,

J. R.

,

Viscous Flow, Cambridge Texts in Applied Mathematics

(

Cambridge University Press

,

1995

).

Google Scholar

55.

Owens

,

R. G.

, and

Phillips

,

T. N.

,

Computational Rheology

(

Imperial College Press

,

2002

).

Google Scholar

Crossref

56.

Phan-Thien

,

N.

, and

Tanner

,

R. I.

, “

A new constitutive equation derived from network theory

,”

J. Non-Newtonian Fluid Mech.

2

(

4

),

353

–

365

(

1977

).

https://doi.org/10.1016/0377-0257(77)80021-9

Google Scholar

Crossref

57.

Phillips

,

T. N.

, “

Singular matched eigenfunction expansions for Stokes flow around a corner

,”

IMA J. Appl. Math.

42

(

1

),

13

–

26

(

1989

).

https://doi.org/10.1093/imamat/42.1.13

Google Scholar

Crossref

58.

Poole

,

R. J.

,

Davoodi

,

M.

, and

Zografos

,

K.

, “

On the similarities between the simplified Phan-Thien Tanner (sPTT) and FENE-P models

,”

Br. Soc. Rheol. Bull.

60

(

2

),

29

–

34

(

2019

).

Google Scholar

59.

Rallison

,

J. M.

, and

Hinch

,

E. J.

, “

The flow of an Oldroyd fluid past a reentrant corner: The downstream boundary layer

,”

J. Non-Newtonian Fluid Mech.

116

(

2–3

),

141

–

162

(

2004

).

https://doi.org/10.1016/j.jnnfm.2003.10.001

Google Scholar

Crossref

60.

Renardy

,

M.

, “

Mathematical analysis of viscoelastic flows

,”

Ann. Rev. Fluid Mech.

21

,

21

(

1989

).

https://doi.org/10.1146/annurev.fl.21.010189.000321

Google Scholar

Crossref

61.

Renardy

,

M.

, “

A matched solution for corner flow of the upper convected Maxwell fluid

,”

J. Non-Newtonian Fluid Mech.

58

,

83

(

1995

).

Google Scholar

62.

Renardy

,

M.

, “

Re-entrant corner behavior of the PTT fluid

,”

J. Non-Newtonian Fluid Mech.

69

(

1

),

99

–

104

(

1997

).

https://doi.org/10.1016/S0377-0257(96)01506-6

Google Scholar

Crossref

63.

Renardy

,

M.

, “

Current issues in non-Newtonian flows: A mathematical perspective

,”

J. Non-Newtonian Fluid Mech.

90

(

2–3

),

243

–

259

(

2000a

).

https://doi.org/10.1016/S0377-0257(99)00081-6

Google Scholar

Crossref

64.

Renardy

,

M.

,

Mathematical Analysis of Viscoelastic Flows

(

SIAM

,

2000b

).

Google Scholar

Crossref

65.

Schlichting

,

H.

, and

Gersten

,

K.

,

Boundary-Layer Theory

(

Springer

,

Berlin, Heidelberg

,

2000

).

Google Scholar

Crossref

66.

Shapira

,

M.

,

Degani

,

D.

, and

Weihs

,

D.

, “

Stability and existence of multiple solutions for viscous flow in suddenly enlarged channels

,”

Comput. Fluids

18

(

3

),

239

–

258

(

1990

).

https://doi.org/10.1016/0045-7930(90)90009-M

Google Scholar

Crossref

67.

van Dyke

,

M.

,

Perturbation Methods in Fluid Mechanics

(

Academic Press

,

1964

).

Google Scholar

68.

Walters

,

K.

, and

Webster

,

M. F.

, “

The distinctive CFD challenges of computational rheology

,”

Int. J. Numer. Methods Fluids

43

(

5

),

577

–

596

(

2003

).

https://doi.org/10.1002/fld.522

Google Scholar

Crossref

69.

Wang

,

C. Y.

, “

Flow in a channel with longitudinal ribs

,”

J. Fluids Eng.

116

(

2

),

233

–

237

(

1994

).

https://doi.org/10.1115/1.2910260

Google Scholar

Crossref

2022

Author(s)

Numerical verification of sharp corner behavior for Giesekus and Phan-Thien–Tanner fluids

I. INTRODUCTION

II. GOVERNING EQUATIONS

III. REVISITING THE ASYMPTOTIC RESULTS

IV. UNSTEADY PLANAR 4:1 CONTRACTION FLOW

A. Numerical verification: Corner vortex analysis and Couette correction

B. Numerical verification: Sharp corner asymptotics for $W i = O (1)$

C. Numerical study at high Weissenberg numbers

V. DISCUSSION

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

References

Citing articles via

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

CSF	$\overset{▽}{T^{p}}$	$\frac{κ}{(1 - β)} g (T^{p})$	$T^{p} / Wi$	$2 (1 - β) D / Wi$
PTT	$r^{- 0.7841}$	$r^{- 0.6572}$	$r^{- 0.3286}$	$r^{- 0.4555}$
Giesekus	$r^{- 0.7351}$	$r^{- 0.5592}$	$r^{- 0.2796}$	$r^{- 0.4555}$
NSF (2.18)	$v \cdot \nabla λ$	$2 μ \nabla \cdot w$	$\frac{κ}{Wi} g_{λ}$	$\frac{1}{Wi} (λ - \frac{1}{\| v \|^{2}})$
PTT	$r^{- 1.8743}$	$r^{- 1.6713}$	$r^{- 1.7462}$	$r^{- 1.4176}$
Giesekus	$r^{- 1.8241}$	$r^{- 1.5445}$	$r^{- 1.6482}$	$r^{- 1.3686}$
NSF (2.19)	$v \cdot \nabla μ$	$ν \nabla \cdot w$	$\frac{κ}{Wi} g_{μ}$	$\frac{1}{Wi} μ$
PTT	$r^{- 0.5823}$	$r^{- 0.3807}$	$r^{- 0.4554}$	$r^{- 0.1268}$
Giesekus	$r^{- 0.4555}$	$r^{- 0.1759}$	$r^{- 0.2796}$	r⁰
NSF (2.20)	$v \cdot \nabla ν$	⋯	$\frac{κ}{Wi} g_{ν}$	$\frac{1}{Wi} (ν - \| v \|^{2})$
PTT	$r^{0.7083}$		$r^{0.7604}$	$r^{1.089}$
Giesekus	$r^{0.9131}$		$r^{1.089}$	$r^{1.089}$

Numerical verification of sharp corner behavior for Giesekus and Phan-Thien–Tanner fluids

I. INTRODUCTION

II. GOVERNING EQUATIONS

III. REVISITING THE ASYMPTOTIC RESULTS

IV. UNSTEADY PLANAR 4:1 CONTRACTION FLOW

A. Numerical verification: Corner vortex analysis and Couette correction

B. Numerical verification: Sharp corner asymptotics for Wi=O(1)

C. Numerical study at high Weissenberg numbers

V. DISCUSSION

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

References

Citing articles via

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

B. Numerical verification: Sharp corner asymptotics for $W i = O (1)$