Home | History | Annotate | Line # | Download | only in ev67 
      1 dnl  Alpha ev67 mpn_gcd_11 -- Nx1 greatest common divisor.
      2 
      3 dnl  Copyright 2003, 2004 Free Software Foundation, Inc.
      4 
      5 dnl  This file is part of the GNU MP Library.
      6 dnl
      7 dnl  The GNU MP Library is free software; you can redistribute it and/or modify
      8 dnl  it under the terms of either:
      9 dnl
     10 dnl    * the GNU Lesser General Public License as published by the Free
     11 dnl      Software Foundation; either version 3 of the License, or (at your
     12 dnl      option) any later version.
     13 dnl
     14 dnl  or
     15 dnl
     16 dnl    * the GNU General Public License as published by the Free Software
     17 dnl      Foundation; either version 2 of the License, or (at your option) any
     18 dnl      later version.
     19 dnl
     20 dnl  or both in parallel, as here.
     21 dnl
     22 dnl  The GNU MP Library is distributed in the hope that it will be useful, but
     23 dnl  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
     24 dnl  or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
     25 dnl  for more details.
     26 dnl
     27 dnl  You should have received copies of the GNU General Public License and the
     28 dnl  GNU Lesser General Public License along with the GNU MP Library.  If not,
     29 dnl  see https://www.gnu.org/licenses/.
     30 
     31 include(`../config.m4')
     32 
     33 
     34 C ev67: 3.4 cycles/bitpair for 1x1 part
     35 
     36 
     37 C mp_limb_t mpn_gcd_1 (mp_srcptr xp, mp_size_t xsize, mp_limb_t y);
     38 C
     39 C In the 1x1 part, the algorithm is to change x,y to abs(x-y),min(x,y) and
     40 C strip trailing zeros from abs(x-y) to maintain x and y both odd.
     41 C
     42 C The trailing zeros are calculated from just x-y, since in twos-complement
     43 C there's the same number of trailing zeros on d or -d.  This means the cttz
     44 C runs in parallel with abs(x-y).
     45 C
     46 C The loop takes 5 cycles, and at 0.68 iterations per bit for two N-bit
     47 C operands with this algorithm gives the measured 3.4 c/l.
     48 C
     49 C The slottings shown are for SVR4 style systems, Unicos differs in the
     50 C initial gp setup and the LEA.
     51 
     52 
     53 ASM_START()
     54 PROLOGUE(mpn_gcd_11)
     55 	mov	r16, r0
     56 	mov	r17, r1
     57 
     58 	ALIGN(16)
     59 L(top):	subq	r0, r1, r7		C l0  d = x - y
     60 	cmpult	r0, r1, r16		C u0  test x >= y
     61 
     62 	subq	r1, r0, r4		C l0  new_x = y - x
     63 	cttz	r7, r8			C U0  d twos
     64 
     65 	cmoveq	r16, r7, r4		C l0  new_x = d if x>=y
     66 	cmovne	r16, r0, r1		C u0  y = x if x<y
     67 	unop				C l   \ force cmoveq into l0
     68 	unop				C u   /
     69 
     70 	C				C cmoveq2 L0, cmovne2 U0
     71 
     72 	srl	r4, r8, r0		C U0  x = new_x >> twos
     73 	bne	r7, L(top)		C U1  stop when d==0
     74 
     75 
     76 L(end):	mov	r1, r0			C U0  return y << common_twos
     77 	ret	r31, (r26), 1		C L0
     78 EPILOGUE()
     79 ASM_END()
     80