A dual-phase-lag porous-thermoelastic problem with microtemperatures

N. Bazarra; J. R. Fernández; R. Quintanilla; N. Bazarra; J. R. Fernández; R. Quintanilla

doi:10.3934/era.2022065

Electronic Research Archive

2022, Volume 30, Issue 4: 1236-1262. doi: 10.3934/era.2022065

Previous Article Next Article

Research article

A dual-phase-lag porous-thermoelastic problem with microtemperatures

1.
Departamento de Matemática Aplicada I, Universidade de Vigo, Escola de Enxeñería de Telecomunicación, Campus As Lagoas Marcosende s/n, 36310 Vigo, Spain
2.
Departamento de Matemáticas, E.S.E.I.A.A.T.-U.P.C., Colom 11, 08222 Terrassa, Barcelona, Spain

Received: 21 September 2021 Revised: 10 December 2021 Accepted: 16 December 2021 Published: 15 March 2022

In this work, we consider a multi-dimensional dual-phase-lag problem arising in porous-thermoelasticity with microtemperatures. An existence and uniqueness result is proved by applying the semigroup of linear operators theory. Then, by using the finite element method and the Euler scheme, a fully discrete approximation is numerically studied, proving a discrete stability property and a priori error estimates. Finally, we perform some numerical simulations to demonstrate the accuracy of the approximation and the behavior of the solution in one- and two-dimensional problems.

Keywords:

dual-phase-lag,
porous-thermoelasticity with microtemperatures,
existence and uniqueness,
finite elements,
a priori error estimates,
numerical simulations

Citation: N. Bazarra, J. R. Fernández, R. Quintanilla. A dual-phase-lag porous-thermoelastic problem with microtemperatures[J]. Electronic Research Archive, 2022, 30(4): 1236-1262. doi: 10.3934/era.2022065

Related Papers:

[1]	Dorin Ieşan, Ramón Quintanilla . Second gradient thermoelasticity with microtemperatures. Electronic Research Archive, 2025, 33(2): 537-555. doi: 10.3934/era.2025025
[2]	Antonio Magaña, Alain Miranville, Ramón Quintanilla . On the time decay in phase–lag thermoelasticity with two temperatures. Electronic Research Archive, 2019, 27(0): 7-19. doi: 10.3934/era.2019007
[3]	Nisachon Kumankat, Kanognudge Wuttanachamsri . Well-posedness of generalized Stokes-Brinkman equations modeling moving solid phases. Electronic Research Archive, 2023, 31(3): 1641-1661. doi: 10.3934/era.2023085
[4]	Noelia Bazarra, José R. Fernández, Ramón Quintanilla . On the mixtures of MGT viscoelastic solids. Electronic Research Archive, 2022, 30(12): 4318-4340. doi: 10.3934/era.2022219
[5]	Wenya Qi, Padmanabhan Seshaiyer, Junping Wang . A four-field mixed finite element method for Biot's consolidation problems. Electronic Research Archive, 2021, 29(3): 2517-2532. doi: 10.3934/era.2020127
[6]	Noelia Bazarra, José R. Fernández, Ramón Quintanilla . Numerical analysis of a problem in micropolar thermoviscoelasticity. Electronic Research Archive, 2022, 30(2): 683-700. doi: 10.3934/era.2022036
[7]	Hsueh-Chen Lee, Hyesuk Lee . An a posteriori error estimator based on least-squares finite element solutions for viscoelastic fluid flows. Electronic Research Archive, 2021, 29(4): 2755-2770. doi: 10.3934/era.2021012
[8]	Shan Jiang, Li Liang, Meiling Sun, Fang Su . Uniform high-order convergence of multiscale finite element computation on a graded recursion for singular perturbation. Electronic Research Archive, 2020, 28(2): 935-949. doi: 10.3934/era.2020049
[9]	Yiyuan Qian, Haiming Song, Xiaoshen Wang, Kai Zhang . Primal-dual active-set method for solving the unilateral pricing problem of American better-of options on two assets. Electronic Research Archive, 2022, 30(1): 90-115. doi: 10.3934/era.2022005
[10]	Jiwei Jia, Young-Ju Lee, Yue Feng, Zichan Wang, Zhongshu Zhao . Hybridized weak Galerkin finite element methods for Brinkman equations. Electronic Research Archive, 2021, 29(3): 2489-2516. doi: 10.3934/era.2020126

Abstract

1. Introduction

In the second part of the past century, there were developed different theories where the microstructure of elastic materials was taken into account. Two good books concerning these aspects were written by Eringen ^[1] and Ieşan ^[2]. In particular, Goodman and Cowin ^[3] established the theory of continuum for granular materials with intersticial voids. They proposed the bulk density of the materials as a product of the matrix density and the volume fraction. Later, Cowin and Nunziato ^[4,5,6] stated a theory of porous elastic solids modeling the deformations of materials with small voids distributed within them. Of course, this theory was extended to include the thermal effects ^[7]. Contributions in this theory became huge and we can recall a few of them ^{[8,9,10,11,12,13,14,15,16,17,18,19]}. In general, the microstructure components of elastic materials have deserved much attention. One of the possibilities could be the so-called "microtemperatures" ^[20,21].

If we revise the literature regarding the heat conduction we should agree that most of the contributions consider the Fourier's law. That is, the heat flux vector can be obtained as a linear expression of the gradient of temperature. However, it is also well-known that this assumption brings us to the fact that the thermal waves propagate instantaneously. Therefore, the causality principle is not satisfied. This paradox has motivated many scientists to look for alternative constitutive equations replacing Fourier's law in the heat flux vector. Perhaps the most known is the proposition of Cattaneo and Maxwell ^[22], who suggested to introduce a relaxation parameter to the Fourier law, bringing to a damped hyperbolic equation describing the heat conduction. Some other authors have proposed other alternatives, but in this paper we focus our attention in the Tzou theory ^[23]. In this case, the author introduced two relaxation parameters and the theory of Cattaneo and Maxwell can be seen as a particular case. Many contributions have been obtained in the context of the dual-phase-lag thermoelasticity ^[24,25,26]. In a very recent contribution, it was suggested a way to extend the dual-phase-lag theory to include the microtemperatures ^[27,28]. In fact, several contributions have been obtained concerning the stability of the problem.

Here, we center our attention to the case of porous-thermo-elastic solids with microtemperatures within the context of the dual-phase-lag theory. For this reason, we consider the case where the heat and the microheat are determined by means of the dual-phase-lag. Therefore, in this work we first show the existence of solutions in the multi-dimensional setting and then, we propose a numerical study of the problem, proving a discrete stability property and a priori error estimates, and performing some numerical simulations.

The paper is structured in the following form. In the next section, we describe the model and the basic assumptions. Then, in Section 3 we prove that the thermomechanical problem has a unique solution in the general multi-dimensional setting by means of the semigroup of linear operators theory. The numerical approximation of this problem is presented in Section 4, by using the finite element method and the implicit Euler scheme to approximate the spatial variable and to discretize the time derivatives, respectively. A discrete stability property and a priori error estimates are shown. Finally, some numerical simulations are described in Section 5 to demonstrate the accuracy of the algorithm and the behavior of the solution.

2. The thermomechanical model

We recall in this section the basic equations and assumptions under we will work in this paper. We will consider a domain $B\subset\mathbb{R}^d$ , $d = 1, 2, 3$ , with a boundary $\Gamma = \partial B$ assumed to be smooth enough.

The basic equations for the dual-phase-lag thermoelastic problem with microtemperatures are given by the evolution equations (see ^[28]):

$\begin{eqnarray*} &&\rho \ddot {u}_i = t_{ij,j},\\ && J \ddot{\phi} = h_{i,i}+g,\\ && \rho T_0 \dot {\eta} = q_{i,i},\\ &&\rho \dot{\epsilon}_i = P_{ij,j}+q_i-Q_i, \end{eqnarray*}$

and the constitutive equations^*:

^*Though we could consider different relaxation parameters for the macrotemperatures and the microtemperatures, we assume here that they agree. The reason is basically the big difficulty to treat with the general case. We do recognize that we are not able to deal with it.

$\begin{eqnarray*} &&t_{ij} = \lambda e_{rr}\delta_{ij}+2\mu e_{ij}+\mu_0\phi \delta_{ij}-\beta_{0}\theta \delta_{ij}, \end{eqnarray*}$

$\begin{eqnarray*} &&h_i = a\phi_{,i}-\mu_2 T_i,\\ &&g = -\mu_0 e_{ii}-\zeta \phi+\beta_1\theta,\\ &&\rho \eta = \beta_0e_{ii}+\beta_1\phi+a\theta,\\ && \rho \epsilon_i = -\mu_2\phi_{,i}-b T_i,\\ && q_i+\tau_1 \dot{q}_i+\frac{\tau_1^2}{2}\ddot{q}_i = (k\theta_{,i}+\kappa_1 T_i)+\tau_2 (k\dot{\theta}_{,i}+\kappa_1 \dot{T}_i),\\ && P_{ij}+\tau_1 \dot{P}_{ij}+\frac{\tau_1^2}{2}\ddot{P}_{ij} = -\kappa_4 T_{r,r}\delta_{ij}-\kappa_5T_{j,i}-\kappa_6 T_{i,j}\\ &&\qquad -\tau_2(\kappa_4 \dot{T}_{r,r}\delta_{ij}+\kappa_5\dot{T}_{j,i}+\kappa_6 \dot{T}_{i,j}),\\ && Q_i+\tau_1 \dot{Q}_{i}+\frac{\tau_1^2}{2}\ddot{Q}_{i} = (\kappa-\kappa_1)\theta_{,i}+(\kappa_1-\kappa_2) T_{i}\\ &&\qquad +\tau_2((\kappa-\kappa_1) \dot{\theta}_{,i}+(\kappa_1-\kappa_2) \dot{T}_{i}), \end{eqnarray*}$

where $\mathit{\boldsymbol{u}} = (u_i)$ is the displacement field, $\phi$ is the porosity function, $\theta$ is the temperature and $\mathit{\boldsymbol{T}} = (T_i)$ is the microtemperatures vector. Moreover, $\bf{e} = (e_{ij})$ denotes the linearized strain tensor given by

$e_{ij} = \frac{1}{2}(u_{i,j}+u_{j,i}),$

$t_{ij}$ is the stress tensor, $h_i$ is the equilibrated stress, $g$ is the equilibrated body force, $\eta$ is the entropy, $\epsilon_i$ is the first moment of the energy, $P_{ij}$ is the first heat flux moment, $Q_i$ is the microheat flux average, $\tau_1$ and $\tau_2$ are the time relaxation parameters, and $\rho,$ $J$ , $\lambda$ , $\mu$ , $\mu_0$ , $\beta_0$ , $a$ , $\mu_2$ , $\zeta$ , $\beta_1$ , $b$ , $\kappa_1$ , $\kappa_4$ , $\kappa_5$ , $\kappa_6$ , $\kappa_2$ and $\kappa$ are constitutive parameters. From now on, we assume that $T_0 = 1$ to simplify the calculations. Anyway, this assumption is not restrictive.

After substitution of the constitutive equations into the evolution equations we obtain the following linear system:

$\begin{eqnarray*} &&\rho \ddot {u}_i = \mu u_{i,jj}+(\lambda+\mu)u_{j,ji}+\mu_0\phi_{,i}-\beta_0\theta_{,i},\\ && J\ddot{\phi} = a_0\phi_{,ii}-\mu_2 T_{i,i}-\mu_0 u_{i,i}-\zeta\phi+\beta_1\theta,\\ && a(\dot{\theta}+\tau_1\ddot{\theta}+\frac{\tau_1^2}{2}\dddot{\theta}) = \kappa(\theta_{,ii}+\tau_2\dot{\theta}_{,ii})+ \kappa_1(T_{i,i}+\tau_2\dot{T}_{i,i})\\ && \qquad -\beta_0(\dot{u}_{i,i}+\tau_1\ddot{u}_{i,i}+\frac{\tau_1^2}{2}\dddot{u}_{i,i}) -\beta_1(\dot{\phi}+\tau_1 \ddot{\phi}+\frac{\tau_1^2}{2}\dddot{\phi}), \\ &&b(\dot{T}_i+\tau_1\ddot{T}_i+\frac{\tau_1^2}{2}\dddot{T}_i) = \kappa_6(T_{i,jj}+\tau_2\dot{T}_{i,jj})+ (\kappa_4+\kappa_5)(T_{j,ji}+\tau_2\dot{T}_{j,ji})\\ && \qquad -\kappa_2(T_i+\tau_2\dot{T}_i)-\kappa_1(\theta_{,i}+\tau_2\dot{\theta}_{,i}) -\mu_2(\dot{\phi}_{,i}+\tau_1 \ddot{\phi}_{,i}+\frac{\tau_1^2}{2}\dddot{\phi}_{,i}). \end{eqnarray*}$

Proceeding as in ^[28], defining $\hat{f} = f+\tau_1\dot{f}+\frac{\tau_1^2}{2}\ddot{f}$ we can write the previous system as follows:

$\begin{eqnarray} &&\rho \ddot {\hat{u}}_i = \mu \hat{u}_{i,jj}+(\lambda+\mu)\hat{u}_{j,ji}+\mu_0\hat{\phi}_{,i}-\beta_0(\theta_{,i}+\tau_1 \dot{\theta}_{,i}+\frac{\tau_1^2}{2}\ddot{\theta}_{,i}), \end{eqnarray}$

(2.1)

$\begin{eqnarray} &&J\ddot{\hat{\phi}} = a_0\hat{\phi}_{,ii}-\mu_2(T_{i,i}+\tau_1 \dot{T}_{i,i}+\frac{\tau_1^2}{2}\ddot{T}_{i,i})-\mu_0 \hat{u}_{i,i}-\zeta\hat{\phi}+\beta_1(\theta+\tau_1 \dot{\theta}+\frac{\tau_1^2}{2}\ddot{\theta}), \end{eqnarray}$

(2.2)

$\begin{eqnarray} &&a(\dot{\theta}+\tau_1\ddot{\theta}+\frac{\tau_1^2}{2}\dddot{\theta}) = \kappa(\theta_{,ii}+\tau_2\dot{\theta}_{,ii})+ \kappa_1(T_{i,i}+\tau_2\dot{T}_{i,i}) -\beta_0\dot{\hat{u}}_{i,i} -\beta_1\dot{\hat{\phi}}, \end{eqnarray}$

(2.3)

$\begin{eqnarray} &&b(\dot{T}_i+\tau_1\ddot{T}_i+\frac{\tau_1^2}{2}\dddot{T}_i) = \kappa_6(T_{i,jj}+\tau_2\dot{T}_{i,jj})+ (\kappa_4+\kappa_5)(T_{j,ji}+\tau_2\dot{T}_{j,ji})\\ &&\qquad -\kappa_2(T_i+\tau_2\dot{T}_i)-\kappa_1(\theta_{,i}+\tau_2\dot{\theta}_{,i}) -\mu_2\dot{\hat{\phi}}_{,i}. \end{eqnarray}$

(2.4)

However, in order to simplify the writing, we remove the superscript hat over all the variables.

Since we assume homogeneous Dirichlet boundary conditions, it follows that

$\begin{equation} u_i( \mathit{\boldsymbol{x}},t) = \phi ( \mathit{\boldsymbol{x}},t) = \theta( \mathit{\boldsymbol{x}},t) = T_i( \mathit{\boldsymbol{x}},t) = 0, \quad t\in [0,\infty), \quad \mathit{\boldsymbol{x}} \in \partial B, \end{equation}$

(2.5)

and, to define a well posed problem, we will need to prescribe some initial conditions, for a.e. $\mathit{\boldsymbol{x}}\in B$ ,

$\begin{equation} \begin{array}{l} u_i( \mathit{\boldsymbol{x}},0) = u_i^0( \mathit{\boldsymbol{x}}),\quad \dot{u}_i( \mathit{\boldsymbol{x}},0) = v_i^0( \mathit{\boldsymbol{x}}),\quad \phi( \mathit{\boldsymbol{x}},0) = \phi^0( \mathit{\boldsymbol{x}}), \quad \dot{\phi}( \mathit{\boldsymbol{x}},0) = \psi^0( \mathit{\boldsymbol{x}}),\\ \theta( \mathit{\boldsymbol{x}},0) = \theta^0( \mathit{\boldsymbol{x}}),\quad T_i( \mathit{\boldsymbol{x}},0) = T_i^0( \mathit{\boldsymbol{x}}),\quad \dot{\theta}( \mathit{\boldsymbol{x}},0) = \vartheta^0( \mathit{\boldsymbol{x}}),\quad \dot{T}_i( \mathit{\boldsymbol{x}},0) = S_i^0( \mathit{\boldsymbol{x}}), \\ \ddot{\theta}( \mathit{\boldsymbol{x}},0) = \xi^0( \mathit{\boldsymbol{x}}),\quad \ddot{T}_i( \mathit{\boldsymbol{x}},0) = R_i^0( \mathit{\boldsymbol{x}}). \end{array} \end{equation}$

(2.6)

In this paper, we are going to assume that

$\begin{equation} \rho,\, J,\, a,\, b,\, \mu,\, a_0,\, \kappa,\, \kappa_6,\, \kappa_4+\kappa_5,\, \tau_1,\, \lambda+\mu,\, (\lambda+\mu)\zeta-\mu_0^2,\, \kappa\kappa_2-\kappa_1^2, \, \tau_2-\tau_1/2 \end{equation}$

(2.7)

are positive numbers. These assumptions are the natural ones for this theory from the mechanical and thermal points of view.

Let us define the internal energy $\mathcal{E}$ and the dissipation $\bf{D}$ as follows:

