Random noncommutative geometry can be seen as a Euclidean path-integral quantization approach to the theory defined by the Spectral Action in noncommutative geometry (NCG). With the aim of investigating phase transitions in random NCG of arbitrary dimension, we study the nonperturbative Functional Renormalization Group for multimatrix models whose action consists of noncommutative polynomials in Hermitian and anti-Hermitian matrices. Such structure is dictated by the Spectral Action for the Dirac operator in Barrett’s spectral triple formulation of fuzzy spaces. The present mathematically rigorous treatment puts forward “coordinate-free” language that might be useful also elsewhere, all the more so because our approach holds for general multimatrix models. The toolkit is a noncommutative calculus on the free algebra that allows to describe the generator of the renormalization group flow—a noncommutative Laplacian introduced here—in terms of Voiculescu’s cyclic gradient and Rota–Sagan–Stein noncommutative derivative. We explore the algebraic structure of the Functional Renormalization Group equation and, as an application of this formalism, we find the beta -functions, identify the fixed points in the large-N limit and obtain the critical exponents of two-dimensional geometries in two different signatures.