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.