mMax = 50;
Ffuncs[r_] = Table[F[n, m][r], {n, mMax}, {m, mMax}];
(nsols = NSolve[
Flatten@{Thread[
Flatten[D[Ffuncs[r], r]] == Flatten[Ffuncs[r].Ffuncs[r]]]},
Flatten@Thread[Flatten[D[Ffuncs[r], r]]]
] /. 1. -> 1) // {ListQ[#], Dimensions[#]} & // AbsoluteTiming
nf = Experimental`CreateNumericalFunction[
Join[{r}, Flatten@Ffuncs[r] /. subs],
Flatten[D[Ffuncs[r], r]] /. First@nsols /.
rhs_ :> (Hold[#] &@subs /. Rule -> Set /.
Hold[v_] :> Block[{F}, v; rhs]),