Metisse integrate merge

.github/workflows/build_wheels_and_publish.yml

@@ -10,8 +10,10 @@ env:
   CIBW_ARCHS_LINUX: "x86_64"
   CIBW_SKIP: "*-win32 *musllinux*"
   CIBW_MANYLINUX_X86_64_IMAGE: manylinux2014
+  CIBW_ARCHS_MACOS: "x86_64 arm64"
   CIBW_BEFORE_BUILD: pip install numpy --config-settings=setup-args="-Dallow-noblas=true"
   CIBW_BUILD_VERBOSITY: "1"
+  CIBW_ENVIRONMENT_MACOS: "FC=gfortran"
 
 jobs:
   build-wheels-and-dist:
@@ -26,11 +28,15 @@ jobs:
     - uses: actions/checkout@v4
       with:
         fetch-depth: 0
+        submodules: recursive
+    - name: Switch METISSE submodule to track develop branch
+      run: |
+        git submodule set-branch --branch develop src/cosmic/src/METISSE
+        git submodule update --remote src/cosmic/src/METISSE
     - uses: actions/setup-python@v5
       name: Install Python
       with:
         python-version: "3.10"
-
     - name: link gfortran and hdf5
       if: runner.os == 'macOS'
       run: |

.github/workflows/python-package.yml

@@ -15,7 +15,15 @@ jobs:
         python-version: ['3.9', '3.10']
 
     steps:
-    - uses: actions/checkout@v4
+    - name: Checkout repository with submodules
+      uses: actions/checkout@v4
+      with:
+        submodules: recursive
+        fetch-depth: 0  # needed for full submodule history
+    - name: Switch METISSE submodule to track develop branch
+      run: |
+        git submodule set-branch --branch develop src/cosmic/src/METISSE
+        git submodule update --remote src/cosmic/src/METISSE
     - name: Set up Python ${{ matrix.python-version }}
       uses: actions/setup-python@v4
       with:

.gitmodules

@@ -0,0 +1,4 @@
+[submodule "src/cosmic/src/METISSE"]
+	path = src/cosmic/src/METISSE
+	url = https://github.com/TeamMETISSE/METISSE.git
+	branch = develop

bin/cosmic-pop

@@ -70,7 +70,7 @@ def parse_commandline():
 
     defaults = {}
     if not (args.inifile is None and (('-h' in remaining_argv) or ('--help' in remaining_argv))):
-        BSEDict, seed_int, filters, convergence, sampling = utils.parse_inifile(args.inifile)
+        BSEDict, SSEDict, seed_int, filters, convergence, sampling = utils.parse_inifile(args.inifile)
         defaults.update(sampling)
         defaults.update(filters)
         defaults.update(convergence)
@@ -167,8 +167,9 @@ if __name__ == '__main__':
     pool = schwimmbad.choose_pool(mpi=args.mpi, processes=args.nproc)
     if isinstance(pool, MPIPool):
         if not pool.is_master():
+            from sys import exit
             pool.wait()
-            sys.exit(0)
+            exit(0)
         nproc = len(pool.workers)
     else:
         nproc = args.nproc
@@ -179,7 +180,7 @@ if __name__ == '__main__':
 
     # READ AND PARSE INIFILE
     ###########################################################################
-    BSEDict, seed_int, filters, convergence, sampling = utils.parse_inifile(args.inifile)
+    BSEDict, SSEDict, seed_int, filters, convergence, sampling = utils.parse_inifile(args.inifile)
 
     # we now overwrite the inifile values with what was specified from the command line
     # (which could mean not overwriting anything at all because they are populated
@@ -215,8 +216,8 @@ if __name__ == '__main__':
                               "with {0}={2} from the commandline".format(argument, seed_int, getattr(args, argument)))
                 seed_int = getattr(args, argument)
 
-    # Check that the values in BSEDict, filters, and convergence are valid
-    utils.error_check(BSEDict, filters, convergence, sampling)
+    # Check that the values in BSEDict, SSEDict, filters, and convergence are valid
+    utils.error_check(BSEDict, SSEDict, filters, convergence, sampling)
 
     if seed_int != 0:
         np.random.seed(seed_int)
@@ -267,7 +268,7 @@ if __name__ == '__main__':
         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']), 'w')
 
     # save configuration settings to output file
