\section{The refinement and coarsening implementation}\label{S:rc-impl}

\subsection{The refinement routines}%
\label{S:refinement_routines}%
\idx{mesh refinement!implementation|(}

For the refinement of a mesh the following symbolic constant is defined
and the refinement is done by the functions
\fdx{refine()@{\code{refine()}}}%
\fdx{global_refine()@{\code{global\_refine()}}}
\cdx{MESH_REFINED@{\code{MESH\_REFINED}}}
\bv\begin{verbatim}
#define MESH_REFINED   1

U_CHAR refine(MESH *mesh, FLAGS fill_flags);
U_CHAR global_refine(MESH *mesh, int n_bisections, FLAGS fill_flags);
\end{verbatim}\ev
\begin{compatibility}
  In previous versions, the last parameter \code{fill\_flags} was
  missing. To obtain the old behaviour, \code{FILL\_NOTHING} should be
  passed for the parameter \code{fill\_flags}.
\end{compatibility}
Description:
\begin{descr}
\kitem{refine(mesh, fill\_flags)} refines all leaf elements with a
  \emph{positive} element marker \code{mark} times (this mark is
  usually set by some adaptive procedure); the routine loops over all
  leaf elements and refines the elements with a positive marker until
  there is no element left with a positive marker; the return value is
  \code{MESH\_REFINED}, if at least one element was refined, and
  \code{0} otherwise. Every refinement has to be done via this
  routine. The basic steps of this routine are described below.
  \begin{description}
  \item[Parameters]\hfill
    \begin{descr}
    \kitem{mesh} The mesh which will be refined, possibly.
    \kitem{FLAGS fill\_flags} Request additional data filled in during
      the mesh-traversal, useful for custom mesh-adaptation call-back
      in DOF-vectors and -matrices. Additionally -- if any
      mesh-adaptation call-backs have been registered -- then the set
      of fill-flags will be augmented by the requirements of the
      related basis-function sets.
    \end{descr}
  \item[Return Value] Either \code{MESH\_REFINED}, if at least one
    element has been sub-divided, or \code{0} otherwise.
  \end{description}

  \hrulefill

\kitem{global\_refine(mesh, num\_bisections, fill\_flags)} sets all element
  markers for leaf elements of \code{mesh} to \code{mark}; the mesh is
  then refined by \code{refine()} which results in a \code{mark}
  global refinement of the mesh; the return value is
  \code{MESH\_REFINED}, if \code{mark} is positive, and \code{0}
  otherwise.
  \begin{description}
  \item[Parameters]\hfill
    \begin{descr}
    \kitem{MESH *mesh} The mesh which will be refined, possibly.
    \kitem{num\_bisections} The number of bisections to perform on
      each element. This is an upper limit: the resulting mesh will
      have no ``hanging-nodes'', this conformal closure may require
      more refinement steps than requested by \code{num\_bisections}.
    \kitem{FLAGS *fill\_flags} Request additional data filled in
      during the mesh-traversal, useful for custom mesh-adaptation
      call-back in DOF-vectors and -matrices.
    \end{descr}
  \item[Return Value] Either \code{MESH\_REFINED}, if at least one
    element has been sub-divided, or \code{0} otherwise.
  \end{description}       
\end{descr}

\subsubsection{Basic steps of the refinement algorithm}%
\label{S:basic_steps_rc}

The refinement of a mesh is principally done in two steps --- each 
step corresponding to one mesh traversal. In the
first step no coordinate information is necessary, only a
topological refinement is performed. If new nodes are created that 
belong to a a then these can be projected in the second
step where coordinate information is calculated.

Again using the notion of ``refinement edge'' for the element itself
in 1d, the algorithm performs the following steps:
\begin{descr}
\item[1.] 
The whole mesh is refined only topologically. This part consists of
\begin{descr}
\item the collection of a compatible refinement patch; this includes
      the recursive refinement of adjacent elements with an incompatible 
      refinement edge;
\item the topological bisection of the patch elements;
\item the transformation of leaf data from parent to child, if
      such a function is available in the \code{leaf\_data\_info}
      structure;
\item allocation of new DOFs;
\item handing on of DOFs from parent to the children;
\item interpolation of DOF vectors from the coarse grid to 
      the fine one on the whole refinement patch, if the function
      \code{refine\_interpol()} is available for these DOF vectors
      (compare Section \ref{S:DOF_INTERPOL}); these routines must not
      use coordinate information;
\item a deallocation of DOFs on the parent when 
      \code{preserve\_coarse\_dofs == 0}, see \secref{S:access_fe_space}.
\end{descr}
This process is described in detail below.
\item[2.]  New nodes which belong to the curved part of the boundary
  are now projected onto the curved boundary via the
  \code{active\_projection()}%
  \idx{active_projection()@{\code{active\_projection()}}}
  function in the \code{EL\_INFO} structure. This entry is passed down
  from the corresponding macro element during the mesh traversal of this 
  step.
  The coordinates of the projected node are stored in a
  \code{REAL\_D}--vector and the pointers \code{el->new\_coord} of all
  parents \code{el} which belong to the refinement patch are set to
  this vector.
\end{descr}

The topological refinement is done by the recursive refinement
\algorithmref{book:A:recursive_refine}. In 1d, no recursion is needed.  In
2d and 3d, all elements at the refinement edge of a marked element are
collected. If a neighbour with an incompatible refinement edge is
found, this neighbour is refined first by a recursive call of the
refinement function. Thus, after looping around the refinement edge,
the patch of simplices at this edge is always a compatible
refinement patch.  The elements of this patch are stored in a vector
\code{ref\_list} with elements of type
\code{RC\_LIST\_EL}\ddx{RC_LIST_EL@{\code{RC\_LIST\_EL}}}, compare
\secref{S:rc_list_el}. This vector is an argument for the functions
for interpolation of DOF vectors during refinement, compare
\secref{S:DOF_INTERPOL}.

In 1d the vector has length \code{1}. In 2d the length is \code{2} if
the refinement edge is an interior edge; for a boundary edge the
length is \code{1} since only the element itself has to be refined.
For 1d and 2d, only the \code{el} entry of the components is set and
used.

In 3d this vector is allocated with length \code{mesh->max\_edge\_neigh}.
As mentioned in \secref{S:rc_list_el} we can define an orientation
of the edge and by this orientation we can define the right and left
neighbors (inside the patch) of an element at this edge.

The patch is bisected by first inserting a new vertex at the midpoint
of the refinement edge. Then all elements of the refinement patch are
bisected. This includes the allocation of new DOFs, the adjustment of
DOF pointers, and the memory allocation for leaf data (if initialized by
the user) and transformation of leaf data from parent to
child (if a pointer to a function \code{refine\_leaf\_data()}
\idx{leaf data!transformation during refinement} is provided by the
user in the \code{init\_leaf\_data()} call). Then memory for
parents' leaf data is freed and information stored there is definitely
lost.

In the case of higher order elements we also have to add new DOFs on the
patch and if we do not need information about the higher order DOFs
on coarser levels they are removed from the parents. 
There are some basic rules for adding and removing DOFs
which are important for the prolongation and restriction of data
(see \secref{S:DOF_INTERPOL}):%
\idx{DOFs!adding and removing of DOFs|(}
\begin{descr}
\item[1.] Only DOFs of the same kind (i.e. \code{VERTEX}, \code{EDGE}, or
\code{FACE}) and whose nodes have the same geometrical position on parent
and child are handed on to this child from the parent;
\item[2.] DOFs at a vertex, an edge or a face belong to all elements sharing
this vertex, edge, face, respectively;
\item[3.] DOFs on the parent are only removed if the entry 
\code{preserve\_coarse\_dofs} in the corresponding \code{DOF\_ADMIN} data 
structure is \false; in that case only DOFs which are not handed on to a 
child are removed on the parent.
\end{descr}

A direct consequence of 1. is that only DOFs inside the
patch are added or removed; DOFs on the patch boundary stay untouched.
\code{CENTER} DOFs can not be handed from parent to child since
the centers of the parent and the children are always at different
positions.

Using standard Lagrange finite elements, only DOFs that are not
handed from parent to child have to be set while interpolating
a finite element function to the finer grid; all values of the 
other DOFs stay the same (the same holds during coarsening and
interpolating to the coarser grid).

Due to 2. it is clear that DOFs shared by more than one
element have to be allocated only once and pointers to these DOFs are
set correctly for all elements sharing it.

Now, we take a closer look at DOFs that are handed on by the parents
and those that have to be allocated: In 1d we have
\idx{refinement!DOFs!handed from parent to children}
%
\bv\begin{verbatim}
child[0]->dof[0] = el->dof[0];
child[1]->dof[1] = el->dof[1];
\end{verbatim}\ev
%
in 2d
%
\bv\begin{verbatim}
child[0]->dof[0] = el->dof[2];
child[0]->dof[1] = el->dof[0];
child[1]->dof[0] = el->dof[1];
child[1]->dof[1] = el->dof[2];
\end{verbatim}\ev
%
In 3d for \code{child[1]} this passing of DOFs additionally depends
on the element type \code{el\_type} of the parent. For
\code{child[0]} we always have
%
\bv\begin{verbatim}
child[0]->dof[0] = el->dof[0];
child[0]->dof[1] = el->dof[2];
child[0]->dof[2] = el->dof[3];
\end{verbatim}\ev
%
For \code{child[1]} and a parent of type \code{0} we have
\bv\begin{verbatim}
child[1]->dof[0] = el->dof[1];
child[1]->dof[1] = el->dof[3];
child[1]->dof[2] = el->dof[2];
\end{verbatim}\ev
%
and for a parent of type \code{1} or \code{2}
%
\bv\begin{verbatim}
child[1]->dof[0] = el->dof[1];
child[1]->dof[1] = el->dof[2];
child[1]->dof[2] = el->dof[3];
\end{verbatim}\ev
%
In 1d
\bv\begin{verbatim}
child[0]->dof[1] = child[1]->dof[0]
\end{verbatim}\ev
and in 3d and 3d
\bv\begin{verbatim}
child[0]->dof[DIM] = child[1]->dof[DIM]
\end{verbatim}\ev
%
is the newly allocated DOF at the midpoint of the refinement edge
(compare Figure~\ref{book:F:node_numbering_2d} on
page~\pageref{book:F:node_numbering_2d} for the 1d and 2d situation
and Figure~\ref{book:F:node_numbering_3d} on
page~\pageref{book:F:node_numbering_3d} for the 3d situation).

In the case that we have DOFs at the midpoint of edges (only 2d and
3d) the following DOFs are passed on 
(let \code{enode = mesh->node[EDGE]} be the offset for DOFs at edges):
for 2d
%
\bv\begin{verbatim}
child[0]->dof[enode+2] = el->dof[enode+1];
child[1]->dof[enode+2] = el->dof[enode+0];
\end{verbatim}\ev
and for 3d
\bv\begin{verbatim}
child[0]->dof[enode+0] = el->dof[enode+1];
child[0]->dof[enode+1] = el->dof[enode+2];
child[0]->dof[enode+3] = el->dof[enode+5];
\end{verbatim}\ev
for \code{child[0]} a for \code{child[1]} of a parent of type \code{0}
\bv\begin{verbatim}
child[1]->dof[enode+0] = el->dof[enode+4];
child[1]->dof[enode+1] = el->dof[enode+3];
child[1]->dof[enode+3] = el->dof[enode+5];
\end{verbatim}\ev
and finally for \code{child[1]} of a parent of type \code{1} or \code{2}
\bv\begin{verbatim}
child[1]->dof[enode+0] = el->dof[enode+3];
child[1]->dof[enode+1] = el->dof[enode+4];
child[1]->dof[enode+3] = el->dof[enode+5];
\end{verbatim}\ev

\begin{figure}[htbp]
\centerline{\includegraphics[width=0.9\hsize]{EPS/ref_tetra_e}}
\caption[Mesh-refinement, maintenance of DOFs]{Edge DOFs that are
  freed $\bullet$, passed on $\circ$, and newly alloc ated $\Box$}
\label{F:edge_numbering_3d}
\end{figure}

We also have to create new DOFs (compare Figure~\ref{F:edge_numbering_3d}).
Two additional DOFs are created in the refinement
edge which are shared by all patch elements. Pointers to these
DOFs are adjusted for 
\idx{refinement!DOFs!newly created}
\bv\begin{verbatim}
child[0]->dof[enode+0],
child[1]->dof[enode+1]
\end{verbatim}\ev
in 2d and 
\bv\begin{verbatim}
child[0]->dof[enode+2],
child[1]->dof[enode+2]
\end{verbatim}\ev
in 3d for all patch elements.

In 3d, for each interior face of the refinement patch
there is a new edge where we have to add a new DOF vector.
These DOFs are shared by two children in the case of a boundary
face; otherwise it is shared by four children and pointers of 
\bv\begin{verbatim}
child[0]->dof[enode+4] = child[1]->dof[enode+{5,4,4}],
child[0]->dof[enode+5] = child[1]->dof[enode+{4,4,5}]
\end{verbatim}\ev
are adjusted for those elements.

In 3d, there may be also DOFs at faces; the face DOFs in the boundary of the
patch are passed on (let \code{fnode = mesh->node[FACE]} be the offset
for DOFs at faces):
\bv\begin{verbatim}
child[0]->dof[fnode+3] = el->dof[fnode+1];
child[1]->dof[fnode+3] = el->dof[fnode+0];
\end{verbatim}\ev
For the common face of \code{child[0]} and \code{child[1]} we have
to allocate a new face DOF vector which is located at 
\bv\begin{verbatim}
child[0]->dof[fnode+0] = child[1]->dof[fnode+0]
\end{verbatim}\ev
and finally for each interior face of the patch two new face DOF vectors
are created and pointers for adjacent children are adjusted:
\bv\begin{verbatim}
child[0]->dof[fnode+1],
child[0]->dof[fnode+2],
child[1]->dof[fnode+1],
child[1]->dof[fnode+2]
\end{verbatim}\ev
Each of these DOF vectors may be shared with another child of a patch element.

If DOFs are located at the barycenter they have to be allocated
for both children in 2d and 3d (let \code{cnode = mesh->node[CENTER]}
be the offset for DOFs at the center)
\bv\begin{verbatim}
child[0]->dof[cnode],
child[1]->dof[cnode].
\end{verbatim}\ev
  
After adding and passing on of DOFs on the patch we can interpolate
data from the coarse to the fine grid on the whole patch. This is an
operation on the whole patch since new DOFs can be shared by more than
one patch element and usually the value(s) of such a DOF should only
be calculated once.

\idx{refinement!interpolation of DOF vectors}
All DOF vectors and matrices having a pointer to a function
\code{refine\_interpol()}\idx{refine_interpol()@{\code{refine\_interpol()}}}
in the corresponding data structure are
interpolated to the fine grid. Such a function essentially depends
on the described passing on and new allocation of DOFs. 
An abstract description of such functions can be found in 
\secref{book:S:inter_restrict} and a more detailed one for Lagrange
elements in \secref{S:Lagrange_elements}.

After such an interpolation, DOFs of higher degree on parent elements
may no longer be of interest (when not using a higher order multigrid
method).

In such a case \code{DOF\_ADMIN.flags \& ADM\_PRESERVE\_COARSE\_DOFS}
should evaluate to $0$ and in this case all DOFs on the parent that
are not handed over to the children will be removed. The following
DOFs are removed on the parent for all patch elements (some DOFs are
shared by several elements): The DOFs at the center
%
\idx{refinement!DOFs!removed on the parent} 
\bv\begin{verbatim}
el->dof[mesh->node[CENTER]]
\end{verbatim}\ev
 are removed in all dimensions. In 2d, additionally DOFs in the 
refinement edge
\bv\begin{verbatim} 
el->dof[mesh->node[EDGE]+2]
\end{verbatim}\ev
are removed and in 3d the DOFs in the refinement edge and the DOFs in the
two faces adjacent to the refinement edge
\bv\begin{verbatim}
el->dof[mesh->node[EDGE]+0],
el->dof[mesh->node[FACE]+2],
el->dof[mesh->node[FACE]+3],
el->dof[mesh->node[CENTER]]
\end{verbatim}\ev
are deleted on the parent. Note that only the information about DOF
indices is deleted, the pointers \code{el->dof[n]},
\code{n}$\in\{$\code{0}$,\dots,$\code{ mesh->n\_node\_el-1}$\}$,
themselves remain available after refinement. This setting of DOF pointers and
pointers to children is the main part of the refinement module.
\idx{mesh refinement!implementation|)}%
\idx{DOFs!adding and removing of DOFs|)}