$\begin{array}{l} \mathcal{E}(t) = \frac{1}{2}\int_B \Big[ \rho \dot{\hat{u}}_i\dot{\hat{u}}_i+J|\dot{\hat{\phi}}|^2+2\mu_0\hat{u}_{i,i}\hat{\phi}+\zeta |\hat{\phi}|^2+a\hat{\phi}_{,i}\hat{\phi}_{,i}+\mu \hat{u}_{i,j}\hat{u}_{i,j}+(\lambda+\mu)\hat{u}_{i,i}\hat{u}_{j,j}\Big]\,d \mathit{\boldsymbol{x}}\\ \qquad +\int_B\Big[a|\hat{\theta}|^2+b\hat{T}_i\hat{T}_i+\kappa(\tau_1+\tau_2)|\nabla\theta|^2 + \frac{\kappa\tau_1^2\tau_2}{2}|\nabla\dot{\theta}|^2+\kappa_2(\tau_1+\tau_2)T_iT_i\\ \qquad + \frac{\kappa_2\tau_1^2\tau_2}{2}\dot{T}_i\dot{T}_i+\kappa\tau_1^2\theta_{,i}\dot{\theta}_{,i} +\kappa_2\tau_1^2T_i\dot{T}_i+2(\tau_1+\tau_2)\kappa_1\theta_{,i}T_i +\kappa_1\tau_1^2\tau_2\dot{T}_i\dot{\theta}_{,i}\\ \qquad +\kappa_1\tau_1^2(\theta_{,i}\dot{T}_i+\dot{\theta}_{,i}T_i) +\kappa_6(\tau_1+\tau_2)T_{i,j}T_{i,j}+\kappa_6\frac{\tau_1^2\tau_2}{2}\dot{T}_{i,j}\dot{T}_{i,j} +\kappa_6\tau_1^2T_{i,j}\dot{T}_{i,j}\\ \qquad +(\kappa_4+\kappa_5)(\tau_1+\tau_2)T_{i,i}T_{j,j} +(\kappa_4+\kappa_5)\frac{\tau_1^2\tau_2}{2}\dot{T}_{i,i}\dot{T}_{j,j}+(\kappa_4+\kappa_5)\tau_1^2T_{i,i}\dot{T}_{j,j}\Big]\, d \mathit{\boldsymbol{x}},\\ \bf{D}(t) = \int_B\Big[\kappa( |\nabla\theta|^2+(\tau_1\tau_2-\frac{\tau_1^2}{2})|\nabla \dot{\theta}|^2)+\kappa_2(T_iT_i+(\tau_1\tau_2-\frac{\tau_1^2}{2})\dot{T}_i\dot{T}_i)\\ \qquad +2\kappa_1(\theta_{,i}T_i+(\tau_1\tau_2-\frac{\tau_1^2}{2})\dot{T}_i\dot{\theta}_{,i}) +\kappa_6(T_{i,j}T_{i,j}+(\tau_1\tau_2-\frac{\tau_1^2}{2})\dot{T}_{i,j}\dot{T}_{i,j})\\ \qquad +(\kappa_4+\kappa_5)(T_{i,i}T_{j,j}+(\tau_1\tau_2-\frac{\tau_1^2}{2}) \dot{T}_{i,i}\dot{T}_{j,j})\Big]\, d \mathit{\boldsymbol{x}}. \end{array}$

It can be proved that, for a.e. $t\geq 0$ ,

$\mathcal{E}(t)+\int_0^t \bf{D}(s)\,ds = \mathcal{E}(0).$

In view of the previous assumptions, $\mathcal{E}(t)$ and $\bf{D}(t)$ are positive definite functions, and the last equality shows the stability of the solutions.

3. Existence and uniqueness

In this section, we give an existence and uniqueness result for the problem determined by system (2.1)-(2.4), boundary conditions (2.5) and initial conditions (2.6).

To this end we will work in the Hilbert space:

$\begin{array}{l} \mathcal{H} = \mathit{\boldsymbol{W}}_0^{1,2}(B)\times \mathit{\boldsymbol{L}}^2(B)\times W_0^{1,2}(B)\times L^2(B)\times W_0^{1,2}(B)\times W_0^{1,2}(B)\times L^2(B)\\ \qquad \qquad \times \mathit{\boldsymbol{W}}^{1,2}_0(B)\times \mathit{\boldsymbol{W}}_0^{1,2}(B)\times \mathit{\boldsymbol{L}}^2(B). \end{array}$

If we denote by

$U = ( \mathit{\boldsymbol{u}}, \mathit{\boldsymbol{v}},\phi,\psi,\theta,\vartheta,\xi, \mathit{\boldsymbol{T}}, \mathit{\boldsymbol{S}}, \mathit{\boldsymbol{R}}),\quad U^* = ( \mathit{\boldsymbol{u}}^*, \mathit{\boldsymbol{v}}^*,\phi^*,\psi^*,\theta^*,\vartheta^*,\xi^*, \mathit{\boldsymbol{T}}^*, \mathit{\boldsymbol{S}}^*, \mathit{\boldsymbol{R}}^*),$

we define the inner product:

$\begin{array}{l} < U,U^* > = \frac{1}{2}\int_B\Big[ \rho v_iv_i^*+J\psi\psi^*+\mu_0(u_{i,i}\phi^*+u_{i,i}^*\phi)+a\phi_{,i}\phi_{,i}^*+\mu u_{i,j}u_{i,j}^*\\ \qquad +(\lambda+\mu)u_{i,i}u_{j,j}^*+a\hat{\theta}\hat{\theta}^*+b\hat{T}_i\hat{T}_i^* +\kappa (\tau_1+\tau_2)\theta_{,i}\theta_{,i}^*+\frac{\kappa\tau_1^2\tau_2}{2}\vartheta_{,i}\vartheta_{,i}^*\\ \qquad + \kappa_2(\tau_1+\tau_2)T_iT_i^*+\frac{\kappa\tau_1^2\tau_2}{2}S_iS_i^*+\frac{\kappa \tau_1^2}{2}(\theta_{,i}\vartheta_{,i}^*+\theta_{,i}^*\vartheta_{,i})+\frac{\kappa_2 \tau_1^2}{2}(T_iS_i^*+T_iS_i^*)\\ \qquad +(\tau_1+\tau_2)\kappa_1(\theta_{,i}T_i^*+\theta_{,i}^*T_i) + \frac{\kappa_1 \tau_1^2}{2}(\theta_{,i}S_i^*+\theta_{,i}^*S_i+\vartheta_{,i}T_i^*+\vartheta_{,i}^*T_i)\\ \qquad +\frac{\kappa_1 \tau_1^2\tau_2}{2}(\vartheta_{,i}S_i^*+\vartheta_{,i}^* S_i)+\kappa_6(\tau_1+\tau_2)T_{i,j}T_{i,j}^*+\kappa_6\frac{\tau_1^2}{2}(T_{i,j}S_{i,j}^*+T_{i,j}^*S_{i,j})\\ \qquad +\kappa_6 \frac{\tau_1^2\tau_2}{2}S_{ij}S_{ij}^*+(\kappa_4+\kappa_5)(\tau_1+\tau_2)(T_{i,i}S_{j,j}^*+T_{i,i}^*S_{j,j}) +(\kappa_4+\kappa_5)\frac{\tau_1^2\tau_2}{2}S_{i,i}S_{j,j}^*\\ \qquad +(\kappa_4+\kappa_5)(\tau_1+\tau_2)T_{i,i}T_{j,j}^* +(\kappa_4+\kappa_5)\frac{\tau_1^2}{2}(T_{i,i}S_{j,j}^*+T_{i,i}^*S_{j,j})\Big]\, d \mathit{\boldsymbol{x}}, \end{array}$

where $\hat{\theta} = \theta+\tau_1\vartheta+\frac{\tau_1^2}{2}\xi$ and $\hat{ \mathit{\boldsymbol{T}}} = \mathit{\boldsymbol{T}}+\tau_1 \mathit{\boldsymbol{S}}+\frac{\tau_1^2}{2} \mathit{\boldsymbol{R}}$ .

It is worth noting that this inner product defines a norm which is equivalent to the usual one in the Hilbert space because of the previous assumptions.

Our problem can be written as

$\frac{d U}{dt}(t) = \mathcal{A} U(t),\quad U(0) = U^0,$

where $U^0 = (\mathit{\boldsymbol{u}}^0, \mathit{\boldsymbol{v}}^0, \phi^0, \psi^0, \theta^0, \vartheta^0, \xi^0, \mathit{\boldsymbol{T}}^0, \mathit{\boldsymbol{S}}^0, \mathit{\boldsymbol{R}}^0)$ and

$\mathcal{A} = \left(\begin{array}{c c c c c c c c c c } 0 & \mathit{\boldsymbol{I}} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ \mathit{\boldsymbol{A}} & 0 & \mathit{\boldsymbol{B}} & 0 & \mathit{\boldsymbol{C}} & \mathit{\boldsymbol{D}} & \mathit{\boldsymbol{E}} & 0 & 0 & 0\\ 0 & 0 & 0 & I & 0 &0 & 0 & 0 & 0 & 0 \\ G & 0 & F & 0 & H & J & K & L & M & N\\ 0 & 0 & 0 & 0 & 0 & I & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & I & 0 & 0 & 0 \\ 0 & O & 0 & P & Q & U & V & W & X & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & I & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & I\\ 0 & 0 & 0 & \mathit{\boldsymbol{Y}} & \mathit{\boldsymbol{Z}} & \mathit{\boldsymbol{A}}^* & 0 & \mathit{\boldsymbol{B}}^* & \mathit{\boldsymbol{C}}^* & \mathit{\boldsymbol{D}}^* \end{array}\right),$

where $I$ is the identity operator and

$\begin{array}{l} A_i \mathit{\boldsymbol{u}} = \frac{1}{\rho}\Big[\mu u_{i,jj}+(\lambda+\mu) u_{j,ji}\Big],\quad \mathit{\boldsymbol{A}} = (A_i),\\ B_i\phi = \frac{1}{\rho}\mu_0\phi_{,i},\quad \mathit{\boldsymbol{B}} = (B_i),\quad C_i\theta = -\frac{1}{\rho} \beta_0\theta_{,i},\quad \mathit{\boldsymbol{C}} = (C_i),\\ D_i \vartheta = -\frac{\tau_1}{\rho}\beta_0\vartheta_{,i},\quad \mathit{\boldsymbol{D}} = (D_i),\quad E_i \xi = -\frac{\tau_1^2}{2}\beta_0\xi_{,i},\quad \mathit{\boldsymbol{E}} = (E_i),\\ F\phi = \frac{1}{J}[a_0\phi_{,ii}-\zeta\phi],\quad G \mathit{\boldsymbol{u}} = \frac{-1}{J} \mu_0 u_{i,i},\\ H\theta = \frac{\beta_1}{J}\theta, \quad J \vartheta = \frac{\beta_1 \tau_1}{J} \vartheta,\quad K\xi = \frac{\beta_1\tau_1^2}{2J}\xi,\\ L \mathit{\boldsymbol{T}} = -\frac{\mu_2}{J} T_{i,i},\quad M \mathit{\boldsymbol{S}} = \frac{-\mu_2\tau_1}{J}S_{i,i},\quad N \mathit{\boldsymbol{R}} = -\frac{\mu_2\tau_1^2}{2J}R_{i,i},\\ O \mathit{\boldsymbol{v}} = -\frac{2\beta_0}{a\tau_1^2}v_{i,i},\quad P\psi = -\frac{2\beta_1}{a\tau_1^2}\psi,\\ Q\theta = -\frac{2\kappa}{a\tau_1^2}\theta_{,ii},\quad U\vartheta = -\frac{2}{a\tau_1^2}[\kappa\tau_2 \vartheta_{,ii}-\frac{\tau_1^2}{2}\vartheta],\\ V\xi = -\frac{2}{\tau_1}\xi,\quad W \mathit{\boldsymbol{T}} = \frac{2}{a\tau_1^2}\kappa_1 T_{i,i},\\ X \mathit{\boldsymbol{S}} = \frac{2\tau_2}{a\tau_1^2}S_{i,i},\quad Y_i\psi = -\frac{2\mu_2}{b\tau_1^2}\psi_{,i}, \quad \mathit{\boldsymbol{Y}} = (Y_i),\\ Z_i\theta = -\frac{2\kappa_1}{b\tau_1^2}\theta_{,i}, \quad \mathit{\boldsymbol{Z}} = (Z_i),\quad A_i^*\vartheta = -\frac{2\kappa_1\tau_2}{b\tau_1^2}\vartheta_{,i},\quad \mathit{\boldsymbol{A}}^* = (A_i^*),\\ B_i^* \mathit{\boldsymbol{T}} = \frac{2}{b\tau_1^2}[\kappa_6 T_{i,jj}+ (\kappa_4+\kappa_5)T_{j,ji}], \quad \mathit{\boldsymbol{B}}^* = (B_i^*),\\ C_i^* \mathit{\boldsymbol{S}} = \frac{2}{b\tau_1^2}[\kappa_6\tau_2S_{i,jj}+(\kappa_4+\kappa_5)\tau_2T_{j,ji}-b S_i], \quad \mathit{\boldsymbol{C}}^* = (C_i^*),\\ D_i^* \mathit{\boldsymbol{R}} = -\frac{2}{\tau_1}R_i,\quad \mathit{\boldsymbol{D}}^* = (D_i^*). \end{array}$

We note that the domain of our operator is the subset $(\mathit{\boldsymbol{u}}, \mathit{\boldsymbol{v}}, \phi, \psi, \theta, \vartheta, \xi, \mathit{\boldsymbol{T}}, \mathit{\boldsymbol{S}}, \mathit{\boldsymbol{R}})\in \mathcal{H}$ such that

$\begin{array}{l} \mathit{\boldsymbol{A}} \mathit{\boldsymbol{u}}\in \mathit{\boldsymbol{L}}^2(B),\quad \xi \in W_0^{1,2}(B),\quad F\phi\in L^2(B),\quad \mathit{\boldsymbol{R}}\in \mathit{\boldsymbol{W}}_0^{1,2}(B),\quad \mathit{\boldsymbol{v}}\in \mathit{\boldsymbol{W}}_0^{1,2}(B),\\ P\theta+Q\vartheta\in L^2(B),\quad \mathit{\boldsymbol{B}}^* \mathit{\boldsymbol{T}}+ \mathit{\boldsymbol{C}}^* \mathit{\boldsymbol{S}}\in \mathit{\boldsymbol{L}}^2(B), \quad \psi \in W_0^{1,2}(B). \end{array}$

It is clear that this is a dense subspace of our Hilbert space. On the other side, we have that

$\begin{array}{l} Re < \mathcal{A} U,U > = -\frac{1}{2}\int_B \Big[ \kappa(\theta_{,i}\theta_{,i}+(\tau_1\tau_2-\frac{\tau_1^2}{2})\vartheta_{,i}\vartheta_{,i})\\ \qquad +\kappa_2(T_i T_i+(\tau_1\tau_2 - \frac{\tau_1^2}{2})S_iS_i)+2\kappa_1(\theta_{,i} T_i+(\tau_1\tau_2 - \frac{\tau_1^2}{2})\vartheta_{,i}S_i)\\ \qquad +\kappa_6(T_{i,j}T_{i,j} +(\tau_1\tau_2 - \frac{\tau_1^2}{2})S_{i,j}S_{i,j})\\ \qquad +(\kappa_4+\kappa_5)(T_{i,i}T_{j,j}+ (\tau_1\tau_2 - \frac{\tau_1^2}{2})S_{i,i}S_{j,j})\Big]\, d \mathit{\boldsymbol{x}}\leq 0 \end{array}$

for every $U$ in the domain of the operator denoted by $\mathcal{D(A)}$ .

To show that $\mathcal{A}$ generates a contractive semigroup, it is sufficient to prove that zero belongs to the resolvent of the operator.

To this end, we consider $F = (\mathit{\boldsymbol{f}}_1, \mathit{\boldsymbol{f}}_2, f_3, f_4, f_5, f_6, f_7, \mathit{\boldsymbol{f}}_8, \mathit{\boldsymbol{f}}_9, \mathit{\boldsymbol{f}}_{10})\in \mathcal{H}$ . We need to prove that the equation

$\mathcal{A} U = F$

has a solution in the domain of the operator. That is, we must solve the following system:

$\begin{array}{l} \mathit{\boldsymbol{v}} = \mathit{\boldsymbol{f}}_1,\\ \mathit{\boldsymbol{A}} \mathit{\boldsymbol{u}}+ \mathit{\boldsymbol{B}}\phi+ \mathit{\boldsymbol{C}}\theta+ \mathit{\boldsymbol{D}}\vartheta+ \mathit{\boldsymbol{E}}\xi = \mathit{\boldsymbol{f}}_2,\\ \psi = f_3,\\ G \mathit{\boldsymbol{u}}+F\phi +H\theta+J\vartheta+K\xi+L \mathit{\boldsymbol{T}}+M \mathit{\boldsymbol{S}}+N \mathit{\boldsymbol{R}} = f_4,\\ \vartheta = f_5,\\ \xi = f_6,\\ O \mathit{\boldsymbol{v}}+P\psi+Q\theta+U\vartheta+V\xi+W \mathit{\boldsymbol{T}}+X \mathit{\boldsymbol{S}} = f_7,\\ \mathit{\boldsymbol{S}} = \mathit{\boldsymbol{f}}_8,\\ \mathit{\boldsymbol{R}} = \mathit{\boldsymbol{f}}_9,\\ \mathit{\boldsymbol{Y}}\psi+ \mathit{\boldsymbol{Z}} \theta+ \mathit{\boldsymbol{A}}^* \vartheta+ \mathit{\boldsymbol{B}}^* \mathit{\boldsymbol{T}}+ \mathit{\boldsymbol{C}}^* \mathit{\boldsymbol{S}}+ \mathit{\boldsymbol{D}}^* \mathit{\boldsymbol{R}} = \mathit{\boldsymbol{f}}_{10}. \end{array}$

We obtain the solution for $\mathit{\boldsymbol{v}}, \, \psi, \mathit{\boldsymbol{R}}$ and $\xi$ satisfying the regularity conditions and we also find that

$\begin{array}{l} \mathit{\boldsymbol{A}} \mathit{\boldsymbol{u}}+ \mathit{\boldsymbol{B}}\phi+ \mathit{\boldsymbol{C}}\theta = \mathit{\boldsymbol{f}}_2- \mathit{\boldsymbol{D}} f_5- \mathit{\boldsymbol{E}} f_6,\\ G \mathit{\boldsymbol{u}}+F\phi +H\theta+L \mathit{\boldsymbol{T}} = f_4-Jf_5-Kf_6-M \mathit{\boldsymbol{f}}_8-N \mathit{\boldsymbol{f}}_9,\\ Q\theta+W \mathit{\boldsymbol{T}} = f_7-O \mathit{\boldsymbol{f}}_1-P f_3-Uf_5-V f_6-X \mathit{\boldsymbol{f}}_8,\\ \mathit{\boldsymbol{Z}}\theta+ \mathit{\boldsymbol{B}}^* \mathit{\boldsymbol{T}} = \mathit{\boldsymbol{f}}_{10}- \mathit{\boldsymbol{Y}} f_3- \mathit{\boldsymbol{A}}^* f_5- \mathit{\boldsymbol{C}}^* \mathit{\boldsymbol{f}}_8- \mathit{\boldsymbol{D}}^* \mathit{\boldsymbol{f}}_9. \end{array}$

If we consider the last two equations we see that the right-hand side belongs to $W^{-1, 2}(B)\times \mathit{\boldsymbol{W}}^{-1, 2}(B)$ . On the other side, the bilinear form

