DivU constraint with arbitrary sources

[dmontgomeryNREL]

This PR updates the divu constraint to be compatible with arbitrary source terms. This is necessary for system closure in any model that introduces external sources on density (e.g. Spray). This should be consistent with the user defined sources proposed in PR #427.

[dmontgomeryNREL]

After further review, it appears that the current divu constraint works for arbitrary source terms. I restored the code from this PR to the most recent commit in development, but added some additional clarification in the Docs derivation.

[baperry2]

@dmontgomeryNREL can you add a comment in the compute_divu source code with a brief explanation as well?

@drummerdoc - here is the derivation from @dmontgomeryNREL and I that supports this, can you take a quick look and make sure you also agree because it's an important point.

[dmontgomeryNREL]

@dmontgomeryNREL can you add a comment in the compute_divu source code with a brief explanation as well? Done!

[drummerdoc]

In (23) above, Ym is the mass fraction of the local gas (same as the Ym in (32)), correct? Note that this says that mass is added at the local composition. Where do you account for mass being added at a different composition (like an evaporating fuel) - if that's even possible.... taking a minute to wrap my head around this.

[baperry2]

$Y_m$ is the local composition. Note Equation 23 is just a mathematical rearrangement using the product rule.

The physical meaning is more apparent in Equation 24, which is the complete equation for $DY_m/Dt$. That equation has two source terms that give appropriate behavior for the affect on $Ym$: `\dot{S}{\rho Y_m} - Ym\dot{S}{\rho}`

Consider an H2O drop evaporating into pure N2, so $Y{N2} =1$, $Y{H2O} = 0$, and \dot{S}_{\rho Y_N2} = 0 and \dot{S}_{\rho Y_H2O} = \dot{S}_{\rho}. For N2, the first term is 0 so the mass fraction is decreasing due to H2O evaporating. For H2O, the second term is 0 so the mass fraction is increasing.

[drummerdoc]

The subtlety is the blending that occurs between the ambient composition and that of the source, yet the sum of the mass fractions has to remain 1. In chemical reactions, a preferential loss of a given species occurs while respecting a conservation of atoms, so nothing really vanishes or appears, just switches bonds. Here, like at an inflow boundary, we are adding specific atoms, and then we re-normalize so that the sum of the new Ym's again equal 1. The head-hurting part is how to do that consistently as a volumetric source and make sure that when we get to the end of the time step that you have properly accounted for the injection of all the added atoms. Maybe there's nothing weird here, but I can certainly imagine some sanity checks to ensure detailed conservation is being satisfied properly. Let's try to add something like that in zero or one dimension as a regression test - maybe as a later PR.

[baperry2]

I think the bookkeeping should all be fine because we do things in terms sources of rhoYs for the actual effect on the transported species. In testing this PR out I didn't see any issues with species conservation, just stuff relating to divu as described in #429.

So I think this PR is good to go, although we may need to revisit some aspects.