When using the SIMPLE method on a mesh with a collocated variable arrangement, the following interpolation is used for the advecting velocities:
\begin{equation} u_f = \overline{u}_f - \overline{D}_f\left(\left(\nabla P\right)_f - \overline{\left(\nabla P\right)}_f\right) \end{equation}
where the over bar denotes a geometrically interpolated quantity and
\begin{equation} D_P = \frac{V_P}{a_P} \end{equation}
where $V_P$ is the volume of cell $P$ and $a_P$ is the central coefficient arising from the discretized momentum equation,
\begin{equation} a_Pu_P + \sum_{F}a_Fu_F = -V_P\nabla P \end{equation}
Anyways, my question is, how do you correctly determine $\overline{D}_f$ at boundaries? I recently wrote a SIMPLE solver for unstructured meshes, but have noticed that the largest continuity errors are always occurring at the boundary cells, leading me to think I may have accounted for this term incorrectly (I simply extrapolate it to the boundary face). This term can be very important in the computations as it also shows up in the pressure correction equation. I also have some difficulty getting the solver to converge when using inlet/outlet boundary conditions.