(*^ ::[ Information = "This is a Mathematica Notebook file. It contains ASCII text, and can be transferred by email, ftp, or other text-file transfer utility. It should be read or edited using a copy of Mathematica or MathReader. If you received this as email, use your mail application or copy/paste to save everything from the line containing (*^ down to the line containing ^*) into a plain text file. On some systems you may have to give the file a name ending with ".ma" to allow Mathematica to recognize it as a Notebook. The line below identifies what version of Mathematica created this file, but it can be opened using any other version as well."; FrontEndVersion = "Macintosh Mathematica Notebook Front End Version 2.2"; MacintoshStandardFontEncoding; fontset = title, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, e8, 24, "Times"; fontset = subtitle, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, e6, 18, "Times"; fontset = subsubtitle, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeTitle, center, M7, italic, e6, 14, "Times"; fontset = section, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, grayBox, M22, bold, a20, 18, "Times"; fontset = subsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, blackBox, M19, bold, a15, 14, "Times"; fontset = subsubsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, whiteBox, M18, bold, a12, 12, "Times"; fontset = text, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, B65535, 14, "Times"; fontset = smalltext, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, 10, "Times"; fontset = input, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeInput, M42, N23, bold, L-5, 12, "Courier"; fontset = output, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L-5, 12, "Courier"; fontset = message, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, R65535, L-5, 12, "Courier"; fontset = print, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L-5, 12, "Courier"; fontset = info, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, B65535, L-5, 12, "Courier"; fontset = postscript, PostScript, formatAsPostScript, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeGraphics, M7, l34, o1, w167, h167, 12, "Courier"; fontset = name, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, italic, R7192, G22316, B34762, 10, "Geneva"; fontset = header, inactive, noKeepOnOnePage, preserveAspect, M7, 12, "Times"; fontset = leftheader, inactive, L2, 12, "Times"; fontset = footer, inactive, noKeepOnOnePage, preserveAspect, center, M7, 12, "Times"; fontset = leftfooter, inactive, L2, 12, "Times"; fontset = help, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, 10, "Times"; fontset = clipboard, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, 12, "Times"; fontset = completions, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, 12, "Times"; fontset = special1, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, underline, B65535, r0, b0, 60, "Times"; fontset = special2, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, R65535, r58981, g58981, b58981, 144, "Times"; fontset = special3, inactive, nohscroll, noKeepOnOnePage, preserveAspect, center, M7, bold, R65535, G65535, r43689, g43689, b43689, 18, "Apple Chancery"; fontset = special4, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, R65535, B65535, 14, "Times"; fontset = special5, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, R65535, G65533, B12818, r58187, g20165, b20850, 18, "B Courier Bold"; paletteColors = 128; automaticGrouping; currentKernel; ] :[font = input; initialization; preserveAspect] *) (**************************************************************************) (* *) (* Algorithm DecomposeBinomial *) (* *) (* Version of April 17 2000 *) (* programmed by David Callan *) (* *) (* This file contains a Mathematica implementation of the algorithm *) (* Decomposebinomial which finds an integer k (when one exists) *) (* such that k times (an+b)choose(cn+d) divided by the *) (* polynomials specified by list1 and list2, yields a *) (* sequence of integers. Here list1 is a sequence (possibly empty) *) (* of positive integers specifying terms in {an+b,an+b-1,an+b-2,...} *) (* and list2 is likewise for {an+b-(cn+d)+1,an+b-(cn+d)+2,...}. *) (* *) (**************************************************************************) BeginPackage["DecomposeBinomial`"] DecomposeBinomial::usage = "DecomposeBinomial[{a,b,c,d},list1,list2] or decomp[{a,b,c,d},list1,list2] finds an integer k (if one exists) such that k times (an+b)choose(cn+d) divided by the polynomials specified by list1 and list2 yields a sequence of integers. Here list1 is a sequence (possibly empty) of positive integers specifying terms in {an+b,an+b-1,an+b-2,...} and list2 is likewise for {an+b-(cn+d)+1,an+b-(cn+d)+2,...}. It also returns a Certificate of Integrality in the form of an integral linear combination of binomial coefficients that yields the sequence. For example, DecomposeBinomial[{2,0,1,-3},{1},{}] returns: A suitable integer multiplier is 6. Your input sequence, (2 n + 0)choose( n - 3)/(2 n), when multiplied by 6, begins (starting at n = 3) {1, 6, 27, 110, 429, 1638, 6188} and is equal to the following integral linear combination of binomial coefficients: -1 (2n-1)choose(n-4) +1 (2n-1)choose(n-3) Note the multiplier returned need not be minimal, although it often is."; Integralizer::usage = "Integralizer[{a,b,c,d},list1,list2] finds an integer k (if one exists) such that k times (an+b)choose(cn+d) divided by the polynomials specified by list1 and list2 yields a sequence of integers. Here list1 is a sequence (possibly empty) of positive integers specifying terms in {an+b,an+b-1,an+b-2,...} and list2 is likewise for {an+b-(cn+d)+1,an+b-(cn+d)+2,...}." decomp::usage = "decomp is an abbreviation for DecomposeBinomial."; Begin["`Private`"] DecomposeBinomial[list0_,list1_,list2_]:= Module[{a,b,c,d,e,f,u,v,k,n,p,inputfactors,ip,top,bot,numfactors, denfactors,cancelindices,dropindices,lowerusedindices, upperusedindices,degree,cancelfactors,lhs,lhstest,rhs,cc, variables,equations,solution,fractions,magicnumber, integermults,aa,bb,ck,dd,s1,start}, ( max[list_]:=If[list=={},0,Max[list]]; prod[list_]:=Product[list[[i]],{i,Length[list]}]; prod[list_,j_]:=Product[list[[i]],{i,j}]; coeff[poly_, i_] := If[i == 0, poly /. n -> 0, Coefficient[poly, n^i]]; lcm[list_] := Fold[LCM, 1, list]; nu[j_]:=NumberForm[j,NumberSigns->{"-","+"}]; nu0[j_]:=If[j==0,"",NumberForm[j,NumberSigns->{"-","+"}]]; drop1[j_]:=If[j==1,"",j]; {a,b,c,d}=list0; e=a-c;f=b-d; u=max[list1];v=max[list2]; inputfactors=Join[a n+b+1-list1,e n +f +list2]; top=Table[e n+f+v+1-j,{j,u+v}]; bot=Table[c n+d-u-v+i,{i,u+v}]; numfactors=Join[Table[a n+b+1-i,{i,u}],Table[e n+f+j,{j,v}]]; denfactors=bot; cancelindices=dropindices={}; (* Find factors to be dropped *) k=a/c-1; Do[(test=Flatten[Position[top,Expand[k bot[[i]]],1]]; If[test!={} && test[[1]]<=i,(PrependTo[dropindices,i]; AppendTo[cancelindices,i])]), {i,u+v}]; lowerusedindices=Complement[Range[u+v],cancelindices]; upperusedindices=Complement[Range[u+v],dropindices]; degree=Length[upperusedindices]; cancelfactors=denfactors[[cancelindices]]; lhstest=Simplify[prod[numfactors]/ (prod[inputfactors]prod[cancelfactors])]; If[PolynomialQ[lhstest,n]==False||PolynomialQ[rhs,n]==False, (Print["No such k exists."];Abort[])]; lhs=Collect[lhstest,n]; lowerusedfactors=bot[[lowerusedindices]]; p[0]=prod[lowerusedfactors]; p[j_]:=cc[j] Simplify[p[0]prod[top,j]/prod[bot,j]]; rhs=Collect[cc[0] p[0]+Apply[Plus, Map[p,upperusedindices]],{n}]; variables=Join[{cc[0]},Map[cc,upperusedindices]]; equations=Table[coeff[lhs,i]==coeff[rhs,i],{i,0,degree}]; solution=Flatten[Solve[equations,variables]]; If[solution=={},(Print["No such k exists!"];Abort[])]; fractions=variables/.solution; magicnumber=lcm[Denominator[fractions]]; integermults=magicnumber fractions; ip=inputfactors; s1=If[ip=={},-Infinity,1+Floor[Max[Table[(-ip[[i]]/.n->0)/ ((ip[[i]]/.n->1)-(ip[[i]]/.n->0)), {i,Length[ip]}]]]]; start=Ceiling[Max[s1,-d/c]]; Print["A suitable integer multiplier is ",magicnumber,". "]; Print["Your input sequence, "]; Print["(",HoldForm[aa n+bb]/.{aa->a,bb->b,n->"n"},")choose(", HoldForm[ck n +dd]/.{ck->drop1[c],dd->d,n->"n"},")/(", Block[{i},Product[inputfactors[[i]]/.n->"n", {i,Length[inputfactors]}]],"),"]; Print["when multiplied by ",magicnumber, ", begins (starting at n = ",start,")"]; Print[Block[{n},Table[magicnumber Binomial[a n+b,c n+d]/ prod[inputfactors],{n,start,start+6}]]]; Print["and is equal to the following integral linear"]; Print["combination of binomial coefficients: "]; Print[""]; Print[nu[integermults[[1]]]," (",a,"n",nu0[b-u],")choose(", drop1[c],"n",nu0[d-u-v],")"]; Do[If[integermults[[i+1]]!=0, Print[nu[integermults[[i+1]]]," (",a,"n",nu0[b-u], ")choose(",drop1[c],"n",nu0[d-u-v+upperusedindices[[i]]],")"]], {i,degree}] ) ]; decomp[list0_,list1_,list2_]:=DecomposeBinomial[list0,list1,list2]; Integralizer[list0_,list1_,list2_]:= Module[{a,b,c,d,e,f,u,v,k,n,p,inputfactors,top,bot,numfactors, denfactors,cancelindices,dropindices,lowerusedindices, upperusedindices,degree,cancelfactors,lhs,lhstest,rhs,cc, variables,equations,solution,fractions,magicnumber}, ( max[list_]:=If[list=={},0,Max[list]]; prod[list_]:=Product[list[[i]],{i,Length[list]}]; prod[list_,j_]:=Product[list[[i]],{i,j}]; coeff[poly_, i_] := If[i == 0, poly /. n -> 0, Coefficient[poly, n^i]]; lcm[list_] := Fold[LCM, 1, list]; nu[j_]:=NumberForm[j,NumberSigns->{"-","+"}]; nu0[j_]:=If[j==0,"",NumberForm[j,NumberSigns->{"-","+"}]]; drop1[j_]:=If[j==1,"",j]; {a,b,c,d}=list0; e=a-c;f=b-d; u=max[list1];v=max[list2]; inputfactors=Join[a n+b+1-list1,e n +f +list2]; top=Table[e n+f+v+1-j,{j,u+v}]; bot=Table[c n+d-u-v+i,{i,u+v}]; numfactors=Join[Table[a n+b+1-i,{i,u}],Table[e n+f+j,{j,v}]]; denfactors=bot; cancelindices=dropindices={}; (* Find factors to be dropped *) k=a/c-1; Do[(test=Flatten[Position[top,Expand[k bot[[i]]],1]]; If[test!={} && test[[1]]<=i,(PrependTo[dropindices,i]; AppendTo[cancelindices,i])]), {i,u+v}]; lowerusedindices=Complement[Range[u+v],cancelindices]; upperusedindices=Complement[Range[u+v],dropindices]; degree=Length[upperusedindices]; cancelfactors=denfactors[[cancelindices]]; lhstest=Simplify[prod[numfactors]/ (prod[inputfactors]prod[cancelfactors])]; If[PolynomialQ[lhstest,n]==False||PolynomialQ[rhs,n]==False, (Print["No such k exists."];Abort[])]; lhs=Collect[lhstest,n]; lowerusedfactors=bot[[lowerusedindices]]; p[0]=prod[lowerusedfactors]; p[j_]:=cc[j] Simplify[p[0]prod[top,j]/prod[bot,j]]; rhs=Collect[cc[0] p[0]+Apply[Plus, Map[p,upperusedindices]],{n}]; variables=Join[{cc[0]},Map[cc,upperusedindices]]; equations=Table[coeff[lhs,i]==coeff[rhs,i],{i,0,degree}]; solution=Flatten[Solve[equations,variables]]; If[solution=={},(Print["No such k exists!"];Abort[])]; fractions=variables/.solution; magicnumber=lcm[Denominator[fractions]] )] End[] EndPackage[] (* ;[s] 3:0,0;6030,1;8065,0;8129,-1; 2:2,13,10,Courier,1,12,0,0,0;1,13,10,Courier,1,12,0,0,65535; ^*)