CliMA / LandHydrology.jl

The CliMA land hydrology model, including soil, snow, and surface processes
6 stars 0 forks source link

surface fluxes boundary conditions #54

Closed kmdeck closed 2 years ago

kmdeck commented 2 years ago

PR adding in surface fluxes (prescribed atmos/reanalysis driving soil model)

This PR only supports this type of boundary condition with the full soil model - prognostic energy and prognostic water.

This does not yet include sublimation of ice or precipitation, the flux associated with a change in gravitational potential, or the resistance due to the soil dry layer.

codecov[bot] commented 2 years ago

Codecov Report

Merging #54 (cccd858) into main (c574172) will increase coverage by 0.11%. The diff coverage is 100.00%.

:exclamation: Current head cccd858 differs from pull request most recent head bd08710. Consider uploading reports for the commit bd08710 to get more accurate results Impacted file tree graph

@@            Coverage Diff             @@
##             main      #54      +/-   ##
==========================================
+ Coverage   98.61%   98.72%   +0.11%     
==========================================
  Files          10       10              
  Lines         432      472      +40     
==========================================
+ Hits          426      466      +40     
  Misses          6        6              
Impacted Files Coverage Δ
src/SoilModel/models.jl 100.00% <ø> (ø)
src/SoilModel/parameters.jl 100.00% <ø> (ø)
src/SoilModel/SoilWaterParameterizations.jl 100.00% <100.00%> (ø)
src/SoilModel/boundary_conditions.jl 97.72% <100.00%> (+0.98%) :arrow_up:

Continue to review full report at Codecov.

Legend - Click here to learn more Δ = absolute <relative> (impact), ø = not affected, ? = missing data Powered by Codecov. Last update c574172...bd08710. Read the comment docs.

kmdeck commented 2 years ago

"Hiding" comment b/c it may be outdated. TBD.

As evaporation dries out the surface, we expect the matric potential at the surface to diverge to -inf. According to e.g. Eqn 7.26 of Bonan's text, this implies the specific humidity in the soil at the surface goes to zero, which is what we'd expect...it's harder to evaporate if the liquid water is in a deep potential in the soil compared with e.g. outside of that potential. However, the factor (Eqn 7.27) does not necessarily approach zero quickly. To decrease q from the saturation specific humidity by a factor of e, the matric potential must be -14,000. For some soils, (e.g. sand), this won't happen unless the effective saturation is ~1e-9. This seems to cause numerical issues - evaporation continues to happen at a rate close to the potential evaporation rate, while the soil water content approaches the residual, and I have to take very tiny steps to correctly remove water at the surface due to evaporation without overshooting the residual amount. (E*dt >> theta - theta_residual, ignoring the flux from below towards the surface).

[ I tried "softening": I used theta = max(theta, theta_residual + epsilon). But Im not sure this will help in all cases. in the sand case above, if epsilon > 1e-9, we would continue evaporation at close to the potential rate (or around 1/e smaller), b/c the factor 7.27 would be >~1/e. ] evap_1e-6

evap_1e-8

evap_1e-10

The soil resistance term will help (it isnt included yet; it's meant to account for diffusion within a dry surface layer forming between the evaporative front and the surface), as it will decrease evaporation as the soil grows drier. But even in the limit of zero water content, this resistance does not diverge to infinity (possibly a flaw in construction?) - it goes towards a constant, I'm not yet sure if it will end up affecting the evaporative flux that much (have to check how it compares to aerodynamic resistance values).

Is there something else missing? should we be estimating the evaporating area of the soil surface as proportional to volumetric water content^(2/3), for example? Otherwise evaporating from very dry sand (without the soil resistance) is essentially the same as evaporating from a body of water!

kmdeck commented 2 years ago

heat_fluxes moisture_fluxes profiles

Latest iteration seems to make sense. Currently computing the SHF and evaporation with SurfaceFluxes.jl. The total heat flux boundary condition at the top is a combination of the SHF and the volumetric_internal_energy of liquid water in the soil*evaporation. Note that if I didnt include the term, I got funny results - see below evaporation results.pdf

I have not added in latent heat fluxes yet.

kmdeck commented 2 years ago

bors try

bors[bot] commented 2 years ago

try

Build succeeded:

kmdeck commented 2 years ago

In the end, I added exactly what is in the design doc. I think we don't need to include the volumetric energy loss from liquid water as well - only the contribution from the vapor loss. The results look reasonable for a case where T_atm > T_surf (so negative SHF), but positive evaporation -> positive LHF dominating at the beginning. Things are moving (very slowly) towards equilibrium, which looks to be consistent with T_surf = T_atm , q_atm = q_surf, E = 0, SHF = 0.

heat_fluxes moisture_fluxes profiles

kmdeck commented 2 years ago

@charleskawczynski @yairchn @akshaysridhar @Yujie-W

does this look OK to merge? thanks!

kmdeck commented 2 years ago

bors try

bors[bot] commented 2 years ago

try

Build succeeded:

kmdeck commented 2 years ago

bors try

bors[bot] commented 2 years ago

try

Build failed:

kmdeck commented 2 years ago

bors r+

bors[bot] commented 2 years ago

Build succeeded: