Submarine groundwater discharge from subterranean estuaries is affected by tides, which are represented in computational models as time-dependent boundary conditions on the seaward boundary. Conventionally, a small time step is used in the numerical model to phase-resolve the tidal signal so as to ensure accurate results, although at the cost of excessive computation times for long-term simulations. This study proposes a highly efficient alternative method for modeling the tidal signal, in which a phase-averaged pressure is assigned to the seawater boundary with a much larger time step. The assigned pressure condition is first determined from an analytical solution of the time-independent pressure boundary condition. Along with the analytical solution, a single calibration factor is introduced at the beach face to account for the conductance at the beach. This results in good agreement between the results for phase-averaged and phase-resolved simulations. The new method is verified by comparison of the results for a wide range of physical cases determined using TOUGHREACT, a model for simulating coupled hydrodynamic, thermodynamic, and geochemical processes. This comparison shows that the phase-averaged results give good agreement except for a small underestimation of the mixing zone over the saltwater wedge region. These results confirm that the new boundary condition is suitable for efficient, long-term simulations of coastal aquifers subjected to tidal forcing.