BETAS, a spectral code for three-dimensional magnetohydrodynamic equilibrium and nonlinear stability calculations. A fully three-dimensional spectral code has been developed and implemented. It computes three-dimensional equilibria using the variational approach, by minimizing the potential energy subject to appropriate constraints. A second minimization, with an additional constraint, is used to examine the stability of solutions. The magnetic field representation allows for non-nested flux surfaces. Thus it can be used to study solutions with islands and variation of the potential energy with respect to topology changes. A spectral representation in the toroidal and poloidal angles, coupled with a special choice of collocation points in the radial direction results in greatly enhanced resolution. A fast iterative method was developed, with the number of iterations required for convergence independent of the mesh size. Residuals converge exponentially to the round-off error, allowing the potential energy to be computed to the eight-digit accuracy required for nonlinear stability analysis. A small amount of artificial viscosity may be needed for convergence in cases where low-order resonances are present in the plasma region. Numerical results show convergence to known axially symmetric equilibrium solutions, as well as close agreement with eigenvalue calculations in helically symmetric stellarator configurations. Solutions exhibiting island formations are found to have lower potential energy than that of the nearby nested case. Internal modes are found to localize inside the resonant surface, with the eigenfunction having sharp gradients at the resonant surface.