#367 closed defect - FPE (fixed)

crash in prt_comment.cpp

Reported by: peter Owned by: nobody
Priority: major Milestone: C17_branch
Component: output Version: trunk
Keywords: Cc:

Description

Nick Abel found that the attached sim crashes on an FP exception with the trunk@11275. Closer inspection revealed that the culprit is this line in prt_comment.cpp:

differ = fabs(1.-iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauIn()* rt.DoubleTau/iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauTot())*100.;

since iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauTot() is zero.

Replacing the "case B" command with "database H-like Lyman pumping off" seems to make the bug go away, so this command maybe the culprit. But then, this bug seems to be quite sensitive to details in the input script, so that may be a coincidence too.

The sim crashes at the end of the first grid point. It takes a few minutes to get there.

Attachments (1)

nick.tar.xz (206.3 KB) - added by peter 15 months ago.
input files

Download all attachments as: .zip

Change History (3)

Changed 15 months ago by peter

Attachment: nick.tar.xz added

input files

comment:1 Changed 15 months ago by peter

This looks like it is some late fallout from setting opac.taumin to zero. On the first iteration TauTot() is zero and opac.taumin used to protect against that.

One way of dealing with this would be to invert the calculation of differ:

differ = fabs(1.-iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauTot()/ (iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauIn()*rt.DoubleTau))*100.;

and similar for the Lya test further down. This would be unconditionally safe as TauIn() is guaranteed to be non-zero. It would slightly alter the test, but not in a significant way I would think.

PS - going over the code in rt_tau_init.cpp, I noted that there are still explicit constants 1e-20 present in the code. These should be replaced by opac.taumin.

comment:2 Changed 15 months ago by peter

Resolution: fixed
Status: newclosed

Fixed in r11292.

Note: See TracTickets for help on using tickets.