Acessibilidade / Reportar erro

High-order finite volume method for solving viscoelastic fluid lows

Abstract

Computational Fluid Dynamics (CFD) is widely used by polymer processing industries in order to evaluate polymeric fluid flows. A successful computational code must provide reliable predictions (modeling) in a fast and efficient way (simulation). In this work, a new approach to solve the governing equations of viscoelastic fluid flows is proposed. It is based on the finite-volume method with collocated arrangement of the variables, using high-order approximations for the linear and nonlinear average fluxes in the interfaces and for the nonlinear terms obtained from the discretization of the constitutive equations. The approximations are coupled to the Weighted Essentially Non-Oscillatory (WENO) scheme to avoid oscillations in the solution. The Oldroyd-B model is used to describe the rheological behavior of the viscoelastic fluid. The average values of the variables in the volumes are used during the resolution, and the point values are recovered in the post-processing step by deconvolution of the average values. The nonlinear system, resulting from the discretization of the equations, is solved simultaneously using a Newton-like method. The obtained solutions are oscillation-free and accurate, demonstrated by the application on a classic problem in computational fluid dynamics, the slip-stick flow.

Non-Newtonian Fluids; Viscoelasticity; Finite-volume method; High-order interpolation schemes; WENO scheme


PROCESS SYSTEMS ENGINEERING

High-order finite volume method for solving viscoelastic fluid lows

A. R. Muniz; A. R. Secchi* * To whom correspondence should be addressed ; N. S. M. Cardozo

Escola de Engenharia - Departamento de Engenharia Química, Universidade Federal do Rio Grande do Sul, Phone: +(55) (51) 3316-3528; Fax: +(55) (51) 3316-3277, Rua Sarmento Leite, 288/24 CEP: 90050-170, Porto Alegre - RS, Brazil E-mail: arge@enq.ufrgs.br, E-mail: muniz@enq.ufrgs.br, E-mail: nilo@enq.ufrgs.br

ABSTRACT

Computational Fluid Dynamics (CFD) is widely used by polymer processing industries in order to evaluate polymeric fluid flows. A successful computational code must provide reliable predictions (modeling) in a fast and efficient way (simulation). In this work, a new approach to solve the governing equations of viscoelastic fluid flows is proposed. It is based on the finite-volume method with collocated arrangement of the variables, using high-order approximations for the linear and nonlinear average fluxes in the interfaces and for the nonlinear terms obtained from the discretization of the constitutive equations. The approximations are coupled to the Weighted Essentially Non-Oscillatory (WENO) scheme to avoid oscillations in the solution. The Oldroyd-B model is used to describe the rheological behavior of the viscoelastic fluid. The average values of the variables in the volumes are used during the resolution, and the point values are recovered in the post-processing step by deconvolution of the average values. The nonlinear system, resulting from the discretization of the equations, is solved simultaneously using a Newton-like method. The obtained solutions are oscillation-free and accurate, demonstrated by the application on a classic problem in computational fluid dynamics, the slip-stick flow.

Keywords: Non-Newtonian Fluids; Viscoelasticity; Finite-volume method; High-order interpolation schemes; WENO scheme.

INTRODUCTION

The simulation of polymeric fluid flows can be used for many purposes in engineering, like equipment design, process optimization, developing of new processes, etc. This has been an area of intensive research in last two decades, and many challenges still remain.

Polymeric fluids are classified as non-Newtonian fluids due to its shear-thinning viscosity, but their rheological behavior is more complex compared to this single characteristic. These fluids combine properties of viscous liquids and elastic solids, being called viscoelastic fluids. This behavior leads to occurrence of curious phenomena in some flows, as can be seen in Bird et al. (1987) and Tanner (1999) for example, making their study very interesting.

Two topics must be focused in the simulation of such flows: the mathematical model, consisting of the conservation equations, a constitutive equation, and adequate boundary conditions; and the numerical technique, used to solve the partial differential equations (PDE) system.

