|
1 | 1 | #include "COWD.h" |
2 | 2 |
|
3 | | -/* For COWDs, calculate: |
4 | | - * |
5 | | - * (a) the maximum mass acceptance rate of this star, as the accretor, during mass transfer, and |
6 | | - * (b) the retention efficiency parameter |
7 | | - * |
8 | | - * |
9 | | - * For a given mass transfer rate, this function computes the amount of mass a WD would retain after |
10 | | - * flashes, as given by appendix B of Claeys+ 2014. |
11 | | - * https://ui.adsabs.harvard.edu/abs/2014A%26A...563A..83C/abstract |
12 | | - * |
13 | | - * |
14 | | - * DBL_DBL CalculateMassAcceptanceRate(const double p_DonorMassRate, const bool p_IsHeRich) |
15 | | - * |
16 | | - * @param [IN] p_DonorMassRate Mass transfer rate from the donor |
17 | | - * @param [IN] p_IsHeRich Material is He-rich or not |
18 | | - * @return Tuple containing the Maximum Mass Acceptance Rate (Msun/yr) and Retention Efficiency Parameter |
19 | | - */ |
20 | | -DBL_DBL COWD::CalculateMassAcceptanceRate(const double p_DonorMassRate, const bool p_IsHeRich) { |
21 | | - |
22 | | - m_AccretionRegime = DetermineAccretionRegime(p_IsHeRich, p_DonorMassRate); |
23 | | - |
24 | | - double acceptanceRate = 0.0; // acceptance mass rate - default = 0.0 |
25 | | - double fractionAccreted = 0.0; // accretion fraction - default = 0.0 |
26 | | - |
27 | | - acceptanceRate = p_DonorMassRate * CalculateEtaHe(p_DonorMassRate); |
28 | | - if (!p_IsHeRich) acceptanceRate *= CalculateEtaH(p_DonorMassRate); |
29 | | - |
30 | | - fractionAccreted = acceptanceRate / p_DonorMassRate; |
31 | | - |
32 | | - return std::make_tuple(acceptanceRate, fractionAccreted); |
33 | | -} |
34 | | - |
35 | | - |
36 | | -/* |
37 | | - * Determine the WD accretion regime based on the MT rate and whether the donor is He rich. Also, |
38 | | - * initialize He-Shell detonation or Off-center ignition when necessary, by changing the value |
39 | | - * of m_HeShellDetonation or m_OffCenterIgnition (respectively). |
40 | | - * |
41 | | - * The accretion regime is one of the options listed in enum ACCRETION_REGIME (constants.h) |
42 | | - * |
43 | | - * Note that we have merged the different flashes regimes from Piersanti+ 2014 into a single regime. |
44 | | - * |
45 | | - * ACCRETION_REGIME DetermineAccretionRegime(const bool p_HeRich, const double p_DonorMassLossRate) |
46 | | - * |
47 | | - * @param [IN] p_HeRich Whether the accreted material is helium-rich or not |
48 | | - * @param [IN] p_DonorMassLossRate Donor mass loss rate, in units of Msol / Myr |
49 | | - * @return Current WD accretion regime |
50 | | - */ |
51 | | -ACCRETION_REGIME COWD::DetermineAccretionRegime(const bool p_HeRich, const double p_DonorMassLossRate) { |
52 | | - |
53 | | - double logMdot = log10(p_DonorMassLossRate / MYR_TO_YEAR); // logarithm of the accreted mass (M_sun/yr) |
54 | | - ACCRETION_REGIME regime = ACCRETION_REGIME::ZERO; |
55 | | - |
56 | | - if (p_HeRich) { |
57 | | - // The following coefficients in logMassTransfer limits come from table A1 in Piersanti+ 2014. |
58 | | - double logMassTransferCrit = WD_LOG_MT_LIMIT_PIERSANTI_RG_SS_0 + WD_LOG_MT_LIMIT_PIERSANTI_RG_SS_1 * m_Mass; |
59 | | - double logMassTransferStable = WD_LOG_MT_LIMIT_PIERSANTI_SS_MF_0 + WD_LOG_MT_LIMIT_PIERSANTI_SS_MF_1 * m_Mass; // Piersanti+2014 has several Flashes regimes. Here we group them into one. |
60 | | - double logMassTransferDetonation = WD_LOG_MT_LIMIT_PIERSANTI_SF_Dt_0 + WD_LOG_MT_LIMIT_PIERSANTI_SF_Dt_1 * m_Mass; // critical value for double detonation regime in Piersanti+ 2014 |
61 | | - if (utils::Compare(logMdot, logMassTransferStable) < 0) { |
62 | | - if (utils::Compare(logMdot, logMassTransferDetonation) > 0) { |
63 | | - regime = ACCRETION_REGIME::HELIUM_FLASHES; |
64 | | - } |
65 | | - else { |
66 | | - regime = ACCRETION_REGIME::HELIUM_ACCUMULATION; |
67 | | - if ((utils::Compare(m_Mass, MASS_DOUBLE_DETONATION_CO) >= 0) && (utils::Compare(m_HeShell, WD_HE_SHELL_MCRIT_DETONATION) >= 0)) { |
68 | | - m_HeShellDetonation = true; |
69 | | - } |
70 | | - } |
71 | | - } |
72 | | - else if (utils::Compare(logMdot, logMassTransferCrit) > 0) { |
73 | | - regime = ACCRETION_REGIME::HELIUM_OPT_THICK_WINDS; |
74 | | - } |
75 | | - else { |
76 | | - regime = ACCRETION_REGIME::HELIUM_STABLE_BURNING; |
77 | | - if ((utils::Compare(logMdot, COWD_LOG_MDOT_MIN_OFF_CENTER_IGNITION) > 0) && (utils::Compare(m_Mass, COWD_MASS_MIN_OFF_CENTER_IGNITION) > 0)) { |
78 | | - m_OffCenterIgnition = true; |
79 | | - } |
80 | | - } |
81 | | - } |
82 | | - else { |
83 | | - // The following coefficients in logMassTransfer limits come from quadratic fits to Nomoto+ 2007 results (table 5) in Mass vs log10 Mdot space, to cover the low-mass end. |
84 | | - double m_Mass_2 = m_Mass * m_Mass; |
85 | | - double logMassTransferCrit = WD_LOG_MT_LIMIT_NOMOTO_REDGIANT_0 + WD_LOG_MT_LIMIT_NOMOTO_REDGIANT_1 * m_Mass + WD_LOG_MT_LIMIT_NOMOTO_REDGIANT_2 * m_Mass_2; |
86 | | - double logMassTransferStable = WD_LOG_MT_LIMIT_NOMOTO_STABLE_0 + WD_LOG_MT_LIMIT_NOMOTO_STABLE_1 * m_Mass + WD_LOG_MT_LIMIT_NOMOTO_STABLE_2 * m_Mass_2; |
87 | | - |
88 | | - if (utils::Compare(logMdot, logMassTransferStable) < 0) { |
89 | | - regime = ACCRETION_REGIME::HYDROGEN_FLASHES; |
90 | | - } |
91 | | - else if (utils::Compare(logMdot, logMassTransferCrit) > 0) { |
92 | | - regime = ACCRETION_REGIME::HYDROGEN_OPT_THICK_WINDS; |
93 | | - } |
94 | | - else { |
95 | | - regime = ACCRETION_REGIME::HYDROGEN_STABLE_BURNING; |
96 | | - } |
97 | | - } |
98 | | - |
99 | | - return regime; |
100 | | -} |
101 | | - |
102 | | - |
103 | 3 | /* |
104 | 4 | * Specifies next stage, if the star changes its phase. |
105 | 5 | * |
|
0 commit comments