Some updates based on recent changes.

git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/EHFinder/trunk@156 2a26948c-0e4f-0410-aee8-f1d3e353619c

1 files changed, 67 insertions, 31 deletions

diff --git a/doc/documentation.tex b/doc/documentation.tex
index ba33704..f8e7431 100644
--- a/doc/documentation.tex
+++ b/doc/documentation.tex
@@ -169,7 +169,7 @@ More details about the actual use of the thorn in
  section~\ref{AEIDevelopment_EHFinder_UseThorn}
 
 \section{Re-initialization}
-\label{AEIDevelopment_EHFinder_reparam}
+\label{AEIDevelopment_EHFinder_re_init}
 The evolution equation for $f$, 
 equation~(\ref{AEIDevelopment_EHFinder_evolve}), causes steepening of
 the gradient of $f$, which is undesireble numerically. For that reason, $f$
@@ -179,7 +179,7 @@ $f$ are changed so that the the value of $f$ in a grid point is equal to the
 evolving $f$ according to the following evolution equation (in the parameter
 $\lambda$)
 \begin{equation}
-\frac{df}{d\lambda} = -\frac{f}{\sqrt{f^{2}+1}}\left (|\nabla f|-1\right ).
+\frac{df}{d\lambda} = -\frac{f}{\sqrt{f^{2}+1}}\left (|\nabla f|-1\right )
 \label{AEIDevelopment_EHFinder_reinit}
 \end{equation}
 until a steady state is achieved. This method is called the {\tt pde}-method
@@ -258,9 +258,11 @@ Changes to the excision regions are only done after re-initialization,
 since it is only at this time that $f$ is a distance function. The excision
 regions can move across the grid, tracking the surface $f=0$.
 
-Also if the numerical was run done with excision using {\tt SpaceMask},
+Also if the numerical run was done with excision using {\tt SpaceMask},
 the {\tt EHFinder} excision region is guaranteed to completely cover the 
-numerical excision region.
+numerical excision region. The {\tt EHFinder} currently only supports the
+old style excision mask but the support of the new style excision mask
+should be trivial and fast to implement.
 \section{Upwinding}
 \label{AEIDevelopment_EHFinder_upwind}
 All finite differences used in the evolution of the null surface are second
@@ -268,7 +270,7 @@ order one sided differences. For that reason a {\tt ghost\_size} larger or
 equal to 2 should always be used. It is possible to choose between three
 different upwinding schemes. This is done by setting the parameter
 {\tt ehfinder::upwind\_type} to either {\tt intrinsic}, {\tt shift} or
-{\tt characteristic} {\em [Note: get rid of the shift upwinding]}.
+{\tt characteristic}.
 
 The {\tt intrinsic} scheme, looks at the values of $f$ itself, to determine
 the direction of the stencil. This is basically in order to be able to handle
@@ -309,6 +311,11 @@ Figure~\ref{AEIDevelopment_EHFinder_upwindfig}. However, if the
 re-initialization is done often enough, this turns out not to cause
 any problems.
 
+The {\tt shift} scheme was implemented as an improvement to the {\tt intrinsic}
+scheme but it turned out that the {\tt characteristic} scheme was superior.
+Therefore {\tt shift} upwinding should not be used. It might be completely
+removed in future revisions.
+
 For the re-initialization the default is to use the {\tt intrinsic} 
 second order scheme (the re-initialization doesn't depend on where the
 surface is moving), so {\tt characteristic} upwinding is not applicable.
@@ -322,15 +329,14 @@ Here the most important parameters are described.
 \begin{itemize}
 \item {\tt ehfinder::mode} \\  
   The mode can either be set to {\tt normal} (normal event horizon finder
-  mode), {\tt test\_reparam} (mode to test the re-initialization routine),
-  {\tt analysis} (compare a previously calculated level set function to
+  mode), {\tt analysis} (compare a previously calculated level set function to
   a small number of analytic spacetimes) and {\tt generator} (only evolve
   the generators while keeping the level set function fixed). The
   default is {\tt normal} and should normally not be changed. The other modes
   are only useful for debugging and testing purposes.
 \item{\tt ehfinder::eh\_number\_level\_sets} \\
   An integer parameter specifying how many individual level set functions to
-  evolve at a time.
+  evolve at a time. Currently it has to be between 1 and 10.
 \item {\tt ehfinder::eh\_metric\_type} \\
   The metric type can either be set to {\tt numerical} or {\tt analytic}. If
   it is set to {\tt numerical} the {\tt EHFinder} will attempt to read in the
@@ -365,7 +371,7 @@ Here the most important parameters are described.
   A vector parameter specifying how much to translate the initial surface in
   the $y$-direction ($y_{0}$ in equation~\ref{AEIDevelopment_EHFinder_sphere}).
   Also used for the initial ellipsoid. The default is 0.
-\item {\tt translate\_z} \\
+\item {\tt translate\_z[i]} \\
   A vector parameter specifying how much to translate the initial surface in
   the $z$-direction ($z_{0}$ in equation~\ref{AEIDevelopment_EHFinder_sphere}).
   Also used for the initial ellipsoid. The default is 0.
@@ -393,28 +399,33 @@ Here the most important parameters are described.
   as inactive and are not evolved as described in
   section~\ref{AEIDevelopment_EHFinder_excise}. The default is $7$
   gridspacings.
+\item {\tt use\_inner\_excision} \\
+  A boolean parameter specifying whether the interior excision should be
+  used or not.
+\item {\tt use\_outer\_excision} \\
+  A boolean parameter specifying whether the exterior excision should be
+  used or not.
 \item {\tt upwind\_type} \\
