Log-conformation representation of hiperbolic conservation laws with source term

Log-conformation representation of hiperbolic conservation laws with source term

L.A. Bielça^I; M. Meneguette^{II, *}

^IMSc in Computacional and Applied Mathematics, FCT/UNESP, Universidade Estadual Paulista, 19060-900 Presidente Prudente, SP, Brazil. E-mail: [e-mail]

^IIDepartment of Mathematics and Computing, FCT/UNESP, Universidade Estadual Paulista, 19060-900 Presidente Prudente, SP, Brazil. E-mail: [e-mail]

ABSTRACT

The objective of this work is to study, through a simple equation, the statement that the numerical instability associated to the high Weissenberg number in equations with source term can be resolved by the use of the so called logarithmic conformation representation. We will focus on hyperbolic conservation laws, but more specifically on the advection equation with a source term. The source term imposes a necessity of an elastic balance, as well as the CFL convective balance for stability. Will be seen that the representation of such equation by the log-conformation removes the restriction of stability inherent to the elastic balance pointed out by Fattal & Kupferman [3] as the cause of the high Weissenberg number problem (HWNP).

Keywords: source term, log-conformation representation (LCR), high Weissenberg number problem (HWNP).

RESUMO

O objetivo deste trabalho é estudar, por meio de uma equação simples, a afirmação de que a instabilidade numérica associada ao elevado número Weissenberg em equações com termo fonte, pode ser resolvida pela utilização da denominada representação por conformação logarítmica. Vamos nos concentrar nas leis de conservação hiperbólicas, mas mais especificamente na equação de advecção com termo fonte. O termo fonte impõe a necessidade de um equilíbrio elástico, tanto quanto um equilíbrio convectivo (CFL), para a estabilidade. Será visto que a representação de tal equação pela conformação logarítmica remove a restrição de estabilidade inerente para o equilíbrio elástico apontado por Fattal & Kupferman [3] como causa do problema do alto número de Weissenberg.

Palavras-chave: termo fonte, representação por conformação logarítmica, problema do alto número de Weissenberg.

1 INTRODUCTION

This work explores some important aspects of the numerical treatment, as well as the comparisons with the exact solution, of a simple partial differential equation with a source term, in analogy to what happens in the simulation of viscoelastic fluid flows [1,3]. The search for robust techniques is still an important object of research, mainly because the most complete equations lead to new developments; this is the case of an equation with a source term. A naive approach could lead one to believe that the equations without source and with a source are alike and, therefore, every numerical method designed to the first would easily be applied to the second. In fact, this is not always the case. In a wider context, one case is the synthesized situation by the Weissenberg (Wi) number. In short, the higher is Wi, the more intense the interactions of elastic nature and the more decisive is the source term contribution for the flow. Classical methods usually bump into restrictions regarding on the magnitude of Wi, about which is still to be discovered the cause: if it is of the numerical nature or, even, if intrinsically to the viscoelastic model itself. This situation has troubled, but also arose the interest of many researchers in the last decade, constituting the so-called High Weissenberg Number Problem (HWNP). Recently, through the logarithmic conformation representation (LCR), there has been a better indication that the cause might be numerical [3]. In the same way that the numerical solution without source requires the correct convection balance, the solution with source needs to take into account a correct balance between the terms of convective nature (that generates the hyperbolic character) as well as those of elastic source (that introduces a stiff character in the equation).

In this work, we aim to effectively analyse and reproduce the mechanisms associated to the HWNP through the study of a simplified equation. A simple equation that allows the broad study of these aspects is the advection equation with a source term.

Giving a hyperbolic equation in conservative form, one can consider a source term so that the equation is still hyperbolic. However, not all extensions of the numerical methods designed for the first case are adequate to the new equation.

The advection equation is a hyperbolic equation, whose numerical solution must respect the CFL condition for the correct balance of convection [2,5].

Let us consider the advection equation with a source term

where

u = u(x, t), x ∈[0, L], t > 0, b = b(x)> 0, a = a(x)

is the speed of advection and Wi is the Weissenberg number. This hyperbolic equation when approximated by finite difference methods, as much as the balance imposed by the convective CFL condition, requires an analysis for the elastic balance resulting from the addition of the source term.

The analytic solution of equation (1.1), for unitary initial condition, can be deducted as

where c = b(x) –. For at < x < L the solution is actually constant.

We will see that the restriction of stability inherent to the elastic balance, pointed by [3] as the cause of the HWNP, when solving (1.1) can be removed by the log-conformation representation.

2 LOG-CONFORMATION REPRESENTATION

In order to clarify the instability problem mentioned above, we detailed and tested the proposed methodology in [3]. Without loss of generality, in (1.1), we will consider a(x) =a > 0. Using in (1.1) the first-order Upwind Method (for the convective term) and Euler explicit (for the transient term), we obtain the scheme

and then

We consider U_{i, j} = U(x_i, t_j) as being the solution provided by the numerical scheme in (x_i, t_j), with δx the spatial step, δt the time step and b_i = b(x_i).

Then the numerical method (2.2) will be stable if

and

Observe that, for (2.4), it is enough

or

