cooling sim convergence

[Marios Chatzikos]

Hello,

I have been looking at cooling sims that set the time step as a fraction
of gas cooling time.  By default that fraction is taken to be 4%, but
I've been looking at whether converged results are obtained by using 1%
fractions.

The quantities I have used as a metric for convergence are line emission
integrals (the usual Γ factors), defined in the context of the cooling
flow model.  Many Γ's are discrepant between these two runs.  Oddly,
some of them are from ions that are abundant at temperatures above the
peak of the cooling function (10^5K), where the gas is essentially in
equilibrium.  The poster child of this is O VII λ22.1012A, a transition
that occurs between the two lowest levels in the ion (1^1S - 2^3S).  See
first figure.

After some investigation, it appears that the problem comes from not
using level-resolved advective terms in the linear system, but rather
splitting the total creation rate for the ion among levels using the
partition function.  This guarantees that the correct limit is obtained
at LTE, but it produces incorrect results at low temperatures.

The second plot illustrates the evolution of the level 2 population
under different integration schemes.  Those labeled 'level' are the
resolved runs, while those labeled 'EQ' have their populations forced to
equilibrium (no advective terms are used in the linear system).  The
most discrepant curves are obtained with the partitioned source terms. 
The rest are in very close agreement with each other.

I think this proves that use of level-resolved advection terms produces
the desired result at low densities and high temperatures.  I am now
looking at the behavior of the system at low densities and temperatures
(below 10^5K).  Later I'll look if this formalism recovers the LTE
limit.  I'll let you know when I have something worth sharing.

Let me know what you think.

Thanks,

Marios

[Robin Williams]

Hi Marios --

The obvious interpretation of this -- /assuming/ that there's no bug in the level-resolved physics -- is that non-equilibrium excitation of the levels is important.  Is there a suitable form of critical density for the system which concurs with this?

What is the relative cost of the calculations with levels set by LTE vs. integrated explicitly?  I think the reason why the partition function version is the default was the concern that always doing the detailed calculation would be too expensive.

Yes, checking everything agrees in LTE would obviously be very useful...

Many thanks

  Robin 

--
--
http://groups.google.com/group/cloudy-dev
---
You received this message because you are subscribed to the Google Groups "cloudy-dev" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cloudy-dev+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/cloudy-dev/c64ca96a-33e2-6900-252b-6d8df1914f81%40gmail.com.

[Gary J. Ferland]

Hi Robin,

On Wed, Feb 17, 2021 at 1:45 PM 'Robin Williams' via cloudy-dev <cloud...@googlegroups.com> wrote:

Hi Marios --

The obvious interpretation of this -- /assuming/ that there's no bug in the level-resolved physics -- is that non-equilibrium excitation of the levels is important.  Is there a suitable form of critical density for the system which concurs with this?

For atoms, the levels that do not have E1 decay routes are highly metastable.  But that is not true of such highly charged systems - non-E1 transitions are surprisingly fast - the transition rate goes up as a power of nuclear charge Z.  The environment has a very low density so the radiative decays should be the fastest way out of a level.  The rates levels are populated are far slower - so source terms have to be compared with the recombination rate.  Here I believe the advective terms dominate.

thanks,

Gary

To view this discussion on the web visit https://groups.google.com/d/msgid/cloudy-dev/CAJBtPxe_Qq5tanFdOWXxH3VPkyv0xKeuarj2BfGBkYGhYBsy6A%40mail.gmail.com.

--

Gary J. Ferland
Physics, Univ of Kentucky
Lexington KY 40506 USA
Tel: 859 257-8795
https://pa.as.uky.edu/users/gary

[Marios Chatzikos]

Hi Robin,

My focus was on the impact of the dynamical terms on the matrix system, so I did not look for an error in the level-resolved physics of the iso-solver.  If there is one, it may not show up in the plots I shared, as they simply illustrate that the equilibrium solution is recovered correctly, for this ion.  If there was a problem with the level-resolved dynamical terms, I don't think the equilibrium solution would have been recovered.

