\[ \renewcommand{\vec}[1]{\boldsymbol{#1}} \renewcommand{\tilde}[1]{\widetilde{#1}} \newcommand{\N}{\mathcal{N}} \newcommand{\R}{\mathbb{R}} \newcommand{\T}{\mathcal{T}} \newcommand{\dd}{\,\mathrm{d}} \newcommand{\dt}{\frac{\mathrm{d}}{\mathrm{d}t}} \newcommand{\md}{\partial^\bullet} \newcommand{\abs}[1]{\left|#1\right|} \newcommand{\norm}[1]{\|#1\|} \newcommand{\diam}{\mathop{\rm diam}\nolimits} \]
Tom Ranner
joint work with Charles M. Elliott and Pierre Stepanov
Slides available at tomranner.org/fbp2021
Today: error analysis
Given initial condition \(u_0\), (possibly discontinuous) diffusion coefficient \(A\), and data \(f_1, f_2, g\), find \(u \colon \Omega \to \R\) such that
\[ \begin{aligned} u_t - \nabla \cdot( A|_{\Omega_i(t)} \nabla u ) & = f_i && \mbox{ in } \Omega_i(t) \\ [A \partial_\nu u ] = g \quad \mbox{ and } \quad [u] & = 0 && \mbox{ on } \Gamma(t) \\ u & = 0 && \mbox{ on } \partial \Omega \\ u(\cdot, 0) & = u_0 && \mbox{ in } \Omega. \end{aligned} \]
\([\eta]\) denotes jump in \(\eta\) across \(\Gamma(t)\).
We assume the (nice, smooth) flow of the domain given by the flow \(\Phi^t \colon \Omega \to \Omega\).
We assume the flow preserves partition of the domain so that \(\Phi^t(\Gamma_0) = \Gamma(t)\).
We define a global velocity field \(\vec{w}\) by \[ \dt \Phi_t(x) = w(\Phi_t(x), t) \qquad \mbox{ for all } x \in \Omega. \]
We identify a pair \(v = (v_1, v_2)\) for a function \(v\) defined on \(\Omega\) by \(v_i = v|_{\Omega_i}\).
We define the spaces: \[ \begin{aligned} H(t) & = L^2(\Omega_1(t)) \times L^2(\Omega_2(t)) \\ V(t) & = \{ v \in H^1(\Omega_1(t)) \times H^1(\Omega_2(t)) \,|\, u_1|_{\Gamma(t)} = u_2|_{\Gamma(t)} \mbox{ and } u_2|_{\partial \Omega} = 0 \}. \end{aligned} \]
\(V(t) \subset H(t) \subset V^*(t)\) form an evolving Hilbert triple.
Each spaces forms a compatible pair with the push-forward (or pull-back) map given by \[ \begin{aligned} \phi_t v &= (v_1(\Phi_{-t}), v_2(\Phi_{-t})) && \mbox{ for } v \in H(0) \\ \phi_{-t} v &= (v_1(\Phi_{t}), v_2(\Phi_{t})) && \mbox{ for } v \in H(t). \end{aligned} \]
This allows us to define moving Sobolev spaces (for \(X = V, H\)) \[ L^2_X := \{ v \colon [0,T] \to \bigcup_{t \in [0,T]} X(t) \times \{ t \} \,|\, t \mapsto (\hat{v}(t), t), \phi_{-t} \hat{v}(t) \in L^2(0, T; X_0) \} \]
and a material derivative for smooth function \(v\) by \[ \md_t v := \phi_t \left( \dt \phi_{-t} v \right) = \partial_t v + \vec{w} \cdot \nabla \eta. \]
Let \(f \in L^2_{V^*}\) and \(g \in L^2_{H^{-1/2}(\Gamma)}\), find \(u \in \{ v \in L^2_V : \md_t v \in L^2_{V^*}\}\) such that \[ \begin{aligned} \int_0^T \langle \md_t u, v \rangle_{V(t)} \dd t + \underbrace{\sum_i \int_0^T \int_{\Omega_i(t)} A_i \nabla u_i \cdot \nabla v_i - w \cdot \nabla u_i v_i \dd x}_{ =: a(t; u, v)} \dd t \\ = \int_0^T \underbrace{\langle f, v \rangle_{V(t)} + \langle g, v \rangle_{H^{1/2}(\Gamma(t))}}_{ =: l(t; v)} \dd t, \end{aligned} \] with the initial condition \(u(0) = u_0 \in H_0\).
The variational form has a unique solution \(u\) that satisfies \[ \norm{u}_{L^2_V}^2 + \norm{\md_t u}_{L^2_{V^*}}^2 \le C \left( \norm{f}_{L^2_{V^*}}^2 + \norm{g}_{L^2_{H^{-1/2}(\Gamma)}}^2 + \norm{u_0}_{H_0}^2 \right). \]
Suppose further that \(f \in L^2_H\), \(g \in L^2_{H^{1/2}(\Gamma)}\) and \(u_0 \in V(0)\) then the solution \(u \in \{v \in L^2_V, \md_t v \in L^2_H\}\) and satisfies \[ \norm{u}_{L^2_V}^2 + \norm{\md_t u}_{L^2_{H}}^2 \le C \left( \norm{f}_{L^2_{H}}^2 + \norm{g}_{L^2_{H^{1/2}(\Gamma)}}^2 + \norm{u_0}_{V_0}^2 \right). \]
We construct a moving isoparametric finite element method of order \(k\).
We denote by \(\tilde{\Gamma}_h\) the initial approximation of \(\Gamma_0\).
We construct a bijection \(\Psi_h\colon \Omega \to \Omega\) such that
\[ \Psi_h(x) = x + (\mu^*(x))^{\textcolor{green}{k+2}} (p(y(x), 0) - y(x)) \quad \mbox{ if } \mu^*(x) \neq 0 \]
The initial mesh \(F_{K}(\hat{x}, 0)\) is defined by the interpolation of \(\Psi_h\) of order \(k\)
We move each Lagrange node by the Lagrange interpolation (order \(k\)) of the smooth velocity \(w\).
We assume that this gives a good mesh.
On going work - how to generalise/relax this condition.
We construct a bijection between the computational domain and the exact domain. We call this the lift denoted by \(\Lambda_h \colon \Omega \to \Omega\).
We can also lift functions defined on the triangulation: \(v_h \mapsto v_h^\ell\).
The process also defines a lifted material derivative. (Dziuk & Elliott, 2012)
The procedure is similar to construction of the initial domain.
Important: The lift is the identity mapping for all elements away from the interface.
\[ \begin{aligned} S_h^i(t) := \Big\{ & \chi_h = (\chi_K)_{K \in \T_h^i(t)} \in \prod_{K \in \T_h^i(t)} \{ \hat\chi \circ F^{-1}_{K} : \hat\chi \in P_k(\hat{K}) \} : \\ & \chi_K(a(t)) = \chi_{K'}(a(t)) \mbox{ for all } K, K' \in \T^i(a(t)), \\ & \qquad \mbox{ for all } a(t) \in \N_h^i(t) \Big\} \end{aligned} \]
\[ \begin{aligned} S_h(t) := \big\{ & \chi_h = (\chi^1_h, \chi^2_h) \in S_h^1(t) \times S_h^2(t) : \\ & \chi_h^1(a(t)) = \chi_h^2(a(t)) \mbox{ for all } a(t) \in \N_h^\Gamma(t), \\ & \chi_h^2(a) = 0 \mbox{ for all } a \in \N_h^{\partial \Omega} \big\} \end{aligned} \]
Given \(f_h\), \(g_h\), diffusion coefficients \(A_h^i\) and initial data \(U_{h,0}\), find \(U_h \in \{ v_h \in L^2_{S_h} : \md_h v_h \in L^2_{S_h}\}\) such that \[ \begin{aligned} \dt \int_{\Omega(t)} U_h v_h \dd x + \underbrace{\sum_i \int_{\Omega_{h}^i(t)} A_{h}^i \nabla U_{h}^i \cdot \nabla v_h - W_h \cdot \nabla U_h^i v_h \dd x}_{ =: a_h(t; U_h, v_h)} \\ = \underbrace{\int_{\Omega_h(t)} f_h v_h \dd x + \int_{\Gamma_h(t)} g_h v_h \dd \sigma_h}_{ =: l_h(t; v_h)} + \int_{\Omega_h(t)} U_h \md_h v_h \end{aligned} \] with the initial condition \(U_h(0) = U_{h,0}\).
The Discrete problem has a unique solution \(U_h \in C^1_{S_h}\) and \[ \sup_{t \in [0,T]} \norm{U_h}_{H_h(t)}^2 + \int_0^T \norm{U_h}_{V_h(t)}^2 \dd t \le C(T) \norm{U_{h,0}}_{H_h(0)}^2. \]
\[ \begin{aligned} & \sup_{t \in [0,T]} \norm{u - U_h^\ell}_{H(t)}^2 + h^2 \int_0^T \norm{u - U_h^\ell}_{V(t)}^2 \\ & \qquad \le \norm{u_0 - U_{h,0}^\ell}_{H(t)}^2 + c h^{2k+2} C(u). \end{aligned} \]
We apply the ideas from (Elliott & TR, 2020).
We use the Ritz projection \(\Pi_h \colon V(t) \to S^h(t)\) given by \[ a_h(t; \Pi_h z, v_h) = a(t; z, v_h^\ell) \qquad \mbox{ for all } v_h \in S^h(t). \]
Split errors: \[ u - U_h^\ell = \rho + \theta = (u - (\Pi_h u)^\ell) + (\Phi_h u - U_h)^\ell. \]
Estimates for \(\rho\) and \(\theta\) both require geometric estimates, interpolation estimates, duality arguments.
The geometric errors compare a lifted problem (with lifted material derivative) to the continuous problem.
We want to show: \[ \tag{P1} \abs{m(t; \eta_h^\ell, v_h^\ell) - m_h(t; \eta_h, v_h)} \le c h^{k+1} \norm{\eta_h^\ell}_{V(t)} \norm{v_h^\ell}_{V(t)}. \]
where
\[ \begin{aligned} m(t; \eta, v) &= \sum_i \int_{\Omega^i} \eta^i v^i \dd x \\ m_h(t; \eta_h, v_h) &= \sum_i \int_{\Omega_h^i} \eta_h^i v_h^i \dd x. \end{aligned} \]
\[ \tag{P1} \abs{m(t; \eta_h^\ell, v_h^\ell) - m_h(t; \eta_h, v_h)} \le c h^{k+1} \norm{\eta_h^\ell}_{V(t)} \norm{v_h^\ell}_{V(t)}. \]
Let \(J_h\) be the Jacobian \(\sqrt{\det( \nabla \Lambda_h^T \nabla \Lambda h)}\) then \[ \sup_{t \in [0, T]} \norm{J_h - 1}_{L^\infty(\Omega_{h,i}(t))} \le c h^k. \]
\[ \begin{aligned} & \abs{m(t; \eta_h^\ell, v_h^\ell) - m_h(t; \eta_h, v_h)} = \abs{\int_{\Omega} \eta_h^\ell v_h^\ell (J_h - 1)} \\ & \qquad \le c h^k \norm{ \eta_h }_{L^2(\Omega)} \norm{ v_h }_{L^2(\Omega)}. \end{aligned} \]
Denote by \([J_j \neq 1]\) be the elements where \(J_h \neq 1\) then \[ \begin{aligned} & \abs{m(t; \eta_h^\ell, v_h^\ell) - m_h(t; \eta_h, v_h)} = \abs{\int_{\Omega} \eta_h^\ell v_h^\ell (J_h - 1)} \\ & \le c h^k \left(\int_{[J_h \neq 1]} (\eta_h^\ell)^2 \right)^{1/2} \left(\int_{[J_h \neq 1]} (v_h^\ell)^2 \right)^{1/2}. \end{aligned} \]
Applying the narrow band trace inequality(Elliott & TR, 2012) we get the result: For any \(\eta \in H^1(\Omega)\) \[ \norm{ \eta }_{L^2([J_h \neq 1])} \le c h^{1/2} \norm{ \eta }_{H^1(\Omega)}. \]
\[ \abs{m(t; \eta_h^\ell, v_h^\ell) - m_h(t; \eta_h, v_h)} \le c h^{k+1} \norm{ \eta_h }_{V(t)} \norm{ v_h }_{V(t)}. \]
\[ \tag{P2} \abs{a(t; \eta_h^\ell, v_h^\ell) - a_h(t; \eta_h, v_h)} \le c h^{k+1} \norm{\eta_h^\ell}_{Z_0(t)} \norm{v_h^\ell}_{Z_0(t)}. \]
We can apply the same ideas as for \(m\) bilinear forms to get order \(k+1\) estimate.
But this needs careful application to ensure that we can apply this for smooth functions.
We do this by using our Ritz projection so that (P2) is only applied to smooth functions.
Define the bilinear form \(b \colon V(t) \times V(t) \to \R\) by \[ \dt a(t; \eta, v) = a(t; \md_t \eta, v) + a(t; \eta, \md_t v) + b(t; \eta, v). \]
Let \(z \in Z_k(t)\) and \(\eta = z - (\Pi_h z)^\ell\). We need to show \[ \abs{b(t; \eta, \zeta)} \le c \left( \norm{\eta}_{H(t)} + h \norm{\eta}_{V(t)} + h^{k+1} \norm{z}_{Z_k(t)} \right) \norm{\zeta}_{Z_0(t)}. \tag{B3} \]
To show this we use:
Implemented using the high-level finite element package firedrake
.
Allows quick implementation using ufl
but with efficient, parallel implementation using PETSc
.
Use appropriate order BDF scheme for time stepping (including mesh movement).
We have presented a numerical method and accompanying analysis for a moving coupled bulk-surface problem.
The method and analysis is based on the fundamental and abstract frameworks from (Elliott & TR, 2020).
Numerical results demonstrate the convergence.
Questions remain over the quality of the moving mesh.
Thank you for your attention!
C.M. Elliott & TR. A unified theory for continuous-in-time evolving finite element space approximations to partial differential equations in evolving domains, IMA Journal of Numerical Analysis, Volume 41, Issue 3, July 2021, Pages 1696–1845, https://doi.org/10.1093/imanum/draa062