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

数学计算

开发平台:

Unix_Linux

  1. dnl  Intel Pentium mpn_modexact_1_odd -- exact division style remainder.
  2. dnl  Copyright 2000, 2001, 2002 Free Software Foundation, Inc.
  3. dnl
  4. dnl  This file is part of the GNU MP Library.
  5. dnl
  6. dnl  The GNU MP Library is free software; you can redistribute it and/or
  7. dnl  modify it under the terms of the GNU Lesser General Public License as
  8. dnl  published by the Free Software Foundation; either version 3 of the
  9. dnl  License, or (at your option) any later version.
  10. dnl
  11. dnl  The GNU MP Library is distributed in the hope that it will be useful,
  12. dnl  but WITHOUT ANY WARRANTY; without even the implied warranty of
  13. dnl  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14. dnl  Lesser General Public License for more details.
  15. dnl
  16. dnl  You should have received a copy of the GNU Lesser General Public License
  17. dnl  along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.
  18. include(`../config.m4')
  19. C P5: 23.0 cycles/limb
  20. C mp_limb_t mpn_modexact_1_odd (mp_srcptr src, mp_size_t size,
  21. C                               mp_limb_t divisor);
  22. C mp_limb_t mpn_modexact_1c_odd (mp_srcptr src, mp_size_t size,
  23. C                                mp_limb_t divisor, mp_limb_t carry);
  24. C
  25. C There seems no way to pair up the two lone instructions in the main loop.
  26. C
  27. C The special case for size==1 saves about 20 cycles (non-PIC), making it
  28. C the same as mpn_mod_1, and in fact making modexact faster than mod_1 at
  29. C all sizes.
  30. C
  31. C Alternatives:
  32. C
  33. C Using mmx for the multiplies might be possible, with pmullw and pmulhw
  34. C having just 3 cycle latencies, but carry bit handling would probably be
  35. C complicated.
  36. defframe(PARAM_CARRY,  16)
  37. defframe(PARAM_DIVISOR,12)
  38. defframe(PARAM_SIZE,   8)
  39. defframe(PARAM_SRC,    4)
  40. dnl  re-using parameter space
  41. define(VAR_INVERSE,`PARAM_SIZE')
  42. TEXT
  43. ALIGN(16)
  44. PROLOGUE(mpn_modexact_1c_odd)
  45. deflit(`FRAME',0)
  46. movl PARAM_DIVISOR, %eax
  47. movl PARAM_CARRY, %edx
  48. jmp L(start_1c)
  49. EPILOGUE()
  50. ALIGN(16)
  51. PROLOGUE(mpn_modexact_1_odd)
  52. deflit(`FRAME',0)
  53. movl PARAM_DIVISOR, %eax
  54. xorl %edx, %edx C carry
  55. L(start_1c):
  56. ifdef(`PIC',`
  57. call L(here) FRAME_pushl()
  58. L(here):
  59. shrl %eax C d/2
  60. movl (%esp), %ecx C eip
  61. addl $_GLOBAL_OFFSET_TABLE_+[.-L(here)], %ecx
  62. movl %ebx, (%esp) C push ebx
  63. andl $127, %eax
  64. movl PARAM_SIZE, %ebx
  65. movl binvert_limb_table@GOT(%ecx), %ecx
  66. subl $2, %ebx
  67. movb (%eax,%ecx), %cl C inv 8 bits
  68. jc L(one_limb)
  69. ',`
  70. dnl non-PIC
  71. shrl %eax C d/2
  72. pushl %ebx FRAME_pushl()
  73. movl PARAM_SIZE, %ebx
  74. andl $127, %eax
  75. subl $2, %ebx
  76. jc L(one_limb)
  77. movb binvert_limb_table(%eax), %cl C inv 8 bits
  78. ')
  79. movl %ecx, %eax
  80. addl %ecx, %ecx C 2*inv
  81. imull %eax, %eax C inv*inv
  82. imull PARAM_DIVISOR, %eax C inv*inv*d
  83. subl %eax, %ecx C inv = 2*inv - inv*inv*d
  84. movl %ecx, %eax
  85. addl %ecx, %ecx C 2*inv
  86. imull %eax, %eax C inv*inv
  87. imull PARAM_DIVISOR, %eax C inv*inv*d
  88. subl %eax, %ecx C inv = 2*inv - inv*inv*d
  89. pushl %esi FRAME_pushl()
  90. ASSERT(e,` C d*inv == 1 mod 2^GMP_LIMB_BITS
  91. movl %ecx, %eax
  92. imull PARAM_DIVISOR, %eax
  93. cmpl $1, %eax')
  94. movl PARAM_SRC, %esi
  95. movl %ecx, VAR_INVERSE
  96. movl (%esi), %eax C src[0]
  97. leal 4(%esi,%ebx,4), %esi C &src[size-1]
  98. xorl $-1, %ebx C -(size-1)
  99. ASSERT(nz)
  100. jmp L(entry)
  101. C The use of VAR_INVERSE means only a store is needed for that value, rather
  102. C than a push and pop of say %edi.
  103. ALIGN(16)
  104. L(top):
  105. C eax scratch, low product
  106. C ebx counter, limbs, negative
  107. C ecx carry bit
  108. C edx scratch, high product
  109. C esi &src[size-1]
  110. C edi
  111. C ebp
  112. mull PARAM_DIVISOR C h:dummy = q*d
  113. movl (%esi,%ebx,4), %eax C src[i]
  114. subl %ecx, %edx C h -= -c
  115. L(entry):
  116. subl %edx, %eax C s = src[i] - h
  117. sbbl %ecx, %ecx C new -c (0 or -1)
  118. imull VAR_INVERSE, %eax C q = s*i
  119. incl %ebx
  120. jnz L(top)
  121. mull PARAM_DIVISOR
  122. movl (%esi), %eax C src high
  123. subl %ecx, %edx C h -= -c
  124. cmpl PARAM_DIVISOR, %eax
  125. jbe L(skip_last)
  126. deflit(FRAME_LAST,FRAME)
  127. subl %edx, %eax C s = src[i] - h
  128. popl %esi FRAME_popl()
  129. sbbl %ecx, %ecx C c (0 or -1)
  130. popl %ebx FRAME_popl()
  131. imull VAR_INVERSE, %eax C q = s*i
  132. mull PARAM_DIVISOR C h:dummy = q*d
  133. movl %edx, %eax
  134. subl %ecx, %eax
  135. ret
  136. C When high<divisor can skip last step.
  137. L(skip_last):
  138. deflit(`FRAME',FRAME_LAST)
  139. C eax src high
  140. C ebx
  141. C ecx
  142. C edx r
  143. C esi
  144. subl %eax, %edx C r-s
  145. popl %esi FRAME_popl()
  146. sbbl %eax, %eax C -1 if underflow
  147. movl PARAM_DIVISOR, %ebx
  148. andl %ebx, %eax C divisor if underflow
  149. popl %ebx FRAME_popl()
  150. addl %edx, %eax C addback if underflow
  151. ret
  152. C Special case for size==1 using a division for r = c-a mod d.
  153. C Could look for a-c<d and save a division sometimes, but that doesn't seem
  154. C worth bothering about.
  155. L(one_limb):
  156. deflit(`FRAME',4)
  157. C eax
  158. C ebx size-2 (==-1)
  159. C ecx
  160. C edx carry
  161. C esi src end
  162. C edi
  163. C ebp
  164. movl %edx, %eax
  165. movl PARAM_SRC, %edx
  166. movl PARAM_DIVISOR, %ecx
  167. popl %ebx FRAME_popl()
  168. subl (%edx), %eax C c-a
  169. sbbl %edx, %edx
  170. decl %ecx C d-1
  171. andl %ecx, %edx C b*d+c-a if c<a, or c-a if c>=a
  172. divl PARAM_DIVISOR
  173. movl %edx, %eax
  174. ret
  175. EPILOGUE()