%I #16 Nov 27 2023 13:05:26
%S 1,4,12,84,3948,88320,815424,14652150
%N Least m such that A081212(m) = n.
%p gsqfree := proc(n) local resul ; resul := n ; while not numtheory[issqrfree](resul) do resul := resul-1 ; od ; RETURN(resul) ; end: A081210 := proc(n) option remember ; local pfd,resul,p ; if n = 1 then RETURN(1) ; else pfd := ifactors(n)[2] ; resul := 1 ; for p from 1 to nops(pfd) do resul := resul*gsqfree( op(1,op(p,pfd))^op(2,op(p,pfd))) ; od ; fi ; end: r := proc(n,k) option remember ; if k= 0 then n ; else r(A081210(n),k-1) ; fi ; end: A081212 := proc(n) local i; i := 0 ; while r(n,i) <> r(n,i+1) do i := i+1 ; od ; RETURN(i) ; end: A081214 := proc() local a,m,h ; a :=[seq(-1,i=1..40)] ; for m from 1 to 8000000 do h := A081212(m) ; if h+1 <= nops(a) then if op(h+1,a) = -1 then a := subsop(h+1=m,a) ; print(a) ; fi ; fi ; od ; RETURN(a) ; end: A081214() ; # _R. J. Mathar_, Apr 04 2007
%t gsf[n_] := gsf[n] = For[k = n, True, k--, If[SquareFreeQ[k], Return[k]]];
%t A081210[n_] := A081210[n] = Times @@ gsf /@ Power @@@ FactorInteger[n];
%t A081212[n_] := A081212[n] = Module[{cnt = 0}, FixedPoint[(cnt++; A081210[#])&, n]; cnt - 1];
%t a[n_] := a[n] = For[m = 1, True, m++, If[A081212[m] == n, Return[m]]];
%t Table[Print[n, " ", a[n]]; a[n], {n, 0, 7}] (* _Jean-François Alcover_, Nov 27 2023 *)
%Y Cf. A081210, A081212.
%K nonn,more
%O 0,2
%A _Reinhard Zumkeller_, Mar 10 2003
%E More terms from _R. J. Mathar_, Apr 04 2007
%E One more term from _Jean-François Alcover_, Nov 27 2023
|