]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update the discussion about broken RT space in step-61. 8984/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 1 Nov 2019 00:38:48 +0000 (18:38 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 1 Nov 2019 04:01:25 +0000 (22:01 -0600)
examples/step-61/doc/intro.dox
examples/step-61/step-61.cc

index 491b44f65df8c7055dd95fc5fba929e2a72bbef9..65e017a7f5e0dea0b9f3825474a4823a1cfe50e3 100644 (file)
@@ -49,7 +49,7 @@ as easy to implement. However, as we will see in the following, this
 additional effort is not prohibitive.
  
 
-<h3>  Weak Galerkin finite element methods</h3>
+<h3> Weak Galerkin finite element methods </h3>
 
 Weak Galerkin Finite Element Methods (WGFEMs) use discrete weak functions 
 to approximate scalar unknowns, and discrete weak gradients to 
@@ -132,7 +132,7 @@ form above. This works for the case where we have to evaluate the
 test function $q_h$ on the boundary (where we would simply take its interface
 part $q_h^\partial$) but we have to be careful with the gradient because
 that is only defined in cell interiors. Consequently,
- the weak Galerkin scheme for the Poisson equation is defined by
+the weak Galerkin scheme for the Poisson equation is defined by
 @f{equation*}
 \mathcal{A}_h\left(p_h,q \right) = \mathcal{F} \left(q_h \right),
 @f}
@@ -150,7 +150,7 @@ and
 @f}
 The key point is that here, we have replaced the gradient $\nabla p_h$ by the
 <i>discrete weak gradient</i> operator
- $ \nabla_{w,d} p_h $ that makes sense for our peculiarly defined approximation $p_h$.
+$\nabla_{w,d} p_h$ that makes sense for our peculiarly defined approximation $p_h$.
  
 The question is then how that operator works. For this, let us first say how we
 think of the discrete approximation $p_h$ of the pressure. As mentioned above,
@@ -173,7 +173,8 @@ $\nabla_{w,d} p_h$ would simply be the exact gradient $\nabla p_h$. But, since
 $p_h|_K$ is not continuous between interior and boundary of $K$, we need a more
 general definition; furthermore, we can not deal with arbitrary functions, and
 so require that $\nabla_{w,d} p_h$ is also in a finite element space (which, since
-the gradient is a vector, has to be vector-valued).
+the gradient is a vector, has to be vector-valued, and because the weak gradient
+is defined on each cell separately, will also be discontinuous between cells).
 
 The way this is done is to define this weak gradient operator $\nabla_{w,d}|_K :
 DGQ_k(K) \times DGQ_r(\partial K) \rightarrow RT_s(K)$ (where $RT_s(K)$ is the
@@ -202,7 +203,17 @@ are equal for all test functions.
   Rather, it is in a "broken" Raviart-Thomas space that below we will
   represent by the symbol $DGRT_s$. (The term "broken" here refers to
   the process of "breaking something apart", and not to the synonym to
-  the expression "not functional".)
+  the expression "not functional".) One might therefore (rightfully) argue that
+  the notation used in the weak Galerkin literature is a bit misleading,
+  but as so often it all depends on the context in which a certain
+  notation is used -- in the current context, references to the
+  Raviart-Thomas space or element are always understood to be to the
+  "broken" spaces.
+
+@note deal.II happens to have an implementation of this broken Raviart-Thomas
+  space: The FE_DGRT class. As a consequence, in this tutorial we will simply
+  always use the FE_DGRT class, even though in all of those places where
+  we have to compute cell-local matrices and vectors, it makes no difference.
 
 
 <h3> Representing the weak gradient </h3>
index 169af388f9e6512ebeb46da2b86e68b97edad3ba..32bf8dcdbfe9702a404074d14d3dfb35b310ab0c 100644 (file)
@@ -362,22 +362,24 @@ namespace Step61
   // exactly the kind of information and operation provided by the
   // DoFHandler class.
   //
-  // On the other hand, we don't have such a DoFHandler object for the
-  // Raviart-Thomas space in this program. In fact, we don't even have
-  // an element that can represent the "broken" Raviart-Thomas space
-  // we really want to use here (i.e., the restriction of the
-  // Raviart-Thomas shape functions to individual cells, without the
-  // need for any kind of continuity across cell interfaces). We solve
-  // this conundrum by using the fact that one can call
+  // We could create a DoFHandler object for the "broken" Raviart-Thomas space
+  // (using the FE_DGRT class), but we really don't want to here: At
+  // least in the current function, we have no need for any globally defined
+  // degrees of freedom associated with this broken space, but really only
+  // need to reference the shape functions of such a space on the current
+  // cell. As a consequence, we use the fact that one can call
   // FEValues::reinit() also with cell iterators into Triangulation
   // objects (rather than DoFHandler objects). In this case, FEValues
   // can of course only provide us with information that only
-  // references information of cells, rather than degrees of freedom
+  // references information about cells, rather than degrees of freedom
   // enumerated on these cells. So we can't use
   // FEValuesBase::get_function_values(), but we can use
   // FEValues::shape_value() to obtain the values of shape functions
   // at quadrature points on the current cell. It is this kind of
-  // functionality we will make use of below.
+  // functionality we will make use of below. The variable that will
+  // give us this information about the Raviart-Thomas functions below
+  // is then the `fe_values_rt` (and corresponding `fe_face_values_rt`)
+  // object.
   //
   // Given this introduction, the following declarations should be
   // pretty obvious:

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.