It is shown that it is possible to construct, within the framework of a basis set expansion method, the full Foldy-Wouthuysen transformation (i.e., to all orders in the inverse velocity of light) for an arbitrary potential once the Dirac equation has been solved. On this basis an iterative procedure to solve the Dirac equation is suggested that involves only the large component, obviating the time-consuming (at least in molecular calculations) introduction of large basis sets for a proper description of just the small components. The methods are used to compare the expectation value of the radial distance operator in the Dirac picture and in the Schrodinger picture for the orbitals of the Uranium atom.