mpn/x86/divrem_1.asm

   1 dnl  x86 mpn_divrem_1 -- mpn by limb division extending to fractional quotient.
   2
   3 dnl  Copyright 1999, 2000, 2001, 2002, 2003, 2007 Free Software Foundation,
   4 dnl  Inc.
   5 dnl
   6 dnl  This file is part of the GNU MP Library.
   7 dnl
   8 dnl  The GNU MP Library is free software; you can redistribute it and/or
   9 dnl  modify it under the terms of the GNU Lesser General Public License as
  10 dnl  published by the Free Software Foundation; either version 3 of the
  11 dnl  License, or (at your option) any later version.
  12 dnl
  13 dnl  The GNU MP Library is distributed in the hope that it will be useful,
  14 dnl  but WITHOUT ANY WARRANTY; without even the implied warranty of
  15 dnl  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  16 dnl  Lesser General Public License for more details.
  17 dnl
  18 dnl  You should have received a copy of the GNU Lesser General Public License
  19 dnl  along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.
  20
  21 include(`../config.m4')
  22
  23
  24 C       cycles/limb
  25 C 486   approx 43 maybe
  26 C P5        44
  27 C P6        39
  28 C P6MMX     39
  29 C K6        22
  30 C K7        42
  31 C P4        58
  32
  33
  34 C mp_limb_t mpn_divrem_1 (mp_ptr dst, mp_size_t xsize,
  35 C                         mp_srcptr src, mp_size_t size, mp_limb_t divisor);
  36 C mp_limb_t mpn_divrem_1c (mp_ptr dst, mp_size_t xsize,
  37 C                          mp_srcptr src, mp_size_t size, mp_limb_t divisor,
  38 C                          mp_limb_t carry);
  39 C
  40 C Divide src,size by divisor and store the quotient in dst+xsize,size.
  41 C Extend the division to fractional quotient limbs in dst,xsize.  Return the
  42 C remainder.  Either or both xsize and size can be 0.
  43 C
  44 C mpn_divrem_1c takes a carry parameter which is an initial high limb,
  45 C effectively one extra limb at the top of src,size.  Must have
  46 C carry<divisor.
  47 C
  48 C
  49 C Essentially the code is the same as the division based part of
  50 C mpn/generic/divrem_1.c, but has the advantage that we get the desired divl
  51 C instruction even when gcc is not being used (when longlong.h only has the
  52 C rather slow generic C udiv_qrnnd().
  53 C
  54 C A test is done to see if the high limb is less than the divisor, and if so
  55 C one less div is done.  A div is between 20 and 40 cycles on the various
  56 C x86s, so assuming high<divisor about half the time, then this test saves
  57 C half that amount.  The branch misprediction penalty on each chip is less
  58 C than half a div.
  59 C
  60 C
  61 C Notes for P5:
  62 C
  63 C It might be thought that moving the load down to pair with the store would
  64 C save 1 cycle, but that doesn't seem to happen in practice, and in any case
  65 C would be a mere 2.2% saving, so it's hardly worth bothering about.
  66 C
  67 C A mul-by-inverse might be a possibility for P5, as done in
  68 C mpn/x86/pentium/mod_1.asm.  The number of auxiliary instructions required
  69 C is a hinderance, but there could be a 10-15% speedup available.
  70 C
  71 C
  72 C Notes for K6:
  73 C
  74 C K6 has its own version of this code, using loop and paying attention to
  75 C cache line boundary crossings.  The target 20 c/l can be had with the
  76 C decl+jnz of the present code by pairing up the load and store in the
  77 C loops.  But it's considered easier not to introduce complexity just for
  78 C that, but instead let k6 have its own code.
  79 C
  80
  81 defframe(PARAM_CARRY,  24)
  82 defframe(PARAM_DIVISOR,20)
  83 defframe(PARAM_SIZE,   16)
  84 defframe(PARAM_SRC,    12)
  85 defframe(PARAM_XSIZE,  8)
  86 defframe(PARAM_DST,    4)
  87
  88         TEXT
  89         ALIGN(16)
  90
  91 PROLOGUE(mpn_divrem_1c)
  92 deflit(`FRAME',0)
  93
  94         movl    PARAM_SIZE, %ecx
  95         pushl   %edi            FRAME_pushl()
  96
  97         movl    PARAM_SRC, %edi
  98         pushl   %esi            FRAME_pushl()
  99
 100         movl    PARAM_DIVISOR, %esi
 101         pushl   %ebx            FRAME_pushl()
 102
 103         movl    PARAM_DST, %ebx
 104         pushl   %ebp            FRAME_pushl()
 105
 106         movl    PARAM_XSIZE, %ebp
 107         orl     %ecx, %ecx
 108
 109         movl    PARAM_CARRY, %edx
 110         jz      L(fraction)
 111
 112         leal    -4(%ebx,%ebp,4), %ebx   C dst one limb below integer part
 113         jmp     L(integer_top)
 114
 115 EPILOGUE()
 116
 117
 118 PROLOGUE(mpn_divrem_1)
 119 deflit(`FRAME',0)
 120
 121         movl    PARAM_SIZE, %ecx
 122         pushl   %edi            FRAME_pushl()
 123
 124         movl    PARAM_SRC, %edi
 125         pushl   %esi            FRAME_pushl()
 126
 127         movl    PARAM_DIVISOR, %esi
 128         orl     %ecx,%ecx
 129
 130         jz      L(size_zero)
 131         pushl   %ebx            FRAME_pushl()
 132
 133         movl    -4(%edi,%ecx,4), %eax   C src high limb
 134         xorl    %edx, %edx
 135
 136         movl    PARAM_DST, %ebx
 137         pushl   %ebp            FRAME_pushl()
 138
 139         movl    PARAM_XSIZE, %ebp
 140         cmpl    %esi, %eax
 141
 142         leal    -4(%ebx,%ebp,4), %ebx   C dst one limb below integer part
 143         jae     L(integer_entry)
 144
 145
 146         C high<divisor, so high of dst is zero, and avoid one div
 147
 148         movl    %edx, (%ebx,%ecx,4)
 149         decl    %ecx
 150
 151         movl    %eax, %edx
 152         jz      L(fraction)
 153
 154
 155 L(integer_top):
 156         C eax   scratch (quotient)
 157         C ebx   dst+4*xsize-4
 158         C ecx   counter
 159         C edx   scratch (remainder)
 160         C esi   divisor
 161         C edi   src
 162         C ebp   xsize
 163
 164         movl    -4(%edi,%ecx,4), %eax
 165 L(integer_entry):
 166
 167         divl    %esi
 168
 169         movl    %eax, (%ebx,%ecx,4)
 170         decl    %ecx
 171         jnz     L(integer_top)
 172
 173
 174 L(fraction):
 175         orl     %ebp, %ecx
 176         jz      L(done)
 177
 178         movl    PARAM_DST, %ebx
 179
 180
 181 L(fraction_top):
 182         C eax   scratch (quotient)
 183         C ebx   dst
 184         C ecx   counter
 185         C edx   scratch (remainder)
 186         C esi   divisor
 187         C edi
 188         C ebp
 189
 190         xorl    %eax, %eax
 191
 192         divl    %esi
 193
 194         movl    %eax, -4(%ebx,%ecx,4)
 195         decl    %ecx
 196         jnz     L(fraction_top)
 197
 198
 199 L(done):
 200         popl    %ebp
 201         movl    %edx, %eax
 202         popl    %ebx
 203         popl    %esi
 204         popl    %edi
 205         ret
 206
 207
 208 L(size_zero):
 209 deflit(`FRAME',8)
 210         movl    PARAM_XSIZE, %ecx
 211         xorl    %eax, %eax
 212
 213         movl    PARAM_DST, %edi
 214
 215         cld     C better safe than sorry, see mpn/x86/README
 216
 217         rep
 218         stosl
 219
 220         popl    %esi
 221         popl    %edi
 222         ret
 223 EPILOGUE()