# HG changeset patch
# User Francesco Potortì  <pot@gnu.org>
# Date 1228738706 -3600
# Node ID a85bde34b8b9d7521b388b0f8b90f4282f11161e
# Parent  af0adfbd3d167ee4e5ec9a1255df1b1bc299c088
Set max_recursion_depth and use a subfunction in nchoosek

diff -r af0adfbd3d16 -r a85bde34b8b9 scripts/ChangeLog
--- a/scripts/ChangeLog	Mon Dec 08 07:47:35 2008 +0100
+++ b/scripts/ChangeLog	Mon Dec 08 13:18:26 2008 +0100
@@ -1,3 +1,7 @@
+2008-11-26  Francesco Potortì  <pot@gnu.org>
+
+	* specfun/nchoosek.m: Set max_recursion_depth and use a subfunction.
+
 2008-11-04  Thorsten Meyer <thorsten.meyier@gmx.de>
 
         * miscellaneous/unpack.m: return directly after recursive handling
diff -r af0adfbd3d16 -r a85bde34b8b9 scripts/specfun/nchoosek.m
--- a/scripts/specfun/nchoosek.m	Mon Dec 08 07:47:35 2008 +0100
+++ b/scripts/specfun/nchoosek.m	Mon Dec 08 13:18:26 2008 +0100
@@ -75,9 +75,26 @@
   elseif (k == n)
      A = v(:).';
   else
-    m = round (exp (sum (log (k:n-1)) - sum (log (2:n-k))));
-    A = [v(1)*ones(m,1), nchoosek(v(2:n),k-1);
-	 nchoosek(v(2:n),k)];
+    oldmax = max_recursion_depth ();
+    unwind_protect
+      max_recursion_depth (n);
+      A = nck (v, k);
+    unwind_protect_cleanup
+      max_recursion_depth (oldmax);
+    end_unwind_protect
   endif
+endfunction
 
+function A = nck (v, k)
+  n = length (v);
+  if (n == 1 || k < 2 || k == n)
+    A = nchoosek (v, k);
+  else
+    m = nchoosek (n-1, k-1);
+    A = [v(1)*ones(m,1), nck(v(2:n),k-1);
+	 nck(v(2:n), k)];
+  endif
 endfunction
+
+%!assert (nchoosek(100,45), bincoeff(100,45))
+%!assert (nchoosek(1:5,3),[1:3;1,2,4;1,2,5;1,3,4;1,3,5;1,4,5;2:4;2,3,5;2,4,5;3:5])