orangesquid (os) wrote,
orangesquid
os

pseudo-base-prime part 6 (some updates to code from part 5)

Here we go: bc code:
/* example usage: */
/* j=num2p(200,f[]) --> 3 0 2 in f[] */
/* n=p2num(f[],j) returns 200 */

define num2p(n, *f[]) {
  auto i, j, l, s[], z;
  z = scale;
  scale = 0;
  l = n;
  for(i = 2; i <= l; ++i) {
    s[i - 2] = 1;
  }
  for(i = 2; i <= l; ++i) {
    for(j = 2; i*j <= l; ++j) {
      s[i*j - 2] = 0;
    }
  }
  for(i = 2; i <= l; ++i) {
    print i, ": ", s[i - 2], "\n";
  }
  j = 0;
  for(i = 2; i <= l; ++i) {
    if(s[i - 2]) {
      f[j] = 0;
      print n, " % ", i, " == ", n % i;
      while(n % i == 0) {
        ++f[j];
        n = n / i;
        print "for ", i, " (p#", j, "): n -> ", n, "; f -> ", f[j], "\n";
      }
      ++j;
    }
  }
  for(i = 0; i < j; ++i) {
    print f[i], " "
  }
  print "\n"
  scale = z;
  return j;
}

define void basep(n, *f[]) {
  auto j, g[], i;
  j = num2p(n,g[]);
  for(i = 0; i < j; ++i) {
    f[i+1] = g[i];
  }
  f[0] = j;
  bpnorm(f[]);
}

define bp2num(f[]) {
  auto i, p[], r, n;
  r=f[0];
  for(i = 0; i < r; ++i) {
    p[i] = f[i + 1];
  }
  n = p2num(p[],r);
  return n;
}

define p2num(p[], r) {
  auto i, j, l, s[], z, o, k;
  z = scale;
  scale = 0;
  l = r;
  k = 0;
  while(k < r) {
    l = l * 2;
    for(i = 2; i <= l; ++i) {
      s[i - 2] = 1;
    }
    for(i = 2; i <= l; ++i) {
      for(j = 2; i*j <= l; ++j) {
        s[i*j - 2] = 0;
      }
    }
    k = 0;
    for(i = 2; i <= l; ++i) {
      if(s[i - 2]) {
        ++k;
      }
    }
  }
  print "l=", l, " k=", k, "\n";
  for(i = 2; i <= l; ++i) {
    print i, ": ", s[i - 2], "\n";
  }
  o = 1;
  j = 0;
  for(i = 2; i <= l; ++i) {
    if(s[i - 2]) {
      if(p[j]) {
        o = o * i^p[j];
        print "after ", i, " (p#", j, "): o -> ", o, " from p:", p[j], "\n";
      }
      ++j;
    }
  }
  scale = z;
  return o;
}

define squidn(a,b) {
  auto m, n, p[], q[], r[], j, k, l;
  j = num2p(a,p[]);
  k = num2p(b,q[]);
  l = j + k;
  print "j:", j, " k:", k, " -> l:", l, "\n";
  for(m = 0; m < l; ++m) {
    r[m] = 0;
  }
  for(m = 0; m < j; ++m) {
    for(n = 0; n < k; ++n) {
      r[m + n] = r[m + n] + p[m]*q[n];
      print "@", m, ",", n, " r->", r[m+n], " from ", p[m], ",", q[n], "\n";
    }
  }
  return p2num(r[],l);
}

define void squid(a[], b[], *c[]) {
  auto m, n, p[], q[], r[], j, k, l;
  j = a[0];
  k = b[0];
  for(n = 0; n < j; ++n) {
    p[n] = a[n + 1];
  }
  for(n = 0; n < k; ++n) {
    q[n] = b[n + 1];
  }
  l = j + k;
  print "j:", j, " k:", k, " -> l:", l, "\n";
  for(m = 0; m < l; ++m) {
    r[m] = 0;
  }
  for(m = 0; m < j; ++m) {
    for(n = 0; n < k; ++n) {
      r[m + n] = r[m + n] + p[m]*q[n];
      print "@", m, ",", n, " r->", r[m+n], " from ", p[m], ",", q[n], "\n";
    }
  }
  for(n = 0; n < l; ++n) {
    c[n + 1] = r[n];
  }
  c[0] = l;
}

