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
  • Post a new comment

    Error

    default userpic

    Your IP address will be recorded 

  • 7 comments