Solution approaches for evaporation-driven density instabilities in a slab of saturated porous media

TL;DR

This paper analyzes gravitational instabilities in saline porous media caused by evaporation, deriving analytical solutions, proposing a novel eigenvalue computation method, and validating results with numerical simulations.

Contribution

It introduces a new fundamental matrix method for eigenvalue problems and provides a comprehensive analysis of evaporation-driven density instabilities in saturated porous media.

Findings

Ground-state salinity profile derived analytically.

Critical Rayleigh number computed using the new method.

Numerical simulations confirm linear stability predictions.

Abstract

This work considers the gravitational instability of a saline boundary layer formed by an evaporation-induced flow through a fully-saturated porous slab. Evaporation of saline waters can eventually result in the formation of salt lakes as salt accumulates. As natural convection can impede the accumulation of salt, establishing a relation between its occurrence and the value of physical parameters such as evaporation rate, height of the slab or porosity is crucial. One step towards determining when gravitational instabilities can arise is to compute the ground-state salinity, that evolves due to the uniform upwards flow caused by evaporation. The resulting salt concentration profile exhibits a sharply increasing salt concentration near the surface, which can lead to a gravitationally unstable setting. In this work, this ground state is analytically derived within the framework of Sturm-Liouville theory. Then, the method of linear stability in conjunction with the quasi-steady state approach is employed to investigate the occurrence of instabilities. These instabilities can develop and grow over time depending on the Rayleigh number and the dimensionless height of the porous medium. To calculate the critical Rayleigh number, which can determine the stability of a particular system, the eigenvalues of the linear perturbation equations have to be computed. Here, a novel fundamental matrix method is proposed to solve this eigenvalue problem and shown to coincide with an established Chebyshev-Galerkin method in their shared range of applicability. Finally, a 2-dimensional direct numerical simulation of the full equation system via the finite volume method is employed to validate the time of onset of convective instabilities predicted by the linear theory. Moreover, the fully nonlinear convection patterns are analyzed.