SpECTRE
v2024.12.16
|
Classes | |
class | Spectral::Swsh::CoefficientsMetadata |
A container for libsharp metadata for the spin-weighted spherical harmonics modal representation. More... | |
class | Spectral::Swsh::SpinWeightedSphericalHarmonic |
A utility for evaluating a particular spin-weighted spherical harmonic function at arbitrary points. More... | |
class | Spectral::Swsh::SwshInterpolator |
Performs interpolation for spin-weighted spherical harmonics by taking advantage of the Clenshaw method of expanding recurrence relations. More... | |
struct | Spectral::Swsh::Tags::Eth |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::Ethbar |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::EthEth |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::EthbarEth |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::EthEthbar |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::EthbarEthbar |
Struct for labeling the | |
struct | Spectral::Swsh::Tags::InverseEth |
Struct for labeling the inverse | |
struct | Spectral::Swsh::Tags::InverseEthbar |
Struct for labeling the inverse | |
struct | Spectral::Swsh::Tags::Derivative< Tag, DerivativeKind > |
Prefix tag representing the spin-weighted derivative of a spin-weighted scalar. More... | |
struct | Spectral::Swsh::Tags::SwshTransform< Tag > |
Prefix tag representing the spin-weighted spherical harmonic transform of a spin-weighted scalar. More... | |
struct | Spectral::Swsh::Tags::LMaxBase |
Base Tag for the maximum spin-weighted spherical harmonic l; sets angular resolution. More... | |
struct | Spectral::Swsh::Tags::LMax |
Tag for the maximum spin-weighted spherical harmonic l; sets angular resolution. More... | |
struct | Spectral::Swsh::Tags::NumberOfRadialPointsBase |
Base Tag for the number of radial grid points in the three-dimensional representation of radially concentric spherical shells. More... | |
struct | Spectral::Swsh::Tags::NumberOfRadialPoints |
Tag for the number of radial grid points in the three-dimensional representation of radially concentric spherical shells. More... | |
struct | Spectral::Swsh::Tags::SwshInterpolator< Tag > |
Tag for a SwshInterpolator associated with a particular set of angular coordinates. More... | |
struct | Spectral::Swsh::SwshTransform< TagList, Representation > |
A DataBox mutate-compatible computational struct for performing several spin-weighted spherical harmonic transforms. Internally dispatches to libsharp. More... | |
struct | Spectral::Swsh::InverseSwshTransform< TagList, Representation > |
A DataBox mutate-compatible computational struct for performing several spin-weighted inverse spherical harmonic transforms. Internally dispatches to libsharp. More... | |
Typedefs | |
template<typename DerivativeTag > | |
using | Spectral::Swsh::coefficient_buffer_tags_for_derivative_tag = implementation defined |
A metafunction for determining the coefficient buffers needed by angular_derivatives() to avoid repeatedly allocating space for modal data each time a derivative is taken. More... | |
template<int Spin, typename TagList > | |
using | Spectral::Swsh::get_tags_with_spin = implementation defined |
Extract from TagList the subset of those tags that have a static int member spin equal to the template parameter Spin . More... | |
template<int Spin, typename TagList > | |
using | Spectral::Swsh::get_prefix_tags_that_wrap_tags_with_spin = implementation defined |
Extract from TagList the subset of those tags that wrap a tag that has a static int member spin equal to the template parameter Spin . More... | |
template<ComplexRepresentation Representation, typename DerivativeTagList > | |
using | Spectral::Swsh::make_transform_list_from_derivative_tags = implementation defined > >, decltype(std::make_integer_sequence< int, 5 >{})>::type |
Assemble a tmpl::list of SwshTransform s given a list of Spectral::Swsh::Tags::Derivative<Tag, Derivative> that need to be computed. The SwshTransform s constructed by this type alias correspond to the Tag s in the list. More... | |
Functions | |
constexpr size_t | Spectral::Swsh::size_of_libsharp_coefficient_vector (const size_t l_max) |
Convenience function for determining the number of spin-weighted spherical harmonics coefficients that are stored for a given l_max More... | |
constexpr double | Spectral::Swsh::sharp_swsh_sign_change (const int from_spin_weight, const int to_spin_weight, const bool real) |
Compute the relative sign change necessary to convert between the libsharp basis for spin weight from_spin_weight to the basis for spin weight to_spin_weight , for the real component coefficients if real is true, otherwise for the imaginary component coefficients. The sign change for a given coefficient is equivalent to the product of sharp_swsh_sign(from_spin, m, real) * sharp_swsh_sign(to_spin, m,
real) . Due to the form of the signs, it does not end up depending on m (the m's in the power of | |
constexpr double | Spectral::Swsh::sharp_swsh_sign (const int spin_weight, const int m, const bool real) |
Compute the sign change between the libsharp convention and the set of spin-weighted spherical harmonics given by the relation to the Wigner rotation matrices. More... | |
const CoefficientsMetadata & | Spectral::Swsh::cached_coefficients_metadata (size_t l_max) |
Generation function for obtaining a CoefficientsMetadata object which is computed by the libsharp calls only once, then lazily cached as a singleton via a static member of a function template. This is the preferred method for obtaining a CoefficientsMetadata when the l_max is not very large. More... | |
template<int Spin> | |
std::complex< double > | Spectral::Swsh::libsharp_mode_to_goldberg (size_t l, int m, size_t l_max, const SpinWeighted< ComplexModalVector, Spin > &libsharp_modes, size_t radial_offset) |
Compute the mode coefficient for the convention of [81]. See the documentation for TransformJob for complete details on the libsharp and Goldberg coefficient representations. More... | |
template<int Spin> | |
void | Spectral::Swsh::goldberg_modes_to_libsharp_modes_single_pair (const LibsharpCoefficientInfo &coefficient_info, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > libsharp_modes, size_t radial_offset, std::complex< double > goldberg_plus_m_mode_value, std::complex< double > goldberg_minus_m_mode_value) |
Set modes of a libsharp-compatible data structure by specifying the modes in the [81] representation. More... | |
template<int Spin> | |
void | Spectral::Swsh::goldberg_modes_to_libsharp_modes_single_pair (size_t l, int m, size_t l_max, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > libsharp_modes, size_t radial_offset, std::complex< double > goldberg_plus_m_mode_value, std::complex< double > goldberg_minus_m_mode_value) |
Set modes of a libsharp-compatible data structure by specifying the modes in the [81] representation. More... | |
constexpr size_t | Spectral::Swsh::goldberg_mode_index (const size_t l_max, const size_t l, const int m, const size_t radial_offset=0) |
Returns the index into a vector of modes consistent with [81]. More... | |
constexpr size_t | Spectral::Swsh::number_of_swsh_collocation_points (const size_t l_max) |
Convenience function for determining the number of spin-weighted spherical harmonic collocation values that are stored for a given l_max for a libsharp-compatible set of collocation points. | |
constexpr size_t | Spectral::Swsh::number_of_swsh_theta_collocation_points (const size_t l_max) |
Returns the number of spin-weighted spherical harmonic collocation values in | |
constexpr size_t | Spectral::Swsh::number_of_swsh_phi_collocation_points (const size_t l_max) |
Returns the number of spin-weighted spherical harmonic collocation values in | |
Mesh< 3 > | Spectral::Swsh::swsh_volume_mesh_for_radial_operations (const size_t l_max, const size_t number_of_radial_points) |
Obtain the three-dimensional mesh associated with a libsharp-compatible sequence of spherical nodal shells. More... | |
template<ComplexRepresentation Representation = ComplexRepresentation::Interleaved, int Spin, typename... ModalThenNodalTypes> | |
void | Spectral::Swsh::swsh_transform (const size_t l_max, const size_t number_of_radial_points, const gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > first_coefficient, const ModalThenNodalTypes &... coefficients_then_collocations) |
Perform a forward libsharp spin-weighted spherical harmonic transform on any number of supplied SpinWeighted<ComplexDataVector, Spin> . More... | |
template<ComplexRepresentation Representation = ComplexRepresentation::Interleaved, int Spin> | |
SpinWeighted< ComplexModalVector, Spin > | Spectral::Swsh::swsh_transform (size_t l_max, size_t number_of_radial_points, const SpinWeighted< ComplexDataVector, Spin > &collocation) |
Perform a forward libsharp spin-weighted spherical harmonic transform on a single supplied SpinWeighted<ComplexDataVector, Spin> . More... | |
template<ComplexRepresentation Representation = ComplexRepresentation::Interleaved, int Spin, typename... NodalThenModalTypes> | |
void | Spectral::Swsh::inverse_swsh_transform (const size_t l_max, const size_t number_of_radial_points, const gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > first_collocation, const NodalThenModalTypes &... collocations_then_coefficients) |
Perform an inverse libsharp spin-weighted spherical harmonic transform on any number of supplied SpinWeighted<ComplexModalVector, Spin> . More... | |
template<ComplexRepresentation Representation = ComplexRepresentation::Interleaved, int Spin> | |
SpinWeighted< ComplexDataVector, Spin > | Spectral::Swsh::inverse_swsh_transform (size_t l_max, size_t number_of_radial_points, const SpinWeighted< ComplexModalVector, Spin > &libsharp_coefficients) |
Perform an inverse libsharp spin-weighted spherical harmonic transform on a single supplied SpinWeighted<ComplexModalVector, Spin> . More... | |
template<int Spin> | |
void | Spectral::Swsh::interpolate_to_collocation (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > target, const SpinWeighted< ComplexDataVector, Spin > &source, size_t target_l_max, size_t source_l_max, size_t number_of_radial_points) |
Convert spin-weighted spherical harmonic data to a new set of collocation points (either downsampling or upsampling) | |
template<ComplexRepresentation Representation, typename TagList > | |
using | Spectral::Swsh::make_transform_list = typename detail::make_transform_list_impl< -2, Representation, TagList, decltype(std::make_integer_sequence< int, 5 >{})>::type |
Assemble a tmpl::list of SwshTransform s or InverseSwshTransform s given a list of tags TagList that need to be transformed. The Representation is the Spectral::Swsh::ComplexRepresentation to use for the transformations. More... | |
template<ComplexRepresentation Representation, typename TagList > | |
using | Spectral::Swsh::make_inverse_transform_list = typename detail::make_inverse_transform_list_impl< -2, Representation, TagList, decltype(std::make_integer_sequence< int, 5 >{})>::type |
Assemble a tmpl::list of SwshTransform s or InverseSwshTransform s given a list of tags TagList that need to be transformed. The Representation is the Spectral::Swsh::ComplexRepresentation to use for the transformations. More... | |
template<int Spin> | |
std::complex< double > | Spectral::Swsh::libsharp_mode_to_goldberg_plus_m (const LibsharpCoefficientInfo &coefficient_info, const SpinWeighted< ComplexModalVector, Spin > &libsharp_modes, size_t radial_offset) |
Compute the mode coefficient for the convention of [81]. More... | |
template<int Spin> | |
std::complex< double > | Spectral::Swsh::libsharp_mode_to_goldberg_minus_m (const LibsharpCoefficientInfo &coefficient_info, const SpinWeighted< ComplexModalVector, Spin > &libsharp_modes, size_t radial_offset) |
Compute the mode coefficient for the convention of [81]. More... | |
template<int Spin> | |
void | Spectral::Swsh::libsharp_to_goldberg_modes (gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > goldberg_modes, const SpinWeighted< ComplexModalVector, Spin > &libsharp_modes, size_t l_max) |
Compute the set of Goldberg Spin-weighted spherical harmonic modes (in the convention of [81]) from a libsharp-compatible series of modes. More... | |
template<int Spin> | |
SpinWeighted< ComplexModalVector, Spin > | Spectral::Swsh::libsharp_to_goldberg_modes (const SpinWeighted< ComplexModalVector, Spin > &libsharp_modes, size_t l_max) |
Compute the set of Goldberg Spin-weighted spherical harmonic modes (in the convention of [81]) from a libsharp-compatible series of modes. More... | |
template<int Spin> | |
void | Spectral::Swsh::goldberg_to_libsharp_modes (gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > libsharp_modes, const SpinWeighted< ComplexModalVector, Spin > &goldberg_modes, size_t l_max) |
Compute the set of libsharp-compatible spin-weighted spherical harmonic modes from a set of Goldberg modes (following the convention of [81]) More... | |
template<int Spin> | |
SpinWeighted< ComplexModalVector, Spin > | Spectral::Swsh::goldberg_to_libsharp_modes (const SpinWeighted< ComplexModalVector, Spin > &goldberg_modes, size_t l_max) |
Compute the set of libsharp-compatible spin-weighted spherical harmonic modes from a set of Goldberg modes (following the convention of [81]) More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_volume_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_max_l, double exponential_alpha, size_t exponential_half_power, gsl::not_null< ComplexDataVector * > buffer, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > transform_buffer) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_volume_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_max_l, double exponential_alpha, size_t exponential_half_power) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_volume_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_min_l, size_t filter_max_l, double exponential_alpha, size_t exponential_half_power, gsl::not_null< ComplexDataVector * > buffer, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > transform_buffer) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_volume_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_min_l, size_t filter_max_l, double exponential_alpha, size_t exponential_half_power) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_boundary_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_max_l, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > transform_buffer) |
Filter a libsharp-compatible set of collocation points on a spherical surface. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_boundary_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_max_l) |
Filter a libsharp-compatible set of collocation points on a spherical surface. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_boundary_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_min_l, size_t filter_max_l, gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > transform_buffer) |
Filter a libsharp-compatible set of collocation points on a spherical surface. More... | |
template<int Spin> | |
void | Spectral::Swsh::filter_swsh_boundary_quantity (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > to_filter, size_t l_max, size_t filter_min_l, size_t filter_max_l) |
Filter a libsharp-compatible set of collocation points on a spherical surface. More... | |
Utilities, tags, and metafunctions for using and manipulating spin-weighted spherical harmonics
using Spectral::Swsh::coefficient_buffer_tags_for_derivative_tag = typedef tmpl::list< Spectral::Swsh::Tags::SwshTransform<typename DerivativeTag::derivative_of>, Spectral::Swsh::Tags::SwshTransform<DerivativeTag> > |
A metafunction for determining the coefficient buffers needed by angular_derivatives()
to avoid repeatedly allocating space for modal data each time a derivative is taken.
angular_derivatives()
rather than the individual angular_derivative()
. using Spectral::Swsh::get_prefix_tags_that_wrap_tags_with_spin = typedef tmpl::filter<TagList, tmpl::bind<detail::wrapped_has_spin, tmpl::_1, std::integral_constant<int, Spin> >> |
Extract from TagList
the subset of those tags that wrap a tag that has a static int member spin
equal to the template parameter Spin
.
using Spectral::Swsh::get_tags_with_spin = typedef tmpl::remove_duplicates<tmpl::filter< TagList, detail::has_spin<tmpl::_1, std::integral_constant<int, Spin> >> > |
Extract from TagList
the subset of those tags that have a static int member spin
equal to the template parameter Spin
.
using Spectral::Swsh::make_inverse_transform_list = typedef typename detail::make_inverse_transform_list_impl< -2, Representation, TagList, decltype(std::make_integer_sequence<int, 5>{})>::type |
Assemble a tmpl::list
of SwshTransform
s or InverseSwshTransform
s given a list of tags TagList
that need to be transformed. The Representation
is the Spectral::Swsh::ComplexRepresentation
to use for the transformations.
Up to five SwshTransform
s or InverseSwshTransform
s will be returned, corresponding to the possible spin values. Any number of transformations are aggregated into that set of SwshTransform
s (or InverseSwshTransform
s). The number of transforms is up to five because the libsharp utility only has capability to perform spin-weighted spherical harmonic transforms for integer spin-weights from -2 to 2.
using Spectral::Swsh::make_transform_list = typedef typename detail::make_transform_list_impl< -2, Representation, TagList, decltype(std::make_integer_sequence<int, 5>{})>::type |
Assemble a tmpl::list
of SwshTransform
s or InverseSwshTransform
s given a list of tags TagList
that need to be transformed. The Representation
is the Spectral::Swsh::ComplexRepresentation
to use for the transformations.
Up to five SwshTransform
s or InverseSwshTransform
s will be returned, corresponding to the possible spin values. Any number of transformations are aggregated into that set of SwshTransform
s (or InverseSwshTransform
s). The number of transforms is up to five because the libsharp utility only has capability to perform spin-weighted spherical harmonic transforms for integer spin-weights from -2 to 2.
using Spectral::Swsh::make_transform_list_from_derivative_tags = implementation defined> >, decltype(std::make_integer_sequence<int, 5>{})>::type |
Assemble a tmpl::list
of SwshTransform
s given a list of Spectral::Swsh::Tags::Derivative<Tag, Derivative>
that need to be computed. The SwshTransform
s constructed by this type alias correspond to the Tag
s in the list.
This is intended as a convenience utility for the first step of a derivative routine, where one transforms the set of fields for which derivatives are required.
const CoefficientsMetadata & Spectral::Swsh::cached_coefficients_metadata | ( | size_t | l_max | ) |
Generation function for obtaining a CoefficientsMetadata
object which is computed by the libsharp calls only once, then lazily cached as a singleton via a static member of a function template. This is the preferred method for obtaining a CoefficientsMetadata
when the l_max
is not very large.
See the comments in the similar implementation found in SwshCollocation.hpp
for more details on the lazy cache.
void Spectral::Swsh::filter_swsh_boundary_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_max_l | ||
) |
Filter a libsharp-compatible set of collocation points on a spherical surface.
A modal Heaviside angular filter is applied which simply sets to zero all l > filter_max_l
modes.
filter_max_l
of l_max - 3
should be used. void Spectral::Swsh::filter_swsh_boundary_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_max_l, | ||
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | transform_buffer | ||
) |
Filter a libsharp-compatible set of collocation points on a spherical surface.
A modal Heaviside angular filter is applied which simply sets to zero all l > filter_max_l
modes.
filter_max_l
of l_max - 3
should be used. void Spectral::Swsh::filter_swsh_boundary_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_min_l, | ||
size_t | filter_max_l | ||
) |
Filter a libsharp-compatible set of collocation points on a spherical surface.
A modal Heaviside angular filter is applied which simply sets to zero all l > max_l
and l < min_l
modes.
void Spectral::Swsh::filter_swsh_boundary_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_min_l, | ||
size_t | filter_max_l, | ||
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | transform_buffer | ||
) |
Filter a libsharp-compatible set of collocation points on a spherical surface.
A modal Heaviside angular filter is applied which simply sets to zero all l > max_l
and l < min_l
modes.
void Spectral::Swsh::filter_swsh_volume_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_max_l, | ||
double | exponential_alpha, | ||
size_t | exponential_half_power | ||
) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells.
Two separate filters are applied. First, an exponential radial filter is applied to each radial ray, with parameters exponential_alpha
and exponential_half_power
(see Spectral::filtering::exponential_filter
for details on these parameters). Next, a modal Heaviside angular filter is applied which simply sets to zero all l > filter_max_l
modes.
dg::Actions::ExponentialFilter
may be necessary. exponential_half_power
of 8, exponential_alpha
of 108, and filter_max_l
of l_max - 3
should be used. This gives a highly aggressive radial filter, though, and for runs not attempting to compare with SpEC it is recommended to use smaller parameters to preserve more of the radial modes. void Spectral::Swsh::filter_swsh_volume_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_max_l, | ||
double | exponential_alpha, | ||
size_t | exponential_half_power, | ||
gsl::not_null< ComplexDataVector * > | buffer, | ||
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | transform_buffer | ||
) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells.
Two separate filters are applied. First, an exponential radial filter is applied to each radial ray, with parameters exponential_alpha
and exponential_half_power
(see Spectral::filtering::exponential_filter
for details on these parameters). Next, a modal Heaviside angular filter is applied which simply sets to zero all l > filter_max_l
modes.
dg::Actions::ExponentialFilter
may be necessary. exponential_half_power
of 8, exponential_alpha
of 108, and filter_max_l
of l_max - 3
should be used. This gives a highly aggressive radial filter, though, and for runs not attempting to compare with SpEC it is recommended to use smaller parameters to preserve more of the radial modes. void Spectral::Swsh::filter_swsh_volume_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_min_l, | ||
size_t | filter_max_l, | ||
double | exponential_alpha, | ||
size_t | exponential_half_power | ||
) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells.
Two separate filters are applied. First, an exponential radial filter is applied to each radial ray, with parameters exponential_alpha
and exponential_half_power
(see Spectral::filtering::exponential_filter
for details on these parameters). Next, a modal Heaviside angular filter is applied which simply sets to zero all l > max_l
or l < min_l
modes.
dg::Actions::ExponentialFilter
may be necessary. void Spectral::Swsh::filter_swsh_volume_quantity | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | to_filter, |
size_t | l_max, | ||
size_t | filter_min_l, | ||
size_t | filter_max_l, | ||
double | exponential_alpha, | ||
size_t | exponential_half_power, | ||
gsl::not_null< ComplexDataVector * > | buffer, | ||
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | transform_buffer | ||
) |
Filter a volume collocation set in the form of consecutive libsharp-compatible spherical shells.
Two separate filters are applied. First, an exponential radial filter is applied to each radial ray, with parameters exponential_alpha
and exponential_half_power
(see Spectral::filtering::exponential_filter
for details on these parameters). Next, a modal Heaviside angular filter is applied which simply sets to zero all l > max_l
or l < min_l
modes.
dg::Actions::ExponentialFilter
may be necessary.
|
constexpr |
Returns the index into a vector of modes consistent with [81].
radial_offset
is used to index into three-dimensional data in which concentric spherical shells are stored as consecutive blocks of angular modes. Goldberg modes are stored in an m varies fastest, and ascending
[(0, 0), (1, -1), (1, 0), (1, 1), (2, -2), ...]
void Spectral::Swsh::goldberg_modes_to_libsharp_modes_single_pair | ( | const LibsharpCoefficientInfo & | coefficient_info, |
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | libsharp_modes, | ||
size_t | radial_offset, | ||
std::complex< double > | goldberg_plus_m_mode_value, | ||
std::complex< double > | goldberg_minus_m_mode_value | ||
) |
Set modes of a libsharp-compatible data structure by specifying the modes in the [81] representation.
This interface takes a LibsharpCoefficientInfo
, so that iterating over the libsharp data structure and performing edits based on the corresponding Goldberg modes is convenient.
coefficient_info | An iterator which stores an TransformJob for details. |
libsharp_modes | The libsharp-compatible modal storage to be altered. |
radial_offset | The index of which angular slice is desired. Modes for concentric spherical shells are stored as consecutive blocks of angular modes. |
goldberg_plus_m_mode_value | The coefficient in the Goldberg representation ( |
goldberg_minus_m_mode_value | The coefficient in the Goldberg representation ( |
void Spectral::Swsh::goldberg_modes_to_libsharp_modes_single_pair | ( | size_t | l, |
int | m, | ||
size_t | l_max, | ||
gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | libsharp_modes, | ||
size_t | radial_offset, | ||
std::complex< double > | goldberg_plus_m_mode_value, | ||
std::complex< double > | goldberg_minus_m_mode_value | ||
) |
Set modes of a libsharp-compatible data structure by specifying the modes in the [81] representation.
This interface sets the ( TransformJob
for details.
libsharp_modes | The libsharp-compatible modal storage to be altered. |
radial_offset | The index of which angular slice is desired. Modes for concentric spherical shells are stored as consecutive blocks of angular modes. |
goldberg_plus_m_mode_value | The coefficient in the Goldberg representation ( |
goldberg_minus_m_mode_value | The coefficient in the Goldberg representation ( |
l | the spherical harmonic |
m | the spherical harmonic |
l_max | the l_max for the provided coefficient set |
SpinWeighted< ComplexModalVector, Spin > Spectral::Swsh::goldberg_to_libsharp_modes | ( | const SpinWeighted< ComplexModalVector, Spin > & | goldberg_modes, |
size_t | l_max | ||
) |
Compute the set of libsharp-compatible spin-weighted spherical harmonic modes from a set of Goldberg modes (following the convention of [81])
Internally iterates and uses the goldberg_modes_to_libsharp_modes_single_pair()
. In many applications where it is not necessary to store the full Goldberg representation, memory allocation can be saved by manually performing the iteration and calls to goldberg_modes_to_libsharp_modes_single_pair()
.
void Spectral::Swsh::goldberg_to_libsharp_modes | ( | gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | libsharp_modes, |
const SpinWeighted< ComplexModalVector, Spin > & | goldberg_modes, | ||
size_t | l_max | ||
) |
Compute the set of libsharp-compatible spin-weighted spherical harmonic modes from a set of Goldberg modes (following the convention of [81])
Internally iterates and uses the goldberg_modes_to_libsharp_modes_single_pair()
. In many applications where it is not necessary to store the full Goldberg representation, memory allocation can be saved by manually performing the iteration and calls to goldberg_modes_to_libsharp_modes_single_pair()
.
void Spectral::Swsh::inverse_swsh_transform | ( | const size_t | l_max, |
const size_t | number_of_radial_points, | ||
const gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | first_collocation, | ||
const NodalThenModalTypes &... | collocations_then_coefficients | ||
) |
Perform an inverse libsharp spin-weighted spherical harmonic transform on any number of supplied SpinWeighted<ComplexModalVector, Spin>
.
This function places the result in one or more SpinWeighted<ComplexDataVector, Spin>
returned via provided pointer. This is a simpler interface to the same functionality as the DataBox mutation compatible InverseSwshTransform
.
The following parameter ordering for the multiple input interface is enforced by the interior function calls, but is not obvious from the explicit parameter packs in this function signature:
size_t l_max
: angular resolution for the transformsize_t number_of_radial_points
: radial resolution (number of consecutive transforms in each of the vectors)gsl::not_null<SpinWeighted<ComplexDataVector, Spin>*>
, all sharing the same Spin
: The return-by-pointer containers for the transformed nodal quantitiesconst SpinWeighted<ComplexModalVector, Spin>&
, with the same Spin
as the previous function argument set : The input containers of modal spin-weighted spherical harmonic data.template parameters:
Representation
: Either ComplexRepresentation::Interleaved
or ComplexRepresentation::RealsThenImags
, indicating the representation for intermediate steps of the transformation. The two representations will give identical results but may help or hurt performance depending on the task. If unspecified, defaults to ComplexRepresentation::Interleaved
.The input is a set of libsharp-compatible coefficients.
SpinWeighted< ComplexDataVector, Spin > Spectral::Swsh::inverse_swsh_transform | ( | size_t | l_max, |
size_t | number_of_radial_points, | ||
const SpinWeighted< ComplexModalVector, Spin > & | libsharp_coefficients | ||
) |
Perform an inverse libsharp spin-weighted spherical harmonic transform on a single supplied SpinWeighted<ComplexModalVector, Spin>
.
This function returns a SpinWeighted<ComplexDataVector, Spin>
by value (causing an allocation). This is a simpler interface to the same functionality as the DataBox mutation compatible InverseSwshTransform
. If you have two or more vectors to transform at once, consider the pass-by-pointer version of Spectral::Swsh::inverse_swsh_transform
or the DataBox interface SwshTransform
, as they are more efficient for performing several transforms at once.
The input is a set of libsharp-compatible coefficients.
std::complex< double > Spectral::Swsh::libsharp_mode_to_goldberg | ( | size_t | l, |
int | m, | ||
size_t | l_max, | ||
const SpinWeighted< ComplexModalVector, Spin > & | libsharp_modes, | ||
size_t | radial_offset | ||
) |
Compute the mode coefficient for the convention of [81]. See the documentation for TransformJob
for complete details on the libsharp and Goldberg coefficient representations.
This interface extracts the (l
, m
) Goldberg mode from the provided libsharp-compatible libsharp_modes
, at angular slice radial_offset
(Modes for concentric spherical shells are stored as consecutive blocks of angular modes). l_max
must also be specified to determine the representation details.
std::complex< double > Spectral::Swsh::libsharp_mode_to_goldberg_minus_m | ( | const LibsharpCoefficientInfo & | coefficient_info, |
const SpinWeighted< ComplexModalVector, Spin > & | libsharp_modes, | ||
size_t | radial_offset | ||
) |
Compute the mode coefficient for the convention of [81].
See the documentation for TransformJob
for complete details on the libsharp and Goldberg coefficient representations. This interface takes a LibsharpCoefficientInfo
, so that iterating over the libsharp data structure and performing computations based on the corresponding Goldberg modes is convenient. Two functions are provided, one for the positive m modes in the Goldberg representation, and one for the negative m modes, as the distinction is not made by the coefficient iterator.
coefficient_info | An iterator which stores an |
libsharp_modes | The libsharp-compatible modal storage to be accessed |
radial_offset | The index of which angular slice is desired. Modes for concentric spherical shells are stored as consecutive blocks of angular modes. |
std::complex< double > Spectral::Swsh::libsharp_mode_to_goldberg_plus_m | ( | const LibsharpCoefficientInfo & | coefficient_info, |
const SpinWeighted< ComplexModalVector, Spin > & | libsharp_modes, | ||
size_t | radial_offset | ||
) |
Compute the mode coefficient for the convention of [81].
See the documentation for TransformJob
for complete details on the libsharp and Goldberg coefficient representations. This interface takes a LibsharpCoefficientInfo
, so that iterating over the libsharp data structure and performing computations based on the corresponding Goldberg modes is convenient. Two functions are provided, one for the positive m modes in the Goldberg representation, and one for the negative m modes, as the distinction is not made by the coefficient iterator.
coefficient_info | An iterator which stores an |
libsharp_modes | The libsharp-compatible modal storage to be accessed |
radial_offset | The index of which angular slice is desired. Modes for concentric spherical shells are stored as consecutive blocks of angular modes. |
SpinWeighted< ComplexModalVector, Spin > Spectral::Swsh::libsharp_to_goldberg_modes | ( | const SpinWeighted< ComplexModalVector, Spin > & | libsharp_modes, |
size_t | l_max | ||
) |
Compute the set of Goldberg Spin-weighted spherical harmonic modes (in the convention of [81]) from a libsharp-compatible series of modes.
Modes for concentric spherical shells are stored as consecutive blocks of angular modes in both representations.
void Spectral::Swsh::libsharp_to_goldberg_modes | ( | gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | goldberg_modes, |
const SpinWeighted< ComplexModalVector, Spin > & | libsharp_modes, | ||
size_t | l_max | ||
) |
Compute the set of Goldberg Spin-weighted spherical harmonic modes (in the convention of [81]) from a libsharp-compatible series of modes.
Modes for concentric spherical shells are stored as consecutive blocks of angular modes in both representations.
|
constexpr |
Returns the number of spin-weighted spherical harmonic collocation values in
The full number of collocation points is the product of the number of
|
constexpr |
Returns the number of spin-weighted spherical harmonic collocation values in
The full number of collocation points is the product of the number of
|
constexpr |
Compute the sign change between the libsharp convention and the set of spin-weighted spherical harmonics given by the relation to the Wigner rotation matrices.
The sign change is obtained via the difference between the libsharp convention and the convention which uses:
See [81]. The sign change is computed for the real component coefficients if real
is true, otherwise for the imaginary component coefficients. For full details on the sign convention used in libsharp, see the documentation for TransformJob. This function outputs the
Note that in this equation, the "real" and "imag" superscripts refer to the set of basis functions used for the decomposition of the real and imaginary part of the spin-weighted collocation points, not real or imaginary parts of the basis functions themselves.
|
constexpr |
Compute the relative sign change necessary to convert between the libsharp basis for spin weight from_spin_weight
to the basis for spin weight to_spin_weight
, for the real component coefficients if real
is true, otherwise for the imaginary component coefficients. The sign change for a given coefficient is equivalent to the product of sharp_swsh_sign(from_spin, m, real) * sharp_swsh_sign(to_spin, m,
real)
. Due to the form of the signs, it does not end up depending on m (the m's in the power of
The sign change is obtained by the difference between the libsharp convention and the convention which uses:
See [81]
|
constexpr |
Convenience function for determining the number of spin-weighted spherical harmonics coefficients that are stored for a given l_max
This includes the factor of 2 associated with needing to store both the transform of the real and imaginary parts, so is the full size of the result of a libsharp swsh transform.
void Spectral::Swsh::swsh_transform | ( | const size_t | l_max, |
const size_t | number_of_radial_points, | ||
const gsl::not_null< SpinWeighted< ComplexModalVector, Spin > * > | first_coefficient, | ||
const ModalThenNodalTypes &... | coefficients_then_collocations | ||
) |
Perform a forward libsharp spin-weighted spherical harmonic transform on any number of supplied SpinWeighted<ComplexDataVector, Spin>
.
This function places the result in one or more SpinWeighted<ComplexModalVector, Spin>
returned via provided pointer. This is a simpler interface to the same functionality as the DataBox mutation compatible SwshTransform
.
The following parameter ordering for the multiple input interface is enforced by the interior function calls, but is not obvious from the explicit parameter packs in this function signature:
size_t l_max
: angular resolution for the transformsize_t number_of_radial_points
: radial resolution (number of consecutive transforms in each of the vectors)gsl::not_null<SpinWeighted<ComplexModalVector, Spin>*>
, all sharing the same Spin
: The return-by-pointer containers for the transformed modal quantitiesconst SpinWeighted<ComplexDataVector, Spin>&
, with the same Spin
as the previous function argument set : The input containers of nodal spin-weighted spherical harmonic data.template parameters:
Representation
: Either ComplexRepresentation::Interleaved
or ComplexRepresentation::RealsThenImags
, indicating the representation for intermediate steps of the transformation. The two representations will give identical results but may help or hurt performance depending on the task. If unspecified, defaults to ComplexRepresentation::Interleaved
.The result is a set of libsharp-compatible coefficients.
SpinWeighted< ComplexModalVector, Spin > Spectral::Swsh::swsh_transform | ( | size_t | l_max, |
size_t | number_of_radial_points, | ||
const SpinWeighted< ComplexDataVector, Spin > & | collocation | ||
) |
Perform a forward libsharp spin-weighted spherical harmonic transform on a single supplied SpinWeighted<ComplexDataVector, Spin>
.
This function returns a SpinWeighted<ComplexModalVector, Spin>
by value (causing an allocation). This is a simpler interface to the same functionality as the DataBox mutation compatible SwshTransform
. If you have two or more vectors to transform at once, consider the pass-by-pointer version of Spectral::Swsh::swsh_transform
or the DataBox interface InverseSwshTransform
, as they are more efficient for performing several transforms at once.
The result is a set of libsharp-compatible coefficients.
Mesh< 3 > Spectral::Swsh::swsh_volume_mesh_for_radial_operations | ( | const size_t | l_max, |
const size_t | number_of_radial_points | ||
) |
Obtain the three-dimensional mesh associated with a libsharp-compatible sequence of spherical nodal shells.