{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "#Each calculation is associated to some portion of the preprint by\n", "#Daniel J. Katz and Courtney M. van der Linden entitled\n", "#\"Peak Sidelobe Level and Peak Crosscorrelation of Golay-Rudin-Shapiro Sequences\"\n", "#arXiv: 2108.07318v3 [cs.IT]" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-982286915957927039/62500000000000000000000\n", "6308285076880354641/1000000000000000000000000\n" ] } ], "source": [ "#Abstract: We speak of \"the positive real root of x^4-3x-6\", which we approximate as 1.842626...\n", "#We show that f(x)=x^4-3x-6 has only one real root, and that it is approximated by this.\n", "#The derivative is 4x^3-3, which is negative for x < (3/4)^(1/3) and positive for x > (3/4)^(1/3).\n", "#So x^4-3x-6 is strictly decreasing, and then strictly increasing, so it has at most 2 real roots.\n", "#Furthermore, f(0)=-6 is negative, and since f is monic and of fourth degree, this means that the\n", "#graph of y=f(x) crosses the negative x-axis once and crosses the posiitve x-axis once, so there \n", "#is precisely one positive real root. \n", "#Then we use rational arithmetic to bracket the positive real root.\n", "\n", "def hjj_poly(r):\n", " return r^4-3*r-6\n", "print(hjj_poly(1842626/(10^6)))\n", "print(hjj_poly(1842627/(10^6)))" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-148715039143273238387761/15625000000000000000000000000000000\n", "56686995311154407122213/1000000000000000000000000000000000000\n" ] } ], "source": [ "#Abstract: We speak of \"the real root of x^3+x^2-2x-4\", which we approximate as 1.658967...\n", "#We show that m(x)=x^3+x^2-2x-4 has only one real root, and that it is approximated by this.\n", "#The discriminant is -236, so there is only one real root.\n", "\n", "def our_poly(r):\n", " return r^3+r^2-2*r-4\n", "print(our_poly(1658967081916/(10^12)))\n", "print(our_poly(1658967081917/(10^12)))" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "#Theorem 1.3 is asserting that 5 alpha_0^(-4)=0.660113..., but since Theorem 1.3 is a part of Theorem 5.2,\n", "#we defer the calculation until then." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "#Theorem 1.4 is asserting that 5 alpha_0^(-3)=1.095107..., but since Theorem 1.4 is a part of Theorem 5.2,\n", "#we defer the calculation until then." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n= 0:[0, 1]\n", "\n", "n= 1:[-1, -1][1, 1]\n", "\n", "n= 2:[-3, 1][-1, 1][1, 3][3, -1]\n", "\n", "n= 3:[-5, -1][-3, -5][-1, 3][1, 1][3, 1]\n", " [5, 1]\n", "\n", "n= 4:[-11, 5][-7, 1][-5, 1][-3, 5][3, 7]\n", " [5, 3][7, 3][11, -5]\n", "\n", "n= 5:[-21, -1][-13, -9][-11, -13][-9, 3][-5, 7]\n", " [5, -3][9, 9][11, -7][13, 5][21, 1]\n", " \n", "\n", "n= 6:[-43, 13][-41, -3][-27, 13][-23, -7][-21, 9]\n", " [-11, 5][11, 7][19, 15][21, -5][27, 7]\n", " [41, 3][43, -13]\n", "\n", "n= 7:[-85, -9][-53, -9][-45, -33][-43, -21][-23, 15]\n", " [-21, -1][21, -27][23, 21][37, 21][43, -31]\n", " [53, 5][85, 9]\n", "\n", "n= 8:[-107, 53][-105, -27][-91, 5][-85, 49][-43, -19]\n", " [43, -1][75, 15][85, -13][105, 15][107, -1]\n", " [171, -21]\n", "\n", "n= 9:[-181, -33][-171, -29][85, -83][149, -3][151, 45]\n", " [171, -55][213, -19]\n", "\n", "n= 10:[-363, 109][-361, -99][-341, 153][299, -57]\n", "\n" ] } ], "source": [ "#Table 1. We calculate the crosscorrelation C_{n,s} of the nth Rudin-Shapiro sequence with its companion at\n", "#for selected values of n and s. First we do a brute-force verification of the table entries by directly\n", "#constructing Rudin-Shapiro pairs and computing their crosscorrelations directly\n", "\n", "def NegateList(L):\n", " return [-a for a in L]\n", "def RSPair(n):\n", " if(0==n):\n", " return [[1],[1]]\n", " small=RSPair(n-1)\n", " return [small[0]+small[1],small[0]+NegateList(small[1])]\n", "def CrossCorr(a,b,s):\n", " if(len(a)!=len(b)):\n", " print(\"error: mismatched sequence lengths\")\n", " return\n", " if(abs(s)>=len(a)):\n", " return 0\n", " if(s >= 0):\n", " return sum(a[j+s]*b[j] for j in range(len(a)-s))\n", " if(s < 0):\n", " return sum(a[j+s]*b[j] for j in range(-s,len(a)))\n", "def RSCrossCorr(n,s):\n", " pair=RSPair(n)\n", " return CrossCorr(pair[0],pair[1],s)\n", "\n", "table_one_shifts=[[0],\n", " [-1,1],\n", " [-3,-1,1,3],\n", " [-5,-3,-1,1,3,5],\n", " [-11,-7,-5,-3,3,5,7,11],\n", " [-21,-13,-11,-9,-5,5,9,11,13,21],\n", " [-43,-41,-27,-23,-21,-11,11,19,21,27,41,43],\n", " [-85,-53,-45,-43,-23,-21,21,23,37,43,53,85],\n", " [-107,-105,-91,-85,-43,43,75,85,105,107,171],\n", " [-181,-171,85,149,151,171,213],\n", " [-363,-361,-341,299]]\n", "\n", "table_one_data_brute=[[[s,RSCrossCorr(n,s)] for s in table_one_shifts[n]] for n in range(len(table_one_shifts))]\n", "\n", "def PresentTableOne(table):\n", " for n in range(len(table)):\n", " print(\"n=\",n,end=\":\")\n", " thingsinrow=0\n", " for pair in table[n]:\n", " print(pair,end=\"\")\n", " thingsinrow+=1\n", " if(5==thingsinrow):\n", " print(\"\\n\",\" \"*3,end=\" \")\n", " thingsinrow=0\n", " print(\"\\n\")\n", "\n", "PresentTableOne(table_one_data_brute)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n= 0:[0, 1]\n", "\n", "n= 1:[-1, -1][1, 1]\n", "\n", "n= 2:[-3, 1][-1, 1][1, 3][3, -1]\n", "\n", "n= 3:[-5, -1][-3, -5][-1, 3][1, 1][3, 1]\n", " [5, 1]\n", "\n", "n= 4:[-11, 5][-7, 1][-5, 1][-3, 5][3, 7]\n", " [5, 3][7, 3][11, -5]\n", "\n", "n= 5:[-21, -1][-13, -9][-11, -13][-9, 3][-5, 7]\n", " [5, -3][9, 9][11, -7][13, 5][21, 1]\n", " \n", "\n", "n= 6:[-43, 13][-41, -3][-27, 13][-23, -7][-21, 9]\n", " [-11, 5][11, 7][19, 15][21, -5][27, 7]\n", " [41, 3][43, -13]\n", "\n", "n= 7:[-85, -9][-53, -9][-45, -33][-43, -21][-23, 15]\n", " [-21, -1][21, -27][23, 21][37, 21][43, -31]\n", " [53, 5][85, 9]\n", "\n", "n= 8:[-107, 53][-105, -27][-91, 5][-85, 49][-43, -19]\n", " [43, -1][75, 15][85, -13][105, 15][107, -1]\n", " [171, -21]\n", "\n", "n= 9:[-181, -33][-171, -29][85, -83][149, -3][151, 45]\n", " [171, -55][213, -19]\n", "\n", "n= 10:[-363, 109][-361, -99][-341, 153][299, -57]\n", "\n", "\n", "It is True that our two versions of the table match\n" ] } ], "source": [ "#Table I: the claim is made that every table entry can be deduced from the following rules:\n", "#direct computation: C_{0,0}=1, C_{1,-1}=-1, C_{1,1}=1\n", "#known vanishing: C_{n,s}=0 if |s|>=2^n or if n>0 and s is even\n", "#Lemma 2.11 (where we drop the conjugates and use l_n=2^n for R-S sequences): \n", "#C_{n,s}= C_{n-1,2^{n-1}-s}+2C_{n-2}(2^{n-2}-s} for s > 0\n", "#C_{n,s}=-C_{n-1,2^{n-1}+s}+2C_{n-2}(2^{n-2}+s} for s < 0\n", "#Clearly the n=0 and n=1 rows of the table come from the direct computation\n", "#For n>1, we check that the s-values for which we must compute table entries are really deducible via Lemma 2.11\n", "#from previous table entries (along with known vanishing items)\n", "\n", "\n", "#This fetches a C_{n,s} value. It first tries to apply the known vanishing rules, and if that does not work,\n", "#Then it searches the furnished table_row (a list of pairs [s,C_{n,s}]). If it cannot find the required shift,\n", "#it quits.\n", "\n", "def FetchRSCrossCorr(n,s,table_row):\n", " if(n<0):\n", " print(\"n must be at least 0\")\n", " exit()\n", " if abs(s)>=2^n:\n", " return 0\n", " if (n>0) and (0==(s%2)):\n", " return 0\n", " for pair in table_row:\n", " if(pair[0]==s):\n", " return pair[1]\n", " print(\"cannot find crosscorrelation value\")\n", " exit()\n", " \n", "#DeducdeRSCrossCorr computes C_{n,s} using the above rules, where we input two previous \"rows\" of the table, \n", "#where the format #of a row is a list of pairs, each pair being [s,C_{k,s}] (k=n-1 for row_one_below, \n", "#and k=n-2 for row_two_below). This only works for n>=2.\n", "\n", "def DeduceRSCrossCorr(n,s,row_one_below,row_two_below):\n", " if(n<2):\n", " print(\"n must be at least two\")\n", " exit()\n", " if abs(s)>=(2^n):\n", " return 0\n", " if 0==(s%2):\n", " return 0\n", " C_one_below=FetchRSCrossCorr(n-1,2^(n-1)-abs(s),row_one_below)\n", " C_two_below=FetchRSCrossCorr(n-2,2^(n-2)-abs(s),row_two_below)\n", " if(s>0):\n", " return C_one_below+2*C_two_below\n", " return -C_one_below+2*C_two_below\n", "\n", "#This calculates the entries that should be on Table I by applying the deduction rules as promised:\n", "\n", "def ProduceTableOne(table_of_shifts):\n", " #these are the direct computations\n", " table=[[[0,1]],[[-1,-1],[1,1]]]\n", " for n in range(2,len(table_of_shifts)):\n", " row=[[s,DeduceRSCrossCorr(n,s,table[n-1],table[n-2])] for s in table_of_shifts[n]]\n", " table+=[row]\n", " return table\n", "\n", "table_one_data_deduced=ProduceTableOne(table_one_shifts)\n", "\n", "PresentTableOne(table_one_data_deduced)\n", "\n", "print(\"\\nIt is\", table_one_data_brute==table_one_data_deduced, \"that our two versions of the table match\")" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "#This is the reduction algorithm from Section 2.4\n", "R.=PolynomialRing(QQ)\n", "\n", "#These functions effect the reduction of a(alpha0,alpha1,alpha2) to the standard reduced form b(alpha0,alpha1) \n", "#described in the second paragraph.\n", "\n", "def Eliminatealpha2(a):\n", " return a.subs(alpha2=-alpha0-alpha1-1)\n", "def Loweralpha1(a):\n", " r=a%(alpha1^2)\n", " q=a//(alpha1^2)\n", " if R(0)==q:\n", " return a\n", " return Loweralpha1(r+q*(-alpha0*alpha1-alpha1-alpha0^2-alpha0+2))\n", "def Loweralpha0(a):\n", " return a%(alpha0^3+alpha0^2-2*alpha0-4) \n", "def ReducePoly(a):\n", " return Loweralpha0(Loweralpha1(Eliminatealpha2(a)))\n", "\n", "#EuclideanOneStep takes a pair [[a,b,c],[d,e,f]]\n", "#where it is supposed that there are some initial inputs S,T to a Euclidean algorithm\n", "#and a and d are consecutive remainders in a Euclidean algorithm with a=bS+cT and d=eS+fT\n", "#then it calculates quotient and remainder q,r such that a=qd+r and then \n", "#we have r=a-qd=(bS+cT)-q(eS+fT)=(b-qe)S+(c-qf)T.\n", "#so EuclideanOneStep outputs [[d,e,f],[r,b-qe,c-qf]]\n", "\n", "def EuclideanOneStep(pair_of_triples):\n", " dividend=pair_of_triples[0][0]\n", " divisor=pair_of_triples[1][0]\n", " quotient=dividend//divisor\n", " remainder=dividend%divisor\n", " new_triple=[remainder,pair_of_triples[0][1]-quotient*pair_of_triples[1][1],pair_of_triples[0][2]-quotient*pair_of_triples[1][2]]\n", " return [pair_of_triples[1],new_triple]\n", "\n", "#Given a and b, this returns the triple with [gcd(a,b),u,v] where gcd(a,b)=u*a+v*b\n", "\n", "def EuclideanAlgorithm(a,b):\n", " for c in [a,b]:\n", " if (R(0)!=(c//alpha1)) or (R(0)!=(c//alpha2)):\n", " print(\"We only run Euclidean Algorithm on polynomials in alpha0\")\n", " exit()\n", " if R(0)==a:\n", " if R(0)==b:\n", " return[R(0),R(0),R(0)]\n", " lead_coeff=b.coefficient(alpha0^b.degree())\n", " return [(1/lead_coeff)*b,R(0),R(1/lead_coeff)]\n", " if R(0)==b:\n", " lead_coeff=a.coefficient(alpha0^a.degree())\n", " return [(1/lead_coeff)*a,R(1/lead_coeff),R(0)]\n", " pair_of_triples=[[a,R(1),R(0)],[b,R(0),R(1)]]\n", " while R(0)!=pair_of_triples[1][0]:\n", " pair_of_triples=EuclideanOneStep(pair_of_triples)\n", " triple=pair_of_triples[0]\n", " gcd=triple[0]\n", " lead_coeff=gcd.coefficient(alpha0^gcd.degree())\n", " triple=[(1/lead_coeff)*u for u in triple]\n", " return triple\n", "\n", "#If e is a polynomial representing a real number, then this gives a polynomial in alpha0 that represents\n", "#the multiplicative inverse of e via the polynomial Euclidean algorithm\n", "\n", "def InvertRealPoly(a):\n", " a=ReducePoly(a)\n", " if R(0)!=(a//alpha1):\n", " print(\"Trying to invert non-real\")\n", " exit()\n", " if R(0)==a:\n", " print(\"Trying to invert zero\")\n", " exit()\n", " return EuclideanAlgorithm(a,alpha0^3+alpha0^2-2*alpha0-4)[1] \n", "\n", "#This reduces the expression c(alpha0,alpha1,alpha2)/d(alpha0,alpha1,alpha2) to a \n", "#standard reduced form as in the third #paragraph of Section 2.5\n", "\n", "def ReduceRF(c,d):\n", " f=InvertRealPoly(d*d.subs(alpha1=alpha2,alpha2=alpha1))\n", " return ReducePoly(c*d.subs(alpha1=alpha2,alpha2=alpha1)*f)\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "#Definition 2.17: This is the signifier for a real quantity represented by a polynomial \n", "\n", "def sig(a):\n", " a=ReducePoly(a)\n", " if R(0)!=(a//alpha1):\n", " print(\"Trying to obtain a signifier of a non-real quantity\")\n", " exit()\n", " p=a.constant_coefficient()\n", " q=a.coefficient(alpha0)\n", " r=a.coefficient(alpha0^2)\n", " return p^3-p^2*q-2*p*q^2+4*q^3+5*p^2*r-10*p*q*r-4*q^2*r+12*p*r^2-8*q*r^2+16*r^3" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "remainder: 0\n", "quotient: r^3*y^3 + 3*q*r^2*y^2 - r^3*y^2 + 3*q^2*r*y - 2*q*r^2*y - 2*r^3*y + q^3 - q^2*r - 2*q*r^2 + 4*r^3\n", "\n", "it is True that the quotient came out as expected\n" ] } ], "source": [ "#Lemma 2.18: check the claim about factorization of n(w(y))\n", "Special_Ring.=PolynomialRing(QQ)\n", "s = -3*p + q - 5*r\n", "t = 3*p^2 - 2*p*q - 2*q^2 + 10*p*r - 10*q*r + 12*r^2\n", "u = -p^3 + p^2*q + 2*p*q^2 - 4*q^3 - 5*p^2*r + 10*p*q*r + 4*q^2*r - 12*p*r^2 + 8*q*r^2 - 16*r^3\n", "n=y^3+s*y^2+t*y+u\n", "w=p+q*y+r*y^2\n", "expected_quotient=y^3*r^3 + 3*y^2*q*r^2 -y^2*r^3 + 3*y*q^2*r - 2*y*q*r^2 -2*y*r^3 + q^3 - q^2*r - 2*q*r^2 + 4*r^3\n", "poly=n.subs(y=w)\n", "remainder=poly%(y^3+y^2-2*y-4)\n", "quotient=poly//(y^3+y^2-2*y-4)\n", "print(\"remainder:\", remainder)\n", "print(\"quotient:\", quotient)\n", "print(\"\\nit is\", quotient==expected_quotient, \"that the quotient came out as expected\")" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "784310082937/1000000000000000000\n", "-17168250811/1953125000000000\n" ] } ], "source": [ "#Paragraph of Section 2.4 before Lemma 2.20: this is another check that alpha_0=1.658967...\n", "print(sig(alpha0-1658967/10^6))\n", "\n", "print(sig(alpha0-1658968/10^6))" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1/2*alpha0^2 - 1/2\n", "31546657703741739416522343271/15625000000000000000000000000000000\n", "-86486793323783612696476801/64000000000000000000000000000000\n" ] } ], "source": [ "#Lemma 2.20: proof that |alpha1/alpha0|=sqrt((alpha0^2-1)/2), which is\n", "#0.935994...., using the fact that alpha1 and alpha2 are conjugates\n", "rho1=ReduceRF(alpha1,alpha0)\n", "rho2=ReduceRF(alpha2,alpha0)\n", "rho1rho2=ReducePoly(rho1*rho2)\n", "print(rho1rho2)\n", "print(sig(rho1rho2-(935994/10^6)^2))\n", "print(sig(rho1rho2-(935995/10^6)^2))" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "It is True that each new row of the table can be obtained from the previous one\n", "\n", "t= 1:[-1, -1, 0, 2, 0] [0, 0, 1, 0, 2] \n", " \n", "t= 2:[-1, 2, -1, 0, -2] [0, 1, 2, 2, 0] \n", " \n", "t= 3:[-2, -3, -2, 2, 0] [-1, 0, 3, 0, 2] \n", " [0, 1, 0, 6, 0] [1, 2, -1, 0, 6] \n", " \n", "t= 4:[-4, 1, 0, -10, 0] [-3, 2, -1, 0, -10] \n", " [1, 6, 1, 0, 2] [2, -3, 6, 2, 0] \n", " \n", "t= 5:[-6, -3, -10, 2, 0] [2, -5, 2, 14, 0] \n", " [4, 9, 0, 6, 0] [5, 2, -9, 0, 6] \n", " \n", "t= 6:[-12, -7, 0, -26, 0] [-11, 2, 7, 0, -26] \n", " [5, 14, -7, 0, -6] [9, 6, 9, 0, 18] \n", " [10, -11, 6, -14, 0] \n", "t= 7:[-23, -26, -7, 0, -14] [-22, 5, -26, 18, 0] \n", " [10, -21, -6, 14, 0] [11, 0, 21, 0, 14] \n", " [18, 3, 18, 30, 0] [21, -14, -17, 0, -10] \n", " \n", "t= 8:[-46, 19, -14, -66, 0] [-43, 18, 31, 0, -42] \n", " [21, 14, -15, 0, -54] [37, 30, -15, 0, 42] \n", " [42, -3, -10, -62, 0] \n", "t= 9:[-91, -66, 33, 0, 10] [-86, 13, -42, 98, 0] \n", " [42, -29, -54, -2, 0] [74, -45, 42, 30, 0] \n", " \n", "t= 10:[-182, 99, 10, -66, 0] [-181, 0, -99, 0, -66] \n", " [-171, 98, 55, 0, -58] [149, 30, -87, 0, -6] \n", " " ] } ], "source": [ "#Table 2\n", "#Code that uses Remark 3.3 to obtain the values recursively\n", "#Note that if t>1, then to get any one of these at (t,j), it suffices to know the ones at (t-1,floor(j/2))\n", "#So the list of js for the (t-1)th row of the table should always contain the floored halves of the js in the\n", "#tth row\n", "\n", "def BigAlpha(t,j):\n", " if(1==t):\n", " if((-1)==j):\n", " return -1\n", " return 0\n", " if(0==(j%2)):\n", " return -BigAlpha(t-1,j/2)+BigBeta(t-1,j/2)\n", " return BigGamma(t-1,(j-1)/2)\n", "\n", "def BigBeta(t,j):\n", " if(1==t):\n", " if(0==j):\n", " return 1\n", " return 0\n", " if(0==(j%2)):\n", " return BigDelta(t-1,j/2)\n", " return BigAlpha(t-1,(j-1)/2)-BigBeta(t-1,(j-1)/2)\n", "def BigGamma(t,j):\n", " if(1==t):\n", " if((-1)==j):\n", " return 2\n", " return 0\n", " if(0==(j%2)):\n", " return 2*BigAlpha(t-1,j/2)+2*BigBeta(t-1,j/2)\n", " return 0\n", "def BigDelta(t,j):\n", " if(1==t):\n", " if(0==j):\n", " return 2\n", " return 0\n", " if(0==(j%2)):\n", " return 0\n", " return 2*BigAlpha(t-1,(j-1)/2)+2*BigBeta(t-1,(j-1)/2)\n", "\n", "#js[t] is the list of j values such that we want BigAlpha_{t,j}, BigBeta_{t,j}, BigGamma_{t,j}, and BigDelta_{t,j}\n", "#in our table\n", "\n", "js=[[],\n", " [-1, 0],\n", " [-1, 0],\n", " [-2, -1, 0, 1],\n", " [-4, -3, 1, 2],\n", " [-6, 2, 4, 5],\n", " [-12, -11, 5, 9, 10],\n", " [-23, -22, 10, 11, 18, 21],\n", " [-46, -43, 21, 37, 42],\n", " [-91, -86, 42, 74],\n", " [-182, -181, -171, 149]]\n", "\n", "#first we check that for t>1, one can obtain the tth row of the table from the (t-1)th row, by checking that\n", "#the floored halves of the items in the tth row do occur in the (t-1)th row.\n", "\n", "table_supports_self=True\n", "for t in range(2,len(js)):\n", " for j in js[t]:\n", " if not(floor(j/2) in js[t-1]):\n", " print(\"in row t=\", t, \"the entry j=\", j, \"requires the presence of floor(j/2)=\", floor(j/2), \"in the previous row, but it is missing\")\n", " table_supports_self=False\n", "print(\"It is\", table_supports_self, \"that each new row of the table can be obtained from the previous one\")\n", "\n", "#now we print the table\n", "for t in range(1,len(js)):\n", " print(\"\\nt=\",t,end=\":\")\n", " thingsinrow=0\n", " for j in js[t]:\n", " print([j,BigAlpha(t,j),BigBeta(t,j),BigGamma(t,j),BigDelta(t,j)],end=\" \")\n", " thingsinrow+=1\n", " if(2==thingsinrow):\n", " print(\"\\n\"+\" \"*4,end=\" \")\n", " thingsinrow=0" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n: 0 s: 0 corr here: 1 corr CPP: 1 difference: 0\n", "n: 1 s: -1 corr here: -1 corr CPP: -1 difference: 0\n", "n: 1 s: 1 corr here: 1 corr CPP: 1 difference: 0\n", "n: 2 s: 1 corr here: 3 corr CPP: 3 difference: 0\n", "n: 3 s: -3 corr here: -5 corr CPP: -5 difference: 0\n", "n: 4 s: 3 corr here: 7 corr CPP: 7 difference: 0\n", "n: 5 s: -11 corr here: -13 corr CPP: -13 difference: 0\n", "n: 6 s: 13 corr here: 19 corr CPP: 19 difference: 0\n", "n: 7 s: -45 corr here: -33 corr CPP: -33 difference: 0\n", "n: 8 s: -107 corr here: 53 corr CPP: 53 difference: 0\n", "n: 9 s: -179 corr here: -85 corr CPP: -85 difference: 0\n", "n: 10 s: -341 corr here: 153 corr CPP: 153 difference: 0\n", "n: 11 s: -717 corr here: -217 corr CPP: -217 difference: 0\n", "n: 12 s: -1451 corr here: 373 corr CPP: 373 difference: 0\n", "n: 13 s: -2867 corr here: -557 corr CPP: -557 difference: 0\n", "n: 14 s: -5453 corr here: 961 corr CPP: 961 difference: 0\n", "n: 15 s: -10923 corr here: -1717 corr CPP: -1717 difference: 0\n", "n: 16 s: -22955 corr here: 2445 corr CPP: 2445 difference: 0\n", "n: 17 s: -43691 corr here: -4285 corr CPP: -4285 difference: 0\n", "n: 18 s: -91733 corr here: 6257 corr CPP: 6257 difference: 0\n", "n: 19 s: -174765 corr here: -11153 corr CPP: -11153 difference: 0\n", "n: 20 s: -349525 corr here: 19041 corr CPP: 19041 difference: 0\n", "n: 21 s: -699059 corr here: -28293 corr CPP: -28293 difference: 0\n", "n: 22 s: -1398101 corr here: 53321 corr CPP: 53321 difference: 0\n", "n: 23 s: -2796237 corr here: -72905 corr CPP: -72905 difference: 0\n", "n: 24 s: -5592403 corr here: 129485 corr CPP: 129485 difference: 0\n", "n: 25 s: -11184811 corr here: -214365 corr CPP: -214365 difference: 0\n" ] } ], "source": [ "#data from the C++ program that calculates PCCs and PSLs for Golay-Rudin-Shapiro sequences\n", "PCCSfromCPP=[1,1,3,5,7,13,19,33,53,85,153,217,373,557,961,1717,2445,4285,6257,11153,19041,28293,53321,72905,129485,214365,342769,640933,860709,1624877,2490985,4188609,7618449,10688117,20617465,29999429,51521697,90947021,133991557,255886741,372089521,668060317,1099665689,1724813029,3146759617,4701529197,8491242153,13498854709,22289746385,38672931645,59901979961]\n", "PCClocsfromCPP=[[0],[-1,1],[1],[-3],[3],[-11],[13],[-45],[-107],[-179],[-341],[-717],[-1451],[-2867],[-5453],[-10923],[-22955],[-43691],[-91733],[-174765],[-349525],[-699059],[-1398101],[-2796237],[-5592403],[-11184811],[-22369613],[-44739243],[-89478451],[-178956971],[-357913941],[-715827885],[-1431655765],[-2863311539],[-5726623061],[-11453246123],[-22906492245],[-45812984491],[-91625968979],[-183251937963],[-366503875925],[-733007751851],[-1466015503701],[-2932031007403],[-5864062014805],[-11728124029611],[-23456248059221],[-46912496118443],[-93824992236885],[-187649984473771],[-375299968947541]]\n", "PCCvalsfromCPP=[[1],[-1,1],[3],[-5],[7],[-13],[19],[-33],[53],[-85],[153],[-217],[373],[-557],[961],[-1717],[2445],[-4285],[6257],[-11153],[19041],[-28293],[53321],[-72905],[129485],[-214365],[342769],[-640933],[860709],[-1624877],[2490985],[-4188609],[7618449],[-10688117],[20617465],[-29999429],[51521697],[-90947021],[133991557],[-255886741],[372089521],[-668060317],[1099665689],[-1724813029],[3146759617],[-4701529197],[8491242153],[-13498854709],[22289746385],[-38672931645],[59901979961]]\n", "\n", "#Some basic checking that the C++ program agress with calculations here for the shorter sequences\n", "for n in range(26):\n", " for j in range(len(PCClocsfromCPP[n])):\n", " s=PCClocsfromCPP[n][j]\n", " corrcpp=PCCvalsfromCPP[n][j]\n", " corrhere=RSCrossCorr(n,s)\n", " diff=corrhere-corrcpp\n", " print(\"n:\",n,\"s:\",s,\"corr here:\",corrhere,\"corr CPP:\",corrcpp,\"difference:\",diff)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "t=1: -16\n", "t=2: nothing to do\n", "t=3: -32\n", "t=4: -8 \n", " -296\n", "t=5: -136\n", "t=6: -14752\n", "t=7: -2704 \n", " -1864\n", "t=8: -158656 \n", " -78424\n", "t=9: -389464 \n", " -477632\n", "t=10: -266744 \n", " -267712 \n", " -1139224\n" ] } ], "source": [ "#Propositon 3.7: this checks the inequalities in the final list\n", "print(\"t=1:\", sig(alpha0-alpha0^2))\n", "print(\"t=2: nothing to do\")\n", "print(\"t=3:\", sig(3*alpha0+2-alpha0^4))\n", "print(\"t=4:\", sig(alpha0+10-alpha0^5), \"\\n \", sig(7*alpha0-alpha0^5))\n", "print(\"t=5:\", sig(11*alpha0-alpha0^6))\n", "print(\"t=6:\", sig(17*alpha0-alpha0^7))\n", "print(\"t=7:\", sig(33*alpha0-alpha0^8), \"\\n \", sig(21*alpha0+14-alpha0^8))\n", "print(\"t=8:\", sig(49*alpha0-alpha0^9), \"\\n \", sig(13*alpha0+62-alpha0^9))\n", "print(\"t=9:\", sig(87*alpha0-alpha0^10), \"\\n \", sig(83*alpha0+2-alpha0^10))\n", "print(\"t=10:\", sig(109*alpha0+66-alpha0^11), \"\\n \", sig(153*alpha0-alpha0^11), \"\\n \", sig(117*alpha0+6-alpha0^11))" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "4/59\n", "-16/59\n", "4/59\n", "-22/59\n" ] } ], "source": [ "#Proposition 4.11: claim that 0