4237779830

Build Type

push

github

Committed by GitHub

Commit Message

Merge pull request #4755 from nilsvu/test_bbh_domain

Run Details

63865 of 66590 relevant lines covered (95.91%)

411287.72 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

86.42

/src/PointwiseFunctions/AnalyticData/GrMhd/MagneticFieldLoop.cpp

// Distributed under the MIT License.
// See LICENSE.txt for details.

#include "PointwiseFunctions/AnalyticData/GrMhd/MagneticFieldLoop.hpp"

#include <cmath>  // IWYU pragma: keep
#include <cstddef>
#include <ostream>
#include <pup.h>

#include "DataStructures/DataVector.hpp"  // IWYU pragma: keep
#include "DataStructures/Tensor/EagerMath/DotProduct.hpp"
#include "DataStructures/Tensor/Tensor.hpp"
#include "PointwiseFunctions/Hydro/LorentzFactor.hpp"
#include "PointwiseFunctions/Hydro/SpecificEnthalpy.hpp"
#include "PointwiseFunctions/Hydro/Tags.hpp"
#include "Utilities/ConstantExpressions.hpp"  // IWYU pragma: keep
#include "Utilities/ContainerHelpers.hpp"
#include "Utilities/GenerateInstantiations.hpp"
#include "Utilities/MakeWithValue.hpp"
#include "Utilities/Math.hpp"  // IWYU pragma: keep
#include "Utilities/StdArrayHelpers.hpp"
#include "Utilities/StdHelpers.hpp"

// IWYU pragma: no_include <complex>

