Finite Element Approximation of the FENE-P Model

TL;DR

This paper develops finite element schemes for the FENE-P model of dilute polymeric fluids, proving stability and convergence in 2D, and establishing the existence of global weak solutions with stress diffusion.

Contribution

It introduces new finite element schemes for the FENE-P model, demonstrating their stability and proving the existence of global weak solutions in two dimensions.

Findings

Both schemes satisfy a free energy bound without time step constraints.

Convergence to global weak solutions is shown for the FENE-P model with stress diffusion in 2D.

Existence of global-in-time weak solutions is established in two spatial dimensions.

Abstract

We extend our analysis on the Oldroyd-B model in Barrett and Boyaval [1] to consider the finite element approximation of the FENE-P system of equations, which models a dilute polymeric fluid, in a bounded domain $D $\subset$ R d , d = 2 or 3$, subject to no flow boundary conditions. Our schemes are based on approximating the pressure and the symmetric conforma-tion tensor by either (a) piecewise constants or (b) continuous piecewise linears. In case (a) the velocity field is approximated by continuous piecewise quadratics ($d = 2$) or a reduced version, where the tangential component on each simplicial edge ($d = 2$) or face ($d = 3$) is linear. In case (b) the velocity field is approximated by continuous piecewise quadratics or the mini-element. We show that both of these types of schemes, based on the backward Euler type time discretiza-tion, satisfy a free energy bound, which involves the logarithm of both the conformation tensor and a linear function of its trace, without any constraint on the time step. Furthermore, for our approximation (b) in the presence of an additional dissipative term in the stress equation, the so-called FENE-P model with stress diffusion, we show (subsequence) convergence in the case $d = 2$, as the spatial and temporal discretization parameters tend to zero, towards global-in-time weak solutions of this FENE-P system. Hence, we prove existence of global-in-time weak solutions to the FENE-P model with stress diffusion in two spatial dimensions.

Equations

$$\text{Re}\left(\frac{\partial \mathbf{u}}{\partial t}+(\mathbf{u}\cdot \nabla )\mathbf{u}\right)$$

$$\mathrm{div}\mathbf{u}$$

$$\frac{\partial \mathbf{\sigma }}{\partial t}+(\mathbf{u}\cdot \nabla )\mathbf{\sigma }$$

$$\mathbf{u}(0,\mathbf{x})$$

$$\mathbf{\sigma }(0,\mathbf{x})$$

$$\mathbf{u}$$

$$A(\mathbf{\varphi }):={\left(1-\frac{\mathrm{tr}(\mathbf{\varphi })}{b}\right)}^{-1}\mathbf{I}-{\mathbf{\varphi }}^{-1}\forall \mathbf{\varphi }\in {\mathbb{R}}_{\mathrm{S},>0,b}^{d\times d}:=\left\{\mathbf{\psi }\in {\mathbb{R}}_{\mathrm{S},>0}^{d\times d}:\mathrm{tr}(\mathbf{\psi })<b\right\}.$$

$$\mathbf{v}\cdot \mathbf{w}$$

$$\mathbf{\varphi }:\mathbf{\psi }$$

$$\nabla \mathbf{\varphi }::\nabla \mathbf{\psi }$$

$$\Vert \mathbf{v}\Vert$$

$$\Vert \mathbf{\varphi }\Vert$$

$$\mathbf{\varphi }{\mathbf{\chi }}^T:\mathbf{\psi }$$

$$\Vert \mathbf{\psi }\mathbf{\varphi }\Vert$$

$$\Vert \mathbf{\varphi }\mathbf{v}\Vert$$

$$\mathbf{\varphi }={\mathbf{O}}^T\mathbf{DO}\Rightarrow \mathrm{tr}\left(\mathbf{\varphi }\right)=\mathrm{tr}\left(\mathbf{D}\right),$$

$$g(\mathbf{\varphi }):={\mathbf{O}}^{T}g(\mathbf{D})\mathbf{O}\Rightarrow \mathrm{tr}\left(g(\mathbf{\varphi })\right)=\mathrm{tr}\left(g(\mathbf{D})\right),$$

$$d^{-1}(\mathrm{tr}(|\mathbf{\varphi }|){)}^2\le {\Vert \mathbf{\varphi }\Vert }^2\le (\mathrm{tr}(|\mathbf{\varphi }|){)}^2$$

$$\left(\mathbf{\varphi }-\mathbf{\psi }\right):g^{\prime }(\mathbf{\psi })\ge \mathrm{tr}\left(g(\mathbf{\varphi })-g(\mathbf{\psi })\right)\ge \left(\mathbf{\varphi }-\mathbf{\psi }\right):g^{\prime }(\mathbf{\varphi })\forall \mathbf{\varphi },\mathbf{\psi }\in {\mathbb{R}}_{\mathrm{S}}^{d\times d},$$

$$\frac{\mathrm{d}}{\mathrm{d}t}\mathrm{tr}\left(g(\mathbf{\varphi })\right)$$

$$\Vert g(\mathbf{\varphi })-g(\mathbf{\psi })\Vert \le g_{\mathrm{Lip}}\Vert \mathbf{\varphi }-\mathbf{\psi }\Vert \forall \mathbf{\varphi },\mathbf{\psi }\in {\mathbb{R}}_{\mathrm{S}}^{d\times d}.$$

$$\mathrm{W}$$

$$\mathrm{S}$$

$${\int }_{\mathcal{D}}\Vert \mathbf{v}{\Vert }^2\mathrm{d}\mathbf{x}\le C_P{\int }_{\mathcal{D}}\Vert \nabla \mathbf{v}{\Vert }^2\mathrm{d}\mathbf{x}\forall \mathbf{v}\in \mathrm{W},$$

$$F(\mathbf{v},\mathbf{\varphi })$$

$$\forall (\mathbf{v},\mathbf{\varphi })\in [L^2(\mathcal{D}){]}^d\times {\mathrm{S}}^{\star }$$

$$F(\mathbf{v},\mathbf{\varphi })\ge \frac{\text{Re}}{2}{\int }_{\mathcal{D}}\Vert \mathbf{v}{\Vert }^2\mathrm{d}\mathbf{x}+\frac{\epsilon }{2\text{Wi}}{\int }_{\mathcal{D}}\mathrm{tr}\left(\mathbf{\varphi }-\ln (\mathbf{\varphi })-\mathbf{I}\right)\mathrm{d}\mathbf{x}\forall (\mathbf{v},\mathbf{\varphi })\in [L^2(\mathcal{D}){]}^d\times {\mathrm{S}}^{\star },$$

$$\frac{\mathrm{d}}{\mathrm{d}t}F(\mathbf{u},\mathbf{\sigma })+(1-\epsilon ){\int }_{\mathcal{D}}\Vert \nabla \mathbf{u}{\Vert }^2\mathrm{d}\mathbf{x}+\frac{\epsilon }{2{\mathrm{Wi}}^2}{\int }_{\mathcal{D}}\mathrm{tr}\left((A(\mathbf{\sigma }){)}^2\mathbf{\sigma }\right)\mathrm{d}\mathbf{x}=⟨\mathbf{f},\mathbf{u}{⟩}_{H_0^1(\mathcal{D})},$$

$${\int }_{\mathcal{D}}\left[\frac{\text{Re}}{2}\frac{\partial \Vert \mathbf{u}{\Vert }^2}{\partial t}+(1-\epsilon )\Vert \nabla \mathbf{u}{\Vert }^2+\frac{\epsilon }{\text{Wi}}{\left(1-\frac{\mathrm{tr}(\mathbf{\sigma })}{b}\right)}^{-1}\mathbf{\sigma }:\nabla \mathbf{u}\right]\mathrm{d}\mathbf{x}$$

$$\displaystyle\hskip 14.45377pt+\frac{\varepsilon}{2{\text{Wi}}}\int_{\mathcal{D}}\Bigg{[}\left(\frac{\partial{\boldsymbol{\sigma}}}{\partial t}+({\boldsymbol{u}}\cdot\boldsymbol{\nabla}){\boldsymbol{\sigma}}\right)+\frac{A({\boldsymbol{\sigma}})\,{\boldsymbol{\sigma}}}{{\text{Wi}}}\Bigg{]}:A({\boldsymbol{\sigma}})\,{\rm d}\boldsymbol{x}$$

$$-\frac{\epsilon }{2\text{Wi}}{\int }_{\mathcal{D}}\left(\left(\nabla \mathbf{u}\right)\mathbf{\sigma }+\mathbf{\sigma }{\left(\nabla \mathbf{u}\right)}^T\right):A(\mathbf{\sigma })\mathrm{d}\mathbf{x}=⟨\mathbf{f},\mathbf{u}{⟩}_{H_0^1(\mathcal{D})}.$$

$$\left(\frac{\partial \mathbf{\sigma }}{\partial t}+(\mathbf{u}\cdot \nabla )\mathbf{\sigma }\right):A(\mathbf{\sigma })=\left(\frac{\partial }{\partial t}+(\mathbf{u}\cdot \nabla )\right)\left(-b\ln \left(1-\frac{\mathrm{tr}(\mathbf{\sigma })}{b}\right)-\mathrm{tr}\left(\ln (\mathbf{\sigma })\right)\right).$$

$$\left(\mathbf{\sigma }{\left(\nabla \mathbf{u}\right)}^T+\left(\nabla \mathbf{u}\right)\mathbf{\sigma }\right):A(\mathbf{\sigma })$$

$$\frac{\partial \mathrm{tr}(\mathbf{\sigma })}{\partial t}+(\mathbf{u}\cdot \nabla )\mathrm{tr}(\mathbf{\sigma })$$

$${\sup F(\mathbf{u}(t,\cdot ),\mathbf{\sigma }(t,\cdot ))}_{t\in (0,T)}+\frac{1-\epsilon }{2}{\int }_{{\mathcal{D}}_T}\Vert \nabla \mathbf{u}{\Vert }^2\mathrm{d}\mathbf{x}\mathrm{d}t+\frac{\epsilon }{2{\mathrm{Wi}}^2}{\int }_{{\mathcal{D}}_T}\mathrm{tr}\left((A(\mathbf{\sigma }){)}^2\mathbf{\sigma }\right)\mathrm{d}\mathbf{x}\mathrm{d}t$$

$$\le 2\left(F({\mathbf{u}}^0,{\mathbf{\sigma }}^0)+\frac{1+C_P}{2(1-\epsilon )}{\Vert \mathbf{f}\Vert }_{L^2(0,T;H^{-1}(\mathcal{D}))}^2\right).$$

$$⟨\mathbf{f},\mathbf{u}{⟩}_{H_0^1(\mathcal{D})}$$

$$\le \frac{1}{2{\nu }^2}{\Vert \mathbf{f}\Vert }_{H^{-1}(\mathcal{D})}^2+\frac{{\nu }^2}{2}(1+C_P){\Vert \nabla \mathbf{u}\Vert }_{L^2(\mathcal{D})}^2.$$

Full text

Keywords: FENE-P model, entropy, finite element method, convergence analysis, stress diffusion, existence of weak solutions

AMS Subject Classification: 35Q30, 65M12, 65M60, 76A10, 76M10, 82D60

1. Introduction

1.1. The FENE-P model

We consider the standard FENE-P model for a dilute polymeric fluid. The fluid, confined to an open bounded domain $\mathcal{D}\subset\mathbb{R}^{d}$ ($d=2$ or $3$) with a Lipschitz boundary $\partial\mathcal{D}$, is governed by the following non-dimensionalized system for a given $b\in\mathbb{R}_{>0}$:

(P) Find ${\boldsymbol{u}}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{u}}(t,\boldsymbol{x})\in\mathbb{R}^{d}$, $p:(t,\boldsymbol{x})\in\mathcal{D}_{T}:=(0,T)\times\mathcal{D}\mapsto p(t,\boldsymbol{x})\in\mathbb{R}$ and ${\boldsymbol{\sigma}}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{\sigma}}(t,\boldsymbol{x})\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}$ such that

[TABLE]

where

[TABLE]

Here $\mathbb{R}^{d\times d}_{\rm S}$ denotes the set of symmetric $\mathbb{R}^{d\times d}$ matrices, and $\mathbb{R}^{d\times d}_{{\rm S},>0}$ the set of symmetric positive definite $\mathbb{R}^{d\times d}$ matrices. In addition, $\boldsymbol{I}\in\mathbb{R}^{d\times d}_{{\rm S},>0}$ is the identity, and $\operatorname{tr}(\cdot)$ denotes trace. The unknowns in (P) are the velocity of the fluid, ${\boldsymbol{u}}$, the hydrostatic pressure, $p$, and the symmetric conformation tensor of the polymer molecules, ${\boldsymbol{\sigma}}$. The latter is linked to the symmetric polymeric extra-stress tensor ${\boldsymbol{\tau}}$ through the relation ${\boldsymbol{\tau}}=\frac{\varepsilon}{{\text{Wi}}}\,A({\boldsymbol{\sigma}})\,{\boldsymbol{\sigma}}$. In addition, $\boldsymbol{f}:(t,\boldsymbol{x})\in\mathcal{D}_{T}\mapsto\boldsymbol{f}(t,\boldsymbol{x})\in\mathbb{R}^{d}$ is the given density of body forces acting on the fluid; and the following given parameters are dimensionless: the Reynolds number ${\text{Re}}\in\mathbb{R}_{>0}$, the Weissenberg number ${\text{Wi}}\in\mathbb{R}_{>0}$, the elastic-to-viscous viscosity fraction $\varepsilon\in(0,1)$, and the FENE-P parameter $b>0$ (related to a maximal admissible extensibility of the polymer molecules within the fluid). For the sake of simplicity, we will limit ourselves to the no flow boundary conditions (1.1f). Finally, we denote $\boldsymbol{\nabla}{\boldsymbol{u}}(t,\boldsymbol{x})\in\mathbb{R}^{d\times d}$ the velocity gradient tensor field with $[\boldsymbol{\nabla}{\boldsymbol{u}}]_{ij}=\frac{\partial{\boldsymbol{u}}_{i}}{\partial\boldsymbol{x}_{j}}$, and $({\rm div}\,{\boldsymbol{\sigma}})(t,\boldsymbol{x})\in\mathbb{R}^{d}$ the vector field with $[{\rm div}\,{\boldsymbol{\sigma}}]_{i}=\sum_{j=1}^{d}\frac{\partial{\boldsymbol{\sigma}}_{ij}}{\partial\boldsymbol{x}_{j}}$.

For data $\boldsymbol{f}\equiv 0$, divergence free ${\boldsymbol{u}}^{0}\in[L^{2}(\mathcal{D})]^{d}$, and ${\boldsymbol{\sigma}}^{0}$, which is symmetric positive definite for a.e. $\boldsymbol{x}\in\mathcal{D}$, satisfying $\ln(1-\frac{\operatorname{tr}({\boldsymbol{\sigma}}^{0})}{b})\in L^{1}(\mathcal{D})$, then the existence of a global-in-time weak solution ${\boldsymbol{u}}\in L^{\infty}(0,T;[L^{2}(\mathcal{D})]^{d})\cap L^{2}(0,T,[H^{1}_{0}(\mathcal{D})]^{d})$, ${\boldsymbol{\sigma}}\in L^{\infty}(0,T;[L^{\infty}(\mathcal{D})]^{d\times d})$ and ${\boldsymbol{\tau}}\in L^{2}(0,T;[L^{2}(\mathcal{D})]^{d\times d})$ to (P), (1.1a–f), was proved in Masmoudi [20].

In this work, we consider finite element approximations of the FENE-P system (P) and the corresponding model with stress diffusion, (Pα), which is obtained by adding the dissipative term $\alpha\,\Delta{\boldsymbol{\sigma}}$ for a given $\alpha\in\mathbb{R}_{>0}$ to the right-hand side of (1.1c) with an additional no flux boundary condition for ${\boldsymbol{\sigma}}$ on $\partial\mathcal{D}$. This paper extends the results in Barrett and Boyaval [1], where finite element approximations of the corresponding Oldroyd-B models, where $A({\boldsymbol{\sigma}})=\boldsymbol{I}-{\boldsymbol{\sigma}}^{-1}$, were introduced and analysed. In fact, the convergence proof of the finite element approximation of the Oldroyd-B model with stress diffusion for $d=2$ in [1] provided the first existence proof of global-in-time weak solutions for this system. Note that $A({\boldsymbol{\sigma}})=\boldsymbol{I}-{\boldsymbol{\sigma}}^{-1}$ is the formal limit of (1.2) for infinite extensibility; that is, $b\rightarrow\infty$.

The model (Pα) has been considered computationally in Sureshkumar and Beris [24]. We recall also that El-Kareh and Leal [14] showed the existence of a weak solution to a modified stationary FENE-P system of equations, which included stress diffusion, but there an additional regularization was also present in their modified system and played an essential role in their proof. We stress that the dissipative term $\alpha\,\Delta{\boldsymbol{\sigma}}$ in (Pα) is not a regularization, but can be physically motivated through the centre-of-mass diffusion in the related microscopic-macroscopic polymer model, though with a positive $\alpha\ll 1$, see Barrett and Süli [3], [5], Schieber [23] and Degond and Liu [13].

