We introduce a model for describing the dynamics of large numbers of interacting cells. The fundamental dynamical variables in the model are subcellular elements, which interact with each other through phenomenological intra- and intercellular potentials. Advantages of the model include i) adaptive cell-shape dynamics, ii) flexible accommodation of additional intracellular biology, and iii) the absence of an underlying grid. We present here a detailed description of the model, and use successive mean-field approximations to connect it to more coarse-grained approaches, such as discrete cell-based algorithms and coupled partial differential equations. We also discuss efficient algorithms for encoding the model, and give an example of a simulation of an epithelial sheet. Given the biological flexibility of the model, we propose that it can be used effectively for modeling a range of multicellular processes, such as tumor dynamics and embryogenesis.