The crystalline swelling properties and interlayer structure of a cesium montmorillonite clay were investigated using molecular computer simulations. Two classes of dry-clay structures, proposed previously to explain X-ray diffraction and NMR experiments, were identified using Monte Carlo annealing calculations. Hydrated clays with water contents ranging from 0.044 to 0.440 g(H2O)/g(clay) were investigated using molecular dynamics simulations. Layer spacings calculated as a function of water content were found to be similar to experimental swelling curves, showing a distinct plateau at the monolayer-hydrate spacing. Hydration energies were calculated as a function of water content and expressed in three complementary forms. The immersion energy form was found to be most useful, revealing an apparently global minimum in the swelling-coordinate energy that corresponds to the monolayer hydrate. This is in agreement with experimental measurements and may help clarify the energetic origins of discrete, crystalline swelling processes in clay minerals. For two- and four-layer hydrates, cesium ions readily formed two different types of inner-sphere complexes with the clay surface. Ions associated with negatively charged tetrahedral substitution sites formed exclusively inner-sphere complexes and occupied hexagonal cavities adjacent to the substitutions. Other cesium ions occupied both inner-sphere and outer-sphere configurations with roughly equal probability. The ease with which cesium associates with the clay surface may be responsible for the formation of monolayer hydrates in cesium-substituted clays and for selective binding of cesium to many clay minerals.