$\mathcal{B}[(\theta, \mathit{\boldsymbol{T}}),(\theta^*, \mathit{\boldsymbol{T}}^*)] = < (\theta, \mathit{\boldsymbol{T}}),(Q\theta^*+W \mathit{\boldsymbol{T}}^*, \mathit{\boldsymbol{Z}}\theta^*+ \mathit{\boldsymbol{B}}^* \mathit{\boldsymbol{T}}^*) >$

defines a bounded, and coercive, bilinear form in $W_0^{1, 2}(B)\times \mathit{\boldsymbol{W}}_0^{1, 2}(B)$ . In view of the Lax-Milgram lemma, we conclude the existence of $(\theta, \mathit{\boldsymbol{T}})\in W_0^{1, 2}(B)\times \mathit{\boldsymbol{W}}_0^{1, 2}(B)$ satisfying the last two equations. Now, if we introduce this solution into the first two ones we obtain

$\begin{array}{l} \mathit{\boldsymbol{A}} \mathit{\boldsymbol{u}}+ \mathit{\boldsymbol{B}}\phi = \mathit{\boldsymbol{f}}_2- \mathit{\boldsymbol{D}} f_5- \mathit{\boldsymbol{E}} f_6- \mathit{\boldsymbol{C}}\theta,\\ G \mathit{\boldsymbol{u}}+F\phi = f_4-Jf_5-Kf_6-M \mathit{\boldsymbol{f}}_8-N \mathit{\boldsymbol{f}}_9-H\theta. \end{array}$

A similar argument shows the existence of $(\mathit{\boldsymbol{u}}, \phi)\in \mathit{\boldsymbol{W}}_0^{1, 2}(B)\times W_0^{1, 2}(B)$ .

In fact, we can also prove an inequality of the type

$\|U\|\leq C^* \|F\|$

for a given positive constant $C^* > 0$ , and so, we have proved the following.

Theorem 3.1. The operator $\mathcal{A}$ generates a contractive semigroup.

Therefore, we can conclude the existence and uniqueness of solution to problem (2.1)-(2.6), that we state as follows.

Theorem 3.2. Assume that constitutive coefficients satisfy conditions(2.7) and that $U^0\in \mathcal{D(A)}$ . Then, there exists aunique solution $U\in C^1([0, \infty); \mathcal{H})\cap C([0, \infty); \mathcal{D(A)})$ to system (2.1)-(2.4), boundary conditions(2.5) and initial conditions (2.6).

4. Fully discrete approximations: a priori error analysis

In this section, we consider a fully discrete aproximation to a variational problem of the above thermomechanical problem by using the finite element method and the implicit Euler scheme.

Let $Y = L^2(B)$ , $H = [L^2(B)]^d$ and $Q = [L^2(B)]^{d\times d}$ , and denote by $(\cdot, \cdot)_Y$ , $(\cdot, \cdot)_H$ and $(\cdot, \cdot)_Q$ the respective scalar products in these spaces, with corresponding norms $\|\cdot\|_Y$ , $\|\cdot\|_H$ and $\|\cdot\|_Q$ . Moreover, let the variational spaces $E$ and $V$ be given as

$\begin{array}{l} E = \{z \in H^1(B)\, ; \, z = 0 \quad \hbox{on} \quad \Gamma\},\\ V = \{ \mathit{\boldsymbol{w}} \in [H^1(B)]^d\, ; \, \mathit{\boldsymbol{w}} = \mathit{\boldsymbol{0}} \quad \hbox{on} \quad \Gamma\}, \end{array}$

with respective scalar products $(\cdot, \cdot)_E$ and $(\cdot, \cdot)_V$ , and norms $\|\cdot\|_E$ and $\|\cdot\|_V$ .

By using Green's formula and taking into account the above boundary conditions, we write the variational formulation of the thermomechanical problem in terms of variables $\mathit{\boldsymbol{v}} = (v_i) = \dot{ \mathit{\boldsymbol{u}}} = (\dot{u}_i)$ and $\psi = \dot{\phi}$ , the temperature speed $\vartheta = \dot{\theta}$ , the temperature acceleration $\xi = \dot{\vartheta}$ , the microtemperature speed $\mathit{\boldsymbol{S}} = (S_i) = \dot{ \mathit{\boldsymbol{T}}} = (\dot{T}_i)$ and the microtemperature acceleration $\mathit{\boldsymbol{R}} = (R_i) = \dot{ \mathit{\boldsymbol{S}}} = (\dot{S}_i)$ .

Problem VP. Find the function $\mathit{\boldsymbol{v}}:[0, T_f]\rightarrow V$ , the function $\psi:[0, T_f]\rightarrow E$ , the temperature acceleration $\xi:[0, T_f]\rightarrow E$ and the microtemperatures acceleration $\mathit{\boldsymbol{R}}:[0, T_f]\rightarrow V$ such that $\mathit{\boldsymbol{v}}(0) = \mathit{\boldsymbol{v}}^0$ , $\psi(0) = \psi^0$ , $\xi(0) = \xi^0$ and $\mathit{\boldsymbol{R}} (0) = \mathit{\boldsymbol{R}}^0$ , and, for a.e. $t\in (0, T_f)$ and $\mathit{\boldsymbol{w}}, \, \mathit{\boldsymbol{\eta}}\in V,$ $r, \, z\in E$ ,

$\begin{eqnarray} &&\rho (\dot{ \mathit{\boldsymbol{v}}}(t) , \mathit{\boldsymbol{w}})_H +\mu (\nabla \mathit{\boldsymbol{u}} (t), \nabla \mathit{\boldsymbol{w}})_Q+(\lambda+\mu) ( \mbox{div}\, \mathit{\boldsymbol{u}}(t), \mbox{div}\, \mathit{\boldsymbol{w}})_Y \\ && \qquad -\mu_0(\nabla \phi(t), \mathit{\boldsymbol{w}})_H+\beta_0(\nabla\Big(\frac{\tau_1^2}{2}\xi(t)+\tau_1\vartheta(t)+\theta(t)\Big), \mathit{\boldsymbol{w}})_H = 0 , \end{eqnarray}$

(4.1)

$\begin{eqnarray} && J(\dot{\psi}(t),r)_Y+a_0(\nabla \phi(t),\nabla r)_H+\zeta (\phi(t),r)_Y = -\mu_0( \mbox{div}\, \mathit{\boldsymbol{u}}(t),r)_Y\\ && \qquad -\mu_2\Big(\frac{\tau_1^2}{2} \mbox{div}\, \mathit{\boldsymbol{R}}(t)+\tau_1 \mbox{div}\, \mathit{\boldsymbol{S}}(t)+ \mbox{div}\, \mathit{\boldsymbol{T}}(t),r\Big)_Y\\ && \qquad +\beta_1\Big(\frac{\tau_1^2}{2}\xi(t)+\tau_1\vartheta(t)+\theta(t),r\Big)_Y, \end{eqnarray}$

(4.2)

$\begin{eqnarray} && a\Big(\frac{\tau_1^2}{2}\dot{\xi}(t)+\tau_1\xi(t)+\vartheta(t),z\Big)_Y+\kappa(\nabla (\theta(t)+\tau_2 \vartheta(t)),\nabla z)_H = -\beta_0( \mbox{div}\, \mathit{\boldsymbol{v}}(t),z)_Y\\ && \qquad -\beta_1 (\psi(t),z)_Y +\kappa_1 ( \mbox{div}\, \mathit{\boldsymbol{T}}(t)+\tau_2 \mbox{div}\, \mathit{\boldsymbol{S}}(t),z)_Y, \end{eqnarray}$

(4.3)

$\begin{eqnarray} && b\Big(\frac{\tau_1^2}{2} \dot{ \mathit{\boldsymbol{R}}}(t)+\tau_1 \mathit{\boldsymbol{R}}(t)+ \mathit{\boldsymbol{S}}(t), \mathit{\boldsymbol{\eta}}\Big)_H+\kappa_6(\nabla ( \mathit{\boldsymbol{T}}(t)+\tau_2 \mathit{\boldsymbol{S}}(t)),\nabla \mathit{\boldsymbol{\eta}})_Q\\ && \qquad +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}(t)+\tau_2 \mathit{\boldsymbol{S}}(t)), \mbox{div}\, \mathit{\boldsymbol{\eta}})_Y +\kappa_2( \mathit{\boldsymbol{T}}(t)+\tau_2 \mathit{\boldsymbol{S}}(t), \mathit{\boldsymbol{\eta}})_H\\ && \quad = -\kappa_1(\nabla \theta(t)+\tau_1\vartheta(t), \mathit{\boldsymbol{\eta}})_H-\mu_2(\nabla \psi(t), \mathit{\boldsymbol{\eta}})_H , \end{eqnarray}$

(4.4)

where functions $\mathit{\boldsymbol{u}}$ , $\phi$ , $\vartheta$ , $\theta$ , $\mathit{\boldsymbol{S}}$ and $\mathit{\boldsymbol{T}}$ are then recovered from the relations:

$\begin{equation} \begin{array}{l} \mathit{\boldsymbol{u}} (t) = \int_0^t \mathit{\boldsymbol{v}}(s)\, ds+ \mathit{\boldsymbol{u}}^{0},\quad \phi(t) = \int_0^t \psi(s)\, ds+\phi^{0},\\ \vartheta (t) = \int_0^t\xi(s)\, ds+\vartheta^{0},\quad \mathit{\boldsymbol{S}}(t) = \int_0^t \mathit{\boldsymbol{R}}(s)\, ds+ \mathit{\boldsymbol{S}}^{0},\\ \theta (t) = \int_0^t\vartheta(s)\, ds+\theta^{0},\quad \mathit{\boldsymbol{T}}(t) = \int_0^t \mathit{\boldsymbol{S}}(s)\, ds+ \mathit{\boldsymbol{T}}^{0}, \end{array} \end{equation}$

(4.5)

and $T_f > 0$ denotes the final time of interest.

Thus, we construct the finite element spaces $V^h\subset V$ and $E^h\subset E$ given by

$\begin{eqnarray} &&V^h = \{ \mathit{\boldsymbol{w}}^h\in [C(\overline{B})]^d \; ; \; \mathit{\boldsymbol{w}}^h_{|_{Tr}}\in [P_1(Tr)]^d \quad \forall Tr\in { \mathcal{T}}^h,\quad \mathit{\boldsymbol{w}}^h = \mathit{\boldsymbol{0}}\, \hbox{ on } \,\Gamma\}, \end{eqnarray}$

(4.6)

$\begin{eqnarray} [2pt] &&E^h = \{r^h\in C(\overline{B}) \; ; \; r^h_{|_{Tr}}\in P_1(Tr) \quad \forall Tr\in { \mathcal{T}}^h,\quad r^h = 0\, \hbox{ on } \,\Gamma\}, \end{eqnarray}$

(4.7)

where $B$ is assumed to be a polyhedral domain, ${ \mathcal{T}}^h$ denotes a regular triangulation, in the sense of ^[29], of $\overline{B}$ , and $P_1(Tr)$ represents the space of polynomials of global degree less or equal to $1$ in element $Tr$ . Here, $h > 0$ denotes the spatial discretization parameter.

In order to discretize the time derivatives, we use a uniform partition of the time interval $[0, T_f]$ denoted by $0 = t_0 < t_1 < \ldots < t_N = T_f$ , with time step $k = T_f/N$ . Moreover, for a continuous function $f(t)$ let $f_n = f(t_n)$ and, for the sequence $\{z_n\}_{n = 0}^N$ , we denote by $\delta z_n = (z_n-z_{n-1})/k$ its corresponding divided differences.

Using the classical implicit Euler scheme, the fully discrete approximation of Problem VP is the following.

Problem VP $^{hk}$ . Find the discrete function $\mathit{\boldsymbol{v}}^{hk} = \{ \mathit{\boldsymbol{v}}_n^{hk}\}_{n = 0}^N \subset V^h$ , the discrete function $\psi^{hk} = \{\psi_n^{hk}\}_{n = 0}^N \subset E^h$ , the discrete temperature acceleration $\xi^{hk} = \{\xi_n^{hk}\}_{n = 0}^N \subset E^h$ and the discrete microtemperature acceleration $\mathit{\boldsymbol{R}}^{hk} = \{ \mathit{\boldsymbol{R}}_n^{hk}\}_{n = 0}^N \subset V^h$ such that $\mathit{\boldsymbol{v}}_0^{hk} = \mathit{\boldsymbol{v}}_0^h$ , $\psi_0^{hk} = \psi_0^h$ , $\xi_0^{hk} = \xi_0^h$ and $\mathit{\boldsymbol{R}}_0^{hk} = \mathit{\boldsymbol{R}}_0^h$ , and, for all $n = 1, \ldots, N$ and $\mathit{\boldsymbol{w}}^h, \, \mathit{\boldsymbol{\eta}}^h\in V^h,$ $r^h, \, z^h\in E^h$ ,

$\begin{eqnarray} &&\rho (\delta{ \mathit{\boldsymbol{v}}_n^{hk}} , \mathit{\boldsymbol{w}}^h)_H +\mu (\nabla \mathit{\boldsymbol{u}}_n^{hk}, \nabla \mathit{\boldsymbol{w}}^h)_Q+(\lambda+\mu) ( \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{w}}^h)_Y \\ && \qquad -\mu_0(\nabla \phi_n^{hk}, \mathit{\boldsymbol{w}}^h)_H+\beta_0(\nabla\Big(\frac{\tau_1^2}{2}\xi_n^{hk}+\tau_1\vartheta_n^{hk}+\theta_n^{hk}\Big), \mathit{\boldsymbol{w}}^h)_H = 0 , \end{eqnarray}$

(4.8)

$\begin{eqnarray} && J(\delta{\psi_n^{hk}},r^h)_Y+a_0(\nabla \phi_n^{hk},\nabla r^h)_H+\zeta (\phi_n^{hk},r^h)_Y = -\mu_0( \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk},r^h)_Y\\ && \qquad -\mu_2\Big(\frac{\tau_1^2}{2} \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk}+\tau_1 \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}+ \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk},r^h\Big)_Y\\ && \qquad +\beta_1\Big(\frac{\tau_1^2}{2}\xi_n^{hk}+\tau_1\vartheta_n^{hk}+\theta_n^{hk},r^h\Big)_Y, \end{eqnarray}$

(4.9)

$\begin{eqnarray} && a\Big(\frac{\tau_1^2}{2}\delta{\xi_n^{hk}}+\tau_1\xi_n^{hk}+\vartheta_n^{hk},z^h\Big)_Y+\kappa(\nabla (\theta_n^{hk}+\tau_2 \vartheta_n^{hk}),\nabla z^h)_H = -\beta_0( \mbox{div}\, \mathit{\boldsymbol{v}}_n^{hk},z^h)_Y\\ && \qquad -\beta_1 (\psi_n^{hk},z^h )_Y +\kappa_1 ( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk},z^h)_Y, \end{eqnarray}$

(4.10)

$\begin{eqnarray} && b\Big(\frac{\tau_1^2}{2} \delta{ \mathit{\boldsymbol{R}}_n^{hk}}+\tau_1 \mathit{\boldsymbol{R}}_n^{hk}+ \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{\eta}}^h\Big)_H+\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}),\nabla \mathit{\boldsymbol{\eta}}^h)_Q\\ && \qquad +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}), \mbox{div}\, \mathit{\boldsymbol{\eta}}^h)_Y +\kappa_2( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{\eta}}^h)_H\\ && \quad = -\kappa_1(\nabla( \theta_n^{hk}+\tau_1\vartheta_n^{hk}), \mathit{\boldsymbol{\eta}}^h)_H-\mu_2(\nabla \psi_n^{hk}, \mathit{\boldsymbol{\eta}}^h)_H , \end{eqnarray}$

(4.11)

where discrete functions $\mathit{\boldsymbol{u}}_n^{hk}$ , $\phi_n^{hk}$ , $\vartheta_n^{hk}$ $\theta_n^{hk}$ , $\mathit{\boldsymbol{S}}_n^{hk}$ and $\mathit{\boldsymbol{T}}_n^{hk}$ are then recovered from the relations:

$\begin{equation} \begin{array}{l} \mathit{\boldsymbol{u}}^{hk}_n = k\sum\limits_{j = 1}^n \mathit{\boldsymbol{v}}^{hk}_j+ \mathit{\boldsymbol{u}}_{0}^h,\quad \phi_n^{hk} = k\sum\limits_{j = 1}^n\psi_j^{hk}+\phi_{0}^h,\\ \vartheta^{hk}_n = k\sum\limits_{j = 1}^n\theta^{hk}_j+\vartheta_{0}^h,\quad \mathit{\boldsymbol{S}}_n^{hk} = k\sum\limits_{j = 1}^n \mathit{\boldsymbol{R}}_j^{hk}+ \mathit{\boldsymbol{S}}_{0}^h,\\ \theta^{hk}_n = k\sum\limits_{j = 1}^n\vartheta^{hk}_j+\theta_{0}^h,\quad \mathit{\boldsymbol{T}}_n^{hk} = k\sum\limits_{j = 1}^n \mathit{\boldsymbol{S}}_j^{hk}+ \mathit{\boldsymbol{T}}_{0}^h, \end{array} \end{equation}$

(4.12)

and the discrete initial conditions $\mathit{\boldsymbol{u}}_0^h$ , $\mathit{\boldsymbol{v}}_0^h$ , $\phi_0^h$ , $\psi_0^h$ , $\theta_0^h$ , $\vartheta_0^h$ , $\xi_0^h$ , $\mathit{\boldsymbol{T}}_0^h$ , $\mathit{\boldsymbol{S}}_0^h$ and $\mathit{\boldsymbol{R}}_0^h$ are defined as follows:

$\begin{equation} \begin{array}{l} \mathit{\boldsymbol{u}}_0^h = \mathcal{P}^{1h} \mathit{\boldsymbol{u}}^0,\quad \mathit{\boldsymbol{v}}_0^h = \mathcal{P}^{1h} \mathit{\boldsymbol{v}}^0,\quad \phi_0^h = \mathcal{P}^{2h} \phi^0,\quad \psi_0^h = \mathcal{P}^{2h} \psi^0,\\ \theta_0^h = \mathcal{P}^{2h} \theta^0,\quad \vartheta_0^h = \mathcal{P}^{2h} \vartheta^0,\quad \xi_0^h = \mathcal{P}^{2h} \xi^0,\quad \mathit{\boldsymbol{T}}_0^h = \mathcal{P}^{1h} \mathit{\boldsymbol{T}}^0,\\ \mathit{\boldsymbol{S}}_0^h = \mathcal{P}^{1h} \mathit{\boldsymbol{S}}^0,\quad \mathit{\boldsymbol{R}}_0^h = \mathcal{P}^{1h} \mathit{\boldsymbol{R}}^0. \end{array} \end{equation}$

