The vertical velocity, w, in three-dimensional circulation models is typically computed from the three-dimensional continuity equation assuming that the depth-varying horizontal velocity field was calculated earlier in the solution sequence. Computing w in this way appears to require the solution of an over-determined system since the continuity equation is first order, yet w must satisfy two boundary conditions (one at the free surface and one at the bottom). At least three methods have been previously proposed to compute w: (i) the “traditional” method that solves the continuity equation with the bottom boundary condition and ignores the free surface boundary condition, (ii) a “vertical derivative” method that solves the vertical derivative of the continuity equation using both boundary conditions and (iii) an “adjoint” approach that minimizes a cost functional comprised of residuals in the continuity equation and in both boundary conditions. The latter solution is equivalent to the "traditional" solution plus a correction that increases linearly over the depth and is proportional to the misfit between the "traditional" solution at the surface and the surface boundary condition. In this paper we show that the "vertical derivative" method yields inaccurate and physically inconsistent results if it is discretized as has been previously proposed. However, if properly discretized the "vertical derivative" method is equivalent to the “adjoint” method if the cost function is weighted to exactly satisfy the boundary conditions. Furthermore, if the horizontal flow field satisfies the depth-integrated continuity equation locally, one of the boundary conditions is redundant and w obtained from the "traditional" method should match the free surface boundary condition. In this case, the “traditional,” “adjoint” and properly discretized “vertical derivative” approaches yield the same results for w. If the horizontal flow field is not locally mass conserving, the mass conservation error is transferred into the solution for w. This is particularly important for models that do not guarantee local mass conservation, such as some finite element models