Barrett and Süli have introduced, and proved the existence of global-in-time weak solutions for $d=2$ and $3$ to, microscopic-macroscopic dumbbell models of dilute polymers with center-of-mass diffusion in the corresponding Fokker–Planck equation for a finitely extensible nonlinear elastic (FENE) spring law or a Hookean-type spring law, see [4] and [6]. Recently, Barrett and Süli [8] have proved rigorously that the macroscopic Oldroyd-B model with stress diffusion is the exact closure of the microscopic-macroscopic Hookean dumbbell model with center-of-mass diffusion for $d=2$, when the existence of global-in-time weak solutions to both models can be proved. In addition, Barrett and Süli [7] have introduced and analysed a finite element approximation for the FENE microscopic-macroscopic dumbbell model with center-of-mass diffusion.

From a physical viewpoint, the FENE-P model is more realistic than the Oldroyd-B model because it accounts for the finite-extensibility of the polymer molecules in the fluid through the non-dimensional parameter $b>0$. From a mathematical viewpoint, compared to the Oldroyd-B model where the nonlinear terms are only the material derivative terms (like $({\boldsymbol{\nabla}{\boldsymbol{u}}}){\boldsymbol{\sigma}}$), the FENE-P model has an additional singular nonlinearity due to the factor $\left(1-\frac{\operatorname{tr}(\cdot)}{b}\right)^{-1}$ in the definition of $A(\cdot)$, which necessitates a careful mathematical treatment. Hence, this paper is not a trivial extension of [1]. In fact, the latter additional nonlinearity is exactly what makes the FENE-P model closer to the physics of polymers than the Oldroyd-B model, and thus also to many other macroscopic models based on different constitutive relations that have been developed by physicists for polymers. We note the FENE-P system is the approximate macroscopic closure of the FENE microscopic-macroscopic dumbbell model, whereas the Oldroyd-B system is the exact macroscopic closure of the Hookean microscopic-macroscopic dumbbell model. Hence, the microscopic-macroscopic dumbbell models corresponding to Oldroyd-B and FENE-P, only the spring laws differ; see e.g. Bird et al. [9] and Renardy [22] for a more complete review of the differences between the Oldroyd-B and the FENE-P models from the physical viewpoint, and for other macroscopic models with more nonlinear effects than the Oldroyd-B model, e.g. the Giesekus model and the Phan–Thien Tanner model.

Similarly to [1], our analysis in the present paper exploits the underlying free energy of the system, see Wapperom and Hulsen [26] and Hu and Lelièvre [19]. In particular, the finite element approximation of (Pα) has to be constructed extremely carefully to inherit this free energy structure, and requires the approximation of $\operatorname{tr}({\boldsymbol{\sigma}})$ as a new unknown. It is definitely not our goal to review all the macroscopic models used in rheology, although similar studies could probably be pursued for other macroscopic models endowed with a free energy. We will point out the main differences with Barrett and Boyaval [1], and we thus hope to sufficiently suggest how our technique could be adapted to any nonlinear model with a free energy. We believe that our approach contributes to a better understanding of the numerical stability of the models used in computational rheology, where numerical instabilities sometimes termed “High-Weissenberg Number Problems”, see HWNP in Owens and Phillips [21], still persist. Indeed, as exposed in Boyaval et al. [11], our point is that to make progress in this area one should identify sufficiently general rules for the derivation of good discretizations of macroscopic models such that they retain the dissipative structure of weak solutions to the system, at least in some benchmark flows.

The outline of this paper is as follows. First, we end this section by introducing our notation and some auxiliary results. In Section 2 we review the formal free energy bound for the FENE-P system (P). In Section 3 we introduce our regularization $G_{\delta}$ of of $G\equiv\ln$, which appears in the definition of the free energy of the FENE-P system (P). We then introduce a regularized problem (Pδ), and show a formal free energy bound for it. In Section 4, on assuming that $\mathcal{D}$ is a polytope for ease of exposition, we introduce our finite element approximation of (Pδ), namely (P${}^{\Delta t}_{\delta,h}$), based on approximating the pressure and the symmetric conformation tensor by piecewise constants; and the velocity field with continuous piecewise quadratics or a reduced version, where the tangential component on each simplicial edge ($d=2$) or face ($d=3$) is linear. Using the Brouwer fixed point theorem, we prove existence of a solution to (P${}^{\Delta t}_{\delta,h}$) and show that it satisfies a discrete regularized free energy bound for any choice of time step; see Theorem 4.1. We conclude by showing that, in the limit $\delta\rightarrow 0_{+}$, these solutions of (P${}^{\Delta t}_{\delta,h})$ converge to a solution of (P${}^{\Delta t}_{h}$) with the approximation of the conformation tensor being positive definite and having a trace strictly less than $b$. Moreover, this solution of (P${}^{\Delta t}_{h}$) satisfies a discrete free energy bound; see Theorem 4.2. Next, in Section 5 we introduce the FENE-P system with stress diffusion, (Pα), where the dissipative term $\alpha\,\Delta{\boldsymbol{\sigma}}$ has been added to the right-hand side of (1.1c). We then introduce the corresponding regularized version (Pα,δ), and show a formal free energy bound for it. In Section 6 we introduce our finite element approximation of (Pα,δ), namely (P${}^{\Delta t}_{\alpha,\delta,h}$), based on approximating the velocity field with continuous piecewise quadratics or the mini element, and the pressure, the symmetric conformation tensor and its trace by continuous piecewise linears. Here we assume that $\mathcal{D}$ is a convex polytope and that the finite element mesh consists of quasi-uniform non-obtuse simplices. Using the Brouwer fixed point theorem, we prove existence of a solution to (P${}^{\Delta t}_{\alpha,\delta,h}$) and show that it satisfies a discrete regularized free energy bound for any choice of time step; see Theorem 6.1. In Section 7 we prove, in the case $d=2$, (subsequence) convergence of the solutions of (P${}_{\alpha,\delta,h}^{\Delta t}$), as the regularization parameter, $\delta$, and the spatial, $h$, and temporal, $\Delta t$, discretization parameters tend to zero, to global-in-time weak solutions of (Pα); see Theorem 7.3. This existence result for (Pα) is new to the literature.

1.2. Notation and auxiliary results

The absolute value and the negative part of a real number $s\in\mathbb{R}$ are denoted by $\left|s\right|:=\max\{s,-s\}$ and $[s]_{-}=\min\{s,0\}$, respectively. We adopt the following notation for inner products

[TABLE]

where $\cdot^{T}$ and $\operatorname{tr}\left(\cdot\right)$ denote transposition and trace, respectively. The corresponding norms are

[TABLE]

We will use on several occasions that $\operatorname{tr}(\boldsymbol{\phi})=\operatorname{tr}(\boldsymbol{\phi}^{T})$ and $\operatorname{tr}(\boldsymbol{\phi}\boldsymbol{\psi})=\operatorname{tr}(\boldsymbol{\psi}\boldsymbol{\phi})$ for all $\boldsymbol{\phi},\boldsymbol{\psi}\in\mathbb{R}^{d\times d}$, and

[TABLE]

For any $\boldsymbol{\phi}\in\mathbb{R}^{d\times d}_{\rm S}$, there exists a decomposition

[TABLE]

where $\mathbf{O}\in\mathbb{R}^{d\times d}$ is an orthogonal matrix and $\mathbf{D}\in\mathbb{R}^{d\times d}$ a diagonal matrix. Hence, for any $g:\mathbb{R}\rightarrow\mathbb{R}$, one can define $g(\boldsymbol{\phi})\in\mathbb{R}^{d\times d}_{\rm S}$ as

[TABLE]

where $g({\bf D})\in\mathbb{R}^{d\times d}_{\rm S}$ is the diagonal matrix with entries $[g({\bf D})]_{ii}=g({\bf D}_{ii})$, $i=1,\ldots,d$. Although the diagonal decomposition (1.6) is not unique, (1.7) uniquely defines $g(\boldsymbol{\phi})$. We note for later purposes that

[TABLE]

One can show via diagonalization, see e.g. [1] for details, that for all concave function $g\in C^{1}(\mathbb{R})$, it holds

[TABLE]

where $g^{\prime}$ denotes the first derivative of $g$. If $g\in C^{1}(\mathbb{R})$ is convex, the inequalities in (1.9) are reversed. It follows from (1.9) and (1.3b) that for any $\boldsymbol{\phi}\in C^{1}([0,T];\mathbb{R}^{d\times d}_{\rm S})$ and any concave or convex $g\in C^{1}(\mathbb{R})$

[TABLE]

Of course, a similar result holds for spatial derivatives. Furthermore, the results (1.9) and (1.10) hold true when $C^{1}(\mathbb{R})$ and $\mathbb{R}^{d\times d}_{\rm S}$ are replaced by $C^{1}(\mathbb{R}_{>0})$ and $\mathbb{R}^{d\times d}_{{\rm S},>0}$ or $C^{1}(0,b)$ and $\mathbb{R}^{d\times d}_{{\rm S},>0,b}$. Finally, one can show that if $g\in C^{0,1}(\mathbb{R})$ with Lipschitz constant $g_{\rm Lip}$, then

[TABLE]

We adopt the standard notation for Sobolev spaces, e.g. $H^{1}(\mathcal{D}):=\{\eta:\mathcal{D}\mapsto\mathbb{R}\,:\,\int_{\mathcal{D}}[\,|\eta|^{2}+\|\nabla\eta\|^{2}\,]\,{\rm d}\boldsymbol{x}<\infty\}$ with $H^{1}_{0}(\mathcal{D})$ being the closure of $C^{\infty}_{0}(\mathcal{D})$ for the corresponding norm $\|\cdot\|_{H^{1}(\mathcal{D})}$. We denote the associated semi-norm as $|\cdot|_{H^{1}(\mathcal{D})}$. The topological dual of the Hilbert space $H^{1}_{0}(\mathcal{D})$, with pivot space $L^{2}(\mathcal{D})$, will be denoted by $H^{-1}(\mathcal{D})$. Such function spaces are naturally extended when the range $\mathbb{R}$ is replaced by $\mathbb{R}^{d}$, $\mathbb{R}^{d\times d}$ and $\mathbb{R}^{d\times d}_{\rm S}$; e.g. $H^{1}(\mathcal{D})$ becomes $[H^{1}(\mathcal{D})]^{d}$, $[H^{1}(\mathcal{D})]^{d\times d}$ and $[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}$ , respectively. For ease of notation, we write the corresponding norms and semi-norms as $\|\cdot\|_{H^{1}(\mathcal{D})}$ and $|\cdot|_{H^{1}(\mathcal{D})}$, respectively, as opposed to e.g. $\|\cdot\|_{[H^{1}(\mathcal{D})]^{d}}$ and $|\cdot|_{[H^{1}(\mathcal{D})]^{d}}$, respectively. We denote the duality pairing between $H^{-1}(\mathcal{D})$ and $H^{1}_{0}(\mathcal{D})$ as $\langle\cdot,\cdot\rangle_{H^{1}_{0}(\mathcal{D})}$, and we similarly write $\langle\cdot,\cdot\rangle_{H^{1}_{0}(\mathcal{D})}$ for the duality pairing between e.g. $[H^{-1}(\mathcal{D})]^{d}$ and $[H^{1}_{0}(\mathcal{D})]^{d}$. For notational convenience, we introduce also convex sets such as $[H^{1}(\mathcal{D})]^{d\times d}_{{\rm S},>0}:=\{\boldsymbol{\phi}\in[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}:\boldsymbol{\phi}\in\mathbb{R}^{d\times d}_{{\rm S},>0}\mbox{ a.e. in }\mathcal{D}\}$, and $[H^{1}(\mathcal{D})]^{d\times d}_{{\rm S},>0,b}:=\{\boldsymbol{\phi}\in[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}:\boldsymbol{\phi}\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}\mbox{ a.e. in }\mathcal{D}\}$.

In order to analyse (P), we adopt the notation

[TABLE]

Throughout the paper $C$ will denote a generic positive constant independent of the regularization parameter $\delta$ and the mesh parameters $h$ and $\Delta t$. Finally, we recall the Poincaré inequality

[TABLE]

where $C_{P}\in\mathbb{R}_{>0}$ depends only on $\mathcal{D}$.

2. Formal free energy bound for the problem (P)

In this section we recall from Hu and Lelièvre [19] the free energy structure of problem $(\bf P)$. Let $F({\boldsymbol{u}},{\boldsymbol{\sigma}})$ denote the free energy associated with a solution $({\boldsymbol{u}},p,{\boldsymbol{\sigma}})$ to problem $(\bf P)$, where we define

[TABLE]

with $\mathrm{S}^{\star}\subset\mathrm{S}_{>0,b}$ such that $F(\cdot,\cdot)$ is well-defined. Here the first term $\frac{{\text{Re}}}{2}\int_{\mathcal{D}}\|\boldsymbol{v}\|^{2}$ corresponds to the usual kinetic energy term, and the second term, which is nonnegative, is a relative entropy term. Moreover, on noting that $\ln$ is a concave function on $\mathbb{R}_{>0}$, we observe

[TABLE]

where the right-hand side is the free energy of the Oldroyd-B model under the same no flow boundary conditions, see e.g. [19] and [1]. Clearly, diagonalization yields that the relative entropy term of this Oldroyd-B model is nonnegative. Of course, the $\boldsymbol{I}$ term in the relative entropy for FENE-P and Oldroyd-B plays no real role, and just means that the minimum relative entropy for Oldroyd-B is zero and is obtained by $\boldsymbol{\phi}=\boldsymbol{I}$. Finally, we note that $\operatorname{tr}\left(\ln(\boldsymbol{\phi})\right)$ is rewritten as $\ln\left(\det(\boldsymbol{\phi})\right)$ in [19], which once again is easily deduced from diagonalization.

Proposition 2.1**.**

With $\boldsymbol{f}\in L^{2}(0,T;$ $[H^{-1}(\mathcal{D})]^{d})$ let $({\boldsymbol{u}},p,{\boldsymbol{\sigma}})$ be a sufficiently smooth solution to problem $(\bf P)$, (1.1a–f), such that ${\boldsymbol{\sigma}}(t,\cdot)\in\mathrm{S}^{\star}$ for $t\in(0,T)$. Then the free energy $F({\boldsymbol{u}},{\boldsymbol{\sigma}})$ satisfies for a.a. $t\in(0,T)$

[TABLE]

where the third term on the left-hand side is positive, via diagonalization, on recalling (1.2).

Proof.

Multiplying the Navier-Stokes equation (1.1a) with ${\boldsymbol{u}}$ and the stress equation (1.1c) with $\frac{\varepsilon}{2{\text{Wi}}}\,A({\boldsymbol{\sigma}})$, summing and integrating over $\mathcal{D}$ yields, after using integrations by parts, the boundary condition (1.1f) and the incompressibility property (1.1b) in the standard way, that

[TABLE]

It follows from the chain rule and (1.10) that

[TABLE]

On integrating (2.5) over $\mathcal{D}$, the $({\boldsymbol{u}}\cdot\boldsymbol{\nabla})$ term on the right-hand side vanishes as ${\boldsymbol{u}}(t,\cdot)\in\mathrm{V}$. On noting (1.2), (1.5a), (1.3b) and (1.1b), we obtain that

[TABLE]

Hence, on combining (2.4)–(2.6) and noting a trace property, we obtain the desired free energy equality (2.3). ∎

For later purposes, we note the following.

Remark 2.1**.**

The step in the above proof of testing (1.1c) with \frac{\varepsilon}{\mbox{\footnotesize2\rm{\text{Wi}}}}\,A({\boldsymbol{\sigma}}) is equivalent to testing (1.1c) with -\frac{\varepsilon}{\mbox{\footnotesize2\rm{\text{Wi}}}}{\boldsymbol{\sigma}}^{-1} and testing the corresponding trace equation

[TABLE]

with \frac{\varepsilon}{\mbox{\footnotesize2\rm Wi}}\left(1-\frac{\operatorname{tr}({\boldsymbol{\sigma}})}{b}\right)^{-1}, and adding.

Recall that in the limit $b\rightarrow\infty$ the FENE-P model formally converges to the Oldroyd-B model. It is thus interesting to note that when $b\rightarrow\infty$, the free energy equality (2.3) formally converges to the corresponding free energy equality for the Oldroyd-B model, on recalling (2.2). Finally, we note the following result.

Corollary 2.1**.**

Under the assumptions of Proposition 2.1 it follows that

[TABLE]

Proof.

One can bound the term $\langle\boldsymbol{f},{\boldsymbol{u}}\rangle_{H^{1}_{0}(\mathcal{D})}$ in (2.3), using the Cauchy-Schwarz and Young inequalities for $\nu\in\mathbb{R}_{>0}$, and the Poincaré inequality (1.13), by

[TABLE]

Combining (2.9) and (2.3) with $\nu^{2}=(1-\varepsilon)/(1+C_{P})$, and integrating in time yields the result (2.8). ∎

