From 049efbd973577b0abbc863db8c322f9957cd4402 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 3 Apr 2024 07:51:17 -0400 Subject: [PATCH 01/18] Added compression options for pandas --- bin/cosmic-pop | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/bin/cosmic-pop b/bin/cosmic-pop index 4176c40cf..e131448c7 100755 --- a/bin/cosmic-pop +++ b/bin/cosmic-pop @@ -119,6 +119,10 @@ def parse_commandline(): parser.add_argument("--seed", type=int) parser.add_argument("--verbose", action="store_true", default=False, help="Run in Verbose Mode") + parser.add_argument("--complib",type=str,default="zlib", + help="HDFStore compression library") + parser.add_argument("--complevel",type=int,default=0, + help="HDFStore compression level") group = parser.add_mutually_exclusive_group() group.add_argument("-n", "--nproc", @@ -232,7 +236,7 @@ if __name__ == '__main__': # Open the hdf5 file to store the fixed population data try: - dat_store = pd.HDFStore('dat_kstar1_{0}_kstar2_{1}_SFstart_{2}_SFduration_{3}_metallicity_{4}.h5'.format(kstar1_range_string, kstar2_range_string, sampling['SF_start'], sampling['SF_duration'], sampling['metallicity'])) + dat_store = pd.HDFStore('dat_kstar1_{0}_kstar2_{1}_SFstart_{2}_SFduration_{3}_metallicity_{4}.h5'.format(kstar1_range_string, kstar2_range_string, sampling['SF_start'], sampling['SF_duration'], sampling['metallicity']),complib=args.complib,complevel=args.complevel) conv_save = pd.read_hdf(dat_store, 'conv') log_file = open('log_kstar1_{0}_kstar2_{1}_SFstart_{2}_SFduration_{3}_metallicity_{4}.txt'.format(kstar1_range_string, kstar2_range_string, sampling['SF_start'], sampling['SF_duration'], sampling['metallicity']), 'a') log_file.write('There are already: '+str(conv_save.shape[0])+' '+kstar1_range_string+'_'+kstar2_range_string+' binaries evolved\n') @@ -246,7 +250,7 @@ if __name__ == '__main__': idx = int(np.max(pd.read_hdf(dat_store, 'idx'))[0]) except: conv_save = pd.DataFrame() - dat_store = pd.HDFStore('dat_kstar1_{0}_kstar2_{1}_SFstart_{2}_SFduration_{3}_metallicity_{4}.h5'.format(kstar1_range_string, kstar2_range_string, sampling['SF_start'], sampling['SF_duration'], sampling['metallicity'])) + dat_store = pd.HDFStore('dat_kstar1_{0}_kstar2_{1}_SFstart_{2}_SFduration_{3}_metallicity_{4}.h5'.format(kstar1_range_string, kstar2_range_string, sampling['SF_start'], sampling['SF_duration'], sampling['metallicity']),complib=args.complib,complevel=args.complevel) total_mass_singles = 0 total_mass_binaries = 0 total_mass_stars = 0 From 17e2fdc7dc8db51fdca763495fc2d298f898a5cb Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 3 Apr 2024 17:42:59 -0400 Subject: [PATCH 02/18] Created maximum wall time argument --- bin/cosmic-pop | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/bin/cosmic-pop b/bin/cosmic-pop index e131448c7..af3149cf8 100755 --- a/bin/cosmic-pop +++ b/bin/cosmic-pop @@ -98,6 +98,8 @@ def parse_commandline(): help="Number of binaries to try before checking for " "convergence, it will check ever Nstep binaries until " "it reach Niter binaries", type=int, default=10000) + parser.add_argument("--max-wall-time", type=int, default=3155760, + help="Maximum wall time (seconds) for sampling binaries") parser.add_argument("--binary_state", nargs='+', type=int) parser.add_argument("--sampling_method") parser.add_argument("--primary_model", help="Chooses the initial primary mass function from: salpeter55, kroupa93, kroupa01", type=str) @@ -288,7 +290,7 @@ if __name__ == '__main__': log_file.write("You have specified both qmin and m2_min.\n") log_file.write("COSMIC will use qmin={} to determine the secondary masses in the initial sample.\n".format(args.qmin)) - while (Nstep < args.Niter) & (np.max(match) > convergence['match']): + while (Nstep < args.Niter) & (np.max(match) > convergence['match']) & ((time.time() - start_time) < args.max_wall_time): # Set random seed such that each iteration gets a unique, determinable seed rand_seed = seed_int + Nstep np.random.seed(rand_seed) From 498fbd31db81743e8d42ea2c3f5133d852c486b8 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Fri, 15 Nov 2024 11:55:10 -0500 Subject: [PATCH 03/18] added new physics --- src/cosmic/src/assign_remnant.f | 56 +++++++++++++++++++++++++++++ src/cosmic/src/mlwind.f | 64 +++++++++++++++++++++++---------- src/cosmic/utils.py | 4 +-- 3 files changed, 104 insertions(+), 20 deletions(-) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 858cde738..50ec530ca 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -13,6 +13,7 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) real*8 mc,mcbagb,mass,mt real*8 frac,kappa,sappa,alphap,polyfit real*8 mcx, bhspin,mrem,mch + real*8 fmix, mcritnsbh, mtemp1, mtemp2 integer kw,kidx @@ -231,6 +232,61 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) endif endif mc = mt + elseif(remnantflag.eq.5)then +* +* Use the Fryer et al. 2022 SN Prescription +* +* For this, we just set the proto-core mass to one + if(mc.le.3.5d0)then + mcx = 1.2d0 + elseif(mc.le.6.d0)then + mcx = 1.3d0 + elseif(mc.le.11.d0)then + mcx = 1.4d0 + elseif(mc.gt.11.d0)then + mcx = 1.6d0 + endif + + if(ecsn.gt.0.d0.and.mcbagb.le.ecsn.and. + & mcbagb.ge.ecsn_mlow)then + mt = 1.38d0 ! ECSN fixed mass, no fallback + else +* Parameters of Fryer2022 model + fmix=1.0 + mcritnsbh=5.75 +* We need mt in multiple places, so temp1 will be the working mt + mtemp1=mt +* mtemp2 is the calculated value of the remnant mass + mtemp2=1.2 + (0.05*fmix) + + & (0.01*((mc/fmix)**2)) + + & EXP(fmix*(mc-mcritnsbh)) +* We don't care about mtemp2 if it's less than zero + if(mtemp2.lt.0.)then + mtemp1 = 0. + kw=15 +* We only care about mtemp2 if it is less than the total +* mass of the star + elseif(mtemp2.lt.mt)then + mtemp1 = mtemp2 +* If mtemp2 is less, we also want to estimate the fallback fraction + fallback=(mtemp1-mcx)/(mt-mcx) + mt = mcx + fallback*(mtemp1 - mcx) + endif + endif + if(bhspinflag.eq.0)then + bhspin = bhspinmag + elseif(bhspinflag.eq.1)then + bhspin = ran3(idum1) * bhspinmag + elseif(bhspinflag.eq.2)then + if(mc.le.13.d0)then + bhspin = 0.9d0 + elseif(mc.lt.27.d0)then + bhspin = -0.064d0*mc + 1.736d0 + else + bhspin = 0.0d0 + endif + endif + mc = mt endif * Specify the baryonic to gravitational remnant mass prescription diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index 2f514747e..64d4588a3 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -123,7 +123,8 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) endif * mlwind = dms - elseif(windflag.eq.2.or.windflag.eq.3)then + elseif(windflag.eq.2.or.windflag.eq.3.or.windflag.eq.5 + & .or.windflag.eq.6.or.windflag.eq.7)then * Vink winds etc according to as implemented following * Belczynski, Bulik, Fryer, Ruiter, Valsecchi, Vink & Hurley 2010. * @@ -167,24 +168,51 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) endif * Apply Vink, de Koter & Lamers (2001) OB star winds. * Next check if hot massive H-rich O/B star in appropriate temperature ranges. - if(teff.ge.12500.and.teff.le.25000)then - if(eddlimflag.eq.0) alpha = 0.85d0 - if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) - dms = -6.688d0 + 2.210d0*LOG10(lum/1.0d+05) - - & 1.339d0*LOG10(mt/30.d0) - 1.601d0*LOG10(1.3d0/2.d0) + - & alpha*LOG10(z/zsun) + 1.07d0*LOG10(teff/2.0d+04) - dms = 10.d0**dms - testflag = 2 - elseif(teff.gt.25000.)then -* Although Vink et al. formulae are only defined until Teff=50000K, -* we follow the Dutch prescription of MESA, and extend to higher Teff - dms = -6.697d0 + 2.194d0*LOG10(lum/1.0d+05) - - & 1.313d0*LOG10(mt/30.d0) - 1.226d0*LOG10(2.6d0/2.d0) + - & alpha*LOG10(z/zsun) +0.933d0*LOG10(teff/4.0d+04) - - & 10.92d0*(LOG10(teff/4.0d+04)**2) - dms = 10.d0**dms - testflag = 2 + if (windflag.eq.2.or.windflag.eq.3.or.windflag.eq.5) then + if(teff.ge.12500.and.teff.le.25000)then + if(eddlimflag.eq.0) alpha = 0.85d0 + if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) + dms = -6.688d0 + 2.210d0*LOG10(lum/1.0d+05) - + & 1.339d0*LOG10(mt/30.d0) - + & 1.601d0*LOG10(1.3d0/2.d0) + + & alpha*LOG10(z/zsun) + 1.07d0*LOG10(teff/2.0d+04) + dms = 10.d0**dms + testflag = 2 + elseif(teff.gt.25000.)then +* Although Vink et al. formulae are only defined until Teff=50000K, +* we follow the Dutch prescription of MESA, and extend to higher Teff + dms = -6.697d0 + 2.194d0*LOG10(lum/1.0d+05) - + & 1.313d0*LOG10(mt/30.d0) - + & 1.226d0*LOG10(2.6d0/2.d0) + + & alpha*LOG10(z/zsun) +0.933d0*LOG10(teff/4.0d+04) - + & 10.92d0*(LOG10(teff/4.0d+04)**2) + dms = 10.d0**dms + testflag = 2 + endif +* Apply fwind cut + if (windflag.eq.5) then +* TODO This should not be 3.0. It should be specified in the .ini file +* It should also be inverted (0.333 rather than 3.0) + dms = dms*0.33 + endif + endif +* Bjorklund 2023 https://doi.org/10.1051/0004-6361/202141948 + if (windflag.eq.6) then + dms = -5.52d0 + 2.39d0*LOG10(lum/1.0d+06) - + & 1.48*LOG10(mt/45.0d0) + 2.12d0 *LOG10(teff/4.5d+04) - + & (0.75d+0 - 1.87d0*LOG10(teff/4.5d+04))*LOG10(z/zsun) + dms = 10.0d0**dms + endif +* Kticka, Kubat, Kritckova 2023 10.48550/arXiv.2311.01257 + if (windflag.eq.7) then + dms = -13.82d0 + 0.358d0*LOG10(z/zsun) + + & (1.52d0 - 0.11d0*LOG10(z/zsun))*LOG10(lum/1.0d+06) + + & 13.82d0*LOG10((1.0d0 + 0.73*LOG10(z/zsun))*EXP( + & -1.0d0 * (teff - 1.416d+4)**2 / (3.58d+3)**2) + + & 3.84d0*EXP(-1.0d0*(teff-3.790d+04)**2 / (5.65d+04)**2)) + dms = 10.0d0**dms endif +* if((windflag.eq.3.or.kw.ge.2).and.kw.le.6)then * LBV-like mass loss beyond the Humphreys-Davidson limit. diff --git a/src/cosmic/utils.py b/src/cosmic/utils.py index 6fa985774..b184e94a6 100644 --- a/src/cosmic/utils.py +++ b/src/cosmic/utils.py @@ -1139,7 +1139,7 @@ def error_check(BSEDict, filters=None, convergence=None, sampling=None): flag = "windflag" if flag in BSEDict.keys(): - if BSEDict[flag] not in [0, 1, 2, 3, 4]: + if BSEDict[flag] not in [0, 1, 2, 3, 4, 5, 6, 7]: raise ValueError( "'{0:s}' needs to be set to either 0, 1, 2, or 3, 4 (you set it to '{1:d}')".format( flag, BSEDict[flag] @@ -1393,7 +1393,7 @@ def error_check(BSEDict, filters=None, convergence=None, sampling=None): flag = "remnantflag" if flag in BSEDict.keys(): - if BSEDict[flag] not in [0, 1, 2, 3, 4]: + if BSEDict[flag] not in [0, 1, 2, 3, 4, 5]: raise ValueError( "'{0:s}' needs to be set to either 0, 1, 2, 3, or 4 (you set it to '{1:d}')".format( flag, BSEDict[flag] From 041a0985ef7a31a532b95e7638c490957f08a9c4 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Tue, 3 Dec 2024 20:00:50 -0500 Subject: [PATCH 04/18] updated solar wind assumptions different from Vink --- src/cosmic/src/mlwind.f | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index 64d4588a3..d974295f0 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -200,14 +200,14 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) if (windflag.eq.6) then dms = -5.52d0 + 2.39d0*LOG10(lum/1.0d+06) - & 1.48*LOG10(mt/45.0d0) + 2.12d0 *LOG10(teff/4.5d+04) - - & (0.75d+0 - 1.87d0*LOG10(teff/4.5d+04))*LOG10(z/zsun) + & (0.75d+0 - 1.87d0*LOG10(teff/4.5d+04))*LOG10(z/0.014d0) dms = 10.0d0**dms endif * Kticka, Kubat, Kritckova 2023 10.48550/arXiv.2311.01257 if (windflag.eq.7) then - dms = -13.82d0 + 0.358d0*LOG10(z/zsun) + - & (1.52d0 - 0.11d0*LOG10(z/zsun))*LOG10(lum/1.0d+06) + - & 13.82d0*LOG10((1.0d0 + 0.73*LOG10(z/zsun))*EXP( + dms = -13.82d0 + 0.358d0*LOG10(z/0.0134d0) + + & (1.52d0 - 0.11d0*LOG10(z/0.0134d0))*LOG10(lum/1.0d+06) + + & 13.82d0*LOG10((1.0d0 + 0.73*LOG10(z/0.0134d0))*EXP( & -1.0d0 * (teff - 1.416d+4)**2 / (3.58d+3)**2) + & 3.84d0*EXP(-1.0d0*(teff-3.790d+04)**2 / (5.65d+04)**2)) dms = 10.0d0**dms From 3a38d08d6ae9750fb99ceaba58a3f70a7d31217c Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Sat, 25 Jan 2025 22:51:54 -0500 Subject: [PATCH 05/18] fixed wind --- src/cosmic/src/mlwind.f | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index 64d4588a3..fdb4fa44f 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -198,12 +198,12 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) endif * Bjorklund 2023 https://doi.org/10.1051/0004-6361/202141948 if (windflag.eq.6) then - dms = -5.52d0 + 2.39d0*LOG10(lum/1.0d+06) - + dms = -5.52d0 + 2.39d0*LOG10(lum/1.0d+06) + & 1.48*LOG10(mt/45.0d0) + 2.12d0 *LOG10(teff/4.5d+04) - & (0.75d+0 - 1.87d0*LOG10(teff/4.5d+04))*LOG10(z/zsun) dms = 10.0d0**dms endif -* Kticka, Kubat, Kritckova 2023 10.48550/arXiv.2311.01257 +* Kritcka, Kubat, Kritckova 2023 10.48550/arXiv.2311.01257 if (windflag.eq.7) then dms = -13.82d0 + 0.358d0*LOG10(z/zsun) + & (1.52d0 - 0.11d0*LOG10(z/zsun))*LOG10(lum/1.0d+06) + From 0bdcdeabe91e1409c42346221e96b8afff4691e2 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Sat, 25 Jan 2025 22:55:07 -0500 Subject: [PATCH 06/18] changed wind to be fixed --- src/cosmic/src/mlwind.f | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index b9c07dd7a..25df1aa68 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -199,7 +199,7 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) * Bjorklund 2023 https://doi.org/10.1051/0004-6361/202141948 if (windflag.eq.6) then dms = -5.52d0 + 2.39d0*LOG10(lum/1.0d+06) + - & 1.48*LOG10(mt/45.0d0) + 2.12d0 *LOG10(teff/4.5d+04) - + & 1.48*LOG10(mt/45.0d0) + 2.12d0 *LOG10(teff/4.5d+04) + & (0.75d+0 - 1.87d0*LOG10(teff/4.5d+04))*LOG10(z/0.014d0) dms = 10.0d0**dms endif From 7dc6074dff1205704230f9d7f47e39b8db89201a Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Thu, 6 Mar 2025 11:11:35 -0500 Subject: [PATCH 07/18] added new WR winds --- src/cosmic/src/mlwind.f | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index 25df1aa68..5589bf8b2 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -170,7 +170,7 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) * Next check if hot massive H-rich O/B star in appropriate temperature ranges. if (windflag.eq.2.or.windflag.eq.3.or.windflag.eq.5) then if(teff.ge.12500.and.teff.le.25000)then - if(eddlimflag.eq.0) alpha = 0.85d0 + if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.85d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = -6.688d0 + 2.210d0*LOG10(lum/1.0d+05) - & 1.339d0*LOG10(mt/30.d0) - @@ -220,7 +220,7 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) * past the limit, rather than just for giant, evolved stars x = 1.0d-5*r*sqrt(lum) if(lum.gt.6.0d+05.and.x.gt.1.d0)then - if(eddlimflag.eq.0) alpha = 0.d0 + if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = 1.5d0*1.0d-04*((z/zsun)**alpha) testflag = 3 @@ -228,9 +228,15 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) elseif(kw.ge.7.and.kw.le.9)then !WR (naked helium stars) * If naked helium use Hamann & Koesterke (1998) WR winds reduced by factor of * 10 (Yoon & Langer 2005), with Vink & de Koter (2005) metallicity dependence - if(eddlimflag.eq.0) alpha = 0.86d0 + if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.86d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = 1.0d-13*(lum**1.5d0)*((z/zsun)**alpha) +* Yang et al (2023) 10.1051/0004-6361/202244770 + if(eddlimflag.eq.2) then + alpha = LOG10(lum) + dms = 10**(0.45d0*alpha**3 - 5.26d0*alpha**2 + 20.93d0*alpha - + & 34.56d0) + endif testflag = 4 endif * From 1b67e1c20e2ae19d327cee75c05930943e1ce656 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Thu, 6 Mar 2025 11:16:13 -0500 Subject: [PATCH 08/18] undid bad decisions --- src/cosmic/src/mlwind.f | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/src/cosmic/src/mlwind.f b/src/cosmic/src/mlwind.f index 5589bf8b2..25df1aa68 100644 --- a/src/cosmic/src/mlwind.f +++ b/src/cosmic/src/mlwind.f @@ -170,7 +170,7 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) * Next check if hot massive H-rich O/B star in appropriate temperature ranges. if (windflag.eq.2.or.windflag.eq.3.or.windflag.eq.5) then if(teff.ge.12500.and.teff.le.25000)then - if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.85d0 + if(eddlimflag.eq.0) alpha = 0.85d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = -6.688d0 + 2.210d0*LOG10(lum/1.0d+05) - & 1.339d0*LOG10(mt/30.d0) - @@ -220,7 +220,7 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) * past the limit, rather than just for giant, evolved stars x = 1.0d-5*r*sqrt(lum) if(lum.gt.6.0d+05.and.x.gt.1.d0)then - if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.d0 + if(eddlimflag.eq.0) alpha = 0.d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = 1.5d0*1.0d-04*((z/zsun)**alpha) testflag = 3 @@ -228,15 +228,9 @@ real*8 FUNCTION mlwind(kw,lum,r,mt,mc,rl,z) elseif(kw.ge.7.and.kw.le.9)then !WR (naked helium stars) * If naked helium use Hamann & Koesterke (1998) WR winds reduced by factor of * 10 (Yoon & Langer 2005), with Vink & de Koter (2005) metallicity dependence - if(eddlimflag.eq.0.or.eddlimflag.eq.2) alpha = 0.86d0 + if(eddlimflag.eq.0) alpha = 0.86d0 if(eddlimflag.eq.1) alpha = MLalpha(mt,lum,kw) dms = 1.0d-13*(lum**1.5d0)*((z/zsun)**alpha) -* Yang et al (2023) 10.1051/0004-6361/202244770 - if(eddlimflag.eq.2) then - alpha = LOG10(lum) - dms = 10**(0.45d0*alpha**3 - 5.26d0*alpha**2 + 20.93d0*alpha - - & 34.56d0) - endif testflag = 4 endif * From e56d53f91fe46786ed65599a07cc4672f6c10e0f Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Thu, 6 Mar 2025 12:46:11 -0500 Subject: [PATCH 09/18] added Chris's ppisn --- src/cosmic/src/assign_remnant.f | 28 ++++++++++++++++++++++++++++ 1 file changed, 28 insertions(+) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 50ec530ca..7056595b5 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -84,6 +84,34 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) * kw = 15 else +* Chris Belczynski Evolutionary Roads Weak PPISN +* This has to happen before the SNa, because it modifies +* the properties of the star during explosion + if(pisn.eq.-4.and.mt.ge.45d0)then + if(mcbagb.ge.65d0) then + mt = 0.d0 + kw = 15 + else +* PPISN + if(mcbagb.ge.60d0) then + mtemp1 = 938d0 - (14.3d0*mcbagb) + elseif(mcbagb.ge.40d0) then + mtemp1 = 55.6d0 + else + mtemp1 = 6.0d0 + (0.83d0*mcbagb) + endif +* Update mass + if(mt.gt.mtemp1) then + mt = mtemp1 + endif + if(mcbagb.gt.mtemp1) then + mcbagb = mtemp1 + endif + if(mc.gt.mtemp1) then + mc = mtemp1 + endif + endif +* Carry on with the Supernovae if(remnantflag.eq.0)then mt = 1.17d0 + 0.09d0*mc elseif(remnantflag.eq.1)then From dcaf7c331485900cc3690844de7d148c3f11d0c7 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Mon, 10 Mar 2025 10:13:40 -0400 Subject: [PATCH 10/18] added weak pisn --- src/cosmic/src/assign_remnant.f | 1 + src/cosmic/utils.py | 3 ++- 2 files changed, 3 insertions(+), 1 deletion(-) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 7056595b5..745103655 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -109,6 +109,7 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) endif if(mc.gt.mtemp1) then mc = mtemp1 + endif endif endif * Carry on with the Supernovae diff --git a/src/cosmic/utils.py b/src/cosmic/utils.py index b184e94a6..7c6a9db2d 100644 --- a/src/cosmic/utils.py +++ b/src/cosmic/utils.py @@ -1359,9 +1359,10 @@ def error_check(BSEDict, filters=None, convergence=None, sampling=None): or (BSEDict[flag] == -1) or (BSEDict[flag] == -2) or (BSEDict[flag] == -3) + or (BSEDict[flag] == -4) ): raise ValueError( - "'{0:s}' needs to be set to either 0, greater than 0 or equal to -1, -2, or -3 " + "'{0:s}' needs to be set to either 0, greater than 0 or equal to -1, -2, -3, or -4 " "(you set it to '{1:0.2f}')".format( flag, BSEDict[flag] ) From c5a24d6681e42e6dd66621663fb296ddbe96ea08 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Thu, 27 Mar 2025 16:55:19 -0400 Subject: [PATCH 11/18] created common block variables --- src/cosmic/src/assign_remnant.f | 8 ++++++++ src/cosmic/src/evolv2.f | 9 +++++++++ 2 files changed, 17 insertions(+) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 745103655..c25b7bed5 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -4,6 +4,8 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) INCLUDE 'const_bse.h' common /fall/fallback + REAL*8 preSNmass,preSNmenv,preSNmassc + COMMON preSNmass,preSNmenv,preSNmassc REAL*8 fallback REAL ran3 EXTERNAL ran3 @@ -84,6 +86,12 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) * kw = 15 else +* Beginning of supernova block +* +* Store values in common block + preSNmass = mt + preSNmenv = mcbagb + preSNmassc = mc * Chris Belczynski Evolutionary Roads Weak PPISN * This has to happen before the SNa, because it modifies * the properties of the star during explosion diff --git a/src/cosmic/src/evolv2.f b/src/cosmic/src/evolv2.f index 3b6597154..a31dfefec 100644 --- a/src/cosmic/src/evolv2.f +++ b/src/cosmic/src/evolv2.f @@ -199,6 +199,8 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, EXTERNAL ran3 * + REAL*8 preSNmass, preSNmenv, preSNmassc + COMMON preSNmass, preSNmenv, preSNmassc * REAL*8 z,tm,tn,m0,mt,rm,lum,mc,rc,me,re,k2,age,dtm,dtr REAL*8 tscls(20),lums(10),GB(10),zpars(20) @@ -1368,6 +1370,13 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, b02_bcm = B(2) endif +* print *, evolve_type +* print *, mass0(k), menv(k), massc(k) +* print *, preSNmass, preSNmenv, preSNmassc + mass0(k) = preSNmass + menv(k) = preSNmenv + massc(k) = preSNmassc +* print *, mass0(k), menv(k), massc(k) CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, From 74925603f09cba0542c057e2b0f53115d35d541a Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 9 Apr 2025 13:59:58 -0400 Subject: [PATCH 12/18] worked on block variables --- src/cosmic/src/assign_remnant.f | 8 +-- src/cosmic/src/comenv.f | 90 ++++++++++++++++++++++++--------- src/cosmic/src/evolv2.f | 8 +++ 3 files changed, 78 insertions(+), 28 deletions(-) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index c25b7bed5..207a86929 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -73,6 +73,10 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) mass = mt * else +* Store values in common block + preSNmass = mt + preSNmenv = mcbagb - mc + preSNmassc = mc if(ecsn.gt.0.d0.and.mcbagb.lt.ecsn_mlow)then * * Star is not massive enough to ignite C burning. @@ -88,10 +92,6 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) else * Beginning of supernova block * -* Store values in common block - preSNmass = mt - preSNmenv = mcbagb - preSNmassc = mc * Chris Belczynski Evolutionary Roads Weak PPISN * This has to happen before the SNa, because it modifies * the properties of the star during explosion diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index 6062b7365..c312875d8 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -41,6 +41,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, REAL*8 bhspin1,bhspin2 REAL*8 deltam_1,deltam_2 common /fall/fallback + REAL*8 preSNmass, preSNmenv, preSNmassc + COMMON preSNmass, preSNmenv, preSNmassc INTEGER formation1,formation2 REAL*8 sigmahold REAL*8 AURSUN,K3 @@ -299,37 +301,47 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), + & massc1_bpp,preSNmassc,rad1_bpp,rad2_bpp, + & M02,preSNmass,lumin(2),lumin(1),teff2,teff1, + & RC2,RC1,menv_bpp(2),preSNmenv,renv_bpp(2), & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), & epoch(1),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') + print *, "case 1a" else teff1 = 1000.d0*((1130.d0*lumin(1)/ & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & preSNmassc,massc2_bpp,rad1_bpp,rad2_bpp, + & preSNmass,M02,lumin(1),lumin(2),teff1,teff2, + & RC1,RC2,preSNmenv,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2, & binstate,mergertype,'bpp') + print *, "case 1b" endif endif CALL kick(KW1,M_postCE,M1,M2,ECC,SEP_postCE, @@ -604,37 +616,47 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), + & massc1_bpp,preSNmassc,rad1_bpp,rad2_bpp, + & M02,preSNmass,lumin(2),lumin(1),teff2,teff1, + & RC2,RC1,menv_bpp(2),preSNmenv,renv_bpp(2), & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), & epoch(1),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') + print *, "case 2a" else teff1 = 1000.d0*((1130.d0*lumin(1)/ & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & preSNmassc,massc2_bpp,rad1_bpp,rad2_bpp, + & preSNmass,M02,lumin(1),lumin(2),teff1,teff2, + & RC1,RC2,preSNmenv,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2, & binstate,mergertype,'bpp') + print *, "case 2b" endif endif * USSN: if ussn flag is set, have reduced kicks for stripped He stars (SN=8) @@ -773,37 +795,47 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), + & preSNmassc,massc2_bpp,rad1_bpp,rad2_bpp, + & preSNmass,M01,lumin(2),lumin(1),teff2,teff1, + & RC2,RC1,preSNmenv,menv_bpp(1),renv_bpp(2), & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), & epoch(1),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') + print *, "case 3a" else teff1 = 1000.d0*((1130.d0*lumin(1)/ & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc1_bpp,preSNmassc,rad1_bpp,rad2_bpp, + & M01,preSNmass,lumin(1),lumin(2),teff1,teff2, + & RC1,RC2,menv_bpp(1),preSNmenv,renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2, & binstate,mergertype,'bpp') + print *, "case 3b" endif endif CALL kick(KW2,M_postCE,M2,M1,ECC,SEP_postCE, @@ -1005,37 +1037,47 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,-1.d0,TB,0.d0, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), + & massc1_bpp,preSNmassc,rad1_bpp,rad2_bpp, + & M02,preSNmass,lumin(2),lumin(1),teff2,teff1, + & RC2,RC1,menv_bpp(2),preSNmenv,renv_bpp(2), & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), & epoch(1),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') + print *, "case 4a" else teff1 = 1000.d0*((1130.d0*lumin(1)/ & (rad1_bpp**2.d0))**(1.d0/4.d0)) teff2 = 1000.d0*((1130.d0*lumin(2)/ & (rad2_bpp**2.d0))**(1.d0/4.d0)) +* Load preSN values for the SN writetab + print *, preSNmass, preSNmenv, preSNmassc + print *, M02, menv_bpp(2), massc2_bpp + print *, M01, menv_bpp(1), massc1_bpp CALL writetab(jp,tphys,evolve_type, & mass1_bpp,mass2_bpp,kstar1_bpp, & kstar2_bpp,-1.d0,TB,0.d0, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & preSNmassc,massc2_bpp,rad1_bpp,rad2_bpp, + & preSNmass,M02,lumin(1),lumin(2),teff1,teff2, + & RC1,RC2,preSNmenv,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2, & binstate,mergertype,'bpp') + print *, "case 4b" endif endif CALL kick(KW,MF,M1,0.d0,0.d0,-1.d0,0.d0,vk,star1, diff --git a/src/cosmic/src/evolv2.f b/src/cosmic/src/evolv2.f index a31dfefec..a01781831 100644 --- a/src/cosmic/src/evolv2.f +++ b/src/cosmic/src/evolv2.f @@ -1330,6 +1330,10 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, else b02_bcm = B(2) endif +* Load preSN values for the SN writetab + mass0(k) = preSNmass + menv(k) = preSNmenv + massc(k) = preSNmassc CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, @@ -3598,6 +3602,10 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, else b02_bcm = B(2) endif +* Load preSN values for the SN writetab + mass0(k) = preSNmass + menv(k) = preSNmenv + massc(k) = preSNmassc CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, From 3a2b23d8884018ef049e2e45dcbf795547508855 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Fri, 18 Apr 2025 14:54:07 -0400 Subject: [PATCH 13/18] worked on documents --- src/cosmic/src/comenv.f | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index c312875d8..02ea4b449 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -21,6 +21,23 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, * Update : P. D. Kiel (for ECSN, fallback and bugs) * Date : cmc version mid 2010 * +* +* Note on indexing binary components: +* +* M01, M1, MC1, AJ1, JSPIN1, KW1, formation1, bhspin1 are selected with star1 +* M02, M2, MC2, AJ2, JSPIN2, KW2, formation2, bhspin1 are selected with star2 +* +* star1 is the donor. star2 is the accretor +* switchedCE is .true. if j1 is 2; .false. otherwise +* Nope! +* +* deltam_1 and deltam_2 are calculated before being passed into comenv +* +* teff, radc, and ospin are calculated inside comenv (or a function) +* +* Other quantities are passed as vectors. +* These are lumin, menv_bpp, tms, rad, renv, bacc, tacc, epoch, and B_0. +* tms and rad are immediately separated into tms1/2_bpp and rad1/2_bpp * INTEGER KW1,KW2,KW,KW1i,KW2i,snp INTEGER star1,star2 From 7bcc8bb00ad2dbb834ec9d55767620b349970e14 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Sun, 20 Apr 2025 14:16:36 -0400 Subject: [PATCH 14/18] bugfix (issue 698) --- src/cosmic/src/comenv.f | 55 ++++++++++++++++++++++++++--------------- 1 file changed, 35 insertions(+), 20 deletions(-) diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index 6062b7365..80cdab0db 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -21,6 +21,21 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, * Update : P. D. Kiel (for ECSN, fallback and bugs) * Date : cmc version mid 2010 * +* Update : V. E. Delfavero (bugs in bpp reporting) +* Date : 20th, April 2025 +* +* Note on indexing binary components: +* +* M01, M1, MC1, AJ1, JSPIN1, KW1, formation1, bhspin1 are selected with star1 +* M02, M2, MC2, AJ2, JSPIN2, KW2, formation2, bhspin1 are selected with star2 +* +* deltam_1 and deltam_2 are calculated before being passed into comenv +* +* teff, radc, and ospin are calculated inside comenv (or a function) +* +* Other quantities are passed as vectors. +* These are lumin, menv_bpp, tms, rad, renv, bacc, tacc, epoch, and B_0. +* tms and rad are immediately separated into tms1/2_bpp and rad1/2_bpp * INTEGER KW1,KW2,KW,KW1i,KW2i,snp INTEGER star1,star2 @@ -305,11 +320,11 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), - & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), - & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), - & epoch(1),bhspin2,bhspin1, + & M02,M01,lumin(1),lumin(2),teff1,teff2, + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), + & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), + & epoch(2),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') else @@ -610,11 +625,11 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), - & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), - & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), - & epoch(1),bhspin2,bhspin1, + & M02,M01,lumin(1),lumin(2),teff1,teff2, + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), + & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), + & epoch(2),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') else @@ -779,11 +794,11 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), - & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), - & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), - & epoch(1),bhspin2,bhspin1, + & M02,M01,lumin(1),lumin(2),teff1,teff2, + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), + & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), + & epoch(2),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') else @@ -1011,11 +1026,11 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(2),lumin(1),teff2,teff1, - & RC2,RC1,menv_bpp(2),menv_bpp(1),renv_bpp(2), - & renv_bpp(1),OSPIN2,OSPIN1,B_0(2),B_0(1), - & bacc(2),bacc(1),tacc(2),tacc(1),epoch(2), - & epoch(1),bhspin2,bhspin1, + & M02,M01,lumin(1),lumin(2),teff1,teff2, + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), + & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), + & epoch(2),bhspin2,bhspin1, & deltam_2,deltam_1,formation2,formation1, & binstate,mergertype,'bpp') else From 4318963fe15eba6f7640c856fb592ea6f8035bcf Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Mon, 21 Apr 2025 10:21:29 -0400 Subject: [PATCH 15/18] Added preSN mass tracking --- src/cosmic/src/assign_remnant.f | 6 ++++ src/cosmic/src/comenv.f | 58 +++++++++++++++++++-------------- src/cosmic/src/evolv2.f | 14 ++++++++ 3 files changed, 54 insertions(+), 24 deletions(-) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 858cde738..68a3464d9 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -4,6 +4,8 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) INCLUDE 'const_bse.h' common /fall/fallback + REAL*8 mass_preSN, mHe_preSN, massc_preSN + COMMON mass_preSN, mHe_preSN, massc_preSN REAL*8 fallback REAL ran3 EXTERNAL ran3 @@ -70,6 +72,10 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,mt,kw,bhspin,kidx) mass = mt * else +* Store values in common block + mass_preSN = mt + mHe_preSN = mcbagb - mc + massc_preSN = mc if(ecsn.gt.0.d0.and.mcbagb.lt.ecsn_mlow)then * * Star is not massive enough to ignite C burning. diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index 80cdab0db..25b17eb59 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -56,6 +56,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, REAL*8 bhspin1,bhspin2 REAL*8 deltam_1,deltam_2 common /fall/fallback + REAL*8 mass_preSN, mHe_preSN, massc_preSN + COMMON mass_preSN, mHe_preSN, massc_preSN INTEGER formation1,formation2 REAL*8 sigmahold REAL*8 AURSUN,K3 @@ -319,9 +321,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(1),lumin(2),teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc1_bpp,massc_preSN,rad1_bpp,rad2_bpp, + & M02,mass_preSN,lumin(1),lumin(2), + & teff1,teff2, + & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -337,9 +340,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc_preSN,massc2_bpp,rad1_bpp,rad2_bpp, + & mass_preSN,M02,lumin(1),lumin(2), + & teff1,teff2, + & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -624,9 +628,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(1),lumin(2),teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc1_bpp,massc_preSN,rad1_bpp,rad2_bpp, + & M02,mass_preSN,lumin(1),lumin(2), + & teff1,teff2, + & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -642,9 +647,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc_preSN,massc2_bpp,rad1_bpp,rad2_bpp, + & mass_preSN,M02,lumin(1),lumin(2), + & teff1,teff2, + & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -793,9 +799,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(1),lumin(2),teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc_preSN,massc2_bpp,rad1_bpp,rad2_bpp, + & mass_preSN,M01,lumin(1),lumin(2), + & teff1,teff2, + & RC2,RC1,mHe_preSN,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -811,9 +818,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,SEP_postCE,TB,ECC, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc1_bpp,massc_preSN,rad1_bpp,rad2_bpp, + & M01,mass_preSN,lumin(1),lumin(2), + & teff1,teff2, + & RC1,RC2,menv_bpp(1),mHe_preSN,renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -1025,9 +1033,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,-1.d0,TB,0.d0, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M02,M01,lumin(1),lumin(2),teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc1_bpp,massc_preSN,rad1_bpp,rad2_bpp, + & M02,mass_preSN,lumin(1),lumin(2), + & teff1,teff2, + & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -1043,9 +1052,10 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & kstar2_bpp,-1.d0,TB,0.d0, & rrl1_bpp,rrl2_bpp, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, - & massc1_bpp,massc2_bpp,rad1_bpp,rad2_bpp, - & M01,M02,lumin(1),lumin(2),teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & massc_preSN,massc2_bpp,rad1_bpp,rad2_bpp, + & mass_preSN,M02,lumin(1),lumin(2), + & teff1,teff2, + & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, diff --git a/src/cosmic/src/evolv2.f b/src/cosmic/src/evolv2.f index 083b6c7b8..d022a36da 100644 --- a/src/cosmic/src/evolv2.f +++ b/src/cosmic/src/evolv2.f @@ -195,6 +195,8 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, REAL*8 deltam1_bcm,deltam2_bcm,b01_bcm,b02_bcm REAL*8 B(2),Bbot,omdot,b_mdot,b_mdot_lim,evolve_type COMMON /fall/fallback + REAL*8 mass_preSN, mHe_preSN, massc_preSN + COMMON mass_preSN, mHe_preSN, massc_preSN REAL ran3 EXTERNAL ran3 * @@ -1328,6 +1330,10 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, else b02_bcm = B(2) endif +* Load preSN values for the SN writetab + mass0(k) = mass_preSN + menv(k) = mHe_preSN + massc(k) = massc_preSN CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, @@ -1368,6 +1374,10 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, b02_bcm = B(2) endif +* Load preSN values for the SN writetab + mass0(k) = mass_preSN + menv(k) = mHe_preSN + massc(k) = massc_preSN CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, @@ -3629,6 +3639,10 @@ SUBROUTINE evolv2(kstar,mass,tb,ecc,z,tphysf, else b02_bcm = B(2) endif +* Load preSN values for the SN writetab + mass0(k) = mass_preSN + menv(k) = mHe_preSN + massc(k) = massc_preSN CALL writetab(jp,tphys,evolve_type, & mass(1),mass(2),kstar(1),kstar(2), & sep,tb,ecc,rrl1,rrl2, From bdbe275969207daba1057fa83bebd2a015fb4fa5 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 11 Mar 2026 16:34:45 -0400 Subject: [PATCH 16/18] removed preSN from assign_remnant --- src/cosmic/src/assign_remnant.f | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/cosmic/src/assign_remnant.f b/src/cosmic/src/assign_remnant.f index 1086517b8..fedf83e7d 100644 --- a/src/cosmic/src/assign_remnant.f +++ b/src/cosmic/src/assign_remnant.f @@ -4,8 +4,6 @@ SUBROUTINE assign_remnant(zpars,mc,mcbagb,mass,kidx,mt,kw,bhspin) INCLUDE 'const_bse.h' common /fall/fallback - REAL*8 mass_preSN, mHe_preSN, massc_preSN - COMMON mass_preSN, mHe_preSN, massc_preSN REAL*8 fallback real*8 zpars(20) From c165703dafee2fb2402d6221f4c73c5e493c5ab9 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 11 Mar 2026 16:44:23 -0400 Subject: [PATCH 17/18] Removed preSN from comenv --- src/cosmic/src/comenv.f | 36 +++++++++++++++++------------------- 1 file changed, 17 insertions(+), 19 deletions(-) diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index 83e7ee744..a0bff86c2 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -56,8 +56,6 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, REAL*8 bhspin1,bhspin2 REAL*8 deltam_1,deltam_2 common /fall/fallback - REAL*8 mass_preSN, mHe_preSN, massc_preSN - COMMON mass_preSN, mHe_preSN, massc_preSN INTEGER formation1,formation2 REAL*8 sigmahold REAL*8 AURSUN,K3 @@ -323,9 +321,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & M02,mass_preSN,lumin(1),lumin(2), + & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -343,9 +341,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & mass_preSN,M02,lumin(1),lumin(2), + & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -632,9 +630,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & M02,mass_preSN,lumin(1),lumin(2), + & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -652,9 +650,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & mass_preSN,M02,lumin(1),lumin(2), + & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -805,9 +803,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & mass_preSN,M01,lumin(1),lumin(2), + & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,mHe_preSN,menv_bpp(2),renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -825,9 +823,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & M01,mass_preSN,lumin(1),lumin(2), + & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,menv_bpp(1),mHe_preSN,renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -1041,9 +1039,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & M02,mass_preSN,lumin(1),lumin(2), + & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),mHe_preSN,renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -1061,9 +1059,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp, & mc_he(1),mc_he(2),mc_co(1),mc_co(2), & rad1_bpp,rad2_bpp, - & mass_preSN,M02,lumin(1),lumin(2), + & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,mHe_preSN,menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -1115,7 +1113,7 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & jp,tphys,11.d0,M1,M2,KW1,KW2,-1.d0,-1.d0,-1.d0,0.d0, & 0.d0,aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp,mc_he(1), & mc_he(2),mc_co(1),mc_co(2),rad(1),rad(2),M01,M02,lumin(1), - & lumin(2),teff1,teff2,RC1,RC2,MENV,mHe_preSN,renv_bpp(1), + & lumin(2),teff1,teff2,RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2),bacc(1),bacc(2), & tacc(1),tacc(2),epoch(1),epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2,2,-1, From 83442e9a5d93ddf3c26336cf32e0d622410fcf25 Mon Sep 17 00:00:00 2001 From: Vera Delfavero Date: Wed, 11 Mar 2026 17:06:39 -0400 Subject: [PATCH 18/18] updated comenv --- src/cosmic/_commit_hash.py | 2 +- src/cosmic/src/comenv.f | 29 +++++++++++++++++++---------- 2 files changed, 20 insertions(+), 11 deletions(-) diff --git a/src/cosmic/_commit_hash.py b/src/cosmic/_commit_hash.py index 78bd0bccc..765e68a4e 100644 --- a/src/cosmic/_commit_hash.py +++ b/src/cosmic/_commit_hash.py @@ -1 +1 @@ -COMMIT_HASH = "ffb858c3c70bc4225b1aed6d57ed3872bf4413bb" +COMMIT_HASH = "c165703dafee2fb2402d6221f4c73c5e493c5ab9" diff --git a/src/cosmic/src/comenv.f b/src/cosmic/src/comenv.f index a0bff86c2..1acfe3e43 100644 --- a/src/cosmic/src/comenv.f +++ b/src/cosmic/src/comenv.f @@ -323,7 +323,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -343,7 +344,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -632,7 +634,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -652,7 +655,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -805,7 +809,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -825,7 +830,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -1041,7 +1047,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M02,M01,lumin(1),lumin(2), & teff1,teff2, - & RC2,RC1,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC2,RC1,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN2,OSPIN1,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin2,bhspin1, @@ -1061,7 +1068,8 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & rad1_bpp,rad2_bpp, & M01,M02,lumin(1),lumin(2), & teff1,teff2, - & RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), + & RC1,RC2,menv_bpp(1),menv_bpp(2), + & renv_bpp(1), & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2), & bacc(1),bacc(2),tacc(1),tacc(2),epoch(1), & epoch(2),bhspin1,bhspin2, @@ -1113,8 +1121,9 @@ SUBROUTINE COMENV(M01,M1,MC1,AJ1,JSPIN1,KW1, & jp,tphys,11.d0,M1,M2,KW1,KW2,-1.d0,-1.d0,-1.d0,0.d0, & 0.d0,aj1_bpp,aj2_bpp,tms1_bpp,tms2_bpp,mc_he(1), & mc_he(2),mc_co(1),mc_co(2),rad(1),rad(2),M01,M02,lumin(1), - & lumin(2),teff1,teff2,RC1,RC2,menv_bpp(1),menv_bpp(2),renv_bpp(1), - & renv_bpp(2),OSPIN1,OSPIN2,B_0(1),B_0(2),bacc(1),bacc(2), + & lumin(2),teff1,teff2,RC1,RC2,menv_bpp(1),menv_bpp(2), + & renv_bpp(1),renv_bpp(2), + & OSPIN1,OSPIN2,B_0(1),B_0(2),bacc(1),bacc(2), & tacc(1),tacc(2),epoch(1),epoch(2),bhspin1,bhspin2, & deltam_1,deltam_2,formation1,formation2,2,-1, & zpars(14)**2.d5,'bpp')