#OK to post homework
#Jason Saied, April 14, 2019, Assignment 21

#Problem 1

#input a functional chain called chain consisting of expressions in a variable x 
#outputs a list of functional chains corresponding to the different (multiplicative) terms in the derivative of the composition 
ListDiff:=proc(chain, x) local i, n, chainList:

	n:=nops(chain):
	chainList:=[]:
	for i from 1 to n do
		chainList:=[op(chainList),[op(1..(i-1), chain),diff(chain[i],x)]]:
	od:

end:

#input a list L of functional chains, corresponding to terms in a product
#outputs a list SumList of all the (additive) terms in the product rule, where the terms are given in the same format as L
ProductRule:=proc(L,x) local i, SumList:

	SumList:=[]:

	for i from 1 to nops(L) do
	
		#multiply the derivative of each term by the other terms, and prepare to sum them up
		SumList:=[op(SumList),[op(ListDiff(L[i],x)), op(1..(i-1), L), op((i+1)..nops(L), L)]]:
	
	od:
	
	return SumList:

end:

#input a list S in the format of SumList above and a variable x
#interpret S as a sum, and for each term use ProductRule to take the derivative, then put the results back into a larger list of terms, the derivative
TakeDer:=proc(S,x) local i, newS:

	newS:=[]:

	for i from 1 to nops(S) do
	
		newS:=[op(newS), op(ProductRule(S[i],x))]:
	
	od:
	
	return newS:

end:

#input a functional chain L, a variable x, a value x0, and a number of derivatives k
#take k derivatives, then plug in x_0
EvalFCdiffk:=proc(L,x,x0,k) local i,j,S, currentProd, currentSum, temp:

	S:=[[L]]:
	
	for i from 1 to k do
	
		S:=TakeDer(S,x):
	
	od:
	
	currentSum:=0:
	
	for i from 1 to nops(S) do
	
		currentProd:=1:
	
		for j from 1 to nops(S[i]) do 
		
			temp:=EvalFC(S[i][j],x,x0):
			temp:=temp[nops(temp)]:
		
			currentProd:=currentProd*temp:
		
		od:
		
		currentSum:= currentSum+currentProd:
	
	od:
	
	return currentSum:

end:




#Problem 2
CountDiff1:=proc(n) local i, termList, coeffSumList:
	termList:=[]:
	coeffSumList:=[]:
	for i from 2 to n do
		termList:=[op(termList), nops(diff(f(g(x)),x$i))]:
		coeffSumList:=[op(coeffSumList), SumCoeff(diff(f(g(x)),x$i))]:
	od:
	termList, coeffSumList:
end:

#CountDiff1(10) gave output [2, 3, 5, 7, 11, 15, 22, 30, 42], [2, 5, 15, 52, 203, 877, 4140, 21147, 115975].
#The number of terms appears to be the number of partitions of i. 
#As we saw in class, the sum of the coefficients appears to be the Bell numbers (the number of ways to partition a set of i labeled elements). 

#Problem 3
CountDiff2:=proc(n) local i, termList, coeffSumList:
	termList:=[]:
	coeffSumList:=[]:
	for i from 2 to n do
		termList:=[op(termList), nops(diff(f(g(h(x))),x$i))]:
		coeffSumList:=[op(coeffSumList), SumCoeff(diff(f(g(x)),x$i))]:
	od:
	termList, coeffSumList:
end:

#CountDiff2(10) gave output [3, 6, 13, 23, 44, 74, 129, 210, 345], [2, 5, 15, 52, 203, 877, 4140, 21147, 115975].
#The number of terms here is in the OEIS, but only as the number of terms in the ith derivative of a composition of 3 functions!
#Oddly, the sum of the coefficients appears to be the same as for the composition of two functions... (The Bell numbers)
#I will briefly see what happens for four functions: 

CountDiff3:=proc(n) local i, termList, coeffSumList:
	termList:=[]:
	coeffSumList:=[]:
	for i from 2 to n do
		termList:=[op(termList), nops(diff(f(g(h(l(x)))),x$i))]:
		coeffSumList:=[op(coeffSumList), SumCoeff(diff(f(g(x)),x$i))]:
	od:
	termList, coeffSumList:
