  Building on the established methods for superconducting circuit quantization,
we present a new theoretical framework for approximate numerical simulation of
Josephson quantum circuits. Simulations based on this framework provide access
to a degree of complexity and circuit size heretofore inaccessible to
quantitative analysis, including fundamentally new kinds of superconducting
quantum devices. This capability is made possible by two improvements over
previous methods: first, physically-motivated choices for the canonical circuit
modes and physical basis states which allow a highly-efficient matrix
representation; and second, an iterative method in which subsystems are
diagonalized separately and then coupled together, at increasing size scales
with each iteration, allowing diagonalization of Hamiltonians in extremely
large Hilbert spaces to be approximated using a sequence of diagonalizations in
much smaller spaces.