define void bpconj(p[],*q[]) {
  auto i;
  if(p[0] == 0) {
    q[0] = 1;
    q[1] = 1;
    return;
  }
  if(p[1] == 0) {
    for(i = 2; i <= p[0]; ++i) {
      q[i - 1] = p[i];
    }
    q[0] = p[0] - 1;
    q[1] = q[1] + 1;
  } else {
    for(i = 1; i <= p[0]; ++i) {
      q[i + 1] = p[i];
    }
    q[1] = 0;
    q[0] = p[0] + 1;
    q[2] = q[2] - 1;
  }
}

define void bpnorm(*p[]) {
  auto i, c, b;
  b=0;
  for(i = p[0]; i > 0; --i) {
    if(b == 0) {
      if(p[i] == 0) {
        c++;
      } else {
        b++;
      }
    }
  }
  p[0] = p[0] - c;
}

define void dispnonl(p[]) {
  auto i;
  print "{";
  for(i = 1; i <= p[0]; ++i) {
    print p[i];
    if(i < p[0]) {
      print ",";
    }
  }
  print "}";
}

define void disp(p[]) {
  dispnonl(p[]);
  print "\n";
}

define void bptree(r,c,*p[]) {
  auto i, j, s, m, o[], z;
  s = scale;
  scale = 0;
  p[0] = 1;
  p[1] = 0;
  for(i = 0; i < r; ++i) {
    print "@r", i, ":";
    m = c / 2;
    m = m * 2;
    m = c - m;
    if(m > 0) {
      print "R";
      o[r - i - 1] = 0;
      c = (c - 1) / 2;
    } else {
      print "L";
      o[r - i - 1] = 1;
      c = c / 2;
    }
    print ",c:=", c, "\n";
    disp(p[]);
  }
  for(i = 0; i < r; ++i) {
    print "op", i, ":", o[i], "\n";
    if(o[i] > 0) {
      p[1] = p[1] + 1;
    } else {
      print "n=", p[0], ":";
      disp(p[]);
      for(j = p[0]; j >= 1; --j) {
        print "j.", j, ".", p[j], "\n";
        p[j + 1] = p[j];
      }
      p[1] = 0;
      p[0] = p[0] + 1;
    }
    disp(p[]);
  }
  scale = s;
}

define void inbptree(p[],*rc[]) {
  auto g, r, c, i, o[], n;
  g = 1;
  n = 0;
  while(g > 0) {
    disp(p[]);
    if(p[0] == 1) {
      if(p[1] == 0) {
        p[0] = 0;
      }
    }
    if(p[0] > 0) {
      if(p[1] > 0) {
        p[1] = p[1] - 1;
        o[n++] = 1;
        print "L\n";
      } else {
        for(i = 1; i < p[0]; ++i) {
          p[i] = p[i + 1];
        }
        p[0] = p[0] - 1;
        o[n++] = 0;
        print "R\n";
      }
    } else {
      g = 0;
    }
  }
  r = 0;
  c = 0;
  for(i = n - 1; i >= 0; --i) {
    if(o[i] > 0) {
      r++;
      c = c * 2;
      print "op1\n";
    } else {
      r++;
      c = c * 2 + 1;
      print "op0\n";
    }
  }
  rc[0] = r;
  rc[1] = c;
  print "(",r,",",c,")\n";
}

define void bpseq(n,*p[]) {
  auto s, r, c, w;
  s = scale;
  scale = 0;
  w = 1;
  r = 0;
  c = 0;
  if(n == 0) {
    p[0] = 0;
    scale = s;
    return;
  }
  while(n >= w) {
    w = w * 2;
    ++r;
  }
  print "2^",w,">n@r=",r,"\n";
  c = n - w / 2;
  bptree(r,c,p[]);
  scale = s;
}

define void nbpord(n,*s[]) {
  auto i, p[];
  for(i = 0; i < n; ++i) {
    bpseq(i,p[]);
    print "seq#",i,":";
    disp(p[]);
    s[i + 1] = bp2num(p[]);
  }
  s[0] = n;
}

define void bpordnth(n,*s[]) {
  auto i, m, j, p[], x;
  k = 1;
  s[0] = n;
  print "bpordnth: 0 to ",n,"\n";
  for(i = 0; i < n; ++i) {
    print "bpordnth: bpseq(",i,")-1=";
    bpseq(i,p[]);
    s[i + 1] = bp2num(p[])-1;
    print s[i+1],"\n";
  }
}

define bporddif(n) {
  auto i, o[], s;
  nbpord(2^n, o[]);
  s = 0;
  if(n > 0) {
    for(i = 1; i <= o[0]; ++i) {
      s = s + o[i];
    }
    s = s - 2^(2*n-1) - 2^(n-1);
  }
  return s;
}

