Computational framework for resolving boundary layers in electrochemical systems using weak imposition of Dirichlet boundary conditions

Abstract

We present a finite element based computational framework to model electrochemical systems. The electrochemical system is represented by the coupled Poisson–Nernst–Planck (PNP) and Navier–Stokes (NS) equations. The key quantity of interest in such simulations is the current (flux) at the system boundaries. Accurately computing the current flux is challenging due to the small critical dimension of the boundary layers (small Debye layer) that require fine mesh resolution at the boundaries. We present a numerical framework which resolves this challenge by utilizing a weak imposition of Dirichlet boundary conditions for Poisson–Nernst–Planck equations. In this numerical framework we utilize a block iterative strategy to solve NS and PNP equations. This allows us to efficiently and easily implement the weak imposition of Dirichlet boundary conditions. The results from our numerical framework shows excellent agreement when compared to strong imposition of boundary conditions (strong imposition requires a much finer mesh). Furthermore, we show that the weak imposition of the boundary conditions allows us to resolve the fluxes in the boundary layers with much coarser meshes compared to strong imposition. We also show that the method converges as we refine the mesh near the boundaries at a much faster rate compared to strong imposition of the boundary layer. We present multiple test cases with varying boundary layer thickness to illustrate the utility of the numerical framework. We illustrate the approach on canonical 3D problems that otherwise would have been computationally intractable to solve accurately. Lastly, we simulate electrokinetic instabilities near a perm-selective membrane with weakly imposed boundary conditions on the membrane. This approach substantially reduces the computational cost of modeling thin boundary layers in electrochemical systems.