Execute eddy parameterization only once per timestep?

[milankl]

@swilliamson7 One other thing I'd try is to not add your S_u, S_v to momentum_u, momentum_v (which are executed on every RK substep) but to here and also only evaluate the Zanna Bolton forcing once with the diffusion operator.

https://github.com/swilliamson7/ShallowWaters.jl/blob/be72c43fb59e227fc5784e9de321b6aeced6531e/src/diffusion.jl#L223-L226

[swilliamson7]

I tried this here: https://github.com/swilliamson7/ShallowWaters.jl/blob/zb2/src/diffusion.jl, unfortunately it seemed to produce the same output as before when I added the forcing with the momentum equations (assuming I successfully placed S_u and S_v in the correct locations)

[milankl]

That looks correct to me (assuming you do you have compensated=false which should be default). But can you clarify what you mean by "same output". Like it should not be bitwise identical, but, e.g. in terms of energy levels reached, do you get the same? If things look the same can you ramp up your kappa and see what it does in a case where it might be unrealistically strong?

[swilliamson7]

It is indeed the default! And yeah, I don't mean identically the same numerically, more that on average they have the same value. I'm attaching a picture I just made with the energy spinups (I need to learn to create better diagnostic plots than this but haven't written that code yet). But yeah increasing the magnitude of the coefficient decreases the energy by quite a bit :/

[milankl]

Okay this looks like a good metric to use to assess what's happening. You could even cut the length of the integration in half.

Good to see that red and pink are initially very similar before branching off. I'd therefore suggest that you leave it in the diffusion time step as it's faster and there's no need to stability reasons to do so. (Rather the opposite actually.)

So just to double check,

a) half the coefficient would be between green and red?

b) flipping the sign of the coefficient also doesn't increase the energy levels? Could you try half the negative so like *-1/2?

I'm asking that because I'm wondering whether there's a ~ linear relationship between the coefficient and the energy levels which I'd expect. If that's not the case there's either still something wrong in the implementation (which I doubt by now) or there's some more complicated interaction we don't understand yet. For example, it could be that whatever the backscatter is trying to do is compensated/overcompensated by the diffusion term.

Strategies I suggest: try to find the coefficient at the edge of stability. Play around with it and see where things blow up and then reduce it by a factor of two or so and then we can see what it's actually doing (spectra for example or just the term itself).

But this energy metric is good. Can you post whatever you try in term of that graph here?

[swilliamson7]

Commenting just to apologize for not responding! I haven't forgotten about this and will get to it, just been working on other stuff momentarily (trying to advance to candidacy...)

[milankl]

Of course, whenever you're ready just throw stuff in here. I get notifications anyway!