To approximate the problem of natural convection in porous media, we utilize the mixed finite element method and interior penalty enriched Galerkin approximation for spatial discretization and the 4th-order backward differentiation formula for time-stepping as our full order model (FOM) as previously implemented by the authors [1]}. This combination is selected to guarantee local mass conservation and accurately capture the gravity-driven flow in the advection-dominated regime [2]}, [3]}. However, this FOM requires substantial computational resources because it is highly nonlinear and needs a very small time-step to satisfy Courant–Friedrichs–Lewy (CFL) condition [4]}, [5]}, [6]}, [7]}. Hence, it is not directly suitable to handle large-scale inverse problems, optimization, or control, in which an extensive set of simulations must be explored [8]}, [9]}, [10]}, [11]}, [12]}, [13]}, [14]}.
