From: Ignacio Tomas (-EXP) Date: Wed, 23 Sep 2020 03:51:29 +0000 (-0600) Subject: Modified dox and cc. Now they are consistent. X-Git-Tag: v9.3.0-rc1~1080^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8c4eb3ee6af75fc6c1c214d30e12e70f3db18bc3;p=dealii.git Modified dox and cc. Now they are consistent. --- diff --git a/examples/step-69/doc/intro.dox b/examples/step-69/doc/intro.dox index 4cd6c2fb9f..2032eb9297 100644 --- a/examples/step-69/doc/intro.dox +++ b/examples/step-69/doc/intro.dox @@ -163,9 +163,9 @@ $\mathbf{u}(\mathbf{x},t)$ remains in $\mathcal{B}$.

Variational versus collocation-type discretizations

-Following Step-9, Step-12, and Step-33, at this point it might look tempting -to base a discretization of Euler's equations on a (semi-discrete) variational -formulation: +Following Step-9, Step-12, Step-33, and Step-67, at this point it might look +tempting to base a discretization of Euler's equations on a (semi-discrete) +variational formulation: @f{align*} (\partial_t\mathbf{u}_{h},\textbf{v}_h)_{L^2(\Omega)} - ( \mathbb{f}(\mathbf{u}_{h}) ,\text{grad} \, \textbf{v}_{h})_{L^2(\Omega)} @@ -193,7 +193,10 @@ techniques (see the early work of @cite Brooks1982 and @cite Johnson1986). They have proven to be some of the best approaches for simulations in the subsonic shockless regime and similarly benign situations. -However, in the transonic and supersonic regime, and shock-hydrodynamics + + +However, in the transonic and supersonic regimes, and shock-hydrodynamics applications the use of variational schemes might be questionable. In fact, at the time of this writing, most shock-hydrodynamics codes are still firmly grounded on finite volume methods. The main reason for failure of @@ -244,7 +247,7 @@ vector-valued DoFHandler (i.e. initialized with an FESystem) it is possible to compute the block-matrices (required in order to advance the solution) with relative ease. However, the interactions that have to be computed in the context of time-explicit collocation-type schemes (such as finite -differences and/or the scheme presented in this tutorial) can be are +differences and/or the scheme presented in this tutorial) can be better described as interactions between nodes (not between DOFs). In addition, in our case we do not solve a linear equation in order to advance the solution. This leaves very little reason to use vector-valued @@ -282,7 +285,7 @@ where (\mathbf{U}_i^{n},\mathbf{U}_j^{n}, \textbf{n}_{ij}), \lambda_{\text{max}} (\mathbf{U}_j^{n}, \mathbf{U}_i^{n}, \textbf{n}_{ji}) \} \|\mathbf{c}_{ij}\|$ if $i \not = j$ is the so - called graph viscosity. The graph viscosity serves as a + called graph-viscosity. The graph-viscosity serves as a stabilization term, it is somewhat the discrete counterpart of $\epsilon \Delta \mathbf{u}$ that appears in the notion of viscosity solution described above. We will base our construction of $d_{ij}$ on @@ -363,7 +366,7 @@ compute all $d_{ij}$ in a separate step prior to actually performing above update. The core principle remains unchanged, though: we do not loop over cells but rather over all edges of the sparsity graph. -@note It is not uncommon to encounter such fully algebraic schemes (i.e. +@note It is not uncommon to encounter such fully-algebraic schemes (i.e. no bilinear forms, no cell loops, and no quadrature) outside of the finite element community in the wider CFD community. There is a rich history of application of this kind of schemes, also called edge-based or @@ -371,7 +374,7 @@ application of this kind of schemes, also called edge-based or @cite Rainald2008 for a historical overview). However, it is important to highlight that the algebraic structure of the scheme (presented in this tutorial) and the node-loops are not just a performance gimmick. Actually, the -structure of this scheme was borne out of theoretical necessity: the proof of +structure of this scheme was born out of theoretical necessity: the proof of pointwise stability of the scheme hinges on the specific algebraic structure of the scheme. In addition, it is not possible to compute the algebraic viscosities $d_{ij}$ using cell-loops since they depend nonlinearly on @@ -413,12 +416,17 @@ For the case of the reflecting boundary conditions we will proceed as follows: @f{align*} \mathbf{m}_i \dealcoloneq \mathbf{m}_i - (\widehat{\boldsymbol{\nu}}_i \cdot \mathbf{m}_i) \widehat{\boldsymbol{\nu}}_i \ \ - \text{for all }\mathbf{x}_i \in \partial\Omega^r + \text{where} \ \ + \widehat{\boldsymbol{\nu}}_i \dealcoloneq + \frac{\int_{\partial\Omega} \phi_i \widehat{\boldsymbol{\nu}} \, + \, \mathrm{d}\mathbf{s}}{\big|\int_{\partial\Omega} \phi_i + \widehat{\boldsymbol{\nu}} \, \mathrm{d}\mathbf{s}\big|} + \ \ \text{for all }\mathbf{x}_i \in \partial\Omega^r + \ \ \ \ \boldsymbol{(1)} @f} - that removes the normal component of $\mathbf{m}$. Here the definition of - nodal normal $\widehat{\boldsymbol{\nu}}_i$ is very much arbitrary (there is - no unique definition) but it should be consistent upon refinement with the - underlying geometry. + that removes the normal component of $\mathbf{m}$. This is a somewhat + naive idea that preserves a few fundamental properties of the PDE as we + explain below. This is approach is usually called "explicit treatment of boundary conditions". The well seasoned finite element person might find this approach questionable. @@ -452,26 +460,37 @@ integrate in space and time $\int_{\Omega}\int_{t_1}^{t_2}$ we would obtain @f{align*} \int_{\Omega} \rho(\mathbf{x},t_2) \, \mathrm{d}\mathbf{x} = \int_{\Omega} \rho(\mathbf{x},t_1) \, \mathrm{d}\mathbf{x} \ , \ \ -\int_{\Omega} E(\mathbf{x},t_2) \, \mathrm{d}\mathbf{x} = -\int_{\Omega} E(\mathbf{x},t_1) \, \mathrm{d}\mathbf{x} \ , \ \ -\int_{\Omega} \mathbf{m}(\mathbf{x},t_2) \, \mathrm{d}\mathbf{x} = -\int_{\Omega} \mathbf{m}(\mathbf{x},t_1) \, \mathrm{d}\mathbf{x} +\int_{\Omega} \mathbf{m}(\mathbf{x},t_2) \, \mathrm{d}\mathbf{x} + \int_{t_1}^{t_2} \! \int_{\partial\Omega} p \boldsymbol{\nu} \, -\mathrm{d}\mathbf{s} \mathrm{d}t\, . +\mathrm{d}\mathbf{s} \mathrm{d}t = +\int_{\Omega} \mathbf{m}(\mathbf{x},t_1) \, +\mathrm{d}\mathbf{x} \ , \ \ +\int_{\Omega} E(\mathbf{x},t_2) \, \mathrm{d}\mathbf{x} = +\int_{\Omega} E(\mathbf{x},t_1) \, \mathrm{d}\mathbf{x} \ \ \ \ +\boldsymbol{(2)} @f} -Note that momentum is not a conserved quantity (interaction with walls leads to -momentum gain/loss). Even though we will not use reflecting boundary conditions -in the entirety of the domain, we would like to know that our implementation of +Note that momentum is NOT a conserved quantity (interaction with walls leads to +momentum gain/loss): however $\mathbf{m}$ has to satisfy a momentum balance. +Even though we will not use reflecting boundary conditions in the entirety of +the domain, we would like to know that our implementation of reflecting boundary conditions is consistent with the conservation properties mentioned -above. In order to guarantee such conservation property it is necessary to -modify the values of the vectors $\mathbf{c}_{ij}$ as follows +above. In particular, if we use the projection $\boldsymbol{(1)}$ in the +entirety of the domain the following discrete mass-balance can be guaranteed: @f{align*} -\mathbf{c}_{ij} \dealcoloneq \mathbf{c}_{ij} + \int_{\partial \Omega^r} -(\boldsymbol{\nu}_j - \boldsymbol{\nu}(\mathbf{s})) \phi_i \phi_j \, -\mathrm{d}\mathbf{s} -\ \ \text{whenever both }\mathbf{x}_i\text{ and }\mathbf{x}_j \text{ lie on the -boundary} +\sum_{i \in \mathcal{V}} m_i \rho_i^{n+1} = +\sum_{i \in \mathcal{V}} m_i \rho_i^{n} \ , \ \ +\sum_{i \in \mathcal{V}} m_i \mathbf{m}_i^{n+1} ++ \tau_n \int_{\partial\Omega} \Big(\sum_{i \in \mathcal{V}} p_i^{n} \phi_i\Big) +\widehat{\boldsymbol{\nu}} \mathrm{d}\mathbf{s} = +\sum_{i \in \mathcal{V}} m_i \mathbf{m}_i^{n} \ , \ \ +\sum_{i \in \mathcal{V}} m_i E_i^{n+1} = \sum_{i \in \mathcal{V}} m_i +E_i^{n} \ \ \ \ +\boldsymbol{(3)} @f} -This consistent modification of the $\mathbf{c}_{ij}$ is a direct consequence -of simple integration by parts arguments, see page 12 of @cite GuermondEtAl2018 -for more details. +where $p_i$ is the pressure at the nodes that lie at the boundary. Clearly +$\boldsymbol{(3)}$ is the discrete counterpart of $\boldsymbol{(2)}$. The +proof of identity $\boldsymbol{(3)}$ is omitted, but we briefly mention that +it hinges on the definition of the nodal normal +$\widehat{\boldsymbol{\nu}}_i$ provided in $\boldsymbol{(1)}$. We also note that +this enforcement of reflecting boundary conditions is different from the one +originally advanced in @cite GuermondEtAl2018. diff --git a/examples/step-69/step-69.cc b/examples/step-69/step-69.cc index f513ab165f..bc7ad26450 100644 --- a/examples/step-69/step-69.cc +++ b/examples/step-69/step-69.cc @@ -992,20 +992,14 @@ namespace Step69 // // @f{align*} // \widehat{\boldsymbol{\nu}}_i \dealcoloneq - // \frac{\boldsymbol{\nu}_i}{|\boldsymbol{\nu}_i|} \ \text{ where } - // \boldsymbol{\nu}_i \dealcoloneq \sum_{T \subset \text{supp}(\phi_i)} - // \sum_{F \subset \partial T \cap \partial \Omega} - // \sum_{\mathbf{x}_{q,F}} \nu(\mathbf{x}_{q,F}) - // \phi_i(\mathbf{x}_{q,F}). + // \frac{\int_{\partial\Omega} \phi_i \widehat{\boldsymbol{\nu}} \, + // \, \mathrm{d}\mathbf{s}}{\big|\int_{\partial\Omega} \phi_i + // \widehat{\boldsymbol{\nu}} \, \mathrm{d}\mathbf{s}\big|} // @f} // - // Here $T$ denotes elements, - // $\text{supp}(\phi_i)$ the support of the shape function $\phi_i$, - // $F$ are faces of the element $T$, and $\mathbf{x}_{q,F}$ - // are quadrature points on such face. Note that this formula for - // $\widehat{\boldsymbol{\nu}}_i$ is nothing else than some form of - // weighted averaging. Other more sophisticated definitions for $\nu_i$ - // are possible but none of them have much influence in theory or practice. + // We will compute the numerator of this expression first and store it in + // OfflineData::BoundaryNormalMap. We will normalize these + // vectors in a posterior loop. template void OfflineData::assemble()