Efficient numerical computations of yield stress fluid flows using second-order cone programming

https://doi.org/10.1016/j.cma.2014.10.008Get rights and content

Highlights

  • A novel optimization method for the simulation of visco-plastic fluid flows is presented.

  • The minimization problem arising in the variational formulation is formulated as a second-order cone program.

  • It is solved using very efficient interior point solvers compared to currently available methods.

  • The method is validated on benchmark examples and confronted to some experimental results.

Abstract

This work addresses the numerical computation of the two-dimensional flow of yield stress fluids (with Bingham and Herschel–Bulkley models) based on a variational approach and a finite element discretization. The main goal of this paper is to propose an alternative optimization method to existing procedures such as penalization and augmented Lagrangian techniques. It is shown that the minimum principle for Bingham and Herschel–Bulkley yield stress fluid steady flows can, indeed, be formulated as a second-order cone programming (SOCP) problem, for which very efficient primal–dual interior point solvers are available. In particular, the formulation does not require any regularization of the visco-plastic model as is usually the case for existing techniques, avoiding therefore the difficult choice of the regularization parameter. Besides, it is also unnecessary to adopt a mixed stress–velocity approach or discretize explicitly auxiliary variables as frequently proposed in existing methods. Finally, the performance of dedicated SOCP solvers, like the Mosek software package, enables us to solve large-scale problems on a personal computer within seconds only. The proposed method will be validated on classical benchmark examples and used to simulate the flow generated around a plate during its withdrawal from a bath of yield stress fluid.

Introduction

Yield stress fluids are encountered in a wide range of applications: toothpastes, cement, mortar, foams, muds, mayonnaise, etc. The fundamental character of these fluids is that they are able to flow (i.e. deform indefinitely) only if they are submitted to a stress larger than a critical value; otherwise they deform in a finite way like solids. Here, we will focus on the flow characteristics of non-thixotropic yield stress fluids which can, as a good first approximation, be described as simple yield stress fluids, i.e. for which the yield stress and more generally the apparent viscosity are independent of the flow history. Even in that case the flow characteristics of such materials are difficult to predict as they involve permanent or transient solid and liquid regions whose location cannot generally be determined a priori.

Uniform flow of these fluids in simple geometries can easily be described analytically. However there exists a wide set of more complex situations for which analytical description remains strongly approximate or even impossible, so that numerical simulations only can provide useful information or predictions. This concerns stationary flows in channels of varying cross-section such as extrusion, expansion, flow through porous medium, or transient flows such as flows around obstacles, spreading, spin-coating, squeeze flow, elongation, etc.

Various numerical methods have been proposed to simulate the flow of visco-plastic fluids. The main difficulty arising in such simulations is the non-regular constitutive relation between stresses and strain rates which distinguish between a rigid solid domain if the stress state is below the yield stress and a liquid flow otherwise. One major approach consists in regularizing the constitutive behavior by introducing an artificial parameter (which can be viewed as a high viscosity below yield stress e.g. the bi-viscous model) to obtain a smoother problem  [1], [2], [3], [4]. However, the distinction between solid and liquid regions is then much harder to detect and the solution process can deteriorate when changing the regularization parameter. The other main approach avoids any regularization and aims at solving the associated non-smooth variational problem with different optimization techniques, the main example being the Augmented Lagrangian approach  [5], [6], [7]. Mesh adaptive strategies have also been proposed to enhance the prediction of the fluid–solid boundaries  [6]. An interesting review of such techniques can be found in  [8].

The simulation of yield stress fluid flows is actually strongly related to the limit analysis (or yield design in a more general manner) of mechanical structures which aims at computing the limit loads and associated failure mechanisms of rigid perfectly-plastic structures  [9], [10]. More precisely, when viscous effects are dominated by yield effects (Bi) both problems become coincident. In the past decades, computational limit analysis also received an important attention, especially due to the development of mathematical programming tools used to solve the arising optimization problems. Interestingly, regularized models have also first been proposed for limit analysis problems  [11], [12], leading to the same difficulties as those previously mentioned. Linear programming solvers have then been proposed in the case of piecewise linear yield surfaces  [13], [14], [15], [16]. This approach presented the advantage of avoiding any regularization parameter and enabled us to solve a linear optimization problem quite efficiently for medium scale problems. However, its computational efficiency was limited by the important number of additional constraints introduced by the linearization of the non-linear yield surface. Efficient algorithms, namely interior point solvers, have then been developed for linear problems. Their extension to non-linear optimization problems and especially second-order cone programming (SOCP) resulted in an important breakthrough as regards the development of computational capabilities for limit analysis  [17], [18], [19]. These algorithms have also been implemented in commercial codes such as the Mosek software package  [20]. Remarkably, a large number of traditional yield surfaces can be expressed using conic constraints so that limit analysis problems can be formulated within SOCP  [21]. Most of the recent publications in the field of limit analysis are now considering such approaches as the most efficient ones and applying them to a large number of problems: 2D plane strain problems  [22], [23], frame structures  [24], and thin plates in bending  [25], [26].

