import

ChangeLog

+2014-07-28
+	* Fix %/% operator
+	* Remove doxygen auto-generated files.
+
+2013-12-03
+	* Fix c++ code to work with new gcc.
+2013-10-29
+	* Fix c++ code to work with new gcc.
+2013-09-12  Antoine Lucas.
+	* Fix Warn / Note.
+
+2013-09-10  Antoine Lucas.
+	* fix man page.
+
+2013-06-09  Antoine Lucas.
+        * fix factorize('33162879029270137')
+
+2013-06-07  Martin Maechler  <maechler@stat.math.ethz.ch>
+	big changes related to  "a/b" and "b^(-1)";
+	also: now give [NA + warning()} rather than an error when
+	there is no inverse, a division to 0, or a modulus mismatch
+
+2013-02-14  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* R/biginteger.R (duplicated.bigz, unique.bigz, all.equal.bigz):
+	new (cheap imp.).
+
+	* R/bigq.R (duplicated.bigq, unique.bigq, all.equal.bigq, mean.bigq):
+	new (cheap imp.)
+
+	* R/biginteger.R (rep.bigz): and
+	* R/bigq.R (rep.bigq): allow ('length.out' and 'each') arguments
+
+2013-02-13  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* src/bigrationalR.cc (bigrational_num): numerator(<NA>)
+
+	* R/matrixz.R (.dimsetZQ): dim(.) <- NULL  now works
+
+2013-02-04  Martin Maechler
+
+	* DESCRIPTION (Version): 0.5-4   (0.5-3 failed in R 2.14.2)
+
+2013-02-03  Antoine Lucas + Martin Maechler
+	* man/sizeinbase.Rd updated
+	* new function asNumeric()
+
+2012-05-07  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* src/bigrational.cc (^): Q^n  and Z^n now also work for negative
+	'n' and with NA's.
+
+	* R/bigq.R (as.bigz.bigq): as.bigz(as.bigq(NA)) now works
+
+2012-05-02  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* src/bigrationalR.cc (bigrational_pow, mpqz_pow):
+	* src/bigrational.cc (operator^): add   <bigq> ^ n   functionality.
+	* TODO: (one less :-)
+
+2012-05-01  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* DESCRIPTION (Version): 0.5-2
+
+	* src/bigintegerR.cc (bigI_choose, bigI_factorial): now vectorized.
+	* man/factorialZ.Rd: updated, including examples
+
+2012-01-20 Antoine Lucas
+
+	* tests/gmp-test.R: update
+
+2012-01-12  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* NAMESPACE,
+	* src/init.cc: new file, providing proper registration
+	* R/zzz.R: new gmpVersion() [to get GMP C library version].
+
+2012-01-09  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* R/biginteger.R (as.bigz): mod=NULL now is equivalent to NA;
+	 important for  'modulus(x) <- NULL' to work.
+
+2012-01-07  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* man/modulus.Rd: example, (and FIXME)
+
+	* TODO: as.bigz(integer()) and as.bigq(integer()) now both work and
+	give a 0-length object; analogously, the basic functions now work
+	with the 0-length big[qz] vectors and return 0-length.
+
+	* src/bigintegerR.cc (create_SEXP): ensure the modulus also has
+	class "bigz".
+
+2012-01-06  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* R/matrix-prods.R: new file, containing all %*%, and new
+	crossprod, tcrossprod methods.
+	<bigz> %*% <double>  now works, too.
+
+2011-12-31  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* src/bigintegerR.cc (biginteger_as_integer): new
+	(biginteger_log, biginteger_log2): new --> allowing proper log().
+
+	* R/biginteger.R (as.integer.bigz): new -- gave nonsense!
+	* R/biginteger.R (log.bigz): now work correctly also for large numbers.
+
+2011-12-29  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* src/bigintegerR.cc (bigI_frexp): new, with R interface
+	* man/frexpZ.Rd:
+
+	* R/bigq.R (Math.bigq): define a "fall-back" group method
+	* R/biginteger.R (Math.bigz): ditto
+
+	* R/biginteger.R (print.bigz, print.bigq): Mention class and length,
+	and then do *not* quote, by default.
+
+2011-12-28  Martin Maechler  <maechler@stat.math.ethz.ch>
+
+	* DESCRIPTION (Version): 0.5-0
+
+	* src/bigintegerR.cc (biginteger_max): na_remove now correct from R code
+	* src/bigrationalR.cc (bigrational_max): ditto
+
+	* R/biginteger.R (factorialZ, chooseZ): new
+
+2010-04-10  Antoine Lucas <antoinelucas@gmail.com>
+
+	* ./Changes : see for 0.4-12 and everything older.

