We analyze the impact of using a penalty method to enforce Dirichlet boundary conditions in the Finite Element method, focusing our attention on the approximation of boundary fluxes. In particular, we demonstrate through both theory and numerical experiments, that the use of the penalty method may lead to large errors in the estimate of the flux in addition to the usual discretization errors. We propose here a new estimator of the boundary flux, which involves a term defined in terms of the derivative of the flux with respect to the penalty parameter. It is shown that this estimator allows one to substantially reduce the error made due to the introduction of the penalty method. The errors in the flux are then dominated by discretization errors, which can be controlled using adjoint- based techniques. An adjoint problem associated with the boundary flux is also proposed based on the analysis of the penalty method. Finally, we present some numerical experiments that demonstrate that the solution of this adjoint problem can be used for goal-oriented adaptive refinement, along with the use of the improved estimator for the flux, allowing us to control both the discretization and penalty errors.