We propose two methods for the unsupervised detection of communities in undirected multiplex networks. These networks consist of multiple layers that record different relationships between the same entities or incorporate data from different sources. Both methods are formulated as gradient flows of suitable energy functionals: the first (MPBTV) builds on the minimization of a balanced total variation functional, which we show to be equivalent to multiplex modularity maximization, while the second (DGFM3) directly maximizes multiplex modularity. The resulting non-linear matrix-valued ordinary differential equations (ODEs) are solved efficiently by a graph Merriman--Bence--Osher (MBO) scheme. Key to the efficiency is the approximate integration of the discrete linear differential operators by truncated eigendecompositions in the matrix exponential function. Numerical experiments on several real-world multiplex networks show that our methods are competitive with the state of the art with respect to various metrics. Their major benefit is a significant reduction of computational complexity leading to runtimes that are orders of magnitude faster for large multiplex networks.