For the first topic, the crucial point is to select a reliable constitutive equation. In the literature, there exist a great number of constitutive equations, split into groups according to their mathematical form and predictability of the rheological properties. A good description of the most known constitutive equations for viscoelastic fluids can be found in Bird et al. (1987) and Larson (1988). In the present work, the well-known Oldroyd-B model (Bird et al., 1987), a nonlinear differential model, is used as the constitutive equation. With the nonlinear differential models, it is possible to predict many rheological characteristics of viscoelastic fluids, but they have some limitations with respect to quantitative predictions (Bird et al., 1987; Macosko, 1994). Their wide use is due to the relative facility of implementation, compared to the integral models, and the ability of making many predictions, at least, qualitatively.

In relation to the second topic, there still are many difficulties to get around, mainly when the elasticity of fluid is more prominent, in other words, in flows under high Weissenberg numbers (Tanner, 1999; Keunings, 1986). To solve the governing PDE system, the finite-difference, finite-element, and finite-volume methods are the most used methods. The first papers appeared in the second half of the 70’s, using finite-difference and finite-element methods (Crochet and Pilate, 1976; Perera and Walters, 1977; Kawahara and Takeuchi, 1977; Caswell, 1979). A good review of several finite-element methods for viscoelastic flows can be found in Baaijens (1998). The finite-volume method started to be used in the beginning of the 90’s (Hu and Joseph, 1990; Yoo and Na, 1991), after a great development of solving Newtonian flows during the 70’s and 80’s.

In the finite-volume method, the domain is divided in small control volumes (or cells), creating a mesh, and the governing equations are integrated over each volume, assuring the conservation of all properties in each cell of the mesh. The resulting discretized equation system is solved for the variables in the center of each volume. The differences between the varieties of finite-volume approaches found in the literature arise in relation to the arrangement of the variables, the interpolation schemes, the solution methods of the discretized equations, the methods for treating the pressure-velocity coupling, etc.

In the beginning, the staggered arrangement of variables was exclusively applied (Yoo and Na, 1991; Darwish et al., 1992), but, in more recent works, the collocated arrangement is adopted (Missirlis et al., 1998; Oliveira et al., 1998), due to its advantages relative to the staggered arrangement (Peric et al., 1988).

In most of the approaches, the discretized nonlinear equation systems are solved by a segregated technique, requiring an explicit equation for pressure, not present in the original system of equations for incompressible flows. The usual way is to obtain a Poisson equation for pressure from the continuity equation. There are many forms to create this equation, and consequently, different iterative procedures. The most used are semi-implicit methods based on pressure correction equations, like SIMPLE (Patankar and Spalding, 1972), SIMPLER (Patankar, 1980), and SIMPLEC (van Doormal and Raithby, 1984). For example, the SIMPLE method was used by Luo (1996), while Yoo and Na (1991) used the SIMPLER method. A more robust method than the previous ones is the SIMPLEC, used in more recent works (Oliveira et al., 1998; Alves et al., 2000).

The finite-volume method requires two levels of approximation in the discretization of the equations: the approximation of the integral of the fluxes over the cell faces and the approximation of cell face values in terms of the cell center values. In the classical approach, used in all the papers cited above, the approximations are, at most, of second order. The integrals of the fluxes are expressed in terms of point values over the interface, the second-order approximation being the most used, where the value of the fluxes along the interfaces is equal to the flux in the midpoint of the face. The average value over the volume is equal to the value in the central point, another approximation of 2nd order. The values of the fluxes in the middle point (or more points) of the face are expressed in terms of the cell center values, using interpolation schemes.

The most used interpolation schemes are the upwind differencing scheme (UDS) and hybrid differencing scheme (HDS). The upwind scheme is a first-order approximation widely used due to its great stability, but it requires a very fine grid to obtain an accurate solution, because the solution is "contaminated" by the presence of strong numerical diffusion (Patankar, 1980). The hybrid scheme is nothing more than the use of UDS in regions where the convective fluxes have more influence, and the central differencing scheme (CDS) in the rest, so the same problem still remains.

With the use of higher-order schemes, coarser grids can be employed. The CDS, power-law, and LUDS (Harten, 1983) are examples of 2nd order approximation schemes, and the QUICK (Leonard, 1979) is formally a 3rd order approximation scheme. However, as a 2nd order approximation is generally used for the fluxes along the interfaces, the resulting order of approximation is 2. The use of these schemes improves the accuracy of the results, but can cause some undesirable problems. In flows involving discontinuities, more specifically, in regions with large gradients, the solutions obtained are characterized by the presence of unreal oscillations (wiggles), besides the possibility of convergence loss in the numerical procedure. The solution becomes unbounded, and the value at the interfaces can be larger than the neighboring center volume value.