3. Formal free energy bound for a regularized problem (Pδ)

3.1. A regularization

Let $G:s\in\mathbb{R}_{>0}\mapsto\ln s\in\mathbb{R}$ denote the logarithm function, whose domain of definition can be straightforwardly extended to the set of symmetric positive definite matrices using (1.6) and (1.7). We define the following concave $C^{1,1}(\mathbb{R})$ regularization of $G$ based on a given parameter $\delta\in(0,1)$:

[TABLE]

We define also the following scalar functions

[TABLE]

Hence, we have that

[TABLE]

For later purposes, we note the following results concerning these functions.

Lemma 3.1**.**

For any $\boldsymbol{\phi},\boldsymbol{\psi}\in\mathbb{R}^{d\times d}_{\rm S}$, $\eta\in\mathbb{R}$ and for any $\delta\in(0,1)$ we have that

[TABLE]

In addition, if $\delta\in(0,\frac{1}{2}]$ we have that

[TABLE]

Proof.

All the results are proved in Lemma 2.1 in [1], except (3.4b) and this can be easily proved via diagonalization. ∎

We introduce the following regularization of $A$, (1.2), for any $\delta\in(0,\frac{1}{2}]$:

[TABLE]

In addition to Lemma 3.1 we will also make use of the following result, which is similar to (3.7).

Lemma 3.2**.**

For any ${s\in\mathbb{R}}$, $b\in\mathbb{R}_{>0}$ and $\delta\in(0,\frac{1}{2}]$, we have that

[TABLE]

Proof.

On recalling (3.1), we first note from the concavity of $G_{\delta}$ that

[TABLE]

From the scalar version of (3.7), we have that

[TABLE]

We note that

[TABLE]

Combining (3.10)–(3.14) yields the desired result (3.9a).

We now consider (3.9b). If $1-\frac{s}{b}\leq\delta$, i.e. $s\geq b(1-\delta)$, then

[TABLE]

If $1-\frac{s}{b}\geq\delta$, i.e. $s\leq b(1-\delta)$, then

[TABLE]

Combining (3.15) and (3.16) yields the desired result (3.9b). ∎

3.2. The regularized problem (Pδ)

Using the regularizations $G_{\delta},\,\beta_{\delta}$ and $A_{\delta}$ introduced above, we consider the following regularization of (P) for a given $\delta\in(0,\frac{1}{2}]$:

(Pδ) Find ${\boldsymbol{u}}_{\delta}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{u}}_{\delta}(t,\boldsymbol{x})\in\mathbb{R}^{d}$, $p_{\delta}:(t,\boldsymbol{x})\in\mathcal{D}_{T}\mapsto p_{\delta}(t,\boldsymbol{x})\in\mathbb{R}$ and ${\boldsymbol{\sigma}}_{\delta}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{\sigma}}_{\delta}(t,\boldsymbol{x})\in\mathbb{R}^{d\times d}_{\rm S}$ such that

[TABLE]

3.3. Formal free energy bound for (Pδ)

In this section, we extend the formal energy results (2.3) and (2.8) for (P) to problem (Pδ). We will assume throughout that

[TABLE]

where $\boldsymbol{n}_{\partial\mathcal{D}}$ is normal to $\partial\mathcal{D}$, $b^{\star},\,\sigma_{\rm min}^{0},\,\sigma_{\rm max}^{0}\in\mathbb{R}_{>0}$ with $b^{\star}<b$. Let $F_{\delta}({\boldsymbol{u}}_{\delta},{\boldsymbol{\sigma}}_{\delta},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta}))$ denote the free energy associated with a solution $({\boldsymbol{u}}_{\delta},p_{\delta},{\boldsymbol{\sigma}}_{\delta})$ to problem $\bf(P_{\delta})$, where we define

[TABLE]

Note that the second term in $F_{\delta}$ has been regularized in comparison with $F$ in (2.1). Similarly to (2.2), we have, on noting (3.10), the inequality

[TABLE]

where the right-hand side in (3.20) is the free energy of the corresponding regularized Oldroyd-B model, see [1] and note (3.4c). It also follows from (3.1) and (3.18) that

[TABLE]

Proposition 3.1**.**

Let $\delta\in(0,\frac{1}{2}]$ and $({\boldsymbol{u}}_{\delta},p_{\delta},{\boldsymbol{\sigma}}_{\delta})$ be a sufficiently smooth solution to problem $(\bf P_{\delta})$, (3.17a–f). Then the free energy $F_{\delta}({\boldsymbol{u}}_{\delta},{\boldsymbol{\sigma}}_{\delta},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta}))$ satisfies for a.a. $t\in(0,T)$

[TABLE]

where the third term on the left-hand side is nonnegative from (3.8) and (3.4b).

Proof.

Similarly to the proof of Proposition 2.1, we multiply the regularized Navier-Stokes equation (3.17a) by ${\boldsymbol{u}}_{\delta}$ and the regularized stress equation (3.17c) with $\frac{\varepsilon}{2{\text{Wi}}}A_{\delta}({\boldsymbol{\sigma}}_{\delta},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta}))$, sum and integrate over $\mathcal{D}$, use integrations by parts, the boundary condition (3.17f) and the incompressibility property (3.17b). This yields the desired result (3.22) on noting the following analogues of (2.5) and (2.6)

[TABLE]

and

[TABLE]

Here we have recalled (3.8) and (1.10) for (3.23), and (1.5a), (1.3b), (3.4a) and (3.17b) for (3.24). ∎

Similarly to Remark 2.1, we note the following.

Remark 3.1**.**

The step in the above proof of testing the regularized stress equation (3.17c) with \frac{\varepsilon}{\mbox{ \footnotesize2\rm{\text{Wi}}}}A_{\delta}({\boldsymbol{\sigma}}_{\delta},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta})) is equivalent to testing (3.17c) with -\frac{\varepsilon}{\mbox{\footnotesize2\rm{\text{Wi}}}}G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\delta}) and testing the corresponding regularized trace equation

[TABLE]

with \frac{\varepsilon}{\mbox{\footnotesize2\rm Wi}}G_{\delta}^{\prime}\left(1-\frac{\operatorname{tr}({\boldsymbol{\sigma}}_{\delta})}{b}\right), and adding.

Corollary 3.1**.**

Under the assumptions of Proposition 3.1 it follows that

[TABLE]

Proof.

The proof of (3.26) follows from (3.22) in the same way as (2.8) follows from (2.3), and in addition noting (3.21). ∎

4. Finite element approximation of (Pδ) and (P)

4.1. Finite element discretization

We now introduce a finite element discretization of the problem (Pδ), which satisfies a discrete analogue of (3.22).

The time interval $[0,T)$ is split into intervals $[t^{n-1},t^{n})$ with $\Delta t_{n}=t^{n}-t^{n-1}$, $n=1,\ldots,N_{T}$. We set $\Delta t:=\max_{n=1,\ldots,N_{T}}\Delta t_{n}$. We will assume throughout that the domain $\mathcal{D}$ is a polytope. We define a regular family of meshes $\{\mathcal{T}_{h}\}_{h>0}$ with discretization parameter $h>0$, which is built from partitionings of the domain $\mathcal{D}$ into regular open simplices so that

[TABLE]

Here $\rho_{k}$ is the diameter of the largest inscribed ball contained in the simplex $K_{k}$ and $h_{k}$ is the diameter of $K_{k}$, so that $h=\max_{k=1,\ldots,N_{K}}h_{k}$. For each element $K_{k}$, $k=1,\ldots,N_{K}$, of the mesh $\mathcal{T}_{h}$ let $\{P^{k}_{i}\}_{i=0}^{d}$ denotes its vertices, and $\{\boldsymbol{n}^{k}_{i}\}_{i=0}^{d}$ the outward unit normals of the edges $(d=2)$ or faces $(d=3)$ with $\boldsymbol{n}^{k}_{i}$ being that of the edge/face opposite vertex $P^{k}_{i}$, $i=0,\ldots,d$. In addition, let $\{\eta^{k}_{i}(\boldsymbol{x})\}_{i=0}^{d}$ denote the barycentric coordinates of $\boldsymbol{x}\in K_{k}$ with respect to the vertices $\{P^{k}_{i}\}_{i=0}^{d}$; that is, $\eta^{k}_{i}\in\mathbb{P}_{1}$ and $\eta^{k}_{i}(P^{k}_{j})=\delta_{ij}$, $i,\,j=0,\ldots,d$. Here $\mathbb{P}_{m}$ denote polynomials of maximal degree $m$ in $\boldsymbol{x}$, and $\delta_{ij}$ the Kronecker delta notation. Finally, we introduce $\partial\mathcal{T}_{h}:=\{E_{j}\}_{j=1}^{N_{E}}$ as the set of internal edges $E_{j}$ of triangles in the mesh $\mathcal{T}_{h}$ when $d=2$, or the set of internal faces $E_{j}$ of tetrahedra when $d=3$.

We approximate the problem (Pδ) by the problem (P${}_{\delta,h}^{\Delta t})$ based on the finite element spaces $\mathrm{W}_{h}^{0}\times\mathrm{Q}_{h}^{0}\times\mathrm{S}_{h}^{0}$. As is standard, we require the discrete velocity-pressure spaces $\mathrm{W}_{h}^{0}\times\mathrm{Q}_{h}^{0}\subset\mathrm{W}\times{\rm Q}$ satisfy the discrete Ladyshenskaya-Babuška-Brezzi (LBB) inf-sup condition

[TABLE]

see e.g. [16, p114]. In the following, we set

[TABLE]

where

[TABLE]

Here, for $k=1,\ldots,N_{K}$ and $i=0,\ldots,d$

[TABLE]

We introduce also

[TABLE]

which approximates $\mathrm{V}$. It is well-known that the choices (4.2a,b) satisfy (4.1), see e.g. [12, p221] for $\mathrm{W}_{h}^{0}=\mathrm{W}_{h}^{2}$ and $d=2$, and Chapter II, Sections 2.1 ($d=2$) and 2.3 ($d=3$) in [16] for $\mathrm{W}_{h}^{0}=\mathrm{W}_{h}^{2,-}$. Moreover, these particular choices of $\mathrm{S}_{h}^{0}$ and $\mathrm{Q}_{h}^{0}$ have the desirable property that

[TABLE]

which makes it a straightforward matter to mimic the free energy inequality (3.22) at a discrete level. Since $\mathrm{S}_{h}^{0}$ is discontinuous, we will use the discontinuous Galerkin method to approximate the advection term $({\boldsymbol{u}}_{\delta}\cdot\nabla){\boldsymbol{\sigma}}_{\delta}$ in the following. Then, for the boundary integrals, we will make use of the following definitions (see e.g. [15, p267]). Given $\boldsymbol{v}\in\mathrm{W}_{h}^{0}$, then for any $\boldsymbol{\phi}\in\mathrm{S}_{h}^{0}$ (or ${\rm Q}_{h}^{0}$) and for any point $\boldsymbol{x}$ that is in the interior of some $E_{j}\in\partial\mathcal{T}_{h}$, we define the downstream and upstream values of $\boldsymbol{\phi}$ at $\boldsymbol{x}$ by

[TABLE]

respectively. In addition, we denote by

[TABLE]

the jump and mean value, respectively, of $\boldsymbol{\phi}$ at the point $\boldsymbol{x}$ of boundary $E_{j}$. From (4.7), it is clear that the values of $\boldsymbol{\phi}^{+\boldsymbol{v}}|_{E_{j}}$ and $\boldsymbol{\phi}^{-\boldsymbol{v}}|_{E_{j}}$ can change along $E_{j}\in\partial\mathcal{T}_{h}$. Finally, it is easily deduced that

[TABLE]

where $\boldsymbol{n}\equiv\boldsymbol{n}(E_{j})$ is a unit normal to $E_{j}$, whose sign is of no importance, and $\boldsymbol{n}_{K_{k}}$ is the outward unit normal vector of boundary $\partial K_{k}$ of $K_{k}$. We note that similar ideas appear in upwind schemes; e.g. see Chapter IV, Section 5 in [16] for the Navier-Stokes equations.

4.2. A free energy preserving approximation

(P${}^{\Delta t}_{\delta,h}$) of (Pδ)

For any source term $\boldsymbol{f}\in L^{2}\left(0,T;[H^{-1}(\mathcal{D})]^{d}\right)$, we define the following piecewise constant function with respect to the time variable

[TABLE]

It is easily deduced that for $n=1,\ldots,N_{T}$

[TABLE]

Throughout this section we choose ${\boldsymbol{u}}_{h}^{0}\in\mathrm{V}_{h}^{0}$ to be the $L^{2}$ projection of ${\boldsymbol{u}}^{0}$ onto $\mathrm{V}_{h}^{0}$ and ${\boldsymbol{\sigma}}_{h}^{0}\in\mathrm{S}_{h}^{0}$ to be the $L^{2}$ projection of ${\boldsymbol{\sigma}}^{0}$ onto $\mathrm{S}_{h}^{0}$. Hence, we have that

[TABLE]

We are now ready to introduce our approximation (P${}_{\delta,h}^{\Delta t}$) of (Pδ) for $\delta\in(0,\frac{1}{2}]$:

(P${}_{\delta,h}^{\Delta t}$) Setting $({\boldsymbol{u}}_{\delta,h}^{0},{\boldsymbol{\sigma}}_{\delta,h}^{0})=({\boldsymbol{u}}_{h}^{0},{\boldsymbol{\sigma}}_{h}^{0})\in\mathrm{V}_{h}^{0}\times(\mathrm{S}_{h}^{0}\cap\mathrm{S}_{>0,b})$ as defined in (4.12a), then for $n=1,\ldots,N_{T}$ find $({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ such that for any test functions $(\boldsymbol{v},\boldsymbol{\phi})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$

[TABLE]

In deriving (P${}_{\delta,h}^{\Delta t}$), we have noted (1.5a) and that

[TABLE]

and we refer to [15, p267] and [11] for the consistency of our approximation of the stress advection term. We note that on replacing $A_{\delta}({\boldsymbol{\sigma}}_{\delta,h}^{n},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta,h}^{n}))$ with $\boldsymbol{I}-G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\delta,h}^{n})$ then (P${}_{\delta,h}^{\Delta t}$), (4.13a,b), collapses to the corresponding finite element approximation of Oldroyd-B studied in [1], see (3.12a,b) there.

Before proving existence of a solution to (P${}^{\Delta t}_{\delta,h}$), we first derive a discrete analogue of the energy bound (3.22) for (P${}^{\Delta t}_{\delta,h}$), which uses the elementary equality

[TABLE]

4.3. Energy bound for (P${}^{\Delta t}_{\delta,h}$)

Proposition 4.1**.**

For $n=1,\ldots,N_{T}$, a solution $\left({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n}\right)\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ to (P${}_{\delta,h}^{\Delta t}$), (4.13a,b), if it exists, satisfies

[TABLE]

Proof.

Similarly to (3.24), we have that

[TABLE]

where we have noted (3.8), (1.3b), (3.4a) and (4.5). Then, similarly to the proof of Proposition 3.1, we choose $\boldsymbol{v}={\boldsymbol{u}}_{\delta,h}^{n}\in\mathrm{V}_{h}^{0}$ in (4.13a) and $\boldsymbol{\phi}=\frac{\varepsilon}{2{\text{Wi}}}A_{\delta}({\boldsymbol{\sigma}}_{\delta,h}^{n},\operatorname{tr}({\boldsymbol{\sigma}}_{\delta,h}^{n}))\in\mathrm{S}_{h}^{0}$ in (4.13b) and obtain, on noting (4.15), (3.4a), (4.5) and (4.17), that

[TABLE]

It follows from (3.8), (3.4d) and the concavity of $G_{\delta}$ that

[TABLE]

Similarly to (4.19), we have, on recalling (3.8), (4.7) and (4.8), that

[TABLE]

Finally, we note from (4.9) and as ${\boldsymbol{u}}_{\delta,h}^{n-1}\in\mathrm{V}_{h}^{0}$ that for all $q\in{\rm Q}_{h}^{0}$

[TABLE]

Combining (4.18)–(4.21) yields the first desired inequality in (4.16). The second inequality in (4.16) follows immediately from (2.9) with $\nu^{2}=(1-\varepsilon)/(1+C_{P})$. ∎

4.4. Existence of a solution to (P${}^{\Delta t}_{\delta,h}$)

Proposition 4.2**.**

Let $\delta\in(0,\frac{1}{2}]$, then, given $({\boldsymbol{u}}_{\delta,h}^{n-1},{\boldsymbol{\sigma}}_{\delta,h}^{n-1})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ and for any time step $\Delta t_{n}>0$, there exists at least one solution $\left({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n}\right)\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ to (P${}_{\delta,h}^{\Delta t}$), (4.13a,b).

Proof.

We introduce the following inner product on the Hilbert space $\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$

[TABLE]

Given $({\boldsymbol{u}}_{\delta,h}^{n-1},{\boldsymbol{\sigma}}_{\delta,h}^{n-1})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$, let $\mathcal{F}:\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}\mapsto\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ be such that for any $(\boldsymbol{w},\boldsymbol{\psi})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$