-    configuration_settings = {'BSEDict' : BSEDict, 'filters' : filters,
+    configuration_settings = {'BSEDict' : BSEDict, 'SSEDict': SSEDict, 'filters' : filters,
                               'convergence' : convergence, 'sampling' : sampling}
 
     with warnings.catch_warnings():
@@ -298,6 +299,9 @@ 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))
 
+    if sampling["sampling_method"] == "multidim" and args.Nstep < args.nproc ** 2:
+            warnings.warn("`Nstep` is less than the square of `nproc`. You will see repeats with multidim sampling. Consider raising `Nstep`.")
+
     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
@@ -322,8 +326,8 @@ if __name__ == '__main__':
                                                             met = sampling['metallicity'],
                                                             size = args.Nstep,
                                                             qmin = args.qmin,
-                                                            params = args.inifile)
-            # if m2_min is specified, use it to sample the initial binary table
+                                                            params = args.inifile,
+                                                            SSEDict = SSEDict)
             elif hasattr(args,'m2_min'):
                 init_samp_list = InitialBinaryTable.sampler(format_ = sampling['sampling_method'],
                                                             final_kstar1 = kstar1_range,
@@ -338,7 +342,8 @@ if __name__ == '__main__':
                                                             met = sampling['metallicity'],
                                                             size = args.Nstep,
                                                             m2_min = args.m2_min,
-                                                            params = args.inifile)
+                                                            params = args.inifile,
+                                                            SSEDict = SSEDict)
             else:
                 raise ValueError("You must specify either qmin or m2_min in the",
                                  " inifile if you are using the independent sampler")
@@ -389,6 +394,7 @@ if __name__ == '__main__':
         bpp, bcm, initCond, kick_info = Evolve.evolve(initialbinarytable=IBT,
                                                       pool=pool,
                                                       BSEDict=BSEDict,
+                                                      SSEDict=SSEDict,
                                                       idx=idx,
                                                       dtp=dtp,
                                                       timestep_conditions=filters['timestep_conditions'])

changelog.md

@@ -58,4 +58,4 @@ See the discussed changes in our previous releases here: https://github.com/COSM
 
 ## 3.6.2
  - Add functions to `cosmic.utils` for initC IO that's more efficient (`save_initC`, `load_initC`) by saving
- identical setting columns separately with only one copy - saves ~1kb per binary
+ identical setting columns separately with only one copy - saves ~1kb per binary

ci/compile_benchmark.sh

@@ -1,2 +1,59 @@
-gfortran -coverage -fprofile-arcs -ftest-coverage -O0 src/cosmic/src/assign_remnant.f src/cosmic/src/benchmarkevolv2.f src/cosmic/src/comenv.f src/cosmic/src/corerd.f src/cosmic/src/deltat.f src/cosmic/src/dgcore.f src/cosmic/src/evolv2.f src/cosmic/src/gntage.f src/cosmic/src/hrdiag.f src/cosmic/src/hrdiag_remnant.f src/cosmic/src/instar.f src/cosmic/src/kick.f src/cosmic/src/mix.f src/cosmic/src/mlwind.f src/cosmic/src/mrenv.f src/cosmic/src/ran3.f src/cosmic/src/rl.f src/cosmic/src/star.f src/cosmic/src/zcnsts.f src/cosmic/src/zfuncs.f src/cosmic/src/concatkstars.f src/cosmic/src/bpp_array.f src/cosmic/src/checkstate.f  -o benchmarkevolv2.exe -I src/cosmic/src -Wl,-rpath,${CONDA_PREFIX}/lib
+#!/usr/bin/env bash
+set -euo pipefail
+
+# Directories
+SRC_DIR="src/cosmic/src"
+METISSE_DIR="$SRC_DIR/METISSE/src"
+
+# Compiler flags (removed rpath to CONDA_PREFIX)
+FFLAGS="-coverage -fprofile-arcs -ftest-coverage -O0 -J$SRC_DIR -I$SRC_DIR"
+
+# Phase 1: Compile METISSE modules in dependency order
+gfortran $FFLAGS -c \
+    $METISSE_DIR/track_support.f90 \
+    $METISSE_DIR/c_m_interface.f90 \
+    $METISSE_DIR/z_support.f90 \
+    $METISSE_DIR/sse_support.f90 \
+    $METISSE_DIR/remnant_support.f90 \
+    $METISSE_DIR/interp_support.f90 \
+    $METISSE_DIR/METISSE_gntage.f90 \
+    $METISSE_DIR/METISSE_deltat.f90 \
+    $METISSE_DIR/METISSE_mlwind.f90 \
+    $METISSE_DIR/METISSE_hrdiag.f90 \
+    $METISSE_DIR/METISSE_star.f90 \
+    $METISSE_DIR/METISSE_zcnsts.f90 \
+    $METISSE_DIR/comenv_lambda.f90 \
+    $METISSE_DIR/METISSE_miscellaneous.f90 \
+    $SRC_DIR/METISSE_utils.f90
+
+# Phase 2: Compile COSMIC and SSE sources + link everything
+gfortran $FFLAGS \
+    $SRC_DIR/hrdiag_remnant.f \
+    $SRC_DIR/assign_remnant.f \
+    $SRC_DIR/benchmarkevolv2.f \
+    $SRC_DIR/corerd.f \
+    $SRC_DIR/comenv.f \
+    $SRC_DIR/dgcore.f \
+    $SRC_DIR/evolv2.f \
+    $SRC_DIR/gntage.f \
+    $SRC_DIR/instar.f \
+    $SRC_DIR/kick.f \
+    $SRC_DIR/mix.f \
+    $SRC_DIR/mrenv.f \
+    $SRC_DIR/ran3.f \
+    $SRC_DIR/rl.f \
+    $SRC_DIR/concatkstars.f \
+    $SRC_DIR/comprad.f \
+    $SRC_DIR/bpp_array.f \
+    $SRC_DIR/checkstate.f \
+    $SRC_DIR/deltat.f \
+    $SRC_DIR/mlwind.f \
+    $SRC_DIR/hrdiag.f \
+    $SRC_DIR/star.f \
+    $SRC_DIR/zcnsts.f \
+    $SRC_DIR/SSE/SSE_*.f \
+    *.o \
+    -o benchmarkevolv2.exe
+
+# Run the benchmark
 ./benchmarkevolv2.exe

