# restart; read `\G:/My Drive/Aek/SchmidtNumbers
# /SchimdtAISol.txt`;

# First version: October 27,2026

with(LinearAlgebra):

# Solution to conjecture4.1 of my paper:
# "q-analog of Schmidt's problem"

####################################
# Section1: Numerical values of
# Schimdt and Conj4.1
#
# Functions:
# Schimdt(r,N), SolSchimdt(j,r),
# GZ(k,r), SolGZ(k,j,r),
# SolveS(k,j,r),
#
# Conj41(d,k,r),
# SolD0(k,r), D0(k,j,r),
# SolD1(k,r), D1(k,j,r),
#
####################################

# Solve for c_k^(r) of Schimdt:
# sum_k binom(n,k)^r*binom(n+k,k)^r = 
# sum_j c[j]*binom(n,j)*binom(n+j,j).
# Remark: refer to my paper:
# "A Simple Proof of Schmidt's Conjecture"
# Try:
# Schimdt(5,15);
Schimdt := proc(r,N) option remember;
local j,n,c,eq,var;

eq := {seq(add(binomial(n,k)^r*binomial(n+k,k)^r,k=0..n) = 
add( c[j]*binomial(n,j)*binomial(n+j,j),j=0..n), n=0..N)};

var := {seq(c[j],j=0..N)};
var := solve(eq,var);
[seq(subs(var,c[j]),j=0..N)];
end:

# Solution to c_k^(r)
# Try:
# Schimdt(3,10);
# seq(SolSchimdt(j,3),j=0..10)];
SolSchimdt := proc(j,r) option remember;
local k;

add(SolGZ(k,j,r),k=ceil(j/r)..j);
end:

# Solve for basis a_j of
# binom(n,k)^r*binom(n+k,k)^r = 
# sum_j a[j]*binom(n,j)*binom(n+j,j).
# Try:
# GZ(2,5);
GZ := proc(k,r) option remember;
local j,n,N,a,eq,var;

N := k*r+k;
eq := {seq(binomial(n,k)^r*binomial(n+k,k)^r = 
add( a[j]*binomial(n,j)*binomial(n+j,j),j=0..n), n=0..N)};

var := {seq(a[j],j=0..N)};
var := solve(eq,var);
[seq(subs(var,a[j]),j=0..N)];
end:

# Output: sequence a[j] which is solution 
# of Guo and Zeng with k and r.
# Try:
# [seq(SolGZ(k,j,r),j=0..k*(r+1))];
SolGZ := proc(k,j,r) option remember;
local i;

if j<k or j>k*r then
	return(0);
elif r=1 then
	return(1);
fi:	

add( binomial(k+i,i)*binomial(k,j-i)*binomial(j,k)*
SolGZ(k,i,r-1),i=0..j);
end:

# Solve for S(k,j,i) in equation 2 of the talk
# Try:
# SolveS(2,3);
SolveS := proc(k,j) option remember;
local i,r,S,N,eq,var;

N:= j+2;
eq := {seq(SolGZ(k,j,r+1)=add(S[i]*SolGZ(k,i,r),i=0..N),r=1..N+3)}; 
var := {seq(S[i],i=0..N)};
print(eq);
var := solve(eq,var);
[seq(subs(var,S[i]),i=0..N)];
end:


##################################################

# Now try my conjecture 4.1 from my paper:
# "q-analog of Schmidt's problem"

# Try:
# Conj41(0,4,5);
# Conj41(1,4,5);
Conj41 := proc(d,k,r) option remember;
local j,n,N,eq,var;

N := k*r+k;
eq := {seq(binomial(n+d*k,k)^r = 
add( a[j]*binomial(n+d*j,j),j=0..N), n=0..N)};

var := {seq(a[j],j=0..N)};
var := solve(eq,var);
[seq(subs(var,a[j]),j=0..N)];

end:

# Output: sequence a[j] which is solution 
# of Conj4.1 with d=0 with k and r.
# Try:
# Conj41(0,4,5);
# SolD0(4,5);
SolD0 := proc(k,r) option remember;
local j,N;

N := k*r+k;
[seq(D0(k,j,r),j=0..N)];
end:

D0 := proc(k,j,r) option remember;
local i;

if r=1 then
	if j=k then
		return(1);
	else
		return(0);
	fi:	
fi:	

add( binomial(k,j-i)*binomial(j,k)*D0(k,i,r-1),i=0..j);
end:

# Output: sequence a[j] which is solution 
# of Conj4.1 with d=1 with k and r.
# Try:
# Conj41(1,4,5);
# SolD1(4,5);
SolD1 := proc(k,r) option remember;
local j,N;

N := k*r+k;
[seq(D1(k,j,r),j=0..N)];
end:

D1 := proc(k,j,r) option remember;
local i;

if r=1 then
	if j=k then
		return(1);
	else
		return(0);
	fi:	
fi:	

