Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Use new common time dependent options in Sphere and BCO domains #6115

Open
wants to merge 2 commits into
base: develop
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
108 changes: 55 additions & 53 deletions src/Domain/Creators/TimeDependentOptions/BinaryCompactObject.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -42,9 +42,9 @@ TimeDependentMapOptions<IsCylindrical>::TimeDependentMapOptions(
std::optional<ShapeMapOptions<domain::ObjectLabel::B>> shape_options_B,
const Options::Context& context)
: initial_time_(initial_time),
expansion_map_options_(expansion_map_options),
rotation_map_options_(rotation_map_options),
translation_map_options_(translation_map_options),
expansion_map_options_(std::move(expansion_map_options)),
rotation_map_options_(std::move(rotation_map_options)),
translation_map_options_(std::move(translation_map_options)),
shape_options_A_(shape_options_A),
shape_options_B_(shape_options_B) {
if (not(expansion_map_options_.has_value() or
Expand Down Expand Up @@ -82,37 +82,51 @@ TimeDependentMapOptions<IsCylindrical>::create_worldtube_functions_of_time()
std::unordered_map<std::string,
std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>>
result{};

// The functions of time need to be valid only for the very first time step,
// after that they need to be updated by the worldtube singleton.
const double initial_expiration_time = initial_time_ + 1e-10;
if (not expansion_map_options_.has_value()) {
ERROR("Initial values for the expansion map need to be provided.");
}
if (not expansion_map_options_->asymptotic_velocity_outer_boundary
.has_value()) {
ERROR(
"The BinaryCompactObject domains don't support using a "
"SettleToConstant function of time for the expansion map. Set the "
"asymptotic velocity to a value and set the decay timescale to "
"'None'");
}
result[expansion_name] =
std::make_unique<FunctionsOfTime::IntegratedFunctionOfTime>(
initial_time_,
std::array<double, 2>{
{{gsl::at(expansion_map_options_.value().initial_values, 0)},
{gsl::at(expansion_map_options_.value().initial_values, 1)}}},
expansion_map_options_.value().initial_values[0][0],
expansion_map_options_.value().initial_values[1][0]},
initial_expiration_time, false);
result[expansion_outer_boundary_name] =
std::make_unique<FunctionsOfTime::FixedSpeedCubic>(
1.0, initial_time_,
expansion_map_options_.value().outer_boundary_velocity,
expansion_map_options_.value().outer_boundary_decay_time);
expansion_map_options_.value()
.asymptotic_velocity_outer_boundary.value(),
expansion_map_options_.value().decay_timescale_outer_boundary);

if (not rotation_map_options_.has_value()) {
ERROR(
"Initial values for the rotation map need to be provided when using "
"the worldtube.");
}
if (rotation_map_options_->decay_timescale.has_value()) {
ERROR(
"The BinaryCompactObject domains don't support using a "
"SettleToConstant function of time for the rotation map. Set the "
"decay timescale to 'None'");
}

result[rotation_name] =
std::make_unique<FunctionsOfTime::IntegratedFunctionOfTime>(
initial_time_,
std::array<double, 2>{
0.,
gsl::at(rotation_map_options_.value().initial_angular_velocity,
2)},
std::array<double, 2>{0., rotation_map_options_.value().angles[1][2]},
initial_expiration_time, true);

// Size and Shape FunctionOfTime for objects A and B. Only spherical excision
Expand Down Expand Up @@ -189,64 +203,52 @@ TimeDependentMapOptions<IsCylindrical>::create_functions_of_time(
// ExpansionMap FunctionOfTime for the function \f$a(t)\f$ in the
// domain::CoordinateMaps::TimeDependent::RotScaleTrans map
if (expansion_map_options_.has_value()) {
if (not expansion_map_options_->asymptotic_velocity_outer_boundary
.has_value()) {
ERROR(
"The BinaryCompactObject domains don't support using a "
"SettleToConstant function of time for the expansion map. Set the "
"asymptotic velocity to a value and set the decay timescale to "
"'None'");
}

result[expansion_name] =
std::make_unique<FunctionsOfTime::PiecewisePolynomial<2>>(
initial_time_,
std::array<DataVector, 3>{
{{gsl::at(expansion_map_options_.value().initial_values, 0)},
{gsl::at(expansion_map_options_.value().initial_values, 1)},
{0.0}}},
initial_time_, expansion_map_options_->initial_values,
expiration_times.at(expansion_name));

// ExpansionMap FunctionOfTime for the function \f$b(t)\f$ in the
// domain::CoordinateMaps::TimeDependent::RotScaleTrans map
result[expansion_outer_boundary_name] =
std::make_unique<FunctionsOfTime::FixedSpeedCubic>(
1.0, initial_time_,
expansion_map_options_.value().outer_boundary_velocity,
expansion_map_options_.value().outer_boundary_decay_time);
expansion_map_options_->initial_values_outer_boundary[0][0],
initial_time_,
expansion_map_options_->asymptotic_velocity_outer_boundary.value(),
expansion_map_options_->decay_timescale_outer_boundary);
}

// RotationMap FunctionOfTime for the rotation angles about each
// axis. The initial rotation angles don't matter as we never
// actually use the angles themselves. We only use their derivatives
// (omega) to determine map parameters. In theory we could determine
// each initial angle from the input axis-angle representation, but
// we don't need to.
// axis.
if (rotation_map_options_.has_value()) {
result[rotation_name] = std::make_unique<
FunctionsOfTime::QuaternionFunctionOfTime<3>>(
initial_time_,
std::array<DataVector, 1>{DataVector{1.0, 0.0, 0.0, 0.0}},
std::array<DataVector, 4>{
{{3, 0.0},
{gsl::at(rotation_map_options_.value().initial_angular_velocity,
0),
gsl::at(rotation_map_options_.value().initial_angular_velocity,
1),
gsl::at(rotation_map_options_.value().initial_angular_velocity,
2)},
{3, 0.0},
{3, 0.0}}},
expiration_times.at(rotation_name));
if (rotation_map_options_->decay_timescale.has_value()) {
ERROR(
"The BinaryCompactObject domains don't support using a "
"SettleToConstant function of time for the rotation map. Set the "
"decay timescale to 'None'");
}

result[rotation_name] =
std::make_unique<FunctionsOfTime::QuaternionFunctionOfTime<3>>(
initial_time_, std::array{rotation_map_options_->quaternions[0]},
rotation_map_options_->angles, expiration_times.at(rotation_name));
}

// TranslationMap FunctionOfTime
if (translation_map_options_.has_value()) {
result[translation_name] = std::make_unique<
FunctionsOfTime::PiecewisePolynomial<2>>(
initial_time_,
std::array<DataVector, 3>{
{{gsl::at(translation_map_options_.value().initial_values, 0)[0],
gsl::at(translation_map_options_.value().initial_values, 0)[1],
gsl::at(translation_map_options_.value().initial_values, 0)[2]},
{gsl::at(translation_map_options_.value().initial_values, 1)[0],
gsl::at(translation_map_options_.value().initial_values, 1)[1],
gsl::at(translation_map_options_.value().initial_values, 1)[2]},
{gsl::at(translation_map_options_.value().initial_values, 2)[0],
gsl::at(translation_map_options_.value().initial_values, 2)[1],
gsl::at(translation_map_options_.value().initial_values, 2)[2]}}},
expiration_times.at(translation_name));
result[translation_name] =
std::make_unique<FunctionsOfTime::PiecewisePolynomial<2>>(
initial_time_, translation_map_options_->initial_values,
expiration_times.at(translation_name));
}

// Size and Shape FunctionOfTime for objects A and B
Expand Down
101 changes: 11 additions & 90 deletions src/Domain/Creators/TimeDependentOptions/BinaryCompactObject.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,8 +17,10 @@
#include "Domain/CoordinateMaps/TimeDependent/Rotation.hpp"
#include "Domain/CoordinateMaps/TimeDependent/Shape.hpp"
#include "Domain/CoordinateMaps/TimeDependent/ShapeMapTransitionFunctions/ShapeMapTransitionFunction.hpp"
#include "Domain/Creators/TimeDependentOptions/ExpansionMap.hpp"
#include "Domain/Creators/TimeDependentOptions/RotationMap.hpp"
#include "Domain/Creators/TimeDependentOptions/ShapeMap.hpp"
#include "Domain/Creators/TimeDependentOptions/Sphere.hpp"
#include "Domain/Creators/TimeDependentOptions/TranslationMap.hpp"
#include "Domain/FunctionsOfTime/FunctionOfTime.hpp"
#include "Domain/Structure/ObjectLabel.hpp"
#include "Options/Auto.hpp"
Expand Down Expand Up @@ -146,101 +148,20 @@ struct TimeDependentMapOptions {
/// The outer boundary radius of the map is always set to
/// the outer boundary of the Domain, so there is no option
/// here to set the outer boundary radius.
struct ExpansionMapOptions {
using type = Options::Auto<ExpansionMapOptions, Options::AutoLabel::None>;
static std::string name() { return "ExpansionMap"; }
static constexpr Options::String help = {
"Options for the expansion map. Specify 'None' to not use this map."};
struct InitialValues {
using type = std::array<double, 2>;
static constexpr Options::String help = {
"Initial value and deriv of expansion."};
};
struct AsymptoticVelocityOuterBoundary {
using type = double;
static constexpr Options::String help = {
"The asymptotic velocity of the outer boundary."};
};
struct DecayTimescaleOuterBoundaryVelocity {
using type = double;
static constexpr Options::String help = {
"The timescale for how fast the outer boundary velocity approaches "
"its asymptotic value."};
};
using options = tmpl::list<InitialValues, AsymptoticVelocityOuterBoundary,
DecayTimescaleOuterBoundaryVelocity>;
ExpansionMapOptions() = default;
ExpansionMapOptions(std::array<double, 2> initial_values_in,
double outer_boundary_velocity_in,
double outer_boundary_decay_time_in)
: initial_values(initial_values_in),
outer_boundary_velocity(outer_boundary_velocity_in),
outer_boundary_decay_time(outer_boundary_decay_time_in) {}

std::array<double, 2> initial_values{
std::numeric_limits<double>::signaling_NaN(),
std::numeric_limits<double>::signaling_NaN()};
double outer_boundary_velocity{
std::numeric_limits<double>::signaling_NaN()};
double outer_boundary_decay_time{
std::numeric_limits<double>::signaling_NaN()};
};

struct RotationMapOptions {
using type = Options::Auto<RotationMapOptions, Options::AutoLabel::None>;
static std::string name() { return "RotationMap"; }
static constexpr Options::String help = {
"Options for a time-dependent rotation map about an arbitrary axis. "
"Specify 'None' to not use this map."};

struct InitialAngularVelocity {
using type = std::array<double, 3>;
static constexpr Options::String help = {"The initial angular velocity."};
};

using options = tmpl::list<InitialAngularVelocity>;
using ExpansionMapOptions =
domain::creators::time_dependent_options::ExpansionMapOptions;

RotationMapOptions() = default;
explicit RotationMapOptions(
std::array<double, 3> initial_angular_velocity_in)
: initial_angular_velocity(initial_angular_velocity_in) {}

std::array<double, 3> initial_angular_velocity{};
};
/// \brief Options for the rotation map
using RotationMapOptions =
domain::creators::time_dependent_options::RotationMapOptions<3>;

/// \brief Options for the Translation Map, the outer radius is always set to
/// the outer boundary of the Domain, so there's no option needed for outer
/// boundary.
struct TranslationMapOptions {
using type = Options::Auto<TranslationMapOptions, Options::AutoLabel::None>;
static std::string name() { return "TranslationMap"; }
static constexpr Options::String help = {
"Options for a time-dependent translation map. Specify 'None' to not "
"use this map."};

struct InitialValues {
using type = std::array<std::array<double, 3>, 3>;
static constexpr Options::String help = {
"Initial position, velocity and acceleration."};
};

using options = tmpl::list<InitialValues>;
TranslationMapOptions() = default;
explicit TranslationMapOptions(
std::array<std::array<double, 3>, 3> initial_values_in)
: initial_values(initial_values_in) {}

std::array<std::array<double, 3>, 3> initial_values{};
};
using TranslationMapOptions =
domain::creators::time_dependent_options::TranslationMapOptions<3>;

// We use a type alias here instead of defining the ShapeMapOptions struct
// because there appears to be a bug in clang-10. If the definition of
// ShapeMapOptions is here inside TimeDependentMapOptions, on clang-10 there
// is a linking error that there is an undefined reference to
// Options::Option::parse_as<TimeDependentMapOptions<A>> (and B). This doesn't
// show up for GCC. If we put the definition of ShapeMapOptions outside of
// TimeDependentMapOptions and just use a type alias here, the linking error
// goes away.
/// \brief Options for the shape map
template <domain::ObjectLabel Object>
using ShapeMapOptions =
domain::creators::time_dependent_options::ShapeMapOptions<
Expand Down
Loading
Loading