%I #31 Jul 15 2024 14:47:56
%S 1,1,1,5,34,535,20848,2120098,572849763,415361983540,815590925440865,
%T 4373589784210012634,64535461714821630421106,
%U 2637732191356603658136444467,300363258297687600380548275359231
%N Number of connected oriented graphs (i.e., connected directed graphs with no bidirected edges) on n nodes.
%H Andrew Howroyd, <a href="/A086345/b086345.txt">Table of n, a(n) for n = 0..50</a>
%H Musa Demirci, Ugur Ana, and Ismail Naci Cangul, <a href="https://doi.org/10.1007/978-981-16-1402-6">Properties of Characteristic Polynomials of Oriented Graphs</a>, Proc. Int'l Conf. Adv. Math. Comp. (ICAMC 2020) Springer, see p. 61.
%H Chathura Kankanamge, <a href="http://hdl.handle.net/10012/14051">Multiple Continuous Subgraph Query Optimization Using Delta Subgraph Queries</a>, Master Maths Thesis, Univ. of Waterloo, 2018.
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/OrientedGraph.html">Oriented Graph</a>
%F Inverse Euler transform of A001174. - _Andrew Howroyd_, Nov 03 2017
%t permcount[v_] := Module[{m = 1, s = 0, k = 0, t}, For[i = 1, i <= Length[v], i++, t = v[[i]]; k = If[i > 1 && t == v[[i - 1]], k + 1, 1]; m *= t*k; s += t]; s!/m];
%t edges[v_] := Sum[GCD[v[[i]], v[[j]]], {i, 2, Length[v]}, {j, 1, i - 1}] + Total @ Quotient[v - 1, 2];
%t a1174[n_] := Module[{s = 0}, Do[s += permcount[p]*3^edges[p], {p, IntegerPartitions[n]}]; s/n!];
%t b = Array[a1174, 15];
%t mob[m_, n_] := If[Mod[m, n] == 0, MoebiusMu[m/n], 0];
%t EULERi[b_] := Module[{a, c, i, d}, c = {}; For[i = 1, i <= Length[b], i++, c = Append[c, i*b[[i]] - Sum[c[[d]]*b[[i - d]], {d, 1, i - 1}]]]; a = {}; For[i = 1, i <= Length[b], i++, a = Append[a, (1/i)*Sum[mob[i, d]*c[[d]], {d, 1, i}]]]; Return[a]];
%t A086345 = EULERi[b] (* _Jean-François Alcover_, Jul 06 2018, after _Andrew Howroyd_ *)
%o (Python)
%o from functools import lru_cache
%o from itertools import combinations
%o from math import prod, factorial, gcd
%o from fractions import Fraction
%o from sympy.utilities.iterables import partitions
%o from sympy import mobius, divisors
%o def A086345(n):
%o @lru_cache(maxsize=None)
%o def b(n): return int(sum(Fraction(3**(sum(p[r]*p[s]*gcd(r,s) for r,s in combinations(p.keys(),2))+sum((q-1>>1)*r+(q*r*(r-1)>>1) for q, r in p.items())),prod(q**r*factorial(r) for q, r in p.items())) for p in partitions(n)))
%o @lru_cache(maxsize=None)
%o def c(n): return n*b(n)-sum(c(k)*b(n-k) for k in range(1,n))
%o return sum(mobius(d)*c(n//d) for d in divisors(n,generator=True))//n if n else 1 # _Chai Wah Wu_, Jul 15 2024
%Y Cf. A054941 (labeled case), A001174, A281446 (multisets).
%K nonn
%O 0,4
%A _Eric W. Weisstein_, Jul 16 2003
%E More terms from _Vladeta Jovovic_, Jul 19 2003
%E a(0)=1 prepended by _Andrew Howroyd_, Sep 09 2018