\subsection{The coarsening routines}%
\label{S:coarsening_routines}%
\idx{mesh coarsening!implementation|(}

For the coarsening of a mesh the following symbolic constant is defined
and the coarsening is done by the functions
\fdx{coarsen()@{\code{coarsen()}}}%
\fdx{global_coarsen()@{\code{global\_coarsen()}}}
\cdx{MESH_COARSENED@{\code{MESH\_COARSENED}}}
\bv\begin{verbatim}
#define MESH_COARSENED 2

U_CHAR coarsen(MESH *mesh, FILL_FLAGS fill_flags);
U_CHAR global_coarsen(MESH *mesh, int num_bisections, FILL_FLAGS flags_flags);
\end{verbatim}\ev
Description:
\begin{descr}
\kitem{coarsen(mesh, fill\_flags)} tries to coarsen all leaf element
  with a \emph{negative} element marker \code{|mark|} times (again,
  this mark is usually set by an adaptive procedure); the return value
  is \code{MESH\_COARSENED} if any element was coarsened, and \code{0}
  otherwise.
\kitem{global\_coarsen(mesh, n\_bisections, fill\_flags)} sets all
  element markers for leaf elements of \code{mesh} to
  \code{n\_bisections}; the mesh is then coarsened by
  \code{coarsen()}; depending on the actual distribution of coarsening
  edges on the mesh, this may not result in a \code{|n\_bisections|}
  global coarsening; the return value is \code{coarsen(mesh)} if
  \code{mark} is negative, and \code{0} otherwise.
\end{descr}

The function \code{coarsen()} implements \algorithmref{book:A:coarsen_mesh}.
For a marked element, the coarsening patch is collected first. This is
done in the same manner as in the refinement procedure. If such a
patch can \emph{definitely} not be coarsened (if one element of the
patch may not be coarsened, e.g.) all coarsening markers for all patch
elements are reset. If we can not coarsen the patch
immediately, because one of the elements has not a common coarsening
edge but is allowed to be coarsened more than once, then nothing is
done in the moment and we try to coarsen this patch later on (compare
\remarkref{book:R:coarsen}).

The coarsening of a patch is the ``inverse'' of the refinement of a
compatible patch. If DOF indices of the parents were deleted during
refinement, then new indices are now allocated. DOF pointers 
on the parents (\code{parent->dof[n]}) do not need to be 
touched, as they remain valid after refinement, see 
\secref{S:refinement_routines}.

If leaf data is stored at the pointer of \code{child[1]}, then
memory for the parent's leaf data is allocated. If a function
\code{coarsen\_leaf\_data} was provided during the call of
\code{init\_leaf\_data()} then leaf data is transformed
from children to parent. Finally, leaf data on both children is
freed.\idx{leaf data!transformation during coarsening}

\idx{coarsening!interpolation of DOF vectors}
\idx{coarsening!restriction of DOF vectors}
Like the interpolation of data during refinement, we now can 
restrict/interpolate data from children to parent. This is done by
the \code{coarse\_restrict()}%
\idx{coarse_restrict()@{\code{coarse\_restrict()}}} functions for all
those DOF vectors and matrices where such a function is available in the
corresponding data structure. Since it does not
make sense to both interpolate and restrict data,
\code{coarse\_restrict()} may be a pointer to a function either for
interpolation or restriction. An abstract description of those
functions can be found in \secref{book:S:inter_restrict} and a more
detailed one for Lagrange elements in \secref{S:Lagrange_elements}.

After these preliminaries the main part of the coarsening can be
performed. DOFs that have been created in the refinement step are now
freed again, and the children of all patch elements are freed and the
pointer to the first child is set to \nil and the pointer to the
second child is adjusted to the \code{leaf\_data} of the parent, or
also set to \nil. Thus, all fine grid information is lost at that
moment, which makes clear that a restriction of data has to be done in
advance.%
\idx{mesh coarsening!implementation|)}

%%% Local Variables: 
%%% mode: latex
%%% TeX-master: "alberta-man"
%%% End: