ASTERIA breaks when running with python 3.11

[jlazar17]

Something has changed with how Python 3.11 handles the EnumMeta, and this breaks an internal iteration with ASTERIA. I attach a working example based on this example that demonstrates this

import os

from astropy.table import Table
import astropy.units as u
import numpy as np

from snewpy import model_path
#from snewpy.neutrino import Flavor
#from snewpy.models.ccsn import Analytic3Species

from asteria.simulation import Simulation
#from asteria import set_rcparams

filename = "AnalyticFluenceExample.dat"
model_folder = f"{model_path}/AnalyticFluence/"

if not os.path.exists(model_folder):
    os.makedirs(model_folder)
file_path = os.path.join(model_folder, filename)

# These numbers _almost_ reproduce the Livermore model included in the SNOwGLoBES repository.
# They are obtained by calculating the total L, <E> and <E^2> from the livermore.dat
# fluence file (which is modelled after a 10kpc supernova).
total_energy = (5.478e+52, 5.485e+52, 4 * 5.55e+52)
mean_energy = (11.5081, 15.4678, 21.0690)
rms_or_pinch = "rms"
rms_energy = (12.8788, 17.8360, 24.3913)

# Make an astropy table with two times, 0s and 1s, with constant neutrino properties
table = Table()
table['TIME'] = np.linspace(0,1,2)
table['L_NU_E'] =  np.linspace(1,1,2)*total_energy[0]
table['L_NU_E_BAR'] = np.linspace(1,1,2)*total_energy[1]
table['L_NU_X'] = np.linspace(1,1,2)*total_energy[2]/4. #Note, L_NU_X is set to 1/4 of the total NU_X energy

table['E_NU_E'] = np.linspace(1,1,2)*mean_energy[0]
table['E_NU_E_BAR'] = np.linspace(1,1,2)*mean_energy[1]
table['E_NU_X'] = np.linspace(1,1,2)*mean_energy[2]

if rms_or_pinch == "rms":
    table['RMS_NU_E'] = np.linspace(1,1,2)*rms_energy[0]
    table['RMS_NU_E_BAR'] = np.linspace(1,1,2)*rms_energy[1]
    table['RMS_NU_X'] = np.linspace(1,1,2)*rms_energy[2]
    table['ALPHA_NU_E'] = (2.0 * table['E_NU_E'] ** 2 - table['RMS_NU_E'] ** 2) / (
        table['RMS_NU_E'] ** 2 - table['E_NU_E'] ** 2)
    table['ALPHA_NU_E_BAR'] = (2.0 * table['E_NU_E_BAR'] ** 2 - table['RMS_NU_E_BAR'] ** 2) / (
        table['RMS_NU_E_BAR'] ** 2 - table['E_NU_E_BAR'] ** 2)
    table['ALPHA_NU_X'] = (2.0 * table['E_NU_X'] ** 2 - table['RMS_NU_X'] ** 2) / (
        table['RMS_NU_X'] ** 2 - table['E_NU_X'] ** 2)
elif rms_or_pinch == "pinch":
    table['ALPHA_NU_E'] = np.linspace(1,1,2)*pinch_values[0]
    table['ALPHA_NU_E_BAR'] = np.linspace(1,1,2)*pinch_values[1]
    table['ALPHA_NU_X'] = np.linspace(1,1,2)*pinch_values[2]
    table['RMS_NU_E'] = np.sqrt((2.0 + table['ALPHA_NU_E'])/(1.0 + table['ALPHA_NU_E'])*table['E_NU_E']**2)
    table['RMS_NU_E_BAR'] =  np.sqrt((2.0 + table['ALPHA_NU_E_BAR'])/(1.0 + table['ALPHA_NU_E_BAR'])*table['E_NU_E_BAR']**2)
    table['RMS_NU_X'] = np.sqrt((2.0 + table['ALPHA_NU_X'])/(1.0 + table['ALPHA_NU_X'])*table['E_NU_X']**2 )
else:
    print("incorrect second moment method: rms or pinch")

table.write(file_path,format='ascii',overwrite=True)

# SNEWPY model dictionary, the format must match the below example for analytic mdoels
model = {
    'name': 'Analytic3Species',
    'param': {
        'filename': file_path
    }
}

sim = Simulation(model=model,
                 distance=10 * u.kpc,
                 Emin=0*u.MeV, Emax=100*u.MeV, dE=1*u.MeV,
                 tmin=-1*u.s,  tmax=10*u.s,    dt=1*u.ms,
                 mixing_scheme='AdiabaticMSW',
                 hierarchy='normal')
sim.run()

[jlazar17]

I have an extremely patchy suggestion for fixing this, which is adding a few lines to simulation.py. Specifically, I think that if you add a few checks at line 246 of the form

                if interaction.name=="requests":
                    continue
                if not self.interactions.requests.value[interaction.name]:
                    continue

it should work as expected with Python 3.11. I have not thought about nor tested how this interacts with other Python versions so I think this is a patch idea at best.

[jlazar17]

I realized I completely forgot to provide traceback. Here it is:

  File "/n/home12/jlazar/ASTERIA_py311_mwe/mwe.py", line 76, in <module>
    sim.run()
  File "/n/home12/jlazar/ASTERIA/python/asteria/simulation.py", line 232, in run
    self.compute_photon_spectra()
  File "/n/home12/jlazar/ASTERIA/python/asteria/simulation.py", line 248, in compute_photon_spectra
    xs = interaction.cross_section(flavor, self.energy).to(u.m ** 2).value
         ^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/n/home12/jlazar/ASTERIA/python/asteria/interactions.py", line 942, in cross_section
    return self.value.cross_section
           ^^^^^^^^^^^^^^^^^^^^^^^^
AttributeError: 'dict' object has no attribute 'cross_section'

[sgriswol]

Thanks for pointing this out, I'm sorry that this issue has become somewhat stale; I should be able to get to it this week. @jakob2508 has run into a similar problem.

The design principle behind using the Enum object was to have an easy-to-use, immutable iterable whose elements could be custom classes (Interaction objects). However, since two implementations for IBD exist, and the user may want to use either one or a different subset of interactions, we enabled interaction selection using EnumMeta. It achieved the desired functionality but made the interactions module obtuse and conflicted with the intended immutability of the enum module, but I digress.

I propose re-implementing the Interactions object as a simple container class. The interaction objects within can be selected at initialization but become immutable (or effectively so) afterward. If done correctly, this should only require changes on the backend, allowing user code to remain unmodified. See this branch for my initial attempt at a fix: https://github.com/icecube/ASTERIA/blob/interactions_rework/python/asteria/interactions.py
I'll revisit this soon, ETA for a fix is the end of the week.