The aim of this paper is then to show how such methods, successfully used in the limit analysis of mechanical structures, can be transposed to the simulation of yield stress fluids. In particular, it will be shown that it is possible to take both Bingham and Herschel–Bulkley models into account and that the computational performances are very interesting as compared to previously mentioned standard techniques. Section  2 is devoted to the variational formulation of visco-plastic flows. Section  3 will deal with the finite element discretization. The formulation of the discrete problem as a second-order cone program, which is the main novelty of this paper, is considered in Section  4. Finally, Section  5 will present different illustrative examples used to validate the procedure by comparing numerical solutions to analytical ones and experimental measures, as well as to assess its numerical efficiency.

Section snippets

Governing equations

Let Ω denote the two-dimensional fluid domain, u¯=uxe¯x+uye¯y be the fluid velocity field and σ¯¯ the stress tensor field. The creeping flow of an incompressible yield stress fluid (characterized by a yield shear stress τ0) is described by the following set of equations and boundary conditions: divσ¯¯+f¯=0¯divu¯=0d¯¯=12(u¯¯+Tu¯¯)in  Ω{s¯¯=2μd¯¯+2τ0d¯¯d¯¯if  12s¯¯:s¯¯τ0d¯¯=0if  12s¯¯:s¯¯<τ0u¯=U¯don  Ω where d¯¯ is the strain rate tensor, f¯ is the volume body force density, U¯d is the

Finite element discretization

In this section, the finite element discretization of the minimum principle (6) will be considered. In particular, only the velocity field v¯ will be discretized and matrices expressing the strain rate–velocity compatibility equation and the divergence-free condition will be written.

As already mentioned, problem (6) is quite similar to those arising in the upper bound kinematical approach of plane strain limit analysis problems with a von Mises strength criterion, the only difference arising

Standard second-order cone programs

Second-order cone programming (SOCP) is a particular class of convex programming which consists of the minimization/maximization of a linear function of the optimization variable vector xRn, under linear equality and/or inequality constraints and particular non-linear constraints, namely second-order cone (SOC) constraints. Hence, in the case when there are no linear inequality constraints, an SOCP problem can be written as:minxcTxs.t.Ax=bxK where ARm×n,bRm,cRn and K is a second-order

Validation on the plane Poiseuille flow

First, the plane Poiseuille flow of a Bingham and a Herschel–Bulkley fluid is considered and numerical velocity fields are compared to analytical solutions.

Consider a rectangular fluid domain Ω=[0;L]×[0;H] subjected to a uniform horizontal pressure-gradient f¯=fe¯x with a no-slip boundary condition on the surfaces y=0 and y=H. Since there is no fluid motion in the y direction for the exact Poiseuille flow, zero vertical velocities uy=0 have been imposed on the surfaces x=0 and x=L where the

Conclusions

A new optimization technique to solve minimization problems arising in the simulation of yield stress fluid flows has been presented. The method relies on a standard finite element discretization, on the formulation of the discretized minimization problem as a second-order cone program and on the use of dedicated efficient interior-point solvers, like the industrial software package Mosek. The advantage of using such techniques over traditional approaches, like the Augmented Lagrangian for

References (31)

  • M. Bercovier et al.

    J. Comput. Phys.

    (1980)
  • D.D. dos Santos et al.

    J. Non-Newton. Fluid Mech.

    (2011)
  • A. Syrakos et al.

    J. Non-Newton. Fluid Mech.

    (2013)
  • P. Saramito et al.

    Comput. Methods Appl. Mech. Engrg.

    (2001)
  • J. Zhang

    Comput. Methods Appl. Mech. Engrg.

    (2010)
  • E.J. Dean et al.

    J. Non-Newton. Fluid Mech.

    (2007)
  • M.A.A. Skordeli et al.

    Eng. Struct.

    (2010)
  • C.V. Le et al.

    Comput. Struct.

    (2010)
  • J. Gondzio

    European J. Oper. Res.

    (2012)
  • T.C. Papanastasiou

    J. Rheol.

    (1987)
  • M. Fortin et al.

    Méthodes de Lagrangien Augmenté: Applications à la Résolution Numérique de Problèmes aux Limites

    (1982)
  • J. Salençon

    Calcul à la Rupture et Analyse Limite

    (1983)
  • J. Salençon

    Yield Design

    (2013)
  • A. Friaâ

    C. R. Acad. Sci., Paris Sér. A–B

    (1978)
  • T. Guennouni et al.

    J. Méc. Théor. Appl.

    (1982)
  • Cited by (0)

    View full text