logcount: new function.

This is in ANSI CL; potentially useful and hard to
implement efficiently in user code.

* arith.c (logcount): New function.

* eval.c (eval_init): Register logcount intrinsic.

* lib.h (logcount): Declared.

* mpi/mi.c (s_mp_count_ones): New static function.
(mp_count_ones): New function.

* mpi/mpi.h (mp_count_ones): Declared.

* txr.1: Documented.

6 files changed, 121 insertions, 0 deletions

diff --git a/arith.c b/arith.c
index e9f8024c..7731558d 100644
--- a/arith.c
+++ b/arith.c
@@ -2697,6 +2697,46 @@ val maskv(struct args *bits)
   return accum;
 }
 
+val logcount(val n)
+{
+  val self = lit("logcount");
+
+  switch (type(n)) {
+  case NUM:
+  case CHR:
+    {
+      int_ptr_t c = c_num(n);
+      uint_ptr_t d = c;
+      if (c < 0)
+        d = ~d;
+#if SIZEOF_PTR == 8
+      d = ((d & 0xAAAAAAAAAAAAAAAA) >>  1) + (d & 0x5555555555555555);
+      d = ((d & 0xCCCCCCCCCCCCCCCC) >>  2) + (d & 0x3333333333333333);
+      d = ((d & 0xF0F0F0F0F0F0F0F0) >>  4) + (d & 0x0F0F0F0F0F0F0F0F);
+      d = ((d & 0xFF00FF00FF00FF00) >>  8) + (d & 0x00FF00FF00FF00FF);
+      d = ((d & 0xFFFF0000FFFF0000) >> 16) + (d & 0x0000FFFF0000FFFF);
+      d = ((d & 0xFFFFFFFF00000000) >> 32) + (d & 0x00000000FFFFFFFF);
+#elif SIZEOF_PTR == 4
+      d = ((d & 0xAAAAAAAA) >>  1) + (d & 0x55555555);
+      d = ((d & 0xCCCCCCCC) >>  2) + (d & 0x33333333);
+      d = ((d & 0xF0F0F0F0) >>  4) + (d & 0x0F0F0F0F);
+      d = ((d & 0xFF00FF00) >>  8) + (d & 0x00FF00FF);
+      d = ((d & 0xFFFF0000) >> 16) + (d & 0x0000FFFF);
+      return unum(d);
+#else
+#error fixme: only 4 or 8 byte pointers supported
+#endif
+    }
+  case BGNUM:
+    {
+      mp_size co = mp_count_ones(mp(n));
+      return unum(co);
+    }
+  default:
+    uw_throwf(error_s, lit("~a: non-integral operand ~s"), self, n, nao);
+  }
+}
+
 /*
  * Source:
  * Better Approximations to Cumulative Normal Functions
diff --git a/eval.c b/eval.c
index 98d8cfdf..e621db2d 100644
--- a/eval.c
+++ b/eval.c
@@ -6481,6 +6481,7 @@ void eval_init(void)
   reg_fun(intern(lit("bit"), user_package), func_n2(bit));
   reg_fun(intern(lit("mask"), user_package), func_n0v(maskv));
   reg_fun(intern(lit("width"), user_package), func_n1(width));
+  reg_fun(intern(lit("logcount"), user_package), func_n1(logcount));
 
   reg_fun(intern(lit("sort-group"), user_package), func_n3o(sort_group, 1));
   reg_fun(intern(lit("unique"), user_package), func_n2ov(unique, 1));
diff --git a/lib.h b/lib.h
index 1c1b9f29..1a81ed83 100644
--- a/lib.h
+++ b/lib.h
@@ -752,6 +752,7 @@ val sign_extend(val num, val nbits);
 val ash(val a, val bits);
 val bit(val a, val bit);
 val maskv(struct args *bits);
+val logcount(val n);
 val string_own(wchar_t *str);
 val string(const wchar_t *str);
 val string_utf8(const char *str);
diff --git a/mpi/mpi.c b/mpi/mpi.c
index 84308e8a..472011d4 100644
--- a/mpi/mpi.c
+++ b/mpi/mpi.c
@@ -2499,6 +2499,65 @@ mp_size mp_count_bits(mp_int *mp)
   return s_highest_bit_mp(mp);
 }
 
+static mp_size s_mp_count_ones(mp_int *mp)
+{
+  mp_size ix;
+  mp_size c;
+  mp_digit *dp = DIGITS(mp);
+
+  for (c = 0, ix = USED(mp) - 1; ix < MP_SIZE_MAX; ix--) {
+    mp_digit d = dp[ix];
+#if MP_DIGIT_SIZE == 8
+    d = ((d & 0xAAAAAAAAAAAAAAAA) >>  1) + (d & 0x5555555555555555);
+    d = ((d & 0xCCCCCCCCCCCCCCCC) >>  2) + (d & 0x3333333333333333);
+    d = ((d & 0xF0F0F0F0F0F0F0F0) >>  4) + (d & 0x0F0F0F0F0F0F0F0F);
+    d = ((d & 0xFF00FF00FF00FF00) >>  8) + (d & 0x00FF00FF00FF00FF);
+    d = ((d & 0xFFFF0000FFFF0000) >> 16) + (d & 0x0000FFFF0000FFFF);
+    d = ((d & 0xFFFFFFFF00000000) >> 32) + (d & 0x00000000FFFFFFFF);
+    c += d;
+#elif MP_DIGIT_SIZE == 4
+    d = ((d & 0xAAAAAAAA) >>  1) + (d & 0x55555555);
+    d = ((d & 0xCCCCCCCC) >>  2) + (d & 0x33333333);
+    d = ((d & 0xF0F0F0F0) >>  4) + (d & 0x0F0F0F0F);
+    d = ((d & 0xFF00FF00) >>  8) + (d & 0x00FF00FF);
+    d = ((d & 0xFFFF0000) >> 16) + (d & 0x0000FFFF);
+    c += d;
+#elif MP_DIGIT_SIZE == 2
+    d = ((d & 0xAAAA) >> 1) + (d & 0x5555);
+    d = ((d & 0xCCCC) >> 2) + (d & 0x3333);
+    d = ((d & 0xF0F0) >> 4) + (d & 0x0F0F);
+    d = ((d & 0xFF00) >> 8) + (d & 0x00FF);
+    c += d;
+#elif MP_DIGIT_SIZE == 1
+    d = ((d & 0xAA) >> 1) + (d & 0x55);
+    d = ((d & 0xCC) >> 2) + (d & 0x33);
+    d = ((d & 0xF0) >> 4) + (d & 0x0F);
+    c += d;
+#else
+#error fixme: unsupported MP_DIGIT_SIZE
+#endif
+  }
+
+  return c;
+}
+
+mp_size mp_count_ones(mp_int *mp)
+{
+  if (SIGN(mp) == MP_NEG) {
+    mp_int tmp;
+    mp_size res;
+    if ((res = mp_init_copy(&tmp, mp)) != MP_OKAY)
+      return res;
+    if ((res = s_mp_sub_d(&tmp, 1) != MP_OKAY))
+      return res;
+    res = s_mp_count_ones(&tmp);
+    mp_clear(&tmp);
+    return res;
+  }
+
+  return s_mp_count_ones(mp);
+}
+
 mp_size mp_is_pow_two(mp_int *mp)
 {
   return s_mp_ispow2(mp) >= 0;
diff --git a/mpi/mpi.h b/mpi/mpi.h
index 06adc1d4..8bd469f7 100644
--- a/mpi/mpi.h
+++ b/mpi/mpi.h
@@ -177,6 +177,7 @@ mp_err mp_to_unsigned_bin(mp_int *mp, unsigned char *str);
 mp_err mp_to_unsigned_buf(mp_int *mp, unsigned char *str, size_t size);
 
 mp_size mp_count_bits(mp_int *mp);
+mp_size mp_count_ones(mp_int *mp);
 mp_size mp_is_pow_two(mp_int *mp);
 
 #if MP_COMPAT_MACROS
diff --git a/txr.1 b/txr.1
index 48a15060..7d6e00bb 100644
--- a/txr.1
+++ b/txr.1
@@ -36263,6 +36263,25 @@ two's complement bitfield 01 denotes 1, and 10 denotes -2.
 
 The argument may be a character.
 
+.coNP Function @ logcount
+.synb
+.mets (logcount << integer )
+.syne
+.desc
+The
+.code logcount
+function considers
+.meta integer
+to have a two's complement representation. If the integer is positive,
+it returns the count of bits in that representation whose value is 1.
+If
+.meta integer
+is negative, it returns the count of zero bits instead. If
+.meta integer
+is zero, the value returned is zero.
+
+The argument may be a character.
+
 .SS* Exception Handling
 
 An