Fig. 1. Location of the Okavango Delta () in northern Botswana, southern Africa.

Fig. 2. Aerial view of islands in the Okavango Delta. Lighter areas represent concentrated patches of salt at the surface.

Fig. 3. Landsat satellite image of the Okavango Delta showing the location of the surveyed Thata Island.

Fig. 4. Layout of borehole-to-surface imaging on Thata Island.

Fig. 5. Results of borehole-to-surface resistivity imaging on Thata Island, in [Ωm].

Fig. 6. Comparison of dispersion treatment in particle tracking methods based on grid cells such as the method of characteristics (MOC) (left) and on particle locations such as the particle strength exchange method (PSE) (right). In this example the sphere of influence in the PSE method in chosen as five times the core size of the particles.

Fig. 7. Thata Island in the Okavango Delta. Evapotranspiration has resulted in the precipitation of salts and the formation of a Trona salt crust. Only the fringe of the island is vegetated.

Fig. 8. Top: Sketch of schematic island in an advanced state. Bottom: Schematic cross-section of the conceptual model. If the main channel is far away from the island, an ‘open’ boundary condition with fixed pressure over the whole depth of the domain is more adequate.

Fig. 9. Test of grid convergence with the d³f model. From top to bottom: grid level 4 (64 × 16 nodes), grid level 5 (128 × 32 nodes), grid level 6 (256 × 64 nodes), grid level 7 (512 × 128 nodes), and grid level 8 (1024 × 256 nodes). Simulation results are shown at an elapsed time of 1000 years. The contour lines represent concentrations from 0.5 to 5 g/l with an interval of 0.5 g/l. For model parameters see Table 1. Domain size is 340 m × 100 m.

Fig. 10. Mass flux through a line parallel to the island surface in 50 m depth at different grid levels.

Fig. 11. Evolution of density fingering for the parameters summarized in Table 1. Left column: d³f results, aquifer is modeled as homogeneous. Middle column: corresponding pse2d results. Right column: pse2d results obtained with a less conductive layer beneath the island. Displayed times are 60, 200, 400, 500, and 1000 years. The contour lines represent concentrations from 0.5 to 5 g/l at intervals of 0.5 g/l. Domain size is 340 m × 100 m.

Fig. 12. Development of the maximum concentration in the center of the island over time for the homogeneous case (—) and the inhomogeneous case (- · - · -).

Fig. 13. Inter-code comparison of the solute distribution for the inhomogeneous case after 2000 years. Numerical results from top to bottom: pse2d, SEAWAT with FD, SEAWAT with TVD, and SEAWAT with MOC. The contour lines represent concentrations between 0.5 and 5 g/l with intervals of 0.5 g/l and concentrations of 6, 8, and 10 g/l. Domain size is 340 m × 100 m.

Fig. 14. Maximum concentrations obtained with the methods: pse2d, TVD, FD and MOC, corresponding to Fig. 13.

Fig. 15. Mass accumulation over time for different numerical methods. Theoretical mass accumulation coincides with FD and TVD results.

Fig. 16. Evolution of density fingering on Thata Island. Results displayed at 100, 150, 200, 250, 260, 300, 400, 500, 1000, and 2000 years. Domain size is 800 m × 100 m.

Fig. 17. Maximum concentration below the island in [g/l] over time corresponding to Fig. 16.

Fig. 18. Ratio between salt mass entering the model domain and salt mass convecting down into deeper aquifer units.

Fig. 19. Pse2d simulation results after 500 years performed with different ET-rates in a homogeneous aquifer. Domain size is 340 m × 100 m.

Fig. 20. Initial velocity fields for evapotranspiration rates of 0.2, 0.4, 0.6 and 1.0 mm/day (from top to bottom). Domain size is 340 m × 100 m.

Fig. 22. Initial velocity fields for the “closed” case and for hydraulic gradients of 0.1‰, 0.5‰, and 1‰ (from top to bottom). Domain size is 340 m × 100 m.

Fig. 21. Pse2d simulation results for the “closed” case after 500 years with different hydraulic gradients in a homogeneous aquifer. Domain size is 340 m × 100 m. The model domain is 340 m long and 100 m deep.

Fig. 23. Pse2d simulation results for the “open” case after 500 years with different hydraulic gradients. Parameters taken for the simulation are summarized in Table 2. The domain displayed is 150 m  x  800 m, and 0 m  z  100 m.

Fig. 24. Initial velocity fields for the “open” case and hydraulic gradients of 0.1‰, 0.5‰ and 1‰ (from top to bottom). The model domain is 800 m long and 100 m deep.

Fig. 25. Maximum concentration in the center of the island over time for different hydraulic gradients.

Fig. 26. 2D simulation results with SEAWAT at t = 500 years with different hydraulic gradients. The domain displayed is 150 m  x  800 m, and 2.5 m  z  100 m.

Fig. 27. 3D simulation results with SEAWAT at t = 500 years with different hydraulic gradients. The island stretches across the whole width of the domain. Displayed are horizontal slices in z = 2.5 m, 12.5 m and 27.5 m depth (from left to right), and 150 m  x  800 m, and 0 m  y  200 m.

Fig. 28. 3D simulation results with SEAWAT at t = 500 years with different hydraulic gradients. The island stretches across the whole width of the domain. Displayed are vertical cross-sections at y = 112.5 m, and 150 m  x  800 m, and 0 m  y  200 m.

Fig. 29. 3D simulation results with SEAWAT at t = 500 years with different hydraulic gradients. The island is elongated with a width of 50 m. Displayed are horizontal slices in z = 2.5 m, 12.5 m and 27.5 m depth (from left to right), and 150 m  x  800 m, and 0 m  y  200 m.

Fig. 30. 3D simulation results with SEAWAT at t = 500 years with different hydraulic gradients. The island is elongated with a width of 50 m. Displayed are vertical cross-sections at y = 112.5 m (top), and 150 m  x  800 m, and 0 m  y  200 m.

Fig. 31. Mass accumulation in the model domain with time for different simulations and different hydraulic gradients. Top: 2D simulation results (Fig. 26). Middle: 3D simulation results with island that stretches across the whole width of the domain (Fig. 27). Bottom: 3D simulation results with island of 50 m width (Fig. 29).

Parameter used for the numerical study with ‘closed’ lateral boundaries

Parameters used for the Thata Island study with ‘open’ lateral boundaries

Parameter list of different simulations performed, where dist. = distribution, fp. = floodplain, isl. = island, hom = homogeneous, het = heterogeneous