diff --git a/ats-src/bench.dats b/ats-src/bench.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/bench.dats
@@ -0,0 +1,43 @@
+#include "share/atspre_staload.hats"
+#include "share/HATS/atspre_staload_prelude.hats"
+#include "share/HATS/atspre_staload_libats_ML.hats"
+#include "share/HATS/atslib_staload_libats_libc.hats"
+#include "ats-src/combinatorics-internal.dats"
+#include "$PATSHOMELOCS/ats-bench-0.2.3/bench.dats"
+
+fun factorial_bench() : void =
+  {
+    val x = fact(160)
+    val _ = intinf_free(x)
+  }
+
+fun double_factorial_bench() : void =
+  {
+    val x = dfact(79)
+    val _ = intinf_free(x)
+  }
+
+fun choose_bench() : void =
+  {
+    val x = choose(322, 16)
+    val _ = intinf_free(x)
+  }
+
+fun catalan_bench() : void =
+  {
+    val x = catalan(300)
+    val _ = intinf_free(x)
+  }
+
+val factorial_delay: io = lam () => factorial_bench()
+val double_factorial_delay: io = lam () => double_factorial_bench()
+val choose_delay: io = lam () => double_factorial_bench()
+val catalan_delay: io = lam () => catalan_bench()
+
+implement main0 () =
+  {
+    val _ = print_slope("factorial", 12, factorial_delay)
+    val _ = print_slope("double factorial", 12, double_factorial_delay)
+    val _ = print_slope("choose", 13, choose_delay)
+    val _ = print_slope("catalan", 9, catalan_delay)
+  }
diff --git a/ats-src/combinatorics-internal.dats b/ats-src/combinatorics-internal.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/combinatorics-internal.dats
@@ -0,0 +1,189 @@
+#include "share/atspre_staload.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mydepies.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mylibies.hats"
+
+staload "$PATSHOMELOCS/atscntrb-hx-intinf/SATS/intinf_vt.sats"
+staload UN = "prelude/SATS/unsafe.sats"
+
+infixr (->) ->>
+
+stadef ->> (b1: bool, b2: bool) = ~b1 || b2
+
+// See [here](http://mathworld.wolfram.com/Derangement.html). I'm not sure how
+// fast this is, but it *seems* to be faster than the Haskell version so that's
+// good.
+fn derangements {n:nat} .<n>. (n : int(n)) : Intinf =
+  let
+    fun loop { n : nat | n > 1 }{ i : nat | i <= n } .<n-i>. (n : int(n), i : int(i), n1 : Intinf, n2 : Intinf) : Intinf =
+      if i < n then
+        let
+          var x = add_intinf0_intinf1(n2, n1)
+          var y = mul_intinf0_int(x, i)
+        in
+          loop(n, i + 1, y, n1)
+        end
+      else
+        let
+          var x = add_intinf0_intinf1(n2, n1)
+          val _ = intinf_free(n1)
+          var y = mul_intinf0_int(x, i)
+        in
+          y
+        end
+  in
+    case+ n of
+      | 0 => int2intinf(1)
+      | 1 =>> int2intinf(0)
+      | 2 =>> int2intinf(1)
+      | n =>> loop(n - 1, 2, int2intinf(1), int2intinf(0))
+  end
+
+dataprop fact_p(int, int) =
+  | fact_p_base(0, 1) of ()
+  | {n:nat}{r:int}{rn:int} fact_p_ind(n + 1, rn) of (fact_p(n, r), MUL(r, n + 1, rn))
+
+stacst fact_b : (int, int) -> bool
+
+stacst mul_b : (int, int, int) -> bool
+
+extern
+praxi fact_b_base : [fact_b(0,1)] unit_p
+
+extern
+praxi mul_b_base0 {n:int} : [mul_b(n,1,n)] unit_p
+
+extern
+praxi mul_b_base1 {n:int} : [mul_b(1,n,n)] unit_p
+
+extern
+praxi mul_b_ind0 {n:int}{m:int}{nm:int} : [mul_b(n,m,nm) ->> mul_b(n+1,m,m+nm)] unit_p
+
+// I have no idea how to actually use this proof
+// I think I need a proof-level function??
+extern
+praxi mul_b_ind1 {n:int}{m:int}{nm:int} : [mul_b(n,m,nm) ->> mul_b(n,m+1,n+nm)] unit_p
+
+extern
+praxi fact_b_ind {n:nat}{r:int}{rn:int} : [fact_b(n,r) && mul_b(r,n+1,rn) ->> fact_b(n+1,rn)] unit_p
+
+extern
+fun fact_v {n:nat} (n : int(n)) : [r:int] (fact_p(n, r) | intinf(r))
+
+extern
+fun imul {m:int}{n:int}{o:int} (x : int(m), y : int(m)) : (MUL(m, n, o) | int(o))
+
+// the fancy proof stuff isn't that useful, but it gets us a tail-recursive (?)
+// implementation which might be good (?)
+// TODO - imul_intinf0_int function
+fun fact {n:nat} .<n>. (k : int(n)) : intinfGte(1) =
+  case+ k of
+    | 0 => int2intinf(1)
+    | 1 => int2intinf(1)
+    | k =>> $UN.castvwtp0(mul_intinf0_int(fact(k - 1), k))
+
+// Double factorial http://mathworld.wolfram.com/DoubleFactorial.html
+fun dfact {n:nat} .<n>. (k : int(n)) : Intinf =
+  case+ k of
+    | 0 => int2intinf(1)
+    | 1 => int2intinf(1)
+    | k =>> let
+      var x = dfact(k - 2)
+      var y = mul_intinf0_int(x, k)
+    in
+      y
+    end
+
+// Number of permutations on n objects using k at a time.
+fn permutations {n:nat}{ k : nat | k <= n }(n : int(n), k : int(k)) : Intinf =
+  let
+    var x = fact(n)
+    var y = fact(n - k)
+    var z = div_intinf0_intinf1(x, y)
+    val _ = intinf_free(y)
+  in
+    z
+  end
+
+// Catalan numbers, indexing starting at zero.
+fn catalan {n:nat}(n : int(n)) : Intinf =
+  let
+    fun numerator_loop { i : nat | i > 1 } .<i>. (i : int(i)) : [ n : nat | n > 0 ] intinf(n) =
+      case+ i of
+        | 2 => int2intinf(n + 2)
+        | i =>> let
+          var x = numerator_loop(i - 1)
+          var y = mul_intinf0_int(x, n + i)
+        in
+          $UN.castvwtp0(y)
+        end
+  in
+    case+ n of
+      | 0 => int2intinf(1)
+      | 1 => int2intinf(1)
+      | k =>> let
+        var x = numerator_loop(k)
+        var y = fact(k)
+        var z = div_intinf0_intinf1(x, y)
+        val _ = intinf_free(y)
+      in
+        $UN.castvwtp0(z)
+      end
+  end
+
+// Number of permutations on n objects using k at a time.
+fn choose {n:nat}{ m : nat | m <= n }(n : int(n), k : int(m)) : Intinf =
+  let
+    fun numerator_loop { m : nat | m > 1 } .<m>. (i : int(m)) : [ n : nat | n > 0 ] intinf(n) =
+      case+ i of
+        | 1 => int2intinf(n)
+        | 2 => $UN.castvwtp0(int2intinf((n - 1) * n))
+        | i =>> let
+          var x = numerator_loop(i - 1)
+          var y = mul_intinf0_int(x, n + 1 - i)
+        in
+          $UN.castvwtp0(y)
+        end
+  in
+    case+ k of
+      | 0 => int2intinf(1)
+      | 1 => int2intinf(n)
+      | k =>> let
+        var x = numerator_loop(k)
+        var y = fact(k)
+        var z = div_intinf0_intinf1(x, y)
+        val _ = intinf_free(y)
+      in
+        $UN.castvwtp0(z)
+      end
+  end
+
+// Sterling numbers of the second kind
+fn sterling {n:nat}{ k : nat | k <= n }(n : int(n), k : int(k)) : Intinf =
+  let
+    fun numerator_loop {j:nat}(j : int(j), acc : Intinf) : Intinf =
+      acc
+  in
+    int2intinf(0)
+  end
+
+// TODO stirling numbers of the second kind.
+// Bell numbers. These can't be called via the FFI because of the mutually
+// recursive functions, so we should probably think of something else.
+fnx bell {n:nat}(n : int(n)) : [ n : nat | n > 0 ] intinf(n) =
+  case- n of
+    | 0 => int2intinf(1)
+    | n when n >= 0 =>> sum_loop(n, n)
+and sum_loop {n:nat}{ m : nat | m >= 1 && m <= n } .<m>. (n : int(n), i : int(m)) : [ n : nat | n > 0 ] intinf(n) =
+  case+ i of
+    | 1 => int2intinf(1)
+    | i =>> let
+      var p = sum_loop(n, i - 1)
+      var b = bell(i)
+      var c = choose(n, i)
+      var pre_ret = mul_intinf0_intinf1(c, b)
+      var ret = add_intinf0_intinf1(pre_ret, p)
+      val _ = intinf_free(b)
+      val _ = intinf_free(p)
+    in
+      $UN.castvwtp0(ret)
+    end
diff --git a/ats-src/combinatorics.dats b/ats-src/combinatorics.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/combinatorics.dats
@@ -0,0 +1,46 @@
+#define ATS_MAINATSFLAG 1
+
+#include "share/atspre_staload.hats"
+#include "ats-src/combinatorics-internal.dats"
+
+extern
+fun choose_ats {n:nat}{ m : nat | m <= n } : (int(n), int(m)) -> Intinf =
+  "mac#"
+
+extern
+fun permutations_ats {n:nat}{ m : nat | m <= n } : (int(n), int(m)) -> Intinf =
+  "mac#"
+
+extern
+fun double_factorial_ats {n:nat} : int(n) -> Intinf =
+  "mac#"
+
+extern
+fun factorial_ats {n:nat} : int(n) -> Intinf =
+  "mac#"
+
+extern
+fun catalan_ats {n:nat} : int(n) -> Intinf =
+  "mac#"
+
+extern
+fun derangements_ats {n:nat} : int(n) -> Intinf =
+  "mac#"
+
+implement choose_ats (n, k) =
+  choose(n, k)
+
+implement permutations_ats (n, k) =
+  choose(n, k)
+
+implement double_factorial_ats (m) =
+  dfact(m)
+
+implement factorial_ats (m) =
+  fact(m)
+
+implement catalan_ats (n) =
+  catalan(n)
+
+implement derangements_ats (n) =
+  derangements(n)
diff --git a/ats-src/number-theory-internal.dats b/ats-src/number-theory-internal.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/number-theory-internal.dats
@@ -0,0 +1,262 @@
+#include "share/atspre_staload.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mydepies.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mylibies.hats"
+#include "ats-src/numerics-internal.dats"
+
+staload "prelude/SATS/integer.sats"
+staload UN = "prelude/SATS/unsafe.sats"
+staload "$PATSHOMELOCS/atscntrb-hx-intinf/SATS/intinf_vt.sats"
+
+#define ATS_MAINATSFLAG 1
+
+// m | n
+fn divides(m : int, n : int) :<> bool =
+  n % m = 0
+
+// Euclid's algorithm
+fnx gcd {k:nat}{l:nat}(m : int(l), n : int(k)) : int =
+  if n > 0 then
+    gcd(n, witness(m % n))
+  else
+    m
+
+fn lcm {k:nat}{l:nat}(m : int(l), n : int(k)) : int =
+  (m / gcd(m, n)) * n
+
+fn is_coprime {k:nat}{l:nat}(m : int(l), n : int(k)) : bool =
+  gcd(m, n) = 1
+
+// stream all divisors of an integer.
+fn divisors(n : intGte(1)) : stream_vt(int) =
+  case+ n of
+    | 1 => $ldelay(stream_vt_cons(1, $ldelay(stream_vt_nil)))
+    | _ => let
+      fun loop { k : nat | k > 0 }{ m : nat | m > 0 }(n : int(k), acc : int(m)) : stream_vt(int) =
+        if acc >= sqrt_int(n) then
+          if n % acc = 0 then
+            if n / acc != acc then
+              let
+                var x: int = n / acc
+              in
+                $ldelay(stream_vt_cons(acc, $ldelay(stream_vt_cons(x, $ldelay(stream_vt_nil)))))
+              end
+            else
+              let
+                
+              in
+                $ldelay(stream_vt_cons(acc, $ldelay(stream_vt_nil)))
+              end
+          else
+            $ldelay(stream_vt_nil)
+        else
+          if n % acc = 0 then
+            let
+              var x: int = n / acc
+            in
+              $ldelay(stream_vt_cons(acc, $ldelay(stream_vt_cons(x, (loop(n, acc + 1))))))
+            end
+          else
+            loop(n, acc + 1)
+    in
+      loop(n, 1)
+    end
+
+// prime divisors of an integer
+fn prime_divisors(n : intGte(1)) : stream_vt(int) =
+  stream_vt_filter_cloptr(divisors(n), lam x => is_prime($UN.cast(x)))
+
+// if n >= 0, p > 1, then n/p >=
+fn div_gt_zero(n : intGte(0), p : intGt(1)) : intGte(0) =
+  $UN.cast(n / p)
+
+// TODO require that p be prime
+fun exp_mod_prime(a : intGte(0), n : intGte(0), p : intGt(1)) : int =
+  let
+    var a1 = a % p
+    var n1 = n % (p - 1)
+  in
+    case+ a of
+      | 0 => 0
+      | x =>> 
+        begin
+          if n > 0 then
+            let
+              var n2: intGte(0) = $UN.cast(half(n1))
+              var i2 = n1 % 2
+              var sq_a: intGte(0) = $UN.cast(a * a % p)
+            in
+              if i2 = 0 then
+                exp_mod_prime(sq_a, n2, p)
+              else
+                let
+                  var y = a * exp_mod_prime(sq_a, n2, p)
+                in
+                  y
+                end
+            end
+          else
+            1
+        end
+  end
+
+// Jacobi symbol for positive integers. See here: http://mathworld.wolfram.com/JacobiSymbol.html
+// I'm pretty sure this is broken in some way, though I'm not really sure why.
+fun jacobi(a : intGte(0), n : Odd) : int =
+  let
+    // TODO make this take p prime only.
+    fun legendre { p : int | p >= 2 }(a : intGte(0), p : int(p)) : intBtwe(~1, 1) =
+      case+ p % a of
+        | 0 => 0
+        | _ => let
+          var i = exp_mod_prime(a, (p - 1) / 2, p)
+        in
+          case+ i of
+            | i when i % (p - 1) = 0 => ~1
+            | i when i % p = 0 => 0
+            | _ => 1
+        end
+    
+    fun get_multiplicity(n : intGte(0), p : intGt(1)) : intGte(0) =
+      case+ n % p of
+        | 0 => 1 + get_multiplicity(div_gt_zero(n, p), p)
+        | _ => 0
+    
+    fun loop { m : int | m > 1 }(acc : int(m)) : int =
+      if acc > n then
+        1
+      else
+        if a % acc = 0 && is_prime(acc) then
+          loop(acc + 1) * exp(legendre(acc, n), get_multiplicity(a, acc))
+        else
+          loop(acc + 1)
+  in
+    loop(2)
+  end
+
+fn count_divisors(n : intGte(1)) : int =
+  stream_vt_length(divisors(n))
+
+vtypedef pair = @{ first = int, second = int }
+
+// aka σ in number theory
+fn sum_divisors(n : intGte(1)) : int =
+  let
+    fun loop { k : nat | k > 0 }{ m : nat | m > 0 }(n : int(k), acc : int(m)) : int =
+      if acc >= sqrt_int(n) then
+        if n % acc = 0 then
+          if n / acc != acc then
+            let
+              var x: int = n / acc
+            in
+              acc + x
+            end
+          else
+            acc
+        else
+          0
+      else
+        if n % acc = 0 then
+          let
+            var x: int = n / acc
+          in
+            acc + x + loop(n, acc + 1)
+          end
+        else
+          loop(n, acc + 1)
+  in
+    loop(n, 1)
+  end
+
+fn is_perfect(n : intGt(1)) : bool =
+  sum_divisors(n) = n
+
+fun rip { n : nat | n > 0 }{ p : nat | p > 0 } .<n>. (n : int(n), p : int(p)) :<> [ r : nat | r <= n && r > 0 ] int(r) =
+  if n % p != 0 then
+    n
+  else
+    if n / p > 0 then
+      let
+        var n1 = n / p
+      in
+        if n1 < n then
+          $UN.cast(rip(n1, p))
+        else
+          1
+      end
+    else
+      1
+
+fun prime_factors(n : intGte(1)) : stream_vt(int) =
+  let
+    fun loop { k : nat | k > 0 }{ m : nat | m > 0 }(n : int(k), acc : int(m)) : stream_vt(int) =
+      if acc >= n then
+        if is_prime(n) then
+          $ldelay(stream_vt_cons(n, $ldelay(stream_vt_nil)))
+        else
+          $ldelay(stream_vt_nil)
+      else
+        if n % acc = 0 && is_prime(acc) then
+          if n / acc > 0 then
+            $ldelay(stream_vt_cons(acc, loop(rip(n, acc), 1)))
+          else
+            $ldelay(stream_vt_cons(acc, $ldelay(stream_vt_nil)))
+        else
+          loop(n, acc + 1)
+  in
+    loop(n, 1)
+  end
+
+// distinct prime divisors
+fn little_omega(n : intGte(1)) : int =
+  let
+    fun loop { k : nat | k > 0 }{ m : nat | m > 0 }(n : int(k), acc : int(m)) :<!ntm> int =
+      if acc >= n then
+        if is_prime(n) then
+          1
+        else
+          0
+      else
+        if n % acc = 0 && is_prime(acc) then
+          if n / acc > 0 then
+            1 + loop(rip(n, acc), 1)
+          else
+            1
+        else
+          loop(n, acc + 1)
+  in
+    loop(n, 1)
+  end
+
+// Euler's totient function.
+fn totient(n : intGte(1)) : int =
+  case+ n of
+    | 1 => 1
+    | n =>> let
+      fn adjust_contents(x : pair, y : int) : pair =
+        @{ first = g0int_mul(x.first, y - 1), second = g0int_mul(x.second, y) }
+      
+      var x: stream_vt(int) = prime_factors(n)
+      var empty_pair = @{ first = 1, second = 1 } : pair
+      var y = stream_vt_foldleft_cloptr(x, empty_pair, lam (acc, next) => adjust_contents(acc, next)) : pair
+    in
+      g0int_div(g0int_mul(n, y.first), y.second)
+    end
+
+// The sum of all φ(m) for m between 1 and n. Note the use of refinement types
+// to prevent 0 from being passed as an argument. This function is actually
+// slower than the Haskell equivalent, as it uses a naïve algorithm.
+fn totient_sum(n : intGte(1)) : Intinf =
+  let
+    fnx loop { n : nat | n >= 1 }{ m : nat | m >= n } .<m-n>. (i : int(n), bound : int(m)) : Intinf =
+      if i < bound then
+        let
+          var x = loop(i + 1, bound)
+          var y = add_intinf0_int(x, witness(totient(i)))
+        in
+          y
+        end
+      else
+        int2intinf(witness(totient(i)))
+  in
+    loop(1, n)
+  end
diff --git a/ats-src/number-theory.dats b/ats-src/number-theory.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/number-theory.dats
@@ -0,0 +1,57 @@
+#include "ats-src/number-theory-internal.dats"
+
+extern
+fun totient_ats { k : nat | k >= 2 }(int(k)) : int =
+  "mac#"
+
+extern
+fun count_divisors_ats { k : nat | k >= 2 }(int(k)) : int =
+  "mac#"
+
+extern
+fun little_omega_ats { n : nat | n > 0 } : int(n) -> int =
+  "mac#"
+
+extern
+fun sum_divisors_ats : { n : nat | n > 1 } int(n) -> int =
+  "mac#"
+
+extern
+fun jacobi_ats : (intGte(0), Odd) -> int =
+  "mac#"
+
+extern
+fun is_perfect_ats : intGt(1) -> bool =
+  "mac#"
+
+extern
+fun totient_sum_ats : intGte(1) -> Intinf =
+  "mac#"
+
+extern
+fun coprime_ats {k:nat}{n:nat} : (int(k), int(n)) -> bool =
+  "mac#"
+
+implement sum_divisors_ats (m) =
+  sum_divisors(m)
+
+implement count_divisors_ats (n) =
+  count_divisors(n)
+
+implement totient_ats (n) =
+  totient(n)
+
+implement little_omega_ats (n) =
+  little_omega(n)
+
+implement is_perfect_ats (n) =
+  is_perfect(n)
+
+implement jacobi_ats (m, n) =
+  jacobi(m, $UN.cast(n))
+
+implement totient_sum_ats (n) =
+  totient_sum(n)
+
+implement coprime_ats (m, n) =
+  is_coprime(m, n)
diff --git a/ats-src/numerics-internal.dats b/ats-src/numerics-internal.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/numerics-internal.dats
@@ -0,0 +1,126 @@
+#include "share/atspre_staload.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mydepies.hats"
+#include "$PATSHOMELOCS/atscntrb-hx-intinf/mylibies.hats"
+
+staload "$PATSHOMELOCS/atscntrb-hx-intinf/SATS/intinf_vt.sats"
+staload "libats/libc/SATS/math.sats"
+staload UN = "prelude/SATS/unsafe.sats"
+
+// Existential types for even and odd numbers. These are only usable with the
+// ATS library.
+typedef Even = [n:nat] int(2*n)
+typedef Odd = [n:nat] int(2*n+1)
+
+// These types work... less well. I'm not sure what the story is with
+// multiplicative constraints in ATS, but in general they're unsolvable due to
+// Gödel's incompleteness theorem.
+typedef gprime(tk: tk, p: int) = { m, n : nat | m < 1 && m <= n && n < p && m*n != p && p > 1 } g1int(tk, p)
+typedef prime(p: int) = gprime(int_kind, p)
+typedef Prime = [p:nat] prime(p)
+
+fn witness(n : int) :<> [m:nat] int(m) =
+  $UN.cast(n)
+
+// Fast computation of Fibonacci numbers via GMP bindings.
+fun fib_gmp(n : intGte(0)) : Intinf =
+  let
+    var z = ptr_alloc()
+    var x = g0int2uint(n + 1)
+    val _ = $GMP.mpz_init(!(z.2))
+    val _ = $GMP.mpz_fib_uint(!(z.2), x)
+  in
+    $UN.castvwtp0(z)
+  end
+
+// Fast integer exponentiation. This performs O(log n) multiplications. This
+// function is mostly useful for exponentiation in modular arithmetic, as
+// it can overflow.
+fun exp {n:nat} .<n>. (x : int, n : int(n)) : int =
+  case+ x of
+    | 0 => 0
+    | x =>> 
+      begin
+        if n > 0 then
+          let
+            var n2 = half(n)
+            var i2 = n % 2
+          in
+            if i2 = 0 then
+              exp(x * x, n2)
+            else
+              let
+                var y = x * exp(x * x, n2)
+              in
+                y
+              end
+          end
+        else
+          1
+      end
+
+// Fast integer exponentiation, that mostly works as we would like.
+fun big_exp {n:nat} .<n>. (x : Intinf, n : int(n)) : Intinf =
+  if compare_intinf_int(x, 0) = 0 then
+    x
+  else
+    if n > 0 then
+      let
+        var n2 = half(n)
+        var i2 = n % 2
+      in
+        if i2 = 0 then
+          let
+            // FIXME copy it
+            var x0 = abs_intinf1(x)
+            var c = mul_intinf0_intinf1(x0, x)
+            val _ = intinf_free(x)
+          in
+            big_exp(c, n2)
+          end
+        else
+          let
+            var x0 = abs_intinf1(x)
+            var c0 = mul_intinf0_intinf1(x0, x)
+            var c1 = big_exp(c0, n2)
+            var c = mul_intinf0_intinf1(c1, x)
+            val _ = intinf_free(x)
+          in
+            c
+          end
+      end
+    else
+      (intinf_free(x) ; int2intinf(1))
+
+// square root is bounded for bounded k.
+fn sqrt_int(k : intGt(0)) :<> [m:nat] int(m) =
+  let
+    var bound: int = g0float2int(sqrt_float(g0int2float(k)))
+  in
+    witness(bound)
+  end
+
+// function to check primality
+fn is_prime(k : intGt(0)) :<> bool =
+  case+ k of
+    | 1 => false
+    | k => 
+      begin
+        let
+          fun loop {n:nat}{m:nat} .<max(0,m-n)>. (i : int(n), bound : int(m)) :<> bool =
+            if i < bound then
+              if k % i = 0 then
+                false
+              else
+                loop(i + 1, bound)
+            else
+              if i = bound then
+                if k % i = 0 then
+                  false
+                else
+                  true
+              else
+                true
+        in
+          loop(2, sqrt_int(k))
+        end
+      end
diff --git a/ats-src/numerics.dats b/ats-src/numerics.dats
new file mode 100644
--- /dev/null
+++ b/ats-src/numerics.dats
@@ -0,0 +1,30 @@
+#define ATS_DYNLOADFLAG 0
+
+#include "share/atspre_staload.hats"
+#include "ats-src/numerics-internal.dats"
+
+staload "$PATSHOMELOCS/atscntrb-hx-intinf/SATS/intinf_vt.sats"
+
+%{^
+#define ATS_MEMALLOC_LIBC
+#include "ccomp/runtime/pats_ccomp_memalloc_libc.h"
+#include "ccomp/runtime/pats_ccomp_runtime_memalloc.c"
+%}
+
+extern
+fun is_prime_ats { n : nat | n > 0 } : int(n) -> bool =
+  "mac#"
+
+extern
+fun exp_ats {m:nat} : ([n:nat] int(n), int(m)) -> int =
+  "mac#"
+
+extern
+fun fib_ats : intGte(0) -> Intinf =
+  "mac#"
+
+implement is_prime_ats (n) =
+  is_prime(n)
+
+implement exp_ats (m, n) =
+  exp(m, n)
diff --git a/cbits/combinatorics.c b/cbits/combinatorics.c
--- a/cbits/combinatorics.c
+++ b/cbits/combinatorics.c
@@ -1,7 +1,7 @@
 /*
 **
 ** The C code is generated by [ATS/Postiats-0-3-10]
-** The starting compilation time is: 2018-5-13: 18h:15m
+** The starting compilation time is: 2018-5-16: 23h:22m
 **
 */
 
