in a [*sufficiently accurate*]{} approximation of the [*quantity of interest*]{} $J(u)\in\IR$, where $$\begin{aligned}
  \label{Ju}
  \J:\Wh\to \IR\end{aligned}$$ is a linear functional. Typically, it is defined as a weighted mean value of $u$ over the computational domain $\Om$ or its boundary $\gom$.

In order to estimate the error $\J(\u) - J(\uh)$, we consider the [*adjoint*]{} (or [*dual*]{}) [*problem*]{} and its discretization. We say that $\zh$ is the [*discrete solution*]{} of the [*dual problem*]{} if $$\begin{aligned}
  \label{eq:DPh}
  \ah(\psi_h, \zh) = \J(\psi_h) \qquad \forall \psi_h\in\Vh,\end{aligned}$$ where $\ah$ and $\J$ are given by and , respectively. Moreover, we assume that the numerical scheme is [*adjoint consistent*]{}, i.e., $$\begin{aligned}
  \label{eq:consD}
  \ah(\psi, \z) = \J(\psi)\qquad \forall \psi\in\Wh,\end{aligned}$$ where $\z$ is the weak solution of the dual problem . This implies the [*Galerkin orthogonality of the error*]{} of the dual problem $$\begin{aligned}
  \label{eq:GOD}
  \ah(\psi_h, \zh - \z) = 0 \qquad\forall \psi_h\in\Vh.\end{aligned}$$ Finally, we define the [*residual of the dual problem*]{} by $$\begin{aligned}
  \label{resD}
  \resD{\zh}{\psi}:= \J(\psi) - \ah(\psi,\zh) = \ah(\psi, \z-\zh),\quad \psi\in\Wh,\end{aligned}$$ where the last equality follows from the adjoint consistency .

The problems and represent two linear algebraic systems whose efficient solution by an iterative solver is developed in §\[sec:algeb\]. Due to iterative and rounding errors, the “exact” 