Appendix A Proof of Theorem 4.1
Thanks to the Scott–Dupont theory (see [15]) we know that for every s with \(0\le s\le k\) and for every \(u\in H^{1+s}(K)\), there exists \(u_\pi \in {\mathbb {P}}_k(K)\), \(k\ge 2\), such that
$$\begin{aligned} \Vert u - u_{\pi } \Vert _{0,K} + h_K |u -u_{\pi } |_{1,K} \lesssim h_K^{1+s} |u|_{1+s,K} \quad \text {for all }K \in {\mathcal {T}}_h. \end{aligned}$$
(A.1)
We can then write the displacement and total pressure error in terms of the poroelastic projector
$$\begin{aligned} ({\varvec{u}}- {\varvec{u}}_h)(t)&= ({\varvec{u}}- I^h_{{\varvec{u}}}{\varvec{u}})(t) + (I^h_{{\varvec{u}}} {\varvec{u}}- {\varvec{u}}_h)(t) := e_{{\varvec{u}}}^I(t) + e_{{\varvec{u}}}^A(t), \\ (\psi - \psi _h)(t)&= (\psi - I^h_{\psi } \psi )(t) + (I^h_{\psi } \psi - \psi _h)(t) := e_{\psi }^I(t) + e_{\psi }^A(t). \end{aligned}$$
Then, a combination of equations (4.1a), (3.4a) and (2.4a) gives
$$\begin{aligned} a_1^h(e_{{\varvec{u}}}^A, \varvec{v}_h)+ b_1(\varvec{v}_h, e_{\psi }^A)&= (a_1({\varvec{u}},\varvec{v}_h) -a_1^h({\varvec{u}}_h, \varvec{v}_h)) + b_1(\varvec{v}_h, \psi - \psi _h) = (F-F^h)(\varvec{v}_h), \end{aligned}$$
and taking as test function \(\varvec{v}_h = \partial _t e_{{\varvec{u}}}^A\), we can write the relation
$$\begin{aligned} a_1^h(e_{{\varvec{u}}}^A, \partial _t e_{{\varvec{u}}}^A) + b_1(\partial _t e_{{\varvec{u}}}^A, e_{\psi }^A) = (F-F^h)(\partial _t e_{{\varvec{u}}}^A). \end{aligned}$$
(A.2)
Now, we write the pressure error in terms of the poroelastic projector as follows
$$\begin{aligned} (p^{\textrm{P}} - p_h^{\textrm{P}})(t) = (p^{\textrm{P}} - I^h_p p^{\textrm{P}})(t) + (I^h_p p^{\textrm{P}} - p_h^{\textrm{P}})(t) := e_{p}^{I,\textrm{P}}(t) + e_{p}^{A, \textrm{P}}(t). \end{aligned}$$
Using (4.1c), (3.4b) and (2.4b), we obtain
$$\begin{aligned}&{\tilde{a}}_2^h(\partial _t e_p^{A,\textrm{P}}, q_h^{\textrm{P}}) + a_2^h(e_p^{A, \textrm{P}}, q_h^{\textrm{P}}) - b_2(q_h^{\textrm{P}}, \partial _t e_{\psi }^{A} ) \\&\quad = ({\tilde{a}}_2^h(\partial _t I^h_p p^{\textrm{P}}, q_h^{\textrm{P}}) - {\tilde{a}}_2(\partial _t p^{\textrm{P}}, q_h^{\textrm{P}})) + b_2(q_h^{\textrm{P}}, \partial _t e_{\psi }^I) + (G- G^h)(q_h^{\textrm{P}}). \end{aligned}$$
We can take \(q_h^{\textrm{P}} = e_p^{A, \textrm{P}}\), which leads to
$$\begin{aligned} \begin{aligned}&{\tilde{a}}_2^h(\partial _t e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) +a_2^h(e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) - b_2(e_p^{A, \textrm{P}}, \partial _t e_{\psi }^A) \\&\quad = ({\tilde{a}}_2^h(\partial _t I^h_p p^{\textrm{P}}, e_p^{A, \textrm{P}}) - {\tilde{a}}_2(\partial _t p^{\textrm{P}}, e_p^{A, \textrm{P}}) ) + b_2(e_p^{A, \textrm{P}}, \partial _t e_{\psi }^I) + (G- G^h)(e_p^{A, \textrm{P}}). \end{aligned} \end{aligned}$$
(A.3)
Next we use (4.1b), (3.4c) and (2.4c), and this implies
$$\begin{aligned}&b_1(e_{{\varvec{u}}}^A, \phi _h)+ b_2(e_p^{A, \textrm{P}} , \phi _h) - a_3(e_{\psi }^A, \phi _h) = -b_2(e_p^{I, \textrm{P}}, \phi _h) + a_3(e_{\psi }^I, \phi _h). \end{aligned}$$
Differentiating the above equation with respect to time and taking \(\phi _h = -e_{\psi }^A\), we can assert that
$$\begin{aligned} - b_1(\partial _t e_{{\varvec{u}}}^A, e_{\psi }^A) - b_2(\partial _t e_p^{A, \textrm{P}}, e_{\psi }^A) + a_3(\partial _t e_{\psi }^A, e_{\psi }^A) = b_2(\partial _t e_p^I, e_{\psi }^A) -a_3(\partial _t e_{\psi }^I, e_{\psi }^A). \end{aligned}$$
(A.4)
Then we simply add (A.2), (A.3) and (A.4), to obtain
$$\begin{aligned} \begin{aligned}&a_1^h(e_{{\varvec{u}}}^A, \partial _t e_{{\varvec{u}}}^A) + {\tilde{a}}_2^h(\partial _t e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) +a_2^h(e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) + a_3(\partial _t e_{\psi }^A, e_{\psi }^A) \\&\qquad - b_2(e_p^{A, \textrm{P}}, \partial _t e_{\psi }^A) - b_2(\partial _t e_p^{A, \textrm{P}}, e_{\psi }^A)\\&\quad = (F- F^h)(\partial _{t} e_{{\varvec{u}}}^A ) + ({\tilde{a}}_2^h(\partial _t I^h_p p^{ \textrm{P}}, e_p^{A, \textrm{P}}) - {\tilde{a}}_2(\partial _{t} p^{\textrm{P}}, e_p^{A, \textrm{P}}))\\&\qquad + b_2(e_p^{A, \textrm{P}}, \partial _{t} e_{\psi }^I) + (G- G^h) (e_p^{A, \textrm{P}}) + b_2(\partial _{t} e_p^{I, \textrm{P}}, e_{\psi }^A) - a_3(\partial _{t} e_{\psi }^I, e_{\psi }^A). \end{aligned} \end{aligned}$$
(A.5)
Regarding the left-hand side of (A.5), repeating arguments to obtain alike to the stability proof, that is,
$$\begin{aligned}&a_1^h(e_{{\varvec{u}}}^A, \partial _t e_{{\varvec{u}}}^A) + {\tilde{a}}_2^h(\partial _t e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) +a_2^h(e_p^{A, \textrm{P}}, e_p^{A, \textrm{P}}) + a_3(\partial _t e_{\psi }^A, e_{\psi }^A)\\&\quad - b_2(e_p^{A, \textrm{P}}, \partial _t e_{\psi }^A) - b_2(\partial _t e_p^{A, \textrm{P}}, e_{\psi }^A) \\&\quad \gtrsim {\mu ^{\min }}\frac{\textrm{d}}{\textrm{d}t} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A)\Vert _0^2 + c_0\frac{\textrm{d}}{\textrm{d}t} \Vert e_p^{A, \textrm{P}}\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{2 \kappa _{\min }}{\eta } \Vert \nabla e_p^A\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{1}{{ \lambda ^{\textrm{E}}}} \frac{\textrm{d}}{\textrm{d}t} \Vert e_{\psi }^{A, \textrm{E}}\Vert _{0, \varOmega ^{\textrm{E}}}^2 \\&\quad + \frac{1}{{ \lambda ^{\textrm{P}}}} \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \biggl ( \alpha ^2\frac{\textrm{d}}{\textrm{d}t}\Vert (I-\varPi ^0_K) e_p^{A,\textrm{P}}\Vert _{0,K}^2 +\frac{\textrm{d}}{\textrm{d}t} \Vert \alpha \varPi ^0_K e_p^{A,\textrm{P}} -e_{\psi }^{A,\textrm{P}}\Vert _{0,K}^2 \biggr ) . \end{aligned}$$
Then integrating equation (A.5) in time and using the consistency of \({\tilde{a}}_2(\cdot , \cdot )\), implies the bound
$$\begin{aligned}&{\mu _{\min }} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(t))\Vert _0^2 + c_0 \Vert e_p^{A, \textrm{P}} (t)\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{1}{{ \lambda ^{\textrm{E}}}} \Vert e_{\psi }^{A, \textrm{E}}(t)\Vert _{0, \varOmega ^{\textrm{E}}}^2 + \frac{\kappa _{\min }}{\eta }\int _0^t \Vert \nabla e_p^{A, \textrm{P}}(s)\Vert _{0, \varOmega ^{\textrm{P}}}^2 \, \textrm{d}s \\&\qquad + \frac{1}{{ \lambda ^{\textrm{P}}}} \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \Bigl (\alpha ^2 \Vert (I-\varPi ^{0,k}_K) e_p^{A, \textrm{P}}(t)\Vert _{0,K}^2 + \Vert (\alpha \varPi ^{0,k}_K e_p^{A, \textrm{P}} -e_{\psi }^{A, \textrm{P}})(t)\Vert _{0,K}^2 \Bigr ) \\&\quad \lesssim {\mu _{\min }} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(0))\Vert _0^2 + c_0 \Vert e_p^{A, \textrm{P}} (0)\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{1}{{ \lambda ^{\textrm{E}}}} \Vert e_{\psi }^{A, \textrm{E}}(0)\Vert _{0, \varOmega ^{\textrm{E}}}^2\\&\qquad + \frac{1}{{ \lambda ^{\textrm{P}}}} \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \Bigl (\alpha ^2 \Vert (I-\varPi ^{0,k}_K) e_p^{A, \textrm{P}}(0)\Vert _{0,K}^2 + \Vert (\alpha \varPi ^{0,k}_K e_p^{A, \textrm{P}} -e_{\psi }^A)(0)\Vert _{0,K}^2 \Bigr ) \\&\qquad + \underbrace{ \int _0^t \bigl ((\varvec{b}- \varvec{b}_h)(s),\partial _{t} e_{{\varvec{u}}}^A(s) \bigr )_{0, \varOmega } \, \textrm{d}s}_{=:D_1} + \underbrace{\int _0^t \bigl ((\ell ^{\textrm{P}} - \ell _h^{ \textrm{P}})(s), e_p^{A, \textrm{P}}(s) \bigr )_{0, \varOmega ^{\textrm{P}}}\, \textrm{d}s}_{=:D_2} \\&\qquad + \underbrace{\int _0^t \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \Bigl ({\tilde{a}}_2^{h,K} \bigl (\partial _t(I^h_p p^{\textrm{P}} - p_{\pi }^{\textrm{P}})(s), e_p^{A, \textrm{P}}(s) \bigr ) - {\tilde{a}}^{K}_2 \bigl (\partial _t(p^{\textrm{P}} - p_{\pi }^{\textrm{P}})(s), e_p^{A, \textrm{P}}(s) \bigr ) \Bigr ) \, \textrm{d}s}_{=:D_3} \\&\qquad + \underbrace{\int _0^t \Bigl (b_2 \bigl (e_p^{A, \textrm{P}}(s), \partial _{t} e_{\psi }^I(s) \bigr ) + b_2 \bigl (\partial _{t} e_p^{I, \textrm{P}} (s), e_{\psi }^A (s)\bigr ) - a_3 \bigl (\partial _{t} e_{\psi }^I(s), e_{\psi }^A(s) \bigr ) \Bigr ) \, \textrm{d}s}_{=:D_4}. \end{aligned}$$
Then we integrate by parts in time and invoke the orthogonality property of the \(L^2\)-projection, which yields
$$\begin{aligned} D_1&= \sum _{K \in \mathcal {T}_h} \Bigl ( \bigl ((\textbf{I}- \varvec{\varPi }^{0,k-2}_K) \varvec{b}(t), (\textbf{I}- \varvec{\varPi }^{0,0}_K) e_{{\varvec{u}}}^A(t) \bigr )_{0, K}\\&\quad - \bigl ((\textbf{I}- \varvec{\varPi }^{0,k-2}_K) \varvec{b}(0), (\textbf{I}- \varvec{\varPi }^{0,0}_K) e_{{\varvec{u}}}^A(0) \bigr )_{0, K} \\&\quad - \int _0^t \bigl (\partial _{t} (\textbf{I}- \varvec{\varPi }^{0,k-2}_K) \varvec{b}(s), (\textbf{I}- \varvec{\varPi }^{0,0}_K) e_{{\varvec{u}}}^A(s) \bigr )_{0, K} \, \textrm{d}s \Bigr ), \end{aligned}$$
and use Cauchy–Schwarz and Young’s inequalities along with Korn’s inequality to arrive at
$$\begin{aligned} D_1&\le \frac{{\mu _{\min }} }{2} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(t))\Vert _0^2 + {C(\mu _{\min })} h^k \biggl (h^k |\varvec{b}(t)|_{k-1}^2 + |\varvec{b}(0)|_{k-1} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(0))\Vert _0 \\&\quad + \int _0^t | \partial _{t} \varvec{b}(s)|_{k-1} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(s))\Vert _{0}\, \textrm{d}s \biggr ), \end{aligned}$$
where we have used standard error estimate for the \(L^2\)-projection \(\varvec{\varPi }_K^{0,k}\) onto piecewise constant functions. Using again Cauchy–Schwarz inequality, standard error estimates for \(\varPi _K^{0,k}\) on the term \(D_2\), Young’s and Poincaré inequalities readily gives
$$\begin{aligned}&D_2 {\le C} h^k \int _0^t |\ell ^{\textrm{P}} (s)|_{k-1, \varOmega ^{\textrm{P}}} \Vert \nabla e_p^{A, \textrm{P}}(s)\Vert _{0, \varOmega ^{\textrm{P}}}\, \,\textrm{d}s\le {C(\kappa _{\min },\eta )} h^{2k} \int _{0}^t |\ell ^{\textrm{P}} (s)|_{k-1, \varOmega ^{\textrm{P}}}^2 \,\textrm{d}s\\&\quad + \frac{\kappa _{\min }}{{4} \eta } \int _{0}^t \Vert \nabla e_p^{A, \textrm{P}}(s)\Vert _{0, \varOmega ^{\textrm{P}}}^2 \,\textrm{d}s. \end{aligned}$$
On the other hand, considering the polynomial approximation \(p_{\pi }^{\textrm{P}}\) (cf. (A.1)) of \(p^{\textrm{P}}\), utilising the triangle inequality, Young’s and Poincaré inequalities yield
$$\begin{aligned} D_3&{\le C(\kappa _{\min },\eta ) } h^{2(k+1)} \bigg (c_0 + \frac{\alpha ^2}{{\lambda ^{\textrm{P}}}}\bigg )^2 \int _0^t |\partial _{t} p^{\textrm{P}}(s)|_{k+1, \varOmega ^{\textrm{P}}}^2 \,\textrm{d}s+ \frac{\kappa _{\min }}{ {4} \eta } \int _{0}^t \Vert \nabla e_p^{A, \textrm{P}}(s)\Vert _{0, \varOmega ^{\textrm{P}}}^2 \,\textrm{d}s. \end{aligned}$$
Also,
$$\begin{aligned} D_4&\le \frac{{C}}{\lambda } h^k \int _0^t \Bigl ( \alpha ( | \partial _{t} \psi ^{\textrm{P}} (s)|_{k, \varOmega ^{\textrm{P}}} + |\partial _t {\varvec{u}}^{\textrm{P}}(s)|_{k+1,\varOmega ^{\textrm{P}}} + |\partial _t {\varvec{u}}^{\textrm{E}}(s)|_{k+1, \varOmega ^{\textrm{E}}} )\Vert e_p^{A, \textrm{P}} (s)\Vert _{0, \varOmega ^{\textrm{P}}} \\&\quad + (\alpha h |\partial _{t} p^{\textrm{P}}(s) |_{k+1, \varOmega ^{\textrm{P}}} + | \partial _{t} \psi ^{\textrm{P}} (s)|_{k, \varOmega ^{\textrm{P}}} + | \partial _{t} \psi ^{\textrm{E}} (s)|_{k, \varOmega ^{\textrm{E}}} + |\partial _t {\varvec{u}}(s)|_{k+1}) \Vert e_{\psi }^A (s)\Vert _0 \Bigr ) \, \textrm{d}s. \end{aligned}$$
Using the discrete inf-sup condition (cf. (3.3)) and a combination of equations (4.1a), (3.4a) and (2.4a), we get
$$\begin{aligned} {{\tilde{\beta }}} \Vert e_{\psi }^A(t) \Vert _0&\le \sup _{\varvec{v}_h \in \textbf{V}_h} \frac{b_1(\varvec{v}_h, e_{\psi }^A(t))}{\Vert \varvec{v}_h \Vert _1} { ~\le C \big (} h^k {\rho } | \varvec{b}(t) |_{k-1} + {\mu _{\max }} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(t))\Vert _0 {\big )}. \end{aligned}$$
(A.6)
Then, with the help of Young’s and Poincaré inequalities, the bound of \(D_4\) becomes
$$\begin{aligned}&D_4 {\le } \frac{{C}}{\lambda } h^k \int _0^t \Bigl ( (\alpha {h} |\partial _{t} p^{\text {P}}(s) |_{k+1, \varOmega ^{\text {P}}} + | \partial _{t} \psi ^{\text {P}} (s)|_{k, \varOmega ^{\text {P}}} + | \partial _{t} \psi ^{\text {E}} (s)|_{k, \varOmega ^{\text {E}}} + |\partial _t {\textbf{u}}^{\text {P}}(s)|_{k+1, \varOmega ^{\text {P}}} \\ {}&\qquad \qquad \quad + |\partial _t {\textbf{u}}^{\text {E}}(s)|_{k+1, \varOmega ^{\text {E}}} ) \times ( h^k {\rho |\textbf{b}(s)|_{k-1}} + {\mu _{\max }} \Vert \mathbf {\varepsilon }(e_{{\textbf{u}}}^A(s))\Vert _0) \\ {}&\qquad \qquad \quad + \alpha \Vert e_p^{A, \text {P}} (s)\Vert _{0, \varOmega ^{\text {P}}} (| \partial _{t} \psi ^{\text {P}} (s)|_{k, \varOmega ^{\text {P}}} + |\partial _t {\textbf{u}}(s)|_{k+1} ) \Bigr ) \,\text {d}s. \end{aligned}$$
Combining the bounds of all \(D_i, i=1,2,3,4\) and proceeding in a similar fashion as for the bounds in the stability proof in [18] (using Lemma 3.1 and (3.6)), we can eventually conclude that
$$\begin{aligned}&{\mu _{\min }} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(t))\Vert _0^2 + c_0 \Vert e_p^{A,\textrm{P}} (t)\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{1}{{\lambda ^{\textrm{E}}}} \Vert e_{\psi }^{A, \textrm{E}}(t)\Vert _{0, \varOmega ^{\textrm{E}}}^2 + \frac{\kappa _{\min }}{\eta } \int _0^t \Vert \nabla e_p^{A, \textrm{P}} (s)\Vert _{0, \varOmega ^{\textrm{P}}}^2 \, \textrm{d} s \\&\quad \lesssim {\mu _{\min }} \Vert \varvec{\varepsilon }(e_{{\varvec{u}}}^A(0))\Vert _0^2 + \Big (c_0 + \frac{\alpha ^2}{{\lambda ^{\textrm{P}}}} \Big ) \Vert e_p^{A, \textrm{P}} (0)\Vert _{0, \varOmega ^{\textrm{P}}}^2 + \frac{1}{{\lambda ^{\textrm{E}}}} \Vert e_{\psi }^{A,\textrm{E}}(0)\Vert _{0, \varOmega ^{\textrm{E}}}^2 \\&\qquad + {C(\mu _{\min }, \kappa _{\min }, \eta ) h^{2k} \Biggl ( \sup _{t \in [0,t_{\text {final}}] } |\varvec{b}(t)|_{k-1}^2 \!+\! \int _0^t \Bigl ( h^2 \big ( c_0 \!+\! \frac{\alpha ^2}{{\lambda ^{\textrm{P}}}}\big )^2 |\partial _{t} p^{\textrm{P}}(s) |_{k+1, \varOmega ^{\textrm{P}}}^2 \!+\! |\ell ^{\textrm{P}} (s)|_{k\!-\!1, \varOmega ^{\textrm{P}}}^2 \Bigr ) \, \textrm{d}s } \\&\qquad { + \biggl ( \int _0^t \Bigl ( |\varvec{b}(s)|_{k-1} + | \partial _{t}\varvec{b}(s)|_{k-1} + h \Big ( c_0 + \frac{\alpha ^2}{{\lambda ^{\textrm{P}}}}\Big ) |\partial _{t} p^{\textrm{P}}(s) |_{k+1, \varOmega ^{\textrm{P}}} } \\&\qquad {+ \frac{1}{{\lambda _{\min }}} \big ( | \partial _{t} \psi ^{\textrm{P}} (s)|_{k, \varOmega ^{\textrm{P}}} +| \partial _{t} \psi ^{\textrm{E}} (s)|_{k, \varOmega ^{\textrm{E}}} +|\partial _t {\varvec{u}}^{\textrm{P}}(s)|_{k+1, \varOmega ^{\textrm{P}}} +|\partial _t {\varvec{u}}^{\textrm{E}}(s)|_{k+1, \varOmega ^{\textrm{E}}} \big ) \Bigr ) \, \textrm{d}s \biggr )^2\Biggr ).} \end{aligned}$$
Then choosing \({\varvec{u}}_h(0): ={\varvec{u}}_I(0)\), \(\psi _h(0): = \varPi ^{0,k-1}\psi (0)\), \(p_h^{\textrm{P}}(0): = p_I^{\textrm{P}}(0)\) and applying the triangle inequality together with (A.6), completes the rest of the proof.
Appendix B: Proof of Theorem 4.2
As in the semidiscrete case we split the individual errors as
$$\begin{aligned} {\varvec{u}}(t_n) - {\varvec{u}}_h^n&= ({\varvec{u}}(t_n) - I^h_{{\varvec{u}}} {\varvec{u}}(t_n)) + (I^h_{{\varvec{u}}} {\varvec{u}}(t_n)- {\varvec{u}}_h^n)=: E_{{\varvec{u}}}^{I,n} + E_{{\varvec{u}}}^{A,n}, \\ \psi (t_n) - \psi _h^n&= (\psi (t_n) - I^h_{\psi } \psi (t_n)) + (I^h_{\psi } \psi (t_n)- \psi _h^n)=: E_{\psi }^{I,n} + E_{\psi }^{A,n}, \\ p^{\textrm{P}}(t_n) - p_h^{n,\textrm{P}}&= (p^{\textrm{P}}(t_n) - I^h_{p} p^{\textrm{P}}(t_n)) + (I^h_{p} p^{\textrm{P}}(t_n)- p_h^{n,\textrm{P}})=: E_{p}^{I,n} + E_{p}^{A,n}, \end{aligned}$$
where the error terms are \(E_{p}^{I,n}:=E_{p}^{I,n} |_{\varOmega ^{\textrm{P}}}, E_{p}^{A,n}:= E_{p}^{A,n} |_{\varOmega ^{\textrm{P}}}\). Then, from estimate (4.2a) and following the steps of the proof of Theorem 4.1 we get the bounds
$$\begin{aligned} \Vert E_{{\varvec{u}}}^{I,n} \Vert _1&\lesssim h^k ( | {\varvec{u}}(0) |_{k+1} + | \psi ^{\textrm{P}}(0) |_{k,\varOmega ^{\textrm{P}}} + | \psi ^{\textrm{E}}(0) |_{k,\varOmega ^{\textrm{E}}} \nonumber \\&\quad + \Vert \partial _t{\varvec{u}}\Vert _{\textbf{L}^1(0,t_n; [H^{k+1}(\varOmega )]^2)} + \Vert \partial _t \psi \Vert _{L^1(0,t_n; k)} ), \end{aligned}$$
(B.1a)
$$\begin{aligned} \Vert E_{\psi }^{I,n} \Vert _0&\lesssim h^k ( | {\varvec{u}}(0) |_{k+1} + | \psi ^{\textrm{P}}(0) |_{k,\varOmega ^{\textrm{P}}} + | \psi ^{\textrm{E}}(0) |_{k,\varOmega ^{\textrm{E}}} \nonumber \\&\quad + \Vert \partial _t{\varvec{u}}\Vert _{\textbf{L}^1(0,t_n; [H^{k+1}(\varOmega )]^2)} + \Vert \partial _t\psi \Vert _{L^1(0,t_n; k)} ), \end{aligned}$$
(B.1b)
$$\begin{aligned} \Vert E_{p}^{I,n} \Vert _{1, \varOmega ^{\textrm{P}}}&\lesssim h^k ( | p^{\textrm{P}}(0) |_{k+1, \varOmega ^{\textrm{P}}} + \Vert \partial _t p^{\textrm{P}}\Vert _{L^1(0,t_n; H^{k+1}(\varOmega ^{\textrm{P}}))}), \end{aligned}$$
(B.1c)
where \(\Vert \partial _t \psi \Vert _{L^1(0,t_n;k)}:= \Vert \partial _t \psi ^{\textrm{P}}\Vert _{L^1(0,t_n; H^k(\varOmega ^{\textrm{P}}))} + \Vert \partial _t \psi ^{\textrm{E}}\Vert _{L^1(0,t_n; H^k(\varOmega ^{\textrm{E}}))}\). From Eqs. (4.1a), (3.10a) and (2.4a), we readily get
$$\begin{aligned} a_1^h(E_{{\varvec{u}}}^{A,n},\varvec{v}_h) + b_1(\varvec{v}_h, E_{\psi }^{A,n}) = F^n(\varvec{v}_h) - F^{h,n}(\varvec{v}_h). \end{aligned}$$
(B.2)
We then use (4.1b) and (3.10c), and proceed to differentiate (2.4c) with respect to time. This implies
$$\begin{aligned} \begin{aligned}&b_1(E_{{\varvec{u}}}^{A,n} - E_{{\varvec{u}}}^{A,n-1}, \phi _h) + b_2(E_p^{A,n}-E_p^{A,n-1}, \phi _h) - a_3(E_{\psi }^{A,n} - E_{\psi }^{A,n-1}, \phi _h) \\&\quad = b_1(({\varvec{u}}(t_n) - {\varvec{u}}(t_{n-1})) - (\varDelta t) \partial _{t} {\varvec{u}}(t_n), \phi _h) + b_2((I_p^h p^{\textrm{P}}(t_n) - I_p^h p^{\textrm{P}}(t_{n-1}))\\&\qquad - (\varDelta t) \partial _{t} p^{\textrm{P}}(t_n), \phi _h) - a_3((I_{\psi }^h \psi (t_n) - I_{\psi }^h \psi (t_{n-1})) - (\varDelta t) \partial _{t} \psi (t_n), \phi _h). \end{aligned} \end{aligned}$$
(B.3)
Choosing \(\varvec{v}_h = E_{{\varvec{u}}}^{A,n} - E_{{\varvec{u}}}^{A,n-1 }\) in (B.2) and \(\phi _h = - E_{\psi }^{A,n}\) in (B.3) and adding the results, gives
$$\begin{aligned}&a_1^h(E_{{\varvec{u}}}^{A,n}, E_{{\varvec{u}}}^{A,n}- E_{{\varvec{u}}}^{A,n-1}) + a_3(E_{\psi }^{A,n} - E_{\psi }^{A,n-1}, E_{\psi }^{A,n}) - b_2(E_{p}^{A,n} - E_{p}^{A,n-1}, E_{\psi }^{A,n}) \nonumber \\&\quad = ( \varvec{b}(t_n) - \varvec{b}^n_h, E_{{\varvec{u}}}^{A,n}- E_{{\varvec{u}}}^{A,n-1 } )_{0, \varOmega } - b_1(({\varvec{u}}(t_n) - {\varvec{u}}(t_{n-1})) - (\varDelta t) \partial _{t} {\varvec{u}}(t_n), E_{\psi }^{A,n}) \nonumber \\&\qquad - b_2((I_p^h p^{\textrm{P}}(t_n) - I_p^h p^{\textrm{P}}(t_{n-1})) - (\varDelta t) \partial _{t} p^{\textrm{P}}(t_n), E_{\psi }^{A,n}) \nonumber \\&\qquad + a_3((I_{\psi }^h \psi (t_n) - I_{\psi }^h \psi (t_{n-1}))- (\varDelta t) \partial _{t} \psi (t_n), E_{\psi }^{A,n}). \end{aligned}$$
(B.4)
Next, and as a consequence of using (4.1c), (3.4b) and (2.4b) with \(q_h^{\textrm{P}} = E_p^{A,n}\), we are left with
$$\begin{aligned}&{\tilde{a}}_2^h(E_{p}^{A,n} - E_{p}^{A,n-1}, E_{p}^{A,n}) + \varDelta t a_2^h(E_{p}^{A,n}, E_{p}^{A,n}) - b_2(E_{p}^{A,n}, E_{\psi }^{A,n} - E_{\psi }^{A,n-1}) \nonumber \\&\quad = \varDelta t ( \ell ^{\textrm{P}}(t_n)- \ell ^{n,\textrm{P}}_h, E_{p}^{A,n} )_{0, \varOmega ^{\textrm{P}}} + {\tilde{a}}_2^h( I^h_p p^{\textrm{P}}(t_n) - I^h_p p^{\textrm{P}}(t_{n-1}), E_{p}^{A,n}) \nonumber \\&\qquad - {\tilde{a}}_2((\varDelta t) \partial _{t} p^{\textrm{P}}(t_n), E_{p}^{A,n}) + b_2(E_{p}^{A,n}, (\varDelta t) \partial _{t} \psi - (I^h_{\psi } \psi (t_n)) - I^h_{\psi } \psi (t_{n-1})). \end{aligned}$$
(B.5)
Adding (B.4)–(B.5) and repeating the arguments used in deriving stability, we can assert that
$$\begin{aligned}&a_3(E_{\psi }^{A,n} - E_{\psi }^{A,n-1}, E_{\psi }^{A,n}) - b_2(E_{p}^{A,n} - E_{p}^{A,n-1}, E_{\psi }^{A,n})\\&\qquad - b_2(E_{p}^{A,n}, E_{\psi }^{A,n} - E_{\psi }^{A,n-1}) + {\tilde{a}}_2^h(E_{p}^{A,n} - E_{p}^{A,n-1}, E_{p}^{A,n}) \\&\quad = (\varDelta t) \bigg ( c_0(\delta _t E_{p}^{A,n} , E_{p}^{A,n})_{0,\varOmega ^{\textrm{P}}} + \frac{1}{\lambda } \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \big (\alpha ^2 (\delta _t (I- \varPi ^{0,k}_K) E_p^{A,n}, (I- \varPi ^{0,k}_K) E_p^{A,n})_{0,K} \\&\qquad -(\delta _t (\alpha \varPi ^{0,k}_K E_p^{A,n} - E_{\psi }^{A,n}), \alpha \varPi ^{0,k}_K E_p^{A,n} - E_{\psi }^{A,n})_{0,K} \big ) + \frac{\varDelta t}{\lambda } (\delta _t E_{\psi }^{A,n} , E_{\psi }^{A,n})_{0,\varOmega ^{\textrm{E}}} \bigg ), \end{aligned}$$
The left-hand side can be bounded using the inequality
$$\begin{aligned} (f_h^n - f_h^{n-1}, f_h^n) \ge \frac{1}{2} \bigl (\Vert f_h^n \Vert _0^2 - \Vert f_h^{n-1}\Vert _0^2 \bigr ), \end{aligned}$$
and then summing over n we get
$$\begin{aligned}&{\mu _{\min }} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,n})\Vert _0^2 + c_0 \Vert E_{p}^{A,n}\Vert _{0, \varOmega ^{\textrm{P}}}^2 + (1/ {\lambda ^{\textrm{E}}}) \Vert E_{\psi }^{A,n}\Vert _{0, \varOmega ^{\textrm{E}}}^2 + (\varDelta t) \frac{\kappa _{\min }}{\eta } \sum _{j=1}^n \Vert \nabla E_{p}^{A,j} \Vert _{0, \varOmega ^{\textrm{P}}}^2 \nonumber \\&\qquad + (1/{\lambda ^{\textrm{P}}}) \sum _{K \in \mathcal {T}_h^{\textrm{P}}}\bigg ( \alpha ^2 \Vert (I- \varPi ^{0,k}_K) E_p^{A,n} \Vert _{0,K}^2 + \Vert \alpha \varPi ^{0,k}_K E_p^{A,n} - E_{\psi }^{A,n}\Vert _{0,K}^2 \bigg ) \\&\quad \le {\mu _{\min }} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,0})\Vert _0^2 + c_0 \Vert E_{p}^{A,0}\Vert _{0, \varOmega ^{\textrm{P}}}^2 + (1/ {\lambda ^{\textrm{E}}}) \Vert E_{\psi }^{A,0}\Vert _{0, \varOmega ^{\textrm{E}}}^2 \\&\qquad + (1/{\lambda ^{\textrm{P}}}) \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \bigg ( \alpha ^2 \Vert (I- \varPi ^{0,k}_K) E_p^{A,0} \Vert _{0,K}^2 + \Vert \alpha \varPi ^{0,k}_K E_p^{A,0} - E_{\psi }^{A,0} \Vert _{0,K}^2 \bigg ) \\&\qquad + \underbrace{ \sum _{j=1}^n ( \varvec{b}(t_j) - \varvec{b}^j_h, E_{{\varvec{u}}}^{A,j}- E_{{\varvec{u}}}^{A,j-1} )_{0,\varOmega }}_{=:L_1} + \underbrace{ \sum _{j=1}^n \varDelta t ( \ell ^{\textrm{P}}(t_j)- \ell ^{j,\textrm{P}}_h, E_{p}^{A,j} )_{0,\varOmega ^{\textrm{P}}}}_{=:L_2} \nonumber \\&\qquad - \underbrace{ \sum _{j=1}^n b_1(({\varvec{u}}(t_j) - {\varvec{u}}(t_{j-1})) - (\varDelta t) \partial _{t} {\varvec{u}}(t_j), E_{\psi }^{A,j})}_{=:L_3} \nonumber \\&\qquad - \underbrace{ \sum _{j=1}^n b_2((I_p^h p^{\textrm{P}}(t_j) - I_p^h p^{\textrm{P}}(t_{j-1})) - (\varDelta t) \partial _{t} p^{\textrm{P}}(t_j), E_{\psi }^{A,j})}_{=:L_4} \nonumber \\&\qquad + \underbrace{ \sum _{j=1}^n a_3((I_{\psi }^h \psi (t_j) - I_{\psi }^h \psi (t_{j-1})) - (\varDelta t) \partial _{t} \psi (t_j), E_{\psi }^{A,j}) }_{:=L_5} \nonumber \\&\qquad + \underbrace{ \sum _{j=1}^n ({\tilde{a}}_2^h( I^h_p p^{\textrm{P}}(t_j) - I^h_p p^{\textrm{P}}(t_{j-1}), E_{p}^{A,j}) - {\tilde{a}}_2((\varDelta t) \partial _{t} p^{\textrm{P}}(t_j), E_{p}^{A,j}) )}_{:=L_6} \nonumber \\&\qquad + \underbrace{ \sum _{j=1}^n b_2(E_{p}^{A,j}, (\varDelta t) \partial _{t} \psi - (I^h_{\psi } \psi (t_j) - I^h_{\psi } \psi (t_{j-1}))}_{:=L_7}. \end{aligned}$$
We bound the term \(L_1\) with the help of the following formula
$$\begin{aligned} \sum _{j=1}^n (f_h^j - f_h^{j-1}, g_h^j )= (f_h^n, g_h^n) - (f_h^0, g_h^0)-\sum _{j=1}^n (f_h^{j-1}, g_h^j- g_h^{j-1}), \end{aligned}$$
and in combination with the orthogonality property of projection \(\varvec{\varPi }^{0,k}_K\) onto piecewise polynomials, and applying Taylor expansion
$$\begin{aligned} f^j - f^{j-1} = (\varDelta t) \partial _{t} f^j - \int _{t_{j-1}}^{t_j} (s-t_{j-1})\partial _{tt} f(s) \,\textrm{d}s, \end{aligned}$$
we get
$$\begin{aligned} L_1&{ = \sum _{K \in \mathcal {T}_h} \bigg ( ((\textbf{I}- \varvec{\varPi }^{0,k}_K)\varvec{b}(t_n), (\textbf{I}- \varvec{\varPi }^{0,0}_K) E_{{\varvec{u}}}^{A,n})_{0,K} - ((\textbf{I}- \varvec{\varPi }^{0,k}_K)\varvec{b}(0), (\textbf{I}- \varvec{\varPi }^{0,0}_K)E_{{\varvec{u}}}^{A,0})_{0,K}} \\&\quad {- \sum _{j=1}^{n} \left( (\textbf{I}- \varvec{\varPi }^{0,k}_K) \bigg ((\varDelta t)\partial _t \varvec{b}^j - \int _{t_{j-1}}^{t_j} (s-t_{j-1}) \partial _{tt}\varvec{b}(s) \,\textrm{d}s\bigg ), E_{{\varvec{u}}}^{A,j-1} \right) _{0,K} \bigg ).} \end{aligned}$$
Then the Cauchy–Schwarz and Korn inequalities, the estimates of projection \(\varvec{\varPi }^{0,k}_K\),and finally using the generalised Young’s inequality, gives the bound
$$\begin{aligned} L_1&{ \le }\frac{{\mu _{\min }}}{2} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,n})\Vert _0^2 +{C(\mu _{\min })} h^{2k} |\varvec{b}(t_n) |_{k-1}^2 + C(\mu _{\min }) h^k |b(t_0)|_{k-1}~ {\mu _{\min }^{{1/2}}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,0})\Vert _0 \\&\qquad +{ (\varDelta t) ~ \sum _{j=1}^n C(\mu _{\min }) \Big (\Vert \partial _{t} \varvec{b}^j \Vert _{0} + (\varDelta t)^{1/2} \Vert \partial _{tt} \varvec{b}\Vert _{\textbf{L}^2(t_{j-1},t_j;[L^2(\varOmega )]^2)} \Big )} {\mu _{\min }^{{1/2}}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j-1}) \Vert _0. \end{aligned}$$
Then the estimate satisfied by the projection \(\varPi ^{0,k}_K\) along with Poincaré and Young’s inequalities yield
$$\begin{aligned} L_2 \,&{ = \sum _{j=1}^n \varDelta t ((I - \varPi ^{0,k}_K)\ell ^{\textrm{P}}(t_j), (I - \varPi ^{0,0}_K)E_{p}^{A,j} )_{0,\varOmega ^{\textrm{P}}} } \\&{\le C(\kappa _{\min },\eta )} {h^{2k} } (\varDelta t) \sum _{j=1}^n | \ell ^{\textrm{P}}(t_j)|_{k-1, \varOmega ^{\textrm{P}}}^2 + (\varDelta t) \frac{\kappa _{\min }}{6 \eta } \sum _{j=1}^n \Vert \nabla E_{p}^{A,j}\Vert _{0, \varOmega ^{\textrm{P}}}^2. \end{aligned}$$
The discrete inf-sup condition (3.3) implies that
$$\begin{aligned} \Vert E_{\psi }^{A,j} \Vert _0 \lesssim { h^k \Vert \varvec{b}^j\Vert _{k-1} + {\mu _{\max }} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j}) \Vert _0 \le h^k \Vert \varvec{b}^j\Vert _{k-1} + C(\mu ) {\mu _{\min }^{1/2}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j}) \Vert _0.} \end{aligned}$$
(B.6)
Applying an expansion in Taylor series, together with (B.6) the Cauchy–Schwarz inequality enable us to write
$$\begin{aligned} L_3 \,&{ \lesssim \sum _{j=1}^n \Vert ({\varvec{u}}(t_j) - {\varvec{u}}(t_{j-1})) - (\varDelta t) \partial _{t} {\varvec{u}}(t_j) \Vert _0 \Vert E_{\psi }^{A,j} \Vert _0 }\\&\lesssim {(\varDelta t)^{3/2}\sum _{j=1}^n \Vert \partial _{tt} {\varvec{u}}\Vert _{\textbf{L}^2(t_{j-1},t_j;[L^2(\varOmega )]^2)} \Big ( h^k \Vert \varvec{b}^j\Vert _{k-1} + C(\mu ) {\mu _{\min }^{1/2}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j}) \Vert _0 \Big ).} \end{aligned}$$
Then, applying again the Cauchy–Schwarz inequality, and after using the interpolation estimates (4.2b), the bound (B.6), we get
$$\begin{aligned} L_4 \,&{\lesssim \frac{\alpha }{\lambda ^{\textrm{P}}} \sum _{j=1}^n (\varDelta t) \Big ( \Vert I_p^h (\delta _t p^{\textrm{P}}(t_j)) - \delta _t p^{\textrm{P}}(t_j) \Vert _{0, \varOmega ^{\textrm{P}}} + \Vert \delta _t p^{\textrm{P}}(t_j)- \partial _{t} p^{\textrm{P}}(t_j) \Vert _{0, \varOmega ^{\textrm{P}}} \Big ) \Vert E_{\psi }^{A,j} \Vert _{0, \varOmega ^{\textrm{P}}}} \\&\lesssim { \frac{\alpha }{{\lambda ^{\textrm{P}}}} (\varDelta t)^{1/2} \sum _{j=1}^n \bigg ( h^{k+1} \Vert \partial _{t}p^{\textrm{P}}\Vert _{L^2(t_{j-1},t_j, H^{k+1}( \varOmega ^{\textrm{P}}))} + (\varDelta t ) \Vert \partial _{tt} p^{\textrm{P}}\Vert _{L^2(t_{j-1},t_j, L^2(\varOmega ^{\textrm{P}}))} \bigg )} \\&\quad \times {( h^k \Vert \varvec{b}^j\Vert _{k-1} + C(\mu ) {\mu _{\min }^{1/2}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j}) \Vert _{0} ).} \end{aligned}$$
On the other hand, the stability of \(a_3(\cdot , \cdot )\), the proof for the bound of \(L_4\) , and the interpolant estimate (4.2a) gives
$$\begin{aligned} L_5&{\lesssim \frac{1}{\lambda } \sum _{j=1}^n (\varDelta t) \Big ( \Vert I_{\psi }^h (\delta _t \psi (t_j) ) - \delta _t \psi (t_j) \Vert _{0} + \Vert \delta _t \psi (t_j) - \partial _{t} {\psi }(t_j) \Vert _{0} \Big ) \Vert E_{\psi }^{A,j} \Vert _{0} } \\&{\lesssim \frac{1}{\lambda _{\min }} (\varDelta t)^{1/2} \sum _{j=1}^n \bigg ( h^k \big (\Vert \partial _t{\varvec{u}}\Vert _{\textbf{L}^2(t_{j-1},t_j; [H^{k+1}(\varOmega )]^2)}^2 + \Vert \partial _t\psi \Vert _{L^2(t_{j-1},t_j; H^k(\varOmega ))}^2 \big )^{1/2} } \\&\quad + {(\varDelta t) \Vert \partial _{tt} \psi \Vert _{L^2(t_{j-1},t_j;L^2(\varOmega ))} \bigg ) ( h^k \Vert \varvec{b}^j\Vert _{k-1} + C(\mu ) {\mu _{\min }^{1/2}} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,j}) \Vert _0 ),} \end{aligned}$$
where \(\Vert \partial _t v\Vert _{\textbf{L}^2(t_{j-1},t_j; H(\varOmega ))}^2:= \sum _{D=\varOmega ^{\textrm{P}}, \varOmega ^{\textrm{E}}} \Vert \partial _t v^D\Vert _{\textbf{L}^2(t_{j-1},t_j; H(D))}^2\) for \(v={\varvec{u}}, \psi \) with space \(H(D)=[H^{k+1}(D)]^2, H^k(D)\) respectively. The polynomial approximation \(\varPi ^{0,k}_K p^{\textrm{P}}\) for fluid pressure, consistency of the bilinear form \({\tilde{a}}_2^h(\cdot , \cdot )\), stability of the bilinear forms \({\tilde{a}}_2(\cdot , \cdot ), {\tilde{a}}_2^h(\cdot , \cdot )\), the usage of the Cauchy–Schwarz, Poincaré and Young’s inequalities gives
$$\begin{aligned} L_6&{\le C \Big ( c_0 +\frac{\alpha ^2}{{\lambda ^{\textrm{P}}}} \Big ) \sum _{j=1}^n \sum _{K \in \mathcal {T}_h^{\textrm{P}}} (\varDelta t) \Big ( \Vert (I^h_p - \varPi ^{0,k}_K) (\delta _t p^{\textrm{P}}(t_j) ) \Vert _{0,K} + \Vert (I- \varPi ^{0,k}_K )(\delta _t p^{\textrm{P}}(t_j) ) \Vert _{0,K} } \\&\qquad \qquad \qquad \qquad \qquad \qquad {+ \Vert \delta _t p^{\textrm{P}}(t_j) - \partial _{t} p^{\textrm{P}}(t_j) \Vert _{0,K} \Big ) \Vert \nabla E_{p}^{A,j} \Vert _{0, K} }\\&{\le C(\kappa _{\min }, \eta ) \Big ( c_0 +\frac{\alpha ^2}{\lambda ^{\textrm{P}}} \Big ) \sum _{j=1}^n (\varDelta t)^{1/2}} \sum _{K \in \mathcal {T}_h^{\textrm{P}}} \! \Big ( h^{k+1} \Vert \partial _t p^{\textrm{P}} \Vert _{L^2(0,t_n;H^{k+1}(\varOmega ^{\textrm{P}}))} \\ {}&\qquad \qquad \qquad \qquad \qquad \qquad + (\varDelta t) \Vert \partial _{tt} p^{\textrm{P}} \Vert _{L^2(0,t_n;L^2(\varOmega ^{\textrm{P}}))} \Big ) \Vert \nabla E_{p}^{A,j} \Vert _{0, K} \\&{\le C(\kappa _{\min }, \eta ) (\varDelta t)} {\Big ( c_0 + \frac{\alpha ^2}{{\lambda ^{\textrm{P}}}} \Big )^2} \Big ( h^{2(k+1)} \Vert \partial _t p^{\textrm{P}} \Vert _{L^2(0,t_n;H^{k+1}(\varOmega ^{\textrm{P}}))}^2 + (\varDelta t)^2 \Vert \partial _{tt} p^{\textrm{P}} \Vert _{L^2(0,t_n;L^2(\varOmega ^{\textrm{P}}))}^2 \Big ) \\&\qquad \qquad \qquad \qquad \qquad \qquad + \varDelta t \frac{\kappa _{\min }}{6 \eta } \sum _{j=1}^n \Vert \nabla E_{p}^{A,j} \Vert _{0, \varOmega ^{\textrm{P}}}^2. \end{aligned}$$
On the other hand, the continuity of \(b_2(\cdot , \cdot )\), the bound derived for the term \(L_5\) and using Young’s inequality implies that
$$\begin{aligned} L_7&{\le \frac{\alpha }{\lambda }\sum _{j=1}^n \bigg ( h^k (\varDelta t)^{1/2} \big (\Vert \partial _t{\varvec{u}}\Vert _{\textbf{L}^2(t_{j-1},t_j; [H^{k+1}(\varOmega )]^2)}^2 + \Vert \partial _t\psi \Vert _{L^2(t_{j-1},t_j; H^k(\varOmega ))}^2 \big )^{1/2} } \\ {}&\qquad \qquad \qquad { + (\varDelta t)^{3/2} \Vert \partial _{tt} \psi \Vert _{L^2(t_{j-1},t_j;L^2(\varOmega ))} \bigg ) \Vert E_{p}^{A,j}\Vert _{0, \varOmega ^{\textrm{P}}} } \\&{\le C(\kappa _{\min }, \eta ) } {\Big ( \frac{\alpha }{{\lambda _{\min }}} \Big )^2}\bigg ( h^{2k} (\Vert \partial _{t} \psi \Vert _{L^2(0,t_n;{H^k(\varOmega )})}^2 + \Vert \partial _{t} {\varvec{u}}\Vert _{\textbf{L}^2(0,t_n;[H^{k+1}(\varOmega )]^2)}^2)\\ {}&\qquad \qquad \qquad + (\varDelta t )^2 \Vert \partial _{tt} \psi \Vert _{L^2(0,t_n;L^2(\varOmega ))}^2 \bigg ) +(\varDelta t) \frac{\kappa _{\min }}{6 \eta } \sum _{j=1}^n \Vert \nabla E_{p}^{A,j} \Vert _{0{, \varOmega ^{\textrm{P}}}}^2. \end{aligned}$$
In turn, putting together the bounds obtained for all \(L_i\)’s, \(i=1, \dots , 7\), using Young’s inequality and Lemma 3.2 allows us to conclude that
$$\begin{aligned}&{\mu _{\min }} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,n})\Vert _0^2 + c_0 \Vert E_{p}^{A,n}\Vert _{0, \varOmega ^{\textrm{P}}}^2 + (1/{\lambda ^{\textrm{E}}}) \Vert E_{\psi }^{A,n}\Vert _{0, \varOmega ^{\textrm{E}}}^2 + (\varDelta t) \frac{\kappa _{\min }}{\eta } \sum _{j=1}^n \Vert \nabla E_{p}^{A,j} \Vert _{0, \varOmega ^{\textrm{P}}}^2 \\&\quad \lesssim {\mu _{\min }} \Vert \varvec{\varepsilon }(E_{{\varvec{u}}}^{A,0})\Vert _0^2 + (c_0 + \alpha ^2/{\lambda ^{\textrm{P}}}) \Vert E_{p}^{A,0}\Vert _{0, \varOmega ^{\textrm{P}}}^2 + (1/{\lambda ^{\textrm{E}}}) \Vert E_{\psi }^{A,0}\Vert _{0, \varOmega ^{\textrm{E}}}^2 \\&\qquad + C(\mu _{\min }) \Big (1 + \varDelta t \Big ) h^{2k} \max _{0 \le j \le n} |\varvec{b}(t_j) |_{k-1}^2 \\&\qquad + { C(\mu _{\min }, \kappa _{\min }, \eta , \alpha , \lambda _{\min }) (\varDelta t)^{1/2} }\\&\qquad { \times \biggl \{ \bigg ( \sum _{j=1}^n \Big ( (\varDelta t)^{1/2} \Vert \partial _{t} \varvec{b}^j\Vert _0 + (\varDelta t) \Vert \partial _{tt} \varvec{b}\Vert _{\textbf{L}^2(t_{j-1},t_j;[L^2(\varOmega )]^2)} } \\&\qquad { + h^k ( h\Vert \partial _{t} p^P \Vert _{L^2(t_{j-1},t_j;H^{k+1}(\varOmega ^{\textrm{P}}))} + (\Vert \partial _{t} {\varvec{u}}\Vert _{\textbf{L}^2(t_{j-1},t_j;[H^{k+1}(\varOmega )]^2)}^2 + \Vert \partial _{t} \psi \Vert _{L^2(t_{j-1},t_j;H^{k}(\varOmega ))}^2)^{1/2}) \Big )} \\&\qquad { + (\varDelta t) \Big ( \Vert \partial _{tt} {\varvec{u}}\Vert _{\textbf{L}^2(t_{j-1},t_j;[L^2(\varOmega )]^2)} + \Vert \partial _{tt} p^P \Vert _{L^2(t_{j-1},t_j;L^2(\varOmega ^{\textrm{P}}))} + \Vert \partial _{tt} \psi \Vert _{L^2(t_{j-1},t_j;L^2(\varOmega ))} \Big ) \bigg )^2}\\&\qquad { + \sum _{j=1}^n h^{2 k} \Big ( (\varDelta t)^{1/2} |\ell ^P(t_j)|_{k-1,\varOmega ^{\textrm{P}}}^2 + |\varvec{b}^j|_{k-1}^2 \Big ) } \\&\qquad { + h^{2k} \Big ( h^2 \Vert \partial _{t} p^P \Vert _{L^2(0,t_n;H^{k+1}(\varOmega ^{\textrm{P}}))}^2 + \Vert \partial _{t} {\varvec{u}}\Vert _{\textbf{L}^2(0,t_n;[H^{k+1}(\varOmega )]^2)}^2 + \Vert \partial _{t} \psi \Vert _{L^2(0,t_n;H^k(\varOmega ))}^2 \Big )} \\&\qquad { + (\varDelta t)^2 \Big ( \Vert \partial _{tt} p^P\Vert _{L^2(0,t_n;L^2(\varOmega ^{\textrm{P}}))} ^2 + \Vert \partial _{tt} {\varvec{u}}\Vert _{\textbf{L}^2(0,t_n; [L^2(\varOmega )]^2)}^2 + \Vert \partial _{tt} \psi \Vert _{L^2(0,t_n;L^2(\varOmega ))}^2 \Big ) \biggr \}. } \end{aligned}$$
And finally, the desired result (4.3) is proven after choosing \({\varvec{u}}_h^0: ={\varvec{u}}_I(0)\), \(\psi _h^0: = \varPi ^{0,k-1}\psi (0)\), \(p_h^{0, \textrm{P}}: = p_I^{\textrm{P}}(0)\) and applying triangle’s inequality together with (B.1a)–(B.1c) and (B.6).
