Journal of Non-Newtonian Fluid Mechanics 249 (2017) 48–52 Contents lists available at ScienceDirect Journal of Non-Newtonian Fluid Mechanics journal homepage: www.elsevier.com/locate/jnnfm Short Communication Transient computations using the natural stress formulation for solving sharp corner flows J.D. Evans a , ∗, C.M. Oishi b a Department of Mathematical Sciences, University of Bath, Bath, BA2 7AY, United Kingdom b Departamento de Matemática e Computação, Faculdade de Ciências e Tecnologia, Universidade Estadual Paulista “Júlio de Mesquita Filho” 19060-900 Presidente Prudente, Sao Paulo, Brazil a r t i c l e i n f o Article history: Received 26 May 2017 Revised 24 August 2017 Accepted 30 August 2017 Available online 19 September 2017 Keywords: Unsteady viscoelastic flows Oldroyd-B model Natural stress formulation Planar contraction Sharp corner a b s t r a c t In this short communication, we analyse the potential of the natural stress formulation (NSF) (i.e. align- ing the stress basis along streamlines) for computing planar flows of an Oldroyd-B fluid around sharp corners. This is the first attempt to combine the NSF into a numerical strategy for solving a transient fluid flow problem considering the momentum equation in Navier–Stokes form (the elastic stress enter- ing as a source term) and using the constitutive equations for natural stress variables. Preliminary results of the NSF are motivating in the sense that accuracy of the numerical solution for the extra stress tensor is improved near to the sharp corner. Comparison studies among the NSF and the Cartesian stress for- mulation (CSF) (i.e. using a fixed Cartesian stress basis) are conducted in a typical benchmark viscoelastic fluid flow involving a sharp corner: the 4 : 1 contraction. The CSF needs a mesh approximately 10 times smaller to capture similar near singularity results to the NSF. © 2017 Elsevier B.V. All rights reserved. c s s s d b c c p f i e c s b c [ [ r s 1. Introduction Flow through a contraction is a benchmark problem in com- putational rheology [9] , where viscoelastic fluids exhibit regions of strong shearing near the walls and uniaxial extension along the centreline. The complex flow patterns that evolve have been the subject of much interest in the literature, with attention fo- cused on: (i) vortex behaviour, both near the re-entrant corner (so- called lip vortices) and salient corners, (ii) variation of the pressure drops across the contraction with strength of fluid elasticity (Weis- senberg number), (iii) particle paths upstream of the contraction and (iv) velocity overshoots along the axis of symmetry. The main numerical approaches to simulate this flow, have been finite-difference [21] , finite-element [8] and finite-volume [2] . Combinations of the methods, e.g. hybrid finite-element finite- volume [1] as well as Lagrangian and semi-Langrangian methods [15,22] are also commonly employed. However, a key feature of all numerical schemes so far employed is that they discretize the vis- coelastic constitutive equations formulated using a fixed Cartesian basis for the stress. We refer to this formulation as the Cartesian stress formulation (CSF) of the constitutive equations. An alterna- tive approach is to exploit the mathematical structure of the upper ∗ Corresponding author. E-mail addresses: masjde@maths.bath.ac.uk , masjde@bath.ac.uk (J.D. Evans), oishi@fct.unesp.br (C.M. Oishi). i t o https://doi.org/10.1016/j.jnnfm.2017.08.012 0377-0257/© 2017 Elsevier B.V. All rights reserved. onvected derivative and align the stress basis with the flow using treamlines. This formulation uses the velocity field to span the tress field, the formulation of the constitutive equations in this etting being termed the natural stress formulation (NSF). The natural stress formulation was first used by Renardy [17] , to emonstrate its ability to eliminate the downstream stress insta- ility encountered during numerical integration around re-entrant orners [18] . Although the idea of transforming the stress tensor omponents to a basis aligned with streamlines for computation urposes had been recognised previously in [3,11] . However, the ull power of the approach has not yet been exploited numerically n a mathematically systematic way for the full contraction geom- try. A key feature of the geometry is the presence of the re-entrant orner at which the velocity gradients and stress are infinite. The ingularity determination for the Oldroyd-B fluid was first given y Hinch [10] , with the asymptotic structure of the local solution ompleted through the upstream wall boundary layer by Renardy 19] and the boundary layer at the downstream wall in [16] and 4,5] . An important aspect of the solution analysis is that the natu- al stress formulation is an efficient way to transmit the necessary tress information from the upstream to downstream regions. This s explicitly calculated in [6,7] for the UCM fluid, illustrating that he necessary stress information is contained in high-order terms f the asymptotic expansions when using the Cartesian stress com- https://doi.org/10.1016/j.jnnfm.2017.08.012 http://www.ScienceDirect.com http://www.elsevier.com/locate/jnnfm http://crossmark.crossref.org/dialog/?doi=10.1016/j.jnnfm.2017.08.012&domain=pdf mailto:masjde@maths.bath.ac.uk mailto:masjde@bath.ac.uk mailto:oishi@fct.unesp.br https://doi.org/10.1016/j.jnnfm.2017.08.012 J.D. Evans, C.M. Oishi / Journal of Non-Newtonian Fluid Mechanics 249 (2017) 48–52 49 p a e c l 2 v O ∇ R T H i r a t R w s a b t i τ n 2 d v a T T C 2 t R s T t[ W v w s A w d g t A w v T n[ [ [ w | T onents, with the natural stress variables being able to uncouple nd extract this information. The success of the natural stress formulation near the re- ntrant corner singularity, both for asymptotics and numerical omputation, are the motivating reasons to investigate the formu- ation for the full contraction geometry. . Mathematical formulation The governing equations for the incompressible flow of a iscoelastic fluid we adopt are the continuity, momentum and ldroyd-B constitutive equations in dimensionless form · v = 0 , (1) e ( ∂v ∂t + v · ∇v ) = −∇ p + β∇ 2 v + ∇ · T , (2) + Wi ( ∂T ∂t + ( v · ∇ ) T − ( ∇v ) T − T ( ∇v ) T ) = 2(1 − β) D . (3) ere v is the velocity field, p is an arbitrary isotropic pressure, T s the polymeric contribution to the extra-stress tensor and D = 1 2 (∇v + ( ∇v ) T ) is the rate-of-strain tensor. The dimensionless pa- ameters are the Reynolds number Re, Weissenberg number Wi nd retardation parameter β ∈ [0, 1] (the dimensionless retardation ime or solvent viscosity ratio) defined by e = ρUL ηs + ηp , Wi = λ1 U L , β = ηs ηs + ηp , (4) here ρ is the density, U and L are characteristic length and flow peeds respectively, λ1 the relaxation time, ηs the solvent viscosity nd ηp the polymer viscosity. The above governing equations have een made non-dimensional using L for the spatial variables, U for he velocity scaling and (ηp + ηs ) U/L the pressure and stress scal- ngs. The total stress is σ = −pI + τ, with the extra-stress tensor = T + 2 βD being rheologically composed of polymer and Newto- ian solvent contributions. .1. Cartesian stress formulation Denoting by i and j the unit vectors in fixed Cartesian x and y irections, we have = u i + v j = ( u, v ) T (5) nd = T 11 ii T + T 12 ( ij T + ji T ) + T 22 jj T . (6) he component form of the polymer constitutive Eq. (3) for the artesian extra-stresses T 11 , T 12 , T 22 is T 11 + Wi ( ∂T 11 ∂t + u ∂T 11 ∂x + v ∂T 11 ∂y − 2 ∂u ∂x T 11 − 2 ∂u ∂y T 12 ) = 2 ( 1 − β) ∂u ∂x , T 22 + Wi ( ∂T 22 ∂t + u ∂T 22 ∂x + v ∂T 22 ∂y − 2 ∂v ∂y T 22 − 2 ∂v ∂x T 12 ) = 2 ( 1 − β) ∂v ∂y , T 12 + Wi ( ∂T 12 ∂t + u ∂T 12 ∂x + v ∂T 12 ∂y − ∂v ∂x T 11 − ∂u ∂y T 22 ) = ( 1 − β) ( ∂u + ∂v ) . (7) ∂y ∂x T .2. Natural stress formulation Aligning the polymer stress basis along streamlines, introduces he so called natural stress variables. We follow the construction of enardy [17] (see also [20,23] ). Introducing the configuration ten- or A by = (1 − β) Wi ( A − I ) , (8) he polymer constitutive Eq. (3) becomes ∂A ∂t + ( v · ∇ ) A − ( ∇v ) A − A ( ∇v ) T ] + 1 Wi ( A − I ) = 0 . (9) e now express A in terms of the dyadic products of the velocity and an orthogonal vector w defined as follows = 1 | v | 2 ( −v , u ) T , o that = λv v T + μ( v w T + w v T ) + νw w T , (10) here λ, μ, ν are termed the natural stress variables. However, a etraction of this construction is that the basis vectors v, w are de- enerate when the velocity field vanishes. As such, it is convenient o use instead, unit vectors in their directions. Hence, we write = ˆ λˆ v ̂ v T + ˆ μ( ̂ v ̂ w T + ˆ w ̂ v T ) + ˆ ν ˆ w ̂ w T , (11) here ˆ = v | v | , ˆ w = | v | w , ˆ λ = | v | 2 λ, ˆ μ = μ, ˆ ν = ν | v | 2 . (12) he scaled natural stress variables ˆ λ, ˆ μ, ̂ ν then satisfy the compo- ent equations ∂ ̂ λ ∂t + 2 ̂ μ | v | 2 ( v ∂u ∂t − u ∂v ∂t ) + | v | 2 ( v · ∇) ( ˆ λ | v | 2 ) + 2 ̂ μ| v | 2 ∇ · w ] + 1 Wi ( ˆ λ − 1 ) = 0 , ∂ ˆ μ ∂t + ( ˆ λ − ˆ ν | v | 2 )( u ∂v ∂t − v ∂u ∂t ) + ( v · ∇) ̂ μ + ˆ ν| v | 2 ∇ · w ] + ˆ μ Wi = 0 , ∂ ̂ ν ∂t + 2 ̂ μ | v | 2 ( u ∂v ∂t − v ∂u ∂t ) + 1 | v | 2 ( v · ∇) ( ˆ ν| v | 2 ) ] + 1 Wi ( ˆ ν − 1 ) = 0 , (13) ith v | 2 ∇ . w = | v | 2 ( ∂ ∂x ( − v | v | 2 ) + ∂ ∂y ( u | v | 2 )) = 1 | v | 2 ( (v 2 − u 2 ) ( ∂v ∂x + ∂u ∂y ) + 4 u v ∂u ∂x ) . he component form of (8) is T 11 = ( 1 − β) Wi ( −1 + 1 | v | 2 ( ˆ λu 2 − 2 ̂ μuv + ˆ νv 2 )) , T 12 = ( 1 − β) Wi | v | 2 ( ˆ λuv + ˆ μ ( u 2 − v 2 ) − ˆ νuv ) , 22 = ( 1 − β) Wi ( −1 + 1 | v | 2 ( ˆ λv 2 + 2 ̂ μuv + ˆ νu 2 )) , (14) 50 J.D. Evans, C.M. Oishi / Journal of Non-Newtonian Fluid Mechanics 249 (2017) 48–52 Fig. 1. a) Sketch of the 4 : 1 contraction geometry. b) Sketch of the mesh around the re-entrant corner. a d u a T S t λ A f n t s i a b 4 i n w t t a with the inverse relationship being ˆ λ − 1 = Wi ( 1 − β) | v | 2 ( u 2 T 11 + 2 uv T 12 + v 2 T 22 ) , ˆ μ= Wi ( 1 − β) | v | 2 ( −uv T 11 + ( u 2 − v 2 ) T 12 + uv T 22 ) , ˆ ν − 1 = Wi ( 1 − β) | v | 2 ( v 2 T 11 − 2 uv T 12 + u 2 T 22 ) . (15) These expressions link the Cartesian and natural stress components together. 3. Numerical method 3.1. CSF and NSF schemes The first step involved in the algorithm of the schemes is the solution of the conservative Eqs. (1) and (2) which can be written as follows : ∂u ∂x + ∂v ∂y = 0 , (16) ∂u ∂t − β Re ( ∂ 2 u ∂x 2 + ∂ 2 u ∂y 2 ) = − ( ∂ ( uu ) ∂x + ∂ ( vu ) ∂y ) − 1 Re ∂ p ∂x + 1 Re ( ∂T 11 ∂x + ∂T 12 ∂y ) , (17) ∂v ∂t − β Re ( ∂ 2 v ∂x 2 + ∂ 2 v ∂y 2 ) = − ( ∂ ( uv ) ∂x + ∂ ( vv ) ∂y ) − 1 Re ∂ p ∂y + 1 Re ( ∂T 12 ∂x + ∂T 22 ∂y ) . (18) For this purpose, in the context of a finite difference Marker-And- Cell (MAC) methodology, a projection scheme for solving viscoelas- tic problems is applied in order to uncouple velocity and pressure fields in Eqs. (16) –(18) , as carefully detailed in [13] . In summary, the momentum Eqs. (17) and (18) are solved via a semi-implicit scheme maintaining the convective terms, pressure gradient and divergence of the stress tensor discretized in the n time-level. After this, the projection step is applied for obtaining the final velocity and pressure fields (u n +1 , v n +1 , p n +1 ) . Once the final velocity field is computed, the constitutive equa- tions of the CSF or NSF can be solved. In the case of the CSF, the final extra stress components (T n +1 11 , T n +1 12 , T n +1 22 ) are directly calcu- lated solving Eq. (7) by an explicit scheme [13] . For the NSF, the fi- nal values for the NS variables ( ̂ λn +1 , ˆ μn +1 , ̂ νn +1 ) are computed by solving Eq. (13) via an explicit discretization scheme. The relation- ships (14) are then used to construct the necessary Cartesian extra stress components for use in the momentum equation. Therefore, according to the choice of the scheme (CSF or NSF), at this stage of the algorithm, all variables are updated; a new computational cycle then begins and the process continues until steady state is reached. 3.2. Initial and boundary conditions In order to solve Eq. (13) numerically, initial conditions for the natural stress variables ˆ λ, ˆ μ and ˆ ν are required. For instance, at t = 0 , based on the relations (15) and considering that the Cartesian extra-stress tensor is initialized as T = 0 , we have adopted ˆ λ = ˆ ν = 1 and ˆ μ = 0 . The boundary conditions used in this work for the application of the NSF are: prescribed inflows, outflows, and rigid walls. At the inlet, the boundary conditions on the NS variables are constructed s function of the extra-stress tensor. For instance, assuming fully eveloped flow, we have imposed : = 3 4 ( 1 − y 2 16 ) , v = 0 , −4 ≤ y ≤ 4 , (19) nd T 11 = 2 ( 1 − β) Wi ( ∂u ∂y )2 , T 12 = ( 1 − β) ( ∂u ∂y ) , 22 = 0 . (20) ubstituting v = 0 in (15) , we can obtain the inlet boundary condi- ions for the natural stress variables: ˆ = Wi ( 1 − β) T 11 + 1 , ˆ μ = Wi ( 1 − β) T 12 , ˆ ν = Wi ( 1 − β) T 22 + 1 . (21) t an outlet, the homogeneous Neumann condition is employed or the velocity field and for the NS variables. For rigid walls, the o-slip condition is used for the velocity field. It is not necessary o impose boundary conditions for the ˆ λ, ˆ μ and ˆ ν at rigid walls ince the Eq. (13) are considered hyperbolic. However, in order to mprove the stability behaviour of the numerical method, we have pplied a linear extrapolation of the NS variables at this type of oundary. . A preliminary result: planar contraction flow The applicability and accuracy of the NSF is now investigated n a typical benchmark problem involving a sharp re-entrant cor- er: the 4: 1 planar contraction flow (see Fig. 1 a)). For this test, e have considered U = 1 . 0 m/s for the velocity scaling while he length scaling is the half-width of the narrower channel of he contraction with L = 1 . 0 m. Moreover, some simplifications are dopted: • Reduction of the channel lengths. The main objective of this work is to present an alternative strategy for improving the computation of the non-Newtonian stress tensor around sharp corners. Therefore, the contraction domain is reduced in order to save CPU time and to give more emphasis on the numeri- cal solution near to the re-entrant corner of the problem (see Fig. 1 b)). The domain description is presented in Fig. 1 a) where the inlet (upstream) and outlet (downstream) channels are con- structed using a length of 8 L . This simplification is acceptable since we are imposing fully developed velocity and stress ten- sor profiles at the inlet (see Section 3.2 ). In addition, a relatively small value for the Weissenberg number ( Wi = 1 ) is applied in J.D. Evans, C.M. Oishi / Journal of Non-Newtonian Fluid Mechanics 249 (2017) 48–52 51 Fig. 2. Comparison between NSF and CSF for the transient behaviour of T 11 near to singularity using uniform (M1) and non-uniform (M2) meshes with Re = 1 , Wi = 1 and β = 0 . 5 . Table 1 Details of the uniform and non-uniform meshes. Mesh Space-step Numbers of cells M1 (uniform) δx = δy = 0 . 05 32,0 0 0 M2 (non-uniform) δx min = δy min = 0 . 004 39,600 t p f a m g o g t t f w N i p t t s s s w u o e f t t f fi t [ m 5 m m p f s s a t all simulations of this Section to avoid the presence and growth of a lip vortex. • Regularization. The basis vectors used in (10) (or in Eq. (11) ) become degenerate where the velocity field vanishes, typically near solid boundaries and stagnation points. Although using the scaled variables (12) is beneficial, the transformation remains singular and care is needed in using (13) . Therefore, terms are regularised in these equations in order to avoid numerical di- vision by zero, where it may be encountered. In particular, a regularization constant tol = 10 −6 is adopted in all simulations. Further, in all simulations, we have set the following parame- ers: Re = 1 , Wi = 1 . 0 and β = 0 . 5 . The influence of varying these arameters as well as the regularization constant tol will be care- ully investigated in further work. Results for the two formulations re presented on a coarse uniform mesh M1 and a non-uniform esh M2. Details about these meshes can be found in Table 1 . Firstly, we consider the transient behaviour of T 11 near to the eometric singularity of the contraction, specifically, in the center f the shaded cell shown in Fig. 1 b). Fig. 2 shows that the NSF ives considerably larger T 11 values than the CSF in the region of he sharp corner for both meshes. The growth of the extra stress ensor when r → 0, evident in Fig. 2 for the mesh M2, is expected rom the asymptotic analysis which indicates that T 11 ∝ r − 2 3 . It is orth remarking in Fig. 2 that the value of T 11 produced by the SF in the coarse mesh is similar to the value reached by the CSF n a non-uniform mesh (approximately 10 times smaller). In Fig. 3 , we have depicted the asymptotic behaviour of the roperties considering results along of the r line (vertical direc- ion) illustrated in Fig. 1 b) i.e. θ = π/ 2 in Fig. 1 a). We have chosen his direction consistent with the literature for ease of compari- on of the results (e.g. [1,2,12,14] ). The results in Fig. 3 for the CSF cheme on the non-uniform mesh agree very well with those pre- ented in [2,14] for comparable mesh spacing. From these figures, e observe a significant improvement of the numerical solutions sing the NSF. As we can see in Fig. 3 a), even with the application f the NSF on a coarse mesh, the values of log (| T 11 |) for the near- st three points of the corner are fitted with the slope expected or the stress. For the other stress components, the numerical solu- ions obtained by the NSF are in good agreement with the asymp- otic analysis for r → 0. In particular, it can be observed that both ormulations produce similar results for the velocity field. In all gures, we have included the expected slope for the extra stress ensor and the velocity field according to asymptotic analyses of 5,10] . It is worth mentioning that away from the singularity, nu- erical results for both formulations are very similar. . Discussion In this work, we have numerically investigated an alternative ethod for solving sharp corner flows involving the Oldroyd-B odel. The mathematical formulation is based on a dyadic decom- osition of the conformation tensor resulting in the natural stress ormulation originally proposed by Renardy [17] . The main contributions of this short communication can be ummarized as follows: • The NSF equations relevant to modelling two dimensional un- steady flows have been presented and used to numerically solve a complex flow; • The numerical results indicate that the NSF gives an effective formulation for capturing the expected behaviour of the stress near to the singularity even in the application of coarse meshes. In the extension of this work, further details regarding the con- truction of the numerical method for the NSF will be provided long with the possibility of devising a hybrid scheme in which he CSF is used at solid boundaries. 52 J.D. Evans, C.M. Oishi / Journal of Non-Newtonian Fluid Mechanics 249 (2017) 48–52 Fig. 3. Asymptotic results for the uniform (M1) and non-uniform (M2) meshes near to the corner using Re = 1 , Wi = 1 and β = 0 . 5 : a) log (| T 11 |), b) log (| T 12 |), c) log (| T 22 |), and d) log (| u |) and log (| v |). [ Acknowledgements The authors would like to thank the financial support given by SPRINT / FAPESP grant no. 2015/50094-7 and The Royal Soci- ety Newton International Exchanges grant no. 2015/NI150225. C.M. Oishi acknowledges CNPq (Conselho Nacional de Desenvolvimento Científico e Tecnologico), Grant 307459/2016-0 , and FAPESP, Grant 2013/07375-0 . The authors would like to thank Hugo L. França for his support and assistance with the non-uniform mesh generation. References [1] M. Aboubacar , H. Matallah , M.F. Webster , Highly elastic solutions for Oldroyd- B and Phan-Thien/Tanner fluids with a finite volume/element method: planar contraction flows, J. Non-Newtonian Fluid Mech. 103 (2002) 65–103 . [2] M.A. Alves , P.J. Oliveira , F.T. Pinho , Benchmark solutions for the flow of Oldroyd-B and PTT fluids in planar contractions, J. Non-Newtonian Fluid Mech. 110 (2003) 45–75 . [3] S. Dupont , J.M. Marchal , M.J. Crochet , Finite element simulation of viscoelastic fluids of the integral type, J. Non-Newtonian Fluid Mech. 17 (1985) 157–183 . [4] J.D. Evans , Re-entrant corner flows of the upper convected Maxwell fluid, Proc. Roy. Soc. A 461 (2005) 117–142 . [5] J.D. Evans , Re-entrant corner flows of Oldroyd-B fluids, Proc. Roy. Soc. A 461 (2005) 2573–2603 . [6] J.D. Evans , Re-entrant corner flows of UCM fluids: the Cartesian stress basis, J. Non-Newtonian Fluid Mech. 150 (2008) 116–138 . [7] J.D. Evans , Re-entrant corner flows of UCM fluids: the natural stress basis, J. Non-Newtonian Fluid Mech. 150 (2008) 139–153 . [8] R. Guènette , M. Fortin , A new mixed finite element method for computing vis- coelastic flows, J. Non-Newtonian Fluid Mech. 60 (1995) 27–52 . [9] O. Hassager , Working group on numerical techniques, in: proceedings of the fifth international workshop on numerical methods in non-Newtonian flows, lake arrowhead, USA, J. Non-Newtonian Fluid Mech. 29 (1988) 2–5 . [10] E.J. Hinch , The flow of an Oldroyd fluid around a sharp corner, J. Non- Newtonian Fluid Mech. 50 (1993) 161–171 . [11] R.A. Keiller , Entry-flow calculations for the Oldroyd-B and FENE equations, J. Non-Newtonian Fluid Mech. 46 (1993) 143–178 . [12] J.M. Kim , C. Kim , J.H. Kim , C. Chung , K.H. Ahn , S.J. Lee , High-resolution fi- nite element simulation of 4:1 planar contraction flow of viscoelastic fluid, J. Non-Newtonian Fluid Mech. 129 (2005) 23–37 . [13] C.M. Oishi , F.P. Martins , M.F. Tomé, J.A. Cuminato , S. McKee , Numerical solu- tion of the extended pom-pom model for viscoelastic free surface flows, J. Non-Newtonian Fluid Mech. 166 (2011) 165–179 . [14] F. Pimenta , M.A. Alves , Stabilization of an open-source finite-volume solver for viscoelastic fluid flows, J. Non-Newtonian Fluid Mech. 239 (2017) 85–104 . [15] T.N. Phillips , A.J. Williams , Viscoelastic flow through a planar contraction us- ing a semi-Lagrangian finite volume method, J. Non-Newtonian Fluid Mech. 87 (1999) 215–246 . [16] J.M. Rallison , E.J. Hinch , The flow of an Oldroyd fluid past a reentrant cor- ner: the downstream boundary layer, J. Non-Newtonian Fluid Mech. 116 (2004) 141–162 . [17] M. Renardy , How to integrate the upper convected Maxwell (UCM) stresses near a singularity (and maybe elsewhere, too), J. Non-Newtonian Fluid Mech. 52 (1994) 91–95 . [18] M. Renardy , The stresses of an upper convected Maxwell fluid in a Newtonian velocity field near a re-entrant corner, J. Non-Newtonian Fluid Mech. 50 (1993) 127–134 . [19] M. Renardy , A matched solution for corner flow of the upper convected Maxwell fluid, J. Non-Newtonian Fluid Mech. 58 (1995) 83–89 . [20] M. Renardy , The high Weissenberg number limit of the UCM model and the Euler equations, J. Non-Newtonian Fluid Mech. 69 (1997) 293–301 . [21] D. Trebotich , P. Colella , G.H. Miller , A stable and convergent scheme for vis- coelastic flow in contraction channels, J. Comput. Phys 205 (2005) 315–342 . [22] P. Wapperom , R. Keunings , V. Legat , The backward-tracking Lagrangian parti- cle method for transient viscoelastic flows, J. Non-Newtonian Fluid Mech. 91 (20 0 0) 273–295 . 23] P. Wapperom , M. Renardy , Numerical prediction of the boundary layers in the flow around a cylinder using a fixed velocity field, J. Non-Newtonian Fluid Mech. 125 (2005) 35–48 . http://dx.doi.org/10.13039/100004730 http://dx.doi.org/10.13039/501100001807 http://dx.doi.org/10.13039/501100003593 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0001 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0001 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0001 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0001 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0002 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0002 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0002 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0002 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0003 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0003 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0003 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0003 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0004 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0004 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0005 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0005 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0006 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0006 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0007 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0007 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0008 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0008 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0008 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0009 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0009 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0010 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0010 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0011 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0011 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0012 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0013 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0014 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0014 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0014 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0015 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0015 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0015 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0016 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0016 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0016 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0017 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0017 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0018 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0018 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0019 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0019 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0020 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0020 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0021 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0021 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0021 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0021 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0022 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0022 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0022 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0022 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0023 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0023 http://refhub.elsevier.com/S0377-0257(17)30249-5/sbref0023 Transient computations using the natural stress formulation for solving sharp corner flows 1 Introduction 2 Mathematical formulation 2.1 Cartesian stress formulation 2.2 Natural stress formulation 3 Numerical method 3.1 CSF and NSF schemes 3.2 Initial and boundary conditions 4 A preliminary result: planar contraction flow 5 Discussion Acknowledgements References