For water tracers: How should we handle bidirectional fluxes?

[billsacks]

@davidnoone - this is largely a question for you:

For a flux that always takes mass out of state variable X and puts it in state variable Y, we know how to compute the tracer flux given the bulk flux:

tracer_flux = bulk_flux * (tracer_state_X / bulk_state_X)

But we're not sure about how to handle bidirectional fluxes – which can be moving mass either from state X to state Y or from state Y to state X, depending on the sign. This can be the case for belowground water movement, for example.

One way I can see to handle this is to compute the bidirectional bulk flux, but then post-process it into two flux variables: one from state X to state Y, the other from state Y to state X. In a given timestep, one of those would always be 0 and the other non-negative. Then we can do the standard tracer calculation on each of those two fluxes.

Would you recommend that we use that approach, or is there some other way this is typically handled?

[billsacks]

@nusbaume - I'm wondering if you have any insight into this?

[nusbaume]

Hi @billsacks,

To me your proposed solution sounds reasonable, assuming we know where the flux is coming from. For example, if state_X and state_Y are 1-D arrays, and we that know for some array index i, e.g.:

if bulk_flux(i) < 0 then
  state_X(i) =  state_X(i) + bulk_flux
  state_Y(i-1) = state_Y(i-1) - bulk_flux
else
  state_X(i+1) =  state_X(i+1) + bulk_flux
  state_Y(i) = state_Y(i) - bulk_flux
end if

Then your method would work fine. However, if we can't easily map the flux back to a specific set of array indices to know which combination of fluxes to use then we might need to figure out a different strategy.

Do you have a particular code routine that you want to apply this technique to first? If you want I could take a quick look to try and find any possible issues.

[billsacks]

Thanks @nusbaume . I think the inspiration for opening this issue was related to #964 (the relevant part of the tech note is here: https://escomp.github.io/ctsm-docs/versions/master/html/tech_note/Hydrology/CLM50_Tech_Note_Hydrology.html?highlight=richards#soil-water). There, if I understand things correctly, the solution of the Richards equation leads to an update in the water content of each soil layer, where the flux between layers i and (i+1) could be in either direction. I was wondering if there are standard ways of solving this problem in the water isotope modeling community.

As I think about this more, I think there may be a more fundamental issue here: We could have situations where the net bulk water flux doesn't fully capture the gross fluxes between two pools (e.g., between two soil layers). As an extreme example that is easier for me to wrap my head around: I could imagine a situation where the net flux between two layers is zero, but this hides the fact that there is actually some mixing of water between these two layers – it's just that the flux from layer i to (i+1) is exactly balanced by the flux from layer (i+1) to i. For the sake of bulk water, we can ignore this, but if the two layers have different isotopic concentrations, then I think we would need to know the gross fluxes in order to correctly model the final state.

This is making me wonder if the solution to #964 may actually not be to try to back out explicit fluxes when we solve the Richards equation on the bulk water, but instead to solve the Richards equation separately for each tracer / isotope, which according to @davidnoone is what was done in CLM4. David said that a better solution would be to just solve the Richards equation on the bulk water and then apply the fluxes from that solution to each tracer, but now I'm wondering if that carries some issues like what I laid out above.

I can see a similar – and possibly more extreme – situation arising for lakes. There, I don't think we have any movement of bulk water between layers: none is needed, since the water mass in each layer is fixed in time. But I'm thinking we may actually need to add a mixing parameterization to handle the redistribution of isotopes or other water tracers between layers.

For now I'm mainly wondering if there are standard methods for dealing with these kinds of problems, so we don't reinvent the wheel.

[ekluzek]

Moving to a discussion.