In order to avoid those undesired problems, the above schemes can be coupled with the UDS, creating a high-resolution scheme (HRS). This class of schemes, like MINMOD (Harten, 1983) or SMART (Gaskell and Lau, 1988), obeys a criterion that avoids the presence of unbounded solutions, called Convective Boundedness Criterion (CBC), of Gaskell and Lau (1988). However, the order of approximation still remains limited to the 2nd order.

The use of higher-order approximations, in the integration of the fluxes over the cell faces and/or in the interpolations schemes, allows the use of even coarser grids, reducing the CPU time of a simulation. However, high-order interpolation schemes are very prone to produce instability and oscillatory solutions. Therefore, a desired high-order approach must lead to stable and accurate solutions, which are the focus of the present work.

A finite-volume approach to solve the governing equations of viscoelastic flows is proposed. The approach is based on the work of Kobayashi (1999), which uses cell average values as primitive variables instead of the point values. The point values, if desired, are obtained at the end of calculations by a deconvolution procedure. For high-order approximations, this approach simplifies the numerical procedure by avoiding the need of integration of the flux from a set of point values in reconstruction procedures (Pereira et al., 2001). The average fluxes at the cell faces are expressed in terms of the cell average variables, using adequate high-order interpolation schemes, obtained through Taylor series expansion. Special treatment must be given to the nonlinear terms, to maintain the global high-order approximation, as it will be seen later.

The resulting equations are solved simultaneously, using a Newton-like method. To solve the large-scale linear systems resulting from the Newton’s iterations, an iterative method, the GMRES method with ILU preconditioning (Golub and Charles, 1996) is used. The matrix sparsity of the linear systems is also taken into account to reduce the amount of memory and CPU time needed.

In the next sections, the developed methodology is presented. The governing equations of a viscoelastic fluid flow are shown in section 2. In section 3, the proposed approach, the interpolation schemes, and the solution method of the discretized equations are presented. In section 4, some numerical results are presented and discussed.

GOVERNING EQUATIONS

he governing equations of an isothermal and incompressible fluid flow are the mass and momentum conservation equations, being given below, in their dimensionless form:

where all the variables are normalized by characteristic values, being the velocity vector, the stress tensor, p the pressure, and Re the Reynolds number. The characteristic values are L for length and V for velocity, defined according to the case studied. For stress and pressure, the characteristic values are p0 = t0 = h0V / L. The Reynolds number is given by Re = r V L / h0, r being the specific mass, and h0 the viscosity.

