Opened 2 years ago

Closed 2 years ago

#331 closed defect - convergence (fixed)

grain assert thrown

Reported by: Gary J. Ferland Owned by: peter
Priority: blocker Milestone: C17_branch
Component: grains Version: trunk
Keywords: Cc:

Description

The following sim throws an assert at the head of the trunk but worked in c13_branch

laser linear 0.6
intensity  -20 range  0.4412 to 1 Ryd
c must provide source of ionization for chemistry 
cosmic rays background 
c
c commands for density & abundances =========
hden 0
grains ism
c
c commands controlling geometry  =========
stop zone 1
set dr 0 
c
constant temperature 2000 
c

The assert:

PROBLEM DISASTER
 An assert has been thrown, this is bad.
 Failed: Ehs[i] > 0. && Ehs[i] <= Ehi[i]
 It happened in the file ../grains.cpp at line number 3898
 This is iteration 1, nzone 0, fzone 0.00, lgSearch=T.

Change History (2)

comment:1 Changed 2 years ago by peter

This is an interesting bug with several facets. The main problem is due to cancellation errors in y2s(). This problem always existed since we implemented WDB06, but was dormant so far. The routine in question is very prone to cancellation errors and it has several series expansions in it already as a result of that...

The problem was flushed out by r10306 which changed the initial estimate for the grain charge. Combined with the fact that the laser command sets the intensity to a very low non-zero value outside the laser range, the code was tricked into believing that the radiation field was extremely hard and it assumed that the grains would be very highly charged as a result. This high charge then flushed out the cancellation errors.

I wonder if we need to set the intensity to a non-zero value outside of the laser. It sounds like a remnant of old code. Cloudy should be able to handle zero flux. So I propose to change line 279 of cont_ffun.cpp from

ffun1_v = SMALL;

to

ffun1_v = 0.;

comment:2 Changed 2 years ago by peter

Resolution: fixed
Status: newclosed

Both issues were fixed in r10481.

Note: See TracTickets for help on using tickets.