NOAA-OWP / LGAR-C

Lumped Arid/Semi-arid Model (LASAM) simulates infiltration and surface runoff (two important components of the hydrologic cycle) based on Layered Green-Ampt with redistribution (LGAR) model
Other
1 stars 4 forks source link

Ptl add adaptive timestep #25

Closed peterlafollette closed 4 months ago

peterlafollette commented 6 months ago

Adaptive time step was added and tested for stability.

Additions

-The model can now run with an adaptive time step. The adaptive time step can be turned on or off in the config file. If it is not specified in the config file, it defaults to false. The time step is equal to the forcing resolution when precip and ponded head are 0, which saves a lot of runtime in arid areas where large periods of time will not have precip. The time step is equal to two times a user specified minimum time step when the precip is less than or equal to 1 cm/h and there is no ponded head. The time step is equal to a user specified minimum time step (which is just taken from "timestep" in the config file when "adaptive_timestep" is set to true) when ponded head is greater than 0, or when precip is greater than 1 cm/h. The time step can not be greater than the forcing data resolution, and this is enforced.

Removals

-

Changes

-Updated README in configs to include adaptive_timestep.

-The giuh convolution integral is now computed at the temporal resolution of the forcing data rather than at the resolution of the sub timestep. The giuh ordinates are still resampled, but to the resolution of the forcing data rather than to the internal time step, which effectively means that the giuh ordinates are not resampled if the forcing data resolution is 1 hour.

-Discovered a bug where theta is very close to theta_r, and convergence might be possible, but the function lgar_theta_mass_balance doesn't converge after a few minutes. The bug comes from the fact that the condition (fabs(delta_mass - delta_mass_prev) < 1E-15) doesn't trigger; the change in mass is indeed small, but it is 1E-11 rather than 1E-15 and changes very slowly, just because of the sensitivity of the theta-psi relationship for some parameter sets in particular. Rather than changing the tolerance to be less strict, I added a new check where the loop will break for extremely large psi values, which have limited physical meaning in the first place.

-Changed the psi value at which calc_Se_from_h returns 1.0 to be 1E-10 rather than 1E-1, making the function accurate for small psi values. There are some cases where we need Se values that are very close to but slightly less than 1, for psi values that are less than 1E-1.

Testing

  1. Stability testing with 40 k parameter sets successful (adaptive and fixed time steps).
  2. Phillipsburg, Bushland, synthetic, and unit tests give expected outputs.
  3. Integrated tests on GitHub.
  4. Tested in ngen, including 3 layer tests, 1 layer tests, and running multiple catchments from the same realization file, including a mix of adaptive and fixed time step config files.

Screenshots

Notes

-

Todos

-

Checklist

Testing checklist

Target Environment support

Accessibility

Other