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

数学计算

开发平台:

Unix_Linux

  1. dnl  AMD K6 mpn_addmul_1/mpn_submul_1 -- add or subtract mpn multiple.
  2. dnl  Copyright 1999, 2000, 2001, 2002, 2003, 2005 Free Software Foundation,
  3. dnl  Inc.
  4. dnl
  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
  8. dnl  modify it under the terms of the GNU Lesser General Public License as
  9. dnl  published by the Free Software Foundation; either version 3 of the
  10. dnl  License, or (at your option) any later version.
  11. dnl
  12. dnl  The GNU MP Library is distributed in the hope that it will be useful,
  13. dnl  but WITHOUT ANY WARRANTY; without even the implied warranty of
  14. dnl  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  15. dnl  Lesser General Public License for more details.
  16. dnl
  17. dnl  You should have received a copy of the GNU Lesser General Public License
  18. dnl  along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.
  19. include(`../config.m4')
  20. C                           cycles/limb
  21. C P5:
  22. C P6 model 0-8,10-12)            5.94
  23. C P6 model 9  (Banias)
  24. C P6 model 13 (Dothan)           5.57
  25. C P4 model 0  (Willamette)
  26. C P4 model 1  (?)
  27. C P4 model 2  (Northwood)
  28. C P4 model 3  (Prescott)
  29. C P4 model 4  (Nocona)
  30. C K6:                           7.65-8.5 (data dependent)
  31. C K7:
  32. C K8:
  33. dnl  K6:           large multipliers  small multipliers
  34. dnl  UNROLL_COUNT    cycles/limb       cycles/limb
  35. dnl        4             9.5              7.78
  36. dnl        8             9.0              7.78
  37. dnl       16             8.4              7.65
  38. dnl       32             8.4              8.2
  39. dnl
  40. dnl  Maximum possible unrolling with the current code is 32.
  41. dnl
  42. dnl  Unrolling to 16 limbs/loop makes the unrolled loop fit exactly in a 256
  43. dnl  byte block, which might explain the good speed at that unrolling.
  44. deflit(UNROLL_COUNT, 16)
  45. ifdef(`OPERATION_addmul_1', `
  46. define(M4_inst,        addl)
  47. define(M4_function_1,  mpn_addmul_1)
  48. define(M4_function_1c, mpn_addmul_1c)
  49. ',`ifdef(`OPERATION_submul_1', `
  50. define(M4_inst,        subl)
  51. define(M4_function_1,  mpn_submul_1)
  52. define(M4_function_1c, mpn_submul_1c)
  53. ',`m4_error(`Need OPERATION_addmul_1 or OPERATION_submul_1
  54. ')')')
  55. MULFUNC_PROLOGUE(mpn_addmul_1 mpn_addmul_1c mpn_submul_1 mpn_submul_1c)
  56. C mp_limb_t mpn_addmul_1 (mp_ptr dst, mp_srcptr src, mp_size_t size,
  57. C                         mp_limb_t mult);
  58. C mp_limb_t mpn_addmul_1c (mp_ptr dst, mp_srcptr src, mp_size_t size,
  59. C                          mp_limb_t mult, mp_limb_t carry);
  60. C mp_limb_t mpn_submul_1 (mp_ptr dst, mp_srcptr src, mp_size_t size,
  61. C                         mp_limb_t mult);
  62. C mp_limb_t mpn_submul_1c (mp_ptr dst, mp_srcptr src, mp_size_t size,
  63. C                          mp_limb_t mult, mp_limb_t carry);
  64. C
  65. C The jadcl0()s in the unrolled loop makes the speed data dependent.  Small
  66. C multipliers (most significant few bits clear) result in few carry bits and
  67. C speeds up to 7.65 cycles/limb are attained.  Large multipliers (most
  68. C significant few bits set) make the carry bits 50/50 and lead to something
  69. C more like 8.4 c/l.  With adcl's both of these would be 9.3 c/l.
  70. C
  71. C It's important that the gains for jadcl0 on small multipliers don't come
  72. C at the cost of slowing down other data.  Tests on uniformly distributed
  73. C random data, designed to confound branch prediction, show about a 7%
  74. C speed-up using jadcl0 over adcl (8.93 versus 9.57 cycles/limb, with all
  75. C overheads included).
  76. C
  77. C In the simple loop, jadcl0() measures slower than adcl (11.9-14.7 versus
  78. C 11.0 cycles/limb), and hence isn't used.
  79. C
  80. C In the simple loop, note that running ecx from negative to zero and using
  81. C it as an index in the two movs wouldn't help.  It would save one
  82. C instruction (2*addl+loop becoming incl+jnz), but there's nothing unpaired
  83. C that would be collapsed by this.
  84. C
  85. C Attempts at a simpler main loop, with less unrolling, haven't yielded much
  86. C success, generally running over 9 c/l.
  87. C
  88. C
  89. C jadcl0
  90. C ------
  91. C
  92. C jadcl0() being faster than adcl $0 seems to be an artifact of two things,
  93. C firstly the instruction decoding and secondly the fact that there's a
  94. C carry bit for the jadcl0 only on average about 1/4 of the time.
  95. C
  96. C The code in the unrolled loop decodes something like the following.
  97. C
  98. C                                         decode cycles
  99. C mull %ebp                    2
  100. C M4_inst %esi, disp(%edi)        1
  101. C adcl %eax, %ecx              2
  102. C movl %edx, %esi             1
  103. C jnc 1f                    /
  104. C incl %esi                   1
  105. C 1: movl disp(%ebx), %eax      /
  106. C                                              ---
  107. C                                               7
  108. C
  109. C In a back-to-back style test this measures 7 with the jnc not taken, or 8
  110. C with it taken (both when correctly predicted).  This is opposite to the
  111. C measurements showing small multipliers running faster than large ones.
  112. C Don't really know why.
  113. C
  114. C It's not clear how much branch misprediction might be costing.  The K6
  115. C doco says it will be 1 to 4 cycles, but presumably it's near the low end
  116. C of that range to get the measured results.
  117. C
  118. C
  119. C In the code the two carries are more or less the preceding mul product and
  120. C the calculation is roughly
  121. C
  122. C x*y + u*b+v
  123. C
  124. C where b=2^32 is the size of a limb, x*y is the two carry limbs, and u and
  125. C v are the two limbs it's added to (being the low of the next mul, and a
  126. C limb from the destination).
  127. C
  128. C To get a carry requires x*y+u*b+v >= b^2, which is u*b+v >= b^2-x*y, and
  129. C there are b^2-(b^2-x*y) = x*y many such values, giving a probability of
  130. C x*y/b^2.  If x, y, u and v are random and uniformly distributed between 0
  131. C and b-1, then the total probability can be summed over x and y,
  132. C
  133. C  1    b-1 b-1 x*y    1    b*(b-1)   b*(b-1)
  134. C --- * sum sum --- = --- * ------- * ------- = 1/4
  135. C       b^2   x=0 y=1 b^2   b^4      2         2
  136. C
  137. C Actually it's a very tiny bit less than 1/4 of course.  If y is fixed,
  138. C then the probability is 1/2*y/b thus varying linearly between 0 and 1/2.
  139. ifdef(`PIC',`
  140. deflit(UNROLL_THRESHOLD, 9)
  141. ',`
  142. deflit(UNROLL_THRESHOLD, 6)
  143. ')
  144. defframe(PARAM_CARRY,     20)
  145. defframe(PARAM_MULTIPLIER,16)
  146. defframe(PARAM_SIZE,      12)
  147. defframe(PARAM_SRC,       8)
  148. defframe(PARAM_DST,       4)
  149. TEXT
  150. ALIGN(32)
  151. PROLOGUE(M4_function_1c)
  152. pushl %esi
  153. deflit(`FRAME',4)
  154. movl PARAM_CARRY, %esi
  155. jmp L(start_nc)
  156. EPILOGUE()
  157. PROLOGUE(M4_function_1)
  158. push %esi
  159. deflit(`FRAME',4)
  160. xorl %esi, %esi C initial carry
  161. L(start_nc):
  162. movl PARAM_SIZE, %ecx
  163. pushl %ebx
  164. deflit(`FRAME',8)
  165. movl PARAM_SRC, %ebx
  166. pushl %edi
  167. deflit(`FRAME',12)
  168. cmpl $UNROLL_THRESHOLD, %ecx
  169. movl PARAM_DST, %edi
  170. pushl %ebp
  171. deflit(`FRAME',16)
  172. jae L(unroll)
  173. C simple loop
  174. movl PARAM_MULTIPLIER, %ebp
  175. L(simple):
  176. C eax scratch
  177. C ebx src
  178. C ecx counter
  179. C edx scratch
  180. C esi carry
  181. C edi dst
  182. C ebp multiplier
  183. movl (%ebx), %eax
  184. addl $4, %ebx
  185. mull %ebp
  186. addl $4, %edi
  187. addl %esi, %eax
  188. adcl $0, %edx
  189. M4_inst %eax, -4(%edi)
  190. adcl $0, %edx
  191. movl %edx, %esi
  192. loop L(simple)
  193. popl %ebp
  194. popl %edi
  195. popl %ebx
  196. movl %esi, %eax
  197. popl %esi
  198. ret
  199. C -----------------------------------------------------------------------------
  200. C The unrolled loop uses a "two carry limbs" scheme.  At the top of the loop
  201. C the carries are ecx=lo, esi=hi, then they swap for each limb processed.
  202. C For the computed jump an odd size means they start one way around, an even
  203. C size the other.
  204. C
  205. C VAR_JUMP holds the computed jump temporarily because there's not enough
  206. C registers at the point of doing the mul for the initial two carry limbs.
  207. C
  208. C The add/adc for the initial carry in %esi is necessary only for the
  209. C mpn_addmul/submul_1c entry points.  Duplicating the startup code to
  210. C eliminate this for the plain mpn_add/submul_1 doesn't seem like a good
  211. C idea.
  212. dnl  overlapping with parameters already fetched
  213. define(VAR_COUNTER, `PARAM_SIZE')
  214. define(VAR_JUMP,    `PARAM_DST')
  215. L(unroll):
  216. C eax
  217. C ebx src
  218. C ecx size
  219. C edx
  220. C esi initial carry
  221. C edi dst
  222. C ebp
  223. movl %ecx, %edx
  224. decl %ecx
  225. subl $2, %edx
  226. negl %ecx
  227. shrl $UNROLL_LOG2, %edx
  228. andl $UNROLL_MASK, %ecx
  229. movl %edx, VAR_COUNTER
  230. movl %ecx, %edx
  231. shll $4, %edx
  232. negl %ecx
  233. C 15 code bytes per limb
  234. ifdef(`PIC',`
  235. call L(pic_calc)
  236. L(here):
  237. ',`
  238. leal L(entry) (%edx,%ecx,1), %edx
  239. ')
  240. movl (%ebx), %eax C src low limb
  241. movl PARAM_MULTIPLIER, %ebp
  242. movl %edx, VAR_JUMP
  243. mull %ebp
  244. addl %esi, %eax C initial carry (from _1c)
  245. jadcl0( %edx)
  246. leal 4(%ebx,%ecx,4), %ebx
  247. movl %edx, %esi C high carry
  248. movl VAR_JUMP, %edx
  249. leal (%edi,%ecx,4), %edi
  250. testl $1, %ecx
  251. movl %eax, %ecx C low carry
  252. jz L(noswap)
  253. movl %esi, %ecx C high,low carry other way around
  254. movl %eax, %esi
  255. L(noswap):
  256. jmp *%edx
  257. ifdef(`PIC',`
  258. L(pic_calc):
  259. C See mpn/x86/README about old gas bugs
  260. leal (%edx,%ecx,1), %edx
  261. addl $L(entry)-L(here), %edx
  262. addl (%esp), %edx
  263. ret_internal
  264. ')
  265. C -----------------------------------------------------------
  266. ALIGN(32)
  267. L(top):
  268. deflit(`FRAME',16)
  269. C eax scratch
  270. C ebx src
  271. C ecx carry lo
  272. C edx scratch
  273. C esi carry hi
  274. C edi dst
  275. C ebp multiplier
  276. C
  277. C 15 code bytes per limb
  278. leal UNROLL_BYTES(%edi), %edi
  279. L(entry):
  280. forloop(`i', 0, UNROLL_COUNT/2-1, `
  281. deflit(`disp0', eval(2*i*4))
  282. deflit(`disp1', eval(disp0 + 4))
  283. Zdisp( movl, disp0,(%ebx), %eax)
  284. mull %ebp
  285. Zdisp( M4_inst,%ecx, disp0,(%edi))
  286. adcl %eax, %esi
  287. movl %edx, %ecx
  288. jadcl0( %ecx)
  289. movl disp1(%ebx), %eax
  290. mull %ebp
  291. M4_inst %esi, disp1(%edi)
  292. adcl %eax, %ecx
  293. movl %edx, %esi
  294. jadcl0( %esi)
  295. ')
  296. decl VAR_COUNTER
  297. leal UNROLL_BYTES(%ebx), %ebx
  298. jns L(top)
  299. popl %ebp
  300. M4_inst %ecx, UNROLL_BYTES(%edi)
  301. popl %edi
  302. movl %esi, %eax
  303. popl %ebx
  304. jadcl0( %eax)
  305. popl %esi
  306. ret
  307. EPILOGUE()