The interaction and propagation of optical pulses in a nonlinear waveguide array is described by the discrete nonlinear Schrödinger equation (1) i ∂ z ψ n = − D ( ψ n + 1 + ψ n − 1 − 2 ψ n ) − γ | ψ n | 2 ψ n , where D is a dispersion (or diffraction) coefficient, and γ is a measure of the nonlinearity. By means of the variational approximation, we study the discrete soliton solutions of this equation. We use a trial function which contains six parameters, corresponding to: position, phase, amplitude, wavevector (velocity), chirp, and width. With this trial function, we can analytically average the Lagrangian, and then by varing the six parameters, obtain the evolution equations for these six parameters, within the variational approximation. Integration of these equations would give, within the variational approximation, the motion of a moving discrete soliton. Requiring all parameters to be stationary, one obtains the conditions for constructing the solution of the stationary discrete soliton. Here we treat the stationary variational solutions. For them, we find for small amplitudes, that there is only one stationary soliton, a doublet solution, which in the continuous limit, becomes the usual nonlinear Schrödinger soliton. At a certain critical amplitude, there is a pitchfork bifurcation, above which the stable singlet soliton apprears, with the doublet soliton becoming unstable. Lastly, using the variational solutions as a starting point, we iterate the full Lagrangian to obtain numerically, the exact discrete soliton solutions. Comparison between the variational and the exact numerical discrete soliton solutions will be made. From these results, we are also able to make some general and important remarks concerning the validity and utility of the found variational soliton solutions.