import numpy as np
import math


#N = 10
N = 20
#N = 45
#N = 50

#Arrays for storing results, A[i] = f_A(i), B[i] = f_AB(i), f[i] = 3(A[i]+B[i])

A = np.zeros(N+1,dtype = float)
B = np.zeros(N+1,dtype = float)
F = np.zeros(N+1,dtype = float)

A[1] = -1
B[1] = 1

#D[i] = d_i = |connected graphs|; G[i] = g_i =|G, where G and complement of G connected|

D = np.zeros(N+1,dtype = float)
G = np.zeros(N+1,dtype = float)

D[1] = 1
G[1] = 1

#Array for storing partition, P[i] = |connected components with i vertices|
P = np.zeros(N,dtype = np.int64)

#Array for storing partition number, S[i] +1 = ith stirling number 
S = np.zeros(N+1,dtype = int)


#calculate D[n], G[n] and write them in the result array
def both_connected(n):
    sum = 0
    for k in range(1,n):
        sum += math.comb(n,k)*k*D[k]*2**(math.comb((n-k),2))
    temp = 2**(math.comb(n,2))
    D[n] = (n*temp - sum)//n
    G[n] = 2*D[n] - temp


def partition(n, m, i, k):

# n: total vertex number needed to be partitioned                        
# m: number of vertices, which haven't been signed to any components 
# i: current iteration will sign number of i-vertices-components      
# k: already k components formed after previous calls of partition
    
    if m>=0 and i > 1:   
        for j in range (m//i, -1, -1):
            P[i] = j
            partition(n, m-j*i, i-1, k+j)
    else:
        if m>=0 and i == 1:
            P[i] = m
            k += m
            S[n] += 1
                  
            #calculating component combi number for specific partition P
            part = 1
            s = n
            
            for l in range (1,n):
                if P[l] != 0:
                    for p in range(P[l]):
                        part *= math.comb(s,l)                        
                        s-=l
                    part = part // math.factorial(P[l])
                               
            #calculating A
            prod_A = 1   
            #1-vertex-components don't matter
            for l in range (2,n):
                if P[l] != 0:
                    prod_A *= (2*A[l]+3*B[l])**(P[l])
            
            A[n] += part * prod_A
                      
            #calculating B
            prod_B = 1            
            #1-vertex-components don't matter
            for l in range (2,n):
                if P[l] != 0:
                    prod_B *= (A[l]+2*B[l])**(P[l])            
            B[n] += part * prod_B * G[k] 
            
    
for n in range(2,N+1):    
    both_connected(n)
    partition(n, n, n-1, 0)

    
F = 3*(A + B)
F[0] = 1
F[1] = 1

for i in range(N+1):
    f = int(F[i])
    f_a = int(A[i])
    f_b = int(B[i])
    d = int(D[i]) #connected
    g = int(G[i]) #both connected
    
    print( f,",")