-  The type of upwinding to be used (either {\tt intrinsic} or 
+  The type of upwinding to be used (either {\tt intrinsic}, {\tt shift} or 
   {\tt characteristic}). See the detailed description of the upwinding types
   in section~\ref{AEIDevelopment_EHFinder_upwind}. The default is 
   {\tt characteristic}.
-\item {\tt reparam\_undo} \\
+\item {\tt surface\_direction} \\
+  The code can track both outgoing and ingoing null surfaces. Choose the 
+  direction by using {\tt outward} or {\tt inward}. The default is
+  {\tt outward}. Note that the code only works as an event horizon finder
+  when evolving outward going null surfaces backwards in time.
+\item {\tt re\_init\_undo} \\
   Should the re-initialization be undone just before pinch-off or not as
-  described in section~\ref{AEIDevelopment_EHFinder_reparam}. The
+  described in section~\ref{AEIDevelopment_EHFinder_re_init}. The
   default is {\tt "no"}.
-\item {\tt re\_param\_method} \\
-  Choose the re-initialization method. The choices are {\tt approx},
-  {\tt pde} or {\tt mixed}. As described in
-  section~\ref{AEIDevelopment_EHFinder_reparam} the recommended method is
-  {\tt pde}. Only use {\tt approx} if you really know what you are doing. The
-  default is {\tt pde}.
-\item {\tt re\_param\_int\_method} \\
+\item {\tt re\_init\_int\_method} \\
   Choose the integration method in the {\tt pde}-re-initialization method.
   Choose either a simple Euler ({\tt euler}) integration scheme or a second
   order Runge-Kutta ({\tt rk2}) scheme. Since a pde is evolved to steady state,
   it seems that the Euler scheme works just fine and is faster than the
   Runge-Kutta scheme. The default is {\tt euler}.
-\item {\tt re\_param\_max\_iter} \\
+\item {\tt re\_init\_max\_iter} \\
   The maximum number of iterations in the {\tt pde}-re-initialization scheme,
   before giving up. Unless you are working at high resolution the default
   should be enough. The default is 800.
@@ -423,21 +434,46 @@ Here the most important parameters are described.
   re-initialization. Don't ever use anything else than second order
   upwinding ({\tt upwind2}). The other choices ({\tt centered} and {\tt upwind})
   are there only for testing purposes and might be removed.
-\item {\tt reparametrize\_every\_pde} \\
+\item {\tt re\_initialize\_every} \\
   How often to re-initialize using the {\tt pde} re-initialization. This
   depends on the problem. For some problems it is necessary to do it more
   often than for other problems. You'll have to experiment to figure out
   what works best. The default is 10.
-\item {\tt reparametrize\_every\_approx} \\
-  How often to re-initialize using the {\tt approx} re-initialization.
-  Again it depends on the problem. Can be used together with
-  {\tt reparametrize\_every\_pde} if {\tt re\_param\_method} is set to
-  {\tt mixed}. The default is 10.
-\item {\tt last\_iteration\_number} The last iteration number of the numerical
+\item {\tt last\_iteration\_number} \\
+  The last iteration number of the numerical
   simulation that produced the metric data. Active when {\tt eh\_metric\_type}
-  is equal to {\tt numerical}. This is used in the code to figure out which
-  data set iteration number to read in from file
-  ({\tt last\_iteration\_number-cctk\_iteration}).
+  is equal to {\tt numerical}.
+\item {\tt saved\_iteration\_every} \\
+  How often was the numerical data saved?
+  This and the above parameter is used in the code to figure out which
+  data set iteration number to read in from file \\
+  ({\tt last\_iteration\_number-saved\_iteration\_every*cctk\_iteration}).
+\item {\tt ntheta} \\
+  How many points in the $\theta$-direction should be used for integrations
+  over the surface.
+\item {\tt nphi} \\
+  How many points in the $\phi$-direction should be used for integrations
+  over the surface.
+\item {\tt maximum\_surface\_number} \\
+  The maximal number of surfaces expected at any given time during the whole
+  evolution.
+\item {\tt surface\_interpolator} \\
+  Which interpolator to use for the location of points on the surface. It
+  should be the name of a valid interpolator provided by one of the available
+  interpolators. The interpolator should be able to return both the
+  interpolated value and derivatives. At present this is provided by
+  {\tt AEILocalInterp}. The default is {\tt Hermite polynomial interpolation}.
+\item {\tt surface\_interpolation\_order} \\
+  The interpolation order used for finding points on the surface. Higher
+  orders require larger number of ghost zones for parallel runs. The default
+  value is 2 (consistent with a ghostsize of 2).
+\item {\tt area\_interpolator} \\
+  Which interpolator should be used to interpolate metric information onto the
+  surface once the surface points have been found. The default is
+  {\tt Lagrange polynomial interpolation}.
+\item {\tt area\_interpolation\_order} \\
+  The interpolation order used for calculating the area of the surface.
+  The default value is 3 (consistent with a ghostsize of 2).
 \end{itemize}
 {\tt EHFinder} also extends the following parameter from {\tt admbase} in order
 to be able to read in initial data.
@@ -461,7 +497,7 @@ these grid functions, as long as storage is turned on and they are set to the
 correct value initially when {\tt EHFinder} is run. If the evolution was done
 with a conformal factor then {\tt StaticConformal::psi} has to be output as
 well, since it is necessary in order to reconstruct the 4-metric. However it
-is only necessary to output is once since it is constant during the evolution.
+is only necessary to output this once since it is constant during the evolution.
 It is not necessary to output the derivatives of the conformal factor. If
 excision was used in the numerical run then it is also necessary to output 
 {\tt SpaceMask::emask}\footnote{When more thorns have been converted to use