[TABLE]

We note that a solution $({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})$ to (4.13a,b), if it exists, corresponds to a zero of $\mathcal{F}$; that is,

[TABLE]

In addition, it is easily deduced that the mapping $\mathcal{F}$ is continuous. For any $(\boldsymbol{w},\boldsymbol{\psi})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$, on choosing $(\boldsymbol{v},\boldsymbol{\phi})=\left(\boldsymbol{w},\frac{\varepsilon}{2{\text{Wi}}}A_{\delta}(\boldsymbol{\psi},\operatorname{tr}(\boldsymbol{\psi}))\right)$, we obtain analogously to (4.16) that

[TABLE]

Let us now assume that for any $\gamma\in\mathbb{R}_{>0}$, the continuous mapping $\mathcal{F}$ has no zero $({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})$ satisfying (4.24), which lies in the ball

[TABLE]

where

[TABLE]

Then for such $\gamma$, we can define the continuous mapping $\mathcal{G}_{\gamma}:\mathcal{B}_{\gamma}\mapsto\mathcal{B}_{\gamma}$ such that for all $(\boldsymbol{v},\boldsymbol{\phi})\in\mathcal{B}_{\gamma}$

[TABLE]

By the Brouwer fixed point theorem, $\mathcal{G}_{\gamma}$ has at least one fixed point $(\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma})$ in $\mathcal{B}_{\gamma}$. Hence it satisfies

[TABLE]

On noting (4.2c) and (4.29), we have that

[TABLE]

where $\mu_{h}:=[1/(\min_{k\in N_{K}}|K_{k}|)]^{\frac{1}{2}}$. Then (3.19), (3.20), (3.7), (4.30) and (4.29) yield that

[TABLE]

Hence for all $\gamma$ sufficiently large, it follows from (4.25), (4.31) and (3.4b) that

[TABLE]

On the other hand as $(\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma})$ is a fixed point of ${\mathcal{G}}_{\gamma}$, we have that

[TABLE]

It follows from (3.8), (3.9b), (3.7) and similarly to (4.31), on noting (4.30) and (4.29), that

[TABLE]

Therefore on combining (4.33) and (4.34), we have for all $\gamma$ sufficiently large that

[TABLE]

which obviously contradicts (4.32). Hence the mapping $\mathcal{F}$ has a zero in $\mathcal{B}_{\gamma}$ for $\gamma$ sufficiently large, and so there exists a solution $({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})$ to (4.13a,b). ∎

Theorem 4.1**.**

For any $\delta\in(0,\frac{1}{2}]$, $N_{T}\geq 1$ and any partitioning of $[0,T]$ into $N_{T}$ time steps, there exists a solution $\{({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}]^{N_{T}}$ to (P${}^{\Delta t}_{\delta,h}$), (4.13a,b). In addition, it follows for $n=1,\ldots,N_{T}$ that

[TABLE]

which yields that

[TABLE]

Moreover, for some $C(h,\Delta t)\in{\mathbb{R}}_{>0}$, but independent of $\delta$, it follows that for $k=1,\ldots,N_{K}$ and $n=1,\ldots,N_{T}$

[TABLE]

Proof.

Existence of a solution to (P${}^{\Delta t}_{\delta,h}$) and the first inequality in (4.36) follow immediately from Propositions 4.2 and 4.1, respectively. Similarly to (3.21), the second inequality in (4.36) is a direct consequence of (3.19), (3.1), (4.12b,c) and (4.11a). Finally, the third inequality in (4.36) follows from (4.12a–c) and (3.18).

It follows from (4.36) and (3.20) that

[TABLE]

The first three bounds in (4.37) then follow immediately from (4.39) and (3.7). Next we note that (3.4c), (3.19) and (4.36) yield that

[TABLE]

The last bound in (4.37) is then simply obtained by using a scalar version of (3.7).

Next, we deduce from (4.36), (3.8) and (3.4b) that for $n=1,\ldots,N_{T}$

[TABLE]

For any $\delta>0$, $\beta_{\delta}({\boldsymbol{\sigma}}_{\delta,h}^{n})\in\mathbb{R}^{d\times d}_{{\rm S},>0}$ and so it follows from (1.3b), (1.4b), (3.8), (3.4a), (1.5b), (4.41), (1.8), (3.3) and (4.37) that

[TABLE]

The desired result (4.38a) follows immediately from (4.42). Similarly to (4.42), we have from (4.41), (1.8), (3.8) and (3.4a) that

[TABLE]

The desired result (4.38b) follows immediately from (4.43) and (4.38a). ∎

4.5. Convergence of (P${}^{\Delta t}_{\delta,h}$) to (P${}^{\Delta t}_{h}$)

We now consider the corresponding direct finite element approximation of (P), i.e. (P${}^{\Delta t}_{h}$) without the regularization $\delta$:

(P${}^{\Delta t}_{h}$) Given initial conditions $({\boldsymbol{u}}_{h}^{0},{\boldsymbol{\sigma}}_{h}^{0})\in\mathrm{V}_{h}^{0}\times(\mathrm{S}_{h}^{0}\cap\mathrm{S}_{>0,b})$ as defined in (4.12a), then for $n=1,\ldots,N_{T}$ find $({\boldsymbol{u}}_{h}^{n},{\boldsymbol{\sigma}}_{h}^{n})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$ such that for any test functions $(\boldsymbol{v},\boldsymbol{\phi})\in\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}$

[TABLE]

We note that (4.44a,b) and $F({\boldsymbol{u}}_{h}^{n},{\boldsymbol{\sigma}}_{h}^{n})$ are only well-defined if ${\boldsymbol{\sigma}}_{h}^{n}\in\mathrm{S}_{h}^{0}\cap\mathrm{S}_{>0,b}$. We also note that on replacing $A({\boldsymbol{\sigma}}_{h}^{n})$ with $\boldsymbol{I}-({\boldsymbol{\sigma}}_{h}^{n})^{-1}$ then (P${}_{h}^{\Delta t}$), (4.44a,b), collapses to the corresponding finite element approximation of Oldroyd-B studied in [1], see (3.35a,b) there.

Theorem 4.2**.**

For all regular partitionings $\mathcal{T}_{h}$ of $\mathcal{D}$ into simplices $\{K_{k}\}_{k=1}^{N_{K}}$ and all partitionings $\{\Delta t_{n}\}_{n=1}^{N_{T}}$ of $[0,T]$, there exists a subsequence $\{\{({\boldsymbol{u}}_{\delta,h}^{n},{\boldsymbol{\sigma}}_{\delta,h}^{n})\}_{n=1}^{N_{T}}\}_{\delta>0}$, where $\{({\boldsymbol{u}}_{\delta,h}^{n},$ ${\boldsymbol{\sigma}}_{\delta,h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}]^{N_{T}}$ solves (P${}^{\Delta t}_{\delta,h}$), (4.13a,b), and $\{({\boldsymbol{u}}_{h}^{n},{\boldsymbol{\sigma}}_{h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}]^{N_{T}}$ such that for the subsequence

[TABLE]

In addition, for all $t\in[0,T]$ $n=1,\ldots,N_{T}$, ${\boldsymbol{\sigma}}_{h}^{n}\mid_{K_{k}}\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}$, $k=1,\ldots,N_{K},$. Moreover, $\{({\boldsymbol{u}}_{h}^{n},{\boldsymbol{\sigma}}_{h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}]^{N_{T}}$ solves (P${}_{h}^{\Delta t}$), (4.44a,b), and for $n=1,\ldots,N_{T}$

[TABLE]

Proof.

For any integer $n\in[1,N_{T}]$, the desired subsequence convergence results (4.45) follow immediately from (4.37), as $({\boldsymbol{u}}_{\delta,h}^{n},$ ${\boldsymbol{\sigma}}_{\delta,h}^{n})$ are finite dimensional for fixed $\mathrm{V}_{h}^{0}\times S_{h}^{0}$. It also follows from (4.37), (4.45) and (1.11) that $[{\boldsymbol{\sigma}}_{h}^{n}]_{-}$ and $[b-\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{n})]_{-}$ vanish on $\mathcal{D}$, so that ${\boldsymbol{\sigma}}_{h}^{n}$ must be non-negative definite and $\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{n})\leq b$ a.e. on $\mathcal{D}$. Moreover, on noting this, (4.45), (3.3) and (1.11), we have the following subsequence convergence results

[TABLE]

If $\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{n})|_{K_{k}}=b$ on some simplex $K_{k}$, then for the subsequence of (4.45) we have that

[TABLE]

In addition, it follows from (1.8), (3.3) and (4.48) for $\delta$ sufficiently small that

[TABLE]

Hence, (4.38a) and (4.49) yield for the subsequence of (4.45) for all $\delta$ sufficiently small that

[TABLE]

but this contradicts (4.48) on recalling (3.2) and (3.3). Therefore, $\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{n})|_{K_{k}}<b$ on all simplices $K_{k}$, and so it follows from (4.45), (3.2) and (3.3) that

[TABLE]

as $\delta\rightarrow 0_{+}$.

The results (4.38b) and (4.51) yield for $\delta$ sufficiently small that

[TABLE]

Furthermore, it follows from (4.52), (4.47), (3.2), (3.3) and as $[\beta_{\delta}({\boldsymbol{\sigma}}_{\delta,h}^{n})]^{-1}\beta_{\delta}({\boldsymbol{\sigma}}_{\delta,h}^{n})=\boldsymbol{I}$ that the following subsequence result

[TABLE]

holds, and so ${\boldsymbol{\sigma}}_{h}^{n}\mid_{K_{k}}\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}$, $k=1,\ldots,N_{K}$. Therefore, we have from (4.45) and (3.1) that

[TABLE]

Since ${\boldsymbol{u}}_{\delta,h}^{n-1},\,{\boldsymbol{u}}_{h}^{n-1}\in C(\overline{\mathcal{D}})$, it follows from (4.45), (4.7) and (4.8) that for $j=1,\ldots,N_{E}$ and all $\boldsymbol{\phi}\in\mathrm{S}_{h}^{0}$

[TABLE]

Hence using (4.45), (4.47), (4.51) and (4.55), we can pass to the limit $\delta\rightarrow 0_{+}$ for the subsequence in (P${}_{\delta,h}^{\Delta t}$), (4.13a,b), to show that $\{({\boldsymbol{u}}_{h}^{n},{\boldsymbol{\sigma}}_{h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{0}\times\mathrm{S}_{h}^{0}]^{N_{T}}$ solves (P${}_{h}^{\Delta t}$), (4.44a,b). Similarly, using (4.45), (4.47), (4.51), (4.53) and (4.54), and noting (3.19) and (2.1), we can pass to the limit $\delta\rightarrow 0_{+}$ in (4.16) to obtain (4.46). ∎

5. Fene-p model with stress diffusion

5.1. Model (Pα), (P) with stress diffusion

In this section we consider the following modified version of (P), (1.1a–f), with a stress diffusion term for a given constant $\alpha\in\mathbb{R}_{>0}$:

(Pα) Find ${\boldsymbol{u}}_{\alpha}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{u}}_{\alpha}(t,\boldsymbol{x})\in\mathbb{R}^{d}$, $p_{\alpha}:(t,\boldsymbol{x})\in\mathcal{D}_{T}\mapsto p_{\alpha}(t,\boldsymbol{x})\in\mathbb{R}$ and ${\boldsymbol{\sigma}}_{\alpha}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{\sigma}}_{\alpha}(t,\boldsymbol{x})\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}$ such that

[TABLE]

Hence problem (Pα) is the same as (P), but with the added diffusion term $\alpha\Delta{\boldsymbol{\sigma}}_{\alpha}$ for the stress equation (5.1c), and the associated Neumann boundary condition (5.1g). Similarly to (Pδ), (3.17a–f), we introduce a regularization of (Pα,δ) of (Pα) mimicking the free energy structure of (Pα). Moreover, we need to be able to construct a finite element approximation of (Pα,δ) that satisfies a discrete analogue of this free energy structure. Apart from the obvious addition of the stress diffusion term, there are three other distinct differences. First, one has to deal with the advective term in the stress equation differently, see (5.11) below, to the approach used in (3.23) for the stability bound, as the approach there cannot be mimicked at a discrete level using continuous piecewise linear functions to approximate ${\boldsymbol{\sigma}}_{\alpha,\delta}$, the regularization of ${\boldsymbol{\sigma}}_{\alpha}$. Note that one cannot use $\mathrm{S}_{h}^{0}$ with the desirable property (4.6) to approximate ${\boldsymbol{\sigma}}_{\alpha,\delta}$ due to the additional stress diffusion term. Second, as a consequence of this stress advective term, one has to introduce another regularization of $\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$, $\varrho_{\alpha,\delta}$, as well as the obvious candidate $\operatorname{tr}(\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta}))$, and solve for this directly, see Remark 5.1 below. Third, it is desirable for the convergence analysis, as $\delta\rightarrow 0$, to have a uniform $L^{2}(\mathcal{D}_{T})$ bound on the extra stress term $A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})$ in the Navier–Stokes equation (5.5a) below. To achieve this we replace $A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})$ there by $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})$, where we define, for $\delta\in(0,\min\{\frac{1}{2},b\}]$,

[TABLE]

with

[TABLE]

It follows from (5.2), (5.3), (1.3b) and (1.4b) that

[TABLE]

Hence a uniform $L^{2}(\mathcal{D}_{T})$ bound on $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\,\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})$ follows from a uniform $L^{1}(\mathcal{D}_{T})$ on $\operatorname{tr}\left(\left(A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\right)^{2}\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})\right)$, which will follow from the free energy bound. Although $\varrho_{\alpha,\delta}\neq\operatorname{tr}(\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta}))$, we will show, in the limit $\delta\rightarrow 0$, that $\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})\rightarrow{\boldsymbol{\sigma}}_{\alpha}$ and $\varrho_{\alpha,\delta}\rightarrow\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$ with ${\boldsymbol{\sigma}}_{\alpha}(\cdot,\cdot)\in\mathbb{R}^{d\times d}_{{\rm S},>0,b}$, and hence implying that $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})\rightarrow 1$. In order to maintain the free energy bound we need to include $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})$ on the right-hand sides of (5.5c,d) below.

Therefore, we consider the following regularization of (Pα) for a given $\delta\in(0,\min\{\frac{1}{2},b\}]$:

(Pα,δ) Find ${\boldsymbol{u}}_{\alpha,\delta}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{u}}_{\alpha,\delta}(t,\boldsymbol{x})\in\mathbb{R}^{d}$, $p_{\alpha,\delta}:(t,\boldsymbol{x})\in\mathcal{D}_{T}\mapsto p_{\alpha,\delta}(t,\boldsymbol{x})\in\mathbb{R}$, ${\boldsymbol{\sigma}}_{\alpha,\delta}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto{\boldsymbol{\sigma}}_{\alpha,\delta}(t,\boldsymbol{x})\in\mathbb{R}^{d\times d}_{\rm S}$ and $\varrho_{\alpha,\delta}:(t,\boldsymbol{x})\in[0,T)\times\mathcal{D}\mapsto\varrho_{\alpha,\delta}(t,\boldsymbol{x})\in\mathbb{R}$ such that

[TABLE]

We note from (3.3) that $-b\,\boldsymbol{\nabla}\beta_{\delta}\left(1-\frac{\varrho_{\alpha,\delta}}{b}\right)=\boldsymbol{\nabla}\varrho_{\alpha,\delta}$ if $\left(1-\frac{\varrho_{\alpha,\delta}}{b}\right)\geq\delta$, i.e. $\varrho_{\alpha,\delta}\leq b(1-\delta)$. We remark again, due the required regularization of the advective terms in (5.5c,d), that $\varrho_{\alpha,\delta}\neq\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta})$. However, we note that on taking the trace of (5.5c), subtracting (5.5d) and integrating over $\mathcal{D}$, yields, on noting (5.5b,f–h) and (1.3b), that for all $t\in[0,T]$

[TABLE]

5.2. Formal free energy bound for (Pα,δ)

First, similarly to (3.1), we introduce, for $\delta\in(0,1)$, the concave $C^{1,1}(\mathbb{R}_{>0})$ function

[TABLE]

We have the following analogue of Proposition 3.1.

Proposition 5.1**.**

Let $\alpha\in{\mathbb{R}}_{>0}$, $\delta\in(0,\min\{\frac{1}{2},b\}]$ and $({\boldsymbol{u}}_{\alpha,\delta},p_{\alpha,\delta},{\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})$ be a sufficiently smooth solution to problem $(\bf P_{\alpha,\delta})$, (5.5a–h). Then the free energy $F_{\delta}({\boldsymbol{u}}_{\alpha,\delta},{\boldsymbol{\sigma}}_{\alpha,\delta},\varrho_{\alpha,\delta})$ satisfies for a.a. $t\in(0,T)$

[TABLE]

Proof.

