chianti-atomic / ChiantiPy

ChiantiPy is a python package to calculate the radiative properties of astrophysical plasmas based on the CHIANTI atomic database
63 stars 32 forks source link

NaNs generated when calculating spectra under certain conditions #435

Closed jslavin closed 1 year ago

jslavin commented 1 year ago

I've run into a situation in which NaNs are generated within Continuum.py. It's strange because it happens when I use a certain temperature range, but when I do essentially the same calculation but broken up into two temperature ranges, I don't get the NaNs. Basic setup to get the NaNs:

import ChiantiPy.core as ch

wvl = np.linspace(5.,170.,700) T = np.logspace(4.5, 6., 16) density = 1. spec = ch.spectrum(T, density, wvl, minAbund=1.E-6)

Result: /export/slavin/anaconda3/envs/py38/lib/python3.8/site-packages/ChiantiPy/core/Continuum.py:1158: RuntimeWarning: invalid value encountered in multiply fbn[0,itemp] = (const.planckconst.light/(1.e-8wvl))**5 \ fb has some Nans 397

The spectra have NaNs in them for certain wavelengths. Doesn't seem to depend on the setting of minAbund - happens also for minAbund = 2.E-5. Also it oddly does not happen when the T range is smaller, e.g. T = np.logspace(4.5,5.3, 9). So I can do the calculations in two segments and get no NaNs but if the T range is over the whole range I get NaNs. I'm using ChiantPy version 0.15.0

kdere commented 1 year ago

Hi Jonathan,

I was not able to reproduce this, although I did use minAbund=1.e-4

the error message that you get is a bit curious:

fbn[0,itemp] = (const.planckconst.light/(1.e-8wvl))**5

there should be a * between 1.e-8 and wvl. maybe your code is corrupted

regards, Ken

jslavin commented 1 year ago

That missing * was just a typo by me, It's there in the code. I tried it with minAbund = 1.E-4 and didn't get the NaNs. Could you try it with minAbund = 2.E-5 and see what happens?

jslavin commented 1 year ago

This seems to derive from expf[0, itemp] becoming infinite. Would it work to just set fbn[0,itemp] to 0 where expf[0,itemp] is infinite? It seems these are missed by the masking that's currently done.

kdere commented 1 year ago

I have tracked the problem down to si_14

Your suggestion might work but I can't guarantee it

kdere commented 1 year ago

it looks to me like the problems occurs at your lowest temperature only

this definitely needs a fix

kdere commented 1 year ago

the freeBound method was producing invalid numbers and a warning was given. I had to get the code to raise and exception and then deal with the problem. The masked arrays did not help at all. I expect to have a revised method of github soon and a version release shortly after

jslavin commented 1 year ago

Great, thanks Ken

On Tue, May 2, 2023 at 9:40 AM Ken Dere @.***> wrote:

the freeBound method was producing invalid numbers and a warning was given. I had to get the code to raise and exception and then deal with the problem. The masked arrays did not help at all. I expect to have a revised method of github soon and a version release shortly after

— Reply to this email directly, view it on GitHub https://github.com/chianti-atomic/ChiantiPy/issues/435#issuecomment-1531499122, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAE5OY3EPXFDSSFUEJ6REYLXEEFF3ANCNFSM6AAAAAAXLPWVLE . You are receiving this because you authored the thread.Message ID: @.***>

kdere commented 1 year ago

a corrected version of continuum.freeBound is now in the github repository. Next I will do a release.

kdere commented 1 year ago

version 0.15.1 has been released