]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Modified dox and cc. Now they are consistent. 10926/head
authorIgnacio Tomas (-EXP) <itomas@sandia.gov>
Wed, 23 Sep 2020 03:51:29 +0000 (21:51 -0600)
committerIgnacio Tomas (-EXP) <itomas@sandia.gov>
Wed, 23 Sep 2020 03:51:29 +0000 (21:51 -0600)
examples/step-69/doc/intro.dox
examples/step-69/step-69.cc

index 4cd6c2fb9f04c368da9f291e92a84db5a5d27d1f..2032eb929761805e95755da3226b213d06b6d3ae 100644 (file)
@@ -163,9 +163,9 @@ $\mathbf{u}(\mathbf{x},t)$ remains in $\mathcal{B}$.
 
 <h4>Variational versus collocation-type discretizations</h4>
 
-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
+<!-- In particular, tutorial Step-67 focuses on Euler's equation of gas
+dynamics in the subsonic regime using dG techniques. -->
+
+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 <i>interactions between nodes</i> (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 <i>graph viscosity</i>. The graph viscosity serves as a
+    called <i>graph-viscosity</i>. 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 <i>edge-based</i> or
@@ -371,7 +374,7 @@ application of this kind of schemes, also called <i>edge-based</i> 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 <i>nodal normal</i>
+$\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.
index f513ab165f550fa6ace81a7e4e0d55518286ca69..bc7ad26450ef48e6e91397cd900f1ced5785003d 100644 (file)
@@ -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
+  // <code>OfflineData<dim>::BoundaryNormalMap</code>. We will normalize these
+  // vectors in a posterior loop.
 
   template <int dim>
   void OfflineData<dim>::assemble()

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.