Hyper terms including in the implicit response but not linearly independent

[David Dickinson]

In order to calculate the response we call the regular time advance after perturbing the fields. This will call hyper_diff which will apply hyper terms based on the current fields. The hyper operators as implemented typically couple together all wavenumbers. This means that separate supercells are no longer independent in simulations including hyper coefficients > 0. As such, the exact response calculated depends on if we are perturbing fields in different supercells in parallel or treating each supercell serially. Only with the fully serial approach would we avoid seeing any contribution to the response coming from the hyper terms

Currently both fields_implicit and fields_local treat supercells in parallel, although the parallel pattern is different in both cases. As a result of this we do get a contribution to the response from hyper terms and this contribution is different between the two field options.

‌

For a simple test case with d_hyper=1 the maximum difference in response matrix was found to be ~ 0.04 whilst the maximum response was ~0.8. The difference for fixed field method with and without hyper terms was ~ 0.052 for local and 0.009 for implicit. This perhaps suggests that the hyper terms can lead to a significant modification of the response matrix.

‌

I think we likely want to exclude the hyper terms from the calculation of the response matrix as I believe this approach relies on the response in a supercell being independent of all other supercells, which is not currently the case when hyper terms are included.

[David Dickinson]

Here’s an example comparison of the response for the two field methods with and without hyper terms.

Note this data is the natural log of the difference. The white spots on the right indicate where the difference is 0 (as log(0) maps to NaN).

‌

[David Dickinson]

Here is a comparison of a small nonlinear simulation with and without hyper terms with fields implicit and fields local. I’ve also included a case where I include hyper terms only for istep>0 (i.e. excluded from the response calculation but included in the real time advance)

It can be seen that the two field approaches agree very well for the cases without hyper and with hyper in the advance only, whilst when the hyper terms are included in the response calculation the two field approaches differ more significantly (although the simulation is not run for long enough, or at high enough quality, to state if they really agree statistically or not).

‌

[David Dickinson]

Exclude hyper terms when istep <= 0

This removes these terms from the response calculation. See issue #133 for some discussion.

→ <<cset 05f7c62f2cc3>>

[David Dickinson]

See branch experimental/try_excluding_hyper_during_response_calculation for the “hyper advance only” code.

[David Dickinson]

Here is a comparison of the heat flux vs time for local with hyper as standard/in advance only over a longer time period. By eye these look in good agreement.

‌

[David Dickinson]

PR #528 proposes to bring in the change to enable hyper in advance only.

This may not be a true “fix” rather something slightly wrong in a different way, but in a way that means the run is reproducible on all core counts and with both field options.

As pointed out by @Joseph Parker we might be able to instead add the hyper filter into the explicit source terms and evolve it as a part of the explicit terms. This would, however, require us to ensure the explicit history arrays (gexp_* ) are allocated when either hyper or nonlinear terms are active (so that we can test the impact of hyper in linear runs). We may then need to consider how/if the cfl estimate needs to be updated.

Let's suppose that we've implemented hyper_filter as being gnew = a*gnew with a constant. This means the stored response becomes a*response_no_hyper as a matrix. We then invert the response matrix to get the useful version - useful_response = (1/a) response^{-1} (where response == response_no_hyper). In the time advance we then have gnew = a*gnew again, so the field equation becomes something like fieldeq = a * antot - field_terms so we get dphi = response^{-1}antot - response^{-1} field_terms/a. If we exclude from the response calculation so that useful_response = response^{-1} then we'll end up with dphi = response^{-1}( a * antot) - response^{-1} field_terms . We could recover the original hyper result by dividing dphi by a but this doesn’t really gain us anything beyond ensuring the response matrix is independent of processor count and field option.