end:

#CountDiff3(10) gave output [4, 10, 26, 55, 121, 237, 468, 867, 1597], [2, 5, 15, 52, 203, 877, 4140, 21147, 115975].
#The same thing happened! The number of terms is only in the OEIS as coming from the composition of four functions, and the sum of the coefficients is STILL the same Bell numbers!! Interesting. 

####from class

#C21.txt, April 11, 2019, Functional chains (and networks)
Help21:=proc(): print(`  EvalFC1(L,x,0), EvalFC(L,x,x0) , RP(d,K,x) , RC(d,K,n,x) `): 
print(`EvalFCdiff(L,x,x0), SumCoeff(A) `):
end:

#y=a*x+b #Data-Set [x1,y1], ..., [xn,yn]  L(a,b):-=Sum(a*xi+b-yi)^2,i=1..n);
#diff(L,a)=0, diff(L,b)=0   a*x_{n+1}+b
# [x1,..., xm; y]: Loss(a1, ..., am,a0):=Sum((a1*x1+...+am*xm+a0-y)^2 , data points in the set)
#[m1,m2,m3, ..., mk] m1xm2 +m2xm3+....
#CNN
#functional chain x->f1(x)->f2(f1(x))
#Def: A functional chain is a list of expressions [f1, ..., fn] in x

#EvalFC1(L,x,x0): inputs L=[f1,..., fn] expressions in x and a number (or symbol) x0
#outputs f1(f2(...fn(x0))))) 
EvalFC1:=proc(L,x,x0) local i,L1:

if nops(L)=1 then
   RETURN(subs(x=x0,L[1])):
fi:

L1:=[op(1..nops(L)-1,L)]:
 expand(subs( x=EvalFC1(L1,x,x0),L[nops(L)])):


end:

#RP(d,K,x): a random polynomial of degree d in x with coeff. between -K and K
RP:=proc(d,K,x) local ra,i:
ra:=rand(-K..K):
add(ra()*x^i,i=0..d):
end:

#RC(d,K,n,x): a random functional chain of length n
RC:=proc(d,K,n,x) local i:
[seq(RP(d,K,x),i=1..n)]:
end:

#f(g(x))'= f'(g(x))*g'(x)
#f1(f2(f3(x))'=f1'(f2(f3(x))*f2'(f3(x))*f3'(x)
#(fk(f(k-1)(f(k-2)))'= fk'(f(k-1)....) [(f(k-1)(f(k-2)....) f1'(x0)]

#x0->f1(x0)->f2(f1(x0)->.... fn(f_{n-1}(x)) ...
#EvalFC(L,x,x0): inputs L=[f1,..., fn] expressions in x and a number (or symbol) x0
#outputs [f1(x0),f2(f1(x0)), f3(f2(f1(x0)), ..., fn(....)]: 
EvalFC:=proc(L,x,x0) local i,L1:

if nops(L)=1 then
   RETURN([subs(x=x0,L[1])]):
fi:

L1:=[op(1..nops(L)-1,L)]:

L1:=EvalFC(L1,x,x0):


[op(L1),subs(x=L1[nops(L1)] , L[nops(L)])]:

end:

#EvalFCdiff(L,x,x0): inputs a functional chain L, variable x, and number or symbol x0
#outputs the derivative of EvalFC(L,x,y) w.r.t y and y=x0
#
EvalFCdiff:=proc(L,x,x0) local i,Ld,M,L0:

L0:=EvalFC(L,x,x0):

Ld:=[seq(diff(L[i],x),i=1..nops(L))]:

M:=[subs(x=x0,Ld[1])]:

for i from 2 to nops(L) do
 M:=[op(M),  subs(x=L0[i-1]  ,Ld[i])  ]:
od:

convert(M,`*`):

end:

#SumCoeff(A): the sum of the coefficient in a + expression
SumCoeff:=proc(A) local su,i, lauren:
su:=0:
if not type(A,`+`) then
  RETURN(FAIL):
fi:
for i from 1 to nops(A) do
   lauren:=op(1,op(i,A)):
    if type(lauren,integer) then
       su:=su+lauren:
       else
       su:=su+1:
     fi:
od:
su:
end: