diff --git a/raster/r.sim/r.sim.water/testsuite/data/depth_complex.pack b/raster/r.sim/r.sim.water/testsuite/data/depth_complex.pack index 0b077ff259f..a2ce951b977 100644 Binary files a/raster/r.sim/r.sim.water/testsuite/data/depth_complex.pack and b/raster/r.sim/r.sim.water/testsuite/data/depth_complex.pack differ diff --git a/raster/r.sim/r.sim.water/testsuite/data/depth_infil.pack b/raster/r.sim/r.sim.water/testsuite/data/depth_infil.pack new file mode 100644 index 00000000000..4dbd33c735d Binary files /dev/null and b/raster/r.sim/r.sim.water/testsuite/data/depth_infil.pack differ diff --git a/raster/r.sim/r.sim.water/testsuite/data/discharge_complex.pack b/raster/r.sim/r.sim.water/testsuite/data/discharge_complex.pack index 78d3646962e..8e5f32be6cb 100644 Binary files a/raster/r.sim/r.sim.water/testsuite/data/discharge_complex.pack and b/raster/r.sim/r.sim.water/testsuite/data/discharge_complex.pack differ diff --git a/raster/r.sim/r.sim.water/testsuite/test_r_sim_water.py b/raster/r.sim/r.sim.water/testsuite/test_r_sim_water.py index 2eb3e2ae0d0..adfb6ce9e12 100644 --- a/raster/r.sim/r.sim.water/testsuite/test_r_sim_water.py +++ b/raster/r.sim/r.sim.water/testsuite/test_r_sim_water.py @@ -19,6 +19,7 @@ class TestRSimWater(TestCase): reference_discharge_default = "discharge_default" reference_depth_complex = "depth_complex" reference_discharge_complex = "discharge_complex" + reference_depth_infil = "depth_infil" @classmethod def setUpClass(cls): @@ -46,7 +47,7 @@ def setUpClass(cls): ) cls.runModule( "r.mapcalc", - expression=f"{cls.infil} = if (x() < 636500 && y() > 223500, 0.001, 0)", + expression=f"{cls.infil} = if (x() > 636500 && y() > 223500, 10, 0)", ) cls.runModule( "r.unpack", @@ -58,6 +59,11 @@ def setUpClass(cls): input="data/discharge_complex.pack", output=cls.reference_discharge_complex, ) + cls.runModule( + "r.unpack", + input="data/depth_infil.pack", + output=cls.reference_depth_infil, + ) @classmethod def tearDownClass(cls): @@ -154,6 +160,28 @@ def test_complex(self): precision="0.000001", ) + def test_infiltration(self): + """Test r.sim.water execution with infiltration""" + # Run the r.sim.water simulation + self.assertModule( + "r.sim.water", + elevation=self.elevation, + dx=self.dx, + dy=self.dy, + rain=self.infil, + infil=self.infil, + depth=self.depth, + niterations=30, + random_seed=1, + ) + + # Assert that the output rasters exist + self.assertRasterExists(self.depth) + # Assert that the output rasters are the same + self.assertRastersEqual( + self.depth, reference=self.reference_depth_infil, precision="0.000001" + ) + @unittest.skip("runs too long") class TestRSimWaterLarge(TestCase): diff --git a/raster/r.sim/simlib/hydro.c b/raster/r.sim/simlib/hydro.c index 5631c2e6d06..c106d6ffb9f 100644 --- a/raster/r.sim/simlib/hydro.c +++ b/raster/r.sim/simlib/hydro.c @@ -180,23 +180,26 @@ void main_loop(const Setup *setup, const Geometry *geometry, if (grids->zz[k][l] != UNDEF) { if (grids->inf[k][l] != UNDEF) { /* infiltration part */ - if (grids->inf[k][l] - grids->si[k][l] > 0.) { - double decr = pow( - addac * sim->w[lw].m, + double decr = + pow(addac * sim->w[lw].m, 3. / 5.); /* decreasing factor in m */ - if (grids->inf[k][l] > decr) { - grids->inf[k][l] -= - decr; /* decrease infilt. in cell - and eliminate the walker */ + if (grids->inf[k][l] > decr) { + grids->inf[k][l] -= + decr; /* decrease infilt. in cell + and eliminate the walker */ + sim->w[lw].m = 0.; + continue; + } + else { + sim->w[lw].m -= sim->rwalk * + grids->inf[k][l] / + setup->sisum; + grids->inf[k][l] = 0.; + // eliminate walker + if (sim->w[lw].m < 0.) { sim->w[lw].m = 0.; - } - else { - sim->w[lw].m -= - pow(grids->inf[k][l], 5. / 3.) / - addac; /* use just proportional part - of the walker weight */ - grids->inf[k][l] = 0.; + continue; } } } diff --git a/raster/r.sim/simlib/input.c b/raster/r.sim/simlib/input.c index faaca7d76b2..2a31d9ac1aa 100644 --- a/raster/r.sim/simlib/input.c +++ b/raster/r.sim/simlib/input.c @@ -291,9 +291,6 @@ int grad_check(Setup *setup, const Geometry *geometry, const Settings *settings, } /* DEFined area */ } } - if (grids->inf != NULL && smax < infmax) - G_warning(_("Infiltration exceeds the rainfall rate everywhere! No " - "overland flow.")); if (n == 0) { G_fatal_error(_("No values in the elevation raster.")); } @@ -372,8 +369,6 @@ int grad_check(Setup *setup, const Geometry *geometry, const Settings *settings, /*if(v1[k][l]*v1[k][l]+v2[k][l]*v2[k][l] > cellsize, warning, *napocitaj ak viac ako 10%a*/ /* THIS IS CORRECT SOLUTION currently commented out */ - if (grids->inf) - grids->inf[k][l] *= settings->timesec; if (inputs->wdepth) grids->gama[k][l] = 0.; if (outputs->et) {