namespace grmhd::AnalyticData {

MagneticFieldLoop::MagneticFieldLoop(
    const double pressure, const double rest_mass_density,
    const double adiabatic_index,
    const std::array<double, 3>& advection_velocity,
    const double magnetic_field_magnitude, const double inner_radius,
    const double outer_radius, const Options::Context& context)
    : pressure_(pressure),
      rest_mass_density_(rest_mass_density),
      adiabatic_index_(adiabatic_index),
      advection_velocity_(advection_velocity),
      magnetic_field_magnitude_(magnetic_field_magnitude),
      inner_radius_(inner_radius),
      outer_radius_(outer_radius),
      equation_of_state_{adiabatic_index_} {
  if (magnitude(advection_velocity) >= 1.0) {
    PARSE_ERROR(context, "MagneticFieldLoop: superluminal AdvectionVelocity = "
                             << advection_velocity);
  }
  if (inner_radius >= outer_radius) {
    PARSE_ERROR(context, "MagneticFieldLoop: InnerRadius of "
                             << inner_radius
                             << " is not less than OuterRadius of "
                             << outer_radius);
  }
}

std::unique_ptr<evolution::initial_data::InitialData>
MagneticFieldLoop::get_clone() const {
  return std::make_unique<MagneticFieldLoop>(*this);
}

MagneticFieldLoop::MagneticFieldLoop(CkMigrateMessage* msg)
    : InitialData(msg) {}

void MagneticFieldLoop::pup(PUP::er& p) {
  InitialData::pup(p);
  p | pressure_;
  p | rest_mass_density_;
  p | adiabatic_index_;
  p | advection_velocity_;
  p | magnetic_field_magnitude_;
  p | inner_radius_;
  p | outer_radius_;
  p | equation_of_state_;
  p | background_spacetime_;
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::RestMassDensity<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::RestMassDensity<DataType>> /*meta*/) const {
  return {make_with_value<Scalar<DataType>>(x, rest_mass_density_)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::ElectronFraction<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::ElectronFraction<DataType>> /*meta*/) const {

  return {make_with_value<Scalar<DataType>>(x, 0.1)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::SpatialVelocity<DataType, 3>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::SpatialVelocity<DataType, 3>> /*meta*/) const {
  auto result = make_with_value<tnsr::I<DataType, 3>>(x, 0.0);
  get<0>(result) = advection_velocity_[0];
  get<1>(result) = advection_velocity_[1];
  get<2>(result) = advection_velocity_[2];
  return {std::move(result)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::SpecificInternalEnergy<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>> /*meta*/) const {
  return equation_of_state_.specific_internal_energy_from_density_and_pressure(
      get<hydro::Tags::RestMassDensity<DataType>>(
          variables(x, tmpl::list<hydro::Tags::RestMassDensity<DataType>>{})),
      get<hydro::Tags::Pressure<DataType>>(
          variables(x, tmpl::list<hydro::Tags::Pressure<DataType>>{})));
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::Pressure<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::Pressure<DataType>> /*meta*/) const {
  return {make_with_value<Scalar<DataType>>(x, pressure_)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::MagneticField<DataType, 3>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::MagneticField<DataType, 3>> /*meta*/) const {
  auto result = make_with_value<tnsr::I<DataType, 3>>(x, 0.0);
  const DataType cylindrical_radius =
      sqrt(square(get<0>(x)) + square(get<1>(x)));
  for (size_t s = 0; s < get_size(cylindrical_radius); ++s) {
    const double current_radius = get_element(cylindrical_radius, s);
    // The data is ill-posed at the origin, so it is not clear what to do...
    if (current_radius <= outer_radius_ and current_radius > inner_radius_) {
      get_element(get<0>(result), s) = -magnetic_field_magnitude_ *
                                       get_element(get<1>(x), s) /
                                       current_radius;
      get_element(get<1>(result), s) = magnetic_field_magnitude_ *
                                       get_element(get<0>(x), s) /
                                       current_radius;
    }
    if (current_radius <= inner_radius_) {
      get_element(get<0>(result), s) = -magnetic_field_magnitude_ *
                                       get_element(get<1>(x), s) /
                                       inner_radius_;
      get_element(get<1>(result), s) =
          magnetic_field_magnitude_ * get_element(get<0>(x), s) / inner_radius_;
    }
  }
  return {std::move(result)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::DivergenceCleaningField<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::DivergenceCleaningField<DataType>> /*meta*/) const {
  return {make_with_value<Scalar<DataType>>(x, 0.0)};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::LorentzFactor<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::LorentzFactor<DataType>> /*meta*/) const {
  using velocity_tag = hydro::Tags::SpatialVelocity<DataType, 3>;
  const auto velocity =
      get<velocity_tag>(variables(x, tmpl::list<velocity_tag>{}));
  return {hydro::lorentz_factor(dot_product(velocity, velocity))};
}

template <typename DataType>
tuples::TaggedTuple<hydro::Tags::SpecificEnthalpy<DataType>>
MagneticFieldLoop::variables(
    const tnsr::I<DataType, 3>& x,
    tmpl::list<hydro::Tags::SpecificEnthalpy<DataType>> /*meta*/) const {
  return hydro::relativistic_specific_enthalpy(
      get<hydro::Tags::RestMassDensity<DataType>>(
          variables(x, tmpl::list<hydro::Tags::RestMassDensity<DataType>>{})),
      get<hydro::Tags::SpecificInternalEnergy<DataType>>(variables(
          x, tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>>{})),
      get<hydro::Tags::Pressure<DataType>>(
          variables(x, tmpl::list<hydro::Tags::Pressure<DataType>>{})));
}

PUP::able::PUP_ID MagneticFieldLoop::my_PUP_ID = 0;

bool operator==(const MagneticFieldLoop& lhs, const MagneticFieldLoop& rhs) {
  // there is no comparison operator for the EoS, but should be okay as
  // the adiabatic_indexs are compared
  return lhs.pressure_ == rhs.pressure_ and
         lhs.rest_mass_density_ == rhs.rest_mass_density_ and
         lhs.adiabatic_index_ == rhs.adiabatic_index_ and
         lhs.advection_velocity_ == rhs.advection_velocity_ and
         lhs.magnetic_field_magnitude_ == rhs.magnetic_field_magnitude_ and
         lhs.inner_radius_ == rhs.inner_radius_ and
         lhs.outer_radius_ == rhs.outer_radius_ and
         lhs.background_spacetime_ == rhs.background_spacetime_;
}

bool operator!=(const MagneticFieldLoop& lhs, const MagneticFieldLoop& rhs) {
  return not(lhs == rhs);
}

#define DTYPE(data) BOOST_PP_TUPLE_ELEM(0, data)
#define TAG(data) BOOST_PP_TUPLE_ELEM(1, data)

#define INSTANTIATE_SCALARS(_, data)                     \
  template tuples::TaggedTuple<TAG(data) < DTYPE(data)>> \
      MagneticFieldLoop::variables(                      \
          const tnsr::I<DTYPE(data), 3>& x,              \
          tmpl::list<TAG(data) < DTYPE(data)>> /*meta*/) const;

GENERATE_INSTANTIATIONS(
    INSTANTIATE_SCALARS, (double, DataVector),
    (hydro::Tags::RestMassDensity, hydro::Tags::ElectronFraction,
     hydro::Tags::SpecificInternalEnergy, hydro::Tags::Pressure,
     hydro::Tags::DivergenceCleaningField, hydro::Tags::LorentzFactor,
     hydro::Tags::SpecificEnthalpy))

#define INSTANTIATE_VECTORS(_, data)                        \
  template tuples::TaggedTuple<TAG(data) < DTYPE(data), 3>> \
      MagneticFieldLoop::variables(                         \
          const tnsr::I<DTYPE(data), 3>& x,                 \
          tmpl::list<TAG(data) < DTYPE(data), 3>> /*meta*/) const;

GENERATE_INSTANTIATIONS(INSTANTIATE_VECTORS, (double, DataVector),
                        (hydro::Tags::SpatialVelocity,
                         hydro::Tags::MagneticField))

#undef DTYPE
#undef TAG
#undef INSTANTIATE_SCALARS
#undef INSTANTIATE_VECTORS
}  // namespace grmhd::AnalyticData

1	// Distributed under the MIT License.
2	// See LICENSE.txt for details.
3
4	#include "PointwiseFunctions/AnalyticData/GrMhd/MagneticFieldLoop.hpp"
5
6	#include <cmath> // IWYU pragma: keep
7	#include <cstddef>
8	#include <ostream>
9	#include <pup.h>
10
11	#include "DataStructures/DataVector.hpp" // IWYU pragma: keep
12	#include "DataStructures/Tensor/EagerMath/DotProduct.hpp"
13	#include "DataStructures/Tensor/Tensor.hpp"
14	#include "PointwiseFunctions/Hydro/LorentzFactor.hpp"
15	#include "PointwiseFunctions/Hydro/SpecificEnthalpy.hpp"
16	#include "PointwiseFunctions/Hydro/Tags.hpp"
17	#include "Utilities/ConstantExpressions.hpp" // IWYU pragma: keep
18	#include "Utilities/ContainerHelpers.hpp"
19	#include "Utilities/GenerateInstantiations.hpp"
20	#include "Utilities/MakeWithValue.hpp"
21	#include "Utilities/Math.hpp" // IWYU pragma: keep
22	#include "Utilities/StdArrayHelpers.hpp"
23	#include "Utilities/StdHelpers.hpp"
24
25	// IWYU pragma: no_include <complex>
26
27	namespace grmhd::AnalyticData {
28
29	MagneticFieldLoop::MagneticFieldLoop(	8✔
30	const double pressure, const double rest_mass_density,
31	const double adiabatic_index,
32	const std::array<double, 3>& advection_velocity,
33	const double magnetic_field_magnitude, const double inner_radius,
34	const double outer_radius, const Options::Context& context)	8✔
35	: pressure_(pressure),
36	rest_mass_density_(rest_mass_density),
37	adiabatic_index_(adiabatic_index),
38	advection_velocity_(advection_velocity),
39	magnetic_field_magnitude_(magnetic_field_magnitude),
40	inner_radius_(inner_radius),
41	outer_radius_(outer_radius),
42	equation_of_state_{adiabatic_index_} {	10✔
43	if (magnitude(advection_velocity) >= 1.0) {	8✔
44	PARSE_ERROR(context, "MagneticFieldLoop: superluminal AdvectionVelocity = "	9✔
45	<< advection_velocity);
46	}
47	if (inner_radius >= outer_radius) {	7✔
48	PARSE_ERROR(context, "MagneticFieldLoop: InnerRadius of "	×
49	<< inner_radius
50	<< " is not less than OuterRadius of "
51	<< outer_radius);
52	}
53	}	7✔
54
55	std::unique_ptr<evolution::initial_data::InitialData>
56	MagneticFieldLoop::get_clone() const {	1✔
57	return std::make_unique<MagneticFieldLoop>(*this);	1✔
58	}
59
60	MagneticFieldLoop::MagneticFieldLoop(CkMigrateMessage* msg)	1✔
61	: InitialData(msg) {}	1✔
62
63	void MagneticFieldLoop::pup(PUP::er& p) {	6✔
64	InitialData::pup(p);	6✔
65	p \| pressure_;	6✔
66	p \| rest_mass_density_;	6✔
67	p \| adiabatic_index_;	6✔
68	p \| advection_velocity_;	6✔
69	p \| magnetic_field_magnitude_;	6✔
70	p \| inner_radius_;	6✔
71	p \| outer_radius_;	6✔
72	p \| equation_of_state_;	6✔
73	p \| background_spacetime_;	6✔
74	}	6✔
75
76	template <typename DataType>
77	tuples::TaggedTuple<hydro::Tags::RestMassDensity<DataType>>
78	MagneticFieldLoop::variables(	128✔
79	const tnsr::I<DataType, 3>& x,
80	tmpl::list<hydro::Tags::RestMassDensity<DataType>> /meta/) const {
81	return {make_with_value<Scalar<DataType>>(x, rest_mass_density_)};	256✔
82	}
83
84	template <typename DataType>
85	tuples::TaggedTuple<hydro::Tags::ElectronFraction<DataType>>
86	MagneticFieldLoop::variables(	32✔
87	const tnsr::I<DataType, 3>& x,
88	tmpl::list<hydro::Tags::ElectronFraction<DataType>> /meta/) const {
89
90	return {make_with_value<Scalar<DataType>>(x, 0.1)};	64✔
91	}
92
93	template <typename DataType>
94	tuples::TaggedTuple<hydro::Tags::SpatialVelocity<DataType, 3>>
95	MagneticFieldLoop::variables(	64✔
96	const tnsr::I<DataType, 3>& x,
97	tmpl::list<hydro::Tags::SpatialVelocity<DataType, 3>> /meta/) const {
98	auto result = make_with_value<tnsr::I<DataType, 3>>(x, 0.0);	96✔
99	get<0>(result) = advection_velocity_[0];	96✔
100	get<1>(result) = advection_velocity_[1];	96✔
101	get<2>(result) = advection_velocity_[2];	96✔
102	return {std::move(result)};	128✔
103	}
104
105	template <typename DataType>
106	tuples::TaggedTuple<hydro::Tags::SpecificInternalEnergy<DataType>>
107	MagneticFieldLoop::variables(	64✔
108	const tnsr::I<DataType, 3>& x,
109	tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>> /meta/) const {
110	return equation_of_state_.specific_internal_energy_from_density_and_pressure(
111	get<hydro::Tags::RestMassDensity<DataType>>(	64✔
112	variables(x, tmpl::list<hydro::Tags::RestMassDensity<DataType>>{})),
113	get<hydro::Tags::Pressure<DataType>>(	96✔
114	variables(x, tmpl::list<hydro::Tags::Pressure<DataType>>{})));	160✔
115	}
116
117	template <typename DataType>
118	tuples::TaggedTuple<hydro::Tags::Pressure<DataType>>
119	MagneticFieldLoop::variables(	128✔
120	const tnsr::I<DataType, 3>& x,
121	tmpl::list<hydro::Tags::Pressure<DataType>> /meta/) const {
122	return {make_with_value<Scalar<DataType>>(x, pressure_)};	256✔
123	}
124
125	template <typename DataType>
126	tuples::TaggedTuple<hydro::Tags::MagneticField<DataType, 3>>
127	MagneticFieldLoop::variables(	18✔
128	const tnsr::I<DataType, 3>& x,
129	tmpl::list<hydro::Tags::MagneticField<DataType, 3>> /meta/) const {
130	auto result = make_with_value<tnsr::I<DataType, 3>>(x, 0.0);	45✔
131	const DataType cylindrical_radius =	27✔
132	sqrt(square(get<0>(x)) + square(get<1>(x)));	63✔
133	for (size_t s = 0; s < get_size(cylindrical_radius); ++s) {	144✔
134	const double current_radius = get_element(cylindrical_radius, s);	54✔
135	// The data is ill-posed at the origin, so it is not clear what to do...
136	if (current_radius <= outer_radius_ and current_radius > inner_radius_) {	54✔
137	get_element(get<0>(result), s) = -magnetic_field_magnitude_ *	×
138	get_element(get<1>(x), s) /	×
139	current_radius;
140	get_element(get<1>(result), s) = magnetic_field_magnitude_ *	×
141	get_element(get<0>(x), s) /	×
142	current_radius;
143	}
144	if (current_radius <= inner_radius_) {	54✔
145	get_element(get<0>(result), s) = -magnetic_field_magnitude_ *	×
146	get_element(get<1>(x), s) /	×
147	inner_radius_;	×
148	get_element(get<1>(result), s) =	×
149	magnetic_field_magnitude_ * get_element(get<0>(x), s) / inner_radius_;	×
150	}
151	}
152	return {std::move(result)};	36✔
153	}
154
155	template <typename DataType>
156	tuples::TaggedTuple<hydro::Tags::DivergenceCleaningField<DataType>>
157	MagneticFieldLoop::variables(	18✔
158	const tnsr::I<DataType, 3>& x,
159	tmpl::list<hydro::Tags::DivergenceCleaningField<DataType>> /meta/) const {
160	return {make_with_value<Scalar<DataType>>(x, 0.0)};	36✔
161	}
162
163	template <typename DataType>
164	tuples::TaggedTuple<hydro::Tags::LorentzFactor<DataType>>
165	MagneticFieldLoop::variables(	32✔
166	const tnsr::I<DataType, 3>& x,
167	tmpl::list<hydro::Tags::LorentzFactor<DataType>> /meta/) const {
168	using velocity_tag = hydro::Tags::SpatialVelocity<DataType, 3>;
169	const auto velocity =	16✔
170	get<velocity_tag>(variables(x, tmpl::list<velocity_tag>{}));	32✔
171	return {hydro::lorentz_factor(dot_product(velocity, velocity))};	48✔
172	}
173
174	template <typename DataType>
175	tuples::TaggedTuple<hydro::Tags::SpecificEnthalpy<DataType>>
176	MagneticFieldLoop::variables(	32✔
177	const tnsr::I<DataType, 3>& x,
178	tmpl::list<hydro::Tags::SpecificEnthalpy<DataType>> /meta/) const {
179	return hydro::relativistic_specific_enthalpy(	×
180	get<hydro::Tags::RestMassDensity<DataType>>(	32✔
181	variables(x, tmpl::list<hydro::Tags::RestMassDensity<DataType>>{})),
182	get<hydro::Tags::SpecificInternalEnergy<DataType>>(variables(	48✔
183	x, tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>>{})),
184	get<hydro::Tags::Pressure<DataType>>(	48✔
185	variables(x, tmpl::list<hydro::Tags::Pressure<DataType>>{})));	48✔
186	}
187
188	PUP::able::PUP_ID MagneticFieldLoop::my_PUP_ID = 0;
189
190	bool operator==(const MagneticFieldLoop& lhs, const MagneticFieldLoop& rhs) {	5✔
191	// there is no comparison operator for the EoS, but should be okay as
192	// the adiabatic_indexs are compared
193	return lhs.pressure_ == rhs.pressure_ and	10✔
194	lhs.rest_mass_density_ == rhs.rest_mass_density_ and	5✔
195	lhs.adiabatic_index_ == rhs.adiabatic_index_ and	5✔
196	lhs.advection_velocity_ == rhs.advection_velocity_ and	5✔
197	lhs.magnetic_field_magnitude_ == rhs.magnetic_field_magnitude_ and	5✔
198	lhs.inner_radius_ == rhs.inner_radius_ and	5✔
199	lhs.outer_radius_ == rhs.outer_radius_ and	15✔
200	lhs.background_spacetime_ == rhs.background_spacetime_;	10✔
201	}
202
203	bool operator!=(const MagneticFieldLoop& lhs, const MagneticFieldLoop& rhs) {	1✔
204	return not(lhs == rhs);	1✔
205	}
206
207	#define DTYPE(data) BOOST_PP_TUPLE_ELEM(0, data)
208	#define TAG(data) BOOST_PP_TUPLE_ELEM(1, data)
209
210	#define INSTANTIATE_SCALARS(_, data) \
211	template tuples::TaggedTuple<TAG(data) < DTYPE(data)>> \
212	MagneticFieldLoop::variables( \
213	const tnsr::I<DTYPE(data), 3>& x, \
214	tmpl::list<TAG(data) < DTYPE(data)>> /meta/) const;
215
216	GENERATE_INSTANTIATIONS(
217	INSTANTIATE_SCALARS, (double, DataVector),
218	(hydro::Tags::RestMassDensity, hydro::Tags::ElectronFraction,
219	hydro::Tags::SpecificInternalEnergy, hydro::Tags::Pressure,
220	hydro::Tags::DivergenceCleaningField, hydro::Tags::LorentzFactor,
221	hydro::Tags::SpecificEnthalpy))
222
223	#define INSTANTIATE_VECTORS(_, data) \
224	template tuples::TaggedTuple<TAG(data) < DTYPE(data), 3>> \
225	MagneticFieldLoop::variables( \
226	const tnsr::I<DTYPE(data), 3>& x, \
227	tmpl::list<TAG(data) < DTYPE(data), 3>> /meta/) const;
228
229	GENERATE_INSTANTIATIONS(INSTANTIATE_VECTORS, (double, DataVector),
230	(hydro::Tags::SpatialVelocity,
231	hydro::Tags::MagneticField))
232
233	#undef DTYPE
234	#undef TAG
235	#undef INSTANTIATE_SCALARS
236	#undef INSTANTIATE_VECTORS
237	} // namespace grmhd::AnalyticData

sxs-collaboration / spectre / 4237779830

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous