gcd_1.asm
上传用户:qaz666999
上传日期:2022-08-06
资源大小:2570k
文件大小:4k
源码类别:

数学计算

开发平台:

Unix_Linux

  1. dnl  Alpha ev67 mpn_gcd_1 -- Nx1 greatest common divisor.
  2. dnl  Copyright 2003, 2004 Free Software Foundation, Inc.
  3. dnl  This file is part of the GNU MP Library.
  4. dnl
  5. dnl  The GNU MP Library is free software; you can redistribute it and/or
  6. dnl  modify it under the terms of the GNU Lesser General Public License as
  7. dnl  published by the Free Software Foundation; either version 3 of the
  8. dnl  License, or (at your option) any later version.
  9. dnl
  10. dnl  The GNU MP Library is distributed in the hope that it will be useful,
  11. dnl  but WITHOUT ANY WARRANTY; without even the implied warranty of
  12. dnl  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13. dnl  Lesser General Public License for more details.
  14. dnl
  15. dnl  You should have received a copy of the GNU Lesser General Public License
  16. dnl  along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.
  17. include(`../config.m4')
  18. C ev67: 3.4 cycles/bitpair for 1x1 part
  19. C mp_limb_t mpn_gcd_1 (mp_srcptr xp, mp_size_t xsize, mp_limb_t y);
  20. C
  21. C In the 1x1 part, the algorithm is to change x,y to abs(x-y),min(x,y) and
  22. C strip trailing zeros from abs(x-y) to maintain x and y both odd.
  23. C
  24. C The trailing zeros are calculated from just x-y, since in twos-complement
  25. C there's the same number of trailing zeros on d or -d.  This means the cttz
  26. C runs in parallel with abs(x-y).
  27. C
  28. C The loop takes 5 cycles, and at 0.68 iterations per bit for two N-bit
  29. C operands with this algorithm gives the measured 3.4 c/l.
  30. C
  31. C The slottings shown are for SVR4 style systems, Unicos differs in the
  32. C initial gp setup and the LEA.
  33. C
  34. C Enhancement:
  35. C
  36. C On the jsr, !lituse_jsr! (when available) would allow the linker to relax
  37. C it to a bsr, but probably only in a static binary.  Plain "jsr foo" gives
  38. C the right object code for relaxation, and ought to be available
  39. C everywhere, but we prefer to schedule the GOT ldq (LEA) back earlier, for
  40. C the usual case of running in a shared library.
  41. C
  42. C bsr could perhaps be used explicitly anyway.  We should be able to assume
  43. C modexact is in the same module as us (ie. shared library or mainline).
  44. C Would there be any worries about the size of the displacement?  Could
  45. C always put modexact and gcd_1 in the same .o to be certain.
  46. ASM_START()
  47. PROLOGUE(mpn_gcd_1, gp)
  48. C r16 xp
  49. C r17 size
  50. C r18 y
  51. C ldah C l
  52. C lda C u
  53. ldq r0, 0(r16) C L   x = xp[0]
  54. lda r30, -32(r30) C u   alloc stack
  55. LEA(  r27, mpn_modexact_1c_odd) C L   modexact addr, ldq (gp)
  56. stq r10, 16(r30) C L   save r10
  57. cttz r18, r10 C U0  y twos
  58. cmpeq r17, 1, r5 C u   test size==1
  59. stq r9, 8(r30) C L   save r9
  60. clr r19 C u   zero c for modexact
  61. unop
  62. unop
  63. cttz r0, r6 C U0  x twos
  64. stq r26, 0(r30) C L   save ra
  65. srl r18, r10, r18 C U   y odd
  66. mov r18, r9 C l   hold y across call
  67. cmpult r6, r10, r2 C u   test x_twos < y_twos
  68. cmovne r2, r6, r10 C l   common_twos = min(x_twos,y_twos)
  69. bne r5, L(one) C U   no modexact if size==1
  70. jsr r26, (r27), mpn_modexact_1c_odd   C L0
  71. LDGP( r29, 0(r26)) C u,l ldah,lda
  72. cttz r0, r6 C U0  new x twos
  73. ldq r26, 0(r30) C L   restore ra
  74. L(one):
  75. mov r9, r1 C u   y
  76. ldq r9, 8(r30) C L   restore r9
  77. mov r10, r2 C u   common twos
  78. ldq r10, 16(r30) C L   restore r10
  79. lda r30, 32(r30) C l   free stack
  80. beq r0, L(done) C U   return y if x%y==0
  81. srl r0, r6, r0 C U   x odd
  82. unop
  83. ALIGN(16)
  84. L(top):
  85. C r0 x
  86. C r1 y
  87. C r2 common twos, for use at end
  88. subq r0, r1, r7 C l0  d = x - y
  89. cmpult r0, r1, r16 C u0  test x >= y
  90. subq r1, r0, r4 C l0  new_x = y - x
  91. cttz r7, r8 C U0  d twos
  92. cmoveq r16, r7, r4 C l0  new_x = d if x>=y
  93. cmovne r16, r0, r1 C u0  y = x if x<y
  94. unop C l    force cmoveq into l0
  95. unop C u   /
  96. C C cmoveq2 L0, cmovne2 U0
  97. srl r4, r8, r0 C U0  x = new_x >> twos
  98. bne r7, L(top) C U1  stop when d==0
  99. L(done):
  100. sll r1, r2, r0 C U0  return y << common_twos
  101. ret r31, (r26), 1 C L0
  102. EPILOGUE()
  103. ASM_END()