(4.13)

In the previous definitions, operators $\mathcal{P}^{1h}$ and $\mathcal{P}^{2h}$ are the projection operators over the finite element spaces $V^h$ and $E^h$ , respectively (see, for instance, ^[30]).

Using the previously given assumptions (2.7), applying the well-known Lax-Milgram lemma we can easily show that fully discrete problem $VP^{hk}$ has a unique solution.

In this section, the objective is to prove a main error estimates result regarding the approximation of the solution to Problem VP by the solution to Problem $VP^{hk}$ . So, first we have the following discrete stability property.

Lemma 4.1. Under the assumptions of Theorem 3.2, it follows that thesequences

$\{ \mathit{\boldsymbol{u}}^{hk}, \mathit{\boldsymbol{v}}^{hk},\phi^{hk},\psi^{hk},\theta^{hk},\vartheta^{hk},\xi^{hk}, \mathit{\boldsymbol{T}}^{hk}, \mathit{\boldsymbol{S}}^{hk}, \mathit{\boldsymbol{R}}^{hk}\},$

generated by Problem $VP^{hk}$ , satisfy the stability estimate:

$\begin{array}{l} \| \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\| \mathit{\boldsymbol{u}}_n^{hk}\|^2_V+\|\psi_n^{hk}\|^2_Y+\| \phi_n^{hk}\|^2_V+\|\theta_n^{hk}\|_E^2+\| \vartheta_n^{hk}\|^2_E +\|\xi_n^{hk}\|^2_Y \\ \qquad +\| \mathit{\boldsymbol{R}}_n^{hk}\|^2_H+\| \mathit{\boldsymbol{S}}_n^{hk}\|^2_V+\| \mathit{\boldsymbol{T}}_n^{hk}\|_V^2\leq C,\end{array}$

where $C$ is a positiveconstant assumed to be independent of the discretization parameters $h$ and $k$ .

Proof. For the sake of simplicity, in this proof we need to assume that $\tau_1^2/2 = 1$ . We note that we can extend the analysis provided below to the general case doing some simple modifications.

First, we estimate the terms on the discrete function $\mathit{\boldsymbol{v}}_n^{hk}$ . Taking $\mathit{\boldsymbol{w}}^h = \mathit{\boldsymbol{v}}_n^{hk}$ as a test function in discrete variational equation (4.8) we obtain

$\begin{array}{l}\rho (\delta{ \mathit{\boldsymbol{v}}_n^{hk}} , \mathit{\boldsymbol{v}}_n^{hk})_H +\mu (\nabla \mathit{\boldsymbol{u}}_n^{hk}, \nabla \mathit{\boldsymbol{v}}_n^{hk})_Q+(\lambda+\mu) ( \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{v}}_n^{hk})_Y \\ \qquad -\mu_0(\nabla \phi_n^{hk}, \mathit{\boldsymbol{v}}_n^{hk})_H+\beta_0(\nabla(\xi_n^{hk}+\tau_1\vartheta_n^{hk}+\theta_n^{hk}), \mathit{\boldsymbol{v}}_n^{hk})_H = 0. \end{array}$

From the estimates

$\begin{array}{l} (\delta \mathit{\boldsymbol{v}}_n^{hk}, \mathit{\boldsymbol{v}}^{hk}_n)_H\geq \frac{1}{2k}\left\{\| \mathit{\boldsymbol{v}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{v}}_{n-1}^{hk}\|^2_H\right\},\\ (\nabla \mathit{\boldsymbol{u}}_n^{hk},\nabla \mathit{\boldsymbol{v}}_n^{hk})_Q\geq \frac{1}{2k}\left\{\|\nabla \mathit{\boldsymbol{u}}_n^{hk}\|^2_Q-\|\nabla \mathit{\boldsymbol{u}}_{n-1}^{hk}\|^2_Q\right\},\\ ( \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{v}}_n^{hk})_Y\geq \frac{1}{2k}\left\{\| \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}\|^2_Y-\| \mbox{div}\, \mathit{\boldsymbol{u}}_{n-1}^{hk}\|^2_Y\right\}, \end{array}$

applying the Cauchy's inequality:

$\begin{equation} ab\leq \epsilon a^2+\frac{1}{4\epsilon} b^2,\quad a,b,\epsilon\in \mathbb{R},\,\epsilon > 0, \end{equation}$

(4.14)

it follows that

Now, we estimate the terms on the discrete function $\psi_n^{hk}$ . Taking $r^h = \psi_n^{hk}$ as a test function in discrete variational equation (4.9) we have

$\begin{array}{l} J(\delta{\psi_n^{hk}},\psi_n^{hk})_Y+a_0(\nabla \phi_n^{hk},\nabla \psi_n^{hk})_H+\zeta (\phi_n^{hk},\psi_n^{hk})_Y = -\mu_0( \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk},\psi_n^{hk})_Y\\ \quad -\mu_2( \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk}+\tau_1 \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}+ \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk},\psi_n^{hk})_Y +\beta_1(\xi_n^{hk}+\tau_1\vartheta_n^{hk}+\theta_n^{hk},\psi_n^{hk})_Y, \end{array}$

and, using the estimates

$\begin{array}{l} (\delta \psi_n^{hk},\psi^{hk}_n)_Y\geq \frac{1}{2k}\left\{\|\psi_n^{hk}\|^2_Y-\|\psi_{n-1}^{hk}\|^2_Y\right\},\\ (\nabla \phi_n^{hk},\nabla\psi_n^{hk})_H\geq \frac{1}{2k}\left\{\|\nabla \phi_n^{hk}\|^2_H-\|\nabla \phi_{n-1}^{hk}\|^2_Y\right\},\\ ( \phi_n^{hk},\psi_n^{hk})_Y\geq \frac{1}{2k}\left\{\| \phi_n^{hk}\|^2_Y-\| \phi_{n-1}^{hk}\|^2_Y\right\}, \end{array}$

we find that

$\begin{array}{l} \frac{1}{2k}\left\{\|\psi_n^{hk}\|^2_Y-\|\psi_{n-1}^{hk}\|^2_Y\right\}+\frac{1}{2k}\left\{\|\nabla \phi_n^{hk}\|^2_H-\|\nabla \phi_{n-1}^{hk}\|^2_Y\right\}\\ \qquad +\frac{1}{2k}\left\{\| \phi_n^{hk}\|^2_Y-\| \phi_{n-1}^{hk}\|^2_Y\right\}+\mu_2( \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk},\psi_n^{hk})_Y\\ \quad \leq C\Big(\|\psi_n^{hk}\|^2_Y+\| \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}\|_Y^2+\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|_Y^2+\| \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}\|_Y^2 +\|\vartheta_n^{hk}\|_Y^2+\|\theta_n^{hk}\|_Y^2+\|\xi_n^{hk}\|_Y^2\Big). \end{array}$

Thirdly, we get the estimates on the discrete temperature acceleration $\xi_n^{hk}$ . Thus, taking $z^h = \xi_n^{hk}$ as a test function in discrete variational equation (4.10) we obtain

$\begin{array}{l} a(\delta{\xi_n^{hk}}+\tau_1\xi_n^{hk}+\vartheta_n^{hk},\xi_n^{hk})_Y+\kappa(\nabla (\theta_n^{hk}+\tau_2 \vartheta_n^{hk}),\nabla \xi_n^{hk})_H = -\beta_0( \mbox{div}\, \mathit{\boldsymbol{v}}_n^{hk},\xi_n^{hk})_Y\\ \qquad -\beta_1 (\psi_n^{hk},\xi_n^{hk} )_Y +\kappa_1 ( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk},\xi_n^{hk})_Y, \end{array}$

and therefore, taking into account that

$\begin{array}{l} (\delta \xi_n^{hk},\xi^{hk}_n)_Y\geq \frac{1}{2k}\left\{\|\xi_n^{hk}\|^2_Y-\|\xi_{n-1}^{hk}\|^2_Y\right\},\\ (\nabla \vartheta_n^{hk},\nabla\xi_n^{hk})_H\geq \frac{1}{2k}\left\{\|\nabla \vartheta_n^{hk}\|^2_H-\|\nabla \vartheta_{n-1}^{hk}\|^2_H\right\},\\ ( \vartheta_n^{hk},\xi_n^{hk})_Y\geq \frac{1}{2k}\left\{\| \vartheta_n^{hk}\|^2_Y-\| \vartheta_{n-1}^{hk}\|^2_Y\right\}, \end{array}$

we find that

$\begin{array}{l} \frac{a}{2k}\left\{\|\xi_n^{hk}\|^2_Y-\|\xi_{n-1}^{hk}\|^2_Y\right\}+\frac{\tau_2}{2k}\left\{\|\nabla \vartheta_n^{hk}\|^2_H-\|\nabla \vartheta_{n-1}^{hk}\|^2_H\right\}+\beta_0( \mbox{div}\, \mathit{\boldsymbol{v}}_n^{hk},\xi_n^{hk})_Y\\ \qquad + \frac{a}{2k}\left\{\| \vartheta_n^{hk}\|^2_Y-\| \vartheta_{n-1}^{hk}\|^2_Y\right\}+\kappa (\nabla \theta_n^{hk},\nabla\xi_n^{hk})_H\\ \quad \leq C\Big(\|\xi_n^{hk}\|^2_Y+\|\nabla\theta_n^{hk}\|_H^2+\|\psi_n^{hk}\|^2_Y+\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|_Y^2+\| \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}\|_Y^2 \Big). \end{array}$

Finally, we get the estimates on the discrete microtemperatures acceleration $\mathit{\boldsymbol{R}}_n^{hk}$ . Then, taking $\mathit{\boldsymbol{\eta}}^h = \mathit{\boldsymbol{R}}_n^{hk}$ as a test function in discrete variational equation (4.11) it follows that

$\begin{array}{l} b(\delta{ \mathit{\boldsymbol{R}}_n^{hk}}+\tau_1 \mathit{\boldsymbol{R}}_n^{hk}+ \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{R}}_n^{hk})_H+\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}),\nabla \mathit{\boldsymbol{R}}_n^{hk})_Q\\ \nonumber \qquad +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}), \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y +\kappa_2( \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{R}}_n^{hk})_H\\ \quad = -\kappa_1(\nabla( \theta_n^{hk}+\tau_1\vartheta_n^{hk}), \mathit{\boldsymbol{R}}_n^{hk})_H-\mu_2(\nabla \psi_n^{hk}, \mathit{\boldsymbol{R}}_n^{hk})_H, \end{array}$

and so, keeping in mind that

$\begin{array}{l} (\delta \mathit{\boldsymbol{R}}_n^{hk}, \mathit{\boldsymbol{R}}^{hk}_n)_H\geq \frac{1}{2k}\left\{\| \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}^{hk}\|^2_H\right\},\\ (\nabla \mathit{\boldsymbol{S}}_n^{hk},\nabla \mathit{\boldsymbol{R}}_n^{hk})_Q\geq \frac{1}{2k}\left\{\|\nabla \mathit{\boldsymbol{S}}_n^{hk}\|^2_Q-\|\nabla \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Q\right\},\\ ( \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y\geq \frac{1}{2k}\left\{\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|^2_Y-\| \mbox{div}\, \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Y\right\},\\ ( \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{R}}^{hk}_n)_H\geq \frac{1}{2k}\left\{\| \mathit{\boldsymbol{S}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_H\right\},\\ (\nabla \psi_n^{hk}, \mathit{\boldsymbol{R}}_n^{hk})_H = -(\psi_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y, \end{array}$

we obtain

$\begin{array}{l} \frac{b}{2k}\left\{\| \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}^{hk}\|^2_H\right\}+ \frac{\kappa_6}{2k}\left\{\|\nabla \mathit{\boldsymbol{S}}_n^{hk}\|^2_Q-\|\nabla \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Q\right\}-\mu_2( \psi_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y\\ \qquad +\frac{\kappa_4+\kappa_5}{2k}\left\{\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|^2_Y-\| \mbox{div}\, \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Y\right\}+\frac{\kappa_2\tau_2}{2k}\left\{\| \mathit{\boldsymbol{S}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_H\right\}\\ \qquad + (\kappa_4+\kappa_5)( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y+\kappa_6(\nabla \mathit{\boldsymbol{T}}_n^{hk},\nabla \mathit{\boldsymbol{R}}_n^{hk})_Q\\ \quad \leq C\Big(\| \mathit{\boldsymbol{R}}_n^{hk}\|_H^2+\|\nabla\theta_n^{hk}\|_H^2 +\|\nabla\vartheta_n^{hk}\|_H^2+\| \mathit{\boldsymbol{T}}_n^{hk}\|_H^2\Big). \end{array}$

Combining all these estimates, it leads to

$\begin{array}{l} \frac{\rho}{2k}\left\{\| \mathit{\boldsymbol{v}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{v}}_{n-1}^{hk}\|^2_H\right\} + \frac{\mu}{2k}\left\{\|\nabla \mathit{\boldsymbol{u}}_n^{hk}\|^2_Q-\|\nabla \mathit{\boldsymbol{u}}_{n-1}^{hk}\|^2_Q\right\}\\ \qquad + \frac{\lambda+\mu}{2k}\left\{\| \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}\|^2_Y-\| \mbox{div}\, \mathit{\boldsymbol{u}}_{n-1}^{hk}\|^2_Y\right\}\\ \qquad + \frac{1}{2k}\left\{\|\psi_n^{hk}\|^2_Y-\|\psi_{n-1}^{hk}\|^2_Y\right\}+\frac{1}{2k}\left\{\|\nabla \phi_n^{hk}\|^2_H-\|\nabla \phi_{n-1}^{hk}\|^2_Y\right\}\\ \qquad +\frac{1}{2k}\left\{\| \phi_n^{hk}\|^2_Y-\| \phi_{n-1}^{hk}\|^2_Y\right\} + \frac{1}{2k}\left\{\| \vartheta_n^{hk}\|^2_Y-\| \vartheta_{n-1}^{hk}\|^2_Y\right\} \\ \qquad +\frac{1}{2k}\left\{\|\xi_n^{hk}\|^2_Y-\|\xi_{n-1}^{hk}\|^2_Y\right\}+\frac{1}{2k}\left\{\|\nabla \vartheta_n^{hk}\|^2_H-\|\nabla \vartheta_{n-1}^{hk}\|^2_H\right\}\\ \qquad +\kappa (\nabla \theta_n^{hk},\nabla\xi_n^{hk})_H+ (\kappa_4+\kappa_5)( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y+\kappa_6(\nabla \mathit{\boldsymbol{T}}_n^{hk},\nabla \mathit{\boldsymbol{R}}_n^{hk})_Q\\ \qquad + \frac{1}{2k}\left\{\| \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}^{hk}\|^2_H\right\}+ \frac{1}{2k}\left\{\|\nabla \mathit{\boldsymbol{S}}_n^{hk}\|^2_Q-\|\nabla \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Q\right\}\\ \qquad +\frac{1}{2k}\left\{\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|^2_Y-\| \mbox{div}\, \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_Y\right\}+\frac{1}{2k}\left\{\| \mathit{\boldsymbol{S}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{S}}_{n-1}^{hk}\|^2_H\right\}\\ \quad \leq C\Big(\| \mathit{\boldsymbol{v}}_n^{hk}\|_H^2+\|\nabla \phi_n^{hk}\|_H^2+\|\nabla \theta_n^{hk}\|_H^2+\|\nabla \vartheta_n^{hk}\|_H^2+\|\psi_n^{hk}\|^2_Y+\| \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}\|_Y^2\\ \qquad +\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|_Y^2+\| \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}\|_Y^2 +\|\vartheta_n^{hk}\|_Y^2+\|\theta_n^{hk}\|_Y^2+\|\xi_n^{hk}\|^2_Y+\| \mathit{\boldsymbol{R}}_n^{hk}\|_H^2\\ \qquad +\| \mathit{\boldsymbol{T}}_n^{hk}\|_H^2+\|\nabla \mathit{\boldsymbol{T}}_n^{hk}\|_Q^2\Big). \end{array}$

Multiplying the above estimates by $k$ and summing up to $n$ we get

$\begin{array}{l} \| \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\|\nabla \mathit{\boldsymbol{u}}_n^{hk}\|^2_Q+\| \mbox{div}\, \mathit{\boldsymbol{u}}_n^{hk}\|^2_Y+\|\psi_n^{hk}\|^2_Y+\|\nabla \phi_n^{hk}\|^2_H+\| \phi_n^{hk}\|^2_Y+\| \vartheta_n^{hk}\|^2_Y \\ \qquad +\|\xi_n^{hk}\|^2_Y+\|\nabla \vartheta_n^{hk}\|^2_H+\| \mathit{\boldsymbol{R}}_n^{hk}\|^2_H+\|\nabla \mathit{\boldsymbol{S}}_n^{hk}\|^2_Q+\| \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk}\|^2_Y+\| \mathit{\boldsymbol{S}}_n^{hk}\|^2_H\\ \qquad +k\sum\limits_{j = 1}^n\Big[\kappa (\nabla \theta_j^{hk},\nabla\xi_j^{hk})_H+ (\kappa_4+\kappa_5)( \mbox{div}\, \mathit{\boldsymbol{T}}_j^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_j^{hk})_Y+\kappa_6(\nabla \mathit{\boldsymbol{T}}_j^{hk},\nabla \mathit{\boldsymbol{R}}_j^{hk})_Q\Big] \end{array}$

$\begin{array}{l} \quad \leq Ck\sum\limits_{j = 1}^n\Big(\| \mathit{\boldsymbol{v}}_j^{hk}\|_H^2+\|\nabla \phi_j^{hk}\|_H^2+\|\nabla \theta_j^{hk}\|_H^2+\|\nabla \vartheta_j^{hk}\|_H^2+\|\psi_j^{hk}\|^2_Y+\| \mbox{div}\, \mathit{\boldsymbol{u}}_j^{hk}\|_Y^2 \\ \qquad +\| \mbox{div}\, \mathit{\boldsymbol{S}}_j^{hk}\|_Y^2+\| \mbox{div}\, \mathit{\boldsymbol{T}}_j^{hk}\|_Y^2 +\|\vartheta_j^{hk}\|_Y^2+\|\theta_j^{hk}\|_Y^2+\|\xi_j^{hk}\|^2_Y+\| \mathit{\boldsymbol{R}}_j^{hk}\|_H^2\\ \qquad +\| \mathit{\boldsymbol{T}}_j^{hk}\|_H^2+\|\nabla \mathit{\boldsymbol{T}}_j^{hk}\|_Q^2\Big)+C\Big( \| \mathit{\boldsymbol{v}}_0^{h}\|^2_H+\| \mathit{\boldsymbol{u}}_n^{hk}\|^2_V+\|\psi_0^{h}\|^2_Y+\|\phi_0^{h}\|^2_E+\| \vartheta_0^{h}\|^2_E\\ \qquad +\|\xi_0^{h}\|^2_Y+\| \mathit{\boldsymbol{R}}_0^{h}\|^2_H+\| \mathit{\boldsymbol{S}}_0^{h}\|^2_V\Big). \end{array}$

