Switch to fix time stepping scheme.

.github/workflows/tests.yaml

@@ -12,9 +12,10 @@ on:
 jobs:
   test:
 
-    runs-on: ubuntu-latest
+    runs-on: ${{ matrix.os }}
     strategy:
       matrix:
+        os: [ubuntu-latest, macos-latest]
         python-version: ['3.10', '3.11', '3.12']
 
     steps:

src/mors/parameters.py

@@ -18,7 +18,7 @@ def SetDefaultParameters(paramsDefault):
     """Sets up the default parameters."""
 
     # Time integration solver parameters
-    paramsDefault['TimeIntegrationMethod'] = 'Rosenbrock'     # options are 'ForwardEuler', 'RungeKutta4', 'RungeKuttaFehlberg', 'Rosenbrock'
+    paramsDefault['TimeIntegrationMethod'] = 'ForwardEuler'   # options are 'ForwardEuler', 'RungeKutta4', 'RungeKuttaFehlberg', 'Rosenbrock'
     paramsDefault['nStepMax'] = 10**6                         # maximum number of timesteps to allow before error
     paramsDefault['DeltaDesired'] = 1.0e-5                    # desired Delta value to use in numerical solvers such as Runge-Kutta-Fehlberg and Rosenbrock
     paramsDefault['AgeMinDefault'] = 1.0                      # Myr - default age to start evolutionary tracks

src/mors/stellarevo.py

@@ -747,17 +747,11 @@ def _CheckMassLimit(MstarArray,Mstar):
 def _Interpolate2D(Z1,Z2,Z,Xarray1,Xarray2,X,Yarray1,Yarray2):
     """Takes two sets of 1D arrays for corresponding X and Y values, returns interpolated Y value corresponding to input X."""
 
-    # Round the input values to decimal places determined by nDecValue at top of file
-    #Z = np.round( Z , decimals=nDecValue )
-
     # Do interpolations to get Y at X for both tracks
     Y1 = _Interpolate1D( Xarray1 , Yarray1 , X )
     Y2 = _Interpolate1D( Xarray2 , Yarray2 , X )
 
     # Do linear interpolation between Y1 and Y2 in Z direction
-    #mInterp = ( Y2 - Y1 ) / ( Z2 - Z1 )
-    #cInterp = Y1 - mInterp * Z1
-    #Y = mInterp * Z + cInterp
     delta = (Z-Z1)/(Z2-Z1)
     Y = Y2*delta + Y1*(1-delta)
 
@@ -768,9 +762,6 @@ def _Interpolate1D(Xarray,Yarray,X):
 
     # Note that it is assumed here that Xarray is in ascending order and this won't work if it is not
 
-    # Round the input values to decimal places determined by nDecValue at top of file
-    #X = np.round( X , decimals=nDecValue )
-
     # Make sure X is in limits
     if not ( Xarray[0] <= X <= Xarray[-1] ):
         raise Exception("input value "+str(X)+" is not within limits of "+str(Xarray[0])+" to "+str(Xarray[-1]))
@@ -787,9 +778,6 @@ def _Interpolate1D(Xarray,Yarray,X):
     else:
 
         # Do linear interpolation
-        #mInterp = ( Yarray[iMax] - Yarray[iMin] ) / ( Xarray[iMax] - Xarray[iMin] )
-        #cInterp = Yarray[iMin] - mInterp * Xarray[iMin]
-        #Y = mInterp * X + cInterp
         delta = (X-Xarray[iMin])/(Xarray[iMax]-Xarray[iMin])
         Y = Yarray[iMax]*delta + Yarray[iMin]*(1-delta)
 

tests/test_spada.py

@@ -3,10 +3,10 @@
 from numpy.testing import assert_allclose
 
 TEST_DATA = (
-    ((0.128, 45.2, 8.5e1),(0.21263373, 1.68612613e+31, 1.74497354e+28)),
-    ((1.113, 17.7, 3.2e3),(1.16581856, 6.81769993e+33, 7.95296357e+28)),
-    ((0.9999, 1.0001, 1.0e4),(1.52546223e+00, 8.12491693e+33, 3.19515309e+28)),
-    ((1.000, 1.00, 1.0e4),(1.52583654e+00, 8.13190335e+33, 3.20176523e+28)),
+    ((0.128, 45.2, 8.5e1),(0.21263373, 1.68612614e+31, 1.73275380e+28)),
+    ((1.113, 17.7, 3.2e3),(1.16581827, 6.81769926e+33, 7.94292986e+28)),
+    ((0.995, 1.005, 1.0e4),(1.48904952e+00, 7.80667546e+33, 3.22933432e+28)),
+    ((1.000, 1.00, 1.0e4),(1.52583342e+00, 8.13191736e+33, 3.38245392e+28)),
 )
 
 @pytest.mark.parametrize("inp,expected", TEST_DATA)