Nonparabolicity effects in two-dimensional electron systems are quantitatively analyzed. A formalism has been developed which allows to incorporate a nonparabolic bulk dispersion relation into the Schrödinger equation. As a consequence of nonparabolicity the wave functions depend on the in-plane momentum. Each subband is parametrized by its energy, effective mass and a subband nonparabolicity coefficient. The formalism is implemented in a one-dimensional Schrödinger-Poisson solver which is applicable both to silicon inversion layers and heterostructures.