Added neutrinos and Tcmb0 read from snapshots

[JBorrow]

Fixes #151.

Can @JegerBroxterman please test this and make sure that it is consistent with their expectations?

[JegerBroxterman]

Fixes #151.

Can @JegerBroxterman please test this and make sure that it is consistent with their expectations?

I don't think it works yet. Although I am not entirely sure I did the correct check. I replaced the old configurations.py file with the adjusted one from above. I did not change anything else. For both the old and new files the T_CMB remains off (2.9966 K). The old version gave neutrino masses = [0. 0. 0.] eV, whereas the new one gives [0.06 0.06 0.06] eV. I believe the correct values for the FLAMINGO runs should be [0.06, 0., 0.] eV.

[JBorrow]

Can you check if the CMB temp is being recalculated by putting a print statement in the except: block?

[JegerBroxterman]

The text in the print statement gets printed.

[JBorrow]

Can you show me the cosmology raw?

[JegerBroxterman]

{'Cosmological run': array([1], dtype=int32), 
'Critical density [internal units]': array([12.87106552]), 
'H [internal units]': array([68.09999997]),
 'H0 [internal units]': array([68.09999997]), 
'Hubble time [internal units]': array([0.01468429]),
 'Lookback time [internal units]': array([9.02056208e-16]),
 'M_nu_eV': array([0.06]), 
'N_eff': array([3.04400163]),
 'N_nu': array([1.]), 
'N_ur': array([2.0308]),
 'Omega_b': array([0.0486]),
 'Omega_cdm': array([0.256011]), 
'Omega_g': array([5.33243487e-05]),
 'Omega_k': array([2.5212783e-09]),
 'Omega_lambda': array([0.693922]),
 'Omega_m': array([0.304611]),
 'Omega_nu': array([0.00138908]), 
'Omega_nu_0': array([0.00138908]), 
'Omega_r': array([7.79180471e-05]), 
'Omega_ur': array([2.45936984e-05]),
 'Redshift': array([0.]),
 'Scale-factor': array([1.]), 
'T_CMB_0 [K]': array([2.7255]), 
'T_CMB_0 [internal units]': array([2.7255]), 
'T_nu_0 [eV]': array([0.00016819]), 
'T_nu_0 [internal units]': array([1.9517578]), 
'Universe age [internal units]': array([0.01407376]), 
'a_beg': array([0.03125]),
 'a_end': array([1.]), 
'deg_nu': array([1.]), 
'deg_nu_tot': array([1.]), 
'h': array([0.681]), 
'time_beg [internal units]': array([9.66296122e-05]), 
'time_end [internal units]': array([0.01407376]),
 'w': array([-1.]), 
'w_0': array([-1.]), 
'w_a': array([0.])}

Still gives the correct values, as was the case before.

[JBorrow]

Ok, that should be fixed now. @MatthieuSchaller can you comment on what the neutrino masses should be? AFAIK they should be uniform, as you model all neutrino species as one particle type.

[JegerBroxterman]

The T_CMB is correct now. As I understand it, FLAMINGO contains one massive species (with mass 0.06eV) and two massless species but @MatthieuSchaller knows this better.

[MatthieuSchaller]

@wullm, could you advise here what is the best interpretation of the neutrinos? And, maybe also what we ought to do in other more complex cases where we have >1 species?

[MatthieuSchaller]

In the flamingo runs, we have one massive nu of 0.06 eV and 2.0308 massless neutrinos. How astropy would like to register that information, I don't know.

[wullm]

Astropy assumes that the (integer) number of species is floor(N_eff), so we need to pass an array of the form [0.06, 0., 0.] * u.eV. Otherwise, it assumes that all neutrinos have mass 0.06 eV.

[JBorrow]

The code should really output an array of masses then as this is really ambiguous at the moment. I'll change this on Monday in swiftsimio.

[wullm]

It does, but only the massive species. The number of massless species is given by the N_ur attribute.

[wullm]

Note that there is also the possibility of degenerate species. For instance, one species with mass 0.05 and degeneracy 2. Presumably, this should be read by astropy as [0.05, 0.05] * u.eV.

[JBorrow]

I mean the mass array should read [0.06, 0.0, 0.0]. This would be clearest, and then also naturally includes degeneracy. For instance the integers we currently store couldn't represent [0.06, 0.02, 0.02, 0.02, 0.0, 0.0] as far as I understand it.

[wullm]

I'm just thinking of non-standard cases, in which a long list is not ideal. For instance, I have simulations with 30 degenerate species. (This is to get a large Omega_nu, but with realistic free streaming lengths.) Moreover, the degeneracies can be non-integer. Same for the massless species. But it's definitely worth spelling things out more clearly.

[wullm]

Yes, [0.06, 0.02, 0.02, 0.02, 0.0, 0.0] is possible. That could be m_nu = [0.06, 0.02], deg_nu = [1, 3], N_ur = 2.

[MatthieuSchaller]

@JegerBroxterman could you also test this on the cosmology variation boxes? Some of them use degeneracies (c.f. flamingo slack).

[JegerBroxterman]

With the current version, the neutrino values in metadata.cosmologyare still off. For Flamingo HYDRO_FIDUCIAL it still gives gives m_nu=[0.06 0.06 0.06] eV whereas cosmology_raw gives 'M_nu_eV': array([0.06]), 'N_eff': array([3.04400163]), 'N_nu': array([1.]), 'N_ur': array([2.0308]), i.e. 1 massive species with mass 0.06eV.

[wullm]

Thanks Jeger, I will propose a fix.

[JBorrow]

Thanks Willem. Can this please be tested?

[MatthieuSchaller]

Thanks both!

[MatthieuSchaller]

@JegerBroxterman could you check whether this now all works?

[JegerBroxterman]

Yes, all looks good now.

[MatthieuSchaller]

It's probably worth making a tag release and push to pip.

[JBorrow]

Thanks all. We need to keep track of the changes to the VR library too.