"""Calculates the number of different types of Semilabelled Binary Trees with n leaves and k labels.
    Answers given in this order. Rooted (R), Rooted using all labels (V),
    Marked (M),  Marked using all labels (VM),
    Unmarked (U), Unmarked using all labels (VU)       
    The number of times each label is used is not specified in first set.  
    Each label used at least once in second answer set
    AUTHOR: Virginia Johnson (2011-07) initial version
    EXAMPLES:  input:  T(n,k) """
def T(n,k):
    #Gets input and will return the number of trees with leaves 0-n on k labels"""
    #first section calculates the rooted binary trees R_k in documentation number of leaves varies, number of labels fixed  
    LL=[]        #stores r_n,0, r_n,1, ...r_n,k
    for p in range(k+1):
        L=[0]*(n+1)     #stores r_0,k,  r_1,k, ...r_n,k
        LL.append(L)
        for i in range(n+1):
            #"0 if no leaves"
            if i==0:
                L[i]=0
            #"p if one leaf"
            elif i==1:
                L[i]=p
            #"if number of leaves is even"
            elif (mod(i,2)==0) and (i!=0):
                L[i]=1/2*L[i/2]
                for j in range(1,i):
                    L[i]+=1/2*L[j]*L[i-j]
            else:
                for j in range(1,i):
                     L[i]+=1/2*L[j]*L[i-j]
            #Calculates Rooted semilabeled binary trees n= number of leaves,
            #"k= number of labels  Each label is used at least once.
    V=[0]*(n+1)
    for i in range(n+1):
         for j in range (0,k):
             V[i]+=(-1)^j*binomial(k,j)*LL[k-j][i]
           
              #this section  calculates the sums needed for a_n;k in documentation""" 
    BA=[]    # this holds values for smaller number of leaves0-k
    for h in range(k+1):
            
        B=[0]*(n+1)
        BA.append(B)
                                   
        for i in range(1,n+1):
            if i==0:
                B[i]=0
            else: 
                B[i]=h*LL[h][i-1]                     #adds in first term           
                for j in [0..floor(i/3)]:         #selects combinations of i,j,k,which sum to n
                    for m in [j..floor((i-j)/2)]:
                        p = i-j-m                 
                        t=[j,m,p]                
                           #t is created to determine how many elements in set to create c_i,j,l documentation""""
                        if (2*j)+p==i and len(set(t))!=1:            #adds in third term first testing for j=m 
                            B[i]+=(1/2)*LL[h][j]*LL[h][p]                    #and eliminating j=m=p which is included  in  
                            #print B[i]                                 #next if statement
                        if j+(2*m)==i:                               #this gets j=m=p and m=p all needed in third term
                            B[i]+=(1/2)*LL[h][j]*LL[h][m]                    # have now added in third term
                            #print B[i]
                        if len(set(t))==1:                          #sets the coeffecient c and adds in second term
                            c=1
                            B[i]+=1/6*c*LL[h][j]*LL[h][m]*LL[h][p]
                        elif len(set(t))==2:
                            c=3
                            B[i]+=1/6*c*LL[h][j]*LL[h][m]*LL[h][p]
                        elif len(set(t))==3:
                            c=6
                      
                            B[i]+=1/6*c*LL[h][j]*LL[h][m]*LL[h][p]               #have now completed adding in 2nd term   
                
               
                
                
    """this section calculates the numbers of Marked trees...M in documentation"""   
    MA=[]           # this holds values for smaller number of leaves0-k
    for h in range(k+1):
                        
        M=[0]*(n+1)
        MA.append(M)            #calculates the final sum
        for i in range(n+1):
            if i==0:
                M[i]=0
            elif i==1:
                M[i]=h
            elif (mod(i,3)==0) and (i!=0):
                M[i]=BA[h][i]+(1/3)*LL[h][i/3]
            else:
                M[i]=BA[h][i]
                
    #This section calculates M^* trees in documentation . Each label is used
    VM=[0]*(n+1)
    for i in range(n+1):
            for j in range (0,k):
                VM[i]+=(-1)^j*binomial(k,j)*MA[k-j][i]
                        #This section calculated unrooted binary trees....U in documentation"""
    AU=[]         # this holds values for smaller number of leaves0-k
    for h in range (k+1): 
        U=[0]*(n+1)
        AU.append(U)
   
        for i in range(n+1):
            if i==0:
                U[i]=0
            elif i==1:
                U[i]=h
            elif (mod(i,2)==0) and (i!=0): 
                U[i]=MA[h][i]-LL[h][i]+LL[h][i/2]
            else:U[i]=MA[h][i]-LL[h][i]
    #This section calculates U^* in documentation : unrooted  binary MUL trees using all k labels
    VU=[0]*(n+1)
    for i in range(n+1):
        for j in range (0,k):
            VU[i]+=(-1)^j*binomial(k,j)*AU[k-j][i]
 #__________________________
  
  
   
             #This section returns the calculated numbers"""
    
    print "Number of leaves,number of labels"
    print n,k
    print "Rooted MUL Binary Trees"    
    print L
    
    print "Rooted MUL Binary Trees using all k labels" 
    print V
    
    print "Marked  MUL Binary Trees"
    print M
    print "Marked MUL Binary Trees using all k labels"
    print VM
    print "Unrooted MUL Binary Trees"
    print U
    print "Unrooted MUL Binary Trees using all k labels"
    print VU