diff --git a/docs/source/notebook_luminosity_function.ipynb b/docs/source/notebook_luminosity_function.ipynb
index 984e21e8ecb890093b9ad5dac1da5b7c2d531e36..7c6163c1a750f199e45daf788bed23f1081b2f86 100644
--- a/docs/source/notebook_luminosity_function.ipynb
+++ b/docs/source/notebook_luminosity_function.ipynb
@@ -36,12 +36,7 @@
"## Setting up the Population object\n",
"To set up and configure the population object we need to make an object instance of the `Population` object, and add configuration via the `.set()` function.\n",
"\n",
- "There are three categories of options that the Population object can set:\n",
- "- BSE (stellar evolution) options: these options will be used for the binary_c calls, and are recognized by comparing the arguments to a known list of available arguments of binary_c. To see which options are available, see section [`binary_c parameters` in the documentation](https://ri0005.pages.surrey.ac.uk/binary_c-python/binary_c_parameters.html). You can access these through `population.bse_options['<bse option name>']` after you have set them. \n",
- "\n",
- "- Grid options: these options will be used to configure the behaviour of the Population object. To see which options are available, see section [`Population grid code options` in the documentation](https://ri0005.pages.surrey.ac.uk/binary_c-python/grid_options_descriptions.html). They can be accessed via `population.grid_options['<grid option name>']` after you have set them. \n",
- "\n",
- "- Custom options: these options are not recognized as either of the above, so they will be stored in the custom_options, and can be accessed via `population.custom_options['<custom option name>']`"
+ "In our case, we only need to set the maximum evolution time to something short, because we care only about zero-age main sequence stars which have, by definition, age zero."
]
},
{
@@ -54,7 +49,7 @@
"name": "stdout",
"output_type": "stream",
"text": [
- "adding: max_evolution_time=15000 to BSE_options\n",
+ "adding: max_evolution_time=0.1 to BSE_options\n",
"verbosity is 1\n"
]
}
@@ -72,7 +67,7 @@
"# All other values that you set are put in a custom_options dict\n",
"population.set(\n",
" # binary_c physics options\n",
- " max_evolution_time=15000, # maximum stellar evolution time in Myr\n",
+ " max_evolution_time=0.1, # maximum stellar evolution time in Myr\n",
" )\n",
"\n",
"# We can access the options through \n",
@@ -295,9 +290,8 @@
"Total starcount for this run will be: 40\n",
"Generating grid code\n",
"Constructing/adding: M_1\n",
- "Population-e9197a7369484dc9ab3e0a46d90ddae4 finished! The total probability was: 1.0000000000000002. It took a total of 2.6588425636291504s to run 40 systems on 2 cores\n",
- "There were no errors found in this run.\n",
- "OrderedDict([('luminosity distribution', OrderedDict([(0.25, 0.025), (3.75, 0.05), (4.25, 0.05), (2.25, 0.025), (3.25, 0.025), (5.25, 0.2), (4.75, 0.1), (5.75, 0.39999999999999997), (6.25, 0.125)]))])\n"
+ "Population-4dfe2a0d894b4c0ca1670b4a43ec1235 finished! The total probability was: 1.0000000000000002. It took a total of 2.68017315864563s to run 40 systems on 2 cores\n",
+ "There were no errors found in this run.\n"
]
}
],
@@ -313,8 +307,8 @@
"# Evolve the population - this is the slow, number-crunching step\n",
"analytics = population.evolve() \n",
"\n",
- "# Show the results\n",
- "print (population.grid_results)"
+ "# Show the results (debugging)\n",
+ "# print (population.grid_results)"
]
},
{
@@ -335,7 +329,7 @@
"name": "stdout",
"output_type": "stream",
"text": [
- "{'population_name': 'e9197a7369484dc9ab3e0a46d90ddae4', 'evolution_type': 'grid', 'failed_count': 0, 'failed_prob': 0, 'failed_systems_error_codes': [], 'errors_exceeded': False, 'errors_found': False, 'total_probability': 1.0000000000000002, 'total_count': 40, 'start_timestamp': 1631110656.7371156, 'end_timestamp': 1631110659.3959582, 'total_mass_run': 2001.4, 'total_probability_weighted_mass_run': 50.035000000000004, 'zero_prob_stars_skipped': 0}\n"
+ "{'population_name': '4dfe2a0d894b4c0ca1670b4a43ec1235', 'evolution_type': 'grid', 'failed_count': 0, 'failed_prob': 0, 'failed_systems_error_codes': [], 'errors_exceeded': False, 'errors_found': False, 'total_probability': 1.0000000000000002, 'total_count': 40, 'start_timestamp': 1631111693.006838, 'end_timestamp': 1631111695.6870112, 'total_mass_run': 2001.3999999999999, 'total_probability_weighted_mass_run': 50.035, 'zero_prob_stars_skipped': 0}\n"
]
}
],
@@ -397,6 +391,14 @@
"\n"
]
},
+ {
+ "cell_type": "markdown",
+ "id": "7d7b275e-be92-4d59-b44d-ef6f24023cc3",
+ "metadata": {},
+ "source": [
+ "Does this look like a reasonable stellar luminosity function to you? The implication is that the most likely stellar luminosity is 10<sup>5.8</sup> L<sub>☉</sub>! Clearly, this is not very realistic... let's see what went wrong."
+ ]
+ },
{
"cell_type": "markdown",
"id": "44586e42-b7cb-4a55-be0a-330b98b20de4",
@@ -412,7 +414,7 @@
"source": [
"## ZAMS Luminosity distribution with the initial mass function\n",
"\n",
- "We now include an initial mass function as a three-part power law based on Kroupa (2001).\n"
+ "In the previous example, all the stars in our grid had an equal weighting. This is very unlikely to be true in reality: indeed, we know that low mass stars are far more likely than high mass stars. So we now include an initial mass function as a three-part power law based on Kroupa (2001). Kroupa's distribution is a three-part power law: we have a function that does this for us (it's very common to use power laws in astrophysics).\n"
]
},
{
@@ -446,9 +448,8 @@
"Total starcount for this run will be: 40\n",
"Generating grid code\n",
"Constructing/adding: M_1\n",
- "Population-f9545fb319e04882816b3125d699dd00 finished! The total probability was: 0.21822161894107872. It took a total of 2.6835670471191406s to run 40 systems on 2 cores\n",
- "There were no errors found in this run.\n",
- "OrderedDict([('luminosity distribution', OrderedDict([(2.25, 0.0164166), (3.75, 0.0037453900000000004), (4.25, 0.0014346559999999999), (0.25, 0.189097), (3.25, 0.00515685), (5.25, 0.0007493004), (4.75, 0.001171479), (5.75, 0.00039801020000000003), (6.25, 5.2369339999999996e-05)]))])\n"
+ "Population-313b45782bbb4191a0ab756bbe893768 finished! The total probability was: 0.21822161894107872. It took a total of 2.811074733734131s to run 40 systems on 2 cores\n",
+ "There were no errors found in this run.\n"
]
}
],
@@ -457,8 +458,8 @@
"population.clean()\n",
"analytics = population.evolve() \n",
"\n",
- "# Show the results\n",
- "print (population.grid_results)"
+ "# Show the results (debugging)\n",
+ "# print (population.grid_results)"
]
},
{
@@ -511,7 +512,11 @@
"id": "0546f2f3-4732-4841-8ef3-565fbf6b9961",
"metadata": {},
"source": [
- "This distribution is peaked at low luminosity, as one expects from observations, but the resolution is clearly not great. We could simply increase the resolution to compensate, but this is very CPU intensive and a waste of time. Instead, let's try sampling the masses of the stars in a more smart way."
+ "This distribution is peaked at low luminosity, as one expects from observations, but the resolution is clearly not great because it's not smooth - it's spiky! \n",
+ "\n",
+ "If you noticed above, the total probability of the grid was about 0.2. Given that the total probability of a probability distribution function should be 1.0, this shows that our sampling is (very) poor. \n",
+ "\n",
+ "We could simply increase the resolution to compensate, but this is very CPU intensive and a complete waste of time and resources. Instead, let's try sampling the masses of the stars in a smarter way."
]
},
{
@@ -533,7 +538,7 @@
"metadata": {},
"outputs": [],
"source": [
- "# Rename the old variable because we want it to be called lnM_1 now\n",
+ "# Rename the old variable (M_1) because we want it to be called lnM_1 now\n",
"population.rename_grid_variable(\"M_1\",\"lnM_1\")"
]
},
@@ -550,15 +555,7 @@
"execution_count": 15,
"id": "108d470a-bb21-40b0-8387-2caa7ab0f923",
"metadata": {},
- "outputs": [
- {
- "name": "stdout",
- "output_type": "stream",
- "text": [
- "{'lnM_1': {'name': 'lnM_1', 'longname': 'Primary mass', 'valuerange': [0.07, 100.0], 'resolution': '40', 'spacingfunc': 'const(math.log(0.07), math.log(100.0), 40)', 'precode': 'M_1=math.exp(lnM_1)', 'probdist': 'three_part_powerlaw(M_1, 0.1, 0.5, 1.0, 150, -1.3, -2.3, -2.3)*M_1', 'dphasevol': 'dlnM_1', 'parameter_name': 'M_1', 'condition': '', 'gridtype': 'centred', 'branchpoint': 0, 'grid_variable_number': 0}}\n"
- ]
- }
- ],
+ "outputs": [],
"source": [
"# update the sampling, note that the IMF is dprob/dM1, and the phase \n",
"# space is now sampled in lnM1, so we multiply by M_1 to \n",
@@ -571,7 +568,7 @@
" parameter_name=\"M_1\",\n",
" precode=\"M_1=math.exp(lnM_1)\",\n",
")\n",
- "print(population.grid_options[\"_grid_variables\"])"
+ "# print(population.grid_options[\"_grid_variables\"]) # debugging"
]
},
{
@@ -591,9 +588,8 @@
"Total starcount for this run will be: 40\n",
"Generating grid code\n",
"Constructing/adding: lnM_1\n",
- "Population-5f7b694e1e06428d803afae080c9fbfe finished! The total probability was: 0.9956307907476223. It took a total of 1.7812578678131104s to run 40 systems on 2 cores\n",
- "There were no errors found in this run.\n",
- "OrderedDict([('luminosity distribution', OrderedDict([(-0.25, 0.0268827), (1.25, 0.0104553), (0.25, 0.0212294), (2.75, 0.00321118), (-3.25, 0.0), (6.25, 7.34708e-05), (-0.75, 0.0771478), (0.75, 0.030004499999999996), (2.25, 0.00921541), (3.25, 0.0045385), (1.75, 0.014776889999999999), (3.75, 0.00283037), (4.25, 0.002380189), (4.75, 0.000869303), (5.25, 0.0007310379999999999), (5.75, 0.00036002859999999996), (-2.75, 0.1961345), (-1.75, 0.2181597), (-2.25, 0.2568974), (-1.25, 0.11973310000000001)]))])\n"
+ "Population-a7b836b9f38f4c97901a027c7206330f finished! The total probability was: 0.9956307907476223. It took a total of 1.7914540767669678s to run 40 systems on 2 cores\n",
+ "There were no errors found in this run.\n"
]
}
],
@@ -602,8 +598,16 @@
"population.clean()\n",
"analytics = population.evolve() \n",
"\n",
- "# Show the results\n",
- "print (population.grid_results)"
+ "# Show the results (debugging)\n",
+ "# print (population.grid_results)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "182b1094-5057-4ccf-bac6-9b0e560ad4f6",
+ "metadata": {},
+ "source": [
+ "You should see that the total probability is very close to 1.0, as you would expect for a well-sampled grid. The total will never be exactly 1.0, but that is because we are running a simulation, not a perfect copy of reality."
]
},
{
@@ -658,9 +662,9 @@
"source": [
"Most stars are low mass red dwarfs, with small luminosities. Without the IMF weighting, our model population would have got this completely wrong! \n",
"\n",
- "As you increase the resolution, you will see this curve becomes even smoother. The wiggles in the curve are (usually) sampling artefacts because the curve should monotonically brighter than about log(*L*/L<sub>☉</sub>)=-2. \n",
+ "As you increase the resolution, you will see this curve becomes even smoother. The wiggles in the curve are (usually) sampling artefacts because the curve should monotonically brighten above about log(*L*/L<sub>☉</sub>)=-2. \n",
" \n",
- "Remember you can play with the binwidth too. If you want a very accurate distribution you need a narrow binwidth, but then you'll also need high resolution so lots of CPU time, hence cost, CO<sub>2</sub>, etc."
+ "Remember you can play with the binwidth too. If you want a very accurate distribution you need a narrow binwidth, but then you'll also need high resolution (lots of stars) so lots of CPU time, hence cost, CO<sub>2</sub>, etc."
]
},
{