(* ::Package:: *)

(*Package: BinomialDiophantinePerturbations
  Author: Luis A. Medina
  email: luis.medina17@upr.edu

  Mathematica Version: 11.0
*)


BeginPackage["BinomialDiophantinePerturbations`"];


(* ::Text:: *)
(*Message.*)


Print["BinomialDiophantinePerturbations: This program was written by 

Luis A. Medina 
University of Puerto Rico, Rio Piedras
San Juan, PR

Written: December 21, 2016.
Last Update: December 21, 2016.

Please report bugs to: luis.medina17@upr.edu

The most current version of this  package is available 
from http://emmy.uprrp.edu/lmedina/

Below you can see the list of functions included in this
package.  For information about a particular function, 
simply click on it or type ?NameOfFunction, e.g. ?LinRec.
Type ?RecExpSumsSP`* if you want to see again the list of 
functions included in the package."]


(* ::Text:: *)
(*Functions visible to the user.*)


S::usage="S[P]: Given a Binary polynomial, returns the exponential sums associated to P."


ANF::usage="ANF[valuedvector]: returns the Algebraic Normal Form of a valued vector."
Perturbation::usage="Perturbation[n,k, F[x]]: returns the exponential sum S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, k\)]\)+F(X)).
Perturbation[n,{\!\(\*SubscriptBox[\(k\), \(1\)]\),...,\!\(\*SubscriptBox[\(k\), \(s\)]\)},F[x]]: returns the exponential sum S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(1\)]\)]\)+...+\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(s\)]\)]\)+F(X))"
PerturbationTuple::usage="PerturbationTuple[n,k,F[x]]: returns a tuple related to S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, k\)]\)+F(X)).
PerturbationTuple[n,{\!\(\*SubscriptBox[\(k\), \(1\)]\),...,\!\(\*SubscriptBox[\(k\), \(s\)]\)},F[x]]: returns a tuple related to S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(1\)]\)]\)+...+\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(s\)]\)]\)+F(X))"
PerturbationSolution::usage="PerturbationSolution[n,k,F[x]]: If S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, k\)]\)+F(X)) is balanced, then it returns the solution to the Diophantine equation related to it.
Otherwise, returns $Failed.  PerturbationTuple[n,{\!\(\*SubscriptBox[\(k\), \(1\)]\),...,\!\(\*SubscriptBox[\(k\), \(s\)]\)},F[x]]: If S(\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(1\)]\)]\)+...+\!\(\*SubscriptBox[\(\[Sigma]\), \(n, \*SubscriptBox[\(k\), \(s\)]\)]\)+F(X)) is balanced, then
it returns returns the solution to the Diophantine equation related to it."
Omega::usage="Omega[n,j]: returns all solutions to the Diophantine equation
\!\(\*UnderoverscriptBox[\(\[Sum]\), \(l = 0\), \(n\)]\)Binomial[n,l] \!\(\*SubscriptBox[\(x\), \(l\)]\)=0,
where \!\(\*SubscriptBox[\(x\), \(l\)]\) is an element of the set \!\(\*SubscriptBox[\(\[CapitalGamma]\), \(j\)]\)={0,\[PlusMinus]1,\[PlusMinus]2,\[PlusMinus]3,...,\!\(\*SuperscriptBox[\(2\), \(j - 1\)]\)}."
OmegaQuotient::usage="OmegaQuotient[n,j]: returns the set of all equivalent classes on the set \[CapitalOmega](n,j)."
EquivalentTuple::usage="EquivalentTuple[tuple]: returns the equivalence class of the tuple."
SporadicPerturbationQ::usage="SporadicPerturbationQ[n,k,F]: If the pertubation \!\(\*SubscriptBox[\(\[Sigma]\), \(n, k\)]\)+F(X) is balanced, it returns True if it is sporadic balanced and False otherwise.
SporadicPerturbationQ[n_,{\!\(\*SubscriptBox[\(k\), \(1\)]\),...,\!\(\*SubscriptBox[\(k\), \(s\)]\)},F]: If the pertubation \!\(\*SubscriptBox[\(\[Sigma]\), \(n, \([\*SubscriptBox[\(k\), \(1\)],  ... , \*SubscriptBox[\(k\), \(s\)]\)\)]\)]+F(X) is balanced, it returns True if it is sporadic balanced and False otherwise."

SymmetricANF::usage="SymmetricANF[tuple]: it assumes that tuple comes from a symmetric Boolean function and that the information on tuple is compressed."
X::usage="Represents a variable to be used by our polynomials."


Begin["Private`"];


S[P_]:=Module[{vars,n},
vars=Sort[Variables[P]];
n=Length[vars];
Total[(-1)^(P/.#&/@(Thread[vars->#]&/@Tuples[{0,1},n]))]
]


ReducedK[list_]:=Module[{ks=list,count=0,L,reducedks={},cuantos,singularks,thek},
 ks=Select[ks,#>=0&];
 singularks=Union[ks];
 L=Length[singularks];
 While[(count=count+1)<=L,
 thek=singularks[[count]];
 cuantos=Count[ks,thek];
  If[OddQ[cuantos],
    reducedks=Append[reducedks,thek]
  ]
 ];
reducedks
]

BS[n_,k_]:=Module[{myk},
 myk=ReducedK[k];
 If[Min[myk]==0,
  -BinomS[n,DeleteCases[myk,0]],
   BinomS[n,myk]
 ]
]

Perturbation[n_Integer,k_Integer,F_?PolynomialQ]/;k>=1:=Module[{CC,vars,tuples,j},
vars=Variables[F];
j=Length[vars];
tuples=Tuples[{0,1},j];
CC=Table[Total[(-1)^(F/.Thread[vars->#])&/@Select[tuples,Total[#]==m&]],{m,0,j}];
Sum[CC[[m+1]]BS[n-j,DeleteCases[Table[
If[Mod[Binomial[m,i],2]!=0,
(k-i),
{}
]
,{i,0,m}],{}]],{m,0,j}]
]

Perturbation[n_Integer,k_List,F_?PolynomialQ]:=Module[{CC,vars,tuples,j},
vars=Variables[F];
j=Length[vars];
tuples=Tuples[{0,1},j];
CC=Table[Total[(-1)^(F/.Thread[vars->#])&/@Select[tuples,Total[#]==m&]],{m,0,j}];
Sum[CC[[m+1]]BS[n-j,DeleteCases[Flatten[Table[
If[Mod[Binomial[m,i],2]!=0,
(k-i),
{}
]
,{i,0,m}]],{}]],{m,0,j}]
]

PerturbationTable[k_Integer,F_?PolynomialQ,min_,max_]/;k>=1:=Module[{CC,vars,tuples,j},
vars=Variables[F];
j=Length[vars];
tuples=Tuples[{0,1},j];
CC=Table[Total[(-1)^(F/.Thread[vars->#])&/@Select[tuples,Total[#]==m&]],{m,0,j}];
Table[Sum[CC[[m+1]]BS[n-j,DeleteCases[Table[
If[Mod[Binomial[m,i],2]!=0,
(k-i),
{}
]
,{i,0,m}],{}]],{m,0,j}],{n,min,max}]
]

PerturbationTable[k_List,F_?PolynomialQ,min_,max_]:=Module[{CC,vars,tuples,j},
vars=Variables[F];
j=Length[vars];
tuples=Tuples[{0,1},j];
CC=Table[Total[(-1)^(F/.Thread[vars->#])&/@Select[tuples,Total[#]==m&]],{m,0,j}];
Table[Sum[CC[[m+1]]BS[n-j,DeleteCases[Flatten[Table[
If[Mod[Binomial[m,i],2]!=0,
(k-i),
{}
]
,{i,0,m}]],{}]],{m,0,j}],{n,min,max}]
]

PerturbationCoeffs[F_]:=Module[{vars,tuples,j},
vars=Variables[F];
j=Length[vars];
tuples=Tuples[{0,1},j];
Table[Total[(-1)^(F/.Thread[vars->#])&/@Select[tuples,Total[#]==m&]],{m,0,j}]
]
signfirstnonzero[list_]:=Module[{bool=True,L=Length[list],c=0},
While[bool&&((c=c+1)<=L),
bool=(list[[c]]==0)
];
If[!bool,list[[c]]/Abs[list[[c]]],1]
]

EquivalentTuple[sol_]:=Module[{L,gcd,misol},
L=Length[sol];
If[EvenQ[L],
misol=Table[sol[[i]]+sol[[L+1-i]],{i,1,L/2}],
misol=Append[Table[sol[[i]]+sol[[L+1-i]],{i,1,(L-1)/2}],sol[[(L+1)/2]]]
];
gcd=GCD@@misol;
If[gcd!=0,misol=misol/gcd];
misol/signfirstnonzero[misol]
]

coef[k_,poly_?PolynomialQ,J_]:=Module[{coeffs=PerturbationCoeffs[poly],j},
j=Length[coeffs]-1;
coeffs.Table[(-1)^Sum[Binomial[i,m]Binomial[J,k-m],{m,0,i}],{i,0,j}]
]

coef[k_List,poly_?PolynomialQ,J_]:=Module[{coeffs=PerturbationCoeffs[poly],j},
j=Length[coeffs]-1;
coeffs.Table[(-1)^Sum[Binomial[i,m]Total[(Binomial[J,#-m]&/@k)],{m,0,i}],{i,0,j}]
]

PerturbationTuple[M_,k_,poly_?PolynomialQ]:=Module[{j=Length[Variables[poly]]},
Table[coef[k,poly,J],{J,0,M-j}]/2
]

PerturbationSolution[M_,k_,poly_?PolynomialQ]:=Module[{j=Length[Variables[poly]],ptuple},
ptuple=PerturbationTuple[M,k,poly];
If[ptuple.Table[Binomial[M-j,i],{i,0,M-j}]==0,
ptuple,
$Failed]
]

OmegaQuotient[n_,j_]:=
If[OddQ[n],
Union[If[#!=ConstantArray[0,(n+1)/2],signfirstnonzero[#] #/GCD@@#,#]&/@Select[Tuples[Range[-2^j,2^j],(n+1)/2],#.Table[Binomial[n,i],{i,0,(n-1)/2}]==0&]],
Union[If[#!=ConstantArray[0,n/2+1],signfirstnonzero[#] #/GCD@@#,#]&/@Select[Select[Tuples[Range[-2^j,2^j],n/2+1],#.Table[Binomial[n,i],{i,0,n/2}]==0&],Abs[Last[#]]<=2^(j-1)&]]
]

Omega[n_,j_]:=Select[Tuples[Range[-2^(j-1),2^(j-1)],n+1],#.Table[Binomial[n,i],{i,0,n}]==0&]
mitabla[k_,poly_,M_]:=Flatten[Position[Table[Perturbation[n,k,poly],{n,k,M}],0]+k-1]
s[m_,l_,k_]:=Module[{r=Floor[Log[2,k]]+1},
1/2^r Sum[(-1)^Sum[Binomial[m,i]Binomial[a,k-i],{i,0,m}] Exp[a (\[Pi] I)/2^(r-1) l],{a,0,2^r-1}]
]

PerturbationClosedForm[n_,k_,F_?PolynomialQ]:=Module[{j=Length[Variables[F]],CC,r=Floor[Log[2,k]]+1},
CC=PerturbationCoeffs[F];
Sum[CC[[m]]Sum[s[m-1,l,k](1+Exp[-((\[Pi] I l)/2^(r-1))])^(n-j),{l,0,2^r-1}],{m,1,j+1}]
]

TransformToTruth[tuple_]:=Module[{n=Length[tuple]-1,count=-1,list={}},
While[(count=count+1)<=n,
list=Join[list,ConstantArray[tuple[[count+1]],Binomial[n,count]]]
];
list
]

ANF[valuedvec_]/;((Length[valuedvec]==2^IntegerExponent[Length[valuedvec],2])&&Union[valuedvec]=={0,1}):=Module[{cuantos,pos1,pos0,tuples,pos,c=0,pol,n,alltuples},
pos=Flatten[Position[valuedvec,1]];
cuantos=Length[pos];
n=Floor[Log[2,Length[valuedvec]]];
alltuples=Sort[Tuples[{0,1},n],Total[#1]<=Total[#2]&];
tuples=alltuples[[pos]];
PolynomialMod[Total[(Times@@(1+X/@Flatten[Position[#,0]])) (Times@@(X/@Flatten[Position[#,1]]))&/@tuples],2]
]

SporadicPerturbationQ[n_Integer,k_Integer,F_?PolynomialQ]/;(k>=0):=Module[{tuple,equivtuple,tupleL},
  tuple=PerturbationSolution[n,k,F];
  If[tuple===$Failed,
    $Failed,
    equivtuple=EquivalentTuple[tuple];
    tupleL=Length[tuple];
    If[EvenQ[tupleL],
      equivtuple!=ConstantArray[0,tupleL/2],
      !(equivtuple==ConstantArray[0,Ceiling[tupleL/2]]||equivtuple==Append[Table[2(-1)^i,{i,0,Floor[tupleL/2]-1}],(-1)^((tupleL-1)/2)])
    ]
  ]
]

SporadicPerturbationQ[n_Integer,k_List,F_?PolynomialQ]:=Module[{tuple,equivtuple,tupleL},
  tuple=PerturbationSolution[n,k,F];
  If[tuple===$Failed,
    $Failed,
    equivtuple=EquivalentTuple[tuple];
    tupleL=Length[tuple];
    If[EvenQ[tupleL],
      equivtuple!=ConstantArray[0,tupleL/2],
      !(equivtuple==ConstantArray[0,Ceiling[tupleL/2]]||equivtuple==Append[Table[2(-1)^i,{i,0,Floor[tupleL/2]-1}],(-1)^((tupleL-1)/2)])
    ]
  ]
]


SymmetricANF[tuple_]:=ANF[TransformToTruth[tuple]];


End[];
EndPackage[];