The critical density is defined for a two-level system, and does not account for advective terms.  I don't think it can easily be extended to many-level systems with advective terms, but if you can recall a paper that discusses a generalized critical density (without advection), please let me know.  Otherwise, one could be sought for empirically, but I'm not sure how much value there is in that.

I did not plan on looking into the efficiency of the code until after I arrived at the correct solution to the problem.  But since you ask, a cooling calculation from 30 MK to 3 K, using timesteps of 1% of the cooling time took 1682 iterations to complete with both the resolved and partitioned runs, and, respectively, 2067 and 2159 seconds.  This amounts to a ~4% speedup in this sim, FWIW.  In any case, this sim was not negatively impacted by the level-resolved physics.

Thanks,

Marios

PS: I attach plots of the instantaneous and accumulated emission in O VII 22.1012A including the two level-resolved runs.  Note that there is a 1.5% offset between the two Γ's that will look into later.

To view this discussion on the web visit https://groups.google.com/d/msgid/cloudy-dev/CAJBtPxe_Qq5tanFdOWXxH3VPkyv0xKeuarj2BfGBkYGhYBsy6A%40mail.gmail.com.

[Marios Chatzikos]

Hello Gary and Robin,

Going in detail over the linear system to be solved, what vary between the equilibrium and dynamical runs are the matrix rate and source rates of the ground state, that is, elements z[0][0] and creation[0], respectively.

At the temperature where the level 2 population peaks, the dynamical terms do indeed dominate over recombination by a factor of about 5.  But, as Gary said, the dominant rate for the ground state (in element z[0][0]) are the regular collisional and radiative rates, which contribute about 60% of the total (with advective 34%, recombination 6%).

By contrast, the creation rate of the ground state is boosted by a factor of ~22 between the equilibrium and dynamical (state-resolved) runs.  These rates are also boosted in the 'partitioned' run, but only by a factor of ~5 in the creation[0] rate.  It seems to me that it is the ground state creation rate that drives the system to the equilibrium solution.

Thanks,

Marios

To view this discussion on the web visit https://groups.google.com/d/msgid/cloudy-dev/CAOxTjFQpAqKOxiRXzzMxqBukQSf9HNBXEtC_FbUV-O2B2bKD7Q%40mail.gmail.com.

[Marios Chatzikos]

Hello again,

In relation to the offset in the accumulated emission of the O VII line mentioned before: I extended the check on Γ's between the 1% and 4% runs to all lines in the sim, and it looks like that discrepancies are now within 3% for iso-sequence lines, and to within 10% for database lines.

It may be useful to extend the use of level-specific advective terms to database species, as well.

Marios

[Robin Williams]

Hi Marios --

The advection terms aren't really physical rates, but a means of representing varying state in the static balance for which Cloudy solves.  I think an accountant might call this an https://en.wikipedia.org/wiki/Accrual

If they dominate, then some other process changing slowly is driving this.  In some cases it may be metastable levels.  From what Gary says, here it sounds more like frozen ionization is leading to the populations in higher levels being controlled by the recombination cascade rather than collisional excitation.  (FWIW, long-timescale under-excitation can also happen in flows which are heating, e.g. the Balmer emission from strong shocks and around supernovae.)

It might be possible to come up with a dynamic equilibrium model of how to distribute the advective sources better than thermal equilibrium balance, based e.g. on how the relative populations of the neighbouring ions compares to LTE.  But in any model there are modelling choices to validate (for which it becomes harder to know if you've missed some crucial detail in some limit, as the situations get more complex), while tracking the levels explicitly is more like "ground truth".

Yes, it would be worth at least having the option of including advective terms in database species.  This will increase the memory footprint of the runs, but will often make it easier to solve the level balance.  It may be possible to store the previous state more economically than the present scheme.

  Robin

P.S. I (re-)sent you some old talk slides on critical densities, which are somewhat relevant to this discussion.

To view this discussion on the web visit https://groups.google.com/d/msgid/cloudy-dev/adae68b3-ccd5-146d-5767-cbc3f1a49f7f%40gmail.com.