Similarly to the proof of Proposition 3.22, on noting Remark 3.1, we multiply the Navier-Stokes equation (5.5a) with ${\boldsymbol{u}}_{\alpha,\delta}$ and the stress equation (5.5c) with $-\frac{\varepsilon}{2{\text{Wi}}}G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta})$, the trace equation (5.5d) with $\frac{\varepsilon}{2{\text{Wi}}}G_{\delta}^{\prime}\left(1-\frac{\varrho_{\alpha,\delta}}{b}\right)$, sum and integrate over $\mathcal{D}$. This yields, after performing integration by parts and noting (5.5b,g,h), (1.5a), (3.4a) and (3.8), that

[TABLE]

Using (1.10), we have that

[TABLE]

As $\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta})\equiv H_{\delta}^{\prime}(G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta}))$, on recalling (5.7), we have, on noting ${\boldsymbol{u}}_{\alpha,\delta}\in\mathrm{V}$ and the spatial version of (1.10), that

[TABLE]

Similarly to (5.11), we have that

[TABLE]

Similarly to (3.4e), we have that

[TABLE]

Combining (5.9)–(5.13a,b) yields the desired result (5.8). ∎

Remark 5.1**.**

We note if we multiply the advection term in (5.5c) by $-G_{\delta}^{\prime}(1-\frac{\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta})}{b})\,\boldsymbol{I}$ and integrate over $\mathcal{D}$, then we obtain, on noting ${\boldsymbol{u}}_{\alpha,\delta}\in\mathrm{V}$ and (5.12), that

[TABLE]

Hence, the need for the new variable $\varrho_{\alpha,\delta}$ in order to mimic the free energy structure of (Pα)**.

The following Corollary follows from (5.8) on noting the proof of Corollary 3.1.

Corollary 5.1**.**

Under the assumptions of Proposition 5.1 it follows that

[TABLE]

6. Finite element approximation of (Pα,δ)

6.1. Finite element discretization

We now introduce a conforming finite element discretization of (Pα,δ), (5.5a–h), which satisfies a discrete analogue of (5.8). As noted in Section 5, we cannot use $\mathrm{S}_{h}^{0}$ with the desirable property (4.6) to approximate ${\boldsymbol{\sigma}}_{\alpha,\delta}$, as we now have the added diffusion term. In the following, we choose

[TABLE]

where $\mathrm{W}_{h}^{2}$ is defined as in (4.3a) and, on recalling the barycentric coordinate notation used in (4.4),

[TABLE]

The velocity-pressure choice, $\mathrm{W}_{h}^{2}\times{\rm Q}_{h}^{1}$, is the lowest order Taylor-Hood element. It satisfies (4.1) with $\mathrm{W}_{h}^{0}$ and ${\rm Q}_{h}^{0}$ replaced by $\mathrm{W}_{h}^{2}$ and ${\rm Q}_{h}^{1}$, respectively, provided, in addition to $\{{\mathcal{T}}_{h}\}_{h>0}$ being a regular family of meshes, that each simplex has at least one vertex in $\mathcal{D}$, see p177 in Girault and Raviart [16] in the case $d=2$ and Boffi [10] in the case $d=3$. Of course, this is a very mild restriction on $\{{\mathcal{T}}_{h}\}_{h>0}$. The velocity-pressure choice, $\mathrm{W}_{h}^{1,+}\times{\rm Q}_{h}^{1}$, is called the mini-element. It satisfies (4.1) with $\mathrm{W}_{h}^{0}$ and ${\rm Q}_{h}^{0}$ replaced by $\mathrm{W}_{h}^{1,+}$ and ${\rm Q}_{h}^{1}$, respectively; see Chapter II, Section 4.1 in Girault and Raviart [16] in the case $d=2$ and Section 4.2.4 in Ern and Guermond [15] in the case $d=3$. Hence for both choices of $\mathrm{W}_{h}^{1}$, it follows that for all $\boldsymbol{v}\in\mathrm{V}$ there exists a sequence $\{\boldsymbol{v}_{h}\}_{h>0}$, with $\boldsymbol{v}_{h}\in\mathrm{V}_{h}^{1}$, such that

[TABLE]

We recall the well-known local inverse inequality for ${\rm Q}^{1}_{h}$

[TABLE]

We recall a similar well-known local inverse inequality for $\mathrm{S}_{h}^{1}$ and $\mathrm{V}_{h}^{1}$

[TABLE]

We introduce the interpolation operator $\pi_{h}:C(\overline{\mathcal{D}})\rightarrow{\rm Q}^{1}_{h}$, and extended naturally to $\pi_{h}:[C(\overline{\mathcal{D}})]^{d\times d}_{\rm S}\rightarrow{\rm S}^{1}_{h}$, such that for all $\eta\in C(\overline{\mathcal{D}})$ and $\boldsymbol{\phi}\in[C(\overline{\mathcal{D}})]^{d\times d}_{\rm S}$

[TABLE]

where $\{P_{p}\}_{p=1}^{N_{P}}$ are the vertices of $\mathcal{T}_{h}$. As $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$ and $q\in\mathrm{Q}_{h}^{1}$ do not imply that $G_{\delta}^{\prime}(\boldsymbol{\phi})\in\mathrm{S}_{h}^{1}$ and $G_{\delta}^{\prime}\left(1-\frac{q}{b}\right)\in\mathrm{Q}_{h}^{1}$, we have to test the finite element approximation of (5.5c,d) with $-\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}[G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})]\in\mathrm{S}_{h}^{1}$ and $\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}\left[G_{\delta}^{\prime}\left(1-\frac{\varrho_{\alpha,\delta,h}^{n}}{b}\right)\right]\in\mathrm{Q}_{h}^{1}$, respectively, where ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n}\in\mathrm{S}_{h}^{1}$ and $\varrho_{\alpha,\delta,h}^{n}\in\mathrm{Q}_{h}^{1}$ are our finite element approximations to ${\boldsymbol{\sigma}}_{\alpha,\delta}$ and $\varrho_{\alpha,\delta}$ at time level $t_{n}$. Therefore the finite element approximation of (5.5c,d) have to be constructed to mimic the results (5.9)–(5.13a,b), when tested with $-\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}[G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})]\in\mathrm{S}_{h}^{1}$ and $\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}\left[G_{\delta}^{\prime}\left(1-\frac{\varrho_{\alpha,\delta,h}^{n}}{b}\right)\right]\in\mathrm{Q}_{h}^{1}$, respectively.

In order to mimic (5.10) and the (Pα,δ) analogue of (3.24), we need to use numerical integration (vertex sampling). We note the following results. As the basis functions associated with ${\rm Q}^{1}_{h}$ and $\mathrm{S}_{h}^{1}$ are nonnegative and sum to unity everywhere, we have, on noting (1.5b), for $k=1,\ldots,N_{K}$ that

[TABLE]

where $r_{1},\,r_{2}\in(1,\infty)$ satisfy $r_{1}^{-1}+r_{2}^{-1}=1$. In addition, we have for $k=1,\ldots,N_{K}$ that

[TABLE]

The first inequality in (6.8) follows immediately from (6.7b), and the second from applying (6.4) and a Cauchy–Schwarz inequality. Of course, scalar versions of (6.7a,b) and (6.8) hold with $[C(\overline{K_{k}})]^{d\times d}_{\rm S}$ and $\mathrm{S}_{h}^{1}$ replaced by $C(\overline{K_{k}})$ and $\mathrm{Q}_{h}^{1}$, respectively.

Furthermore, for later use, we recall the following well-known results concerning the interpolant $\pi_{h}$ for $k=1,\ldots,N_{K}$:

[TABLE]

In order to mimic (5.11) and (5.12), we have to carefully construct our finite element approximation of the advective terms in (5.5c,d). Our construction is a non-trivial extension of an approach that has been used in the finite element approximation of fourth-order degenerate nonlinear parabolic equations, such as the thin film equation; see e.g. Grün and Rumpf [17] and Barrett and Nürnberg [2]. Let $\{\boldsymbol{e}_{i}\}_{i=1}^{d}$ be the orthonormal vectors in $\mathbb{R}^{d}$, such that the $j^{th}$ component of $\boldsymbol{e}_{i}$ is $\delta_{ij}$, $i,j=1,\ldots,d$. Let $\widehat{K}$ be the standard open reference simplex in $\mathbb{R}^{d}$ with vertices $\{\widehat{P}_{i}\}_{i=0}^{d}$, where $\widehat{P}_{0}$ is the origin and $\widehat{P}_{i}=\boldsymbol{e}_{i}$, $i=1,\ldots,d$. Given a simplex $K_{k}\in\mathcal{T}_{h}$ with vertices $\{P^{k}_{i}\}_{i=0}^{d}$, then there exists a non-singular matrix $B_{k}$ such that the linear mapping

[TABLE]

maps vertex $\widehat{P}_{i}$ to vertex $P^{k}_{i}$, $i=0,\ldots,d$. Hence $\mathcal{B}_{k}$ maps $\widehat{K}$ to $K_{k}$. For all $\eta\in{\rm Q}^{1}_{h}$ and $K_{k}\in\mathcal{T}_{h}$, we define

[TABLE]

where for all $\widehat{\boldsymbol{x}}\in\widehat{K}$

[TABLE]

Such notation is easily extended to $\boldsymbol{\phi}\in\mathrm{S}^{1}_{h}$.

Given $q\in\mathrm{Q}_{h}^{1}$ and $K_{k}\in\mathcal{T}_{h}$, then first, for $j=1,\ldots,d$, we find the unique $\widehat{\Lambda}_{\delta,j}^{k}(\widehat{q})\in\mathbb{R}$, which is continuous on $q$, such that

[TABLE]

where $(\widehat{\pi}_{h}\widehat{\eta})(\widehat{\boldsymbol{x}})\equiv(\pi_{h}\eta)(\mathcal{B}_{k}\widehat{\boldsymbol{x}})$ for all $\widehat{\boldsymbol{x}}\in\widehat{K}$ and $\eta\in C(\overline{K_{k}})$. We set

[TABLE]

where we have noted that $G_{\delta}^{\prime}(\beta_{\delta}(s))=G_{\delta}^{\prime}(s)$ for all $s\in\mathbb{R}$ and $G_{\delta}^{\prime}(\cdot)$ is strictly decreasing on $[\delta,\infty)$. Clearly, $\widehat{\Lambda}_{\delta,j}^{k}(\widehat{q})\in\mathbb{R}$, $j=1,\ldots,d$, satisfies (6.13) and depends continuously on $q\mid_{K_{k}}$.

Next, we extend the construction (6.16) for a given $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$ and $K_{k}\in\mathcal{T}_{h}$, to find for $j=1,\ldots,d$ the unique $\widehat{\Lambda}_{\delta,j}^{k}(\widehat{\boldsymbol{\phi}})\in\mathbb{R}^{d\times d}_{\rm S}$, which is continuous on $\boldsymbol{\phi}$, such that

[TABLE]

To construct $\widehat{\Lambda}_{\delta,j}^{k}(\widehat{\boldsymbol{\phi}})$ satisfying (6.17), we note the following. We have from (5.7) and (1.9) that

[TABLE]

Since $G_{\delta}^{\prime}(\beta_{\delta}(s))=G_{\delta}^{\prime}(s)$ for all $s\in\mathbb{R}$, it follows from (3.4e) that

[TABLE]

Therefore the left-hand side of (6.19) is zero if and only if $G_{\delta}^{\prime}(\beta_{\delta}(\boldsymbol{\phi}(P_{j}^{k})))=G_{\delta}^{\prime}(\beta_{\delta}(\boldsymbol{\phi}(P_{0}^{k})))$, which is equivalent to $\beta_{\delta}(\boldsymbol{\phi}(P_{j}^{k}))=\beta_{\delta}(\boldsymbol{\phi}(P_{0}^{k}))$ as $G_{\delta}^{\prime}(\cdot)$ is invertible on $[\delta,\infty)$, the range of $\beta_{\delta}(\cdot)$. Hence, on noting (6.12), (6.18), (6.19) and (1.3b), we have that

[TABLE]

satisfies (6.17) for $j=1,\ldots,d$; where $\lambda_{\delta,j}^{k}\in[0,1]$ is defined as

[TABLE]

Furthermore, $\widehat{\Lambda}_{\delta,j}^{k}(\widehat{\boldsymbol{\phi}})\in\mathbb{R}^{d\times d}_{\rm S}$, $j=1,\ldots,d$, depends continuously on $\boldsymbol{\phi}\mid_{K_{k}}$.

Therefore given $q\in\mathrm{Q}_{h}^{1}$ and $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$, we introduce, for $m,p=1,\ldots,d$,

[TABLE]

It follows from (6.21a,b), (6.13), (6.17) and (6.11) that

[TABLE]

and for $m=1,\ldots,d$

[TABLE]

For a more precise version of (6.22), see Lemma 6.4 below. Of course, for (6.21a) and (6.23a) we can adopt the more compact notation on $K_{k}$, $k=1,\ldots,N_{K}$,

[TABLE]

where $\Lambda_{\delta}^{k}(\widehat{q})\in\mathbb{R}^{d\times d}_{\rm S}$ is diagonal with $[\Lambda_{\delta}^{k}(\widehat{q})]_{jj}=\Lambda_{\delta,j}^{k}(\widehat{q})$, $j=1,\ldots,d$, so that $\Lambda_{\delta}(q)\in\mathbb{R}^{d\times d}_{\rm S}$ with $[\Lambda_{\delta}(q)]_{mp}=\Lambda_{\delta,m,p}(q)$, $m,p=1,\ldots,d$.

Finally, as the partitioning ${\mathcal{T}}_{h}$ consists of regular simplices, we have that

[TABLE]

Hence, it follows from (6.21a,b), (6.25) and (6.20a,b) that for $k=1,\ldots,N_{K}$

[TABLE]

In order to mimic (5.13a,b), we shall assume from now on that the family of meshes, $\{\mathcal{T}_{h}\}_{h>0}$, for the polytope $\mathcal{D}$ consists of non-obtuse simplices only, i.e. all dihedral angles of any simplex in $\mathcal{T}_{h}$ are less than or equal to $\frac{\pi}{2}$. Let $K_{k}$ have vertices $\{P^{k}_{j}\}_{j=0}^{d}$, and let $\eta^{k}_{j}(\boldsymbol{x})$ be the basis functions on $K_{k}$ associated with $\rm{Q}^{1}_{h}$ and $\mathrm{S}_{h}^{1}$, i.e. $\eta^{k}_{j}\mid_{K_{k}}\in\mathbb{P}_{1}$ and $\eta^{k}_{j}(P^{k}_{i})=\delta_{ij}$, $i,j=0,\ldots,d$. As $K_{k}$ is non-obtuse it follows that

[TABLE]

We then have the following result.

Lemma 6.1**.**

Let $g\in C^{0,1}(\mathbb{R})$ be monotonically increasing with Lipschitz constant $g_{\rm Lip}$. As $\mathcal{T}_{h}$ consists of only non-obtuse simplices, then we have for all $q\in{\rm Q}_{h}^{1}$, $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$ that

[TABLE]

Proof.

See the proof of Lemma 5.1 in [1]. ∎

Of course, the construction of a non-obtuse mesh in the case $d=3$ is not straightforward for a general polytope ${\mathcal{D}}$. However, we stress that our numerical method (P${}_{\alpha,\delta,h}^{\Delta t}$), see (6.36a–c) below, does not require this constraint. It is only required to show that (P${}_{\alpha,\delta,h}^{\Delta t}$) mimics the free energy structure of (Pα,δ).

6.2. A free energy preserving approximation, (P${}^{\Delta t}_{\alpha,\delta,h}$),

of (Pα,δ)

In addition to the assumptions on the finite element discretization stated in subsection 6.1, and our definition of $\Delta t$ in subsection 4.1, we shall assume for the convergence analysis, see Section 7, that there exists a $C\in{\mathbb{R}}_{>0}$ such that

[TABLE]

We note that this constraint is not required for the results in this section, in particular Theorem 6.1.

With $\Delta t_{1}$ and $C$ as above, let $\Delta t_{0}\in{\mathbb{R}}_{>0}$ be such that $\Delta t_{1}\leq C\Delta t_{0}$. Given initial data satisfying (3.18), we choose ${\boldsymbol{u}}_{h}^{0}\in\mathrm{V}_{h}^{1}$ and ${\boldsymbol{\sigma}}_{h}^{0}\in\mathrm{S}_{h}^{1}$ throughout the rest of this paper such that

[TABLE]

It follows from (6.30a,b), (6.8) and (3.18) that

[TABLE]

In addition, we note the following result.

Lemma 6.2**.**

For $p=1,\ldots,N_{P}$ we have that

[TABLE]

Proof.

For the proof of the first result in (6.32), see the proof of Lemma 5.2 in [1].

We now prove the second result in (6.32). On choosing $\boldsymbol{\chi}=\boldsymbol{I}\eta$, with $\eta\in{\rm Q}^{1}_{h}$, in (6.30b) yields that $z_{h}:=\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{0})-b^{\star}\in{\rm Q}^{1}_{h}$ satisfies

[TABLE]

where $z:=\operatorname{tr}({\boldsymbol{\sigma}}^{0})-b^{\star}\in L^{\infty}(\mathcal{D})$ and is non-positive on recalling (3.18). Choosing $\eta=\pi_{h}[z_{h}]_{+}\in{\rm Q}^{1}_{h}$, it follows, on noting the ${\rm Q}^{1}_{h}$ version of (6.8) and (6.28) with $g(\cdot)=[\,\cdot\,]_{+}$, that

[TABLE]

Hence $\pi_{h}[z_{h}]_{+}\equiv 0$ and so the second result in (6.32) holds. ∎

Furthermore, it follows from (6.30a,b), (6.31), (3.18), (6.3) and (6.9a,b) that, as $h,\,\Delta t_{0}\rightarrow 0_{+}$,

[TABLE]

Our approximation (P${}_{\alpha,\delta,h}^{\Delta t}$) of (Pα,δ) is then:

(P${}_{\alpha,\delta,h}^{\Delta t}$) Setting $({\boldsymbol{u}}_{\alpha,\delta,h}^{0},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{0},\varrho_{\alpha,\delta,h}^{0})=({\boldsymbol{u}}_{h}^{0},{\boldsymbol{\sigma}}_{h}^{0},\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{0}))\in\mathrm{V}_{h}^{1}\times(\mathrm{S}_{h}^{1}\cap\mathrm{S}_{>0,b})\times\mathrm{Q}_{h}^{1}$, with ${\boldsymbol{u}}_{h}^{0}$ and ${\boldsymbol{\sigma}}_{h}^{0}$ as defined in (6.30a,b), then for $n=1,\ldots,N_{T}$ find $({\boldsymbol{u}}_{\alpha,\delta,h}^{n},$ ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$ such that for any test functions $(\boldsymbol{v},\boldsymbol{\phi},\eta)\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$

[TABLE]

In deriving (P${}^{\Delta t}_{\alpha,\delta,h}$), we have noted (4.14), (1.5a), (6.21a,b) and (6.22). We note that on replacing $A_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\operatorname{tr}(\varrho_{\alpha,\delta,h}^{n}))$ with $\boldsymbol{I}-G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})$ and $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})$ by $1$ then (P${}_{\alpha,\delta,h}^{\Delta t}$), (6.36a,b), collapses to the corresponding finite element approximation of Oldroyd-B with stress diffusion studied in [1], see (5.34a,b) with no $L$ cut-off there.

Before proving existence of a solution to (P${}^{\Delta t}_{\alpha,\delta,h}$), we first derive a discrete analogue of the energy bound (5.8) for (Pα,δ).

6.3. Energy bound

On setting

[TABLE]

we have the following discrete analogue of Proposition 5.1.

Proposition 6.1**.**

For $n=1,\ldots,N_{T}$, a solution $\left({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n}\right)\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$ to (P${}_{\alpha,\delta,h}^{\Delta t}$), (6.36a–c), if it exists, satisfies

[TABLE]

Proof.

The proof is similar to that of Proposition 4.1, we choose as test functions $\boldsymbol{v}={\boldsymbol{u}}_{\alpha,\delta,h}^{n}\in\mathrm{V}_{h}^{1}$, $\boldsymbol{\phi}=-\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}[G_{\delta}^{\prime}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})]\in\mathrm{S}_{h}^{1}$ and $\eta=\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}\left[G_{\delta}^{\prime}\left(1-\frac{\varrho_{\alpha,\delta,h}^{n}}{b}\right)\right]\in\mathrm{Q}_{h}^{1}$ in (6.36a–c), and obtain, on noting (4.15), (3.4a,d), (3.8), (6.28) with $g=-G_{\delta}^{\prime}$ having Lipschitz constant $\delta^{-2}$, (6.23a,b) and (6.37) that

[TABLE]

The first desired inequality in (6.38) follows immediately from (6.39) on noting (6.1a,d), (1.12) and that $\pi_{h}:C(\overline{\mathcal{D}})\rightarrow{\rm Q}^{1}_{h}$. The second inequality in (6.38) follows immediately from (2.9) with $\nu^{2}=(1-\varepsilon)/(1+C_{P})$. ∎

6.4. Existence of discrete solutions

Proposition 6.2**.**

Given $({\boldsymbol{u}}_{\alpha,\delta,h}^{n-1},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n-1},\varrho_{\alpha,\delta,h}^{n-1})\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$ such that \int_{\mathcal{D}}[\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n-1})-\varrho_{\alpha,\delta,h}^{n-1}]\,{\rm d}\boldsymbol{x}$$=0 and for any time step $\Delta t_{n}>0$, then there exists at least one solution $\left({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n}\right)\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$ to (P${}_{\alpha,\delta,h}^{\Delta t}$), (6.36a–c), such that $\int_{\mathcal{D}}[\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})-\varrho_{\alpha,\delta,h}^{n}]\,{\rm d}\boldsymbol{x}=0$.

Proof.

The proof is similar to that of Proposition 4.2. We introduce the following inner product on the Hilbert space $\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$

[TABLE]

Given $({\boldsymbol{u}}_{\alpha,\delta,h}^{n-1},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n-1},\varrho_{\alpha,\delta,h}^{n-1})\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$, let $\mathcal{F}^{h}:\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}\rightarrow\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$ be such that for any $(\boldsymbol{w},\boldsymbol{\psi},\xi)\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$

[TABLE]

A solution $({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})$ to (4.13a–c), if it exists, corresponds to a zero of $\mathcal{F}^{h}$. On recalling (6.21a,b), (6.20a,b) and (6.16), it is easily deduced that the mapping $\mathcal{F}^{h}$ is continuous. For any $(\boldsymbol{w},\boldsymbol{\psi},\xi)\in\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}$, on choosing $(\boldsymbol{v},\boldsymbol{\phi},\eta)=\left(\boldsymbol{w},-\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}[G_{\delta}^{\prime}(\boldsymbol{\psi})],\frac{\varepsilon}{2{\text{Wi}}}\pi_{h}\left[G_{\delta}^{\prime}\left(1-\frac{\xi}{b}\right)\right]\right)$, we obtain analogously to (6.38) that

[TABLE]

Let

[TABLE]

If for any $\gamma\in\mathbb{R}_{>0}$, the continuous mapping $\mathcal{F}^{h}$ has no zero $({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})$, which lies in the ball

[TABLE]

then for such $\gamma$, we can define the continuous mapping $\mathcal{G}_{\gamma}^{h}:\mathcal{B}_{\gamma}^{h}\rightarrow\mathcal{B}_{\gamma}^{h}$ such that for all $(\boldsymbol{v},\boldsymbol{\phi},\eta)\in\mathcal{B}_{\gamma}^{h}$

[TABLE]

By the Brouwer fixed point theorem, $\mathcal{G}_{\gamma}^{h}$ has at least one fixed point $(\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma},\xi_{\gamma})$ in $\mathcal{B}_{\gamma}^{h}$. Hence it satisfies

[TABLE]

In addition, $\left((\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma},\xi_{\gamma}),(\boldsymbol{v},\boldsymbol{\phi},\eta)\right)_{\mathcal{D}}^{h}=\left(\mathcal{G}^{h}(\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma},\xi_{\gamma}),(\boldsymbol{v},\boldsymbol{\phi},\eta)\right)_{\mathcal{D}}^{h}$ with $(\boldsymbol{v},\boldsymbol{\phi},\eta)=({\boldsymbol{0}},\boldsymbol{I},-1)$ yields that

[TABLE]

On noting (6.4), we have that there exists a $\mu_{h}\in\mathbb{R}_{>0}$ such that for all $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$,

[TABLE]

and an equivalent result holding for all $\eta\in\mathrm{Q}_{h}^{1}$. It follows from (6.37), (6.43), (3.7), (3.9a), (6.44) and (6.42) that

[TABLE]

Hence for all $\gamma$ sufficiently large, it follows from (6.41) and (6.45) that

[TABLE]

On the other hand as $(\boldsymbol{w}_{\gamma},\boldsymbol{\psi}_{\gamma},\xi_{\gamma})$ is a fixed point of ${\mathcal{G}}_{\gamma}^{h}$, we have that

[TABLE]

It follows from (6.43), (3.7), (3.9b), and similarly to (6.45), on noting (6.44) and (6.42), that

[TABLE]

Therefore on combining (6.47) and (6.48), we have for all $\gamma$ sufficiently large that

[TABLE]

which obviously contradicts (6.46). Hence the mapping $\mathcal{F}^{h}$ has a zero, $({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})\in\mathcal{B}_{\gamma}^{h}$ for $\gamma$ sufficiently large. Finally, similarly to (6.43), it follows, on choosing $(\boldsymbol{v},\boldsymbol{\phi},\eta)=({\boldsymbol{0}},\boldsymbol{I},-1)$ in $\left(\mathcal{F}^{h}({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n}),(\boldsymbol{v},\boldsymbol{\phi},\eta)\right)_{\mathcal{D}}^{h}=0$, that $\int_{\mathcal{D}}[\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})-\varrho_{\alpha,\delta,h}^{n}]\,{\rm d}\boldsymbol{x}=\int_{\mathcal{D}}[\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n-1})-\varrho_{\alpha,\delta,h}^{n-1}]\,{\rm d}\boldsymbol{x}=0$. ∎

We now have the analogue of Theorem 4.1.

Theorem 6.1**.**

For any $\delta\in(0,\min\{\frac{1}{2},b\}]$, $N_{T}\geq 1$ and any partitioning of $[0,T]$ into $N_{T}$ time steps, there exists a solution $\{({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}]^{N_{T}}$ to (P${}^{\Delta t}_{\alpha,\delta,h}$), (6.36a–c).

In addition, it follows for $n=1,\ldots,N_{T}$ that $\int_{\mathcal{D}}[\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n})-\varrho_{\alpha,\delta,h}^{n}]\,{\rm d}\boldsymbol{x}=0$ and

[TABLE]

which yields that

[TABLE]

where $C$ is independent of $\alpha$, as well as $\delta$, $h$ and $\Delta t$.

Proof.

Existence and the stability result (6.50) follow immediately from Propositions 6.2 and 6.1, respectively, on noting (6.37), (3.1), (6.31), (6.32), (4.11a) and (3.18). The bounds (6.51) follow immediately from (6.50), on noting (3.8), (3.4b), (6.37), (3.7) and that

[TABLE]

∎

Remark 6.1**.**

We recall that we have used $\mathrm{S}_{h}^{1}$ for the approximation of ${\boldsymbol{\sigma}}_{\alpha,\delta}$ in (P${}_{\alpha,\delta,h}^{\Delta t}$), (6.36a–c), due to the presence of the diffusion term in (5.5c). Secondly, due to the advective term in (5.5c), one has to introduce the variable $\varrho_{\alpha,\delta,h}^{n}$ and its equation (6.36c) in (P${}_{\alpha,\delta,h}^{\Delta t}$) in order to obtain the entropy bound (6.38). However, we now have a bound on $\pi_{h}[[b-\varrho_{\alpha,\delta,h}^{n}]_{-}]$ in (6.51), as opposed to $[b-\operatorname{tr}({\boldsymbol{\sigma}}_{\delta,h}^{n})]_{-}$ in (4.37). Now, it does not seem possible to pass to the limit $\delta\rightarrow 0$ in (P${}_{\alpha,\delta,h}^{\Delta t}$) to prove well-posedness of the corresponding direct approximation of (Pα)**, i.e. (P${}_{\alpha,h}^{\Delta t}$) without the regularization $\delta$, as we did for (P${}_{\delta,h}^{\Delta t}$) in subsection 4.5.

Finally, we note the following Lemmas for later purposes.

Lemma 6.3**.**

For all $K_{k}\in{\mathcal{T}}_{h}$ and $\boldsymbol{\phi}\in[C(\overline{K_{k}})]^{d\times d}_{S}$, we have, for $r\in[1,\infty)$, that

[TABLE]

Proof.

It follows immediately from (6.4) that

[TABLE]

and hence the desired result(6.53). ∎

Lemma 6.4**.**

Let $g\in C^{0,1}({\mathbb{R}})$ with Lipschitz constant $g_{\rm Lip}$. For all $K_{k}\in{\mathcal{T}}_{h}$, and for all $q\in\mathrm{Q}_{h}^{1}$, $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$ we have that

[TABLE]

In addition, if $g$ is monotonic then, for all $K_{k}\in{\mathcal{T}}_{h}$ and for all $q\in\mathrm{Q}_{h}^{1}$, we have that

[TABLE]

Proof.

The results (6.54a) are proved in Lemma 5.3 of [1] for the case when $\beta_{\delta}$, $\Lambda_{\delta,m,p}$ and $\mathrm{S}_{h}^{1}$ are replaced by $\beta$, $\Lambda_{m,p}$ and $\mathrm{S}_{h}^{1}\cap\mathrm{S}_{>0}$. The proofs given there are trivially adapted to the present case. In fact, the proof of the first result in (6.54a) in [1] is easily adapted to any function $g\in C^{0,1}({\mathbb{R}})$. Hence, we have the results (6.54b).

The result (6.55) is a simple variation of (6.54b) and follows on noting that

[TABLE]

where $\{P^{k}_{j}\}_{j=0}^{d}$ are the vertices of $K_{k}$. ∎

7. Convergence of (P${}^{\Delta t}_{\alpha,\delta,h}$) to

(Pα) in the case $d=2$

7.1. Stability

Before proving our stability results, we introduce some further notation. We require the $L^{2}$ projector ${\mathcal{R}}_{h}:\mathrm{V}\rightarrow\mathrm{V}_{h}^{1}$ defined by

[TABLE]

In addition, we require ${\mathcal{P}}_{h}:[L^{2}(\mathcal{D})]^{d\times d}_{\rm S}\rightarrow\mathrm{S}_{h}^{1}$ defined by

[TABLE]

It is easily deduced for $p=1,\ldots,N_{P}$ and $i,\,j=1,\ldots,d$ that

[TABLE]

where $\eta_{p}\in{\rm Q}_{h}^{1}$ is such that $\eta_{p}(P_{r})=\delta_{pr}$ for $p,\,r=1,\ldots,N_{P}$. It follows from (7.2) and (6.7b) with $\boldsymbol{\phi}={\mathcal{P}}_{h}\boldsymbol{\chi}$, in both cases, that

[TABLE]

We shall assume from now on that $\mathcal{D}$ is convex and that the family $\{{\mathcal{T}}_{h}\}_{h>0}$ is quasi-uniform, i.e. $h_{k}\geq C\,h$, $k=1,\ldots,N_{K}$. It then follows that

[TABLE]

see Lemma 4.3 in Heywood and Rannacher [18]. Similarly, it is easily established that

[TABLE]

We also require the scalar analogue of ${\mathcal{P}}_{h}$, where $d=1$ and $\mathrm{S}_{h}^{1}$ is replaced by $\mathrm{Q}_{h}^{1}$, satisfying the analogues of (7.2)–(7.4) and (7.6).

Let $([H^{1}(\mathcal{D})]^{d\times d}_{\rm S})^{\prime}$ be the topological dual of $[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}$ with $[L^{2}(\mathcal{D})]^{d\times d}_{\rm S}$ being the pivot space. Let ${\mathcal{E}}:([H^{1}(\mathcal{D})]^{d\times d}_{\rm S})^{\prime}\rightarrow[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}$ be such that ${\mathcal{E}}\boldsymbol{\chi}$ is the unique solution of the Helmholtz problem

[TABLE]

where $\langle\cdot,\cdot\rangle_{H^{1}(\mathcal{D})}$ denotes the duality pairing between $([H^{1}(\mathcal{D})]^{d\times d}_{\rm S})^{\prime}$ and $[H^{1}(\mathcal{D})]^{d\times d}_{\rm S}$. We note that

[TABLE]

and $\|{\mathcal{E}}\cdot\|_{H^{1}(\mathcal{D})}$ is a norm on $([H^{1}(\mathcal{D})]^{d\times d}_{\rm S})^{\prime}$. We also employ this operator in the scalar case, $d=1$, i.e. ${\mathcal{E}}:H^{1}(\mathcal{D})^{\prime}\rightarrow H^{1}(\mathcal{D})$.

Let $\mathrm{V}^{\prime}$ be the topological dual of $\mathrm{V}$ with the space of weakly divergent free functions in $[L^{2}(\mathcal{D})]^{d}$ being the pivot space. Let ${\mathcal{S}}:\mathrm{V}^{\prime}\rightarrow\mathrm{V}$ be such that ${\mathcal{S}}\boldsymbol{w}$ is the unique solution to the Helmholtz-Stokes problem

[TABLE]

where $\langle\cdot,\cdot\rangle_{\mathrm{V}}$ denotes the duality pairing between $\mathrm{V}^{\prime}$ and $\mathrm{V}$. We note that

[TABLE]

and $\|{\mathcal{S}}\cdot\|_{H^{1}(\mathcal{D})}$ is a norm on the reflexive space $\mathrm{V}^{\prime}$.

We recall the following well-known Gagliardo-Nirenberg inequality. Let $r\in[2,\infty)$ if $d=2$, and $r\in[2,6]$ if $d=3$ and $\theta=d(\frac{1}{2}-\frac{1}{r})$. Then, there exists a positive constant $C(\mathcal{D},r,d)$ such that