define void bpordlst(n,*s[]) {
  auto i;
  if(n < 3) {
    n = 3;
  }
  s[0] = n + 1;
  for(i = 0; i <= n; ++i) {
    s[i + 1] = bporddif(i);
  }
}

define void bpordmul(n,*m[]) {
  auto i, s[], d;
  if(n < 4) {
    m[0] = 0;
    return;
  }
  d = scale;
  scale = 10;
  m[0] = n - 3;
  bpordlst(n, s[]);
  for(i = 4; i <= n; ++i) {
    m[i - 3] = s[i + 1] / s[i];
    print "bpordmul@", i - 3, ": ", m[i - 3], "\n";
  }
  scale = d;
}

bpordm10[0] = 7;
bpordm10[1] = 16;
bpordm10[2] = 10.6875;
bpordm10[3] = 9.28070175439;
bpordm10[4] = 8.98361688721;
bpordm10[5] = 9.25657571719;
bpordm10[6] = 9.89187776102;
bpordm10[7] = 10.7296193394;

/* I don't remember what this is... */
bpwtf[0] = 27;
bpwtf[1] = 1;
bpwtf[2] = 2;
bpwtf[3] = 2;
bpwtf[4] = 4;
bpwtf[5] = 2;
bpwtf[6] = 6;
bpwtf[7] = 2;
bpwtf[8] = 8;
bpwtf[9] = 4;
bpwtf[10] = 10;
bpwtf[11] = 2;
bpwtf[12] = 18;
bpwtf[13] = 2;
bpwtf[14] = 14;
bpwtf[15] = 6;
bpwtf[16] = 16;
bpwtf[17] = 2;
bpwtf[18] = 12;
bpwtf[19] = 2;
bpwtf[20] = 50;
bpwtf[21] = 10;
bpwtf[22] = 22;
bpwtf[23] = 2;
bpwtf[24] = 54;
bpwtf[25] = 4;
bpwtf[26] = 26;
bpwtf[27] = 8;

define void lstadd0s(*a[], *b[]) {
  auto i;
  if(a[0] >= b[0]) {
    for(i = b[0] + 1; i <= a[0]; ++i) {
      b[i] = 0;
    }
    b[0] = a[0];
  } else {
    for(i = a[0] + 1; i <= b[0]; ++i) {
      a[i] = 0;
    }
    a[0] = b[0];
  }
}

define void bpmax(a[], b[], *o[]) {
  auto i;
  lstadd0s(a[], b[]);
  o[0] = a[0];
  for(i = 1; i <= a[0]; ++i) {
    if(a[i] >= b[i]) {
      o[i] = a[i];
    } else {
      o[i] = b[i];
    }
  }
}

define void bpmin(a[], b[], *o[]) {
  auto i;
  lstadd0s(a[], b[]);
  o[0] = a[0];
  for(i = 1; i <= a[0]; ++i) {
    if(a[i] >= b[i]) {
      o[i] = b[i];
    } else {
      o[i] = a[i];
    }
  }
  bpnorm(o[]);
}

define void bpfunc(n, *t[]) {
  /* before calling, define a function bpfuncf: */
  /* define bpfuncf(a[], b[], *c[]) {           */
  /*   ...                                      */
  /* }                                          */
  /* for example:                               */
  /* define bpfuncf(a[], b[], *c[]) {           */
  /*   squid(a[], b[], *c[])                    */
  /* }                                          */
  auto r, c, x[], y[], z[];
  t[0] = -n;
  t[1] = -n;
  for(r = 0; r < n; ++r) {
    for(c = 0; c < n; ++c) {
      print "bpfunc: f(",r+1,",",c+1,")\n"
      basep(r + 1, x[]);
      basep(c + 1, y[]);
      print "bpfuncf(";
      dispnonl(x[]);
      print ",";
      dispnonl(y[]);
      print ") = \n";
      bpfuncf(x[], y[], z[]);
      print "bpfuncf result ";
      dispnonl(z[]);
      print "\n";
      t[r * n + c + 2] = bp2num(z[]);
      print "bpfunc: result ",t[r * n + c + 2], " into [",r*n+c+2,"]\n";
    }
  }
}

define void bpfuncf(a[], b[], *c[]) {
  squid(a[], b[], c[]);
}

define void disptbl(t[]) {
  auto i, j;
  for(j = 0; j < -t[0]; ++j) {
    print "[";
    for(i = 0; i < -t[1]; ++i) {
      print t[-t[1] * j + i + 2];
      if(i < -t[1]-1) {
        print ",";
      }
    }
    print "]\n";
  }
}

