SpECTRE
v2024.09.29
|
Items related to spherical harmonics. More...
Namespaces | |
namespace | Tags |
Holds tags and ComputeItems associated with a ylm::Strahlkorper . | |
Classes | |
class | Spherepack |
Defines the C++ interface to SPHEREPACK. More... | |
class | SpherepackIterator |
Iterates over spectral coefficients stored in SPHEREPACK format. More... | |
class | Strahlkorper |
A star-shaped surface expanded in spherical harmonics. More... | |
Enumerations | |
enum class | AngularOrdering { Strahlkorper , Cce } |
Label for the ordering of spherical harmonic points on a sphere. More... | |
Functions | |
std::ostream & | operator<< (std::ostream &os, AngularOrdering ordering) |
template<typename Frame > | |
void | change_expansion_center_of_strahlkorper (gsl::not_null< Strahlkorper< Frame > * > strahlkorper, const std::array< double, 3 > &new_center) |
Changes the expansion center of a Strahlkorper, where the expansion center is defined as the point about which the spectral basis of the Strahlkorper is expanded, which is the quantity returned by Strahlkorper::expansion_center() . | |
template<typename Frame > | |
void | change_expansion_center_of_strahlkorper_to_physical (gsl::not_null< Strahlkorper< Frame > * > strahlkorper) |
Changes the expansion center of a Strahlkorper to the physical center. Because Strahlkorper::physical_center() returns only an approximate quantity, change_expansion_center_of_strahlkorper_to_physical is iterative, and does not return exactly the same result as passing Strahlkorper::physical_center() to change_expansion_center_of_strahlkorper . | |
template<typename Frame > | |
void | fill_ylm_legend_and_data (gsl::not_null< std::vector< std::string > * > legend, gsl::not_null< std::vector< double > * > data, const ylm::Strahlkorper< Frame > &strahlkorper, double time, size_t max_l) |
Fills the legend and row of spherical harmonic data to write to disk. More... | |
template<typename Frame > | |
std::vector< ylm::Strahlkorper< Frame > > | read_surface_ylm (const std::string &file_name, const std::string &surface_subfile_name, size_t requested_number_of_times_from_end) |
Returns a list of ylm::Strahlkorper s constructed from reading in spherical harmonic data for a surface at a requested list of times. More... | |
template<typename Frame > | |
ylm::Strahlkorper< Frame > | read_surface_ylm_single_time (const std::string &file_name, const std::string &surface_subfile_name, double time, double relative_epsilon, bool check_frame=true) |
Similar to ylm::read_surface_ylm , this reads in spherical harmonic data for a surface and constructs a ylm::Strahlkorper . However, this function only does it at a specific time and returns a single ylm::Strahlkorper . More... | |
bool | operator== (const Spherepack &lhs, const Spherepack &rhs) |
bool | operator!= (const Spherepack &lhs, const Spherepack &rhs) |
bool | operator== (const SpherepackIterator &lhs, const SpherepackIterator &rhs) |
bool | operator!= (const SpherepackIterator &lhs, const SpherepackIterator &rhs) |
std::ostream & | operator<< (std::ostream &os, const SpherepackIterator::CoefficientArray &coefficient_array) |
template<typename Frame > | |
bool | operator== (const Strahlkorper< Frame > &lhs, const Strahlkorper< Frame > &rhs) |
template<typename Frame > | |
bool | operator!= (const Strahlkorper< Frame > &lhs, const Strahlkorper< Frame > &rhs) |
template<typename Frame > | |
void | write_strahlkorper_coords_to_text_file (const Strahlkorper< Frame > &strahlkorper, const std::string &output_file_name, AngularOrdering ordering, bool overwrite_file=false) |
Writes the collocation points of a ylm::Strahlkorper to an output text file. More... | |
void | write_strahlkorper_coords_to_text_file (double radius, size_t l_max, const std::array< double, 3 > ¢er, const std::string &output_file_name, AngularOrdering ordering, bool overwrite_file=false) |
Writes the collocation points of a ylm::Strahlkorper to an output text file. More... | |
void | real_spherical_harmonic (gsl::not_null< DataVector * > spherical_harmonic, const DataVector &theta, const DataVector &phi, size_t l, int m) |
Evaluates a real spherical harmonic of order l at the requested angles | |
DataVector | real_spherical_harmonic (const DataVector &theta, const DataVector &phi, size_t l, int m) |
Evaluates a real spherical harmonic of order l at the requested angles | |
template<typename Fr > | |
tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > | theta_phi (const Strahlkorper< Fr > &strahlkorper) |
template<typename Fr > | |
void | theta_phi (const gsl::not_null< tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > * > theta_phi, const Strahlkorper< Fr > &strahlkorper) |
template<typename Fr > | |
tnsr::i< DataVector, 3, Fr > | rhat (const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
void | rhat (const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > r_hat, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
ylm::Tags::aliases::Jacobian< Fr > | jacobian (const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
void | jacobian (const gsl::not_null< ylm::Tags::aliases::Jacobian< Fr > * > jac, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
ylm::Tags::aliases::InvJacobian< Fr > | inv_jacobian (const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
void | inv_jacobian (const gsl::not_null< ylm::Tags::aliases::InvJacobian< Fr > * > inv_jac, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
ylm::Tags::aliases::InvHessian< Fr > | inv_hessian (const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
void | inv_hessian (const gsl::not_null< ylm::Tags::aliases::InvHessian< Fr > * > inv_hess, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
Scalar< DataVector > | radius (const Strahlkorper< Fr > &strahlkorper) |
template<typename Fr > | |
void | radius (const gsl::not_null< Scalar< DataVector > * > result, const Strahlkorper< Fr > &strahlkorper) |
template<typename Fr > | |
tnsr::I< DataVector, 3, Fr > | cartesian_coords (const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius, const tnsr::i< DataVector, 3, Fr > &r_hat) |
template<typename Fr > | |
void | cartesian_coords (const gsl::not_null< tnsr::I< DataVector, 3, Fr > * > coords, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius, const tnsr::i< DataVector, 3, Fr > &r_hat) |
template<typename Fr > | |
tnsr::I< DataVector, 3, Fr > | cartesian_coords (const Strahlkorper< Fr > &strahlkorper) |
template<typename Fr > | |
tnsr::i< DataVector, 3, Fr > | cartesian_derivs_of_scalar (const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius_of_strahlkorper, const ylm::Tags::aliases::InvJacobian< Fr > &inv_jac) |
template<typename Fr > | |
void | cartesian_derivs_of_scalar (const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > dx_scalar, const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius_of_strahlkorper, const ylm::Tags::aliases::InvJacobian< Fr > &inv_jac) |
template<typename Fr > | |
tnsr::ii< DataVector, 3, Fr > | cartesian_second_derivs_of_scalar (const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius_of_strahlkorper, const ylm::Tags::aliases::InvJacobian< Fr > &inv_jac, const ylm::Tags::aliases::InvHessian< Fr > &inv_hess) |
template<typename Fr > | |
void | cartesian_second_derivs_of_scalar (const gsl::not_null< tnsr::ii< DataVector, 3, Fr > * > d2x_scalar, const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius_of_strahlkorper, const ylm::Tags::aliases::InvJacobian< Fr > &inv_jac, const ylm::Tags::aliases::InvHessian< Fr > &inv_hess) |
template<typename Fr > | |
Scalar< DataVector > | laplacian_of_scalar (const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
void | laplacian_of_scalar (const gsl::not_null< Scalar< DataVector > * > laplacian, const Scalar< DataVector > &scalar, const Strahlkorper< Fr > &strahlkorper, const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > &theta_phi) |
template<typename Fr > | |
ylm::Tags::aliases::Jacobian< Fr > | tangents (const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius, const tnsr::i< DataVector, 3, Fr > &r_hat, const ylm::Tags::aliases::Jacobian< Fr > &jac) |
template<typename Fr > | |
void | tangents (const gsl::not_null< ylm::Tags::aliases::Jacobian< Fr > * > result, const Strahlkorper< Fr > &strahlkorper, const Scalar< DataVector > &radius, const tnsr::i< DataVector, 3, Fr > &r_hat, const ylm::Tags::aliases::Jacobian< Fr > &jac) |
template<typename Fr > | |
tnsr::i< DataVector, 3, Fr > | normal_one_form (const tnsr::i< DataVector, 3, Fr > &dx_radius, const tnsr::i< DataVector, 3, Fr > &r_hat) |
template<typename Fr > | |
void | normal_one_form (const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > one_form, const tnsr::i< DataVector, 3, Fr > &dx_radius, const tnsr::i< DataVector, 3, Fr > &r_hat) |
template<typename Fr > | |
std::vector< std::array< double, 4 > > | fit_ylm_coeffs (const DataVector ×, const std::vector< Strahlkorper< Fr > > &strahlkorpers) |
template<typename Frame > | |
void | time_deriv_of_strahlkorper (gsl::not_null< Strahlkorper< Frame > * > time_deriv, const std::deque< std::pair< double, Strahlkorper< Frame > > > &previous_strahlkorpers) |
Compute the time derivative of a Strahlkorper from a number of previous Strahlkorpers. More... | |
Scalar< double > | ylm_to_stf_0 (const ModalVector &l0_coefs) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l. More... | |
template<typename Frame > | |
tnsr::i< double, 3, Frame > | ylm_to_stf_1 (const ModalVector &l1_coefs) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l. More... | |
template<typename Frame > | |
tnsr::ii< double, 3, Frame > | ylm_to_stf_2 (const ModalVector &l2_coefs) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l. More... | |
Items related to spherical harmonics.
Contains functions that depend on a Strahlkorper but not on a metric.
|
strong |
Label for the ordering of spherical harmonic points on a sphere.
Strahlkorper
refers to points on a sphere ordered by SPHEREPACK because Strahlkorper
s hold YlmSpherePacks internally. Cce
refers to points on a sphere ordered by Libsharp because Cce uses Libsharp internally.
void ylm::cartesian_coords | ( | const gsl::not_null< tnsr::I< DataVector, 3, Fr > * > | coords, |
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius, | ||
const tnsr::i< DataVector, 3, Fr > & | r_hat | ||
) |
coords | The returned Cartesian coordinates. |
strahlkorper | The Strahlkorper surface. |
radius | The radius as a function of angle, as returned by ylm::radius . |
r_hat | The Euclidean radial unit vector as returned by ylm::rhat . |
tnsr::I< DataVector, 3, Fr > ylm::cartesian_coords | ( | const Strahlkorper< Fr > & | strahlkorper | ) |
This overload computes radius
, theta_phi
, and r_hat
internally. Use the other overloads if you already have these quantities.
strahlkorper | The Strahlkorper surface. |
tnsr::I< DataVector, 3, Fr > ylm::cartesian_coords | ( | const Strahlkorper< Fr > & | strahlkorper, |
const Scalar< DataVector > & | radius, | ||
const tnsr::i< DataVector, 3, Fr > & | r_hat | ||
) |
cartesian_coords(i)
is
strahlkorper | The Strahlkorper surface. |
radius | The radius as a function of angle, as returned by ylm::radius . |
r_hat | The Euclidean radial unit vector as returned by ylm::rhat . |
void ylm::cartesian_derivs_of_scalar | ( | const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > | dx_scalar, |
const Scalar< DataVector > & | scalar, | ||
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius_of_strahlkorper, | ||
const ylm::Tags::aliases::InvJacobian< Fr > & | inv_jac | ||
) |
dx_scalar | The returned derivatives of the scalar. |
scalar | The scalar to be differentiated. |
strahlkorper | The Strahlkorper surface. |
radius_of_strahlkorper | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
inv_jac | The inverse Jacobian as returned by ylm::inv_jacobian |
tnsr::i< DataVector, 3, Fr > ylm::cartesian_derivs_of_scalar | ( | const Scalar< DataVector > & | scalar, |
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius_of_strahlkorper, | ||
const ylm::Tags::aliases::InvJacobian< Fr > & | inv_jac | ||
) |
dx_scalar(i)
is
scalar | The scalar to be differentiated. |
strahlkorper | The Strahlkorper surface. |
radius_of_strahlkorper | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
inv_jac | The inverse Jacobian as returned by ylm::inv_jacobian |
void ylm::cartesian_second_derivs_of_scalar | ( | const gsl::not_null< tnsr::ii< DataVector, 3, Fr > * > | d2x_scalar, |
const Scalar< DataVector > & | scalar, | ||
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius_of_strahlkorper, | ||
const ylm::Tags::aliases::InvJacobian< Fr > & | inv_jac, | ||
const ylm::Tags::aliases::InvHessian< Fr > & | inv_hess | ||
) |
d2x_scalar | The returned 2nd derivatives of the scalar. |
scalar | The scalar to be differentiated. |
strahlkorper | The Strahlkorper surface. |
radius_of_strahlkorper | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
inv_jac | The inverse Jacobian as returned by ylm::inv_jacobian |
inv_hess | The inverse Hessian as returned by `ylminv_hessian. |
tnsr::ii< DataVector, 3, Fr > ylm::cartesian_second_derivs_of_scalar | ( | const Scalar< DataVector > & | scalar, |
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius_of_strahlkorper, | ||
const ylm::Tags::aliases::InvJacobian< Fr > & | inv_jac, | ||
const ylm::Tags::aliases::InvHessian< Fr > & | inv_hess | ||
) |
d2x_scalar(i,j)
is
scalar | The scalar to be differentiated. |
strahlkorper | The Strahlkorper surface. |
radius_of_strahlkorper | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
inv_jac | The inverse Jacobian as returned by ylm::inv_jacobian |
inv_hess | The inverse Hessian as returned by `ylminv_hessian. |
void ylm::fill_ylm_legend_and_data | ( | gsl::not_null< std::vector< std::string > * > | legend, |
gsl::not_null< std::vector< double > * > | data, | ||
const ylm::Strahlkorper< Frame > & | strahlkorper, | ||
double | time, | ||
size_t | max_l | ||
) |
Fills the legend and row of spherical harmonic data to write to disk.
The number of coefficients to write is based on max_l
, the maximum value that the input strahlkorper
could possibly have. When strahlkorper.l_max() < max_l
, coefficients with strahlkorper.l_max()
will simply be zero. Assuming the same max_l
is always used for a given surface, we will always write the same number of columns for each row, as max_l
sets the number of columns to write
std::vector< std::array< double, 4 > > ylm::fit_ylm_coeffs | ( | const DataVector & | times, |
const std::vector< Strahlkorper< Fr > > & | strahlkorpers | ||
) |
The linear least squares fit of the polynomial of order 3 given a std::vector
of Strahlkorper
s to their Strahlkorper
are the same, and the returned vector consists of std::array<double, 4>
s, each of which consists of the four coefficients that define the best fit cubic to each Strahlkorper
as a function of time.
times | The time corresponding to each Strahlkorper to be fit to. |
strahlkorpers | The Strahlkorper surfaces which consists of a set of Strahlkorper at a particular time. |
void ylm::inv_hessian | ( | const gsl::not_null< ylm::Tags::aliases::InvHessian< Fr > * > | inv_hess, |
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
InvHessian(k,i,j)
is InvHessian
is not symmetric because the Jacobians are Pfaffian. InvHessian
doesn't depend on the shape of the surface.
ylm::Tags::aliases::InvHessian< Fr > ylm::inv_hessian | ( | const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ) |
InvHessian(k,i,j)
is InvHessian
is not symmetric because the Jacobians are Pfaffian. InvHessian
doesn't depend on the shape of the surface.
void ylm::inv_jacobian | ( | const gsl::not_null< ylm::Tags::aliases::InvJacobian< Fr > * > | inv_jac, |
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
InvJacobian(0,i)
is InvJacobian(1,i)
is InvJacobian
doesn't depend on the shape of the surface.
ylm::Tags::aliases::InvJacobian< Fr > ylm::inv_jacobian | ( | const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ) |
InvJacobian(0,i)
is InvJacobian(1,i)
is InvJacobian
doesn't depend on the shape of the surface.
void ylm::jacobian | ( | const gsl::not_null< ylm::Tags::aliases::Jacobian< Fr > * > | jac, |
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
Jacobian(i,0)
is Jacobian(i,1)
is Jacobian
doesn't depend on the shape of the surface.
ylm::Tags::aliases::Jacobian< Fr > ylm::jacobian | ( | const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ) |
Jacobian(i,0)
is Jacobian(i,1)
is Jacobian
doesn't depend on the shape of the surface.
void ylm::laplacian_of_scalar | ( | const gsl::not_null< Scalar< DataVector > * > | laplacian, |
const Scalar< DataVector > & | scalar, | ||
const Strahlkorper< Fr > & | strahlkorper, | ||
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
Scalar< DataVector > ylm::laplacian_of_scalar | ( | const Scalar< DataVector > & | scalar, |
const Strahlkorper< Fr > & | strahlkorper, | ||
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
void ylm::normal_one_form | ( | const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > | one_form, |
const tnsr::i< DataVector, 3, Fr > & | dx_radius, | ||
const tnsr::i< DataVector, 3, Fr > & | r_hat | ||
) |
one_form | The returned normal one form. |
dx_radius | The Cartesian derivatives of the radius, as returned by ylm::cartesian_derivs_of_scalar with ylm::radius passed in as the scalar. |
r_hat | The radial unit vector as returned by ylm::rhat . |
tnsr::i< DataVector, 3, Fr > ylm::normal_one_form | ( | const tnsr::i< DataVector, 3, Fr > & | dx_radius, |
const tnsr::i< DataVector, 3, Fr > & | r_hat | ||
) |
normal_one_form(i)
is Rhat
and DxRadius
. See Eq. (8) of [14]. Note on the word "normal":
dx_radius | The Cartesian derivatives of the radius, as returned by ylm::cartesian_derivs_of_scalar with ylm::radius passed in as the scalar. |
r_hat | The radial unit vector as returned by ylm::rhat . |
void ylm::radius | ( | const gsl::not_null< Scalar< DataVector > * > | result, |
const Strahlkorper< Fr > & | strahlkorper | ||
) |
(Euclidean) distance
Scalar< DataVector > ylm::radius | ( | const Strahlkorper< Fr > & | strahlkorper | ) |
(Euclidean) distance
ylm::Strahlkorper< Frame > ylm::read_surface_ylm_single_time | ( | const std::string & | file_name, |
const std::string & | surface_subfile_name, | ||
double | time, | ||
double | relative_epsilon, | ||
bool | check_frame = true |
||
) |
Similar to ylm::read_surface_ylm
, this reads in spherical harmonic data for a surface and constructs a ylm::Strahlkorper
. However, this function only does it at a specific time and returns a single ylm::Strahlkorper
.
epsilon
of the time
, then an error will occur. Similarly, if no time
is found within the epsilon
, then an error will occur as well.file_name | name of the H5 file containing the surface's spherical harmonic data |
surface_subfile_name | name of the subfile (with no leading slash nor the .dat extension) within file_name that contains the surface's spherical harmonic data to read in |
time | Time to read the coefficients at. |
relative_epsilon | How much error is allowed when looking for a specific time. This is useful so users don't have to know the specific time to machine precision. |
check_frame | Whether to check the frame in the subfile or not. |
void ylm::rhat | ( | const gsl::not_null< tnsr::i< DataVector, 3, Fr > * > | r_hat, |
const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ||
) |
r_hat(i)
is
We need to choose upper vs lower indices for rhat; it doesn't matter because rhat is a quantity defined with a Euclidean metric, so we choose the lower index arbitrarily.
tnsr::i< DataVector, 3, Fr > ylm::rhat | ( | const tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > & | theta_phi | ) |
r_hat(i)
is
We need to choose upper vs lower indices for rhat; it doesn't matter because rhat is a quantity defined with a Euclidean metric, so we choose the lower index arbitrarily.
void ylm::tangents | ( | const gsl::not_null< ylm::Tags::aliases::Jacobian< Fr > * > | result, |
const Strahlkorper< Fr > & | strahlkorper, | ||
const Scalar< DataVector > & | radius, | ||
const tnsr::i< DataVector, 3, Fr > & | r_hat, | ||
const ylm::Tags::aliases::Jacobian< Fr > & | jac | ||
) |
result | The computed tangent vectors. |
strahlkorper | The Strahlkorper surface. |
radius | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
r_hat | The radial unit vector as returned by ylm::rhat . |
jac | The jacobian as returned by ylm::jacobian . |
ylm::Tags::aliases::Jacobian< Fr > ylm::tangents | ( | const Strahlkorper< Fr > & | strahlkorper, |
const Scalar< DataVector > & | radius, | ||
const tnsr::i< DataVector, 3, Fr > & | r_hat, | ||
const ylm::Tags::aliases::Jacobian< Fr > & | jac | ||
) |
tangents(i,j)
is cartesian_coords
) and are considered functions of
tangents(i,0)
and tangents(i,1)
are orthogonal to the normal_one_form
tangents(i,0)
and tangents(i,1)
are not necessarily orthogonal to each other, since orthogonality between 2 vectors (as opposed to a vector and a one-form) is metric-dependent.
strahlkorper | The Strahlkorper surface. |
radius | The radius of the Strahlkorper at each point, as returned by ylm::radius . |
r_hat | The radial unit vector as returned by ylm::rhat . |
jac | The jacobian as returned by ylm::jacobian . |
void ylm::theta_phi | ( | const gsl::not_null< tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > * > | theta_phi, |
const Strahlkorper< Fr > & | strahlkorper | ||
) |
We need to choose upper vs lower indices for theta_phi; it doesn't matter because these are coordinates and not geometric objects, so we choose lower indices arbitrarily.
tnsr::i< DataVector, 2, ::Frame::Spherical< Fr > > ylm::theta_phi | ( | const Strahlkorper< Fr > & | strahlkorper | ) |
We need to choose upper vs lower indices for theta_phi; it doesn't matter because these are coordinates and not geometric objects, so we choose lower indices arbitrarily.
void ylm::time_deriv_of_strahlkorper | ( | gsl::not_null< Strahlkorper< Frame > * > | time_deriv, |
const std::deque< std::pair< double, Strahlkorper< Frame > > > & | previous_strahlkorpers | ||
) |
Compute the time derivative of a Strahlkorper from a number of previous Strahlkorpers.
Does simple 1D FD with non-uniform spacing using fd::non_uniform_1d_weights
.
time_deriv | Strahlkorper whose coefficients are the time derivative of previous_strahlkorpers ' coefficients. |
previous_strahlkorpers | All previous Strahlkorpers and the times they are at. They are expected to have the most recent Strahlkorper in the front and the Strahlkorper furthest in the past in the back of the deque. |
void ylm::write_strahlkorper_coords_to_text_file | ( | const Strahlkorper< Frame > & | strahlkorper, |
const std::string & | output_file_name, | ||
AngularOrdering | ordering, | ||
bool | overwrite_file = false |
||
) |
Writes the collocation points of a ylm::Strahlkorper
to an output text file.
The ordering of the points can be either the typical ylm::Spherepack
ordering or CCE ordering that works with libsharp. Also, an error will occur if the output file already exists, but the output file can be overwritten with the overwrite_file
option.
The second overload will construct a spherical ylm::Strahlkorper
with the given radius
, l_max
, and center
.
The output file format will be as follows with no comment or header lines, a space as the delimiter, and decimals written in scientific notation:
void ylm::write_strahlkorper_coords_to_text_file | ( | double | radius, |
size_t | l_max, | ||
const std::array< double, 3 > & | center, | ||
const std::string & | output_file_name, | ||
AngularOrdering | ordering, | ||
bool | overwrite_file = false |
||
) |
Writes the collocation points of a ylm::Strahlkorper
to an output text file.
The ordering of the points can be either the typical ylm::Spherepack
ordering or CCE ordering that works with libsharp. Also, an error will occur if the output file already exists, but the output file can be overwritten with the overwrite_file
option.
The second overload will construct a spherical ylm::Strahlkorper
with the given radius
, l_max
, and center
.
The output file format will be as follows with no comment or header lines, a space as the delimiter, and decimals written in scientific notation:
Scalar< double > ylm::ylm_to_stf_0 | ( | const ModalVector & | l0_coefs | ) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l.
Spherical harmonics of degree l are equivalent to symmetric trace-free tensors of rank l. This equivalence and the transformation is given e.g. in [191], Eqs. (2.10) - (2.14). The conversion coefficients are hard-coded to numerical precision and implemented up to order l=2.
The spherical harmonic coefficients are expected to be sorted with ascending m, i.e. (-m, -m+1, ... , m)
tnsr::i< double, 3, Frame > ylm::ylm_to_stf_1 | ( | const ModalVector & | l1_coefs | ) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l.
Spherical harmonics of degree l are equivalent to symmetric trace-free tensors of rank l. This equivalence and the transformation is given e.g. in [191], Eqs. (2.10) - (2.14). The conversion coefficients are hard-coded to numerical precision and implemented up to order l=2.
The spherical harmonic coefficients are expected to be sorted with ascending m, i.e. (-m, -m+1, ... , m)
tnsr::ii< double, 3, Frame > ylm::ylm_to_stf_2 | ( | const ModalVector & | l2_coefs | ) |
Converts real spherical harmonic coefficients of degree l into a symmetric trace-free tensor of rank l.
Spherical harmonics of degree l are equivalent to symmetric trace-free tensors of rank l. This equivalence and the transformation is given e.g. in [191], Eqs. (2.10) - (2.14). The conversion coefficients are hard-coded to numerical precision and implemented up to order l=2.
The spherical harmonic coefficients are expected to be sorted with ascending m, i.e. (-m, -m+1, ... , m)