%I #12 Jun 26 2024 04:24:59
%S 1,5,15,87,201,3123,5769,148347,913761,11541123,39975849,2616723387,
%T 6227552241,230557039443,4151870901369,76980002233707,355687471142721,
%U 27886053280896963,121645100796252489,10474674957482235867,135117295282596928401,2811664555920692775603
%N a(n) = n! * Sum_{d|n} 3^(d - 1) / d!.
%H G. C. Greubel, <a href="/A336998/b336998.txt">Table of n, a(n) for n = 1..445</a>
%F E.g.f.: Sum_{k>=1} (exp(3*x^k) - 1) / 3.
%F a(p) = p! + 3^(p - 1), where p is prime.
%t Table[n! Sum[3^(d - 1)/d!, {d, Divisors[n]}], {n, 1, 22}]
%t nmax = 22; CoefficientList[Series[Sum[(Exp[3 x^k] - 1)/3, {k, 1, nmax}], {x, 0, nmax}], x] Range[0, nmax]! // Rest
%o (PARI) a(n) = n! * sumdiv(n, d, 3^(d-1)/d!); \\ _Michel Marcus_, Aug 12 2020
%o (Magma)
%o A336998:= func< n | Factorial(n)*(&+[3^(d-1)/Factorial(d): d in Divisors(n)]) >;
%o [A336998(n): n in [1..40]]; // _G. C. Greubel_, Jun 26 2024
%o (SageMath)
%o def A336998(n): return factorial(n)*sum(3^(k-1)/factorial(k) for k in (1..n) if (k).divides(n))
%o [A336998(n) for n in range(1,41)] # _G. C. Greubel_, Jun 26 2024
%Y Cf. A034730, A053486, A057625, A336997.
%K nonn
%O 1,2
%A _Ilya Gutkovskiy_, Aug 10 2020