docs/cosmic-settings.json

@@ -664,15 +664,15 @@
                 "settings-section-description": "Note: there are cases where a common envelope is forced regardless of the critical mass ratio for unstable mass transfer. In the following cases, a common envelope occurs regardless of the choices below:<ul><li><b>contact</b>: the stellar radii go into contact (common for similar ZAMS systems)</li><li><b>periapse contact</b>: the periapse distance is smaller than either of the stellar radii (common for highly eccentric systems)</li><li><b>core Roche overflow</b>: either of the stellar radii overflow their component's Roche radius (in this case, mass transfer from the convective core is always dynamically unstable)</li></ul>",
                 "name": "alpha1",
                 "description": "Common-envelope efficiency parameter which scales the efficiency of transferring orbital energy to the envelope. See <a href='https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract'>Hurley+2002, Eq. 71</a>.",
-                "type": "number",
-                "options-preface": "This value is equivalent to \\(\\alpha\\) in <a href='https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract'>Hurley+2002, Eq. 71</a>.",
+                "type": "string",
+                "options-preface": "This value is equivalent to \\(\\alpha\\) in <a href='https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract'>Hurley+2002, Eq. 71</a>. A list is supplied such that the primary and secondary can be assigned different values for \\(\\alpha\\).",
                 "options": [
                     {
                         "name": "positive values",
-                        "description": "Sets the common-envelope efficiency parameter, \\(\\alpha\\)"
+                        "description": "Sets the common-envelope efficiency parameter, \\(\\alpha\\), for each binary component"
                     },
                     {
-                        "name": 1.0,
+                        "name": "[1.0, 1.0]",
                         "description": "Default value",
                         "default": true
                     }
@@ -681,7 +681,7 @@
             {
                 "name": "lambdaf",
                 "description": "Binding energy factor for common envelope evolution. The initial binding energy of the stellar envelope is proportional to \\(1 / \\lambda\\). See <a href='https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract'>Hurley+2002, Eq. 69</a>.",
-                "type": "number",
+                "type": "string",
                 "options-preface": "",
                 "options": [
                     {
@@ -851,7 +851,7 @@
                     },
                     {
                         "name": "negative values",
-                        "description": "Same as above settings but using the old Kiel & Hurley 2009 prescription for changing the orbital configuration of the binary, available for reproducibility purposes but not recommended for new work"
+                        "description": "Same as above settings but using the old Kiel & Hurley 2009 prescription for changing the orbital configuration of the binary, available for reproducibility purposes but not recommended for new work. kickflag=-5 is NOT implemented."
                     }
                 ]
             },
@@ -1298,30 +1298,34 @@
             },
             {
                 "name": "acc_lim",
-                "description": "Limits the amount of mass accreted during Roche-lobe overflow",
-                "type": "number",
+                "description": "Limits the amount of mass accreted during Roche-lobe overflow where each binary component must have its own limit specified in a list",
+                "type": "string",
                 "options-preface": "",
                 "options": [
                     {
-                        "name": -1,
-                        "description": "limited to 10x the thermal rate of the accretor for MS/HG/CHeB and unlimited for GB/EAGB/AGB stars",
-                        "default": true
+                        "name": "-1",
+                        "description": "limited to 10x the thermal rate of the accretor for MS/HG/CHeB and unlimited for GB/EAGB/AGB stars"
                     },
                     {
-                        "name": -2,
+                        "name": "-2",
                         "description": "limited to 1x the thermal rate of the accretor for MS/HG/CHeB and unlimited for GB/EAGB/AGB stars"
                     },
                     {
-                        "name": -3,
+                        "name": "-3",
                         "description": "limited to 10x the thermal rate of the accretor for all stars"   
                     },
                     {
-                        "name": -4,
+                        "name": "-4",
                         "description": "limited to 1x the thermal rate of the accretor for all stars"
                     },
                     {
                         "name": ">= 0",
                         "description": "sets overall fraction of donor material that is accreted, with the rest being lost from the system (<code>acc_lim = 0.5</code> assumes 50% accretion efficiency as in <a class='reference external' href='https://ui.adsabs.harvard.edu/abs/2008ApJS..174..223B/abstract'>Belczynski+2008</a>)"
+                    },
+                    {
+                        "name": "[-1, -1]",
+                        "description": "Default choice",
+                        "default": true
                     }
                 ]
             },