From the definition of $\theta_n^{hk}$ and $\mathit{\boldsymbol{T}}_n^{hk}$ we can easily show that

$\begin{array}{l} \|\theta_n^{hk}\|_Y^2\leq Ck\sum\limits_{j = 1}^n \|\vartheta_j^{hk}\|_Y^2 +\|\theta_0^h\|_Y^2,\\ \|\nabla\theta_n^{hk}\|_H^2\leq Ck\sum\limits_{j = 1}^n \|\nabla\vartheta_j^{hk}\|_H^2 +\|\nabla\theta_0^h\|_H^2,\\ \| \mathit{\boldsymbol{T}}_n^{hk}\|_H^2\leq Ck\sum\limits_{j = 1}^n \| \mathit{\boldsymbol{S}}_j^{hk}\|_H^2 +\| \mathit{\boldsymbol{T}}_0^h\|_H^2,\\ \|\nabla \mathit{\boldsymbol{T}}_n^{hk}\|_Q^2\leq Ck\sum\limits_{j = 1}^n \|\nabla \mathit{\boldsymbol{T}}_j^{hk}\|_Q^2 +\|\nabla \mathit{\boldsymbol{T}}_0^h\|_Q^2,\\ \| \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}\|_Y^2\leq Ck\sum\limits_{j = 1}^n \| \mbox{div}\, \mathit{\boldsymbol{T}}_j^{hk}\|_Y^2 +\| \mbox{div}\, \mathit{\boldsymbol{T}}_0^h\|_Y^2. \end{array}$

Observing that

$\begin{array}{l} k\sum\limits_{j = 1}^n(\nabla \theta_n^{hk},\nabla \xi_n^{hk})_H = \sum\limits_{j = 1}^n (\nabla \theta_j^{hk},\nabla (\vartheta_j^{hk}-\vartheta_{j-1}^{hk}))_H\\ \qquad = (\nabla \theta_n^{hk},\nabla \vartheta_n^{hk})_H+\sum\limits_{j = 1}^{n-1}(\nabla(\theta_j^{hk}-\theta_{j+1}^{hk}),\nabla \vartheta_j^{hk})_H+(\nabla \theta_1^{hk},\nabla\vartheta_0^h)_H,\\ \sum\limits_{j = 1}^{n-1}(\nabla(\theta_j^{hk}-\theta_{j+1}^{hk}),\nabla \vartheta_j^{hk})_H\leq Ck\sum\limits_{j = 1}^n\|\nabla \vartheta_j^{hk}\|_H^2+\frac{C}{k}\sum\limits_{j = 1}^{n-1}\|\nabla(\theta_j^{hk}-\theta_{j+1}^{hk})\|_H^2\\ \qquad \leq Ck\sum\limits_{j = 1}^n \|\nabla \vartheta_j^{hk}\|_H^2,\\ k\sum\limits_{j = 1}^n(\nabla \mathit{\boldsymbol{T}}_n^{hk},\nabla \mathit{\boldsymbol{R}}_n^{hk})_Q = \sum\limits_{j = 1}^n (\nabla \mathit{\boldsymbol{T}}_j^{hk},\nabla ( \mathit{\boldsymbol{S}}_j^{hk}- \mathit{\boldsymbol{S}}_{j-1}^{hk}))_Q\\ \qquad = (\nabla \mathit{\boldsymbol{T}}_n^{hk},\nabla \mathit{\boldsymbol{S}}_n^{hk})_Q+\sum\limits_{j = 1}^{n-1}(\nabla( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk}),\nabla \mathit{\boldsymbol{S}}_j^{hk})_Q+(\nabla \mathit{\boldsymbol{T}}_1^{hk},\nabla \mathit{\boldsymbol{S}}_0^h)_Q,\\ \sum\limits_{j = 1}^{n-1}(\nabla( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk}),\nabla \mathit{\boldsymbol{S}}_j^{hk})_Q\leq Ck\sum\limits_{j = 1}^n\|\nabla \mathit{\boldsymbol{S}}_j^{hk}\|_Q^2+\frac{C}{k}\sum\limits_{j = 1}^{n-1}\|\nabla( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk})\|_Q^2\\ \qquad \leq Ck\sum\limits_{j = 1}^n \|\nabla \mathit{\boldsymbol{S}}_j^{hk}\|_Q^2, \end{array}$

$\begin{array}{l} k\sum\limits_{j = 1}^n( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{R}}_n^{hk})_Y = \sum\limits_{j = 1}^n ( \mbox{div}\, \mathit{\boldsymbol{T}}_j^{hk}, \mbox{div}\, ( \mathit{\boldsymbol{S}}_j^{hk}- \mathit{\boldsymbol{S}}_{j-1}^{hk}))_Y\\ \qquad = ( \mbox{div}\, \mathit{\boldsymbol{T}}_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{S}}_n^{hk})_Y+\sum\limits_{j = 1}^{n-1}( \mbox{div}\,( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk}), \mbox{div}\, \mathit{\boldsymbol{S}}_j^{hk})_Y+( \mbox{div}\, \mathit{\boldsymbol{T}}_1^{hk}, \mbox{div}\, \mathit{\boldsymbol{S}}_0^h)_Y,\\ \sum\limits_{j = 1}^{n-1}( \mbox{div}\,( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk}), \mbox{div}\, \mathit{\boldsymbol{S}}_j^{hk})_Y\leq Ck\sum\limits_{j = 1}^n\| \mbox{div}\, \mathit{\boldsymbol{S}}_j^{hk}\|_Y^2+\frac{C}{k}\sum\limits_{j = 1}^{n-1}\| \mbox{div}\,( \mathit{\boldsymbol{T}}_j^{hk}- \mathit{\boldsymbol{T}}_{j+1}^{hk})\|_Y^2\\ \qquad \leq Ck\sum\limits_{j = 1}^n \| \mbox{div}\, \mathit{\boldsymbol{S}}_j^{hk}\|_Y^2, \end{array}$

and, using a discrete version of Gronwall's inequality (see ^[31]), we conclude the discrete stability property.

Now, we will derive a main error estimates result on the numerical errors $\mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}$ , $\psi_n-\psi_n^{hk}$ , $\xi_n-\xi_n^{hk}$ and $\mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}$ . We have the following.

Theorem 4.2. Under the assumptions of Lemma 4.1, if we denote by $(\mathit{\boldsymbol{v}}, \psi, \xi, \mathit{\boldsymbol{R}})$ the solution to Problem $VP$ and by $(\mathit{\boldsymbol{v}}^{hk}, \psi^{hk}, \xi^{hk}, \mathit{\boldsymbol{R}}^{hk})$ the solution to Problem $VP^{hk}$ , then we have the following a priori error estimates, forall $\{ \mathit{\boldsymbol{w}}^{h}_j\}_{j = 0}^N,$ $\{ \mathit{\boldsymbol{\eta}}^{h}_j\}_{j = 0}^N\subset V^h$ and $\{r^{h}_j\}_{j = 0}^N, \{z^{h}_j\}_{j = 0}^N\subset E^h$ ,

$\begin{eqnarray*} &&\nonumber\max\limits_{0\leq n\leq N}\Big\{ \| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\| \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}\|^2_V+\|\phi_n-\phi_n^{hk}\|^2_E +\|\psi_n-\psi_n^{hk}\|^2_Y+\|\xi_n-\xi_n^{hk}\|^2_Y\\ &&\nonumber +\|\vartheta_n-\vartheta_n^{hk}\|^2_E+\|\theta_n-\theta_n^{hk}\|^2_E +\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H+\| \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}\|^2_V+\| \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}\|^2_V\Big\}\\ &&\nonumber \leq Ck\sum\limits_{j = 1}^N\Big(\|\dot{ \mathit{\boldsymbol{v}}}_j-\delta \mathit{\boldsymbol{v}}_j\|^2_H+ \| \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^h\|_V^2+\|\dot{ \mathit{\boldsymbol{u}}}_j-\delta \mathit{\boldsymbol{u}}_j\|^2_V+ \|\dot{\psi}_j-\delta \psi_j\|^2_Y+\|\dot{\phi}_j-\delta \phi_j\|^2_E\\ && +\|\psi_j-r^h_j\|^2_E+ \|\dot{\xi}_j-\delta \xi_j\|^2_Y+\|\dot{\vartheta}_j-\delta \vartheta_j\|^2_E+\|\xi_j-z^h_j\|^2_E +\|\dot{ \mathit{\boldsymbol{R}}}_j-\delta \mathit{\boldsymbol{R}}_j\|^2_H\\ &&\nonumber + \| \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{\eta}}_j^h\|_V^2+\|\dot{ \mathit{\boldsymbol{S}}}_j-\delta \mathit{\boldsymbol{S}}_j\|^2_V +I_j^2+J_j^2\Big)\\ && +C\max\limits_{0\leq n\leq N}\Big\{\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}_n^h\|_H^2+\|\psi_n-r_n^h\|_Y^2+\|\xi_n-z_n^h\|_Y^2+\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}_n^h\|_H^2\Big\}\\ && +\frac{C}{k}\sum\limits_{j = 1}^{N-1}\Big[\| \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^{h}-( \mathit{\boldsymbol{v}}_{j+1}- \mathit{\boldsymbol{w}}_{j+1}^h)\|_H^2+\|\psi_j-r_j^{h}-(\psi_{j+1}-r_{j+1}^h)\|_Y^2\\ && +\| \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{\eta}}_j^{h}-( \mathit{\boldsymbol{R}}_{j+1}- \mathit{\boldsymbol{\eta}}_{j+1}^h)\|_H^2+\|\xi_j-z_j^{h}-(\xi_{j+1}-z_{j+1}^h)\|_Y^2\Big] \\ && +C\Big( \| \mathit{\boldsymbol{v}}^0- \mathit{\boldsymbol{v}}_0^{h}\|^2_H+\|\phi^0-\phi_0^{h}\|^2_E +\| \mathit{\boldsymbol{u}}^0- \mathit{\boldsymbol{u}}_0^{h}\|^2_V+\|\psi^0-\psi_0^{h}\|^2_Y+ \|\xi^0-\xi_0^{h}\|^2_Y\\ &&\nonumber + \|\vartheta^0-\vartheta_0^{h}\|^2_E+ \|\theta^0-\theta_0^{h}\|^2_E +\| \mathit{\boldsymbol{R}}^0- \mathit{\boldsymbol{R}}_0^{h}\|^2_H+\| \mathit{\boldsymbol{S}}^0- \mathit{\boldsymbol{S}}_0^{h}\|^2_V+\| \mathit{\boldsymbol{T}}^0- \mathit{\boldsymbol{T}}_0^{h}\|^2_V\Big), \end{eqnarray*}$

where $C > 0$ is a positive constant which is be independent of thediscretization parameters $h$ and $k$ , but depending on thecontinuous solution, $\delta \xi_j = (\xi_j-\xi_{j-1})/k$ , $\delta\phi_j = (\phi_j-\phi_{j-1})/k$ , $\delta\psi_j = (\psi_j-\psi_{j-1})/k$ , $\delta \mathit{\boldsymbol{v}}_j = (\mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{v}}_{j-1})/k$ , $\delta \mathit{\boldsymbol{u}}_j = (\mathit{\boldsymbol{u}}_j- \mathit{\boldsymbol{u}}_{j-1})/k$ , $\delta \mathit{\boldsymbol{R}}_j = (\mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{R}}_{j-1})/k$ , $\delta \mathit{\boldsymbol{T}}_j = (\mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_{j-1})/k$ , $\delta \mathit{\boldsymbol{S}}_j = (\mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_{j-1})/k$ , $\delta\vartheta_j = (\vartheta_j-\vartheta_{j-1})/k$ and $\delta\theta_j = (\theta_j-\theta_{j-1})/k$ , and $I_n$ and $J_n$ are theintegration errors defined as:

$\begin{eqnarray} && I_n = \Big\|\int_0^{t_n} \vartheta(s)\,ds-k\sum\limits_{j = 1}^n \theta_j\Big\|_E, \end{eqnarray}$

(4.15)

$\begin{eqnarray} && J_n = \Big\|\int_0^{t_n} \mathit{\boldsymbol{S}}(s)\,ds-k\sum\limits_{j = 1}^n \mathit{\boldsymbol{S}}_j\Big\|_V. \end{eqnarray}$

(4.16)

Proof. As we did in the proof of Lemma 4.1, we also assume that the time relaxation parameter $\frac{\tau_1^2}{2} = 1$ .

First, we obtain a priori error estimates on function $\mathit{\boldsymbol{v}}$ . Then, subtracting variational equation (4.1) at time $t = t_n$ for a test function $\mathit{\boldsymbol{w}} = \mathit{\boldsymbol{w}}^h\in V^h\subset V$ and discrete variational equation (4.8), we have, for all $\mathit{\boldsymbol{w}}^h\in V^h$ ,

$\begin{array}{l} \rho (\dot{ \mathit{\boldsymbol{v}}}_n-\delta{ \mathit{\boldsymbol{v}}_n^{hk}} , \mathit{\boldsymbol{w}}^h)_H +\mu (\nabla ( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \nabla \mathit{\boldsymbol{w}}^h)_Q+(\lambda+\mu) ( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \mbox{div}\, \mathit{\boldsymbol{w}}^h)_Y \\ -\mu_0(\nabla (\phi_n- \phi_n^{hk}), \mathit{\boldsymbol{w}}^h)_H+\beta_0(\nabla(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk}), \mathit{\boldsymbol{w}}^h)_H = 0, \end{array}$

and so, we obtain

$\begin{array}{l} \rho (\dot{ \mathit{\boldsymbol{v}}}_n-\delta{ \mathit{\boldsymbol{v}}_n^{hk}} , \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk})_H +\mu (\nabla ( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \nabla ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}))_Q\\ \qquad +(\lambda+\mu) ( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}))_Y-\mu_0(\nabla (\phi_n- \phi_n^{hk}), \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk})_H \\ \qquad +\beta_0(\nabla(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk}), \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk})_H\\ \quad = \rho (\dot{ \mathit{\boldsymbol{v}}}_n-\delta{ \mathit{\boldsymbol{v}}_n^{hk}} , \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h)_H +\mu (\nabla ( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \nabla ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h))_Q\\ \qquad +(\lambda+\mu) ( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h))_Y -\mu_0(\nabla (\phi_n- \phi_n^{hk}), \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h)_H\\ \qquad +\beta_0(\nabla(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk}), \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h)_H \quad \forall \mathit{\boldsymbol{w}}^h\in V^h. \end{array}$

Keeping in mind that

$\begin{array}{l} ( \dot{ \mathit{\boldsymbol{v}}}_n-\delta \mathit{\boldsymbol{v}}_n^{hk}, \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk})_H\geq ( \dot{ \mathit{\boldsymbol{v}}}_n-\delta \mathit{\boldsymbol{v}}_n, \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk})_H + \frac{1}{2k}\left[\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{v}}_{n-1}- \mathit{\boldsymbol{v}}_{n-1}^{hk}\|_H^2\right],\\[0.2cm] (\nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}),\nabla( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}))_Q\geq (\nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}),\nabla(\dot{ \mathit{\boldsymbol{u}}}_n-\delta \mathit{\boldsymbol{u}}_n))_Q\\ \qquad \qquad +\frac{1}{2k}\left[\|\nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q-\|\nabla( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Q\right], \\ ( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}))_Y\geq ( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}), \mbox{div}\,(\dot{ \mathit{\boldsymbol{u}}}_n-\delta \mathit{\boldsymbol{u}}_n))_Y\\ \qquad \qquad +\frac{1}{2k}\left[\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Y\right],\\[0.2cm] \beta_0(\nabla(\xi_n-\xi_n^{hk}), \mathit{\boldsymbol{w}})_H = -\beta_0(\xi_n-\xi_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{w}})_Y, \end{array}$

using inequality (4.14) several times, after some easy calculations we get

$\begin{eqnarray} &&\frac{\rho}{2k}\left[\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{v}}_{n-1}- \mathit{\boldsymbol{v}}_{n-1}^{hk}\|^2_H\right]-\beta_0(\xi_n-\xi_n^{hk}, \mbox{div}\,( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}))_Y\\ && +\frac{\mu}{2k} \left[\| \nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q-\| \nabla( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Q\right]\\ && +\frac{\lambda+\mu}{2k} \left[\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Y\right]\\ && \leq C\Big(\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\|\dot{ \mathit{\boldsymbol{v}}}_n-\delta \mathit{\boldsymbol{v}}_n\|^2_H+\|\nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y\\ && + \| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h\|_V^2+\|\dot{ \mathit{\boldsymbol{u}}}_n-\delta \mathit{\boldsymbol{u}}_n\|^2_V+ (\delta \mathit{\boldsymbol{v}}_n-\delta \mathit{\boldsymbol{v}}_n^{hk}, \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h)_H+\|\xi_n-\xi_n^{hk}\|^2_Y\\ && +\|\nabla(\phi_n-\phi_n^{hk})\|^2_H+\|\nabla(\theta_n-\theta_n^{hk})\|^2_H+\|\nabla(\theta_n-\theta_n^{hk})\|^2_H\Big) \quad \forall \mathit{\boldsymbol{w}}^h\in V^h. \end{eqnarray}$

(4.17)

Secondly, we obtain the error estimates on function $\psi$ . If we subtract variational equation (4.2) at time $t = t_n$ for a test function $r = r^h\in E^h\subset E$ and discrete variational equation (4.9), we find that, for all $r^h\in E^h$ ,

$\begin{array}{l} J(\dot{\psi}_n-\delta{\psi_n^{hk}},r^h)_Y+a_0(\nabla (\phi_n-\phi_n^{hk}),\nabla r^h)_H+\zeta (\phi_n-\phi_n^{hk},r^h)_Y\\ \qquad +\mu_0( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}),r^h)_Y -\beta_1(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk},r^h)_Y\\ \qquad +\mu_2( \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})+\tau_1 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})+ \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),r^h)_Y = 0 . \end{array}$

Then, we have, for all $r^h\in E^h$ ,