add( (-1)^(i+j+k)*binomial(k,j-i)*binomial(j,k)
*D1(k,i,r-1),i=0..j);
end:


####################################
# Section2: Solution to Conj4.1
# by OpenAI
#
# Functions:
# Van(d,t,n), MyW(d,t,j),
# InvB(d,t,n),
#
# BasisC(n,d,T), BasisB(n,d,T),  
#
# Eq6(n,p,q), MyS(d,k,j,i),
# MyA(d,k,j,r), 
#
####################################

# Section 1:

# Vandermonde identity
# Try:
# Van(2,3,11);
# seq(Van(d,3,11),d=0..11);
Van := proc(d,t,n) option remember;
local v,B,R;

B := binomial(n+d*t,t);
R := add( binomial(d*t,t-v)*binomial(n,v),v=0..t);

[B,B-R];
end:

# Try:
# seq(MyW(2,8,j),j=0..10);
MyW := proc(d,t,j) option remember;
local v;

if j=t then
	return(1);
fi:	

-add( binomial(d*t,t-v)*MyW(d,v,j)   ,v=0..t-1);
end:

# Check the last equation on page 1
# Try:
# seq(InvB(2,t,31),t=0..11);
InvB := proc(d,t,n) option remember;
local j,C,R;

C := binomial(n,t);
R := add( MyW(d,t,j)*binomial(n+d*j,j),j=0..t);

[C,C-R];
end:

###################################################

# Write B and C in term of matrix 
# multiplication of each other. 
# So I can see inverse matrix.

# Write forward basis: B = U*C;
# Try: 
# BasisC(111,3,11);
BasisC := proc(n,d,T) option remember;
local t,v,B,U,C,M;

B := [seq(binomial(n+d*t,t),t=0..T)];
C := Matrix([seq(binomial(n,t),t=0..T)]);
U := Matrix([seq([seq(binomial(d*t,t-v),v=0..t)],t=0..T)]);
M := U.Transpose(C);
[B, B-[seq(M[i][1],i=1..T+1)]];
end:


# Write forward basis: C = W*B,
# where W is actually U^(-1).
# Try: 
# BasisB(111,13,11);
BasisB := proc(n,d,T) option remember;
local t,v,B,U,C,N;

C := [seq(binomial(n,t),t=0..T)];
B := Matrix([seq(binomial(n+d*t,t),t=0..T)]);
U := Matrix([seq([seq(binomial(d*t,t-v),v=0..t)],t=0..T)]);
print(U^(-1));
N := U^(-1).Transpose(B);
[C, C-[seq(N[i][1],i=1..T+1)]];
end:





###################################################
# Equation 6 on page2
# seq(Eq6(n,4,7),n=7..20);
Eq6 := proc(n,p,q) option remember;
local u,L,R;

L := binomial(n,p)*binomial(n,q);
R := add( binomial(q+u,p)*binomial(p,u)*binomial(n,q+u)
,u=0..p);

[L,L-R];
end:

# S_d(k,j,i) that is defined in equation7
# Try:
# seq(MyS(d,2,4,3),d=0..10);
# seq(MyS(5,k,4,3),k=0..10);
#
# j:=15:k:=4:
# seq(MyS(0,k,j,i),i=0..20);
# seq(binomial(k,j-i)*binomial(j,k),i=0..20);
#
# j:=15:k:=4:
# seq(MyS(1,k,j,i),i=0..20);
# seq((-1)^(i+j+k)*binomial(k,j-i)*binomial(j,k),i=0..20);
MyS := proc(d,k,j,i) option remember;
local p,q,u;

add(add(add(  
binomial(d*k,k-p)*binomial(d*i,i-q)*binomial(q+u,p)
*binomial(p,u)*MyW(d,q+u,j) 
,u=0..p),q=0..i),p=0..k);
end:


# Definition of MyA as in eq(2).
# Try:
# seq(MyA(3,4,7,r),r=1..10);
# seq(Conj41(3,4,r)[8],r=1..10);
MyA := proc(d,k,j,r) option remember;
local i;

if r=1 then
	if j=k then
		return(1);
	else
		return(0);
	fi:	
fi:	

add( MyS(d,k,j,i)*MyA(d,k,i,r-1),i=0..(r-1)*k);
end:


###########################################
# Section 3: My checking to the solution
# by OpenAI
#
# Functions:
# Check1(n,i,k,d),
#
###########################################

# Check B_i(n)*B_k(n) = sum_j(S(k,j,i)*B_j(n),
# where B_k(n) = binomial(n+d*k,k).
# Try:
# seq(Check1(n,4,3,11),n=0..20);
Check1 := proc(n,i,k,d) option remember;
local j,L,R;

L := binomial(n+d*i,i)*binomial(n+d*k,k);
R := add(  MyS(d,k,j,i)*binomial(n+d*j,j),j=0..i+k);

[L,L-R];
end: