Return to sqr_basecase.c CVS log

Up to [local] / OpenXM_contrib / gmp / mpn / cray / ieee

Annotation of OpenXM_contrib/gmp/mpn/cray/ieee/sqr_basecase.c, Revision 1.1

1.1     ! ohara       1: /* Cray PVP/IEEE mpn_sqr_basecase.
        !             2:
        !             3: Copyright 2000, 2001 Free Software Foundation, Inc.
        !             4:
        !             5: This file is part of the GNU MP Library.
        !             6:
        !             7: The GNU MP Library is free software; you can redistribute it and/or modify
        !             8: it under the terms of the GNU Lesser General Public License as published by
        !             9: the Free Software Foundation; either version 2.1 of the License, or (at your
        !            10: option) any later version.
        !            11:
        !            12: The GNU MP Library is distributed in the hope that it will be useful, but
        !            13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
        !            14: or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
        !            15: License for more details.
        !            16:
        !            17: You should have received a copy of the GNU Lesser General Public License
        !            18: along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
        !            19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
        !            20: MA 02111-1307, USA.  */
        !            21:
        !            22: /* This is just mpn_mul_basecase with trivial modifications.  */
        !            23:
        !            24: #include <intrinsics.h>
        !            25: #include "gmp.h"
        !            26: #include "gmp-impl.h"
        !            27:
        !            28: void
        !            29: mpn_sqr_basecase (mp_ptr rp,
        !            30:                  mp_srcptr up, mp_size_t un)
        !            31: {
        !            32:   mp_limb_t cy[un + un];
        !            33:   mp_limb_t ul;
        !            34:   mp_limb_t a, b, r, s0, s1, c0, c1;
        !            35:   mp_size_t i, j;
        !            36:   int more_carries;
        !            37:
        !            38:   for (i = 0; i < un + un; i++)
        !            39:     {
        !            40:       rp[i] = 0;
        !            41:       cy[i] = 0;
        !            42:     }
        !            43:
        !            44: #pragma _CRI novector
        !            45:   for (j = 0; j < un; j++)
        !            46:     {
        !            47:       ul = up[j];
        !            48:
        !            49:       a = up[0] * ul;
        !            50:       r = rp[j];
        !            51:       s0 = a + r;
        !            52:       rp[j] = s0;
        !            53:       c0 = ((a & r) | ((a | r) & ~s0)) >> 63;
        !            54:       cy[j] += c0;
        !            55:
        !            56: #pragma _CRI ivdep
        !            57:       for (i = 1; i < un; i++)
        !            58:        {
        !            59:          a = up[i] * ul;
        !            60:          b = _int_mult_upper (up[i - 1], ul);
        !            61:          s0 = a + b;
        !            62:          c0 = ((a & b) | ((a | b) & ~s0)) >> 63;
        !            63:          r = rp[j + i];
        !            64:          s1 = s0 + r;
        !            65:          rp[j + i] = s1;
        !            66:          c1 = ((s0 & r) | ((s0 | r) & ~s1)) >> 63;
        !            67:          cy[j + i] += c0 + c1;
        !            68:        }
        !            69:       rp[j + un] = _int_mult_upper (up[un - 1], ul);
        !            70:     }
        !            71:
        !            72:   more_carries = 0;
        !            73: #pragma _CRI ivdep
        !            74:   for (i = 1; i < un + un; i++)
        !            75:     {
        !            76:       r = rp[i];
        !            77:       c0 = cy[i - 1];
        !            78:       s0 = r + c0;
        !            79:       rp[i] = s0;
        !            80:       c0 = (r & ~s0) >> 63;
        !            81:       more_carries += c0;
        !            82:     }
        !            83:   /* If that second loop generated carry, handle that in scalar loop.  */
        !            84:   if (more_carries)
        !            85:     {
        !            86:       mp_limb_t cyrec = 0;
        !            87:       for (i = 1; i < un + un; i++)
        !            88:        {
        !            89:          r = rp[i];
        !            90:          c0 = (r < cy[i - 1]);
        !            91:          s0 = r + cyrec;
        !            92:          rp[i] = s0;
        !            93:          c1 = (r & ~s0) >> 63;
        !            94:          cyrec = c0 | c1;
        !            95:        }
        !            96:     }
        !            97: }