Message Boards Message Boards

Evolution of Stars: Hertzsprung-Russell diagram

GROUPS:

Hi,

I just made a little Hertzsprung-Russell diagram for a class of mine. It's not quite perfect yet, but I wanted to post it anyway.

Mathematica (wikipedia?) knows what that diagram is:

WikipediaSearch["Hertzsprung Russell diagram"]

gives:

<|"Title" -> "Hertzsprung-Russell diagram",

"Snippet" -> "The Hertzsprung-Russell diagram, abbreviated H-R diagram or HRD, is a scatter graph of stars showing the relationship between the stars' absolute magnitudes or luminosities versus their spectral classifications or effective temperatures."|

All the data I need is in Mathematica.

stars = StarData[];

contains data for more than 107,000 stars. I will only download data for 10,000 stars as to not overload my final image.

starChoice = RandomChoice[stars, 10000]; 
starData = StarData[starChoice, {"Luminosity", "EffectiveTemperature", "Mass", "Radius", "BVColorIndex", "Color"}];

I access luminosity, surface temperature, mass, radius, and color data for the stars. I will want to use BubbleChart to generate the diagram. The position of the bubbles will be given by (the logarithms of) the luminosity and the surface temperature; note that the x-asis is reversed. The size of the bubbles will encode the radius and I would like to color the bubbles in a way that relates to their BVColorIndex. The most difficult issue is the colour. As a matter of fact, I will cheat a bit to make the colour distribution clearer. I use some formulas from this website for the conversion of BV to RGB colour.

{\[Lambda], x, y, z} = Import["http://www.cvrl.org/database/data/cmfs/ciexyzjv.csv"]//Transpose;
XYZ[t_] := Module[{h = 6.62607*10^-34, c = 2.998*10^8, k = 1.38065*10^-23}, {x, y, z}.((2 h c^2)/((-1 + E^((h c/k)/(t \[Lambda]*10^-9))) (\[Lambda]*10^-9)^5)) // #/#[[2]] &]

The plot is quite easy to generate from here:

normsun = QuantityMagnitude[StarData["Sun", "Luminosity"]];
BubbleChart[
 Style[{Log[10, QuantityMagnitude[#[[2]]]], QuantityMagnitude[#[[1]]]/normsun, QuantityMagnitude[#[[4]]]}, 
    ColorConvert[XYZColor @@ XYZ[4600*((1/(0.92 #[[5]] + 1.7)^2.4 + 1/(0.92 #[[5]] + 0.62)^2.4))], "RGB"]] & /@ 
  Select[starData, FreeQ[#, Missing] &], ScalingFunctions -> {"Reverse", "Log"}, Background -> Black, 
 FrameLabel -> {"Log(Temperature)", "Log(Luminosity)"}, LabelStyle -> Directive[Bold, 15, White], ImageSize -> Full]

Note that there is a "fudge exponent" of 2.4 which exaggerates the colours. This gives:

enter image description here

If you compare this figure to the corresponding images on the wikipedia page, you will note that is looks quite similar to the second figure. The first figure on that website, shows a fine structure of all the sequences Ia and Ib for the supergiants, II for the bright giants, III for giants, IV for subgiants, and V for the main sequence. Then there are the white dwarfs. It appears that our diagram shows a similar sub-structure, but I might be mistaken.

Our own sun is on the main sequence at the intersection at a horizontal line of normalised Log-luminosity 1. (The axis is normalised to our sun's luminosity.) So our sun is "yellow-white".

As we are already on it, we can also compare where our sun lies in terms of its mass. I therefore compute a SmoothHistogram of the masses of my 10000 stars.

SmoothHistogram[
 Select[QuantityMagnitude[starData[[All, 3]]], NumberQ], 
 PlotRange -> {{0, 1.5*10^31}, All}, PlotTheme -> "Marketing", 
 FrameLabel -> {"Mass", "Density"}, ImageSize -> Full, 
 LabelStyle -> Directive[Bold, Medium], 
 Epilog -> {Red, {PointSize[0.02], 
    Point[{QuantityMagnitude[StarData["Sun", "Mass"]], 0}]}, Green, 
   Line[{{Median[
       Select[QuantityMagnitude[starData[[All, 3]]], 
        NumberQ]], -0.1}, {Median[
       Select[QuantityMagnitude[starData[[All, 3]]], NumberQ]], 
      7.5*10^-31}}], Yellow, 
   Line[{{Mean[
       Select[QuantityMagnitude[starData[[All, 3]]], 
        NumberQ]], -0.1}, {Mean[
       Select[QuantityMagnitude[starData[[All, 3]]], NumberQ]], 
      7.5*10^-31}}]}]

enter image description here

Our sun's mass is marked by the red dot on the x-axis. The green vertical line shows the median of the distribution and the yellow line the mean. I have truncated the diagram at the large masses. The following plot shows mass vs radius:

ListPlot[QuantityMagnitude /@ Select[starData[[All, {3, 4}]], FreeQ[#, Missing] &], PlotTheme -> "Marketing", LabelStyle -> Directive[Bold, 15, Black],FrameLabel -> {"Mass", "Radius"}]

enter image description here

Cheers,

Marco

POSTED BY: Marco Thiel
Answer
2 years ago

enter image description here - another post of yours has been selected for the Staff Picks group, congratulations !

We are happy to see you at the tops of the "Featured Contributor" board. Thank you for your wonderful contributions, and please keep them coming!

POSTED BY: Moderation Team
Answer
1 year ago

Thanks for this post Marco. Unfortunately I was unable to reproduce your beautiful first plot. I get the following error: The front end encountered an error while processing a "NotebookPredictions" packet. I'm running Mathematica 11.0.1 for Linux x86 (64-bit) (September 21, 2016).

One other thing that you may want to mention in your post (or maybe think about lowering the number of stars?): it took over 6 hours to download the data for the 10000 stars (60000 data items). It was not due to a slow connection on my end. I wonder if Wolfram purposely slows down data for large requests or if the request is coming from someone (like me) who only has the Home Edition of Mathematica? The data started coming down fairly fast but got slower and slower. By the end the "downloading" box was counting 1 data item about every 3 seconds!

Thanks again, Joe

POSTED BY: Joe Gilray
Answer
6 months ago

Dear Joe,

I am sorry to hear that you are having trouble downloading the data. I have just run the code again (on MMA 11.1.1) and do not have any indication that anything is going wrong. I do not get an error when I run the code. (There was however, one typo, I just noticed. It was \Transpose in one line instead of //Transpose; fixed now.)

enter image description here

And you are right that downloading can take some time, but it is far below 6 hours for me. I think that it depends a bit on where you are and which server is used to provide the data to you. You are in Oregon (?), and I am not sure which server Wolfram use for you. I am in Scotland and I have to wait a little bit, but after 13 minutes I am done downloading. Here is a screenshot of the results:

enter image description here

I am not sure how to help with the server issue. Have you tried using a proxy-server to another location? The amount of data transferred is rather low.

starData // ByteCount

is just under 10 MB. So, if you get a better server, that might help. Another thing is that there might be a bit of load balancing going on. I am about to post another thing where I use a lot of data from Wolfram servers - not 10000 but rather 100000s of data points. I have noticed that downloading them in smaller batches helps (and probably eats up more points for API requests).

Best wishes,

Marco

POSTED BY: Marco Thiel
Answer
6 months ago

Hi Marco,

Thanks for your thoughts. I fixed the typo and got your code to run on 1000 stars... very nice. Thanks for a thought-provoking article!

Regards, -Joe

POSTED BY: Joe Gilray
Answer
6 months ago

Group Abstract Group Abstract