Opened 8 years ago

Last modified 12 months ago

#102 reopened defect - code aborts

dynamics sim that crashes

Reported by: Gary J. Ferland Owned by: nobody
Priority: major Milestone: C19_branch
Component: dynamics Version: trunk
Keywords: Cc:

Description (last modified by peter)

In r3332 I committed the sim experimental/dynamics_orion_recom.in: a fast recombination front. It crashes in the 3rd iteration.

Remove this sim and document the problems on the known problems page if we can't fix in time for the c10 branch.

Change History (7)

comment:1 Changed 8 years ago by Gary J. Ferland

Version: C08.00trunk

comment:2 Changed 8 years ago by peter

Component: pressure convergencethermal convergence
Description: modified (diff)
Resolution: duplicate
Status: newclosed

This PR is a duplicate of #67 (after the initial fix in r2045). The pressure solver is not involved, this is a failure of the temperature solver. After the heating and cooling curves detach at some point, the solver fails to find a new temperature equilibrium and eventually hits the temperature minimum of the code.

comment:3 Changed 8 years ago by peter

Component: thermal convergencedynamics
Description: modified (diff)
Resolution: duplicate
Status: closedreopened
Summary: two sims that crashdynamics sim that crashes
Type: defect - convergencedefect - code aborts

The previous comment that this PR is a duplicate of #67 is only partially true. The fix r3349 for PR #67 fixed slow/self_gravity.in, but not experimental/dynamics_orion_recom.in.

Removing comments about slow/self_gravity.in from the original report as they are incorrect.

Closer analysis reveals that the problem in experimental/dynamics_orion_recom.in is completely different. At the point where the code crashes there actually is no solution for the temperature. Over the entire temperature range 2.8 - 1e7 K the cooling is larger than the heating. It seems very unlikely that there is a solution at higher temperatures as the previous zone has a temperature of 7631 K. This appears to be a problem in the dynamics code.

comment:4 Changed 8 years ago by peter

Ryan commented the following:

"I think the problem might be missing gravitational terms in phycon.EnthalpyDensity at lines 832 and 843 of pressure_total.cpp. Do we need to add phycon.IntegRhoGravity in as well? And is this also missing from !TotalPressure_v in the same file?"

comment:5 Changed 5 years ago by Ryan Porter

Milestone: C10 branchc13 release

Milestone C10 branch deleted

comment:6 Changed 12 months ago by mchatzikos

Milestone: c13 releaseC18_branch

comment:7 Changed 12 months ago by mchatzikos

The code executes a number of relaxation iterations (currently 3, from 2 when this ticket was first posted), and then attempts to advance the solution forward in time. The failure occurs because the cooling misbehaves at the outermost zone in the first time-dependent iteration (currently 4, from 3). In particular, the following are the partial contents of the 'save cooling' output in the zones leading up to the failure:

2.38650e+17     2.0019e+02      6.7859e-21      1.7580e-19      adve 0.0        2195.7551270    O  1 0.0        0.7054235       C  2 0.0        0.2456607
2.38653e+17     1.5170e+02      6.8068e-21      1.2090e-19      adve 0.0        3192.4633789    O  1 0.0        0.6308184       C  2 0.0        0.3311111
2.38655e+17     1.0797e+02      6.8545e-21      7.2165e-20      adve 0.0        5348.1459961    O  1 0.0        0.4943877       C  2 0.0        0.4833097
2.38657e+17     6.6899e+01      6.9344e-21      3.3057e-20      adve 0.0        11674.2314453   C  2 0.0        0.7513695       O  1 0.0        0.2433844
2.38660e+17     3.7392e+01      7.0947e-21      1.0928e-20      adve 0.0        35312.9375000   C  2 0.0        0.9566543       O  1 0.0        0.0428150
2.38662e+17     5.7704e+00      7.5982e-21      1.9462e-24      adve 0.0        198268048.000   O  1 0.0        0.5644801       hvFB 0.0        0.3244929
#>>>>  Temperature not converged.
2.38665e+17     2.8000e+00      7.0000e-21      1.5992e-24      adve 0.0        241287520.000   O  1 0.0        0.6914103       hvFB 0.0        0.2645316

That is, as the outer region is approached, the cooling remains higher than the heating. It is primarily due to Oxygen, but switches to C+ near the edge of the calculation. When trying to decrease the temperature from ~40K, the C+ cooling vanishes and the cooling becomes 3 dex lower than the heating. The odd behavior of the total and C+ cooling is the first problem.

The second problem is that the temperature solver does not attempt to increase the temperature to bring heating and cooling to closer agreement. Instead, it further decreases the temperature from ~6K to the CMB temperature. This would suggest that the derivative of the total cooling with temperature has the wrong sign.

For the record, the outermost zone in the final relaxation iteration has a depth of ~3e17, so the previous solution is not exceeded when the failure happens. The sim was run with the additional commands

stop temperature off
stop thickness 17.5

to guard against that possibility.

Note: See TracTickets for help on using tickets.