We consider a large class of random geometric graphs constructed from samples $\mathcal{X}_n = \{X_1,X_2,\ldots,X_n\}$ of independent, identically distributed observations of an underlying probability measure $ν$ on a bounded domain $D\subset \mathbb{R}^d$. The popular `modularity' clustering method specifies a partition $\mathcal{U}_n$ of the set $\mathcal{X}_n$ as the solution of an optimization problem. In this paper, under conditions on $ν$ and $D$, we derive scaling limits of the modularity clustering on random geometric graphs. Among other results, we show a geometric form of consistency: When the number of clusters is a priori bounded above, the discrete optimal partitions $\mathcal{U}_n$ converge in a certain sense to a continuum partition $\mathcal{U}$ of the underlying domain $D$, characterized as the solution of a type of Kelvin's shape optimization problem.