SpECTRE  v2024.12.16
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
domain::CoordinateMaps::TimeDependent::Shape Class Reference

Distorts a distribution of points radially according to a spherical harmonic expansion while preserving angles. More...

#include <Shape.hpp>

Public Types

using FunctionsOfTimeMap = std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > >
 

Public Member Functions

 Shape (const std::array< double, 3 > &center, size_t l_max, size_t m_max, std::unique_ptr< ShapeMapTransitionFunctions::ShapeMapTransitionFunction > transition_func, std::string shape_function_of_time_name, std::optional< std::string > size_function_of_time_name=std::nullopt)
 
 Shape (Shape &&)=default
 
Shapeoperator= (Shape &&)=default
 
 Shape (const Shape &rhs)
 
Shapeoperator= (const Shape &rhs)
 
template<typename T >
std::array< tt::remove_cvref_wrap_t< T >, 3 > operator() (const std::array< T, 3 > &source_coords, double time, const FunctionsOfTimeMap &functions_of_time) const
 
std::optional< std::array< double, 3 > > inverse (const std::array< double, 3 > &target_coords, double time, const FunctionsOfTimeMap &functions_of_time) const
 
template<typename T >
std::array< tt::remove_cvref_wrap_t< T >, 3 > frame_velocity (const std::array< T, 3 > &source_coords, double time, const FunctionsOfTimeMap &functions_of_time) const
 
template<typename T >
tnsr::Ij< tt::remove_cvref_wrap_t< T >, 3, Frame::NoFramejacobian (const std::array< T, 3 > &source_coords, double time, const FunctionsOfTimeMap &functions_of_time) const
 
template<typename T >
tnsr::Ij< tt::remove_cvref_wrap_t< T >, 3, Frame::NoFrameinv_jacobian (const std::array< T, 3 > &source_coords, double time, const FunctionsOfTimeMap &functions_of_time) const
 
void coords_frame_velocity_jacobian (gsl::not_null< std::array< DataVector, 3 > * > source_and_target_coords, gsl::not_null< std::array< DataVector, 3 > * > frame_vel, gsl::not_null< tnsr::Ij< DataVector, 3, Frame::NoFrame > * > jac, double time, const FunctionsOfTimeMap &functions_of_time) const
 An optimized call that computes the target coordinates, frame velocity and jacobian at once to avoid duplicate calculations. More...
 
void pup (PUP::er &p)
 
const std::unordered_set< std::string > & function_of_time_names () const
 

Static Public Member Functions

static bool is_identity ()
 

Static Public Attributes

static constexpr size_t dim = 3
 

Friends

bool operator== (const Shape &lhs, const Shape &rhs)
 

Detailed Description

Distorts a distribution of points radially according to a spherical harmonic expansion while preserving angles.

Details

The shape map distorts the distance r between a point and the center while leaving the angles θ, ϕ between them preserved by applying a spherical harmonic expansion with time-dependent coefficients λlm(t). There are two ways to specify the time-dependent coefficients λlm(t):

  1. A single FunctionOfTime which specifies all coefficients. This FunctionOfTime should have ylm::Spherepack::spectral_size() number of components. These are in Spherepack order and should be the Spherepack coefficients, not the spherical harmonic coefficients. See the note below. To use this, set the size_function_of_time_name argument of the constructor to std::nullopt.
  2. Two different FunctionOfTimes. The first is similar to 1.) in that it should have the same number of components, be in Spherepack order, and be the Spherepack coefficients. The only difference is that the l=0 coefficient should be identically 0. The second FunctionOfTime should have a single component which will be the l=0 coefficient. This component should be stored as the spherical harmonic coefficient and not a Spherepack coefficient. See the note below. To use this method, set the size_function_of_time_name argument of the constructor to the name of the FunctionOfTime that's in the cache. This method is useful if we have control systems because we have a separate control system controlling a separate function of time for the l=0 coefficient than we do for the other coefficients.
Note
The quantities stored in the "shape" FunctionOfTime (the shape_function_of_time_name argument in the constructor that must always be specified) are not the complex spherical-harmonic coefficients λlm(t), but instead are the real-valued SPHEREPACK coefficients alm(t) and blm(t) used by Spherepack. This is the same for both methods of specifying FunctionOfTimes above. The relationship between these two sets of coefficients is

(1)al0=2πλl0l0,(2)alm=(1)m2πRe(λlm)l1,m1,(3)blm=(1)m2πIm(λlm)l1,m1.