[TABLE]

We recall also the following compactness results. Let ${\mathcal{Y}}_{0}$, ${\mathcal{Y}}$ and ${\mathcal{Y}}_{1}$ be real Banach spaces, ${\mathcal{Y}}_{i}$, $i=0,1$, reflexive, with a compact embedding ${\mathcal{Y}}_{0}\hookrightarrow{\mathcal{Y}}$ and a continuous embedding ${\mathcal{Y}}\hookrightarrow{\mathcal{Y}}_{1}$. Then, for $\mu_{i}>1$, $i=0,1$, the following embedding is compact:

[TABLE]

see Theorem 2.1 on p184 in Temam [25]. Let ${\mathcal{X}}_{0}$, ${\mathcal{X}}$ and ${\mathcal{X}}_{1}$ be real Hilbert spaces with a compact embedding ${\mathcal{X}}_{0}\hookrightarrow{\mathcal{X}}$ and a continuous embedding ${\mathcal{X}}\hookrightarrow{\mathcal{X}}_{1}$. Then, for $\gamma>0$, the following embedding is compact:

[TABLE]

where $D^{\gamma}_{t}\eta$ is the time derivative of order $\gamma$ of $\eta$, which can be defined in terms of the Fourier transform of $\eta$; see Theorem 2.2 on p186 in Temam [25].

Finally, we recall the discrete Gronwall inequality:

[TABLE]

Theorem 7.1**.**

Under the assumptions of Theorem 6.1, there exists a solution \{({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},$$\varrho_{\alpha,\delta,h}^{n})\}_{n=1}^{N_{T}}\in[\mathrm{V}_{h}^{1}\times\mathrm{S}_{h}^{1}\times\mathrm{Q}_{h}^{1}]^{N_{T}} of (P${}^{\Delta t}_{\alpha,\delta,h}$), (6.36a–c), such that the bounds (6.50) and (6.51) hold.

Moreover, if $d=2$, (6.29) holds and $\Delta t\leq C_{\star}(\zeta^{-1})\,\alpha^{1+\zeta}\,h^{2}$, for a $\zeta>0$ and a $C_{\star}(\zeta^{-1})\in\mathbb{R}_{>0}$ sufficiently small, then the following bounds hold:

[TABLE]

Proof.

Existence and the bounds (6.50) and (6.51) were proved in Theorem 6.1.

On choosing $\boldsymbol{\phi}={\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n}$ in (6.36b), it follows from (3.2), and (4.15), on applying a Young’s inequality, for any $\zeta>0$ that

[TABLE]

We deduce from (3.3), (1.8), (6.7b), (1.5b), (5.2), (5.3), (6.53) and (7.11), as $d=2$, that

[TABLE]

Similarly, as $d=2$, it follows from (6.26b), (3.3), (6.4) and (7.11) that for all $\zeta>0$

[TABLE]

In addition, as $d=2$, we note from (7.11), (1.13) and (6.51) that for all $\zeta>0$

[TABLE]

Combining (7.16)–(7.19), yields, on applying a Young’s inequality, that for all $\zeta>0$

[TABLE]

Hence, summing (7.20) from $n=1,\ldots,m$ for $m=1,\ldots,N_{T}$ yields, on noting (6.7b), that for any $\zeta>0$

[TABLE]

Applying the discrete Gronwall inequality (7.14) to (7.21), and noting (6.29), (6.8), (6.31), (6.5b), (6.50), (6.51), and that $\Delta t\leq C_{\star}(\zeta^{-1})\,\alpha^{1+\zeta}\,h^{2}$, for a $\zeta>0$ where $C_{\star}(\zeta^{-1})$ is sufficiently small, yields the bounds (7.15a).

Similarly to (7.16) on choosing $\eta=\varrho_{\alpha,\delta,h}^{n}$ in (6.36c), it follows from (4.15), (1.3b), (1.4b) and (1.8), on applying a Young’s inequality, for any $\zeta>0$ that

[TABLE]

Similarly to (7.17b), as $d=2$, we deduce from (6.7b), (5.2), (5.3), (3.3), (6.4), (7.11) and (7.15a) that

[TABLE]

Similarly to (7.18), as $d=2$, it follows from (6.26a), (3.3), (6.4) and (7.11) that for all $\zeta>0$

[TABLE]

Combining (7.22)–(7.24) and (7.19), yields, on applying a Young’s inequality, that for all $\zeta>0$

[TABLE]

Hence, summing (7.25) from $n=1,\ldots,m$ for $m=1,\ldots,N_{T}$ yields, on noting (6.7b) and (6.50), that for any $\zeta>0$

[TABLE]

Applying the discrete Gronwall inequality (7.14) to (7.26), and noting (6.29), (6.32), (6.5a,b), (6.50), (6.51), (6.8), (7.15a), and that $\Delta t\leq C_{\star}(\zeta^{-1})\,\alpha^{1+\zeta}\,h^{2}$, for a $\zeta>0$ where $C_{\star}(\zeta^{-1})$ is sufficiently small, yields the bounds (7.15b).

Finally, the desired result (7.15c) follows immediately from (7.18) and (7.24) with $\zeta=2$, on noting (7.15a,b) and (6.8). ∎

Remark 7.1**.**

Our final convergence result will be restricted to $d=2$ for the same reason as why our result for Oldroyd-B in [1] was restricted to $d=2$. For example, the control of the term (7.17b) necessitates the restriction to $d=2$. This could be overcome for FENE-P by replacing the regularization $\beta_{\delta}$ in (Pα,δ)**, (5.5a–h), by ${\beta_{\delta}^{b}}$ and using test functions based on $G_{\delta}^{b}$ in place of $G_{\delta}$, where $(G_{\delta}^{b})^{\prime}(s)={\beta_{\delta}^{b}}(s)$ for all $s\in\mathbb{R}$. On making a similar change to our numerical approximation (P${}_{\alpha,\delta}^{h,\Delta t}$), (6.36a–c), it is then possible to prove the analogues of Theorem 6.1 and Theorem 7.1 for $d=2$ and $3$ with no restriction on $\Delta t$. One can then establish analogues of Lemmas 7.1, 7.2 and Theorem 7.2, below, but now the limits involve the cut-off $b$, i.e. ${\boldsymbol{u}}_{\alpha}^{b}$, ${\boldsymbol{\sigma}}_{\alpha}^{b}$ and $\varrho_{\alpha}^{b}$. Moreover, it does not seem possible to establish that $\varrho_{\alpha}^{b}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha}^{b})$, as we now have $\Lambda_{\delta,m,p}^{b}$ in place of $\Lambda_{\delta,m,p}$ in the analogue of (7.60). Hence, in taking the limit $\delta,\,h,\,\Delta t\rightarrow 0_{+}$, the last term in the analogue of (7.61) would have $b\,\beta^{b}\left(1-\frac{\varrho_{\alpha}^{b}}{b}\right)+\operatorname{tr}(\beta^{b}({\boldsymbol{\sigma}}_{\alpha}^{b}))$ instead of $-(\varrho_{\alpha}^{b}-\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha}^{b}))$, where $\beta^{b}(s)=\min\{s,b\}$. Therefore, the extension of the convergence analysis in this paper to $d=3$ and with a weaker restriction on $\Delta t$ will be a topic of further research.

Lemma 7.1**.**

Under all of the assumptions of Theorem 7.1, the solution \{({\boldsymbol{u}}_{\alpha,\delta,h}^{n},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},$$\varrho_{\alpha,\delta,h}^{n})\}_{n=1}^{N_{T}} of (P${}^{\Delta t}_{\alpha,\delta,h}$), (6.36a–c), satisfies the following bounds:

[TABLE]

where $\vartheta\in(2,4]$ and $C$ in (7.27a,c) is independent of $\alpha$, as well as $\delta$, $h$ and $\Delta t$.

Proof.

On choosing $\boldsymbol{v}={\mathcal{R}}_{h}\left[{\mathcal{S}}\left(\frac{{\boldsymbol{u}}_{\alpha,\delta,h}^{n}-{\boldsymbol{u}}_{\alpha,\delta,h}^{n-1}}{\Delta t_{n}}\right)\right]\in\mathrm{V}_{h}^{1}$ in (6.36a) yields, on noting (7.1), (7.10), (7.5) and Sobolev embedding, that

[TABLE]

where $\theta>0$ as $d=2$. It follows from (6.7b) and (5.4) that

[TABLE]

Applying the Cauchy–Schwarz and the algebraic-geometric mean inequalities, in conjunction with (7.11), for $d=2$, and the Poincaré inequality (1.13) yields that

[TABLE]

Similarly, we have for any $\theta\in(0,1)$, as $d=2$, but now using a Young’s inequality

[TABLE]

On taking the $\frac{2}{\vartheta}$ power of both sides of (7.28), multiplying by $\Delta t_{n}$, summing from $n=1,\ldots,N_{T}$ and noting (7.29), (7.30), (7.31) with $\theta=\frac{\vartheta-2}{6-\vartheta}\Leftrightarrow\vartheta=\frac{2(1+3\theta)}{(1+\theta)}\in(2,4)$, (6.29), (4.11a), (6.50), (6.51) and (6.31) yields that

[TABLE]

and hence the bound (7.27a) with $C$ independent of $\alpha$.

Choosing $\boldsymbol{\phi}={\mathcal{P}}_{h}\left[{\mathcal{E}}\left(\frac{{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n}-{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n-1}}{\Delta t_{n}}\right)\right]\in\mathrm{S}_{h}^{1}$ in (6.36b) yields, on noting (7.2) and (7.8), that

[TABLE]

It follows from (1.5a,b), (1.3b), (1.4b) and (1.8) that for any $\zeta\in\mathbb{R}_{>0}$

[TABLE]

Similarly to (7.17b), it follows from (6.7b), (1.5b), (5.2), (5.3) and (1.8) that

[TABLE]

In addition, (3.3), (6.4), (6.8) and (7.15a) imply that for all $\boldsymbol{\phi}\in\mathrm{S}_{h}^{1}$

[TABLE]

Combining (7.33)–(7.36), yields, on noting (7.6), that

[TABLE]

Similarly to (7.33)–(7.37), choosing $\eta={\mathcal{P}}_{h}\left[{\mathcal{E}}\left(\frac{\varrho_{\alpha,\delta,h}^{n}-\varrho_{\alpha,\delta,h}^{n-1}}{\Delta t_{n}}\right)\right]\in\mathrm{Q}_{h}^{1}$ in (6.36c) yields that

[TABLE]

Multiplying (7.37) and (7.38) by $\Delta t_{n}$, summing from $n=1,...,N_{T}$ and noting (6.50), (6.51), (7.30), (6.29), (6.31) and (7.15a–c) yields the bounds (7.27b).

Multiplying (7.29) by $\Delta t_{n}$, summing from $n=1,...,N_{T}$ and noting (6.50) yields the result (7.27c) with $C$ independent of $\alpha$. Finally, it follows from (6.7a) that

[TABLE]

Multiplying (7.39) by $\Delta t_{n}$, summing from $n=1,...,N_{T}$ and noting (1.3b), (1.4b), (6.50), (1.8), (3.3), (6.53), (7.11) with $d=2$, (6.8) and (7.15a) yields that

[TABLE]

Hence, we have the desired result (7.27d). ∎

Unfortunately, the bound (7.27a) is not useful for obtaining compactness via (7.12), see the discussion in the proof of Theorem 7.2 below. Instead one has to exploit the compactness result (7.13). This we now do, by following the proof of Lemma 5.6 on p237 in Temam [25]. Here the introduction of $\kappa_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{n},\varrho_{\alpha,\delta,h}^{n})$ in the extra stress term, as discussed above (5.2), is crucial, as this yields an $L^{2}(\mathcal{D}_{T})$ bound on this stress term via the bound (5.4). For this purpose, we introduce the following notation in line with (4.10). Let ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t}\in C([0,T];\mathrm{V}_{h}^{1})$ and ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t,\pm}\in L^{\infty}(0,T;\mathrm{V}_{h}^{1})$ be such that for $n=1,\ldots,N_{T}$

[TABLE]

We note that

[TABLE]

where $t^{n}_{+}:=t^{n}$ and $t^{n}_{-}:=t^{n-1}$. We shall adopt ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t(,\pm)}$ as a collective symbol for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t}$, ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t,\pm}$. We also define ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t(,\pm)}$ and $\varrho_{\alpha,\delta,h}^{\Delta t(,\pm)}$ similarly to (7.41a,b).

Using the notation (7.41a,b), (P${}_{\alpha,\delta,h}^{\Delta t}$), i.e. (6.36a–c) multiplied by $\Delta t_{n}$ and summed for $n=1,\ldots,N_{T}$, can be restated as:

[TABLE]

subject to the initial conditions ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t}(0)={\boldsymbol{u}}_{h}^{0}$, ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t}(0)={\boldsymbol{\sigma}}_{h}^{0}$ and $\varrho_{\alpha,\delta,h}^{\Delta t}(0)=\operatorname{tr}({\boldsymbol{\sigma}}_{h}^{0})$.

Lemma 7.2**.**

Under all of the assumptions of Theorem 7.1, the solution $({\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t},\varrho_{\alpha,\delta,h}^{\Delta t})$ of (P${}^{\Delta t}_{\alpha,\delta,h}$), (7.43a–c), satisfies the following bound:

[TABLE]

where $\gamma\in(0,\textstyle\frac{1}{4})$ and $C$ is independent of $\alpha$, as well as $\delta$, $h$ and $\Delta t$.

Proof.

Equation (7.43a) can be reinterpreted as

[TABLE]

where $\boldsymbol{g}^{\Delta t,+}(t)\in\mathrm{V}_{h}^{1}$ is defined by

[TABLE]

Similarly to (7.28) and (7.29) with $\theta\in(0,1)$, it follows from (7.46) that

[TABLE]

On noting (7.30), (7.31) and the bound on ${\boldsymbol{u}}_{\alpha,\delta,h}^{n}$ in (6.51), we deduce from (7.47) that

[TABLE]

Similarly to (7.32), on recalling (6.31), (4.11a) and (6.51), we deduce from (7.48) that

[TABLE]

where $\vartheta\in(2,4]$ and $C$ is independent of $\alpha$, as well as $\delta$, $h$ and $\Delta t$. The rest of the proof follows as on p1825–6 in [1], which is based on the proof of Lemma 5.6 on p237 in [25]. ∎

7.2. Convergence

It follows from (6.50), (6.51), (7.15a–c), (6.31), (6.8), (7.27a–d), (7.44) and (7.41a–c) that

[TABLE]

where $\vartheta\in(2,4]$, $\gamma\in(0,\frac{1}{4})$ and $C$ in (7.50a,d,e,f,h) is independent of $\alpha$, as well as $\delta$, $h$ and $\Delta t$.

We are now in a position to prove the following convergence result concerning (P${}^{\Delta t}_{\alpha,\delta,h}$).

Theorem 7.2**.**

Under all of the assumptions of Theorem 7.1, there exists a subsequence of $\{({\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t},{\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t},\varrho_{\alpha,\delta,h}^{\Delta t})\}_{\delta>0,h>0,\Delta t>0}$, and functions

[TABLE]

such that, as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$,

[TABLE]

[TABLE]

and

[TABLE]

where $\vartheta\in(2,4]$ and $r\in[1,\infty)$.

Proof.

The results (7.52a–c) follow immediately from the bounds (7.50a,f) on noting the notation (7.41a–c). The denseness of $\bigcup_{h>0}{\rm Q}_{h}^{1}$ in $L^{2}(\mathcal{D})$ and (6.1d) yield that ${\boldsymbol{u}}_{\alpha}\in L^{2}(0,T;\mathrm{V})$. Hence, the result (7.51a) holds on noting (7.10) and (6.35), where ${\boldsymbol{u}}_{\alpha}:[0,T]\rightarrow{\rm H}$ is weakly continuous.

The strong convergence result (7.52d) for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t}$ and $r=2$ follows immediately from (7.50a) and the second bound in (7.50f) and (7.13) with ${\mathcal{X}}_{0}=[H^{1}(\mathcal{D})]^{2}$ and ${\mathcal{X}}={\mathcal{X}}_{1}=[L^{2}(\mathcal{D})]^{2}$. Here we note that $H^{1}(\mathcal{D})$ is compactly embedded in $L^{2}(\mathcal{D})$. We note here also that one cannot appeal to (7.12) for this strong convergence result with $\mu_{0}=2$, $\mu_{1}=4/\vartheta$, ${\mathcal{Y}}_{0}=[H^{1}({\mathcal{D}})]^{2}$, ${\mathcal{Y}}_{1}={\mathrm{V}}^{\prime}$ with norm $\|{\mathcal{S}}\cdot\|_{H^{1}(\mathcal{D})}$ and ${\mathcal{Y}}=[L^{2}({\mathcal{D}})]^{2}$ for the stated values of $\vartheta$, since $[L^{2}({\mathcal{D}})]^{2}$ is not continuously embedded in ${\mathrm{V}}^{\prime}$ as ${\mathrm{V}}$ is not dense in $[L^{2}({\mathcal{D}})]^{2}$.