$\begin{array}{l} J(\dot{\psi}_n-\delta{\psi_n^{hk}},\psi_n-\psi_n^{hk})_Y+a_0(\nabla (\phi_n-\phi_n^{hk}),\nabla (\psi_n-\psi_n^{hk}))_H\\ \qquad +\zeta (\phi_n-\phi_n^{hk},\psi_n-\psi_n^{hk})_Y+\mu_0( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}),\psi_n-\psi_n^{hk})_Y\\ \qquad +\mu_2( \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})+\tau_1 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})+ \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),\psi_n-\psi_n^{hk})_Y\\ \qquad -\beta_1(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk},\psi_n-\psi_n^{hk})_Y\\ \quad = J(\dot{\psi}_n-\delta{\psi_n^{hk}},\psi_n-r^h)_Y+a_0(\nabla (\phi_n-\phi_n^{hk}),\nabla (\psi_n-r^h))_H\\ \qquad +\zeta (\phi_n-\phi_n^{hk},\psi_n-r^h)_Y+\mu_0( \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}),\psi_n-r^h)_Y\\ \qquad +\mu_2( \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})+\tau_1 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})+ \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),\psi_n-r^h)_Y\\ \qquad -\beta_1(\xi_n-\xi_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})+\theta_n-\theta_n^{hk},\psi_n-r^h)_Y, \end{array}$

and, taking into account that

$\begin{array}{l} ( \dot{\psi}_n-\delta \psi_n^{hk},\psi_n-\psi_n^{hk})_Y\geq ( \dot{\psi}_n-\delta \psi_n,\psi_n-\psi_n^{hk})_Y + \frac{1}{2k}\left[\|\psi_n-\psi_n^{hk}\|^2_Y-\|\psi_{n-1}-\psi_{n-1}^{hk}\|^2_Y\right],\\ (\nabla(\phi_n-\phi_n^{hk}),\nabla(\dot{\phi}_n-\psi_n^{hk}))_H = (\nabla(\phi_n-\phi_n^{hk}),\nabla(\dot{\phi}_n-\delta \phi_n))_H\\ \qquad +\frac{1}{2k}\Big[\|\nabla(\phi_n-\phi_n^{hk})\|^2_H-\|\nabla(\phi_{n-1}-\phi_{n-1}^{hk})\|_H^2 +\|\nabla(\phi_n-\phi_n^{hk}-(\phi_{n-1}-\phi_{n-1}^{hk}))\|^2_H\Big],\\ (\phi_n-\phi_n^{hk},\dot{\phi}_n-\psi_n^{hk})_Y\geq (\phi_n-\phi_n^{hk},\dot{\phi}_n-\delta \phi_n)_Y +\frac{1}{2k}\Big[\|\phi_n-\phi_n^{hk}\|^2_Y-\|\phi_{n-1}-\phi_{n-1}^{hk}\|^2_Y\Big],\\ ( \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}),\psi_n-r^h)_Y = -( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk},\nabla (\psi_n-r^h))_H, \end{array}$

using again several times inequality (4.14) we get, for all $r^h\in E^h$ ,

$\begin{eqnarray} &&\frac{J}{2k}\left\{\|\psi_n-\psi_n^{hk}\|^2_Y-\|\psi_{n-1}-\psi_{n-1}^{hk}\|^2_Y\right\}+\frac{\zeta}{2k}\Big\{\|\phi_n-\phi_n^{hk}\|^2_Y -\|\phi_{n-1}-\phi_{n-1}^{hk}\|^2_Y\Big\}\\ && \quad + \frac{a_0}{2k}\Big\{\|\nabla(\phi_n-\phi_n^{hk})\|_H^2-\|\nabla(\phi_{n-1}-\phi_{n-1}^{hk})\|^2_H\Big\}\\ && \quad +\mu_2( \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}),\psi_n-\psi_n^{hk})_Y \\ && \leq C\Big(\|\dot{\psi}_n-\delta \psi_n\|^2_Y+\|\dot{\phi}_n-\delta \phi_n\|^2_E+\|\psi_n-r^h\|^2_E+\|\vartheta_n-\vartheta_n^{hk}\|_Y^2\\ &&\quad +\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y+\|\theta_n-\theta_n^{hk}\|^2_Y+\|\xi_n-\xi_n^{hk}\|^2+\| \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|_Y^2\\ && \quad +\|\phi_n-\phi_n^{hk}\|_Y^2+\|\nabla(\phi_n-\phi_n^{hk})\|_H^2+ \| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|_Y^2+\|\psi_n-\psi_n^{hk}\|^2_Y\\ &&\quad+\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|_H^2+(\delta \psi_n-\delta \psi_n^{hk},\psi_n-r^h)_Y\Big). \end{eqnarray}$

(4.18)

Thirdly, we derive the estimates on the temperature acceleration $\xi$ . Therefore, subtracting variational equation (4.3) at time $t = t_n$ for a test function $z = z^h\in E^h\subset E$ and discrete variational equation (4.10), for all $z^h\in E^h$ we find that

$\begin{array}{l} a(\dot{\xi}_n-\delta{\xi_n^{hk}}+\tau_1(\xi_n-\xi_n^{hk})+\vartheta_n-\vartheta_n^{hk},z^h)_Y+\kappa(\nabla (\theta_n-\theta_n^{hk}+\tau_2 (\vartheta_n-\vartheta_n^{hk})),\nabla z^h)_H\\ \qquad +\beta_0( \mbox{div}\, ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}),z^h)_Y +\beta_1 (\psi_n-\psi_n^{hk},z^h )_Y\\ \qquad -\kappa_1 ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})+\tau_2 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}),z^h)_Y = 0. \end{array}$

Thus, it follows that

$\begin{array}{l} a(\dot{\xi}_n-\delta{\xi_n^{hk}}+\tau_1(\xi_n-\xi_n^{hk})+\vartheta_n-\vartheta_n^{hk},\xi_n-\xi_n^{hk})_Y+\beta_0( \mbox{div}\, ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}),\xi_n-\xi_n^{hk})_Y\\ \qquad +\beta_1 (\psi_n-\psi_n^{hk},\xi_n-\xi_n^{hk} )_Y -\kappa_1 ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})+\tau_2 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}),\xi_n-\xi_n^{hk})_Y\\ \qquad +\kappa(\nabla (\theta_n-\theta_n^{hk}+\tau_2 (\vartheta_n-\vartheta_n^{hk})),\nabla (\xi_n-\xi_n^{hk}))_H\\ \quad = a(\dot{\xi}_n-\delta{\xi_n^{hk}}+\tau_1(\xi_n-\xi_n^{hk})+\vartheta_n-\vartheta_n^{hk},\xi_n-z^h)_Y+\beta_0( \mbox{div}\, ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}),\xi_n-z^h)_Y\\ \qquad +\beta_1 (\psi_n-\psi_n^{hk},\xi_n-z^h )_Y -\kappa_1 ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})+\tau_2 \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}),\xi_n-z^h)_Y\\ \qquad +\kappa(\nabla (\theta_n-\theta_n^{hk}+\tau_2 (\vartheta_n-\vartheta_n^{hk})),\nabla (\xi_n-z^h))_H \quad \forall z^h\in E^h. \end{array}$

Using the following estimates

$\begin{array}{l} ( \dot{\xi}_n-\delta \xi_n^{hk},\xi_n-\xi_n^{hk})_Y\geq ( \dot{\xi}_n-\delta \xi_n,\xi_n-\xi_n^{hk})_Y + \frac{1}{2k}\left[\|\xi_n-\xi_n^{hk}\|^2_Y-\|\xi_{n-1}-\xi_{n-1}^{hk}\|^2_Y\right],\\ (\nabla(\vartheta_n-\vartheta_n^{hk}),\nabla(\dot{\vartheta}_n-\xi_n^{hk}))_H\geq (\nabla(\vartheta_n-\vartheta_n^{hk}),\nabla(\dot{\vartheta}_n-\delta \vartheta_n))_H\\ \qquad +\frac{1}{2k}\Big[\|\nabla(\vartheta_n-\vartheta_n^{hk})\|^2_H-\|\nabla(\vartheta_{n-1}-\vartheta_{n-1}^{hk})\|_H^2 \Big],\\ (\vartheta_n-\vartheta_n^{hk},\dot{\vartheta}_n-\xi_n^{hk})_Y\geq (\vartheta_n-\vartheta_n^{hk},\dot{\vartheta}_n-\delta \vartheta_n)_Y +\frac{1}{2k}\Big[\|\vartheta_n-\vartheta_n^{hk}\|^2_Y-\|\vartheta_{n-1}-\vartheta_{n-1}^{hk}\|^2_Y\Big],\\ ( \mbox{div}\, ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}),\xi_n-z^h)_Y = -( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk},\nabla (\xi_n-z^h))_H, \end{array}$

applying several times inequality (4.14) we find that

$\begin{eqnarray} &&\frac{a}{2k}\left\{\|\xi_n-\xi_n^{hk}\|^2_Y-\|\xi_{n-1}-\xi_{n-1}^{hk}\|^2_Y\right\} +\frac{a}{2k}\Big\{\|\vartheta_n-\vartheta_n^{hk}\|^2_Y -\|\vartheta_{n-1}-\vartheta_{n-1}^{hk}\|^2_Y\Big\}\\ && \quad + \frac{\kappa}{2k}\Big\{\|\nabla(\vartheta_n-\vartheta_n^{hk})\|_H^2-\|\nabla(\vartheta_{n-1}-\vartheta_{n-1}^{hk})\|^2_H\Big\}\\ && \quad +\beta_0( \mbox{div}\,( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}),\xi_n-\xi_n^{hk})_Y+(\nabla (\theta_n-\theta_n^{hk}),\nabla (\xi_n-\xi_n^{hk}))_H \\ && \leq C\Big(\|\dot{\xi}_n-\delta \xi_n\|^2_Y+\|\nabla(\theta_n-\theta_n^{hk})\|_H^2+\|\dot{\vartheta}_n-\delta \vartheta_n\|^2_E+\|\xi_n-z^h\|^2_E\\ &&\quad +\|\nabla(\vartheta_n-\vartheta_n^{hk})\|_H^2 +\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\|\xi_n-\xi_n^{hk}\|^2+\| \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|_Y^2\\ &&\quad +\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|_Y^2+\|\psi_n-\psi_n^{hk}\|_Y^2+(\delta \xi_n-\delta \xi_n^{hk},\xi_n-z^h)_Y\Big). \end{eqnarray}$

(4.19)

Finally, we get the error estimates for the microtemperature acceleration $\mathit{\boldsymbol{R}}$ . Then, if we subtract variational equation (4.4) at time $t = t_n$ for a test function $\mathit{\boldsymbol{\eta}} = \mathit{\boldsymbol{\eta}}^h\in V^h\subset V$ and discrete variational equation (4.11), we have, for all $\psi^h\in V^h$ ,

and therefore, we obtain

$\begin{array}{l} b(\dot{ \mathit{\boldsymbol{R}}}_n- \delta{ \mathit{\boldsymbol{R}}_n^{hk}}+\tau_1( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})+ \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H\\ \qquad +\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 ( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})),\nabla ( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Q\\ \qquad +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})), \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Y\\ \qquad +\kappa_2( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 ( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H+\mu_2(\nabla (\psi_n-\psi_n^{hk}), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H\\ \qquad +\kappa_1(\nabla( \theta_n-\theta_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H \\ \quad = b(\dot{ \mathit{\boldsymbol{R}}}_n- \delta{ \mathit{\boldsymbol{R}}_n^{hk}}+\tau_1( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})+ \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H\\ \qquad +\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 ( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})),\nabla ( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h))_Q\\ \qquad +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})), \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h))_Y\\ \qquad +\kappa_2( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}+\tau_2 ( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H+\mu_2(\nabla (\psi_n-\psi_n^{hk}), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H\\ \qquad +\kappa_1(\nabla( \theta_n-\theta_n^{hk}+\tau_1(\vartheta_n-\vartheta_n^{hk})), \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H\quad \forall \mathit{\boldsymbol{\eta}}^h\in V^h. \end{array}$

Keeping in mind that

$\begin{array}{l} ( \dot{ \mathit{\boldsymbol{R}}}_n-\delta \mathit{\boldsymbol{R}}_n^{hk}, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H\geq ( \dot{ \mathit{\boldsymbol{R}}}_n-\delta \mathit{\boldsymbol{R}}_n, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk})_H + \frac{1}{2k}\left[\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}- \mathit{\boldsymbol{R}}_{n-1}^{hk}\|_H^2\right],\\[0.2cm] (\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}),\nabla( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Q\geq (\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}),\nabla(\dot{ \mathit{\boldsymbol{S}}}_n-\delta \mathit{\boldsymbol{S}}_n))_Q\\ \qquad \qquad +\frac{1}{2k}\left[\|\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q-\|\nabla( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Q\right],\\[0.2cm] ( \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Y\geq ( \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}), \mbox{div}\,(\dot{ \mathit{\boldsymbol{S}}}_n-\delta \mathit{\boldsymbol{S}}_n))_Y\\ \qquad \qquad +\frac{1}{2k}\left[\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Y\right],\\ (\nabla (\psi_n-\psi_n^{hk}), \mathit{\boldsymbol{\eta}})_H = -(\psi_n-\psi_n^{hk}, \mbox{div}\, \mathit{\boldsymbol{\eta}})_Y, \end{array}$

using again inequality (4.14), after easy calculations it leads to

$\begin{eqnarray} &&\frac{b}{2k}\left[\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}- \mathit{\boldsymbol{R}}_{n-1}^{hk}\|^2_H\right]-\mu_2(\psi_n-\psi_n^{hk}, \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Y\\ && +\frac{\kappa_6}{2k} \left[\| \nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q-\| \nabla( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Q\right]\\ && +\frac{\kappa_4+\kappa_5}{2k} \left[\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Y\right]\\ && +\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),\nabla ( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Q +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Y\\ && \leq C\Big(\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H+\|\dot{ \mathit{\boldsymbol{R}}}_n-\delta \mathit{\boldsymbol{R}}_n\|^2_H+\|\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Y\\ && + \| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h\|_V^2+\|\dot{ \mathit{\boldsymbol{S}}}_n-\delta \mathit{\boldsymbol{S}}_n\|^2_V+ (\delta \mathit{\boldsymbol{R}}_n-\delta \mathit{\boldsymbol{R}}_n^{hk}, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H+\| \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}\|_H^2\\ && +\|\psi_n-\psi_n^{hk}\|^2_Y+\|\nabla(\theta_n-\theta_n^{hk})\|^2_H+\|\nabla(\vartheta_n-\vartheta_n^{hk})\|^2_H+\| \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}\|_H^2\\ && +\|\nabla( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|^2_Y \Big)\quad \forall \mathit{\boldsymbol{\eta}}^h\in V^h. \end{eqnarray}$

(4.20)

Combining estimates (4.17)–(4.20) it follows that, for all $\mathit{\boldsymbol{w}}^h, \, \mathit{\boldsymbol{\eta}}^h\in V^h$ and $r^h, \, z^h\in E^h$ ,

$\begin{eqnarray*} &&\nonumber\frac{\rho}{2k}\left[\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{v}}_{n-1}- \mathit{\boldsymbol{v}}_{n-1}^{hk}\|^2_H\right]+\frac{\zeta}{2k}\Big\{\|\phi_n-\phi_n^{hk}\|^2_Y -\|\phi_{n-1}-\phi_{n-1}^{hk}\|^2_Y\Big\}\\ &&\nonumber \quad +\frac{\mu}{2k} \left[\| \nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q-\| \nabla( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Q\right]\\ &&\nonumber \quad +\frac{\lambda+\mu}{2k} \left[\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{u}}_{n-1}- \mathit{\boldsymbol{u}}_{n-1}^{hk})\|^2_Y\right] \\ &&\nonumber \quad + \frac{J}{2k}\left\{\|\psi_n-\psi_n^{hk}\|^2_Y-\|\psi_{n-1}-\psi_{n-1}^{hk}\|^2_Y\right\} +(\nabla (\theta_n-\theta_n^{hk}),\nabla (\xi_n-\xi_n^{hk}))_H\\ && \nonumber\quad + \frac{a_0}{2k}\Big\{\|\nabla(\phi_n-\phi_n^{hk})\|_H^2-\|\nabla(\phi_{n-1}-\phi_{n-1}^{hk})\|^2_H\Big\}\\ &&\quad +\nonumber \frac{a}{2k}\left\{\|\xi_n-\xi_n^{hk}\|^2_Y-\|\xi_{n-1}-\xi_{n-1}^{hk}\|^2_Y\right\} +\frac{a}{2k}\Big\{\|\vartheta_n-\vartheta_n^{hk}\|^2_Y -\|\vartheta_{n-1}-\vartheta_{n-1}^{hk}\|^2_Y\Big\}\\ && \nonumber\quad + \frac{\kappa}{2k}\Big\{\|\nabla(\vartheta_n-\vartheta_n^{hk})\|_H^2-\|\nabla(\vartheta_{n-1}-\vartheta_{n-1}^{hk})\|^2_H\Big\} \\ &&\nonumber\quad + \frac{b}{2k}\left[\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H-\| \mathit{\boldsymbol{R}}_{n-1}- \mathit{\boldsymbol{R}}_{n-1}^{hk}\|^2_H\right]\\ &&\nonumber \quad +\frac{\kappa_6}{2k} \left[\| \nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q-\| \nabla( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Q\right]\\ &&\nonumber \quad +\frac{\kappa_4+\kappa_5}{2k} \left[\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Y-\| \mbox{div}\,( \mathit{\boldsymbol{S}}_{n-1}- \mathit{\boldsymbol{S}}_{n-1}^{hk})\|^2_Y\right]\\ &&\nonumber \quad +\kappa_6(\nabla ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),\nabla ( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Q +(\kappa_4+\kappa_5) ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}))_Y\\ &&\nonumber \leq C\Big(\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\|\dot{ \mathit{\boldsymbol{v}}}_n-\delta \mathit{\boldsymbol{v}}_n\|^2_H+\|\nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y\\ &&\nonumber + \| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h\|_V^2+\|\dot{ \mathit{\boldsymbol{u}}}_n-\delta \mathit{\boldsymbol{u}}_n\|^2_V+ (\delta \mathit{\boldsymbol{v}}_n-\delta \mathit{\boldsymbol{v}}_n^{hk}, \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}^h)_H+\|\xi_n-\xi_n^{hk}\|^2_Y\\ && +\|\nabla(\phi_n-\phi_n^{hk})\|^2_H+\|\nabla(\theta_n-\theta_n^{hk})\|^2_H + \|\dot{\psi}_n-\delta \psi_n\|^2_Y+\|\dot{\phi}_n-\delta \phi_n\|^2_E\\ && +\|\psi_n-r^h\|^2_E+\|\vartheta_n-\vartheta_n^{hk}\|_Y^2 +\|\theta_n-\theta_n^{hk}\|^2_Y+\| \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|_Y^2\\ &&\nonumber + \| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|_Y^2+\|\psi_n-\psi_n^{hk}\|^2_Y+\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|_H^2+(\delta \psi_n-\delta \psi_n^{hk},\psi_n-r^h)_Y\\ && + \|\dot{\xi}_n-\delta \xi_n\|^2_Y+\|\dot{\vartheta}_n-\delta \vartheta_n\|^2_E+\|\xi_n-z^h\|^2_E +\|\nabla(\vartheta_n-\vartheta_n^{hk})\|_H^2 \nonumber \\ && +\|\psi_n-\psi_n^{hk}\|_Y^2+(\delta \xi_n-\delta \xi_n^{hk},\xi_n-z^h)_Y +\|\dot{ \mathit{\boldsymbol{R}}}_n-\delta \mathit{\boldsymbol{R}}_n\|^2_H+\|\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q\\ &&\nonumber + \| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h\|_V^2+\|\dot{ \mathit{\boldsymbol{S}}}_n-\delta \mathit{\boldsymbol{S}}_n\|^2_V+ (\delta \mathit{\boldsymbol{R}}_n-\delta \mathit{\boldsymbol{R}}_n^{hk}, \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{\eta}}^h)_H+\| \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}\|_H^2\\ &&\nonumber +\|\psi_n-\psi_n^{hk}\|^2_Y+\| \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}\|_H^2+\|\nabla( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk})\|^2_Q \Big). \end{eqnarray*}$

