AbstractAlthough hydraulic groundwater (GW) theory has been recognized as a promising tool for understanding the role of the aquifer(s) in the surface‐subsurface hydrologic cycle, the integrated modeling community still lacks a proper hydrologic structure to apply the well‐studied theory to large‐scale hydrologic predictions. This study aims to present a novel hydrologic structure that enables the Boussinesq equation‐based depiction of the bidirectional stream‐hillslope processes for applying hydraulic GW theory to large‐scale model configurations. We integrated the BE3S's (Hong et al., 2020, https://doi.org/10.1029/2020wr027571) representation scheme of the catchment‐scale stream‐hillslope continuum into the National Water Model (NWM) and applied the modified NWM (i.e., the NWM‐BE3S) to three major basins in Texas (i.e., the Trinity, Brazos, and Colorado River basins). Since the NWM currently relies on a single reservoir model for baseflow simulation, we used the Boussinesq aquifer as an alternative subsurface hydrology routine and evaluated its predictive skill and efficacy. We identified that the implemented Boussinesq aquifer(s) in the NWM‐BE3S yielded noticeable improvements in predicting streamflow for aquifers that exhibited higher nonlinearities in the observed recessions. The varying degree of improvements in streamflow predictions per the recession nonlinearities demonstrated not only (a) the algorithmic enhancement of subsurface hydrology (physics) but also (b) the applicability of the Boussinesq theory‐based depiction of the stream‐hillslope two‐way continuum. We diagnosed each stream's state based on the bidirectional stream‐hillslope exchanges and identified the dominant processes (i.e., river infiltration or baseflow) that were represented spatially in the NWM‐BE3S.