diff --git a/cbits/number-theory.c b/cbits/number-theory.c
--- a/cbits/number-theory.c
+++ b/cbits/number-theory.c
@@ -1,7 +1,7 @@
 /*
 **
 ** The C code is generated by [ATS/Postiats-0-3-10]
-** The starting compilation time is: 2018-5-13: 18h:15m
+** The starting compilation time is: 2018-5-16: 23h:22m
 **
 */
 
diff --git a/cbits/numerics.c b/cbits/numerics.c
--- a/cbits/numerics.c
+++ b/cbits/numerics.c
@@ -1,7 +1,7 @@
 /*
 **
 ** The C code is generated by [ATS/Postiats-0-3-10]
-** The starting compilation time is: 2018-5-13: 18h:15m
+** The starting compilation time is: 2018-5-16: 23h:22m
 **
 */
 
diff --git a/fast-arithmetic.cabal b/fast-arithmetic.cabal
--- a/fast-arithmetic.cabal
+++ b/fast-arithmetic.cabal
@@ -1,6 +1,6 @@
 cabal-version: 1.18
 name: fast-arithmetic
-version: 0.6.0.5
+version: 0.6.0.6
 license: BSD3
 license-file: LICENSE
 copyright: Copyright: (c) 2018 Vanessa McHale
