WIP: wgtC resolution controls

[AstroBarker]

This is ongoing work to port the wgtC resolution controls in from nubhlight. I have added ComputeTotalEmissivity which (ideally) does just that, and SetWeight which uses that to set wgtC. Currently, ComputeTotalEmissivity is called in UpdateParticleResolution and SetWeight in UpdateTuningParameters'. However, it's crashing when trying to run the thermalization problem (i.e., leptoneq). Would appreciate eyes on this.

• [ ] Adds a test for any bugs fixed. Adds tests for new features.
• [ ] Format your changes by calling scripts/bash/format.sh.
• [ ] Explain what you did.

[Yurlungur]

@brryan might see things I don't here.

[AstroBarker]

Yeah, the crash is usually a root find failure, though issues start to pop up much earlier. Jtot is initially roughly constant at ~ 2e-3 (gotta do the math + units etc) but eventually blows up very fast to nan. Looks like ye somewhere becomes very large and negative -> the emissivity becoming unphysical

[AstroBarker]

Update: the behavior seems to persist even when I turn off the SetWeight bit and only use the other tuning controls (e.g., to update tune_emission). Perhaps there's something there instead?

[Yurlungur]

Update: the behavior seems to persist even when I turn off the SetWeight bit and only use the other tuning controls (e.g., to update tune_emission). Perhaps there's something there instead?

How much is tune_emission that the code arrives at, vs the one you hardcoded previously?

[brryan]

Yeah, the crash is usually a root find failure, though issues start to pop up much earlier. Jtot is initially roughly constant at ~ 2e-3 (gotta do the math + units etc) but eventually blows up very fast to nan. Looks like ye somewhere becomes very large and negative -> the emissivity becoming unphysical

I don't immediately see a problem in the code here but the Ye issue is interesting and I don't have good intuition for that. We don't currently apply any clamps on the value of Ye I don't think; that might need to be done along with the other floor/ceilings enforced in src/fixup/fixup.cpp:ApplyFloorsImpl

[Yurlungur]

I think Ye hitting a floor implies Ye supercooling though, so probably wgtC is too large. Actually what are you setting for the target photon number @AstroBarker ?

[brryan]

I think Ye hitting a floor implies Ye supercooling though, so probably wgtC is too large. Actually what are you setting for the target photon number @AstroBarker ?

I think that's right but even with good particle resolution supercooling may occasionally occur, and even one event might destroy the simulation, so a Ye floor may still be necessary

EDIT: did/do you have such a floor in nubhlight @Yurlungur?

[Yurlungur]

There's no Ye floor in nubhlight.

[Yurlungur]

That said for robustness it still might be a good idea.

[AstroBarker]

Update: the behavior seems to persist even when I turn off the SetWeight bit and only use the other tuning controls (e.g., to update tune_emission). Perhaps there's something there instead?

How much is tune_emission that the code arrives at, vs the one you hardcoded previously?

Previously it was simply 1.0. Now it steadily decreases orders of magnitude and continues till the end/crash.

[AstroBarker]

I think Ye hitting a floor implies Ye supercooling though, so probably wgtC is too large. Actually what are you setting for the target photon number @AstroBarker ?

num_particles is 1e5

[Yurlungur]

Update: the behavior seems to persist even when I turn off the SetWeight bit and only use the other tuning controls (e.g., to update tune_emission). Perhaps there's something there instead?

How much is tune_emission that the code arrives at, vs the one you hardcoded previously?

Previously it was simply 1.0. Now it steadily decreases orders of magnitude and continues till the end/crash.

Hmm... This suggests, I think that it's creating particles until it hits your target number then it's turning off. Does it hit the target number quickly? When running with fixed tune_emission, what number of particles did you see? Did it eventually reach a steady state? Or did it grow without bound?

I wonder if the resolution controls are correctly seeing particle losses due to absorption or outflow out the outer boundary? (I guess this is a periodic problem, so the later case isn't relevant...)

[Yurlungur]

Actually just to double check: Are you running with absorption or outflow boundary conditions? I.e., is there a way for neutrinos to be lost?

[AstroBarker]

Actually just to double check: Are you running with absorption or outflow boundary conditions? I.e., is there a way for neutrinos to be lost?

Absorption is on yes

[Yurlungur]

Possible that if absorption is reduced compared to emission, supercooling is happening early on before equilibrium is achieved. What if you change the time scale for the resolution controls? Maybe multiply it by a big factor, like 10x. Obviously that won't test the controls as quickly. But the controls can be sensitive and problem dependent... so it's possible everything is implemented correctly, but things are going haywire because the controls aren't tuned quite right.

[AstroBarker]

Possible that if absorption is reduced compared to emission, supercooling is happening early on before equilibrium is achieved. What if you change the time scale for the resolution controls? Maybe multiply it by a big factor, like 10x. Obviously that won't test the controls as quickly. But the controls can be sensitive and problem dependent... so it's possible everything is implemented correctly, but things are going haywire because the controls aren't tuned quite right.

Tried t_tune_emission x10, x100 and same sorts of behavior. Ye seems to start approaching the right value and then starts falling again

[AstroBarker]

Seems to be working as intended now? One issue was setting num_particles = 1e5 wasn't okay as it expected an integer and was just set to 1. I also added a check if (correction < 0.0) correction = 4. / 3. to mirror the behavior in nubhlight. This did (positively) affect the results. Next is to queue up a collapsar.

[Yurlungur]

Sounds good. Assuming tests pass, perhaps we should merge into develop?