#OK to post homework
#Karnaa Mistry, 10/25/20, Assignment HW #14

with(combinat):

# 1. 

#CountRuns(L,k): inputs list L of numbers and outputs the number of "runs" of length k in L
CountRuns:=proc(L,k) local i,j,ct,s:

ct:=0:

for i from 1 to nops(L) do:
    if (i+k-1) > nops(L) then break: fi:
    s:={op(i..i+k-1,L)}:
    if nops(s)=1 then ct:=ct+1: fi:
od:

return(ct):
end:




# 2.

#AvePathRuns(m,n,k,N): records the number of runs and finds the sample average and standard deviation
AvePathRuns:=proc(m,n,k,N) local runs,L,sumo,i,average,sqsum:

runs:=[]:

sumo:=0:
sqsum:=0:
for i from 1 to N do:
    L:=RPath(m,n):
    runs:=[op(runs),CountRuns(L,k)]:
od:

for i from 1 to N do:
    sumo:=sumo+runs[i]:
od:

average:=sumo/N:

for i from 1 to N do:
    sqsum:=sqsum+(runs[i]-average)^2:
od:

print(evalf(average)):

print(evalf(sqrt(sqsum/(N-1)))):

end:

# I got values: (note: my procedure takes a long time to run for N=3000)
# 24.04300000,7.380016206
# 24.27300000,7.485529287
# 24.06900000,7.541785719
# 24.29766667,7.452040359
# 24.11566667,7.566419549
# 24.02933333,7.518154501
# 24.14933333,7.547362926
# 24.07700000,7.329641571
# 24.01400000,7.602526122
# 24.21733333,7.302685073

# Yes, they are quite close together




# 3.

#AveGPathRuns(m,n,k,N): records the number of runs and finds the sample average and standard deviation (Gpaths version)
AveGPathRuns:=proc(m,n,k,N) local runs,L,sumo,i,average,sqsum:

runs:=[]:

sumo:=0:
sqsum:=0:
for i from 1 to N do:
    L:=RGPath(m,n):
    runs:=[op(runs),CountRuns(L,k)]:
od:

for i from 1 to N do:
    sumo:=sumo+runs[i]:
od:

average:=sumo/N:

for i from 1 to N do:
    sqsum:=sqsum+(runs[i]-average)^2:
od:

print(evalf(average)):

print(evalf(sqrt(sqsum/(N-1)))):

end:

# I get values like 24.34066667, 7.605517353; 23.90633333, 7.382573672; 24.20933333, 7.524474652; 24.20166667, 7.677107303;
# 24.15566667, 7.545224407... etc. and others that are also similar.  It is quite similar to my previous answer.




# 4.

#CountMaximalRuns(L,k): inputs list L of numbers and outputs the number of "maximal runs" of length k in L
CountMaximalRuns:=proc(L,k) local i,j,ct,s:

ct:=0:

for i from 1 to nops(L) do:
    if (i+k-1) > nops(L) then break: fi:
    if i=1 then 
		s:={op(i..i+k-1,L)}: 
        if nops(s)=1 and L[i+k]<>L[i] then ct:=ct+1: fi:
    elif (i+k) > nops(L) then
        s:={op(i..i+k-1,L)}:
        if nops(s)=1 and L[i-1]<>L[i] then ct:=ct+1: fi:
    else
        s:={op(i..i+k-1,L)}:
        if nops(s)=1 and (L[i-1]<>L[i] and L[i+k]<>L[i]) then ct:=ct+1: fi:
    fi:
od:


return(ct):
end:


#AvePathMaximalRuns(m,n,k,N): records the number of max-runs and finds the sample average and standard deviation
AvePathMaximalRuns:=proc(m,n,k,N) local runs,L,sumo,i,average,sqsum:

runs:=[]:

sumo:=0:
sqsum:=0:
for i from 1 to N do:
    L:=RPath(m,n):
    runs:=[op(runs),CountMaximalRuns(L,k)]:
od:

for i from 1 to N do:
    sumo:=sumo+runs[i]:
od:

average:=sumo/N:

for i from 1 to N do:
    sqsum:=sqsum+(runs[i]-average)^2:
od:

print(evalf(average)):

print(evalf(sqrt(sqsum/(N-1)))):

end:


#AveGPathMaximalRuns(m,n,k,N): records the number of max-runs and finds the sample average and standard deviation (Gpaths ver.)
AveGPathMaximalRuns:=proc(m,n,k,N) local runs,L,sumo,i,average,sqsum:

runs:=[]:

sumo:=0:
sqsum:=0:
for i from 1 to N do:
    L:=RGPath(m,n):
    runs:=[op(runs),CountMaximalRuns(L,k)]:
od:

for i from 1 to N do:
    sumo:=sumo+runs[i]:
od:

average:=sumo/N:

for i from 1 to N do:
    sqsum:=sqsum+(runs[i]-average)^2:
od:

print(evalf(average)):

print(evalf(sqrt(sqsum/(N-1)))):

end:

# Values like 6.176000000, 2.347948385   6.179000000, 2.378552397  6.150666667, 2.367241890  6.204666667, 2.398338777  6.250000000, 2.394871980, etc.
# were very similar to others (Gpaths version) like: 6.132666667, 2.379836290  6.252000000, 2.353934933  6.159000000, 2.422603696, and so on




# 5.

#GoodBrother(L): inputs list of odd length that adds up to 1, outputs unique cycle shift of length 2*n with 1 appended at the end
GoodBrother:=proc(L) local i,j,H,curr,N,gb:

if nops(L) mod 2 = 0 then return("error"): fi:

H:=convert(CycShi(L),list):
N:=nops(L):

gb:=[]:

for i from 1 to N do:
    curr:=0:
    if op(N,H[i])<>1 then next: fi: #check if 1 is appended at the end
    for j from 1 to N do:
        if (curr < 0) then break: fi: #check if H[i] violates Dyck Path definition
        curr := curr + H[i][j]:
        if j = N then gb:=H[i]: fi:
    od:
od:

return(gb):

end:




# 6.

# We use the idea that if two cyclic shifts of a list L (a 'word' in {-1,1} that adds up to 1) are identical, there is some word L1 that L is equal to 
# L1 repeated k times, for some positive integer k.  So, we have something like:

# L = L1 L1 L1 ... (k times)
# Sum(L) = Sum(L1)*k
# Note that, as L is a word in {-1,1} that adds up to 1, Sum(L) = 1.
# A sum of lists, Sum(L1) must be an integer, but how can we express 1 as a product of 2 integers?  It is a fact that for any two integers a and b,
# if a*b = 1, a=1 and b=1.  Therefore, Sum(L) = Sum(L1)*1.  L1 is repeated once only, which means that L is its own unique cyclic shift.

# By trying to find out how two cyclic shifts of L are identical, we were able to conclude that L can only be made up of 1 copies of a smaller list,
# which really means that it can only be made up by itself (unique). (Contradiction after assuming L equals one of its cyclic shifts)