Program Listing for File Specsyn.hpp
↰ Return to documentation for file (src/specsyn/Specsyn.hpp)
#ifndef SPECSYN_HPP
#define SPECSYN_HPP
#include "../interpolation/Interpolator1D.hpp"
#include "../pdfs/PDF.hpp"
#include "../tracks/TrackCommons.hpp"
#include "../utils/Constants.hpp"
#include <array>
#include <cmath>
#include <cstddef>
#include <limits>
#include <memory>
#include <numbers>
#include <stdexcept>
#include <string>
#include <utility>
#include <vector>
namespace io
{
class SimControls;
} // namespace io
namespace specsyn
{
class Specsyn
{
public:
// Shorten type names: Segment is a single isochrone segment
// (an element of the Isochrone returned by
// Tracks2D::getIsochrone -- an isochrone can have more than
// one disjoint segment for non-monotonic tracks), and
// StarData is the type returned by evaluating a Segment at a
// given mass
using Segment = interp::Interpolator1D<static_cast<size_t>(tracks::FieldIdx::nTrackQty)>;
using Isochrone = std::vector<std::unique_ptr<Segment>>;
using StarData = std::array<double,
static_cast<size_t>(tracks::FieldIdx::nTrackQty)>;
explicit Specsyn(const io::SimControls& controls) :
controls_(controls) { }
virtual ~Specsyn() = default;
// Copyable and movable (matching the C++ reference semantics
// of controls_, which just gets rebound to the same referent
// on copy), but never actually copied/assigned in practice --
// every Specsyn is held via unique_ptr (see SimControls's own
// specsyn_ and SpecsynLibChained's libs_)
Specsyn(const Specsyn&) = default;
Specsyn(Specsyn&&) = default;
auto operator=(const Specsyn&) -> Specsyn& = delete;
auto operator=(Specsyn&&) -> Specsyn& = delete;
[[nodiscard]] auto wl() const -> const std::vector<double>& { return wl_; }
[[nodiscard]] auto controls() const -> const io::SimControls& { return controls_; }
[[nodiscard]] auto wlObs() const -> std::vector<double>;
[[nodiscard]] virtual auto fehMin() const -> double
{ return -std::numeric_limits<double>::infinity(); }
[[nodiscard]] virtual auto fehMax() const -> double
{ return std::numeric_limits<double>::infinity(); }
[[nodiscard]] auto requestedFehMin() const -> double { return requestedFehMin_; }
[[nodiscard]] auto requestedFehMax() const -> double { return requestedFehMax_; }
[[nodiscard]] virtual auto loggMin() const -> double
{ return -std::numeric_limits<double>::infinity(); }
[[nodiscard]] virtual auto loggMax() const -> double
{ return std::numeric_limits<double>::infinity(); }
[[nodiscard]] virtual auto spec(const StarData& props, double feh) const
-> std::vector<double> = 0;
[[nodiscard]] virtual auto specForce(const StarData& props, double feh) const
-> std::vector<double>
{
auto result = spec(props, feh);
if (result.empty())
{
throw std::runtime_error(
"Specsyn::specForce: unable to force a spectrum for this star");
}
return result;
}
[[nodiscard]] virtual auto specForIntegration(const StarData& props, double feh) const
-> std::vector<double>
{ return spec(props, feh); }
[[nodiscard]] auto spec(double m, const Segment& segment, double feh) const -> std::vector<double>
{
return spec(segment(m), feh);
}
[[nodiscard]] auto specWl(double m, const Segment& segment, double feh) const -> std::vector<double>
{
auto result = spec(m, segment, feh);
for (std::size_t i = 0; i < result.size(); ++i)
{
result[i] *= wl_[i]; // NOLINT(cppcoreguidelines-pro-bounds-avoid-unchecked-container-access) -- result has the same size as wl_ (one entry per wavelength), and i is bounded by result.size()
}
return result;
}
[[nodiscard]] auto specWlForIntegration(const StarData& props, double feh) const -> std::vector<double>
{
auto result = specForIntegration(props, feh);
for (std::size_t i = 0; i < result.size(); ++i)
{
result[i] *= wl_[i]; // NOLINT(cppcoreguidelines-pro-bounds-avoid-unchecked-container-access) -- result has the same size as wl_ (one entry per wavelength), and i is bounded by result.size()
}
return result;
}
[[nodiscard]] virtual auto specCts(
const Isochrone& isochrone,
const pdfs::PDF& imf,
double mTot,
double mMin,
double mMax,
double feh
) const -> std::vector<double>;
[[nodiscard]] auto specCts(
const pdfs::PDF& sfr,
const pdfs::PDF& imf,
const pdfs::PDF& fehDist,
double curTime,
double fCluster,
double mMin,
double mMax
) const -> std::vector<double>;
[[nodiscard]] auto specAndLbolCts(
const pdfs::PDF& sfr,
const pdfs::PDF& imf,
const pdfs::PDF& fehDist,
double curTime,
double fCluster,
double mMin,
double mMax
) const -> std::pair<std::vector<double>, double>;
[[nodiscard]] auto intRelTol() const -> double;
[[nodiscard]] auto intAbsTol() const -> double;
[[nodiscard]] auto intMaxIter() const -> std::size_t;
protected:
[[nodiscard]] static auto describeStar(const StarData& props, double feh) -> std::string;
[[nodiscard]] static auto getSAandLogg(const StarData& props) -> std::pair<double, double>
{
constexpr double pi = std::numbers::pi_v<double>;
const double logL = props[static_cast<size_t>(tracks::FieldIdx::logL)]; // NOLINT(cppcoreguidelines-pro-bounds-avoid-unchecked-container-access) -- StarData is a fixed-size std::array, and logL is one of its compile-time-known indices
const double logTeff = props[static_cast<size_t>(tracks::FieldIdx::logTe)]; // NOLINT(cppcoreguidelines-pro-bounds-avoid-unchecked-container-access) -- see above
const double mass = props[static_cast<size_t>(tracks::FieldIdx::mass)]; // Msun // NOLINT(cppcoreguidelines-pro-bounds-avoid-unchecked-container-access) -- see above
const double temperature = std::pow(10.0, logTeff); // K
const double luminosity = std::pow(10.0, logL) * utils::Lsun; // erg/s
const double temperature4 = temperature * temperature * temperature * temperature;
const double radius = std::sqrt(luminosity / (4.0 * pi * utils::sigmaSB * temperature4)); // cm
const double area = 4.0 * pi * radius * radius; // cm^2
const double g = utils::G * mass * utils::Msun / (radius * radius); // cm/s^2
const double logg = std::log10(g);
return { area, logg };
}
const io::SimControls& controls_; // NOLINT(cppcoreguidelines-avoid-const-or-ref-data-members) -- deliberately a live reference, not a copy: see this member's own comment for why (controls_'s own tolerances/redshift must be readable live, not snapshotted). This class is only ever used through the same non-copyable, non-movable ownership pattern (unique_ptr in SimControls's own specsyn_) as every other class with a reference member in this codebase (e.g. SpecsynLibNoWind's FeH_/logg_/logTeff_), so the usual objection (disabling implicit copy/move assignment) doesn't apply in practice.
std::vector<double> wl_;
void setRequestedFehRange(const double fehMin, const double fehMax)
{
requestedFehMin_ = fehMin;
requestedFehMax_ = fehMax;
}
private:
double requestedFehMin_ = -std::numeric_limits<double>::infinity();
double requestedFehMax_ = std::numeric_limits<double>::infinity();
[[nodiscard]] auto specCtsImpl(
const Isochrone& isochrone,
const pdfs::PDF& imf,
double mTot,
double mMin,
double mMax,
double feh,
bool forIntegration,
bool computeLbol
) const -> std::vector<double>;
[[nodiscard]] auto specCtsHelper(
const pdfs::PDF& sfr,
const pdfs::PDF& imf,
const pdfs::PDF& fehDist,
double curTime,
double fCluster,
double mMin,
double mMax,
bool computeLbol
) const -> std::vector<double>;
[[nodiscard]] auto continuousSpecIntegrand(
double age, const pdfs::PDF& imf, double mMin, double mMax,
bool computeLbol, double feh) const -> std::vector<double>;
};
} // namespace specsyn
#endif // SPECSYN_HPP