===================================================================
RCS file: /home/cvs/OpenXM_contrib/gmp/mpn/cray/Attic/sub_n.c,v
retrieving revision 1.1
retrieving revision 1.1.1.2
diff -u -p -r1.1 -r1.1.1.2
--- OpenXM_contrib/gmp/mpn/cray/Attic/sub_n.c	2000/09/09 14:12:23	1.1
+++ OpenXM_contrib/gmp/mpn/cray/Attic/sub_n.c	2003/08/25 16:06:19	1.1.1.2
@@ -1,97 +1,82 @@
-/* mpn_sub_n -- Subtract two limb vectors of equal, non-zero length.
-   For Cray vector processors.
+/* Cray PVP mpn_sub_n -- subtract two limb vectors and store their difference
+   in a third limb vector.
 
-   Copyright (C) 1996, 2000 Free Software Foundation, Inc.
+Copyright 1996, 2000, 2001 Free Software Foundation, Inc.
 
-   This file is part of the GNU MP Library.
+This file is part of the GNU MP Library.
 
-   The GNU MP Library is free software; you can redistribute it and/or modify
-   it under the terms of the GNU Lesser General Public License as published by
-   the Free Software Foundation; either version 2.1 of the License, or (at your
-   option) any later version.
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation; either version 2.1 of the License, or (at your
+option) any later version.
 
-   The GNU MP Library is distributed in the hope that it will be useful, but
-   WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
-   or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
-   License for more details.
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+License for more details.
 
-   You should have received a copy of the GNU Lesser General Public License
-   along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
-   the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
-   MA 02111-1307, USA.  */
+You should have received a copy of the GNU Lesser General Public License
+along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA.  */
 
+/* This code runs at 4 cycles/limb.  It may be possible to bring it down
+   to 3 cycles/limb.  */
+
 #include "gmp.h"
 #include "gmp-impl.h"
 
 mp_limb_t
-mpn_sub_n (c, a, b, n)
-     mp_ptr c;
-     mp_srcptr a, b;
-     mp_size_t n;
+mpn_sub_n (mp_ptr rp, mp_srcptr up, mp_srcptr vp, mp_size_t n)
 {
+  mp_limb_t cy[n];
+  mp_limb_t a, b, r, s0, c0, c1;
   mp_size_t i;
-  mp_size_t nm1 = n - 1;
-  int more_carries = 0;
-  int carry_out;
+  int more_carries;
 
-  /* For small operands the non-vector code is faster.  */
-  if (n < 16)
-    goto sequential;
-
-  if (a == c || b == c)
+  /* Main subtract loop.  Generate a raw output difference in rp[] and a
+     borrow vector in cy[].  */
+#pragma _CRI ivdep
+  for (i = 0; i < n; i++)
     {
-      TMP_DECL (marker);
-      TMP_MARK (marker);
-      if (c == a)
-	{
-	  /* allocate temp space for a */
-	  mp_ptr ax = (mp_ptr) TMP_ALLOC (n * BYTES_PER_MP_LIMB);
-	  MPN_COPY (ax, a, n);
-	  a = (mp_srcptr) ax;
-	}
-      if (c == b)
-	{
-	  /* allocate temp space for b */
-	  mp_ptr bx = (mp_ptr) TMP_ALLOC (n * BYTES_PER_MP_LIMB);
-	  MPN_COPY (bx, b, n);
-	  b = (mp_srcptr) bx;
-	}
-      carry_out = mpn_sub_n (c, a, b, n);
-      TMP_FREE (marker);
-      return carry_out;
+      a = up[i];
+      b = vp[i];
+      s0 = a - b;		/* a = s0 + b */
+      rp[i] = s0;
+      c0 = ((s0 & b) | ((s0 | b) & ~a)) >> 63;
+      cy[i] = c0;
     }
-
-  carry_out = a[nm1] < b[nm1];
-
-#pragma _CRI ivdep			/* Cray PVP systems */
-  for (i = nm1; i > 0; i--)
+  /* Borrow subtract loop.  Subtract the borrow vector cy[] from the raw
+     difference rp[] and store the new difference back to rp[0].  If this
+     generates further borrow, set more_carries.  */
+  more_carries = 0;
+#pragma _CRI ivdep
+  for (i = 1; i < n; i++)
     {
-      int cy_in; mp_limb_t t;
-      cy_in = a[i - 1] < b[i - 1];
-      t = a[i] - b[i];
-      more_carries += t < cy_in;
-      c[i] = t - cy_in;
+      r = rp[i];
+      c0 = cy[i - 1];
+      s0 = r - c0;		/* r = s0 + c0 */
+      rp[i] = s0;
+      c0 = (s0 & ~r) >> 63;
+      more_carries += c0;
     }
-  c[0] = a[0] - b[0];
-
+  /* If that second loop generated borrow, handle that in scalar loop.  */
   if (more_carries)
     {
-      /* This won't vectorize, but we should come here rarely.  */
-      int cy;
-    sequential:
-      cy = 0;
-      for (i = 0; i < n; i++)
+      mp_limb_t cyrec = 0;
+      /* Look for places where rp[k] contains just ones and cy[k-1] is
+	 non-zero.  These are where we got a recurrency borrow.  */
+      for (i = 1; i < n; i++)
 	{
-	  mp_limb_t ai, ci, t;
-	  ai = a[i];
-	  t = b[i] + cy;
-	  cy = t < cy;
-	  ci = ai - t;
-	  cy += ci > ai;
-	  c[i] = ci;
+	  r = rp[i];
+	  c0 = (~r == 0 && cy[i - 1] != 0);
+	  s0 = r - cyrec;
+	  rp[i] = s0;
+	  c1 = (s0 & ~r) >> 63;
+	  cyrec = c0 | c1;
 	}
-      carry_out = cy;
+      return cyrec | cy[n - 1];
     }
 
-  return carry_out;
+  return cy[n - 1];
 }