Multiplying the above estimates by $k$ and summing up to $n$ it follows that

$\begin{eqnarray*} &&\nonumber\| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|^2_H+\|\phi_n-\phi_n^{hk}\|^2_Y +\| \nabla( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk})\|^2_Y+\|\psi_n-\psi_n^{hk}\|^2_Y\\ &&\nonumber \quad +\|\nabla(\phi_n-\phi_n^{hk})\|_H^2+ \|\xi_n-\xi_n^{hk}\|^2_Y + Ck\sum\limits_{j = 1}^n(\nabla (\theta_j-\theta_j^{hk}),\nabla (\xi_j-\xi_j^{hk}))_H\\ &&\quad +\|\vartheta_n-\vartheta_n^{hk}\|^2_Y +\|\nabla(\vartheta_n-\vartheta_n^{hk})\|_H^2+\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|^2_H+\| \nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Q\\ &&\nonumber \quad +\| \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk})\|^2_Y+Ck\sum\limits_{j = 1}^n (\nabla ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}),\nabla ( \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{R}}_j^{hk}))_Q\\ &&\nonumber \quad +Ck\sum\limits_{j = 1}^n ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}), \mbox{div}\,( \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{R}}_j^{hk}))_Y \\ &&\nonumber \leq Ck\sum\limits_{j = 1}^n\Big(\| \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{v}}_j^{hk}\|^2_H+\|\dot{ \mathit{\boldsymbol{v}}}_j-\delta \mathit{\boldsymbol{v}}_j\|^2_H+\|\nabla( \mathit{\boldsymbol{u}}_j- \mathit{\boldsymbol{u}}_j^{hk})\|^2_Q+\| \mbox{div}\,( \mathit{\boldsymbol{u}}_j- \mathit{\boldsymbol{u}}_j^{hk})\|^2_Y\\ &&\nonumber + \| \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^h\|_V^2+\|\dot{ \mathit{\boldsymbol{u}}}_j-\delta \mathit{\boldsymbol{u}}_j\|^2_V+ (\delta \mathit{\boldsymbol{v}}_j-\delta \mathit{\boldsymbol{v}}_j^{hk}, \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}^h_j)_H+\|\xi_j-\xi_j^{hk}\|^2_Y\\ && +\|\nabla(\phi_j-\phi_j^{hk})\|^2_H+\|\nabla(\theta_j-\theta_j^{hk})\|^2_H + \|\dot{\psi}_j-\delta \psi_j\|^2_Y+\|\dot{\phi}_j-\delta \phi_j\|^2_E\\ && +\|\psi_j-r^h_j\|^2_E+\|\vartheta_j-\vartheta_j^{hk}\|_Y^2 +\|\theta_j-\theta_j^{hk}\|^2_Y+\| \mbox{div}\,( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk})\|_Y^2\\ &&\nonumber + \| \mbox{div}\,( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk})\|_Y^2+\|\psi_j-\psi_j^{hk}\|^2_Y+\| \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{R}}_j^{hk}\|_H^2+(\delta \psi_j-\delta \psi_j^{hk},\psi_j-r^h_j)_Y\\ && + \|\dot{\xi}_j-\delta \xi_j\|^2_Y+\|\dot{\vartheta}_j-\delta \vartheta_j\|^2_E+\|\xi_j-z^h_j\|^2_E +\|\nabla(\vartheta_j-\vartheta_j^{hk})\|_H^2 \nonumber \\ && +\|\psi_j-\psi_j^{hk}\|_Y^2+(\delta \xi_j-\delta \xi_j^{hk},\xi_j-z_j^h)_Y +\|\dot{ \mathit{\boldsymbol{R}}}_j-\delta \mathit{\boldsymbol{R}}_j\|^2_H+\|\nabla( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk})\|^2_Q\\ &&\nonumber + \| \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{\eta}}_j^h\|_V^2+\|\dot{ \mathit{\boldsymbol{S}}}_j-\delta \mathit{\boldsymbol{S}}_j\|^2_V+ (\delta \mathit{\boldsymbol{R}}_j-\delta \mathit{\boldsymbol{R}}_j^{hk}, \mathit{\boldsymbol{R}}_j- \mathit{\boldsymbol{\eta}}^h_j)_H+\| \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}\|_H^2\\ &&\nonumber +\|\psi_j-\psi_j^{hk}\|^2_Y+\| \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}\|_H^2+\|\nabla( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk})\|^2_Q \Big)\\ && +C\Big( \| \mathit{\boldsymbol{v}}^0- \mathit{\boldsymbol{v}}_0^{h}\|^2_H+\|\phi^0-\phi_0^{h}\|^2_E +\| \mathit{\boldsymbol{u}}^0- \mathit{\boldsymbol{u}}_0^{h}\|^2_V+\|\psi^0-\psi_0^{h}\|^2_Y+ \|\xi^0-\xi_0^{h}\|^2_Y\\ &&\nonumber + \|\vartheta^0-\vartheta_0^{h}\|^2_E +\| \mathit{\boldsymbol{R}}^0- \mathit{\boldsymbol{R}}_0^{h}\|^2_H+\| \mathit{\boldsymbol{S}}^0- \mathit{\boldsymbol{S}}_0^{h}\|^2_V\Big). \end{eqnarray*}$

Taking into account that

$\begin{array}{l} k\sum\limits_{j = 1}^n (\delta \mathit{\boldsymbol{v}}_j-\delta \mathit{\boldsymbol{v}}_j^{hk}, \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^h)_H = \sum\limits_{j = 1}^n ( \mathit{\boldsymbol{v}}_{j}- \mathit{\boldsymbol{v}}_{j}^{hk}-( \mathit{\boldsymbol{v}}_{j-1}- \mathit{\boldsymbol{v}}_{j-1}^{hk}), \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^h)_H\\ \qquad = ( \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}, \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{w}}_n^h)_H+( \mathit{\boldsymbol{v}}_0^h- \mathit{\boldsymbol{v}}_0, \mathit{\boldsymbol{v}}_1- \mathit{\boldsymbol{w}}_1^h)_H\\ \qquad \qquad + \sum\limits_{j = 1}^{n-1} ( \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{v}}_j^{hk}, \mathit{\boldsymbol{v}}_j- \mathit{\boldsymbol{w}}_j^h-( \mathit{\boldsymbol{v}}_{j+1}- \mathit{\boldsymbol{w}}_{j+1}^h))_H, \end{array}$

where similar estimates can be derived for the terms involving functions $\psi$ , $\xi$ and $\mathit{\boldsymbol{R}}$ , and using the estimates:

$\begin{array}{l} k\sum\limits_{j = 1}^n (\nabla (\theta_j-\theta_j^{hk}),\nabla (\delta \vartheta_j-\delta\vartheta_j^{hk}))_H = \sum\limits_{j = 1}^n (\nabla (\theta_j-\theta_j^{hk}),\nabla (\vartheta_j-\vartheta_j^{hk}-(\vartheta_{j-1}-\vartheta_{j-1}^{hk})))_H \\ \quad = (\nabla(\theta_n-\theta_n^{hk}),\nabla(\vartheta_n-\vartheta_n^{hk}))_H+ (\nabla(\theta_1^{hk}-\theta_1),\nabla(\vartheta^0-\vartheta_0^{h}))_H\\ \qquad + \sum\limits_{j = 1}^{n-1}(\nabla(\theta_j-\theta_j^{hk}-(\theta_{j-1}-\theta_{j-1}^{hk})),\nabla(\vartheta_j-\vartheta_j^{hk}))_H, \\ \sum\limits_{j = 1}^{n-1}(\nabla(\theta_j-\theta_j^{hk}-(\theta_{j-1}-\theta_{j-1}^{hk})),\nabla(\vartheta_j-\vartheta_j^{hk}))_H\\ \qquad \leq C\Big(k\sum\limits_{j = 1}^n\|\nabla(\dot{\theta}_j-\delta\theta_j)\|_H^2+k\sum\limits_{j = 1}^n\|\nabla(\vartheta_j-\vartheta_j^{hk})\|_H^2\Big), \\ \|\theta_n-\theta_n^{hk}\|_E^2\leq C\Big(\|\theta^0-\theta_0^h\|_E^2+I_n^2+k\sum\limits_{j = 1}^n\|\vartheta_j-\vartheta_j^{hk}\|_E^2\Big), \\ k\sum\limits_{j = 1}^n (\nabla ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}),\nabla (\delta \mathit{\boldsymbol{S}}_j-\delta \mathit{\boldsymbol{S}}_j^{hk}))_Q = \sum\limits_{j = 1}^n (\nabla ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}),\nabla ( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}-( \mathit{\boldsymbol{S}}_{j-1}- \mathit{\boldsymbol{S}}_{j-1}^{hk})))_Q \\ \quad = (\nabla( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}),\nabla( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}))_Q+ (\nabla( \mathit{\boldsymbol{T}}_1^{hk}- \mathit{\boldsymbol{T}}_1),\nabla( \mathit{\boldsymbol{S}}^0- \mathit{\boldsymbol{S}}_0^{h}))_Q\\ \qquad + \sum\limits_{j = 1}^{n-1}(\nabla( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}-( \mathit{\boldsymbol{T}}_{j-1}- \mathit{\boldsymbol{T}}_{j-1}^{hk})),\nabla( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}))_Q, \end{array}$

$\begin{array}{l} \sum\limits_{j = 1}^{n-1}(\nabla( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}-( \mathit{\boldsymbol{T}}_{j-1}- \mathit{\boldsymbol{T}}_{j-1}^{hk})),\nabla( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}))_Q\\ \qquad \leq C\Big(k\sum\limits_{j = 1}^n\|\nabla(\dot{ \mathit{\boldsymbol{T}}}_j-\delta \mathit{\boldsymbol{T}}_j)\|_Q^2+k\sum\limits_{j = 1}^n\|\nabla( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk})\|_Q^2\Big),\\ \| \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}\|_V^2\leq C\Big(\| \mathit{\boldsymbol{T}}^0- \mathit{\boldsymbol{T}}_0^h\|_V^2+J_n^2+k\sum\limits_{j = 1}^n\| \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}\|_V^2\Big),\\ k\sum\limits_{j = 1}^n ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}), \mbox{div}\, (\delta \mathit{\boldsymbol{S}}_j-\delta \mathit{\boldsymbol{S}}_j^{hk}))_Y = \sum\limits_{j = 1}^n ( \mbox{div}\, ( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}), \mbox{div}\, ( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}-( \mathit{\boldsymbol{S}}_{j-1}- \mathit{\boldsymbol{S}}_{j-1}^{hk})))_Y \\ \quad = ( \mbox{div}\,( \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}), \mbox{div}\,( \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}))_Y+ ( \mbox{div}\,( \mathit{\boldsymbol{T}}_1^{hk}- \mathit{\boldsymbol{T}}_1), \mbox{div}\,( \mathit{\boldsymbol{S}}^0- \mathit{\boldsymbol{S}}_0^{h}))_Y\\ \qquad + \sum\limits_{j = 1}^{n-1}( \mbox{div}\,( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}-( \mathit{\boldsymbol{T}}_{j-1}- \mathit{\boldsymbol{T}}_{j-1}^{hk})), \mbox{div}\,( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}))_Y,\\ \sum\limits_{j = 1}^{n-1}( \mbox{div}\,( \mathit{\boldsymbol{T}}_j- \mathit{\boldsymbol{T}}_j^{hk}-( \mathit{\boldsymbol{T}}_{j-1}- \mathit{\boldsymbol{T}}_{j-1}^{hk})), \mbox{div}\,( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk}))_Y\\ \qquad \leq C\Big(k\sum\limits_{j = 1}^n\| \mbox{div}\,(\dot{ \mathit{\boldsymbol{T}}}_j-\delta \mathit{\boldsymbol{T}}_j)\|_Y^2+k\sum\limits_{j = 1}^n\| \mbox{div}\,( \mathit{\boldsymbol{S}}_j- \mathit{\boldsymbol{S}}_j^{hk})\|_Y^2\Big), \end{array}$

where $I_n$ and $J_n$ are the integration errors defined in (4.15) and (4.16), respectively, applying a discrete version of Gronwall's inequality (see, again, ^[31]) we conclude the desired a priori error estimates.

We note that the error estimates provided in Theorem 4.2 can be used to obtain the convergence order of the approximations provided by the fully discrete problem $VP^{hk}$ . As a particular case, if we assume the additional regularity:

$\begin{equation} \begin{array}{l} \phi \in H^3(0,T;Y)\cap W^{1,\infty}(0,T;H^2(B))\cap H^2(0,T;H^1(B)),\\ \theta\in H^4(0,T;Y)\cap W^{2,\infty}(0,T;H^2(B))\cap H^3(0,T;H^1(B)),\\ \mathit{\boldsymbol{u}}\in H^3(0,T;H)\cap W^{1,\infty}(0,T;[H^2(B)]^d)\cap H^2(0,T;[H^1(B)]^d),\\ \mathit{\boldsymbol{T}}\in H^4(0,T;H)\cap W^{2,\infty}(0,T;[H^2(B)]^d)\cap H^3(0,T;[H^1(B)]^d), \end{array} \end{equation}$

(4.21)

we can conclude the following.

Theorem 4.3. Under the assumptions of Theorem 4.2 and theadditional regularity conditions (4.21), it follows thelinear convergence of the approximations given by Problem VP $^{hk}$ ; that is, there exists a constant $C > 0$ , independent of parameters $h$ and $k$ , such that

$\begin{array}{l} \max\limits_{0\leq n\leq N}\Big\{ \| \mathit{\boldsymbol{v}}_n- \mathit{\boldsymbol{v}}_n^{hk}\|_H+\| \mathit{\boldsymbol{u}}_n- \mathit{\boldsymbol{u}}_n^{hk}\|_V +\|\psi_n-\psi_n^{hk}\|_Y+ \|\phi_n-\phi_n^{hk}\|_E+ \|\theta_n-\theta_n^{hk}\|_E+ \|\vartheta_n-\vartheta_n^{hk}\|_E\\ \qquad \qquad +\|\xi_n-\xi_n^{hk}\|_Y+\| \mathit{\boldsymbol{T}}_n- \mathit{\boldsymbol{T}}_n^{hk}\|_V +\| \mathit{\boldsymbol{S}}_n- \mathit{\boldsymbol{S}}_n^{hk}\|_V+\| \mathit{\boldsymbol{R}}_n- \mathit{\boldsymbol{R}}_n^{hk}\|_H\Big\}\leq C(h+k).\end{array}$

5. Numerical results

In this section, we describe some numerical results obtained solving one- and two-dimensional problems. First, the numerical convergence and the discrete energy decay are shown in the one-dimensional example. Secondly, the numerical dependence on the coupling parameter $\mu_2$ is studied in a two-dimensional example.

5.1. First example: numerical convergence in a one-dimensional example

As an academical example, in order to show the accuracy of the approximations the following one-dimensional problem is considered:

$\begin{eqnarray*} && \rho \ddot {u} = \mu u_{xx}+\mu_0 \phi_x-\beta_0 (\theta_{x}+\tau_1 \dot {\theta}_{x}+\dfrac{\tau_1^2}{2}\ddot{\theta}_{x}) +F_1 \quad \hbox{in } (0,\ell)\times (0,T_f),\\ && J\ddot{\phi} = a0\phi_{xx}-\mu_0 u_x-\mu_2(T_{x}+\tau_1 \dot {T}_{x}+\dfrac{\tau_1^2}{2}\ddot{T}_{x})+\beta_1(\theta+\tau_1 \dot {\theta} +\dfrac{\tau_1^2}{2}\ddot{\theta})-\xi\phi+F_2 \quad \hbox{in }(0,\ell)\times (0,T_f),\\ && a(\dot {\theta}+\tau_1 \ddot {\theta}+\dfrac{\tau_1^2}{2}\dddot{\theta}) = -\beta_0 \dot {u}_{x}-\beta_1\dot{\phi}+\kappa(\theta_{xx}+\tau_2\dot{\theta}_{xx})\\ &&\qquad \qquad +\kappa_1(T_x+\tau_2\dot{T}_{x}) +F_3 \quad \hbox{in } (0,\ell)\times (0,T_f),\\ && b(\dot {T}+\tau_1 \ddot {T}+\dfrac{\tau_1^2}{2}\dddot{T}) = \mu_1\dot{\phi}_{x} +\kappa_4T_{xx}+\tau_2\kappa_6\dot{T}_{xx}-\kappa_2(T+\tau_2\dot{T})\\ &&\qquad \qquad -\kappa_1(\theta_x+\tau_2\dot{\theta}_{x})-\mu_2\dot{\phi} +F_4 \quad \hbox{in } (0,\ell)\times (0,T_f),\\ &&u(0,t) = u(\ell,t) = \phi(0,t) = \phi(\ell,t) = 0\quad \hbox{for a.e. } t\in (0,T_f),\\ && \theta(0,t) = \theta(\ell,t) = R(0,t) = R(\ell,t) = 0\quad \hbox{for a.e. } t\in (0,T_f), \\ &&u^0(x) = v^0(x) = \phi^0(x) = \psi^0(x) = \theta^0(x) = x(x-1)\quad \hbox{for a.e. }x\in (0,\ell),\\ && \vartheta^0(x) = \xi^0(x) = T^0(x) = S^0(x) = R^0(x) = x(x-1)\quad \hbox{for a.e. }x\in (0,\ell). \end{eqnarray*}$

