/*
 * @(#) lmuldiv.c - Long (64-bit result) mul/div operations library (source).
 * (c) 1998 Ivan Maidanski <ivmai@chat.ru> http://ivmai.chat.ru
 * Freeware function library source. All rights reserved.
 **
 * Language: ANSI C
 * Tested with: Microsoft Visual C++ (R) v4.2, Borland C/C++ v3.1
 * Last modified: 1998-07-27 11:10:00 GMT+04:00
 */

#include "lmuldiv.h" /* library header */

/*
 Note:
  declare FORCE_IX86 macro equals non-zero to force to use inline asm anyway,
  declare FORCE_IX86 macro equals 0 not to use inline asm even if possible.
*/

#include <stddef.h> /* NULL */

#if (defined(FORCE_IX86) && FORCE_IX86) || (!defined(FORCE_IX86) && \
     !defined(__STDC__) && defined(_MSC_VER) && \
     defined(_M_IX86) && (_M_IX86>=300))
/*only if forced to use asm or if MSC for platforms on intel80386 or higher*/

unsigned long lu_mul_shr(unsigned long a, unsigned long b)
{
 _asm mov eax,a;
 _asm mul b;
 _asm mov a,edx;
 return a;
}

long l_mul_shr(long a, long b)
{
 _asm mov eax,a;
 _asm imul b;
 _asm mov a,edx;
 return a;
}

unsigned long lu_div_mod(unsigned long *prem,
                         unsigned long hv, unsigned long lv, unsigned long d)
{
 _asm mov eax,hv;
 _asm xor edx,edx;
 _asm div d;
 _asm mov eax,lv;
 _asm div d;
 _asm mov lv,eax;
 _asm mov d,edx;
 if (prem!=NULL)
  *prem=d;
 return lv;
}

long l_div_mod(long *prem, long hv, long lv, long d)
{
 _asm mov eax,hv;
 _asm mov ebx,lv;
 _asm mov ecx,d;
 _asm xor edx,edx;
 _asm or eax,eax;
 _asm js minus;
 _asm or ecx,ecx;
 _asm js plusminus;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm jmp finish;
plusminus:
 _asm neg ecx;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm jmp negfin;
minus:
 _asm neg eax;
 _asm neg ebx;
 _asm sbb eax,0;
 _asm or ecx,ecx;
 _asm jns minusplus;
 _asm neg ecx;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm neg edx;
 _asm jmp finish;
minusplus:
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm neg edx;
negfin:
 _asm neg eax;
finish:
 _asm mov lv,eax;
 _asm mov d,edx;
 if (prem!=NULL)
  *prem=d;
 return lv;
}

unsigned long lu_mul_div(unsigned long *prem,
                         unsigned long a, unsigned long b, unsigned long d)
{
 _asm mov eax,a;
 _asm xor ecx,ecx;
 _asm mul b;
 _asm xchg eax,ecx;
 _asm xchg eax,edx;
 _asm div d;
 _asm xchg eax,ecx;
 _asm div d;
 _asm mov a,eax;
 _asm mov d,edx;
 if (prem!=NULL)
  *prem=d;
 return a;
}

long l_mul_div(long *prem, long a, long b, long d)
{
 _asm mov eax,a;
 _asm imul b;
 _asm xchg eax,ebx;
 _asm mov ecx,d;
 _asm xchg eax,edx;
 _asm xor edx,edx;
 _asm or eax,eax;
 _asm js minus;
 _asm or ecx,ecx;
 _asm js plusminus;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm jmp finish;
plusminus:
 _asm neg ecx;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm jmp negfin;
minus:
 _asm neg eax;
 _asm neg ebx;
 _asm sbb eax,0;
 _asm or ecx,ecx;
 _asm jns minusplus;
 _asm neg ecx;
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm neg edx;
 _asm jmp finish;
minusplus:
 _asm div ecx;
 _asm xchg eax,ebx;
 _asm div ecx;
 _asm neg edx;
negfin:
 _asm neg eax;
finish:
 _asm mov a,eax;
 _asm mov d,edx;
 if (prem!=NULL)
  *prem=d;
 return a;
}

#else /* for other platforms */

#include <limits.h> /* CHAR_BIT */

#define HALFLONGVAL ((unsigned long)1<<(sizeof(long)*CHAR_BIT/2))

unsigned long lu_mul_shr(unsigned long a, unsigned long b)
{
 unsigned long c=a/HALFLONGVAL;
 unsigned long d=b/HALFLONGVAL;
 a%=HALFLONGVAL;
 b%=HALFLONGVAL;
 a=(a*b)/HALFLONGVAL+a*d;
 a+=b*=c;
 return c*d+a/HALFLONGVAL+(a<b? HALFLONGVAL : 0);
}

long l_mul_shr(long a, long b)
{
 return a>=0? (b>=0? lu_mul_shr(a,b) : -(long)lu_mul_shr(a,-b)-!!(a*b)) :
              (b>=0? -(long)lu_mul_shr(-a,b)-!!(a*b) : lu_mul_shr(-a,-b));
}

unsigned long lu_div_mod(unsigned long *prem,
                         unsigned long hv, unsigned long lv, unsigned long d)
{
 int p;
 unsigned long r,m;
 hv%=d;
 if (!hv)
 {
  if (prem!=NULL)
   *prem=lv%d;
  return lv/d;
 }
 if (d<=HALFLONGVAL)
 {
  hv=hv*HALFLONGVAL+lv/HALFLONGVAL;
  lv=(hv%d)*HALFLONGVAL+lv%HALFLONGVAL;
  if (prem!=NULL)
   *prem=lv%d;
  return (hv/d)*HALFLONGVAL+lv/d;
 }
 if (!(d%HALFLONGVAL))
 {
  d/=HALFLONGVAL;
  m=(hv%d)*HALFLONGVAL+lv/HALFLONGVAL;
  if (prem!=NULL)
   *prem=(m%d)*HALFLONGVAL+lv%HALFLONGVAL;
  return (hv/d)*HALFLONGVAL+m/d;
 }
 for (r=0,m=1;m;m*=2)
 {
  r*=2;
  p=(long)hv<0;
  hv*=2;
  if ((long)lv<0)
   hv++;
  lv*=2;
  if (p || hv>=d)
  {
   hv-=d;
   r++;
  }
 }
 if (prem!=NULL)
  *prem=hv;
 return r;
}

long l_div_mod(long *prem, long hv, long lv, long d)
{
 if (hv<0)
 {
  hv=-hv-!!lv;
  lv=-lv;
  d=d>=0? -(long)lu_div_mod((unsigned long *)prem,hv,lv,d) :
          lu_div_mod((unsigned long *)prem,hv,lv,-d);
  if (prem!=NULL)
   *prem=-*prem;
 }
  else d=d>=0? lu_div_mod((unsigned long *)prem,hv,lv,d) :
               -(long)lu_div_mod((unsigned long *)prem,hv,lv,-d);
  return d;
}

unsigned long lu_mul_div(unsigned long *prem,
                         unsigned long a, unsigned long b, unsigned long d)
{
 return lu_div_mod(prem,lu_mul_shr(a,b),a*b,d);
}

long l_mul_div(long *prem, long a, long b, long d)
{
 return l_div_mod(prem,l_mul_shr(a,b),a*b,d);
}

#endif