From 3f33a5857183426fd5a38b89b9838e6f953eee8d Mon Sep 17 00:00:00 2001
From: Ilya Mandel <ilya.mandel@monash.edu>
Date: Thu, 5 Dec 2024 19:22:29 +1100
Subject: [PATCH] M*a should be constant on wind mass loss

---
 src/BaseBinaryStar.cpp | 5 ++++-
 1 file changed, 4 insertions(+), 1 deletion(-)

diff --git a/src/BaseBinaryStar.cpp b/src/BaseBinaryStar.cpp
index 624fa9f93..44aff43a7 100644
--- a/src/BaseBinaryStar.cpp
+++ b/src/BaseBinaryStar.cpp
@@ -1925,8 +1925,11 @@ void BaseBinaryStar::CalculateWindsMassLoss() {
             double mWinds1 = m_Star1->CalculateMassLossValues(true);                                                            // calculate new values assuming mass loss applied
             double mWinds2 = m_Star2->CalculateMassLossValues(true);                                                            // calculate new values assuming mass loss applied
 
-            double aWinds  = m_SemiMajorAxisPrev / (2.0 - ((m_Star1->MassPrev() + m_Star2->MassPrev()) / (mWinds1 + mWinds2))); // new semi-major axis for circularlised orbit
+            double aWinds  = m_SemiMajorAxisPrev * (m_Star1->Mass() + m_Star2->Mass()) / (mWinds1 + mWinds2);                   // new semi-major axis after wind mass loss, integrated to ensure a*M conservation
 
+            std::cout<<"Types"<<(int)m_Star1->StellarType()<<" "<<(int)m_Star2->StellarType()<<"prev masses"<<m_Star1->Mass()<<" "<<m_Star2->Mass()<<std::endl;
+            std::cout<<" new masses"<<mWinds1<<" "<<mWinds2<<"old a "<<m_SemiMajorAxisPrev<<"new a "<<aWinds<<std::endl;
+            
             m_Star1->SetMassLossDiff(mWinds1 - m_Star1->Mass());                                                                // JR: todo: find a better way?
             m_Star2->SetMassLossDiff(mWinds2 - m_Star2->Mass());                                                                // JR: todo: find a better way?