Thus, besides the restriction CFL (2.3), we have (2.6) that is a restriction over the spatial step of the mesh, imposed by the elastic balance. Notice that (2.5) and (2.6) are affected by the Weissenberg number: in (2.6), we can see that, the higher is Wi the smaller must be δx .On the other hand, (2.5) shows the direct relation between Wi e b. So, there is a maximum Wi permitted. These restrictions have consequences when a viscoelastic simulation is being carried out, mainly near stagnation points, for example.

Next, we will see that it is possible to remove the restriction (2.6) using the representation by log-conformation. That is based on the following statement:

The representation of a partial equation by log-conformation consists in replacing a range of unknowns in the partial equations by logarithmic unknowns. Thus, to represent a partial equation whose unknown is u, one must make the change of variables

where

Introducing (2.8) into (1.1) we obtain

that corresponds to LCR version of equation (1.1).

We discretized (2.9) by finite differences, considering Ψ_{i, j}= Ψ(x_i, t_j) as being a solution provided by the numerical scheme at the point (x_i, t_j). So, approximating (2.9) by following (2.1), we obtain the following scheme:

Now, one can see that the numerical solution provided by such scheme will be stable when

Then, the advection equation with a source term under the log-conformation representation does not impose restrictions for stability on the spatial step δx.

It is clear that a numerical break down will happen when the stability restrictions are not observed.

On the other hand, we will present some numerical results, all within the stability range, that ilustrate another aspect still in analysis, ie. given a value for t and fixed Wi, the numerical solution (without LCR), seems to approach the constant part of the analytical solution only asymptotically, even for higher order methods, that is for δt → 0, while for all parameters the simple LCR works well.

To ilustrated the well behaviour of the LCR against other methods we consider three diferent well known methods: upwind (first order), superbee (TVD – 2nd order) and Koren (TVD – 3rd order).

3 NUMERICAL RESULTS

Next, we present numerical results obtained using MATLAB to calculate the exact solution of the equation (1.1), its numerical upwind solution (2.2), its LCR version (2.10) and also two TVD schemes (Koren Limiter and Superbee, cf. [4,2]).

3.1 Influence of δ _t

Taking a = 1,b = 2, Wi = 100, L = 4.5,δx = 0.1and varying δt, Figure 1 represents the perfil, at t = 2, of the numerical solution of the equation (1.1) without LCR and with LCR. The solution provided by the LCR (2.10) always follows the growth of the exact solution. As expected the method with LCR suffers a dissipation at the contact discontinuity point of the solution; the same occurs with the other methods as well. Note that in cases without LCR the numerical solution only follows the adequate increasing of the exact solution as δt decreases.

Table 1 shows the relative error (norm 2) at t = 2. Note that the error by norm 2 , for being global, reflects the dissipation, but for the interval where the solution is constant, the LCR is always precise: this can be seen by the absolute error at x = 4. Also, the TVD methods get better near the contact discontinuity as δt diminishes.

3.2 Influence of Wi

Taking a = 1,b = 2,L = 4.5,δx = 0.1and δt = 0.05, Table 2 presents the relative error, at t = 2, by the numerical schemes for Wi = 10 and Wi = 500. The numerical solutions of the equation (1.1) when we take Wi = 10,500 are given in Figure 2. It is worth mentioning that, for other values of Wi, the pictures obtained are all "similar". It is seen in Table 2 that the absolute error, for a fixed δt gets worse when Wi increases, but the LCR always work well.

4 CONCLUSION

By using first order upwind in the convective term and explicit Euler in the transient term, we note that when including a source term in the advection equation it is also "inserted" a restriction of stability (2.6) to the spatial mesh that is influenced by the Weissenberg number. Thus, the higher Wi is, much more refined must be the spatial mesh. This situation is related to the HWNP.

Also, corroborating with [3] and showing in a very concrete way, through a simple problem with exact solution, we can see that the numerical solution by methods without LCR does not follow the adequate growth of the analytical solution. This will certainly have consequences when a viscoelastic simulation is carried out. As δt → 0 there will be convergence and the higher order TVD schemes will approach the contact discontinuit better; upwind and LCR present strong dissipation near the contact discontinuity as expected. Numerical results verified these assumptions.

ACKNOWLEDGEMENTS

Thanks to FAPESP, for the financial support; Paper presented at CMAC – Sudeste in 2013.

Received on November 30, 2013

Accepted on November 9, 2014

• [1] A.M. Afonso, F.T. Pinho & M.A. Alves. The kernel-conformation constitutive laws. J. Non-Newtonian Fluid Mech., 167168 (2012), 3037.
• [2] J.A. Cuminato & M. Meneguette. Discretização de equações diferenciais parciais: técnicas de diferenças finitas. SBM, (2013).
• [3] R. Fattal & R. Kupferman. Time-dependent simulation of viscoelastic flows at high Weissenberg number using the log-conformation representation. J. Non-Newtonian Fluid Mech., 126 (2005), 2337.
• [4] B. Koren. A robust upwind discretisation method for advection, diffusion and source terms, In: C.B. Vreugdenhil and B. Koren (Eds.), Numerical Methods for Advection-Diffusion Problems (pp. 117 138). Vieweg, (1993).
• [5] R.J. LeVeque. Numerical methods for conservation laws. Birkhauser, (1992).

*

Corresponding author: Messias Meneguette

Publication Dates

• Publication in this collection
  09 Feb 2015
• Date of issue
  Dec 2014

History

• Received
  30 Nov 2013
• Accepted
  09 Nov 2014