diff options
Diffstat (limited to 'src/gr/expansion.cc')
-rw-r--r-- | src/gr/expansion.cc | 27 |
1 files changed, 25 insertions, 2 deletions
diff --git a/src/gr/expansion.cc b/src/gr/expansion.cc index e0b7ba8..4d7ae34 100644 --- a/src/gr/expansion.cc +++ b/src/gr/expansion.cc @@ -256,8 +256,17 @@ return expansion_success; // *** NORMAL RETURN *** //****************************************************************************** // -// This function sets up the global xyz positions of the grid points -// in the gridfns global_[xyz]. These will be used by interplate_geometry(). +// This function sets up the xyz-position gridfns: +// * global_[xyz] (used by interplate_geometry() ) +// * global_{xx,xy,xz,yy,yz,zz} (used by higher-level code to compute +// quadrupole moments of horizons) +// +// Bugs: +// * We initialize the global_{xx,xy,xz,yy,yz,zz} gridfns every time +// this function is called, i.e. at each horizon-finding Newton +// iteration, even though they're only needed at the end of the +// horizon-finding process. In practice the extra cost is small, +// though, so it's probably not worth fixing this... // namespace { void setup_xyz_posns(patch_system& ps, bool print_msg_flag) @@ -290,6 +299,20 @@ if (print_msg_flag) p.gridfn(gfns::gfn__global_x, irho,isigma) = global_x; p.gridfn(gfns::gfn__global_y, irho,isigma) = global_y; p.gridfn(gfns::gfn__global_z, irho,isigma) = global_z; + + const fp global_xx = global_x * global_x; + const fp global_xy = global_x * global_y; + const fp global_xz = global_x * global_z; + const fp global_yy = global_y * global_y; + const fp global_yz = global_y * global_z; + const fp global_zz = global_z * global_z; + + p.gridfn(gfns::gfn__global_xx, irho,isigma) = global_xx; + p.gridfn(gfns::gfn__global_xy, irho,isigma) = global_xy; + p.gridfn(gfns::gfn__global_xz, irho,isigma) = global_xz; + p.gridfn(gfns::gfn__global_yy, irho,isigma) = global_yy; + p.gridfn(gfns::gfn__global_yz, irho,isigma) = global_yz; + p.gridfn(gfns::gfn__global_zz, irho,isigma) = global_zz; } } } |