The following data have been used:

$\begin{array}{l} \ell = 1,\quad T_f = 1,\quad \rho = 1, \quad \mu = 2, \quad \mu_0 = 1,\quad \mu_2 = 1, \quad J = 1, \quad a_0 = 1, \\ \beta_0 = 1, \quad \beta_1 = 1,\quad \xi = 2 \quad \kappa = 3, \quad a = 1, \quad b = 2,\quad \kappa_1 = 1, \quad \kappa_2 = 1,\\ \kappa_4 = 2, \quad \tau_1 = 1, \quad \tau_2 = 2, \end{array}$

and the artificial supply terms $F_i$ , $i = 1, 2, 3, 4$ , are given by, for a.e. $(x, t)\in (0, 1)\times (0, 1)$ ,

$\begin{array}{l} F_1(x,t) = e^t(3x + x(x - 1) - 11/2),\\ F_2(x,t) = e^t(7x + x(x - 1)/2 - 11/2),\\ F_3(x,t) = -e^t(4x - 7x(x-1)/2 + 16),\\ F_4(x,t) = e^t(8x + 8x(x-1) - 16). \end{array}$

We note that the numerical scheme was implemented on a 3.2 Ghz PC using MATLAB, and a typical run $(h = k = 0.001)$ took about 0.13 seconds of CPU time.

In this case, the exact solution to this one-dimensional problem can be easily calculated and it has the following form, for $(x, t)\in [0, 1]\times [0, 1]$ :

$\begin{eqnarray*} u(x,t) = \phi(x,t) = \theta(x,t) = T(x,t) = e^tx(x-1). \end{eqnarray*}$

Thus, the approximation errors estimated by

$\begin{array}{l} \max\limits_{0\leq n\leq N}\Big\{ \|v_n-v_n^{hk}\|_Y+\|u_n-u_n^{hk}\|_E +\|\phi_n-\phi_n^{hk}\|_E+\|\psi_n-\psi_n^{hk}\|_Y+ \|\theta_n-\theta_n^{hk}\|_E\\ \qquad + \|\vartheta_n-\vartheta_n^{hk}\|_E +\|\xi_n-\xi_n^{hk}\|_Y+\|T_n-T_n^{hk}\|_E+\|S_n-S_n^{hk}\|_E+\|R_n-R_n^{hk}\|_Y\Big\}\end{array}$

are presented in for different values of $h$ and $k$ . Moreover, the evolution of the error depending on the parameter $h+k$ is plotted in Figure 1. We notice that the convergence of the algorithm is clearly observed, and the linear convergence, stated in Theorem 4.3, is achieved.

Table 1. Example 1: Numerical errors for some $h$ and $k$ .

$h\downarrow k\rightarrow$	0.01	0.005	0.002	0.001	0.0005	0.0002	0.0001
$1/2^3$	0.574576	0.570780	0.568601	0.567893	0.567542	0.567332	0.567262
$1/2^4$	0.290868	0.286851	0.284646	0.283947	0.283604	0.283400	0.283332
$1/2^5$	0.149884	0.145325	0.142983	0.142270	0.141926	0.141724	0.141657
$1/2^6$	0.080538	0.074898	0.072236	0.071480	0.071127	0.070923	0.070856
$1/2^7$	0.047632	0.040251	0.036965	0.036114	0.035737	0.035528	0.035460
$1/2^8$	0.033167	0.023814	0.019518	0.018481	0.018056	0.017833	0.017764
$1/2^9$	0.027415	0.016591	0.011111	0.009758	0.009240	0.008989	0.008916
$1/2^{10}$	0.025348	0.013719	0.007310	0.005556	0.004879	0.004576	0.004495
$1/2^{11}$	0.024679	0.012688	0.005748	0.003656	0.002778	0.002384	0.002288
$1/2^{12}$	0.024485	0.012356	0.005169	0.002876	0.001829	0.001316	0.001192
$1/2^{13}$	0.024434	0.012260	0.004976	0.002586	0.001439	0.000821	0.000658

| Show Table

DownLoad: CSV

Figure 1. Example 1: Asymptotic constant error.

DownLoad: Full-Size Img PowerPoint

If we assume now that there are not supply terms and we use the final time $T_f = 30$ , the same data than in the previous example and the initial conditions:

$\begin{array}{l} u^0 = v^0 = \phi^0 = \psi^0 = \theta^0 = \vartheta^0 = \xi^0 = 0,\\ T^0(x) = R^0(x) = S^0(x) = x(x-1)\quad \hbox{for a.e. }x\in (0,1),\end{array}$

taking the discretization parameters $h = k = 0.001$ , the evolution in time of the discrete energy

$\begin{array}{l}E_n^{hk} = ||v_n^{hk}||_Y^2+2||u_n^{hk}||_E^2+\|\psi_n^{hk}\|_Y^2+\|\phi_n^{hk}\|_E^2+\|\xi_n^{hk}\|_Y^2+3 \|\vartheta_n^{hk}\|_E^2+3 \|\theta_n^{hk}\|_E^2\\ \qquad \qquad +2\|R_n^{hk}\|_Y^2+2\|S_n^{hk}\|_E^2+2\|T_n^{hk}\|_E^2 \end{array}$

is plotted in Figure 2 (in both natural and semi-log scales). As can be seen, it converges to zero and an exponential decay seems to be achieved.

Figure 2. Example 1: Evolution in time of the discrete energy (natural and semi-log scales).

DownLoad: Full-Size Img PowerPoint

5.2. Second example: dependence on the coupling parameter $\mu_2$

We consider now the two-dimensional domain $B = (0, 1)\times(0, 1)$ . In this case, our aim is to study the dependence on the coupling parameter $\mu_2$ .

The following data have been used in the simulations of this example:

$\begin{array}{l} \rho = 1,\quad T_f = 1, \quad \lambda = 1,\quad \mu = 2, \quad \mu_0 = 1, \quad J = 1, \quad a_0 = 1,\quad \beta_0 = 1, \quad \beta_1 = 1,\\ \kappa = 3, \quad a = 1, \quad b = 2,\quad \kappa_1 = 1, \quad \kappa_2 = 1, \quad \kappa_4 = 1,\quad \kappa_5 = 1,\quad \kappa_6 = 2, \quad \xi = 2, \\ \tau_1 = 1, \quad \tau_2 = 2, \end{array}$

with the initial conditions:

$\mathit{\boldsymbol{u}}^0 = \mathit{\boldsymbol{v}}^0 = \mathit{\boldsymbol{T}}^0 = \mathit{\boldsymbol{R}}^0 = \mathit{\boldsymbol{S}}^0 = \mathit{\boldsymbol{0}},\quad \phi^0 = \psi^0 = \theta^0 = \vartheta^0 = \xi^0 = 0,$

and null Dirichlet boundary conditions. We also add a supply term in the equation for function $\mathit{\boldsymbol{u}}$ given by $\mathit{\boldsymbol{f}}(x, y, t) = (tx(x-1)y(y-1), tx(x-1)y(y-1))$ for $(x, y)\in (0, 1)\times (0, 1)$ and $t\in [0, 1]$ .

Using the time discretization parameter $k = 0.01$ and the uniform finite element mesh obtained dividing the unit square into 2048 triangles, in we plot the norm of function $\mathit{\boldsymbol{u}}$ and the microtemperatures, at middle point $\mathit{\boldsymbol{x}} = (0.5, 0.5)$ , for some values of parameter $\mu_2$ . It seems that, for the highest value $\mu_2 = 20$ , an oscillation is produced.

Figure 3. Example 2: Evolution in time of the norm of $\mathit{\boldsymbol{u}}$ and norm of the microtemperatures at the middle point for different values of parameter $\mu_2$ .

DownLoad: Full-Size Img PowerPoint

Finally, in we plot function $\phi$ and the temperature $\theta$ at final time for the value of parameter $\mu_2 = 20$ . We can note that they have been generated by the resulting deformation.

Figure 4. Example 2: Function $\phi$ and temperature $\theta$ at final time for parameter $\mu_2 = 20$ .

DownLoad: Full-Size Img PowerPoint

Acknowledgments

The work of J.R. Fernández has been partially supported by Ministerio de Ciencia, Innovación y Universidades under the research project PGC2018-096696-B-I00 (FEDER, UE).

The work of R. Quintanilla has been supported by Ministerio de Ciencia, Innovación y Universidades under the research project "Análisis matemático aplicado a la termomecánica" (PID2019-105118GB-I00).

Conflict of interest

The authors declare there is no conflict of interest.

References

[1]	A. C. Eringen, Microcontinuum Field Theories I. Foundations and Solids, Springer-Verlag, New York, 1999. https://doi.org/10.1007/978-1-4612-0555-5
[2]	D. Ieşan, Thermoelastic Models of Continua, Kluwer Academic, Dordrecht, 2004. https://doi.org/10.1007/978-1-4020-2310-1_5
[3]	M. A. Goodman, S. C. Cowin, A continuum theory for granular materials, Arch. Rat. Mech. Anal., 44 (1972), 249–266. https://doi.org/10.1007/BF00284326 doi: 10.1007/BF00284326
[4]	S. C. Cowin, J. W. Nunziato, Linear elastic materials with voids, J. Elast., 13 (1983), 125–147. https://doi.org/10.1007/BF00041230 doi: 10.1007/BF00041230
[5]	S. C. Cowin, The viscoelastic behavior of linear elastic materials with voids, J. Elast., 15 (1985), 185–191. https://doi.org/10.1007/BF00041992 doi: 10.1007/BF00041992
[6]	J. W. Nunziato, S. C. Cowin, A nonlinear theory of elastic materials with voids, Arch. Ration. Mech. Anal., 72 (1979), 175–201. https://doi.org/10.1007/BF00249363 doi: 10.1007/BF00249363
[7]	D. Ieşan, A theory of thermoelastic materials with voids, Acta Mech., 60 (1986), 67–89. https://doi.org/10.1007/BF01302942 doi: 10.1007/BF01302942
[8]	T. A. Apalara, Exponential decay in one-dimensional porous elasticity, Quart. J. Mech. Appl. Math., 70 (2017), 363–372. https://doi.org/10.1093/qjmam/hbx012 doi: 10.1093/qjmam/hbx012
[9]	T. A. Apalara, General decay of solutions in one-dimensional porous-elastic system with memory, J. Math. Anal. Appl., 469 (2019), 457–471 https://doi.org/10.1016/j.jmaa.2017.08.007 doi: 10.1016/j.jmaa.2017.08.007
[10]	N. Bazarra, J. R. Fernández, R. Quintanilla, An a priori error analysis of poro-thermoviscoelastic problems, Appl. Math. Comput., 379 (2020), 125268. https://doi.org/10.1016/j.amc.2020.125268 doi: 10.1016/j.amc.2020.125268
[11]	B. Feng, Uniform decay of energy for porous thermoelastic systems with past history, Appl. Anal., 97 (2018), 210–229. https://doi.org/10.1080/00036811.2016.1258116 doi: 10.1080/00036811.2016.1258116
[12]	B. Feng, T. A. Apalara, Optimal decay for a porous elasticity system with memory, J. Math. Anal. Appl., 470 (2019), 1108–1128. https://doi.org/10.1016/j.jmaa.2018.10.052 doi: 10.1016/j.jmaa.2018.10.052
[13]	B. Feng, L. Yan, D. S. Almeida Junior, Stabilization for an inhomogeneous porous-elastic system with temperatures and microtemperatures, Z. Angew. Math. Mech., 101 (2021), e202000058. https://doi.org/10.1002/zamm.202000058 doi: 10.1002/zamm.202000058
[14]	B. Feng, M. Yin, Decay of solutions for one-dimensional porous elasticity system with memory: the case of nonequal waves speed, Math. Mech. Solids, 24 (2019), 2361–2373. https://doi.org/10.1177/1081286518757299 doi: 10.1177/1081286518757299
[15]	J. R. Fernández, M. Masid, Analysis of a problem arising in porous thermoelasticity of type II, J. Therm. Stress., 39 (2016), 513–531. https://doi.org/10.1080/01495739.2016.1158609 doi: 10.1080/01495739.2016.1158609
[16]	J. R. Fernández, M. Masid, Numerical analysis of a thermoelastic diffusion problem with voids, Int. J. Numer. Anal. Model., 14 (2017), 153–174.
[17]	M. C. Leseduarte, A. Magaña, R. Quintanilla, On the time decay of solutions in porous-thermo-elasticity of type II, Discrete Contin. Dyn. Syst. Ser. B, 13 (2010), 375–391. https://doi.org/10.3934/dcdsb.2010.13.375 doi: 10.3934/dcdsb.2010.13.375
[18]	A. Magaña, R. Quintanilla, On the exponential decay of solutions in one-dimensional generalized porous-thermoelasticity, Asymptot. Anal., 49 (2006), 173–187.
[19]	A. Miranville, R. Quintanilla, Exponential stability in type III thermoelasticity with voids, Appl. Math. Lett., 94 (2019), 30–37. https://doi.org/10.1016/j.aml.2019.02.014 doi: 10.1016/j.aml.2019.02.014
[20]	D. Ieşan, R. Quintanilla, On a theory of thermoelasticity with microtemperatures, J. Therm. Stress., 23 (2000), 195–215. https://doi.org/10.1080/014957300280407 doi: 10.1080/014957300280407
[21]	D. Ieşan, R. Quintanilla, On thermoelastic bodies with inner structure and microtemperatures, J. Math. Anal. Appl., 354 (2009), 12–23. https://doi.org/10.1016/j.jmaa.2008.12.017 doi: 10.1016/j.jmaa.2008.12.017
[22]	C. Cattaneo, On a form of heat equation which eliminates the paradox of instantaneous propagation, C. R. Acad. Sci. Paris, 247 (1958), 431–433.
[23]	D. Y. Tzou, A unified approach for heat conduction from macro to micro-scales, ASME J. Heat Transfer, 117 (1995), 8–16. https://doi.org/10.1115/1.2822329 doi: 10.1115/1.2822329
[24]	Z. Liu, R. Quintanilla, Time decay in dual-phase-lag thermoelasticity: critical case, Comm. Pure Appl. Anal., 17 (2018), 177–190. https://doi.org/10.3934/cpaa.2018011 doi: 10.3934/cpaa.2018011
[25]	F. Maes, K. Van Bockstal, Thermoelastic problem in the setting of dual-phase-lag heat conduction: Existence and uniqueness of a weak solution, J. Math. Anal. Appl., 503 (2021), 125304. https://doi.org/10.1016/j.jmaa.2021.125304 doi: 10.1016/j.jmaa.2021.125304
[26]	R. Quintanilla, R. Racke, Qualitative aspects in dual-phase-lag thermoelasticity, SIAM J. Appl. Math., 66 (2006), 977–1001. https://doi.org/10.1137/05062860X doi: 10.1137/05062860X
[27]	N. Bazarra, M. I. M. Copetti, J. R. Fernández, R.Quintanilla, Numerical analysis of a dual-phase-lag model with microtemperatures, Appl. Num. Math., 166 (2021), 1–25. https://doi.org/10.1016/j.apnum.2021.03.016 doi: 10.1016/j.apnum.2021.03.016
[28]	Z. Liu, R. Quintanilla, Y. Wang, Dual-phase-lag heat conduction with microtemperatures, Z. Angew. Math. Mech., 101 (2021), e202000167. https://doi.org/10.1002/zamm.202000167 doi: 10.1002/zamm.202000167
[29]	P. G. Ciarlet, Basic error estimates for elliptic problems. In: Handbook of Numerical Analysis, P.G. Ciarlet and J.L. Lions eds., vol II (1993), 17–351. https://doi.org/10.1016/S1570-8659(05)80039-0
[30]	Ph. Clement, Approximation by finite element functions using local regularization, RAIRO Math. Model. Numer. Anal., 9 (1975), 77–84. https://doi.org/10.1051/m2an/197509R200771 doi: 10.1051/m2an/197509R200771
[31]	M. Campo, J. R. Fernández, K. L. Kuttler, M. Shillor and J. M. Viaño, Numerical analysis and simulations of a dynamic frictionless contact problem with damage, Comput. Methods Appl. Mech. Engrg., 196 (2006), 476–488. https://doi.org/10.1016/j.cma.2006.05.006 doi: 10.1016/j.cma.2006.05.006

This article has been cited by:

1.	Shengda Zeng, Moncef Aouadi, Temporal behavior in nonlocal strain gradient porous thermoelasticity with relaxed microtemperatures, 2025, 0149-5739, 1, 10.1080/01495739.2025.2458107
2.	Pooja Rani, Rajneesh Kumar, Geeta Partap, Wave Analysis and Deformation in Poro-thermoelastic Medium with Microtemperature Under Nonlocal Elasticity Theory and Dual Phase Lag Model, 2025, 2228-6187, 10.1007/s40997-025-00848-4

Reader Comments

Your name:*

Email:*
© 2022 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)

通讯作者: 陈斌, bchen63@163.com

1.
沈阳化工大学材料科学与工程学院沈阳 110142

Electronic Research Archive

1 1.3

Metrics

Article views(1667) PDF downloads(72) Cited by(2)

Preview PDF

Download XML

Export Citation

Article outline

Show full outline

Figures and Tables

Figures(4) / Tables(1)

Electronic Research Archive

A dual-phase-lag porous-thermoelastic problem with microtemperatures

Related Papers:

Abstract

1. Introduction

2. The thermomechanical model

3. Existence and uniqueness

4. Fully discrete approximations: a priori error analysis

5. Numerical results

5.1. First example: numerical convergence in a one-dimensional example

5.2. Second example: dependence on the coupling parameter $\mu_2$

Acknowledgments

Conflict of interest

References

This article has been cited by:

Reader Comments

通讯作者: 陈斌, bchen63@163.com

Metrics

Figures and Tables

Other Articles By Authors

Catalog

Electronic Research Archive

A dual-phase-lag porous-thermoelastic problem with microtemperatures

Related Papers:

Abstract

1. Introduction

2. The thermomechanical model

3. Existence and uniqueness

4. Fully discrete approximations: a priori error analysis

5. Numerical results

5.1. First example: numerical convergence in a one-dimensional example

5.2. Second example: dependence on the coupling parameter \mu_2

Acknowledgments

Conflict of interest

References

This article has been cited by:

Reader Comments

通讯作者: 陈斌, bchen63@163.com

Metrics

Figures and Tables

Other Articles By Authors

Related pages

Tools

Export File

Citation

Format

Content

Catalog

5.2. Second example: dependence on the coupling parameter $\mu_2$