docs/pages/about.rst

@@ -5,7 +5,13 @@ About
 #####
 COSMIC (Compact Object Synthesis and Monte Carlo Investigation Code) is a rapid binary population synthesis suite with a special focus of generating compact binary populations. 
 
-COSMIC currently implements stellar evolution using SSE (`Hurley, Pols, and Tout 2000 <https://ui.adsabs.harvard.edu/abs/2000MNRAS.315..543H/abstract>`_) and binary interactions using BSE (`Hurley, Tout, and Pols 2002 <https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract>`_). Several modifications have been applied to BSE to account for recent updates to binary evolution especially important to compact binary formation (e.g. metallicity-dependent stellar winds or black hole natal kick strengths). For a detailed discussion of these modifications, see `Breivik et al. 2020 <https://ui.adsabs.harvard.edu/abs/2019arXiv191100903B/abstract>`_.
+COSMIC implements stellar evolution using the fitting formulae described in
+SSE (`Hurley, Pols, and Tout 2000 <https://ui.adsabs.harvard.edu/abs/2000MNRAS.315..543H/abstract>`_) 
+or METISSE which uses pre-computed detailed single star grids that follow the MIST-style equivalent 
+evolutionary phase (EEP) format. For more details on METISSE see the `METISSE documentation <https://cosmic-popsynth.github.io/Metisse/>`_. 
+
+Binary interactions are implemented based on the BSE algorithm (`Hurley, Tout, and Pols 2002 <https://ui.adsabs.harvard.edu/abs/2002MNRAS.329..897H/abstract>`_). 
+Several modifications have been applied to BSE to account for recent updates to binary evolution especially important to compact binary formation (e.g. metallicity-dependent stellar winds or black hole natal kick strengths). For a detailed discussion of these modifications, see `Breivik et al. 2020 <https://ui.adsabs.harvard.edu/abs/2019arXiv191100903B/abstract>`_.
 
 ************
 Using COSMIC

docs/pages/cite.rst

@@ -6,6 +6,9 @@ If you use COSMIC in your work, please cite both:
 
 - The COSMIC release paper: `Breivik et al. 2020 <https://ui.adsabs.harvard.edu/abs/2019arXiv191100903B/abstract>`_
 - The specific Zenodo DOI associated with version of COSMIC you used
+- The METISSE stellar evolution package: JOSS citation TBA
+
+If you used pre-supplied stellar tracks from METISSE, please also cite the relevant papers for those tracks as listed in the `METISSE documentation <https://cosmic-popsynth.github.io/Metisse/steller_models.html>`_.
 
 We recommend using the `Software Citation Station <https://tomwagg.com/software-citation-station>`_ to generate a citation
 for COSMIC (and any other software you use in your work!) - it's a great way to ensure you're citing software correctly and