campo-sirio/gfm/dmbyn.c

/* void _DivUnsArrByUnsArr(c,a,m,b,n,rf)
 *
 * ARGUMENT
 *      unsigned        a[],b[],c[];
 *      int             rf            rounding flag
 *      int             n,m;          number of digits
 *
 * DESCRIPTION
 *      Divides an m-word integer a by a n-word integer b, storing the
 *  result in c. m must be >=n. The quotient is rounded when rf is set,
 *  otherwise it is truncated.
 *
 *
 * SIDE EFFECTS
 *      None.
 *
 * RETURNS
 *      None.
 *
 * AUTHOR
 *      Brugnoli Giugno 1992
 *
 * MODIFICATIONS
 *
 */

#include <stdio.h>
#include <math.h>
#include "gm.h"
#include "gmsystem.h"

void    _DivUnsArrByUnsArr(c,a,m,b,n,rf)
unsigned  SHORT       a[],c[],b[];
int             n,m,rf;
{
  int i,j,k,nc,ft = 0;
  unsigned long dvd,pqt,pp;
  float nqt;
  unsigned SHORT copa[10],copb[10],af[10],ppt[10],ppd[10],pqa[2];

  for(j=0;j<10;j++)
    af[j]=copa[j]=copb[j]=0;
 
	for(j=0;j<n;j++)
    copb[j]=b[j];
  
	for(j=0;j<m;j++)
    copa[j]=a[j];
  
	for(j=0;j<5;j++)
    c[j]=0;
  
  copa[m]=0;
  nc=m;
  af[0]=1;

  for(i=m-n;i>=0;i--)
  {
    dvd=(long)copa[nc]*(MAXUNSINT+1)+(long)copa[nc-1];
    pqt=dvd/copb[n-1];
    pqa[1]=(unsigned SHORT)(pqt>>16);
    pqa[0]=(unsigned SHORT)(pqt&MAXUNSINT);
    for(j=0;j<10;j++)
      ppd[j]=((j>n-1)?0:copb[j]);
  
		for(j=1;j<=i;j++)
    {
      for(k=9;k>=1;k--)
        ppd[k]=ppd[k-1];
    
			ppd[0]=0;
    }
    _MulUnsArrByUnsArr(ppd,pqa,ppt,10,2,ft);
    if (_CompareUnsArr(copa,ppt,nc+1)==-1)
    {
      pp=(long)ppt[nc]*(MAXUNSINT+1)+(long)ppt[nc-1];
      nqt=(float)dvd;
      nqt/=(float)pp;
      nqt*=(float)pqt;
      nqt+=(float)5;
// Brugnoli was HERE
      pqa[1]=(unsigned SHORT)(((long int)nqt) >> BITSPERUI);
      pqa[0]=(unsigned SHORT)(((long int)nqt) & MAXUNSINT);

      _MulUnsArrByUnsArr(ppd,pqa,ppt,10,2,ft);
    }
    
		while (_CompareUnsArr(copa,ppt,nc+1)==-1)
    {
      _SubUnsArr(pqa,af,2);
      _MulUnsArrByUnsArr(ppd,pqa,ppt,10,2,ft);
    }
    c[i]=pqa[0];
    _SubUnsArr(copa,ppt,m);
    nc--;
  }

  /* qui copa <20> = al resto della divisione!!  */

  if (rf)
  {
    _SubUnsArr(copb,af,n);
    _HalveUnsArr(copb,n);
    if (_CompareUnsArr(copa,copb,m)==1)
      _IncrementUnsArr(c);
  }
}