diff --git a/src/BaseBinaryStar.cpp b/src/BaseBinaryStar.cpp index 4ef7a6487..0025303b4 100644 --- a/src/BaseBinaryStar.cpp +++ b/src/BaseBinaryStar.cpp @@ -1201,6 +1201,10 @@ void BaseBinaryStar::ResolveSupernova() { m_Supernova->CalculateSNAnomalies(eccentricityPrev); double cosEccAnomaly = cos(m_Supernova->SN_EccentricAnomaly()); double sinEccAnomaly = sin(m_Supernova->SN_EccentricAnomaly()); + if ((utils::Compare(eccentricityPrev, 0.0) == 0) && m_Companion->IsOneOf(SN_REMNANTS)) { // If circular and first SN, fix eccentric anomaly to 0 + cosEccAnomaly = 1; + sinEccAnomaly = 0; + } // Derived quantities double aPrev = semiMajorAxisPrev_km; diff --git a/src/changelog.h b/src/changelog.h index a658ad9b8..3743c765a 100644 --- a/src/changelog.h +++ b/src/changelog.h @@ -1401,6 +1401,9 @@ // - Fix issue (likely introduced in 03.08.00) with the accretor not gaining mass appropriately // 03.09.01 RTW - Nov 27, 2024 - Enhancement: // - Added systemic velocity components x, y, and z to the output -const std::string VERSION_STRING = "03.09.01"; +// 03.09.02 RTW - Nov 27, 2024 - Defect repair, enhancement: +// - Fixed bugs in vector3d related to indexing and rotation +// - Added tweak for circular systems at first SN, to fix the x-axis along the separation vector +const std::string VERSION_STRING = "03.09.02"; # endif // __changelog_h__ diff --git a/src/typedefs.h b/src/typedefs.h index 23dc6e6cd..a60d38dbc 100755 --- a/src/typedefs.h +++ b/src/typedefs.h @@ -246,6 +246,13 @@ const STELLAR_TYPE_LIST COMPACT_OBJECTS = { STELLAR_TYPE::MASSLESS_REMNANT }; +// (convenience) initializer list for SN REMNANTS +const STELLAR_TYPE_LIST SN_REMNANTS = { + STELLAR_TYPE::NEUTRON_STAR, + STELLAR_TYPE::BLACK_HOLE, + STELLAR_TYPE::MASSLESS_REMNANT +}; + // (convenience) initializer list for GIANTS const STELLAR_TYPE_LIST GIANTS = { STELLAR_TYPE::FIRST_GIANT_BRANCH, diff --git a/src/vector3d.cpp b/src/vector3d.cpp index d7dadd57f..86539b2b5 100644 --- a/src/vector3d.cpp +++ b/src/vector3d.cpp @@ -81,8 +81,6 @@ Vector3d Vector3d::ChangeBasis(const double p_ThetaE, const double p_PhiE, const #define sPhi sin(p_PhiE) #define sPsi sin(p_PsiE) - Vector3d result = *this; // default return is this vector - // Define the Rotation Matrix std::vector rotationMatrix = { { cPhi * cPsi - sPhi * cTheta * sPsi , -cPhi * sPsi - sPhi * cTheta * cPsi , sTheta * sPhi }, @@ -90,14 +88,17 @@ Vector3d Vector3d::ChangeBasis(const double p_ThetaE, const double p_PhiE, const { sTheta * sPsi , sTheta * cPsi , cTheta } }; + Vector3d oldVector = *this; // input + Vector3d newVector = Vector3d(0,0,0); // output + // Apply rotation for (size_t row = 0; row < 3; row++) { for (size_t col = 0; col < 3; col++) { - result[row] += result[col] * rotationMatrix[row][col]; + newVector[row] += oldVector[col] * rotationMatrix[row][col]; } } - return result; + return newVector; #undef cTheta #undef cPhi diff --git a/src/vector3d.h b/src/vector3d.h index ed3d86d5d..d450c4834 100644 --- a/src/vector3d.h +++ b/src/vector3d.h @@ -41,8 +41,8 @@ class Vector3d { THROW_ERROR_IF(p_i < 0 || p_i > 2, ERROR::INDEX_OUT_OF_RANGE); // this is a code defect - if (p_i == 1) return m_x; - else if (p_i == 2) return m_y; + if (p_i == 0) return m_x; + else if (p_i == 1) return m_y; else return m_z; }