define void bpmirror(i[], *o[]) {
  auto n, x, s;
  bpnorm(i[]);
  n = 0;
  for(x = 1; i[x] == 0; ++x) {
    ++n;
  }
  print "chop ",n," zeroes\n";
  for(x = 1; x <= n + 1; ++x) {
    i[x] = i[x + n];
  }
  i[0] = i[0] - n;
  s = scale;
  scale = 0;
  n = i[0] * 2 - i[0] % 2 - n;
  print "expand to ",n," elems\n";
  scale = s;
  for(x = i[0] + 1; x <= n; ++x) {
    i[x] = 0;
  }
  o[0] = n;
  for(x = 1; x <= n; ++x) {
    print "moving ",i[x],"\n";
    o[n - x + 1] = i[x];
  }
  bpnorm(o[]);
}

define bpfact(n) {
  auto z, r;
  r = 2;
  for(z = 2; z <= n; ++z) {
    print "call squidn on ",r,z,"\n";
    r = squidn(r, z);
  }
  return r;
}

define bptfunc(a, b) {
  /* just like bpfunc: define bptfuncf but with scalars: */
  /* define bptfunc(a, b) {                              */
  /*   ...                                               */
  /*   return c;                                         */
  /* }                                                   */
  auto p[], q[], x, y, z, r[], w[];
  bpseq(a,p[]);
  bpseq(b,q[]);
  x = bp2num(p[]);
  y = bp2num(q[]);
  z = bptfuncf(x, y);
  basep(z, r[]);
  inbptree(r[], w[]);
  if(w[0] == 0) {
    return 0;
  }
  return 2^(w[0]-1)+w[1];
}

define bptfuncf(a, b) {
  return a + b;
}

define bp2rle(p[], *r[]) {
  auto i, v, c;
  v = p[1];
  r[0] = 0;
  c = 1;
  for(i = 2; i <= p[0]; ++i) {
    if((p[i] == v)) {
      ++c;
    } else {
      print v,": ",c," @ ",r[0],"\n";
      r[r[0]+1] = v;
      r[r[0]+2] = c;
      r[0] = r[0] + 2;
      c = 1;
      v = p[i];
    }
  }
  print v,": ",c," @ ",r[0],"\n";
  r[r[0]+1] = v;
  r[r[0]+2] = c;
  r[0] = r[0] + 2;
}

define void rle2bp(a[], *b[]) {
  auto v, c, j, m, p, s;
  p = 1;
  s = scale;
  scale = 0;
  print "rle groups: ", a[0]/2, "\n";
  scale = s;
  for(m = 1; m < a[0]; m = m + 2) {
    v = a[m];
    c = a[m + 1];
    print m, "@", p, ": ", c, " ", v, "'s\n"
    for(j = 0; j < c; ++j) {
      b[p++] = v;
    }
  }
  b[0] = p - 1;
}

define bp2binstr(p[]) {
  auto i, j, c, o[], v, q;
  c = 1;
  for(i = 1; i <= p[0]; ++i) {
    v = p[i];
    if(v > 0) {
      for(j = 0; j < v; ++j) {
        print "1@",j+c,"\n";
        o[j + c] = 1;
      }
      c = c + j;
      o[c++] = 0;
    } else {
      o[c++] = 0;
    }
  }
  o[0] = c;
  disp(o[]);
  q = 0;
  for(i = 0; i < c; ++i) {
   q = o[i + 1] * 2 ^ i + q;
  }
  return q;
}

define void binstr2bp(b, *p[]) {
  auto s, c, i;
  s = scale;
  scale = 0;
  for(i = 1; b; b = b / 2) {
    print "i=",i," bit:",b%2," c=",c,"\n";
    if(b%2 == 0) {
      print "hit 0\n"
      p[i++] = c;
      c = 0;
    } else {
      ++c;
    }
  }
  p[i] = c;
  print i,"\n";
  p[0] = i;
  scale = s;
}

---snip---
I had to re-define two functions on the HP49 and define a new one to match
the above (note that squid was renamed to squidn in bc).
BPMIN:
<< { 0 } + SWAP { 0 } +
  LSTADD0S 2 << MIN >> DOLIST
  BPNORM >>

BPMAX:
<< { 0 } + SWAP { 0 } +
  LSTADD0S 2 << MAX >> DOLIST >>

This corrects the errors when passing in {} as an argument, which allows
BPMIN and BPMAX to be passed to BPFUNC.

Tags: baseprime, math, primes
Subscribe
  • Post a new comment

    Error

    default userpic

    Your reply will be screened

    Your IP address will be recorded 

    When you submit the form an invisible reCAPTCHA check will be performed.
    You must follow the Privacy Policy and Google Terms of use.
  • 7 comments