Merge pull request #1326 from jeffriley/ch_spindown

Fixes initialisation in MS_gt_07::Initialise() for CORE_MASS_PRESCRIPTION::SHIKAUCHI
TeamCOMPAS · Jan 30, 2025 · acc2772 · acc2772
2 parents 43961fe + 3b22bc1
commit acc2772
Show file tree

Hide file tree

Showing 3 changed files with 16 additions and 11 deletions.
diff --git a/src/BaseBinaryStar.cpp b/src/BaseBinaryStar.cpp
@@ -2850,19 +2850,16 @@ void BaseBinaryStar::EmitGravitationalWave(const double p_Dt) {
  */
 double BaseBinaryStar::ChooseTimestep(const double p_Multiplier) {
 
-    // Timesteps required by individual stars
-    double dt1 = m_Star1->CalculateTimestep();
-    double dt2 = m_Star2->CalculateTimestep();
-    double dt = std::min(dt1, dt2);
+    double dt = std::min(m_Star1->CalculateTimestep(), m_Star2->CalculateTimestep());       // dt = smaller of timesteps required by individual stars
 
-    if (!IsUnbound()){                                                                      // Check that binary is bound
+    if (!IsUnbound()) {                                                                     // check that binary is bound
 
         if (OPTIONS->EmitGravitationalRadiation()) {                                        // emitting GWs?
             dt = std::min(dt, -1.0E-2 * m_SemiMajorAxis / m_DaDtGW);                        // yes - reduce timestep if necessary to ensure that the orbital separation does not change by more than ~1% per timestep due to GW emission
         }
 
-        if (OPTIONS->TidesPrescription() == TIDES_PRESCRIPTION::KAPIL2024) {                                // tides prescription = KAPIL2024
-                                                                                                            // yes - need to adjust dt
+        if (OPTIONS->TidesPrescription() == TIDES_PRESCRIPTION::KAPIL2024) {                // tides prescription = KAPIL2024
+                                                                                            // yes - need to adjust dt
 
             double omega                   = OrbitalAngularVelocity();
 
@@ -2898,7 +2895,7 @@ double BaseBinaryStar::ChooseTimestep(const double p_Multiplier) {
 
     dt *= p_Multiplier;	
 
-    return std::max(std::round(dt / TIMESTEP_QUANTUM) * TIMESTEP_QUANTUM, TIDES_MNIMUM_FRACTIONAL_NUCLEAR_TIME * NUCLEAR_MINIMUM_TIMESTEP);    // quantised and not less than minimum
+    return std::max(std::round(dt / TIMESTEP_QUANTUM) * TIMESTEP_QUANTUM, TIDES_MNIMUM_FRACTIONAL_NUCLEAR_TIME * NUCLEAR_MINIMUM_TIMESTEP); // quantised and not less than minimum
 }
 
 

diff --git a/src/MS_gt_07.h b/src/MS_gt_07.h
@@ -43,7 +43,11 @@ class MS_gt_07: virtual public BaseStar, public MainSequence {
         // Age for MS_GT_07 is carried over from CH stars switching to MS after spinning down, so not set to 0.0 here
 
         // Initialise core mass, luminosity, radius, and temperature if Shikauchi core mass prescription is used
-        if ((OPTIONS->MainSequenceCoreMassPrescription() == CORE_MASS_PRESCRIPTION::SHIKAUCHI) && (utils::Compare(m_MZAMS, SHIKAUCHI_LOWER_MASS_LIMIT) >= 0)) {
+        // Only do this once - this should not be done if a CH star spins down and becomes a MS star (when using CHE_MODE::PESSIMISTIC)
+        if (OPTIONS->MainSequenceCoreMassPrescription() == CORE_MASS_PRESCRIPTION::SHIKAUCHI &&                                                         // Shikauchi core mass prescription?
+            utils::Compare(m_MZAMS, SHIKAUCHI_LOWER_MASS_LIMIT) >= 0                         &&                                                         // ZAMS mass >= SHIKAUCHI_LOWER_MASS_LIMIT?
+            m_Time <= 0.0) {                                                                                                                            // star not yet aged past creation?
+                                                                                                                                                        // yes - initialise
             m_InitialMainSequenceCoreMass = MainSequence::CalculateInitialMainSequenceCoreMass(m_MZAMS);
             m_MainSequenceCoreMass        = m_InitialMainSequenceCoreMass;
             m_Luminosity                  = MainSequence::CalculateLuminosityShikauchi(m_MainSequenceCoreMass, m_InitialHeliumAbundance, m_Age);

diff --git a/src/changelog.h b/src/changelog.h
@@ -1444,9 +1444,13 @@
 // 03.12.01   JR - Jan 17, 2025     - Defect repair:
 //                                      - (partial?) fix for issue #1149 - remove conditional from TPAGB::IsSupernova().  Whether it fixes issue 1149 completely or not, the conditional shouldn't be there...
 // 03.12.02   SS - Jan 20, 2025     - Defect repair:
-//                                      -  fix for issue #1324 - allow for evolution of pulsars formed from main sequence merger products. Changed CalculateTimestep to ChooseTimestep in MR.h and added check for bound binary to BaseBinaryStar::ChooseTimestep
+//                                      - fix for issue #1324 - allow for evolution of pulsars formed from main sequence merger products.
+//                                        Changed CalculateTimestep to ChooseTimestep in MR.h and added check for bound binary to BaseBinaryStar::ChooseTimestep
+// 03.12.03   JR - Jan 29, 2025     - Defect repair:
+//                                      - fixes initialisation in MS_gt_07::Initialise() for CORE_MASS_PRESCRIPTION::SHIKAUCHI (now allows for CH stars that spin down)
+//                                      - minor code cleanup
 
 
-const std::string VERSION_STRING = "03.12.02";
+const std::string VERSION_STRING = "03.12.03";
 
 # endif // __changelog_h__