Merge pull request #1047 from jeffriley/issue-1034

Requested changes
TeamCOMPAS · Dec 28, 2023 · 12f58ea · 12f58ea
2 parents 2624083 + 175fccb
commit 12f58ea
Show file tree

Hide file tree

Showing 6 changed files with 26 additions and 13 deletions.
diff --git a/src/BH.h b/src/BH.h
@@ -55,7 +55,6 @@ class BH: virtual public BaseStar, public Remnants {
     double  CalculateMassLossRate()                                                 { return 0.0; }                                                     // Ensure that BHs don't lose mass in winds
 
     double  CalculateMomentOfInertia() const                                        { return (2.0 / 5.0) * m_Mass * m_Radius * m_Radius; }
-    double  CalculateMomentOfInertiaAU() const                                      { return CalculateMomentOfInertia() * RSOL_TO_AU * RSOL_TO_AU; }
 
     double  CalculateRadiusOnPhase() const                                          { return CalculateRadiusOnPhase_Static(m_Mass); }                   // Use class member variables - returns radius in Rsol
 };

diff --git a/src/BinaryConstituentStar.cpp b/src/BinaryConstituentStar.cpp
@@ -384,7 +384,7 @@ double BinaryConstituentStar::CalculateCircularisationTimescale(const double p_S
  */
 double BinaryConstituentStar::CalculateSynchronisationTimescale(const double p_SemiMajorAxis) {
 
-    double gyrationRadiusSquared   = CalculateMomentOfInertiaAU() / Mass() / Radius() / Radius();
+    double gyrationRadiusSquared   = CalculateMomentOfInertia() / Mass() / Radius() / Radius();
     double gyrationRadiusSquared_1 = 1.0 / gyrationRadiusSquared;
     double rOverA                  = Radius() / p_SemiMajorAxis;
     double rOverA_6                = rOverA * rOverA * rOverA * rOverA * rOverA * rOverA;

diff --git a/src/GiantBranch.cpp b/src/GiantBranch.cpp
@@ -692,7 +692,7 @@ double GiantBranch::CalculateRemnantRadius() const {
  *
  * @return                                      Radial extent of the star's convective envelope in Rsol
  */
-double GiantBranch::CalculateRadialExtentConvectiveEnvelope() const{
+double GiantBranch::CalculateRadialExtentConvectiveEnvelope() const {
 
 	BaseStar clone = *this;                         // clone this star so can manipulate without changes persisiting
 	clone.ResolveEnvelopeLoss(true);                // update clone's attributes after envelope is lost
@@ -1108,6 +1108,29 @@ double GiantBranch::CalculateLifetimeToHeIgnition(const double p_Mass, const dou
 }
 
 
+///////////////////////////////////////////////////////////////////////////////////////
+//                                                                                   //
+//                               ROTATION CALCULATIONS                               //
+//                                                                                   //
+///////////////////////////////////////////////////////////////////////////////////////
+
+/*
+ * Calculate moment of inertia
+ *
+ * Hurley et al., 2000, paragraph immediately following eq 109 
+ *
+ * 
+ * double GiantBranch::CalculateMomentOfInertia()
+ * 
+ * @return                                      Moment of inertia (Msol AU^2)
+ */
+double GiantBranch::CalculateMomentOfInertia() const {
+    double Rc = CalculateRemnantRadius();
+
+    return (0.1 * (m_Mass - m_CoreMass) * m_Radius * m_Radius) + (0.21 * m_CoreMass * Rc * Rc);
+}
+
+
 ///////////////////////////////////////////////////////////////////////////////////////
 //                                                                                   //
 //                              SUPERNOVA CALCULATIONS                               //
@@ -1673,12 +1696,6 @@ double GiantBranch::CalculateRemnantMassByBelczynski2002(const double p_Mass, co
 }
 
 
-
-
-
-
-
-
 /*
  * Driver function for Core Collapse Supernovas
  *

diff --git a/src/GiantBranch.h b/src/GiantBranch.h
@@ -87,8 +87,7 @@ class GiantBranch: virtual public BaseStar, public MainSequence {
             double          CalculateRemnantMassBySchneider2020(const double p_COCoreMass, const bool p_UseSchneiderAlt = false);
             double          CalculateRemnantMassBySchneider2020Alt(const double p_COCoreMass)               { return CalculateRemnantMassBySchneider2020(p_COCoreMass, true); }
 
-            double          CalculateMomentOfInertia() const                                                { double Rc = CalculateRemnantRadius(); return (0.1 * (m_Mass - m_CoreMass) * m_Radius * m_Radius) + (0.21 * m_CoreMass * Rc * Rc); } // k2 = 0.1 and k3 = 0.21 as defined in Hurley et al. 2000, after eq 109
-            double          CalculateMomentOfInertiaAU() const                                              { return CalculateMomentOfInertia() * RSOL_TO_AU * RSOL_TO_AU; }
+            double          CalculateMomentOfInertia() const;
 
             double          CalculatePerturbationMu() const;
 

diff --git a/src/MainSequence.h b/src/MainSequence.h
@@ -53,7 +53,6 @@ class MainSequence: virtual public BaseStar {
     double          CalculateLuminosityOnPhase() const                                      { return CalculateLuminosityOnPhase(m_Age, m_Mass0, m_LZAMS0); }        // Use class member variables
 
     double          CalculateMomentOfInertia() const                                        { return (0.1 * (m_Mass) * m_Radius * m_Radius); }                      // k2 = 0.1 as defined in Hurley et al. 2000, after eq 109
-    double          CalculateMomentOfInertiaAU() const                                      { return CalculateMomentOfInertia() * RSOL_TO_AU * RSOL_TO_AU; }
 
     double          CalculatePerturbationMu() const                                         { return 5.0; }                                                         // mu(MS) = 5.0 (Hurley et al. 2000, eqs 97 & 98)
 

diff --git a/src/NS.h b/src/NS.h
@@ -73,7 +73,6 @@ class NS: virtual public BaseStar, public Remnants {
 
     double          CalculateMomentOfInertiaCGS() const;                                                                                                            // MoI in CGS
     double          CalculateMomentOfInertia() const                                        { return CalculateMomentOfInertiaCGS() / MSOL_TO_G / RSOL_TO_CM / RSOL_TO_CM; } // MoI (default is solar units)
-    double          CalculateMomentOfInertiaAU() const                                      { return CalculateMomentOfInertia() * RSOL_TO_AU * RSOL_TO_AU; }
 
     double          CalculatePulsarBirthMagneticField();