In[]:=
(*=====================================================================*)(*FASTZETATRANSFORMOVERGF(2)ANDSETSEQUENCEGENERATOR*)(*=====================================================================*)ClearAll[FastZetaTransform,BitwiseORConvolution,ComputeSequence];​​​​(*1.FastZetaTransform(FZT)Modulo2*)​​(*Performsanin-placetransformationinO(NlogN)time*)​​FastZetaTransform[arr_List]:=Module[{a=arr,n=Length[arr],step=1,i,j},While[step<n,i=1;​​While[i<=n,For[j=0,j<step,j++,(*a[i+j+step]=a[i+j+step]\[OPlus]a[i+j]*)a[[i+j+step]]=BitXor[a[[i+j+step]],a[[i+j]]]];​​i+=2*step;];​​step*=2;];​​a];​​​​(*2.BitwiseORConvolutionusingFZT*)​​(*Convolutiontheorem:InverseFZT(FZT(A)·FZT(B))*)​​BitwiseORConvolution[AList_,BList_]:=FastZetaTransform[BitAnd[FastZetaTransform[AList],FastZetaTransform[BList]]];​​​​(*3.CoreSequenceGenerator*)​​(*ComputesthebooleanarrayforS_muptoamaximumintegervalueR*)​​ComputeSequence[m_Integer,R_Integer]:=Module[{k,NSize,A,B,C,nextArr,i},(*Padarraystothenextpowerof2forFZTcompatibility*)k=Ceiling[Log2[R+1]];​​NSize=2^k;​​(*BaseCases:S_1={0},S_2={1}.NoteMathematicais1-indexed.*)A=ConstantArray[0,NSize];A[[1]]=1;​​B=ConstantArray[0,NSize];B[[2]]=1;​​If[m==1,Return[A]];​​If[m==2,Return[B]];​​(*Evaluaterecurrenceuptostepm*)For[i=3,i<=m,i++,(*ComputeA\[OTimes]B*)C=BitwiseORConvolution[A,B];​​(*ComputeA\[OPlus]B\[OPlus]C*)nextArr=BitXor[A,BitXor[B,C]];​​(*ApplyInc()byright-shiftingthebooleanarray*)nextArr=Prepend[Drop[nextArr,-1],0];​​(*Updatehistory*)A=B;​​B=nextArr;];​​B];​​​​(*=====================================================================*)​​(*HELPERFUNCTIONSFORUSAGE*)​​(*=====================================================================*)​​​​(*PredicateFunction:ReturnsTrueifr∈S_m*)​​IsInSetQ[m_Integer,r_Integer]:=Module[{arr},If[r<0,Return[False]];​​arr=ComputeSequence[m,r];​​arr[[r+1]]==1];​​​​(*SetExtractionFunction:ReturnstheactuallistofelementsinS_muptoMaxR*)​​GetSetElements[m_Integer,MaxR_Integer]:=Flatten[Position[ComputeSequence[m,MaxR],1]]-1;​​​​(*ExampleExecution*)​​Print["S_8 up to 20: ",GetSetElements[8,20]];​​Print["Is 16 in S_8? ",IsInSetQ[8,16]];​​Print["Is 13 in S_8? ",IsInSetQ[8,13]];​​Print["Is 500 in S_1000? ",IsInSetQ[1000,500]];​​Print["Is 501 in S_1000? ",IsInSetQ[1000,501]];
S_8 up to 20: {4,6,8,9,12,14,16}
Is 16 in S_8? True
Is 13 in S_8? False
Is 500 in S_1000? True
Is 501 in S_1000? False