Changes

+Sat Apr 10 2010: 0.4.12
+  * fix for sun c++ compiler
+  * add sum() functions.
+
+Tue Oct 27 2009: 0.4.11
+  * minor R & cpp warnings fix
+
+Sat Jan 24 2009: 0.4.9
+  * move printf to Rprintf function
+  * add powm function.
+
+Mon Feb 25 2008: 0.4.8
+  * correct a bug on extract matrix (in affectation x[y]<- value #...)
+  * correct a bug on extract matrix (an invalid memory read, cause invalid result erraticly)
+
+Mon Feb 17 2008: 0.4.7
+  * suppress a bug on rational matrix.
+
+Tue Sep 14 2007: 0.4-6
+  * correction on c++: supress all g++ warnings
+
+Tue Sep 14 2007: 0.4-5
+  * correction of bug ref 698 [abs.bigz]
+  * man pages corrected
+
+Tue May 12 2007: 0.4-4
+  * Other corrections of S3 generic/method consistency (add parameter ... on some functions)
+
+Tue May 10 2007: 0.4-3
+  * correction of S3 generic/method consistency (add parameter ... on some functions)
+
+Sun Jan 21 2007: 0.4-2
+
+Sun Dec 17 2006 
+  * add apply function
+
+Sun Dec 03 2006
+  * add cumsum & prod 
+  * correct bugs on NA & not finite values.
+    move NA_REAL check to R_FINITE check
+    use carefully getValueTemp as it does not contains NA flag.
+
+Sun Nov 05 2006
+  * add max & min (bigintegers)
+
+Mon Oct 09 2006 0.4-1
+  * Correction of bug on function urand.bigz()
+
+Sat Jun 17 2006 0.4-0
+  * extract functions done with matrix (x[i,j] & x[i,j] <- value)
+
+Sun May 28 2006
+  * Matrix integration & solve function in Z/nZ or Q
+
+Thu May 19 2006
+  * rewritten bigrational class with mpq_t value
+  * use of 2 new class bigvec & bigvec_q that includes vector(s) + nrow value.
+  * catch bug in fibnum & lucnum when negative values or NULL.
+
+Wed Mar 15 2006
+  * catch bug when negative integer values (like as.bigz(as.integer(-4)))
+
+Tue Mar 14 2006
+  * catch bug with c.bigz(3,NULL) (and other bugs with NULL argument
+    while add cas NILSXP in create_vector
+  * check for positive values in lucnum & fibnum.
+
+Thu Jan 20 2006 0.3-4
+  * minor changes on data Oakley
+
+Wed Jul 27 2005 0.3-3
+  * minor changes on man page
+
+Fri May 06 2005
+  Minor bug fixed...
+  * error message corrected on function pow
+  * bug fixed (a test to check 0 values with mpz_sgn instead of
+    as_long)
+  * man page modified for as.bigz
+
+Wed Mar 23 2005
+  * catch illegal base argument at as.character.bigz (fix segfault)
+  * added RFC 2409 global unique Diffie-Hellman groups as data constants
+  * typos in DESCRIPTION
+
+Thu Mar 17 2005
+  * add a base 16 output function
+  * Makevars.win
+
+Thu Jan 13 2005 
+  * bug fixed on rational numbers class & use of negatives values
+    allowed for both rational and integer class.
+
+Mon Jan 10 2005
+  * bug fixed on rational numbers class
+
+Wed Jan 5 2005
+  * First try for a class for rational numbers
+
+Fri Dec 10 2004
+  * Function Factorize
+
+Fri Dec 09 2004: 0.2-2
+
+Thu Dec 02 2004:
+   * biginteger class changed to bigz
+   * man page updated.
+
+Sat Nov 27 2004:
+   * New function c.biginteger, rep.biginteger, sizeinbase
+   * former functions nextprime, isprime, gcdex, random rewritten.
+     => no more gmp.cc, gmp.R
+   * Functions for Fibonacci and Lucas number computation !
+   * Some changes in C++ comments to have a doc with doxygen.
+
+Wed Nov 24 2004:
+   * Modifs on doc
+
+Tue Nov 09 2004: 0.2-1
+   * Suppressed macro definitions in Rinternals.h, and replaced macro calls with
+     plain function names. Trying to avoid collisions with stl definitions, i.e. length()
+
+Wed Oct 27 2004: 0.2
+   * Vectorized engine rewritten
+
+Mon Sep 27 2004: 0.1
+   * Initial revision
+

DESCRIPTION

+Package: gmp
+Version: 0.5-13
+Date: 2017-01-02
+Title: Multiple Precision Arithmetic
+Author: Antoine Lucas, Immanuel Scholz, Rainer Boehme <rb-gmp@reflex-studio.de>,
+  Sylvain Jasson <jasson@toulouse.inra.fr>,
+  Martin Maechler <maechler@stat.math.ethz.ch>
+Maintainer: Antoine Lucas <antoinelucas@gmail.com>
+Description: Multiple Precision Arithmetic (big integers and rationals,
+ prime number tests, matrix computation), "arithmetic without limitations"
+ using the C library GMP (GNU Multiple Precision Arithmetic).
+Depends: R (>= 3.0.0)
+Imports: methods
+Suggests: Rmpfr
+SystemRequirements: gmp (>= 4.2.3)
+License: GPL
+BuildResaveData: no
+LazyDataNote: not available, as we use data/*.R *and* our classes
+URL: http://mulcyber.toulouse.inra.fr/projects/gmp

Done

+* TODO's that are done (not yet: see ./TODO !)				-*-org-*-
+
+** Indexing/Subassignment:
+*** Seg.faults in
+ Z <- as.bigz(1:7);     Z[1] <- list(1) # --> seg.fault
+ Q <- as.bigq(1:7, 7);  Q[1] <- list(1) # --> seg.fault
+ Z <- as.bigz(1:7);     Z[1] <- NULL    # --> seg.fault
+ Q <- as.bigq(1:7, 7);  Q[1] <- NULL    # --> seg.fault
+
+** Not working Arithmetic
+*** qq <- as.bigq(c(0:5, NA, 7:8), 12) # shows NA fine, is.na(qq), too, but
+    numerator( qq ) gives 0 instead of NA !!
+
+** Miscellaneous
+*** unique() silently gives nonsense !
+*** duplicated()  ditto
+** Matrix stuff
+***  as.matrix(as.bigq(1:4)),  as.matrix(as.bigz(1:4)) now (~ 2012) work
+**** NB: to go "back", now have asNumeric()  {instead of expecting as.matrix() to give numeric}
+
+
+
+

NAMESPACE

+useDynLib(gmp, .registration=TRUE)
+
+importMethodsFrom("methods", coerce)
+importFrom("methods",
+           setOldClass, setGeneric, setMethod)
+
+## instead of exportPattern() :
+## MM: but *really* should export considerably less !! (i.e. use  S3method(.) more !
+export(
+       "%*%", # unfortunately need to make it into an S3 generic
+       "crossprod", "tcrossprod",#<- *STILL* needed .. unless S4 and we exportMethods() for them
+       ## "abs.bigq", "abs.bigz",
+       "add.bigq", "add.bigz",
+       "apply", "apply.bigq", "apply.bigz", "apply.default",
+       "as.bigq", "as.bigz", "as.bigz.bigq",
+       ## "as.character.bigq", "as.character.bigz",
+       ## "as.double.bigq", "as.double.bigz",
+       ## "as.matrix.bigq", "as.matrix.bigz",
+       ## "as.vector.bigq", "as.vector.bigz",
+       ## "c.bigq", "c.bigz", "cbind.bigq", "cbind.bigz",  "rbind.bigq", "rbind.bigz",
+       ## "cumsum.bigq", "cumsum.bigz",
+       "denominator", "denominator<-",
+       ## "dim<-.bigq", "dim.bigq", "dim<-.bigz", "dim.bigz",
+       "div.bigq", "div.bigz", "divq.bigz",
+       "Eulerian", "Eulerian.all",
+       "Stirling1", "Stirling1.all", "Stirling2", "Stirling2.all",
+       "factorize",
+       "factorialZ", "chooseZ", "frexpZ",
+       "fibnum", "fibnum2", "formatN",
+       "gcd", "gcd.bigz", "gcd.default", "gcdex",
+       "gmpVersion",
+       "inv.bigz", ## "inv",
+       ## "is.na.bigq", "is.na.bigz",
+       "is.bigz", "is.bigq", "is.matrixZQ",
+       "isprime", "is.whole",
+       "lcm.bigz", "lcm.default",
+       ## "length<-.bigq", "length.bigq", "length<-.bigz", "length.bigz",
+       "log10.bigz", "log2.bigz", "log.bigz",
+       "lucnum", "lucnum2",
+       "matrix", "matrix.bigq", "matrix.bigz", "matrix.default",
+       "max.bigq", "max.bigz", "min.bigq", "min.bigz",
+       "mod.bigz",
+       "modulus", "modulus<-", # "modulus<-.bigz", "modulus.bigz",
+       "mul.bigq", "mul.bigz",
+       "ncol.bigq", "ncol.bigz", "nrow.bigq", "nrow.bigz",
+       "nextprime",
+       "numerator", "numerator<-",
+       ## "pow",
+       "pow.bigz", "powm", "pow.bigq",
+       ## "print.bigq", "print.bigz",
+       ## "prod", "prod.default",
+       "prod.bigq", "prod.bigz",
+       "rep.bigq", "rep.bigz",
+       ## "sign.bigq", "sign.bigz",
+       "sizeinbase",
+       "solve.bigq", "solve.bigz",
+       "sub.bigq", "sub.bigz",
+       "sum.bigq", "sum.bigz",
+       ## "t.bigq", "t.bigz",
+       "urand.bigz")
+
+## For Rmpfr (and possibly similar packages)
+export(.as.bigz, ..as.bigz, .as.char.bigz)
+## C routines for other packages to .Call(..) .. hmm, convenient
+## but then also needs mention on a help page:
+## FIXME (2013-09-11): drop these after a while {when dependent packages
+## ----- all use the .as.* and ..as.* substitutes:
+export(biginteger_as, biginteger_as_character)
+
+S3method(gcd,default)
+#S3method(lcm,default)
+S3method(is.whole,default)
+
+S3method("+",bigz)
+S3method("-",bigz)
+S3method("*",bigz)
+S3method("/",bigz)
+S3method("%%",bigz)
+S3method("%/%",bigz)
+S3method("^",bigz)
+S3method("<",bigz)
+S3method("<=",bigz)
+S3method("==",bigz)
+S3method(">=",bigz)
+S3method(">",bigz)
+S3method("!=",bigz)
+S3method("[<-",bigz)
+S3method("[",bigz)
+S3method("[[<-",bigz)
+S3method("[[",bigz)
+##S3method(add,bigz)
+S3method(sum,bigz)
+S3method(sum,bigq)
+#S3method(pow,bigz)
+S3method(gcd,bigz)
+#S3method(lcm,bigz)
+#S3method(inv,bigz)
+S3method(log,bigz)
+S3method(log10,bigz)
+S3method(log2,bigz)
+S3method(c,bigz)
+S3method(cbind,bigz)
+S3method(rbind,bigz)
+S3method(dim,bigz)
+S3method("dim<-", bigz)
+S3method(is.finite,bigz)
+S3method(is.infinite,bigz)
+S3method(is.na,bigz)
+S3method(is.whole,bigz)
+S3method(length,bigz)
+S3method("length<-", bigz)
+S3method(max,bigz)
+S3method(min,bigz)
+S3method(rep,bigz)
+S3method(print,bigz)
+##S3method(urand,bigz)
+S3method(duplicated,bigz)
+S3method(unique,bigz)
+S3method(all.equal,bigz)
+
+S3method(modulus, bigz)
+S3method("modulus<-", bigz)
+
+S3method("+",bigq)
+S3method("-",bigq)
+S3method("*",bigq)
+S3method("/",bigq)
+S3method("^",bigq)
+S3method("<",bigq)
+S3method("<=",bigq)
+S3method("==",bigq)
+S3method(">=",bigq)
+S3method(">",bigq)
+S3method("!=",bigq)
+S3method("[<-",bigq)
+S3method("[",bigq)
+S3method("[[<-",bigq)
+S3method("[[",bigq)
+##S3method(add,bigq)
+S3method(c,bigq)
+S3method(cbind,bigq)
+S3method(rbind,bigq)
+S3method(dim,bigq)
+S3method("dim<-", bigq)
+S3method(is.finite,bigq)
+S3method(is.infinite,bigq)
+S3method(is.na,bigq)
+S3method(is.whole,bigq)
+S3method(length,bigq)
+S3method("length<-", bigq)
+S3method(max,bigq)
+S3method(min,bigq)
+S3method(rep,bigq)
+S3method(print,bigq)
+S3method(duplicated,bigq)
+S3method(unique,bigq)
+S3method(mean,bigq)
+S3method(all.equal,bigq)
+
+S3method(matrix, default)
+S3method(matrix, bigq)
+S3method(matrix, bigz)
+
+## unfortunately, using S4 methods and setOldClass("bigz") is not enough
+## and so we need S3 methods and mask base (see also R/matrix-prods.R):
+S3method("%*%",default)
+S3method("%*%",bigq)
+S3method("%*%",bigz)
+S3method("crossprod",default)
+S3method("crossprod",bigq)
+S3method("crossprod",bigz)
+S3method("tcrossprod",default)
+S3method("tcrossprod",bigq)
+S3method("tcrossprod",bigz)
+
+S3method(apply,default)
+S3method(apply,bigq)
+S3method(apply,bigz)
+
+#S3method(as,bigq)
+#S3method(as,bigz)
+S3method(as.character,bigq)
+S3method(as.character,bigz)
+S3method(as.double,bigq)
+S3method(as.double,bigz)
+S3method(as.integer,bigq)
+S3method(as.integer,bigz)
+S3method(as.matrix,bigq)
+S3method(as.matrix,bigz)
+S3method(as.vector,bigz)
+S3method(as.vector,bigq)
+
+S3method(formatN,default)
+S3method(formatN,integer)
+S3method(formatN,double)
+S3method(formatN,bigq)
+S3method(formatN,bigz)
+
+##----'Math'--------most give "not yet implemented":
+S3method(Math,bigq)
+S3method(Math,bigz)
+
+S3method(abs,bigz)
+S3method(abs,bigq)
+S3method(sign,bigq)
+S3method(sign,bigz)
+
+S3method(floor,bigq)
+S3method(floor,bigz)
+S3method(round,bigq)
+S3method(round,bigz)
+S3method(trunc,bigq)
+S3method(trunc,bigz)
+
+S3method(gamma,bigz)
+
+S3method(cumsum,bigz)
+S3method(cumsum,bigq)
+
+##-- end{'Math'}-------------------
+
+S3method(solve,bigq)
+S3method(solve,bigz)
+
+S3method(t,bigq)
+S3method(t,bigz)
+
+## S4 generic .. also imported in Rmpfr
+export(asNumeric)
+exportMethods(asNumeric)
+
+exportMethods(which.max, which.min)
+
+exportClasses("bigz", "bigq")# setOldClass()'ed, needed in  Rmpfr
+## unfortunately, using S4 methods and setOldClass("bigz") is not enough,
+## see R/matrix-prods.R :
+## exportMethods("%*%", "crossprod", "tcrossprod")#-> also exports S4 generics

R/AllClasses-etc.R

+### Really want  S4 methods for  all the binary operations.
+### Otherwise we "never" make
+###       <bigz> o  <bigq>
+### or    <bigq> o  <bigz>    working --
+##
+## But unfortunately the above seems "impossible", see
+##  see also  setMethod() in ./matrix-prods.R
+
+## OTOH: This *still* helps to define single-dispatch methods for  asNumeric() :
+##       {why does it work there ??}
+
+setOldClass("bigz")#, prototype=as.bigz(integer()))
+setOldClass("bigq")#, prototype=as.bigq(integer()))
+##                cannot use as.bigz() yet which is only defined in ./bigz.R

R/Stirling-n-etc.R

+###--------- Stirling numbers of 1st and 2nd kind ----
+###          ====================================
+### The "double prec" version of this is currently in package 'nacopula'
+### (MM: >>> ../../nacopula/R/special-func.R )
+
+
+##' Compute Stirling numbers of the 1st kind
+##'
+##' s(n,k) = (-1)^{n-k} times
+##' the number of permutations of 1,2,...,n with exactly k cycles
+##'
+##' NIST DLMF 26.8 --> http://dlmf.nist.gov/26.8
+##'
+##' @title  Stirling Numbers of the 1st Kind
+##' @param n
+##' @param k
+##' @return s(n,k)
+##' @author Martin Maechler
+Stirling1 <- function(n,k)
+{
+    ## NOTA BENE: There's no "direct" method available here
+    stopifnot(length(n) == 1, length(k) == 1)
+    if (k < 0 || n < k) stop("'k' must be in 0..n !")
+    if(n == 0) return(as.bigz(1))
+    if(k == 0) return(as.bigz(0))
+    S1 <- function(n,k) {
+        if(k == 0 || n < k) return(as.bigz(0))
+        if(is.na(S <- St[[n]][k])) {
+            ## s(n,k) = s(n-1,k-1) - (n-1) * s(n-1,k) for all n, k >= 0
+            St[[n]][[k]] <<- S <- if(n1 <- n-1L)
+                S1(n1, k-1) - n1* S1(n1, k) else as.bigz(1)
+        }
+        S
+    }
+    if(compute <- (nt <- length(St <- get("S1.tab", .Stirl..env))) < n) {
+        ## extend the "table":
+        length(St) <- n
+        for(i in (nt+1L):n) St[[i]] <- rep(as.bigz(NA), i)
+    }
+    else compute <- is.na(S <- St[[n]][k])
+    if(compute) {
+        S <- S1(n,k)
+        ## store it back:
+        assign("S1.tab", St, envir = .Stirl..env)
+    }
+    S
+}
+
+##' Full Vector of Stirling Numbers of the 1st Kind
+##'
+##' @title  Stirling1(n,k) for all k = 1..n
+##' @param n
+##' @return the same as sapply(1:n, Stirling1, n=n)
+##' @author Martin Maechler
+Stirling1.all <- function(n)
+{
+    stopifnot(length(n) == 1)
+    if(!n) return(as.bigz(numeric(0)))
+    if(get("S1.full.n", .Stirl..env) < n) {
+        assign("S1.full.n", n, envir = .Stirl..env)
+        do.call(c, lapply(seq_len(n), Stirling1, n=n))# which fills "S1.tab"
+    }
+    else get("S1.tab", .Stirl..env)[[n]]
+}
+
+##' Compute Stirling numbers of the 2nd kind
+##'
+##' S^{(k)}_n = number of ways of partitioning a set of $n$ elements into $k$
+##'	non-empty subsets
+##' (Abramowitz/Stegun: 24,1,4 (p. 824-5 ; Table 24.4, p.835)
+##'   Closed Form : p.824 "C."
+##'
+##' @title  Stirling Numbers of the 2nd Kind
+##' @param n
+##' @param k
+##' @param method
+##' @return S(n,k) = S^{(k)}_n
+##' @author Martin Maechler, "direct": May 28 1992
+Stirling2 <- function(n,k, method = c("lookup.or.store","direct"))
+{
+    stopifnot(length(n) == 1, length(k) == 1)
+    if (k < 0 || n < k) stop("'k' must be in 0..n !")
+    method <- match.arg(method)
+    switch(method,
+           "direct" = {
+               sig <- rep(c(1,-1)*(-1)^k, length.out= k+1) # 1 for k=0; -1 1 (k=1)
+               k <- 0:k # (!)
+               ga <- factorialZ(k) # = gamma(k+1)
+               sum( sig * k^n /(ga * rev(ga)))
+           },
+           "lookup.or.store" = {
+               if(n == 0) return(as.bigz(1)) ## else:
+               if(k == 0) return(as.bigz(0))
+               S2 <- function(n,k) {
+                   if(k == 0 || n < k) return(as.bigz(0))