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 non-perturbative 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. As an application of this formalism, we find the $\beta$-functions and identify the fixed points in the large-$N$ limit of $2$-dimensional geometries in two different signatures.