The stress tensor can be split into two parts: a viscous stress tensor ( and an elastic stress tensor (N), as follow:

This split is made because many differential constitutive equations for viscoelastic fluids can be written in this form, and also due to the insertion of diffusive terms in the momentum equations, stabilizing the numerical procedure. This technique, or a similar one, is used in the simulation of viscoelastic flows by many numerical methods, like EVSS - Elastic-Viscous Split Stress method (Rajagopalan et al., 1990) among others.

The viscous stress tensor is given by Newton’s law of viscosity:

where hN is the viscous contribution to the total viscosity h0 and is the shear rate tensor for incompressible flow.

For the elastic stress tensor, an appropriate constitutive equation for viscoelastic fluids must be used. In this work, the Oldroyd-B model (Bird et al., 1987) is used because it can produce high and even infinity values for elongational viscosity at finite elongational rates (Muniz et al., 2004). This is the main reason for the difficulty in obtaining solutions under high We numbers using this model, making it appropriate to test the stability of numerical methods. The dimensionless form for this model is given by:

where We is the Weissenberg number, hE = hP/h0 is a dimensionless parameter of the constitutive equation that represents the contribution of elastic to total viscosity h0 = hN + hP, and P(1)is the advective derivative of the stress tensor, defined in Equation 6. The Weissenberg number is defined as We = lV/L, where l is the relaxation time.

According to the above definitions and knowing that hE + hV = 1, where hV = hN/h0, the momentum conservation equation (Equation 2) becomes:

In this work, for simplicity, only two-dimensional problems in Cartesian coordinates are presented. Obviously, the methodology can be extended to three-dimensional problems and other coordinate systems.

NUMERICAL METHOD

The finite-volume method, using collocated arrangement of variables, is adopted to discretize the governing PDE system. The details of each step of the methodology are given below.

System Discretization

The discretization procedure by the finite-volume method consists of integrating the partial differential equations over the control volumes, leading to algebraic equations (for steady-state problems) or differential-algebraic equations (for transient problems). The discretized equations are written in terms of cell face fluxes and cell average values. The cell face fluxes f in the x and y directions are defined as, respectively:

The cell average value of a quantity f is defined as:

To exemplify this procedure, two discretized equations of the studied system in the steady state are shown: the conservation of momentum in the x direction, Equation 10, and the constitutive equation for the txx component of the stress tensor, Equation 11.

In these equations, there are different types of cell face fluxes: linear and nonlinear convective fluxes and diffusive fluxes, besides some volume averaged products of variables and derivatives. These terms are approximated by the procedure described in the following section.

Approximation of Cell Face Fluxes (Interpolation Schemes)

In the present approach, derived from Kobayashi (1999) and Pereira et al. (2001), the primitive variables in the problem are the cell average values of the former primitive variables. The fluxes over the cell faces and the cell average of nonlinear terms are expressed as a function of the cell average variables in the discretized equations.

Basically, there are three kinds of fluxes to approximate: linear and nonlinear convective fluxes and diffusive fluxes. The schemes used are all derived for uniform spaced grids. The index notation shown in Figure 1 is used to describe the interpolation schemes, valid to any row in the domain.


The interpolation schemes are obtained by proposing an adequate stencil and writing the fluxes in terms of cell average variables multiplied by coefficients to be determined. Taylor series expansion is applied to each term in the scheme, for a given order, according to the definitions given by Equations 8 and 9, for fluxes over the cell faces and cell average values. The coefficients are determined by taking the identity between the two sides of the resulting expression. For more details, see Pereira et al. (2001), that show this procedure to obtain 4th order compact schemes.

To evaluate the convective fluxes, taking the flux for example, a Lagrange 3rd order interpolation scheme is used, in other words, an upwinded 3rd order approximation, which takes two volumes downstream and one volume upstream from the face. Then, the expression depends on the flow direction:

To evaluate the diffusive fluxes , a centered 4th order approximation is used:

Due to its high order, the use of the interpolation scheme given by Equation 12, for the approximation of convective fluxes, can lead to the undesired problems discussed in section 1. Therefore, the 3rd order scheme, using different possible stencils, was used to create a WENO – Weighted Essentially Non-Oscillatory scheme (Liu et al., 1994). This kind of scheme is widely used in problems using finite-difference methods, but not so much explored in the finite-volume context (Liu et al., 1994; Qiu and Shu, 2002).

A WENO scheme consists of a convex combination of different schemes of same order of approximation, but with different stencils, leading to a higher-order scheme. Its main characteristic is to avoid the presence of unreal oscillations in the solution in regions of steep gradients, or discontinuities. When discontinuities are present, the weight of each scheme is dependent on the location of the points of each stencil in relation to the points of the discontinuity, being higher the "smoother" the stencil.

Using 3rd order approximations for the convective fluxes at a generic interface i parallel to y, like Equation 12, considering vx > 0, three (r = 3) possible stencils lead to three different schemes fk:

Then, the flux at interface i can be a convex combination of these schemes, like:

where wk is the weight of each scheme in the resulting scheme (Liu et al., 1994), defined by:

Ck is called optimal weight, ISk is a smoothness measurement of the stencil k, and e is a small positive constant to avoid division by zero, 10-6 being the value used in all calculations in this work. The optimal weights, according to Jiang and Shu (1996), are {C0 = 1/10, C1 = 6/10, C2 = 3/10}. The expressions of ISk for each stencil considered, used in this work, were derived by Jiang and Shu (1996), given below:

The resulting WENO scheme leads to a 5th order approximation for the flux , in smooth regions (Jiang and Shu, 1996).

The approximation is taken line by line in the 2D domain. It is possible to derive more adequate WENO schemes for this domain, where the stencil could take into account volumes in both directions. An example is given by Shi et al. (2002), where it was concluded that the results obtained by the two paths lead to similar solutions, reminding that the two-direction scheme is computationally more expensive compared to the line by line approximation.

To finalize the step of approximations, the nonlinear terms that consist of products of variables at the cell faces, and the products of cell variables, averaged over the volume, must be treated with special attention. Examples for the first case are the inertia terms in the momentum equation, and analogs and, for the second case, the product between velocity and stress terms in the constitutive equations, like .

In the classical approach, the approximation of a product of variables is given by the product of the approximations of each term, resulting in 2nd order approximations:

If these approximations are used, even using high-order approximations for cell face fluxes, the global order of approximation is of 2nd order. To obtain a global high-order scheme, a higher-order approximation for the nonlinear fluxes must be employed too.

The approximation used, based on the work of Pereira et al (2001), is obtained by comparing the Taylor series expansions of the right and left sides of Equations 23 or 24, for a given order. For the nonlinear fluxes at interfaces, using a 4th order approximation, the resulting terms are:

where x means the normal direction to the interface i.

For the nonlinear products averaged over the volume, the terms are:

An approximation of the terms on the RHS of Equations 25 and 26 is needed, and they are expressed in terms of the average cell values in the neighborhood. For example, the resulting expression for Equation 25 is given below.

It must be emphasized that the present approach can be easily extended to more complex and reliable differential constitutive equations, like Phan-Thien-Tanner (Phan-Thien and Tanner, 1977) and Giesekus (1982). The terms that arise in the discretization of these equations are totally analogous to those approximated in the present work..

Boundary Conditions

The boundary conditions are naturally imposed in the discretized equations. When the values of the fluxes are known at the boundary, they are replaced directly in the equation. For example, zero velocity at the walls, zero shear stress at the symmetry lines, among other examples. In situations where the boundary value must be calculated from internal values, interpolation schemes are used, with the same order of approximation of internal cell faces, to maintain the global order. This is valid for all types of fluxes to be approximated.

Solution of the approximated equations

The resulting system of nonlinear equation is solved using a Newton-like method. This method consists of solving a nonlinear system F(x) = 0, through the following iterative procedure:

where xk+1 = xk + dx, xk is the variables vector, at the k iteration, and J(xm) is the Jacobian matrix of the system, evaluated at some previous iterate xm. The Jacobian matrix is calculated by perturbation in finite differences, being kept fixed for a given number of iterations:

The sparsity of the Jacobian matrix is taken into account in the numerical procedure. The linear system is solved using the GMRES method (Saad and Schultz, 1986), which is a Krylov subspace based method, appropriate for large-scale linear systems. The ILU preconditioner is used to improve the convergence of this iterative method (Golub and Charles, 1996).

Deconvolution of Average Values

With the solution of the discretized system of algebraic equations, the cell average values are obtained. Generally, the values of the cell variables at specified points are of interest, and can be obtained by applying a deconvolution procedure to the average values. There are many ways to perform a deconvolution, as for example, to obtain the values at the cell vertices from average values at cell faces using an adequate high-order approximation (Pereira et al., 2001). In this work, the same procedure was adopted. An important point is that the expression used should maintain the order of the approximations used in the numerical procedure.

NUMERICAL RESULTS AND DISCUSSION

The problem solved is the slip-stick flow or entry flow. It is a problem in fluid mechanics frequently employed to test numerical methods and interpolation schemes, due to the presence of a singularity in a boundary condition (slip to no-slip). This singularity causes theoretically an infinite stress at the junction, resulting in profiles of the stress variables with very sharp shape in this region. An illustration of this geometry is shown in Figure 2. Dealing with this type of singularity is important because this is a common feature in many polymer processing flows, like flow in an abrupt contraction or abrupt expansion, exit flow in extrusion, among many others. In the simulations, only the lower half of the shown domain was taken into account, the origin point being located at the slip-stick junction.


The slip-stick problem consists of parallel plate entry flow of a fluid coming from a surface with slip at the walls. The boundary conditions for this problem are a uniform velocity profile and zero stresses at the input, developed profiles at the output, symmetry conditions at the upper and lower boundary of the y-domain before the singularity, and no-slip conditions at the lower boundary after the origin.

In the situation considered in this example, the lengths shown in Figure 1 are H = 1, L1 = 5, and L2 = 7. The values of the dimensionless parameters are Re = 1, We = 0.2, and hE = 0.5.

The problem was solved using two distinct groups of approximations: the first one, referred to in the text as LAG3, uses 3rd order Lagrange interpolation for convective fluxes, 4th order Lagrange for diffusive fluxes, and 4th order approximation for non-linear terms. The second group, referred to as WENO, uses the same approximations mentioned above, together with the WENO scheme for the convective fluxes. Different meshes were used in the calculations, as can be seen in Table 1.

The simulations involving WENO approximations are computationally more expensive, and the simulations using LAG3 present solutions with oscillations near the singularity, bringing instability to the numerical procedure for viscoelastic flows. Far from this critical region, the solutions obtained by both approximations have the same behavior.

In Figure 3, the flow streamlines are shown. Note that the fluid near the entrance wall has the tendency to distribute over the channel height, leading to large variations near the slip-stick junction.


The behavior observed above can be better seen visualizing the variable profiles near to zero on y-axis (near to the singularity). These profiles will be used to observe many characteristics of the solutions obtained by the two groups of approximations previously cited. All the profiles showed hereafter are in the dimensionless values.

In Figure 4, the x-velocity and y-velocity profiles in the horizontal line y = 0.05 obtained by both groups of approximations in mesh M1 are showed.


The x-velocity profile is relatively smooth, and both the approximations lead to similar solutions, oscillation free. However, the use of LAG3 causes undershoots and oscillations in the y-velocity profile near to the singularity, as seen in Figure 4 (b). The use of the WENO scheme eliminates this problem, as seen in the same figure. This behavior also occurs in stress profiles. The shear stress profiles in the horizontal line at y = 0.075, obtained in mesh M5 by both groups of approximations, are shown in Figure 5 (a), and the region near to the singularity is pointed out in Figure 5 (b).


In the other stress profiles, obtained in mesh M5, seen in Figure 6, oscillations are not present, but the solution obtained by the LAG3 presents higher overshoots compared to the WENO solution.


A question that arises when using centered collocated variable arrangement is the possibility of generation of oscillations in the pressure profiles, a common situation in approaches using this kind of arrangement (Peric et al., 1988; Patankar, 1980). Due to the use of high-order interpolation schemes, the pressure drop across a cell depends on the pressure in that cell and on the pressure in the neighboring cells, avoiding the occurrence of this problem. A pressure profile in the horizontal line y = 0.05, obtained in mesh M1, is displayed in Figure 7.


The x-velocity and stress profiles in the horizontal line at y = 0.075, obtained in different meshes, using the WENO scheme, are shown in Figures 8 and 9.



Note that the height of the stress peaks is dependent on the mesh dimension, being higher the more refined the mesh. In regions far from the singularity, the solutions are mesh-independent, after a given mesh size, including horizontal neighbor lines. For the x-velocity profile, the solutions are practically the same, because they do not have large variations like the stress profiles.

Another point to be commented in the obtained solutions is the behavior of the stress profile in the region near to the singularity. In the work of Eggleton et al. (1996), solving the same problem for UCM fluid, the obtained profiles are contaminated by oscillations, even using low-order approximations. Oliveira et al. (1998), solving the same problem, also using the UCM model, point to the possible cause of this behavior as the inadequate approximation of the boundary conditions at the wall, and suggest a correct way to do it, calculating the stress values at the wall by implicit equations derived from the constitutive equations. In the present approach, this is not necessary, due to the high order nature of the approximations, the boundary conditions being easily imposed in an explicit form. The profiles near to the singularity, obtained with the use of the mesh M7 are shown in Figures 10 and 11.



In Figures 10 and 11, a non-oscillatory solution can be observed in the region near to the singularity. The large variation of the variables can be seen again in these profiles, the gradients being larger closer to the singularity, and being more pronounced on the stick side of the flow.

CONCLUSIONS

In this paper, a new approach to simulate viscoelastic flows in Cartesian meshes was presented. The system of equations was discretized by the finite-volume method, using the collocated arrangement for variables, and high-order approximations for the fluxes over cell faces. All the procedure is based on the variable cell averages and, at the end of calculations, the point values are recovered, if needed, by a deconvolution formula. In order to illustrate the proposed approach, an example of application was presented, showing that the solutions obtained are accurate and oscillation-free near to the discontinuities. All the schemes were derived for uniform meshes, it being necessary to use high mesh sizes to capture the solution at the discontinuity.

ACKNOWLEDGEMENTS

The authors gratefully acknowledge the financial support of CNPq.

(Received: July 23, 2006 ; Accepted August 8, 2007)

  • Alves, M.A., Pinho, F.T., Oliveira, P.J.. Effect of a high-resolution differencing scheme on finite-volume predictions of viscoelastic flows, Journal of Non-Newtonian Fluid Mechanics, 93, 287-314 (2000).
  • Baaijens, F.P.T. Mixed finite element methods for viscoelastic flow analysis: a review, Journal of Non-Newtonian Fluid Mechanics, 79, 361-385 (1998).
  • Bird, R.B., Armstrong R.C., Hassager. 0. Dynamics of Polymeric Liquids Vol. 1 - Fluid Mechanics, John Wiley (1987).
  • Caswell, B. Finite element method for fluids with memory, Journal of Non-Newtonian Fluid Mechanics, 5, 199 (1979).
  • Crochet, M.J., Pilate, G. Plane flow of a fluid of second grade through a contraction, Journal of Non-Newtonian Fluid Mechanics, 1, 247258 (1976).
  • Darwish, M.S., Whiteman, J.R., Bevis, M.J. Numerical modelling of viscoelastic liquids using a finite-volume method, Journal of Non-Newtonian Fluid Mechanics, 45, 311337 (1992).
  • Eggleton, C.D., Pulliam, T.H., Ferziger, J.H. Numerical simulation of viscoelastic flow using flux difference splitting at moderate Reynolds numbers, J. Non-Newtonian Fluid Mechanics, 64, 269-298 (1996).
  • Gaskell, P.H., Lau, A.K.C. Curvature-compensated convective transport: SMART, a new boundedness-preserving transport algorithm, Int. J. Num. Methods Fluids, 8, 617-641 (1988).
  • Giesekus, H. A simple constitutive equation for polymer fluids based on the concept of deformation-dependent tensorial mobility, Journal of Non-Newtonian Fluid Mechanics, 11, 69-109 (1982).
  • Golub, G.H, Charles, F.V.L. Matrix Computation, Johns Hopkins, London (1996).
  • Harten, A..J. High Resolution Schemes for hyperbolic conservation laws, J. Computational Physics, 49, 357393 (1983).
  • Hu, H.H., Joseph, D.D. Numerical simulation of viscoelastic flow past a cylinder, J. Non-Newtonian Fluid Mechanics, 37, 347377 (1990).
  • Jiang, G.S., Shu, C.W. Efficient Implementation of Weigthed ENO Schemes, Journal of Computational Physics, 126, 202-228 (1996).
  • Kawahara, M., Takeuchi, N. Mixed Finite element method for analysis of viscoelastic fluid flow, Comput. Fluids, 5, 33-45 (1977).
  • Keunings, R. On the high Weissenberg number problem, J. Non-Newtonian Fluid Mechanics, 20, 209-226 (1986).
  • Kobayashi, M.H. On a class of Padé Finite Volume methods, Journal of Computational Physics, 156, 137-180 (1999).
  • Larson, R.G. Constitutive equations for polymer melts and solutions, Butterworths, Boston (1988).
  • Leonard, B.P. A stable and accurate convective modelling procedure based on quadratic upstream interpolation, Comput. Methods Appl. Mech. Eng, 19, 59-98 (1979).
  • Liu, X., Osher, S., Chan, T. Weighted Essentially Non-oscillatory Schemes, J. Comp. Physics, 115, 200-212 (1994).
  • Luo, X.L. A control volume approach for integral viscoelastic models and its application to contraction flow of polymer melts, Journal of Non-Newtonian Fluid Mechanics, 64, 173-189 (1996).
  • Macosko, C. Rheology: Principles, Measurements and Applications, VCH Publishers (1994).
  • Missirlis, K.A., Assimacopoulos, D., Mitsoulis, E. A finite volume approach in the simulation of viscoelastic expansion flows, J. Non-Newtonian Fluid Mechanics, 78, 91118 (1998).
  • Muniz, A.R., Secchi, A.R. Cardozo, N.S.M. Simulation of Viscoelastic Fluid Flow in Geometry Involving Singularities using Differential Constitutive Equations, Proceedings of the PPS 2004 - Polymer Processing Society - Americas Regional Meeting, Florianópolis, SC, 2p. CDROM #PPS-10-013 (2004).
  • Oliveira, P.J., Pinho, F.T., Pinto, G.A. Numerical simulation of non-linear elastic flows with a general collocated finite-volume method, J. Non-Newtonian Fluid Mechanics, 79, 143 (1998).
  • Patankar, S.V., Spalding, D.B. A calculation procedure for heat, mass and momentum transfer in three-dimensional parabolic flows, Int. J. Heat Mass Transfer, 15, 1787-1806 (1972).
  • Patankar, S.V. Numerical Heat Transfer and Fluid Flow, McGraw-Hill, New York (1980).
  • Pereira, J.M.C., Kobayashi, M.H., Pereira, J.C.F. A fourth-order-accurate finite volume compact method for the incompressible Navier-Stokes solutions, Journal of Computational Physics, 167, 217-243 (2001).
  • Perera, M.G.N., Walters, K. Long-range memory effects in flows involving abrupt changes in geometry. Part I. Flows associated with L-shaped and T-shaped geometries, J. Non-Newtonian Fluid Mechanics, 2, 4981 (1977).
  • Peric, M., Kessler, R., Scheurer, G. Comparison of finite-volume numerical methods with staggered and collocated grids, Computers & Fluids, 16, No. 4, 389-403 (1988).
  • Phan-Thien, N., Tanner, R.I. A new constitutive equation derived from network theory, Journal of Non-Newtonian Fluid Mechanics, 2, 353-365 (1977).
  • Qiu, J., Shu, C.W. On the Construction, Comparison, and Local Characteristic Decomposition for High-Order Central WENO Schemes, Journal of Computational Physics, 183, 187209 (2002).
  • Rajagopalan, D., Armstrong, R.C., Brown, R.A. Finite element methods for calculation of steady, viscoelastic flow using constitutive equations with a Newtonian viscosity, Journal of Non-Newtonian Fluid Mechanics, 36, 159-192 (1990).
  • Saad, Y., Schultz, M.H. GMRES: a generalized minimal residual algorithm for solving non-symmetric linear systems, SIAM J. Sci. Stat. Comp., 7, 856-869 (1986).
  • Shi, J., Hu, C., Shu, C.W. A Technique of Treating Negative Weights in WENO Schemes, Journal of Computational Physics, 175, 108127 (2002).
  • Tanner, R.I. Engineering Rheology, Oxford University Press, Oxford (1999).
  • van Doormal, J.P., Raithby, G.D. Enhancements of the SIMPLE method for predicting incompressible fluid flows, Num. Heat Transfer, 7, 147-163 (1984).
  • Yoo, J.Y., Y. Na. A numerical study of the planar contraction flow of viscoelastic fluids using the SIMPLER algorithm, J. Non-Newtonian Fluid Mechanics, 39, 89 (1991).
  • *
    To whom correspondence should be addressed
  • Publication Dates

    • Publication in this collection
      28 Apr 2008
    • Date of issue
      Mar 2008

    History

    • Accepted
      08 Aug 2007
    • Received
      23 July 2006
    Brazilian Society of Chemical Engineering Rua Líbero Badaró, 152 , 11. and., 01008-903 São Paulo SP Brazil, Tel.: +55 11 3107-8747, Fax.: +55 11 3104-4649, Fax: +55 11 3104-4649 - São Paulo - SP - Brazil
    E-mail: rgiudici@usp.br