The result (7.52d) for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t,\pm}$ and $r=2$ follows immediately from this result for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t}$ and the the bound on the last term on the left-hand side of (7.50a), which yields

[TABLE]

Finally, we note from (7.11), for $d=2$, that, for all $\eta\in L^{2}(0,T;H^{1}(\mathcal{D}))$,

[TABLE]

for all $r\in[2,\infty)$ with $\theta=1-\frac{2}{r}\in(0,1].$ Hence, combining (7.56) and (7.52b,d) for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t(,\pm)}$ with $r=2$ yields (7.52d) for ${\boldsymbol{u}}_{\alpha,\delta,h}^{\Delta t(,\pm)}$ for the stated values of $r$.

Similarly, the results (7.53a–c) follow immediately from (7.50b,g). The strong convergence result (7.53d) for ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t}$ follows immediately from (7.53b,h), (7.8) and (7.12) with $\mu_{0}=\mu_{1}=2$, ${\mathcal{Y}}_{0}=[H^{1}({\mathcal{D}})]^{d\times d}$, ${\mathcal{Y}}_{1}=([H^{1}({\mathcal{D}})]^{d\times d})^{\prime}$ and ${\mathcal{Y}}=[L^{r}({\mathcal{D}})]^{d\times d}$ for the stated values of $r$. Here we note that $H^{1}(\mathcal{D})$ is compactly embedded in $L^{r}(\mathcal{D})$ for the stated values of $r$, and $L^{r}(\mathcal{D})$ is continuously embedded in $(H^{1}(\mathcal{D}))^{\prime}$. Similarly to (7.55) and (7.56), the last bound in (7.50b) then yields that (7.53d) holds for ${\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t(,\pm)}$.

The results (7.54a–d) follow analogously from noting (7.50c,g). Hence, on noting (7.8) and (6.35), the results (7.51b,c) hold, where ${\boldsymbol{\sigma}}_{\alpha}:[0,T]\rightarrow[L^{2}(\mathcal{D})]^{d\times d}_{\rm S}$ and $\varrho_{\alpha}:[0,T]\rightarrow L^{2}(\mathcal{D})$ are weakly continuous, apart from the claims on the non-negative definiteness of ${\boldsymbol{\sigma}}_{\alpha}$ and the bound on $\varrho_{\alpha}$. It remains to prove these, (7.53e,f) and (7.54e,f). It follows from (1.11), (6.54b) and (7.50b,e) that

[TABLE]

The desired non-negative definiteness result on ${\boldsymbol{\sigma}}_{\alpha}$ in (7.51b) then follows from (7.53d). The desired bound on $\varrho_{\alpha}$ in (7.51c) follows similarly from (1.11), (6.54b), (7.50c,e) and (7.54d). The results (7.53e,f) follow immediately from (6.54a), (7.50b), (1.11), (7.53d), (3.3) and the non-negative definiteness result on ${\boldsymbol{\sigma}}_{\alpha}$ in (7.51b). The results (7.54e,f) follow similarly from the scalar version of (6.54a), (7.50c), (1.11), (7.54d), (3.3) and the bound on $\varrho_{\alpha}$ in (7.51c). ∎

Lemma 7.3**.**

Under all of the assumptions of Theorem 7.1, the subsequence of \{({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t},$$\varrho_{\alpha,\delta,h}^{\Delta t})\}_{\delta>0,h>0,\Delta t>0} of Theorem 7.2 and the limiting functions ${\boldsymbol{\sigma}}_{\alpha}$ and $\varrho_{\alpha}$, satisfying (7.51b,c), are such that, as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$,

[TABLE]

In addition, we have with $h=o(\delta)$, as $\delta\rightarrow 0_{+}$, that

[TABLE]

Proof.

Choosing $\boldsymbol{\phi}=\eta\,\boldsymbol{I}$ in (7.43b) and subtracting from (7.43c) yields that

[TABLE]

It follows from (7.52d), (7.53b–d,f), (7.54b–e), (7.7), (7.51a–c) and (6.9a,b) that we may pass to the limit $\delta,\,h,\,\Delta t\rightarrow 0_{+}$ in (7.60) with $\eta\ =\pi_{h}\,\chi$ to obtain

[TABLE]

where $[\varrho_{\alpha}-\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})](0)=0$. For example, in order to pass to the limit on the first term in (7.60), we note that

[TABLE]

Hence the desired first term in (7.61) follows from noting (7.53b,c), (7.54b,c), (7.7) and (6.9a,b). As $C^{\infty}_{0}(0,T;C^{\infty}(\overline{\mathcal{D}}))$ is dense in $L^{2}(0,T;H^{1}(\mathcal{D}))$, we have, on noting (7.51a–c), that (7.61) holds for all $\chi\in L^{2}(0,T;H^{1}(\mathcal{D}))$. It then follows from (7.51a–c) that we can choose $\chi=\varrho_{\alpha}-\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$ in (7.61) to yield that $\varrho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$ as $[\varrho_{\alpha}-\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})](0)=0$. Recalling (7.51b,c), we have that $\varrho_{\alpha}\in[0,b]$ a.e. in $\mathcal{D}_{T}$. The desired convergence result (7.58) then follows from noting this, (5.3), (3.3), (6.54b), (7.50c) and (7.54d).

We now improve on (7.51c) by establishing that $\varrho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})<b$ a.e. in $\mathcal{D}_{T}$. Assuming that $\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})=b$ a.e. in $\mathcal{D}_{T}^{b}\subset\mathcal{D}_{T}$, we have

[TABLE]

We deduce from (3.8) and (3.2) that

[TABLE]

It follows from (1.3b) and (1.4b) that

[TABLE]

Combining (7.64) and (7.65) yields, on noting a scalar version of (6.7a) over $D_{T}^{b}$, (7.50d), (1.8) and (3.3) that

[TABLE]

It follows from (6.53), (7.11), for $d=2$, and (7.50b) that

[TABLE]

Similarly to (6.53), it follows from (6.4) with $K_{k}$ replaced by $(K_{k}\times(t_{n-1},t_{n}))\cap\mathcal{D}^{b}_{T}$, $k=1,\ldots,N_{K}$ and $n=1,\ldots,N_{T}$, that

[TABLE]

where $z:=\pi_{h}\left[\beta_{\delta}\left(1-\frac{\varrho_{\alpha,\delta,h}^{\Delta t,+}}{b}\right)\right]-\left(1-\frac{\varrho_{\alpha}}{b}\right)$. Here we have noted that $\varrho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})=b$ a.e. in $\mathcal{D}_{T}^{b}$. Combining (7.66)–(7.68) yields, on noting (7.11), for $d=2$, that

[TABLE]

It follows from (7.69), (7.51c), (6.28), (3.3), (6.8), (7.50c) and (7.54e) that $T_{1}=0$. In addition, it follows immediately from (7.53e) that $T_{2}=0$. Hence, we conclude from (7.63) that $|\mathcal{D}^{b}_{T}|=0$, and so $\rho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})<b$ a.e. in $\mathcal{D}_{T}$; that is, the second desired result in (7.59).

We now establish the other result in (7.59) that ${\boldsymbol{\sigma}}_{\alpha}$ is symmetric positive definite a.e. in $\mathcal{D}_{T}$, which improves on (7.51b). This result requires the further assumption that $h=o(\delta)$, as $\delta\rightarrow 0_{+}$. Assume that ${\boldsymbol{\sigma}}_{\alpha}$ is not symmetric positive definite a.e. in $\mathcal{D}_{T}^{0}\subset\mathcal{D}_{T}$. Let $\boldsymbol{v}\in L^{\infty}(0,T;[L^{\infty}(\mathcal{D})]^{2})$ be such that ${\boldsymbol{\sigma}}_{\alpha}\,\boldsymbol{v}={\boldsymbol{0}}$ with $\|\boldsymbol{v}\|=1$ a.e. in $\mathcal{D}_{T}^{0}$ and $\boldsymbol{v}={\boldsymbol{0}}$ a.e. in $\mathcal{D}_{T}\setminus\mathcal{D}_{T}^{0}$. It then follows from (3.2) and (3.8) that

[TABLE]

A simple variation of (6.7a), (1.3b), (1.4a,b) and (7.50d) yield that

[TABLE]

Then (7.53e) and the definition of $\boldsymbol{v}$ yield that, as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$,

[TABLE]

so we have that $T_{3}=0$. Similarly to (7.71), on setting $\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}=G_{\delta}^{\prime}\left(1-\frac{\varrho_{\alpha,\delta,h}^{\Delta t,+}}{b}\right)\pi_{h}\left[\beta_{\delta}({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t,+})\right]$, we have from (6.7a) that

[TABLE]

where we have noted from (1.3b), (1.4b), (1.8), (6.53), (3.8) and (7.50i) that

[TABLE]

We will now show, on possibly extracting a further subsequence of \{({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t},$$\varrho_{\alpha,\delta,h}^{\Delta t})\}_{\delta>0,h>0,\Delta t>0}, that

[TABLE]

as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$ with $h=o(\delta)$. It follows immediately from (7.74) and our definition of $\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}$ that $\pi_{h}\left[\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}\right]$ converges weakly in $L^{\frac{8}{5}}(0,T;[L^{\frac{8}{5}}(\mathcal{D})]^{2\times 2})$ to some limit for a subsequence. We just need to show it is the limit stated in (7.75). We have from (7.53e), (7.54d) and (3.1) that $\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}\rightarrow\left(1-\frac{\rho_{\alpha}}{b}\right)^{-1}{\boldsymbol{\sigma}}_{\alpha}$ a.e. on $\mathcal{D}_{T}$, for a subsequence, as we have already established that $\rho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})<b$ a.e. on $\mathcal{D}_{T}$. So it remains to establish that $(I-\pi_{h})\left[\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}\right]$ converges to zero a.e. on $\mathcal{D}_{T}$. As $G_{\delta}^{\prime}\in C^{0,1}(\mathbb{R})$ is monotonic, it follows from (6.9b), (6.55), (6.4), (3.2) and (7.50b,c) that

[TABLE]

Hence, we have for a subsequence that $(I-\pi_{h})\left[\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}\right]$ converges to zero a.e. on $\mathcal{D}_{T}$ as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$ with $h=o(\delta)$. Therefore we have established (7.75).

Similarly to (7.72), we have that (7.73), (7.75) and our definitions of $\boldsymbol{\chi}_{\delta,h}^{\Delta t,+}$ and $\boldsymbol{v}$ yield that $T_{4}=0$. Hence it follows from (7.70) with $T_{3}=T_{4}=0$ that $|\mathcal{D}^{0}_{T}|=0$, and so ${\boldsymbol{\sigma}}_{\alpha}$ is positive definite a.e. in $\mathcal{D}_{T}$; that is, the first desired result in (7.59). ∎

Lemma 7.4**.**

Under all of the assumptions of Lemma 7.3, a further subsequence of the subsequence of $\{({\boldsymbol{\sigma}}_{\alpha,\delta,h}^{\Delta t},\varrho_{\alpha,\delta,h}^{\Delta t})\}_{\delta>0,h>0,\Delta t>0}$ of Lemma 7.3 and the limiting function ${\boldsymbol{\sigma}}_{\alpha}$, satisfying (7.51b) and (7.59), are such that, as $\delta,\,h,\,\Delta t\rightarrow 0_{+}$, with $h=o(\delta)$,

[TABLE]

Proof.

It follows from (6.7a), (1.3b), (1.4b), (5.2), (1.8), (3.3), (6.8) and (7.50b) that

[TABLE]

The desired result (7.77a) then follows immediately from (7.78), (7.53e) and (7.58).

The desired result (7.77b) follows immediately from (3.8), (3.2), (7.75) and (1.2) as $\rho_{\alpha}=\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$.

Similarly to (7.78), it follows from (6.7a), (1.3b), (1.4b), (5.2), (6.8), and (7.50d) that

[TABLE]

It follows immediately from (7.79), (7.53e), (7.58) and (7.50h) that the weak limits in (7.77b,c) are the same. Hence, the desired result (7.77c). ∎

Theorem 7.3**.**

Under all of the assumptions of Lemma 7.3, the limiting functions $({\boldsymbol{u}}_{\alpha},{\boldsymbol{\sigma}}_{\alpha})$ satisfying (7.51a,b) and (7.59) solve the following problem:

(Pα)* Find ${\boldsymbol{u}}_{\alpha}\in L^{\infty}(0,T;{\rm H})\cap L^{2}(0,T;\mathrm{V})\cap W^{1,\frac{4}{\vartheta}}(0,T;\mathrm{V}^{\prime})$ and ${\boldsymbol{\sigma}}_{\alpha}\in L^{\infty}(0,T;[L^{2}(\mathcal{D})]^{2\times 2}_{{\rm S},>0,b})\cap L^{2}(0,T;[H^{1}(\mathcal{D})]^{2\times 2}_{\rm S})\cap H^{1}(0,T;([H^{1}(\mathcal{D})]^{2\times 2}_{\rm S})^{\prime})$, with $A({\boldsymbol{\sigma}}_{\alpha})\,{\boldsymbol{\sigma}}_{\alpha}\in L^{2}(0,T;([L^{2}(\mathcal{D})]^{2\times 2}_{\rm S}))$, such that ${\boldsymbol{u}}_{\alpha}(0)={\boldsymbol{u}}^{0}$, ${\boldsymbol{\sigma}}_{\alpha}(0)={\boldsymbol{\sigma}}^{0}$ and*

[TABLE]

where $\vartheta\in(2,4)$.

Proof.

The function spaces and the initial conditions for $({\boldsymbol{u}}_{\alpha},{\boldsymbol{\sigma}}_{\alpha})$ follow immediately from (7.51a,b), (7.59) and (7.77c). It remains to prove that $({\boldsymbol{u}}_{\alpha},{\boldsymbol{\sigma}}_{\alpha})$ satisfy (7.80a,b). It follows from (6.3), (7.52b–d), (7.77c), (4.11b), (7.9) and (4.14) that we may pass to the limit, $\delta,\,h,\,\Delta t\rightarrow 0_{+}$, with $h=o(\delta)$, in (7.43a) to obtain that $({\boldsymbol{u}}_{\alpha},{\boldsymbol{\sigma}}_{\alpha})$ satisfy (7.80a).

It follows from (7.53b–f), (7.52b,d), (7.77a,b), (7.7), (6.9a,b), (1.5a) and as ${\boldsymbol{u}}_{\alpha}\in L^{2}(0,T;{\rm V})$ that we may pass to the limit $\delta,\,h,\,\Delta t\rightarrow 0_{+}$, with $h=o(\delta)$, in (7.43b) with $\boldsymbol{\chi}=\pi_{h}\boldsymbol{\phi}$ to obtain (7.80b) for any $\boldsymbol{\phi}\in C^{\infty}_{0}(0,T;[C^{\infty}(\overline{\mathcal{D}})]^{2\times 2}_{\rm S})$. For example, similarly to (7.62), in order to pass to the limit on the first term in (7.43b), we note that

[TABLE]

The desired result (7.80b) then follows from noting that $C^{\infty}_{0}(0,T;[C^{\infty}(\overline{\mathcal{D}})]^{2\times 2}_{\rm S})$ is dense in $L^{2}(0,T;$ $[H^{1}(\mathcal{D})]^{2\times 2}_{\rm S})$.

Of course passing to the limit $\delta,\,h,\,\Delta t\rightarrow 0_{+}$, with $h=o(\delta)$, in (7.43c), using in addition (7.54b–d,f), yields the weak formulation for $\operatorname{tr}({\boldsymbol{\sigma}}_{\alpha})$ consistent with (7.80b). ∎

Remark 7.2**.**

Choosing $\frac{1}{2}(\boldsymbol{\phi}+\boldsymbol{\phi}^{T})$ as a test function in (7.80b) for any $\boldsymbol{\phi}\in L^{2}(0,T;[H^{1}(\mathcal{D})]^{2\times 2})$ yields, on noting the symmetry of ${\boldsymbol{\sigma}}_{\alpha}$, (7.80b) with the term $2\,(\boldsymbol{\nabla}{\boldsymbol{u}}_{\alpha}){\boldsymbol{\sigma}}_{\alpha}$ replaced by $(\boldsymbol{\nabla}{\boldsymbol{u}}_{\alpha}){\boldsymbol{\sigma}}_{\alpha}+{\boldsymbol{\sigma}}_{\alpha}(\boldsymbol{\nabla}{\boldsymbol{u}}_{\alpha})^{T}$, which is consistent with (5.1c).

Finally, it follows from (7.50a,f,h), (7.52a–c) and (7.77c) that

[TABLE]

where $\vartheta\in(2,4)$ and $C$ is independent of the stress diffusion coefficient $\alpha$. Of course, in addition, it follows from (7.59) and (1.8) that $\|{\boldsymbol{\sigma}}_{\alpha}\|_{L^{\infty}(\mathcal{D}_{T})}<b$.