We develop an effective mass theory for substitutional donors in silicon in an inhomogeneous environment. Valley-orbit coupling is included perturbatively. We specifically consider the Stark effect in Si:P. In this case, the theory becomes more accurate at high fields, while it is designed to give correct experimental binding energies at zero field. Unexpectedly, the ground state energy for the donor electron is found to increase with electric field as a consequence of spectrum narrowing of the 1s manifold. Our results are of particular importance for the Kane quantum computer.Phosphorus doped silicon is one of the most wellstudied semiconducting systems [1], owing to its importance in the electronics industry. A more recent and exotic application is quantum computing, in which the Kane proposal posits nuclear spin qubits on Si:P, with interactions mediated by donor-bound electrons [2]. By tuning the potentials on aligned electrodes, the electrons are brought to the point of ionization, in order to control nuclear hyperfine interactions and the overlap between neighboring qubits. Spin dependent ionization also provides a means for electrical detection of the qubit state. The precise control of donor-bound electrons in such a complex environment remains an experimental challenge [3], requiring detailed theoretical input [4,5,6].The theory of isolated donors in silicon remains one of the crowning achievements of solid-state physics. One of the most effective treatments for shallow donors is the effective mass approximation (EMA) [7], in which the donor potential is assumed to vary slowly compared to the crystal potential. As a result, the long and shortwavelength physics decouple. An excess electron on the donor can be described by an envelope equation -a Schrödinger equation with an effective mass and a dielectric constant. The theory highlights the most essential feature of silicon's bulk band structure: the six-fold degeneracy of the conduction valleys. Near the impurity core, the assumption of slowly varying potential breaks down, causing valley-orbit interactions to couple envelope functions in different valleys. A careful treatment of the potential very near the impurity (the 'central cell' region) enables estimates for the energy splitting of the six valley states, which are in good agreement with experiments [8].In a more general, inhomogeneous environment, silicon donors can be studied by tight-binding techniques [9]. However, existing EMA theories introduce severe approximations that limit their scope. Many important questions therefore remain open. For example, there is no theory to determine how the weight of a donor wavefunction will redistribute among the six conduction valleys in the presence of a generic (i.e., low symmetry) potential. Since the envelope functions in different valleys have different energies (due to anisotropy of the effective mass), this is an important question.In the context of quantum computing, several recent papers obtain tractible results for donor ionization, by ignoring va...