\documentclass[reqno]{amsart} \AtBeginDocument{{\noindent\small 2004-Fez conference on Differential Equations and Mechanics \newline {\em Electronic Journal of Differential Equations}, Conference 11, 2004, pp. 41--51.\newline ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu \newline ftp ejde.math.txstate.edu (login: ftp)} \thanks{\copyright 2004 Texas State University - San Marcos.} \vspace{9mm}} \setcounter{page}{41} \begin{document} \title[\hfilneg EJDE/Conf/11 \hfil Numerical analysis of Euler-SUPG method] {Numerical analysis of Euler-SUPG modified method for transient viscoelastic flow} \author[M. Bensaada, D. Esselaoui\hfil EJDE/Conf/11 \hfilneg] {Mohammed Bensaada, Driss Esselaoui} % in alphabetical order \address{Mohammed Bensaada \hfill\break Laboratoire des Sciences de l'Ing\'enieur\\ Analyse Num\'erique et Optimisation (SIANO)\\ Facult\'e des Sciences, Universit\'e Ibn Tofail \\ B.P.133, 14000-K\'enitra, Maroc} \address{Driss Esselaoui \hfill\break Laboratoire des Sciences de l'Ing\'enieur \\ Analyse Num\'erique et Optimisation (SIANO)\\ Facult\'e des Sciences, Universit\'e Ibn Tofail \\ B.P.133, 14000-K\'enitra, Maroc} \email{desselaoui@yahoo.fr \; Fax: 212 37 37 27 70} \date{} \thanks{Published October 15, 2004.} \thanks{Supported by CNRST - CNRS, Program STIC01/03, and by CNRST- GRICES Portugal} \subjclass[2000]{65N30, 34K25, 74S05} \keywords{Stabilized method; finite elements; modified SUPG method; \hfill\break\indent transient viscoelastic flow} \begin{abstract} We study a new approximation scheme of transient viscoelastic fluid flow obeying an Oldroyd-B type constitutive law. The approximation stress, velocity, pressure are respectively $P_{1}$-continuous, $P_{2}$-continuous, $P_{1}$-continuous. We use the modified streamline upwinding Petrov-Galerkin method induced by the modified Euler method. We assume that the continuous problem admits a sufficiently smooth and sufficiently small solution. We show that the approximate problem has a solution and we give an error bound. \end{abstract} \maketitle \newtheorem{theorem}{Theorem}[section] \newtheorem{lemma}[theorem]{Lemma} \newtheorem{remark}[theorem]{Remark} \numberwithin{equation}{section} \allowdisplaybreaks \section{Introduction and Presentation of the problem} In the numerical simulation of the viscoelastic fluid flows, the hyperbolic character of the constitutive equation (when using differential models) has to be taken into account (see \cite{EsF,MM}). This hyperbolic character implies that some upwinding is needed to avoid oscillations as in the method of characteristics \cite{EsF,B}, the Lesaint-Raviart discontinuous finite element method \cite{BBS,EsR}, the streamline-upwind method (SU) and the streamline-upwind-Petrov-Galerkin method (SUPG) \cite{EsRZ,MM}. The numerical analysis of the steady case of the viscoelastics fluids flows is abundant. Although the list is not exhaustive, one may see for example \cite{BBS,EsRZ,S}. Moreover the numerical analysis for transient viscoelastic flow remain quite few \cite{K,WB}. For example, some difficulties appear, when we use continuous finite element approximation for $(\sigma,u,p)$ and the standard SUPG method for the convection of the extra stress tensor. To give some response to this difficulties, we develop in this paper the study of continuous finite element(F.E) approximation of a transient viscoelastic fluid flow obeying an Oldroyd-B model. For the convective term of the constitutive equation we use some modified SUPG method linked to a variant of implicit Euler method (see \cite{BEs}). Under this condition we are able to show that the approximate problem is stabilized and has a solution and we give an error bound. The transient viscoelastic fluid flow obeying an Oldroyd-B type constitutive law is considered flowing in a bounded, connected open set $\Omega$ in $IR^{2}$ with lipshitzian boundary $\Gamma$; $n$ is the outward unit normal to $\Gamma$. The basic set of equations of the Oldroyd-B model with a single relaxation time is given by \begin{equation} \label{O} \begin{gathered} \lambda (\frac{\partial\sigma}{\partial t}+(u.\nabla)\sigma+\beta(\sigma,\nabla u) ) +\sigma-2\alpha D(u)=0 \quad \mbox{in } \Omega \times ] 0,T [,\\ Re\;\frac{\partial u}{\partial t}-\nabla.\sigma-2(1-\alpha)\nabla.D(u)+\nabla p =f\quad \mbox{in } \Omega \times] 0,T[, \\ \nabla.u =0 \quad\mbox{in } \Omega\times] 0,T[,\\ u=0 \quad \mbox{on } \Gamma\times] 0,T [,\\ u=u_{0} ,\quad \sigma=\sigma_{0}\quad\mbox{in } \Omega \; t=0. \end{gathered} \end{equation} Where $\lambda > 0$, $Re$ and $0 <\alpha < 1$ are respectively the Weissenberg number, the Reynolds number and the viscosity ratio constant. $f$ is a density of forces. $D(u)=\frac{1}{2}(\nabla u+\nabla u^{\top})$ the rate of strain tensor, and $\beta(\sigma,\nabla u)=-\nabla u\sigma-\sigma\nabla u^{\top}$. \begin{remark} \label{rmk1} \rm The boundary condition $u=0$ on $\Gamma$ can be replaced by $u= u_d$ on $\Gamma$. Regarding $\sigma$ and the hyperbolic character of the constitutive equation, we have to impose $\sigma = \sigma_d$ on $\Gamma^{-}= \{x \in\Gamma;\;u_d.n(x) <0 \}$. \end{remark} \begin{remark} \label{rmk2} \rm The inertia term $(u.\nabla)u$ is neglected in the momentum equation in order to make the analysis simpler. \end{remark} Let us define the following spaces: \begin{gather*} T=\{\tau=(\tau_{ij})_{1\leq i,j\leq 2} : \tau_{ij}=\tau_{ji};\,\tau_{ij} \in L^{2}(\Omega);\,i,j=1,2\}, \quad X=(H_{0}^{1}(\Omega))^{2} \\ Q=\{ q\in L^{2}(\Omega) / {\int_{\Omega}}q\,dx=0\},\quad V=\{v\in X / (q,\nabla.v)=0;\forall q\in Q \}. \end{gather*} The norm and scalar product in $L^{2}(\Omega)$ of functions,\ vectors and tensors are denoted respectively by $\vert \cdot\vert$ and $(\cdot,\cdot)$; ($\vert\cdot\vert_{\Gamma}$ and $(\cdot,\cdot)_{\Gamma}$ in $L^{2}(\Gamma)$); $\langle f,v\rangle $ will denote the duality between $f\in (H^{-1}(\Omega))^{2}$ and $v\in X$. \begin{remark} \label{rmk3} \rm Existence results for problem \eqref{O} are proved in \cite{GS}. In order to make some theoretical analysis of approximate problem of \eqref{O} we use the regularity imposed in \cite{GS}. \end{remark} \section{Description of the approximation scheme} \subsection*{FE approximation} We suppose $\Omega$ polygonal and we consider a triangulation $\Im_{h}$ on $\Omega$ made of triangles $K$ such that $\overline{\Omega}=\{\bigcup K;K\in \Im_{h}\}$ uniformly regular, $\exists \nu_{0},\nu_{1}: \nu_{0}h\leq h_{K}\leq \nu_{1}\varrho_{K}$ where $\varrho_{K}$ is the diameter of the greatest ball included in $K$ and $h_{\rm max}=\max_{K\in\Im_{h}}h_{K}$. We use the Taylor-Hood finite element method for approximations in space of $(u,p)$: $P_{2}$-continuous in velocity, $P_{1}$-continuous in pressure and we consider $P_{1}$-continuous approximation of the stresses. The corresponding $FE$ space are: \begin{gather*} X_{h}=\{v\in X \cap C^{0}(\Omega)^{2}: v_{|K} \in P_{2}(K)^{2}, \forall K\in \Im_{h}\}\\ Q_{h}=\{q\in Q \cap C^{0}(\Omega): q_{|K} \in P_{1}(K), \forall K\in \Im_{h}\}\\ V_{h}=\{ v\in X_{h}: (q,\nabla . v)=0, \,\forall q\in Q_{h}\}. \\ T_{h}=\{\tau\in T\cap C^{0}(\Omega):\tau|_{K}\in P_{1}(K), \forall K\in \Im_{h}\}\,, \end{gather*} where $P_{m}(K)$ denotes the space of polynomials of degrees less or equal to $m$ on $K\in T_{h}$. The term $((u.\nabla)\sigma, \tau)$ is approximated by means of an operator $B$ on $X_{h}\times T_{h}\times T_{h}$ defined by \begin{align*} B(u_{h},\sigma_{h};\tau_{h}) &=((u_{h}(t).\nabla)\sigma_{h}(t),\tau_{h}) +\delta_{0}(h,t)((u_{h}(t)\nabla)\sigma_{h}(t),(u_{h}(t).\nabla)\tau)\\ &\quad +(1/2)((\nabla.u_{h})(t) \sigma,\tau). \end{align*} For the steady case you can see \cite{S}. \subsection*{Numerical method} We propose Euler-SUPG modified scheme, implicit in time, based on the scheme proposed for the transport equation (see \cite{BEs}). We construct an approximation of the solution at each time step $nk, n=0,\dots,N$ in the following way. We start with $u_{h}^{0}=\tilde{u}_{0}$: elliptic projection of $u_{0}$ into $V_{h},\sigma_{h}^{0}=\tilde{\sigma}_{0}$: orthogonal projection of $\sigma_{0}$ into $T_{h}$. Given $u_{h}^{0},\dots,u_{h}^{n}$; $\sigma_{h}^{0},\dots,\sigma_{h}^{n}$, because $(X_{h},Q_{h})$ satisfies the inf sup condition, we look for the solution of the following problem, find $(u_{h}^{n+1},\sigma_{h}^{n+1})\in V_{h}\times T_{h}$, such that \begin{gather} \begin{aligned} &\lambda(\frac{\sigma_{hu_{h}^{n},\delta}^{n+1}-\sigma_{hu_{h}^{n},\delta}^{n}}{k}, \tau_{u_{h}^{n},\lambda})+ (\sigma_{h}^{n+1},\tau_{u_{h}^{n},\lambda})+B(\lambda u_{h}^{n},\sigma_{h}^{n+1};\tau) \\ &-2\alpha (D(u_{h}^{n+1}),\tau_{u_{h}^{n},\lambda})+ \lambda(\beta(\sigma_{h}^{n+},\nabla u_{h}^{n+1}),\tau_{u_{h}^{n},\lambda})=0\quad \forall\tau\in T_{h}; \end{aligned} \label{eq:sc1} \\ \begin{aligned} &Re(\frac{u_{h}^{n+1}-u_{h}^{n}}{k},v)+(\sigma_{h}^{n+1},D(v))+ 2(1-\alpha))(D(u_{h}^{n+1}),D(v)) \\ &=\langle f(t_{n+1},x),v\rangle\quad \forall v\in V_{h}, \end{aligned} \label{eq:sc2} \end{gather} where $\sigma_{hu_{h}^{n},\delta}^{i}=\sigma_{h}^{i}+\delta(h,k)(u_{h}^{n}. \nabla)\sigma_{h}^{i}$ $(i=n,n+1)$; $\tau_{u_{h}^{n},\lambda}=\tau+\lambda \delta_{0}(h,k)(u_{h}^{n}.\nabla)\tau$ for all $\tau\in T_{h}$, and $ \delta (,)$ (resp. $\delta_{0}(,)$) will be specified later. In order to show that equation $(\ref{eq:sc1})-(\ref{eq:sc2})$ defined uniquely $(u_{h}^{n+1},\sigma_{h}^{n+1})$, we multiply equation \eqref{eq:sc2} by $2\alpha$ and add the equation obtained to equation \eqref{eq:sc1}, we get \begin{align*} &\lambda(\frac{\sigma_{hu_{h}^{n},\delta}^{n+1} -\sigma_{hu_{h}^{n},\delta}^{n}}{k},\tau_{u_{h}^{n},\lambda})+2\alpha Re(\frac{u_{h}^{n+1}-u_{h}^{n}}{k},v_{h})+ B(\lambda u_{h}^{n}, \sigma_{h}^{n+1};\tau_{h})\\ &+A(u_{h}^{n};(\sigma_{h}^{n+1},u_{h}^{n+1}),(\tau,v)) +\lambda (\beta(\sigma_{h}^{n+1},\nabla u_{h}^{n+1}),\tau_{u_{h}^{n},\lambda})\\ &=2\alpha\langle f(t_{n+1}),v_{h}\rangle\,,\quad \forall (\tau_{h},v_{h}) \in T_{h} \times V_{h}. \end{align*} where $A(.,.)$ is a bilinear form on $T_{h}\times V_{h}$ defined as $$ A(w;(\sigma,u),(\tau,v))=(\sigma,\tau_{w,\lambda})+2\alpha(\sigma,D(v)) -2\alpha(D(u),\tau_{w,\lambda})+4\alpha(1-\alpha)(D(u),D(v)) $$ From this, we can establish some error bound and then the following existence result \begin{theorem} \label{thm1} There exists $M_{0}$, and $h_{0}$ such that if problem \eqref{O} admits a solution $(\sigma,u,p)$ with, \begin{align*} M=\max\big\{&\Vert\sigma\Vert_{C^{1}([t_{n},t_{n+1}];H^{2})}, \Vert u\Vert_{C^{1}([t_{n},t_{n+1}];H^{3})}, \Vert p\Vert_{C^{0}([t_{n},t_{n+1}];H^{2})}, \\ &\Vert\sigma\Vert_{C^{2}([t_{n},t_{n+1}];L^{2})}, \Vert u\Vert_{C^{2}([t_{n},t_{n+1}],T;L^{2})}\big\}< M_{0}, \end{align*} then if $\delta(,)$ satisfies $\delta(,)=\lambda \delta_{0}(,)$ and $h\leq h_{0}$ there exists a unique solution in $T_{h}\times V_{h}$ of problem \eqref{eq:sc1}- \eqref{eq:sc2}. \end{theorem} \begin{proof} For this purpose we define a mapping $\Phi: T_{h}\times V_{h}\longmapsto T_{h}\times V_{h}$, which to $(\sigma_{1},u_{1})$ associates $(\sigma_{2},u_{2})=\Phi(\sigma_{1},u_{1})$, where $(\sigma_{2},u_{2})\in T_{h}\times V_{h}$ satisfies \begin{align*} &\lambda(\frac{\sigma_{2u_{h}^{n},\delta} -\sigma_{hu_{h}^{n},\delta}^{n}}{k},\tau_{u_{h}^{n},\lambda}) +2\alpha Re(\frac{u_{2}-u_{h}^{n}}{k},v_{h}) + B(\lambda u_{h}^{n},\sigma_{h}^{n+1};\tau_{h})\\ &+A(u_{h}^{n};(\sigma_{2},u_{2}),(\tau,v)) \\ &=-\lambda (\beta(\sigma_{1},\nabla u_{1}),\tau_{u_{h}^{n},\lambda})+2\alpha\langle f(t_{n+1}),v_{h}\rangle\,,\quad \forall (\tau_{h},v_{h}) \in T_{h}\times V_{h}. \end{align*} We define a ball $B_{h}^{(n+1)}$ as follows: let $C^{\ast}$ be given. Then we define \begin{align*} B_{h}^{(n+1)}=\big\{&(\tau_{h},v_{h})\in T_{h}\times V_{h}: [\frac{\alpha Re}{k}\Vert v_{h}-u(t_{n+1})\Vert_{0,\Omega}^{2}\\ &+\frac{\lambda}{2k}\Vert(\tau_{h}-\sigma(t_{n+1}))_{u_{h}^{n},\lambda} \Vert_{0,\Omega}^{2}]^{1/2} \leq C^{\ast}(k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}} +h^{\frac{3}{2}}),\\ & [\frac{1}{4}\{4\alpha(1-\alpha)\Vert D(v_{h}-u(t_{n+1})\Vert_{0,\Omega}^{2} +\Vert\tau_{h}-\sigma(t_{n+1})\Vert_{0,\Omega}^{2}]^{1/2}\\ & \leq C^{\ast}(k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}} +h^{\frac{3}{2}})\big\}. \end{align*} The proof is decomposed in five parts:\\ (a) $\Phi$ is well defined for $h\leq h_{1}=h_{1}(h_{\rm max},\alpha)$ and bounded on bounded sets,\\ (b) let $C_{0}$ be a positive constant independent of $h,k$ and $\lambda,\alpha,Re$. If $h\leq h_{0}=\min\big\{h_{1},\sqrt{k} \big(\frac{C^{*}}{C_{0}M\sqrt{\alpha Re}}\big)^{2/3}, \frac{1}{C_{0}M} \sqrt{k/ \lambda}\big\}$ and $\sqrt{\delta/ k}\leq \frac{1}{C_{0}M}$, we have $B_{h}^{(n+1)}$ is non\-empty. On the other hand, if $M_{0}=C^{*} \frac{\gamma}{\lambda+\gamma}$, $\gamma=\sqrt{\alpha(1-\alpha)}$, we can prove that $\Phi({\overbrace{B}^{\circ}}{_{h}^{(n+1)}}) \subset {\overbrace{B}^{\circ}}{_{h}^{(n+1)}}$, \\ (c) $\Phi$ is continuous on $T_{h}\times V_{h}$. $\Phi(B_{h}^{(n+1)})\subset B_{h}^{(n+1)}$, Brouwer's theorem then gives the existence of fixed point $(\sigma_{h}^{n+1},u_{h}^{n+1})$ of $\Phi$ solution of problem $(1)-(2)$. \\ (d) Furthermore, if $\lambda M$ and $\lambda M\gamma^{-1}$ is sufficiently small, $\Phi$ is a contraction mapping on $B_{h}^{(n+1)}$. Then a result of existence and uniqueness follows from the fixed theorem. \end{proof} \begin{remark} \label{rmk4} When we use only the classical Euler-scheme in time and SUPG method, we can't point out result of (a) in the above proof. \end{remark} \section{Main result and error estimates} Suppose that the continuous problem admits a sufficiently smooth and sufficiently small solution, we can show, the following result. \begin{theorem} \label{thm2} There exists $M_{0}$ and $h_{0}$ such that if problem \eqref{O} admits a solution $(\sigma,u,p)$ with $\sigma\in C^{1}([0,T],(H^{2})^{4})\cap C^{2}([0,T],(L^{2})^{4})$; $u\in C^{1}([0,T],(H^{3})^{2})\cap C^{2}([0,T],(L^{2})^{2})$; $p\in L^{2}([0,T],(H^{2})\cap L_{0}^{2})\cap C^{0}([0,T],H^{2})$, satisfying $$ \max\{\Vert\sigma\Vert_{C^{1}(0,T;H^{2})},\Vert u\Vert_{C^{1}(0,T;H^{3})},\Vert p\Vert_{C^{0}(0,T;H^{2})}, \Vert\sigma\Vert_{C^{2}(0,T;L^{2})},\Vert u\Vert_{C^{2}(0,T;L^{2})}\} \leq M_{0} $$ then for $kh^{-(1+\varepsilon)}(0<\varepsilon\leq 1/2)$ bounded, there exists a constant $C$ independent of $h$ and $k$ such that $$ \max_{0\leq n\leq N}\vert (\sigma_{h}^{n}-\sigma(t_{n}))_{u_{h}^{n-1}}\vert+( \sum_{n=0}^{N}k\vert\sigma_{h}^{n}-\sigma(t_{n})\vert)^{2})^{1/2} \leq C (k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}} +h^{1+\varepsilon}) $$ and $$ \max_{0\leq n\leq N}\vert u_{h}^{n}-u(t_{n})\vert+(\sum_{n=0}^{N}k\vert D(u_{h}^{n}-u(t_{n}))\vert^{2})^{1/2}\leq C (k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}}+h^{1+\varepsilon}) $$ where $N\in N^{\ast}:Nk=T$, $u_{h}^{-1}=u_{h}^{0}$ and $(\sigma_{h}^{i},u_{h}^{i})_{1\leq i\leq N}$ are solutions of $(\eqref{eq:sc1}-\eqref{eq:sc2})_{1\leq i\leq N} $. \end{theorem} \begin{remark} \label{rmk5}\rm The discontinuous stresses approach case with Euler semi-implicit method in time was treated by Baranger and all. \cite{WB}. In the proof of the present result we can fund some similar technical like in \cite{WB}. \end{remark} \begin{proof}[Proof of Theorem \ref{thm2}] For $0\leq n0$, the ball $B_{h,k}^{m}$ for $0\leq m\leq N$, by \begin{align*} &B_{h,k}^{m}\\ &=\big\{(\tau_{i},v_{i})_{i=0,..,m}\in(T_{h}\times V_{h})^{m+1}: \max_{0\leq i\leq m}\{\vert(\tau_{i}-\sigma(t_{i}))_{v_{i-1}}\vert^{2} +\vert v_{i}-u(t_{i})\vert^{2}\}^{1/2}\\ &\quad \leq C_{0}(k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}} +h^{1+\varepsilon})\mbox{ and}\\ &\quad [\sum_{n=0}^{m}k\lbrace\vert\tau_{i}-\sigma(t_{i})\vert^{2}+\vert D(v_{i}-u(t_{i}))\vert^{2}\rbrace]^{1/2} \leq C_{0}(k+\delta+h\sqrt{\delta_{0}}+\frac{h^{2}}{\sqrt{\delta_{0}}} +h^{1+\varepsilon})\big\}. \end{align*} Our aim is to prove that we can choose $M_{0},h_{0},C_{0}$ such that for $M\leq M_{0}$, $h\leq h_{0}$ if $(\sigma_{h}^{n},u_{h}^{n})_{0\leq n\leq m-1}\in B_{h,k}^{m-1}$ for a $C_{0}=C_{0}(M_{0},h_{0},C_{i})$ then $(\sigma_{h}^{n},u_{h}^{n})_{0\leq n\leq m}\in B_{h,k}^{m}$ for the same $C_{0}$, thus for all $m$ such $mk\leq T$. Firstly, by equations (\ref{eq:er1}) to (\ref{eq:er6}) we have $$ \vert (\sigma_{h}^{0}-\sigma(0))_{u_{h}^{0}}\vert+\vert u_{h}^{0} -u(0)\vert\leq Mh^{2}\{C_{3}(1+\lambda M\delta_{0}h^{-1})+C_{1}\} $$ and $$ [k\{\vert\sigma_{h}^{0}-\sigma(0)\vert^{2}+\vert D(u_{h}^{0}-u(0) \vert^{2}\}]^{1/2}\leq \sqrt{2k}(C_{1}+C_{3})Mh^{2}. $$ To ensure that $(\sigma_{h}^{0},u_{h}^{0})=(\tilde{\sigma}_{0},\tilde{u}_{0})\in B_{h,k}^{0}$, it suffices to impose, for $h