SpECTRE Documentation Coverage Report
Current view: top level - __w/spectre/spectre/docs/DevGuide - RedistributingGridpoints.md Hit Total Coverage
Commit: b5ffa4904430ccef0b226f73dcd25c74cb5188f6 Lines: 0 1 0.0 %
Date: 2021-07-28 22:05:18
Legend: Lines: hit not hit

          Line data    Source code
       1           0 : \cond NEVER
       2             : Distributed under the MIT License.
       3             : See LICENSE.txt for details.
       4             : \endcond
       5             : # Redistributing Gridpoints {#redistributing_gridpoints}
       6             : 
       7             : ## Introduction
       8             : The simplest way to construct a volume map from two parameterized surfaces is
       9             : by linearly interpolating between them:
      10             : 
      11             : \f[\vec{x}(\xi,\eta,\zeta) =
      12             : \frac{1-\zeta}{2}\vec{\sigma}_-(\xi,\eta)+
      13             : \frac{1+\zeta}{2}\vec{\sigma}_+(\xi,\eta)\f]
      14             : 
      15             : In the above example, each surface \f$\vec{\sigma}_+\f$ and
      16             : \f$\vec{\sigma}_-\f$ is parameterized using the logical coordinates \f$\xi\f$
      17             : and \f$\eta\f$, and a third coordinate \f$\zeta\in[-1,1]\f$ is used to
      18             : interpolate between them.
      19             : 
      20             : We then distribute gridpoints on this volume by specifying values of the
      21             : coordinates \f$\xi,\eta,\f$ and \f$\zeta\f$ at which the gridpoints are located.
      22             : In SpECTRE these values are the locations of the quadrature nodes. The
      23             : distribution of the gridpoints throughout the volume depends on the
      24             : parameterization used, and the simplest choice of parameterization does not
      25             : necessarily lead to the best gridpoint distribution. In this section we discuss
      26             : situations in which there exist better parameterizations than those obtained by
      27             : linear interpolation.
      28             : 
      29             : ## Generalized Logical Coordinates
      30             : In each of the following examples, we will obtain functions \f$\Xi(\xi),
      31             : \mathrm{H}(\eta),\f$ and \f$\mathrm{Z}(\zeta)\f$ that give better gridpoint
      32             : distributions than using the logical coordinates alone. Where possible, we will
      33             : write the reparameterized map such that the functional form of the map is
      34             : unchanged when replacing \f$\Xi\f$ with \f$\xi\f$, etc. We therefore refer to
      35             : \f$\Xi, \mathrm{H},\f$ and \f$\mathrm{Z}\f$ as the
      36             : *generalized logical coordinates*, as they can also refer to the logical
      37             : coordinates \f$\xi, \eta,\f$ and \f$\zeta\f$ themselves, when the transformation
      38             : is the identity.
      39             : 
      40             : ## Equiangular Maps
      41             : 
      42             : The mapping for a cubed sphere surface can be easily obtained by taking points
      43             : that lie on each face of a cube and normalizing them such that they lie on the
      44             : sphere:
      45             : 
      46             : \f[\vec{\sigma}_{+z}(\xi,\eta) =
      47             : \frac{1}{\sqrt{1 + \xi^2  + \eta^2}}
      48             : \begin{bmatrix}
      49             : \xi\\
      50             : \eta\\
      51             : 1\\
      52             : \end{bmatrix}\f]
      53             : 
      54             : In the above example the parameterization used for the upper \f$+z\f$ surface
      55             : of the cube is linear in \f$\xi\f$ and \f$\eta\f$. However, distances measured
      56             : on the surface of the sphere are not linear in \f$\xi\f$ and \f$\eta\f$. To see
      57             : this, one may compute \f$g_{\xi\xi} = |\frac{\partial\vec{x}}{\partial\xi}|^2\f$
      58             : to see how distances are measured in terms of \f$\xi\f$:
      59             : 
      60             : \f[g_{\xi,\xi}|_{\eta=0} = \frac{1}{(1+\xi^2)^2}\f]
      61             : 
      62             : This metric term demonstrates that a gridpoint distribution uniform in
      63             : \f$\xi\f$ will end up being compressed near \f$\xi=\pm1\f$. Suppose we
      64             : reparameterized the surface using the generalized logical coordinate
      65             : \f$\Xi\in[-1,1]\f$. We would find:
      66             : 
      67             : \f[g_{\xi,\xi}|_{\eta=0} = \frac{\Xi'^2}{(1+\Xi^2)^2}\f]
      68             : 
      69             : Ideally, we would like distances measured along a curvilinear surface to be
      70             : linear in the logical coordinates. We solve the differential equation and
      71             : obtain:
      72             : 
      73             : \f[\Xi = \tan(\xi\pi/4)\f]
      74             : 
      75             : These two parameterizations of the cubed sphere are known as the *equidistant*
      76             : and *equiangular* central projections of the cube onto the sphere. We now
      77             : summarize their usage in SpECTRE CoordinateMaps that have
      78             : `with_equiangular_map` as a specifiable parameter:
      79             : 
      80             : In the case where `with_equiangular_map` is `true`, we have the
      81             : equiangular coordinates
      82             : 
      83             : \f[\textrm{equiangular xi} : \Xi(\xi) = \textrm{tan}(\xi\pi/4)\f]
      84             : 
      85             : \f[\textrm{equiangular eta}  : \mathrm{H}(\eta) = \textrm{tan}(\eta\pi/4)\f]
      86             : 
      87             : with derivatives
      88             : 
      89             : \f[\Xi'(\xi) = \frac{\pi}{4}(1+\Xi^2)\f],
      90             : 
      91             : \f[\mathrm{H}'(\eta) = \frac{\pi}{4}(1+\mathrm{H}^2)\f]
      92             : 
      93             : In the case where `with_equiangular_map` is `false`, we have the equidistant
      94             : coordinates
      95             : 
      96             : \f[ \textrm{equidistant xi}  : \Xi = \xi\f]
      97             : 
      98             : \f[ \textrm{equidistant eta}  : \mathrm{H} = \eta\f]
      99             : 
     100             : with derivatives:
     101             : 
     102             : \f[\Xi'(\xi) = 1\f] \f[\mathrm{H}'(\eta) = 1\f]
     103             : 
     104             : ##  Projective Maps
     105             : 
     106             : The mapping for any convex quadrilateral can be obtained by bilinearly
     107             : interpolating between each vertex \f$\vec{x}_1, \vec{x}_2, \vec{x}_3\f$
     108             : and \f$\vec{x}_4\f$:
     109             : 
     110             : \f[\vec{x}(\xi,\eta) =
     111             : \frac{(1-\xi)(1-\eta)}{4}\vec{x}_1+
     112             : \frac{(1+\xi)(1-\eta)}{4}\vec{x}_2+
     113             : \frac{(1-\xi)(1+\eta)}{4}\vec{x}_3+
     114             : \frac{(1+\xi)(1+\eta)}{4}\vec{x}_4
     115             : \f]
     116             : 
     117             : In the case of a trapezoid where two of the sides are parallel, it is
     118             : appropriate to linearly interpolate along the parallel sides. However,
     119             : linearly interpolating between the two bases results in a less than
     120             : ideal gridpoint distribution. This happens in the case of SpECTRE's Frustum,
     121             : where the logical coordinate \f$\zeta\f$ interpolates between the bases.
     122             : 
     123             : \image html BilinearVProjective.png "Comparison of mappings. (Noah Veltman)"
     124             : 
     125             : As seen in Veltman's [Warp-Off]
     126             : (https://bl.ocks.org/veltman/8f5a157276b1dc18ce2fba1bc06dfb48), linear
     127             : interpolation between the two bases results in a uniformly spaced grid
     128             : between the bases of the frustum. This causes elements near the smaller base
     129             : to be longer in the direction normal to the base, and elements near the larger
     130             : base to be shorter in the direction normal to the base. We desire elements that
     131             : have roughly equal sizes along each of their dimensions.
     132             : 
     133             : We can redistribute the gridpoints in the \f$\zeta\f$ direction using a
     134             : projective map, moving more gridpoints toward the smaller base. We can also see
     135             : in the figure above that a projective map can be applied incorrectly, leaving
     136             : elements distorted at the opposite end. From this we can see that it is
     137             : important to control the degree of projection.
     138             : 
     139             : We adapt a technique from projective geometry to obtain the desired grid
     140             : spacing. The heart of the method lies in the fact that objects arranged in a
     141             : line at equal distances from one another will appear to converge as they
     142             : approach the horizon.
     143             : 
     144             : \image html ProjectionOntoPlane.png "Controlling the degree of projection."
     145             : 
     146             : The above diagram demonstrates how to obtain a nonlinearly parameterized
     147             : object (seen in red) from a linearly parameterized one (seen in purple).
     148             : This is done by lifting the linearly parameterized object into a higher
     149             : spatial dimension \f$w\f$, such that its projection onto the plane remains
     150             : unchanged. As seen above, \f$w_{\delta}\f$ controls the degree of projection
     151             : of one end of the object (purple) into a higher spatial dimension \f$w\f$.
     152             : In projective geometry, these points that exist in the higher dimension are
     153             : labeled with *homogeneous coordinates* \f$\tilde{x}, \tilde{y}, \tilde{z}, w\f$,
     154             : to distinguish them from the Cartesian coordinates that label points that exist
     155             : on the \f$w=1\f$ hyperplane, \f$x,y,z\f$. The resulting grid (seen in red) is
     156             : obtained by projecting back into the \f$w=1\f$ hyperplane. The Cartesian
     157             : coordinates are obtained by dividing each homogeneous coordinate of the
     158             : linearly parameterized object by its respective \f$w\f$ coordinate value.
     159             : 
     160             : The parametric equation for the purple object seen above in homogeneous
     161             : coordinates is:
     162             : \f[\begin{bmatrix}\tilde{x}\\\tilde{y}\\\tilde{z}\\w\\\end{bmatrix}=
     163             : \frac{1-\zeta}{2}\begin{bmatrix}x_1\\y_1\\z_1\\1\\\end{bmatrix}+
     164             : \frac{1+\zeta}{2}\begin{bmatrix}x_2w_{\delta}\\y_2w_{\delta}\\
     165             : z_2w_{\delta}\\w_{\delta}\\\end{bmatrix}\f]
     166             : 
     167             : The equation for the projected red object in Cartesian coordinates is
     168             : obtained by dividing by w:
     169             : \f[\vec{x}(\zeta) = \frac{1}{w(\zeta)}
     170             : \begin{bmatrix}
     171             : \tilde{x}(\zeta)\\
     172             : \tilde{y}(\zeta)\\
     173             : \tilde{z}(\zeta)\\
     174             : \end{bmatrix}\f]
     175             : 
     176             : We wish to cast our parametric equation for the surface into the form:
     177             : \f[\vec{x}(\zeta) =
     178             : \frac{1-\mathrm{Z}}{2}\vec{x}_1 + \frac{1+\mathrm{Z}}{2}\vec{x}_2\f]
     179             : for some appropriate choice of auxiliary variable `projective_zeta`
     180             : \f$ = \mathrm{Z}(\zeta)\f$. We would also like for \f$\mathrm{Z}\f$ to reduce to
     181             : \f$\zeta\f$ when \f$w_{\delta}\ = 1\f$.
     182             : 
     183             : Defining the auxiliary variables \f$w_{\pm} := w_{\delta}\pm 1\f$, the desired
     184             : \f$\mathrm{Z}(\zeta)\f$ is given by:
     185             : \f[\mathrm{Z} = \frac{w_- + \zeta w_+}
     186             : {w_+ + \zeta w_-}\f]
     187             : 
     188             : with derivative:
     189             : \f[\mathrm{Z}' = \frac{\partial\mathrm{Z}}{\partial \zeta} =
     190             : \frac{w_+^2 - w_-^2}{(w_+ + \zeta w_-)^2}\f]
     191             : 
     192             : For SpECTRE CoordinateMaps that have `projective_scale_factor` as a specifiable
     193             : parameter, the value \f$w_{\delta} = 1\f$ should be supplied in case the user
     194             : does not want to use projective scaling. If `auto_projective_scale_factor` is
     195             : set to `true`, the map will compute a value of \f$w_{\delta}\f$ that is
     196             : appropriate.

Generated by: LCOV version 1.14