@@ -16,6 +16,7 @@
 extra-source-files:
     atspkg.dhall
     pkg.dhall
+    ats-src/*.dats
     .atspkg/contrib/ats-includes-0.3.10/ccomp/runtime/*.h
     .atspkg/contrib/ats-includes-0.3.10/ccomp/runtime/*.c
     .atspkg/contrib/ats-includes-0.3.10/prelude/CATS/*.cats
@@ -51,7 +52,7 @@
                   .atspkg/contrib/ats-includes-0.3.10/ .atspkg/contrib
     ghc-options: -Wall -optc-mtune=native -optc-flto -optc-O3
     build-depends:
-        base >=4.7 && <5,
+        base >=4.5 && <5,
         composition-prelude >=1.2.0.0,
         gmpint -any
     
diff --git a/pkg.dhall b/pkg.dhall
--- a/pkg.dhall
+++ b/pkg.dhall
@@ -1,4 +1,4 @@
-let prelude = https://raw.githubusercontent.com/vmchale/atspkg/master/dhall/atspkg-prelude.dhall
+let prelude = https://raw.githubusercontent.com/vmchale/atspkg/master/ats-pkg/dhall/atspkg-prelude.dhall
 
 in λ(x : List Integer) →
   prelude.makeHsPkg { x = x, name = "fast-arithmetic" }