From 2d60d646216c550d6f1822f2393a9ffc6f1fb01e Mon Sep 17 00:00:00 2001 From: Tom Wagg Date: Thu, 19 Feb 2026 13:01:26 -0500 Subject: [PATCH 1/5] call hrdiag with a loop after CE, and actually check for/handle SNe immediately after CEs --- src/cosmic/src/evolv2.f | 157 +++++++++++++++++++++++++++++++++------- 1 file changed, 129 insertions(+), 28 deletions(-) diff --git a/src/cosmic/src/evolv2.f b/src/cosmic/src/evolv2.f index 585c126e..8a5e4fe9 100644 --- a/src/cosmic/src/evolv2.f +++ b/src/cosmic/src/evolv2.f @@ -2495,36 +2495,137 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, formation(j1) = formation(j2) endif +* Calculate stellar properties for both stars + do 117, k = 1,2 + m0 = mass0(k) + mt = mass(k) + mc = massc(k) + kw = kstar(k) + CALL star(kw,m0,mt,tm,tn,tscls,lums,GB,zpars) + CALL hrdiag(m0,aj(k),mt,tm,tn,tscls,lums,GB,zpars, + & rm,lum,kw,mc,rc,me,re,k2,bhspin(k),k) + +* +* Check for a supernova and correct the semi-major axis if so. +* + if((kw.ne.kstar(k).and.kstar(k).le.12.and. + & (kw.eq.13.or.kw.eq.14 + & .or.(kw.eq.15.and.pisn_track(k).eq.7))))then + dms(k) = mass(k) - mt + if(formation(k).ne.11) formation(k) = 1 + if(kw.eq.13.and.ecsn.gt.0.d0)then + if(kstar(k).le.6)then + if(mass0(k).le.zpars(5))then + if(sigma.gt.0.d0.and.sigmadiv.gt.0.d0)then + sigma = -sigmahold/sigmadiv + elseif(sigma.gt.0.d0.and.sigmadiv.lt.0.d0)then + sigma = sigmadiv + endif + formation(k) = 2 + endif + elseif(kstar(k).ge.7.and.kstar(k).le.9)then + if(mass(k).gt.ecsn_mlow.and.mass(k).le.ecsn)then + if(sigma.gt.0.d0.and.sigmadiv.gt.0.d0)then + sigma = -sigmahold/sigmadiv + elseif(sigma.gt.0.d0.and.sigmadiv.lt.0.d0)then + sigma = sigmadiv + endif + formation(k) = 2 + endif + elseif(formation(k).eq.11)then + if(sigma.gt.0.d0.and.sigmadiv.gt.0.d0)then + sigma = -sigmahold/sigmadiv + elseif(sigma.gt.0.d0.and.sigmadiv.lt.0.d0)then + sigma = sigmadiv + endif + formation(k) = 5 + elseif(kstar(k).ge.10.or.kstar(k).eq.12)then + if(sigma.gt.0.d0.and.sigmadiv.gt.0.d0)then + sigma = -sigmahold/sigmadiv + elseif(sigma.gt.0.d0.and.sigmadiv.lt.0.d0)then + sigma = sigmadiv + endif + formation(k) = 4 + endif + endif + + evolve_type = 14.d0 + FLOAT(k) + teff1 = 1000.d0*((1130.d0*lumin(1)/ + & (rad(1)**2.d0))**(1.d0/4.d0)) + teff2 = 1000.d0*((1130.d0*lumin(2)/ + & (rad(2)**2.d0))**(1.d0/4.d0)) + if(B_0(1).eq.0.d0)then !PK. + b01_bcm = 0.d0 + elseif(B_0(1).gt.0.d0.and.B(1).eq.0.d0)then + b01_bcm = B_0(1) + else + b01_bcm = B(1) + endif + if(B_0(2).eq.0.d0)then + b02_bcm = 0.d0 + elseif(B_0(2).gt.0.d0.and.B(2).eq.0.d0)then + b02_bcm = B_0(2) + else + b02_bcm = B(2) + endif + +* Check if PISN occurred, and if so overwrite formation + if(pisn_track(1).ne.0) formation(1) = pisn_track(1) + if(pisn_track(2).ne.0) formation(2) = pisn_track(2) + + CALL writetab(bpp_ind,tphys,evolve_type, + & mass(1),mass(2),kstar(1),kstar(2), + & sep,tb,ecc,rrl1,rrl2, + & aj(1),aj(2),tms(1),tms(2), + & mc_he(1),mc_he(2),mc_co(1),mc_co(2), + & rad(1),rad(2), + & mass0(1),mass0(2),lumin(1),lumin(2), + & teff1,teff2,radc(1),radc(2), + & menv(1),menv(2),renv(1),renv(2), + & ospin(1),ospin(2),b01_bcm,b02_bcm, + & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), + & epoch(2),bhspin(1),bhspin(2), + & deltam1_bcm,deltam2_bcm,formation(1), + & formation(2),binstate,mergertype,z,'bpp') + CALL kick(kw,mass(k),mc_co(k),mt,mass(3-k),ecc,sep,jorb,vk, + & k,rad(3-k),fallback,sigmahold,kick_info, + & disrupt,bkick) + sigma = sigmahold + + ! check if kick caused a merger + if(mass(3-k).lt.0.d0)then + if(kstar(3-k).lt.0.d0) mt = mt-mass(3-k) + if(kw.eq.13.and.mt.gt.mxns) kw = 14 + CALL CONCATKSTARS(kstar(j1), kstar(j2), mergertype) + kstar(k) = kw + mass(k) = mt + epoch(k) = tphys - age + kstar(3-k) = 15 + mass(3-k) = 0.d0 + coel = .true. + binstate = 1 + endif + + tb = (sep/aursun)*SQRT(sep/(aursun*(mt+mass(3-k)))) + oorb = twopi/tb + endif + + mass0(k) = m0 + mass(k) = mt + massc(k) = mc + kstar(k) = kw + rad(k) = rm + lumin(k) = lum + radc(k) = rc + menv(k) = me + renv(k) = re + + +117 continue + * evolve_type = 8.0 - - - mc = massc(1) - rc = radc(1) - CALL star(kstar(1),mass0(1),mass(1),tm,tn,tscls,lums,GB,zpars) - CALL hrdiag(mass0(1),aj(1),mass(1),tm,tn,tscls,lums,GB,zpars, - & rm,lum,kstar(1),mc,rc,me,re,k2,bhspin(1),1) - - rad(1) = rm - lumin(1) = lum - massc(1) = mc - radc(1) = rc - menv(1) = me - renv(1) = re - - - mc = massc(2) - rc = radc(2) - CALL star(kstar(2),mass0(2),mass(2),tm,tn,tscls,lums,GB,zpars) - CALL hrdiag(mass0(2),aj(2),mass(2),tm,tn,tscls,lums,GB,zpars, - & rm,lum,kstar(2),mc,rc,me,re,k2,bhspin(2),2) - - rad(2) = rm - lumin(2) = lum - massc(2) = mc - radc(2) = rc - menv(2) = me - renv(2) = re + mass1_bpp = mass(1) mass2_bpp = mass(2) From 7d05bc062d40ebdbe8572f7aa76fb4e08eb5482d Mon Sep 17 00:00:00 2001 From: Tom Wagg Date: Thu, 19 Feb 2026 13:01:41 -0500 Subject: [PATCH 2/5] update version and changelog --- changelog.md | 7 ++++++- src/cosmic/_version.py | 2 +- 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/changelog.md b/changelog.md index dcf95f2f..2c1d1fda 100644 --- a/changelog.md +++ b/changelog.md @@ -77,7 +77,11 @@ See the discussed changes in our previous releases here: https://github.com/COSM - Bug fix [[#725](https://github.com/COSMIC-PopSynth/COSMIC/issues/725)]: set ``tb=sep=0`` for cases where AIC caused a WD to explode and leave behind no remnant (``evolve_type==9``) - Bug fix [[#724](https://github.com/COSMIC-PopSynth/COSMIC/issues/724)]: remove bug where ``kstar=15`` was assigned its previous epoch mass after merging during a CE and the merger product goes SN -## 3.7.2 +## 3.7.2-3.7.3 + +Issues with PyPi, unreleased. + +## 3.7.4 This release contains _several_ fixes to how CO core masses/remnant masses are handled. It also adds a new PISN prescription, windflag and LBV winds flag. - Fixes: @@ -88,6 +92,7 @@ This release contains _several_ fixes to how CO core masses/remnant masses are h - Update ``mc = mt`` for **ALL** remnants instead of just for certain flags - Only apply conversion from baryonic to gravitational mass once. Previous code **APPLIED IT TWICE!!** for core-collapse BHs! Create subroutine ``baryonic_to_gravitational_mass``, which is called throughout - Ensure stellar type changes during hrdiag are logged to the bpp before SN occurs + - Add check after CE for whether a SN occurs right at the end - Additions/changes: - Linearly grow the core mass of HeMS stars from 0 to the value at the start of HeHG and check for core-collapse for all He stars (not just HeHG and HeGB) diff --git a/src/cosmic/_version.py b/src/cosmic/_version.py index 3bdbbae2..035c71bf 100644 --- a/src/cosmic/_version.py +++ b/src/cosmic/_version.py @@ -1 +1 @@ -__version__ = "3.7.2" +__version__ = "3.7.4" From 81247ed1bc6dfc98d0ce3351755e03ccaaeb15c1 Mon Sep 17 00:00:00 2001 From: Tom Wagg Date: Thu, 19 Feb 2026 13:17:23 -0500 Subject: [PATCH 3/5] fix orbital period etc when disruption occurs --- src/cosmic/src/evolv2.f | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/src/cosmic/src/evolv2.f b/src/cosmic/src/evolv2.f index 8a5e4fe9..e26248b9 100644 --- a/src/cosmic/src/evolv2.f +++ b/src/cosmic/src/evolv2.f @@ -2606,8 +2606,14 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, binstate = 1 endif - tb = (sep/aursun)*SQRT(sep/(aursun*(mt+mass(3-k)))) - oorb = twopi/tb + if(ecc.gt.1.d0)then + tb = -1d0 + sep = -1d0 + oorb = -1d0 + else + tb = (sep/aursun)*SQRT(sep/(aursun*(mt+mass(3-k)))) + oorb = twopi/tb + endif endif mass0(k) = m0 From 3d4213b3570bba879f5679767033fb903d8c5f28 Mon Sep 17 00:00:00 2001 From: Tom Wagg Date: Thu, 19 Feb 2026 13:54:54 -0500 Subject: [PATCH 4/5] clean up issue templates --- .github/ISSUE_TEMPLATE/feature-request.md | 17 ++++++++ .github/ISSUE_TEMPLATE/general_issue.md | 9 ---- .github/ISSUE_TEMPLATE/new-prescription.md | 17 ++++++++ .github/ISSUE_TEMPLATE/question.md | 10 +++++ ...nexpected-stellar-evolutionary-behavior.md | 41 ++++++++++++++++--- 5 files changed, 79 insertions(+), 15 deletions(-) create mode 100644 .github/ISSUE_TEMPLATE/feature-request.md delete mode 100644 .github/ISSUE_TEMPLATE/general_issue.md create mode 100644 .github/ISSUE_TEMPLATE/new-prescription.md create mode 100644 .github/ISSUE_TEMPLATE/question.md diff --git a/.github/ISSUE_TEMPLATE/feature-request.md b/.github/ISSUE_TEMPLATE/feature-request.md new file mode 100644 index 00000000..3b208c38 --- /dev/null +++ b/.github/ISSUE_TEMPLATE/feature-request.md @@ -0,0 +1,17 @@ +--- +name: Feature request +about: I want to request a new feature +title: '' +labels: feature-request +assignees: '' + +--- + +**Describe the new feature you'd like** +A clear and concise description of what you want to happen. E.g. "I love how this package can evolve populations, but I think a great improvement would be for it to also do my taxes." + +**Any suggestions for how to implement this** +Is there a particular way you would implement this? How could this integrate into the existing code? + +**Is your feature request related to a problem? Please describe.** +A clear and concise description of what the problem is. Ex. I'm always frustrated when [...] diff --git a/.github/ISSUE_TEMPLATE/general_issue.md b/.github/ISSUE_TEMPLATE/general_issue.md deleted file mode 100644 index e0973929..00000000 --- a/.github/ISSUE_TEMPLATE/general_issue.md +++ /dev/null @@ -1,9 +0,0 @@ ---- -name: general issue -about: Way to systematically track weird evoltuionary behavior and things tried to - resolve that behavior for better history tracking. -title: "[general issue]" -labels: general ---- - -# Describe your issue in as much detail as possible -- if your issue arises from a cosmic-pop command line call, please include the exact syntax of the command line call in your description. diff --git a/.github/ISSUE_TEMPLATE/new-prescription.md b/.github/ISSUE_TEMPLATE/new-prescription.md new file mode 100644 index 00000000..65f9217d --- /dev/null +++ b/.github/ISSUE_TEMPLATE/new-prescription.md @@ -0,0 +1,17 @@ +--- +name: New physics prescription +about: I want to request/suggest a new physics prescription +title: '' +labels: new-prescription +assignees: '' + +--- + +**Describe the new physics prescription you want to add** +TODO + +**Link to the publication this is based on** +TODO + +**Is this a new option for an existing flag? Or would require a new setting?** +TODO \ No newline at end of file diff --git a/.github/ISSUE_TEMPLATE/question.md b/.github/ISSUE_TEMPLATE/question.md new file mode 100644 index 00000000..02d06d39 --- /dev/null +++ b/.github/ISSUE_TEMPLATE/question.md @@ -0,0 +1,10 @@ +--- +name: Question +about: I want to ask a question about COSMIC +title: '' +labels: question +assignees: '' + +--- + +Use this issue to ask any questions about how functions work, how best to use the package for your project or anything else you need help with! Wherever possible please include code examples to make it easier to help you! diff --git a/.github/ISSUE_TEMPLATE/unexpected-stellar-evolutionary-behavior.md b/.github/ISSUE_TEMPLATE/unexpected-stellar-evolutionary-behavior.md index 55b6bc6c..7187059b 100644 --- a/.github/ISSUE_TEMPLATE/unexpected-stellar-evolutionary-behavior.md +++ b/.github/ISSUE_TEMPLATE/unexpected-stellar-evolutionary-behavior.md @@ -1,19 +1,48 @@ --- name: Unexpected Stellar Evolutionary Behavior -about: Way to systematically track weird evoltuionary behavior and things tried to - resolve that behavior for better history tracking. -title: "[Unexpected Evolutionary Behavior]" +about: I want to report some weird evolutionary behavior and things tried to resolve that behavior for better history tracking. +title: "" labels: Unexpected Evolutionary Behavior assignees: '' --- -# CSV file with Initial Conditions and BSE flags creating this behavior +# Checklist +- [ ] Attached CSV file with initial conditions and BSE flags for reproducing this behavior +- [ ] Included code that produces the BPP/BCM arrays +- [ ] Written description of weird behavior +- [ ] Explained expected outcome +- [ ] Listed flags you have tried to turn on and off that you think relate to the evolutionary behavior +- [ ] Attached a plot if applicable + +# CSV file + +Attach file with initC file here. You can save this like +```python +# create CSV file +initC.to_csv('weird_behavior.csv', index=False) + +# now attach it to this issue! +``` # Code that produces the BPP and BCM arrays -# Description of Weird Behavior +```python +# your code here +``` + +# Description + +Add a description of the weird behavior you've seen -# Flags you have tried to turn on and off that you think relate to the evolutionary behavior +# Expect outcome + +How do you think this system _should_ evolve? + +# Flags + +Which flag have you tried to turn on and off that you think relate to the evolutionary behavior? # Plot (if applicable) + +Plot goes here \ No newline at end of file From c8e47021652621589c31ef5e63dd84761516234e Mon Sep 17 00:00:00 2001 From: Tom Wagg Date: Thu, 19 Feb 2026 16:14:00 -0500 Subject: [PATCH 5/5] don't apply mass modifiers to timesteps provided in initC --- src/cosmic/evolve.py | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/cosmic/evolve.py b/src/cosmic/evolve.py index a27a61fb..173aab7e 100644 --- a/src/cosmic/evolve.py +++ b/src/cosmic/evolve.py @@ -184,6 +184,8 @@ def evolve(self, initialbinarytable, pool=None, bpp_columns=None, bcm_columns=No are logged in the initial conditions table so you can keep track of which systems had their timesteps modified. Avoid overlapping mass ranges for different modifiers as this will result in multiple modifiers being applied in the overlap region. + NOTE: these modifiers are only applied to columns which aren't present in the initialbinarytable + that is passed in (i.e. they only modify values provided by a BSEDict or params.ini file) **kwargs: There are three ways to tell evolve and thus the fortran @@ -250,6 +252,8 @@ def evolve(self, initialbinarytable, pool=None, bpp_columns=None, bcm_columns=No if bcm_columns is None: bcm_columns = BCM_COLUMNS + columns_in_passed_initC = set(initialbinarytable.columns) + # There are three ways to tell evolve and thus the fortran # what you want all the flags and other BSE specific # parameters to be @@ -383,7 +387,8 @@ def evolve(self, initialbinarytable, pool=None, bpp_columns=None, bcm_columns=No initialbinarytable = initialbinarytable.assign(fprimc_array=initialbinarytable[FPRIMC_COLUMNS].values.tolist()) # update timesteps based on mass modifier - if dt_mass_modifiers: + mass_modifier_cols = set(['pts1', 'pts2', 'pts3']).difference(columns_in_passed_initC) + if dt_mass_modifiers and len(mass_modifier_cols) != 0: # warn the user if their mass ranges overlap for i in range(len(dt_mass_modifiers)): for j in range(i + 1, len(dt_mass_modifiers)): @@ -405,9 +410,8 @@ def evolve(self, initialbinarytable, pool=None, bpp_columns=None, bcm_columns=No raise ValueError(f"Timestep modifiers must be positive. You passed {mod} for the " f"mass range {m_low} to {m_high}.") mask = (initialbinarytable['mass_1'] >= m_low) & (initialbinarytable['mass_1'] < m_high) - initialbinarytable.loc[mask, 'pts1'] *= mod - initialbinarytable.loc[mask, 'pts2'] *= mod - initialbinarytable.loc[mask, 'pts3'] *= mod + for col in mass_modifier_cols: + initialbinarytable.loc[mask, col] *= mod # need to ensure that the order of parameters that we pass to BSE # is correct