The "shape" FunctionOfTime stores coefficients only for non-negative m; this is because the function we are expanding is real, so the coefficients for m<0 can be obtained from m>0 coefficients by complex conjugation. If the size_function_of_time_name argument is given to the constructor, then it is asserted that the l=0 coefficient of the "shape" function of time is exactly 0. The l=0 coefficient is then controlled by the "size" FunctionOfTime. Unlike the "shape" FunctionOfTime, the quantity in the "size" FunctionOfTime is the "complex" spherical harmonic coefficient λ00(t), and not the SPHEREPACK coefficient a00(t) ("complex" is in quotes because all m=0 coefficients are always real.) Here and below we write the equations in terms of λlm(t) instead of alm(t) and blm(t), regardless of which FunctionOfTime representation we are using, because the resulting expressions are much shorter.

Note
Also note that the FunctionOfTime coefficients λlm(t) are stored as negative of the coefficients you'd retrieve from a Strahlkorper. This is because you would typically represent the expansion of a strahlkorper as S(r)=+SlmYlm. However, in equation 9 there is a minus sign on the λlmYlm, not a plus sign. Therefore, λlm(t) picks up an extra factor of 1. This is purely a choice of convention.

An additional dimensionless domain-dependent transition function f(r,θ,ϕ) ensures that the distortion falls off correctly to zero at a certain boundary (must be a block boundary). The function f(r,θ,ϕ) is restricted such that

(4)0f(r,θ,ϕ)1

Mapped coordinates

Given a point with cartesian coordinates ξi, let the polar coordinates (r,θ,ϕ) with respect to a center xci be defined in the usual way:

(5)ξ0xc0=rsin(θ)cos(ϕ)(6)ξ1xc1=rsin(θ)sin(ϕ)(7)ξ2xc2=rcos(θ)

The shape map maps the unmapped coordinates ξi to coordinates xi:

(8)xi=ξi(ξixci)f(r,θ,ϕ)rlmλlm(t)Ylm(θ,ϕ).

Or written another way

(9)xi=xci+(ξixci)(1f(r,θ,ϕ)rlmλlm(t)Ylm(θ,ϕ)).

The form in Eq. 9 makes two things clearer

  1. This shape map is just a radial distortion about xci
  2. The coefficients λlm have units of distance because λlm(t)Ylm(θ,ϕ)/r must be dimensionless (because f is dimensionless).

Inverse map

The inverse map is given by:

(10)ξi=xci+(xixci)(r/r~),

where r~ is the radius of ξ, calculated by the transition map. In order to compute r/r~, the following equation must be solved

(11)rr~=11f(r,θ,ϕ)λlm(t)Ylm(θ,ϕ)/r

For more details, see ShapeMapTransitionFunction::original_radius_over_radius .

Frame velocity

The frame velocity vi =dxi/dt is calculated trivially:

(12)vi=(ξixci)f(r,θ,ϕ)rlmλ˙lm(t)Ylm(θ,ϕ).

Jacobian

The Jacobian is given by:

xiξj=δji(1f(r,θ,ϕ)rlmλlm(t)Ylm(θ,ϕ))(13)(ξixci)[(1rξjf(r,θ,ϕ)ξjr3f(r,θ,ϕ))lmλlm(t)Ylm(θ,ϕ)+f(r,θ,ϕ)rlmλlm(t)ξjYlm(θ,ϕ)].

where ξj=ξj.

Inverse Jacobian

The inverse Jacobian is computed by numerically inverting the Jacobian.

For future optimization, the interpolation_info objects calculated in all functions of this class could be cached. Since every element should evaluate the same grid coordinates most time steps, this might greatly decrease computation. Every element has their own clone of the shape map so the caching could be done with member variables. Care must be taken that jacobian currently calculates the interpolation_info with an order higher.

Member Function Documentation

◆ coords_frame_velocity_jacobian()

void domain::CoordinateMaps::TimeDependent::Shape::coords_frame_velocity_jacobian ( gsl::not_null< std::array< DataVector, 3 > * >  source_and_target_coords,
gsl::not_null< std::array< DataVector, 3 > * >  frame_vel,
gsl::not_null< tnsr::Ij< DataVector, 3, Frame::NoFrame > * >  jac,
double  time,
const FunctionsOfTimeMap functions_of_time 
) const

An optimized call that computes the target coordinates, frame velocity and jacobian at once to avoid duplicate calculations.

Details

The first argument source_and_target_coords should contain the source coordinates and will be overwritten in place with the target coordinates.


The documentation for this class was generated from the following file: