Journal of Non-Newtonian Fluid Mechanics 283 (2020) 104338 Contents lists available at ScienceDirect Journal of Non-Newtonian Fluid Mechanics journal homepage: www.elsevier.com/locate/jnnfm An Oldroyd-B solver for vanishingly small values of the viscosity ratio: Application to unsteady free surface flows C. Viezel a , M.F. Tomé a,∗ , F.T. Pinho b , S. McKee c a Department of Applied Mathematics and Statistics, University of São Paulo, 13560-970, São Carlos, São Paulo, Brazil b CEFT, Department of Mechanical Engineering, Faculty of Engineering, University of Porto, Rua Dr. Roberto Frias s/n, 4200-465, Porto, Portugal c Department of Mathematics and Statistics, University of Strathclyde, 26 Richmond Street, G11XH, Glasgow, U.K. a r t i c l e i n f o Keywords: Oldroyd-B model EVSS Free surface Finite difference Drop impact Bouncing drop a b s t r a c t This work is concerned with time-dependent axisymmetric free surface flows of Oldroyd-B fluids for any value of , the ratio of solvent to total viscosities. The Oldroyd-B constitutive equation is dealt with by employing a novel technique to calculate the conformation tensor while an EVSS transformation allows the solution of the momentum equations coupled with the free surface stress conditions: this avoids numerical instabilities that can arise when using small values of . The convergence of this new methodology is verified on pipe flow and also by comparing results from the literature for the time-dependent impacting drop problem. This approach is then used to predict the time-dependent free surface flow after a viscoelastic drop impacts a solid surface for  values in the range [0, 1]. The impacting drop problem is investigated for polymer solutions containing a small solvent contribution ( → 0) or without any solvent viscosity ( = ). In addition, a study of the bouncing drop problem for different values of , Weissenberg and Reynolds numbers is undertaken. 1. Introduction The importance of non-Newtonian free surface flows in industrial processes has attracted the attention of many scientists. Examples of such applications include polymer processing in the plastics industry such as mould filling of complex cavities. Such flows can be modelled by a system of nonlinear equations, but the presence of (multiple) mov- ing free surfaces can make their solution challenging. Furthermore, for a particular polymer, it is not always obvious what the correct constitu- tive equation should be. One choice that must be made is between dif- ferential and integral constitutive models, or indeed a mixture of both. A large number of differential constitutive models have been developed over the past decades: Upper Convected-Maxwell (UCM) [1], Oldroyd- B [2], Phan-Thien-Tanner (PTT) [3], Giesekus [4], Extended Pom-Pom (Pom-Pom) [6,7], among others. On the other hand, integral constitu- tive models have been developed and studied by Papanastasiou et al. [8], Kaye [9] and Luo and Mitsoulis [10], amongst others. Integral con- stitutive equations require more sophisticated approaches to solve the governing equations numerically and require more computational re- sources and, possibly for these reasons, there has been a greater focus on differential constitutive models. In particular, the UCM and Oldroyd- B models have been extensively studied, employing finite element, finite volume and finite difference methods (e.g. [11–16,23–30]). Due to their unbounded elastic normal stresses, these models are arguably the most ∗ Corresponding author. E-mail addresses: murilo@icmc.usp.br (M.F. Tomé), fpinho@fe.up.pt (F.T. Pinho), sean.mckee@strath.ac.uk (S. McKee). challenging viscoelastic constitutive equations from a numerical point of view (see e.g. [33,34]). A decoupling strategy to calculate velocity and pressure has found favour. For instance, Hirt and Nichols [35] introduced the volume of fluid (VOF) method in the early 1980s: this has been used to simulate non-Newtonian flows by many investigators (e.g. [15,32,36–39]). This method, while easy to implement, suffers from numerical diffusion; to overcome this drawback several improved versions have been developed [36,38,40]. Another approach is to represent the free surface by a level set func- tion which is convected with the fluid flow; its evolution in time is ob- tained through the solution of a hyperbolic equation. Osher and Sethian [41] are usually credited with introducing the idea and it has the ca- pability of capturing multiphase flow phenomena. It has been used to simulate filament stretching and jet buckling [29,42,43], mould fill- ing [44] and many other interesting applications (see e.g. [44–49]). A third approach is the front tracking method: unlike VOF or the level set method the front tracking method employs massless markers to describe the fluid interface. In two dimensions the interface between two fluids is described by a set of points (x i , y i ) - the markers - while in three di- mensions it is represented by a set of quadrilaterals and/or triangles. The coordinates of the markers are updated at each time step according to the velocity at the new time step. In two dimensions the interface is visualized by connecting these points by straight lines (i.e. zero or- https://doi.org/10.1016/j.jnnfm.2020.104338 Received 7 February 2020; Received in revised form 20 June 2020; Accepted 24 June 2020 Available online 28 June 2020 0377-0257/© 2020 Elsevier B.V. All rights reserved.