The thin-plate method of modeling neotectonics uses isostasy and vertical integration of lithospheric strength to reduce three-dimensional problems to two dimensions. We introduce new thin-shell continuum elements and fault elements which satisfy the completeness and compatibility requirements, and permit extension of these methods to the lithospheres of spherical planets. Even a coarse grid of these elements with elastic rheology can reproduce low-order toroidal free oscillations. In the program SHELLS, a realistic frictional/dislocation-creep rheology is handled by iteration; this method converges to solutions which can be numerically tested to confirm that they satisfy both local equilibrium and the balance of torques on whole plates. Because this code can incorporate both natural plate shapes (with internal faults) and realistic rheologies, it yields models that are readily tested by their predictions of geodetic velocities, stresses, and fault slip rates. ¿ American Geophysical Union 1995