Calculate explicit fluxes in the course of solving the Richards equation

The solution of the Richards equation in SoilWaterMovementMod is problematic for water tracers / isotopes: This calculates a redistribution of water between layers and adjusts the states without ever calculating explicit fluxes.

From a discussion with @davidnoone in June, 2017:

> Problem: Richards equation: Ideally, this would produce explicit fluxes between components and then apply those fluxes to update the states. But right now that doesn't happen: The solution of the Richards equation updates the states directly, and it's hard (or impossible) to infer the fluxes from there. So for now [referring to the CLM40 isotope implementation], they re-call the Richards equation solving code on the isotopes, even though they'd really like to apply the existing bulk fluxes to the isotopes.
>
> Ideal solution: Have the solution of the Richards equation compute explicit fluxes. Those fluxes would be used to update the state (and if there are any corrections to the fluxes needed to – e.g., – maintain non-negative states, then those corrections would be applied directly to the fluxes). Those same fluxes could then be used to update the isotopic states.

@djk2120 points out some additional issues with respect to transpiration, which I think are related to this:

> Hydraulic redistribution is an emergent behavior of the root water uptake when PHS is on. So there is no special HR code, but rather the PHS code will feature HR under certain conditions (e.g. at night time during a drought).
>
> The specific code that governs this is L313-314 of SoilWaterPlantSinkMod.F90
>
> temp(c) = temp(c) + k_soil_root(p,j) * (smp(c,j) - vegwp(p,4) - grav2)* patch%wtcol(p)
>
> which is accumulating the column root water uptake for a given soil layer ( j ) within a loop over the various PFTs, which will eventually be written to qflx_rootsoi_col(c,j)
>
> When there is HR, the values of qflx_rootsoi_col(c,j) will be positive for certain of the soil layers (and negative in others). This is driven by whether smp(c,j) is greater than or less than vepwp(p,4). (In SMS, this is disallowed and qflx_rootsoi_col values will always be either negative or zero). In the case of HR we are taking water from the places qflx_rootsoi_col is negative and depositing wherever qflx_rootsoi_col is positive, and potentially also sending some of the negative water out through the leaves for transpiration. At this point there is no need to exactly define which water goes where, which I imagine will be an issue with the isotopes. E.g. soil layers 5,6,and 7 are where there is root water uptake and send half the water to layer 2 and half of it out the leaves (non-specific) as compared to layer 5 sends its water to layer 2, and transpiration is supplied from layers 6 and 7 (specific accounting).
>
> We could certainly develop a heuristic for sorting this out (such as mixing the water), but in reality this may venture beyond a software question into being an actual science question.

My very vague understanding of this is that the redistribution within the soil mentioned by @djk2120 is tied in with the general Richards equation issue, but @djk2120 is pointing out some additional potential issues with respect to transpiration: we'll also need to know the explicit source(s) for that flux.