Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <limits>
8 :
9 : #include "DataStructures/Tensor/TypeAliases.hpp"
10 : #include "Options/Options.hpp"
11 : #include "Options/String.hpp"
12 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
13 : #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
14 : #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
15 : #include "Utilities/Gsl.hpp"
16 : #include "Utilities/TMPL.hpp"
17 :
18 : /// \cond
19 : class DataVector;
20 :
21 : namespace PUP {
22 : class er;
23 : } // namespace PUP
24 : /// \endcond
25 :
26 : namespace VariableFixing {
27 :
28 : /*!
29 : * \ingroup VariableFixingGroup
30 : * \brief Adjust the electron fraction (Ye) based on rest mass density (rho).
31 : *
32 : * Based on \cite Liebendorfer2005 spherically symmetric Boltzmann calculations,
33 : * during the collapse phase just before a core-collapse supernova, the electron
34 : * fraction naturally follows the rest mass density. Intuitively, the higher
35 : * the density, the more electrons are captured onto protons, leading to
36 : * neutrino production. As these neutrinos diffuse out of the center of the
37 : * collapsing star, the total number of leptons near center of the CCSN
38 : * drops---the process of deleptonization.
39 : *
40 : * https://iopscience.iop.org/article/10.1086/466517
41 : */
42 1 : class ParameterizedDeleptonization {
43 : public:
44 : /// \brief Use an analytic expression vs rest mass density profile
45 : ///
46 : /// The option to choose between analytic and tabulated will be added in the
47 : /// future. If analytic is chosen (currently the only option), use the below
48 : /// parameters. If tabulated is chosen, a Ye (rho) profile will need to be
49 : /// provided.
50 : // struct DeleptonizationFormat {
51 : // using type = double;
52 : // static type lower_bound() { return 0.0; }
53 : // static constexpr Options::String help = {"Choose an 'analytic' or
54 : // 'tabulated' expression to express electron fraction as a function of
55 : // rest mass density."};
56 : // };
57 :
58 : /// \brief Enable parameterized deleptonizations
59 1 : struct Enable {
60 0 : using type = bool;
61 0 : static constexpr Options::String help = {
62 : "Whether or not to activate parameterized deleptonization for "
63 : "supernovae."};
64 : };
65 :
66 : /// \brief Density near the center of the supernova at bounce, above which
67 : /// the central Ye is assumed
68 : ///
69 : /// In practice: rho(high) ~ 2x10^13 g / cm^3
70 1 : struct HighDensityScale {
71 0 : using type = double;
72 0 : static type lower_bound() { return 0.0; }
73 0 : static constexpr Options::String help = {
74 : "High end of density scale for parameterized deleptonization."};
75 : };
76 : /// \brief Density near the Silicon-Oxygen interface, below which the lower
77 : /// Ye is assumed
78 : ///
79 : /// In practice: rho(low) ~ 2x10^7 g / cm^3
80 1 : struct LowDensityScale {
81 0 : using type = double;
82 0 : static type lower_bound() { return 0.0; }
83 0 : static constexpr Options::String help = {
84 : "Low end of density scale for parameterized deleptonization."};
85 : };
86 : /// \brief Electron fraction of material when the rest mass density is above
87 : /// HighDensityScale
88 : ///
89 : /// In practice: Y_e (high density) ~ 0.28.
90 1 : struct ElectronFractionAtHighDensity {
91 0 : using type = double;
92 0 : static type lower_bound() { return 0.0; }
93 0 : static type upper_bound() { return 0.5; }
94 0 : static constexpr Options::String help = {
95 : "For densities above HighDensityScale, the electron fraction will "
96 : "take this value."};
97 : };
98 : /// \brief Electron fraction of material when the rest mass density is below
99 : /// LowDensityScale
100 : ///
101 : /// In practice: Y_e (low density) ~ 0.5
102 1 : struct ElectronFractionAtLowDensity {
103 0 : using type = double;
104 0 : static type lower_bound() { return 0.0; }
105 0 : static type upper_bound() { return 0.5; }
106 0 : static constexpr Options::String help = {
107 : "For densities below LowDensityScale, the electron fraction will "
108 : "take this value."};
109 : };
110 :
111 : /// \brief Electron fraction correction term. The larger this value, the
112 : /// higher the Ye of matter at densities between LowDensityScale and
113 : /// HighDensityScale.
114 : ///
115 : /// In practice: Y_e (correction) ~ 0.035
116 1 : struct ElectronFractionCorrectionScale {
117 0 : using type = double;
118 0 : static type lower_bound() { return 0.0; }
119 0 : static constexpr Options::String help = {
120 : "For densities between low and high limits, a higher value of "
121 : "ElectronFractionCorrectionScale will increase the value of Ye."};
122 : };
123 :
124 0 : using options =
125 : tmpl::list<Enable, HighDensityScale, LowDensityScale,
126 : ElectronFractionAtHighDensity, ElectronFractionAtLowDensity,
127 : ElectronFractionCorrectionScale>;
128 0 : static constexpr Options::String help = {
129 : "Set electron fraction based on rest mass density. "
130 : "(Low/High)DensityScale sets the limits of the density, beyond which "
131 : "the ElectronFractionAt(Low/High)Density is assumed. At intermediate "
132 : "densities, the higher the ElectronFractionCorrectionScale, the higher "
133 : "the Ye."};
134 :
135 0 : ParameterizedDeleptonization(bool enable, double high_density_scale,
136 : double low_density_scale,
137 : double electron_fraction_at_high_density,
138 : double electron_fraction_at_low_density,
139 : double electron_fraction_correction_scale,
140 : const Options::Context& context = {});
141 :
142 0 : ParameterizedDeleptonization() = default;
143 0 : ParameterizedDeleptonization(const ParameterizedDeleptonization& /*rhs*/) =
144 : default;
145 0 : ParameterizedDeleptonization& operator=(
146 : const ParameterizedDeleptonization& /*rhs*/) = default;
147 0 : ParameterizedDeleptonization(ParameterizedDeleptonization&& /*rhs*/) =
148 : default;
149 0 : ParameterizedDeleptonization& operator=(
150 : ParameterizedDeleptonization&& /*rhs*/) = default;
151 0 : ~ParameterizedDeleptonization() = default;
152 :
153 : // NOLINTNEXTLINE(google-runtime-references)
154 0 : void pup(PUP::er& p);
155 :
156 : // Quantities being mutated
157 0 : using return_tags =
158 : tmpl::list<hydro::Tags::SpecificInternalEnergy<DataVector>,
159 : hydro::Tags::ElectronFraction<DataVector>,
160 : hydro::Tags::Pressure<DataVector>,
161 : hydro::Tags::Temperature<DataVector>>;
162 :
163 : // Things you want from DataBox that won't be change and are passed in as
164 : // const-refs
165 0 : using argument_tags = tmpl::list<hydro::Tags::RestMassDensity<DataVector>,
166 : hydro::Tags::GrmhdEquationOfState>;
167 :
168 : // for use in `db::mutate_apply`
169 : template <size_t ThermodynamicDim>
170 0 : void operator()(
171 : gsl::not_null<Scalar<DataVector>*> specific_internal_energy,
172 : gsl::not_null<Scalar<DataVector>*> electron_fraction,
173 : gsl::not_null<Scalar<DataVector>*> pressure,
174 : gsl::not_null<Scalar<DataVector>*> temperature,
175 : const Scalar<DataVector>& rest_mass_density,
176 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
177 : equation_of_state) const;
178 :
179 : private:
180 : template <size_t ThermodynamicDim>
181 0 : void correct_electron_fraction(
182 : gsl::not_null<Scalar<DataVector>*> specific_internal_energy,
183 : gsl::not_null<Scalar<DataVector>*> electron_fraction,
184 : gsl::not_null<Scalar<DataVector>*> pressure,
185 : gsl::not_null<Scalar<DataVector>*> temperature,
186 : const Scalar<DataVector>& rest_mass_density,
187 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
188 : equation_of_state,
189 : size_t grid_index) const;
190 :
191 : // NOLINTNEXTLINE(readability-redundant-declaration)
192 0 : friend bool operator==(const ParameterizedDeleptonization& lhs,
193 : const ParameterizedDeleptonization& rhs);
194 :
195 0 : bool enable_{false};
196 0 : double high_density_scale_{std::numeric_limits<double>::signaling_NaN()};
197 0 : double low_density_scale_{std::numeric_limits<double>::signaling_NaN()};
198 0 : double electron_fraction_at_high_density_{
199 : std::numeric_limits<double>::signaling_NaN()};
200 0 : double electron_fraction_at_low_density_{
201 : std::numeric_limits<double>::signaling_NaN()};
202 0 : double electron_fraction_half_sum_{
203 : std::numeric_limits<double>::signaling_NaN()};
204 0 : double electron_fraction_half_difference_{
205 : std::numeric_limits<double>::signaling_NaN()};
206 0 : double electron_fraction_correction_scale_{
207 : std::numeric_limits<double>::signaling_NaN()};
208 : };
209 0 : bool operator!=(const ParameterizedDeleptonization& lhs,
210 : const ParameterizedDeleptonization& rhs);
211 :
212 : } // namespace VariableFixing
|
