Sunday, June 09, 2013

An iterative method to compute the integer linear combination of Greatest Common Divisor

Euclid's iterative method to compute the GCD is very well known. You might have also heard that $GCD(a,b)$ can be expressed as integer linear combination of $a,b$ (i.e. $GCD(a,b) = \alpha\times a + \beta \times b \,\, \alpha,\beta \in \mathcal{I}$). The following is an iterative method -- I wrote on a flight between AUS -> DFW just for fun -- which computes this integer linear combination. Runs in $O(\log(\max\{a,b\}))$ time.


template
 static inline Unit gcd_integer_linear_combination(const Unit& a, const Unit& b, Unit& alpha, Unit& beta ){
  Unit divd = a, divs = b, remd=0, quot;
  Unit a_divd[2] = {1,0}; //linear combination of the divident//
  Unit a_divs[2] = {0,1}, temp[2] = {0,0}; //linear combination of the divisor//
  bool b_flip = false;
  if(divd < divs){
   divd = b; divs = a;
   b_flip = true;
  }
  //Invariant's: divd = a_divd[0]*a + a_divd[1]*b //
  //      divs = a_divs[0]*a + a_divs[1]*b //
  while(divs){
   remd = divd%divs; quot = (divd-remd)/divs;
   divd = divs; divs = remd;
   temp[0] = a_divs[0]; temp[1] = a_divs[1];
   a_divs[0] = a_divd[0] - (quot * a_divs[0]); 
   a_divs[1] = a_divd[1] - (quot * a_divs[1]);
   a_divd[0] = temp[0]; a_divd[1] = temp[1];
  }
  if(!b_flip){
   alpha = a_divd[0]; beta = a_divd[1];
  }else{
   alpha = a_divd[1]; beta = a_divd[0];
  }
  return divd;
 }

The following are some examples $GCD(a,b) = (\alpha)*a + (\beta)*b$.
PASSED: 2 = (1613)*42 + (-8)*8468 
PASSED: 1 = (-1684)*6335 + (1641)*6501 
PASSED: 5 = (344)*9170 + (-551)*5725 
PASSED: 1 = (-1949)*1479 + (308)*9359 
PASSED: 1 = (967)*6963 + (-1508)*4465 
PASSED: 2 = (858)*5706 + (-601)*8146 
PASSED: 6 = (181)*3282 + (-87)*6828 
PASSED: 2 = (121)*9962 + (-2450)*492 
PASSED: 1 = (-596)*2996 + (919)*1943 
PASSED: 1 = (-2348)*4828 + (2085)*5437 
PASSED: 1 = (283)*2392 + (-147)*4605 
PASSED: 1 = (-61)*3903 + (1546)*154 
PASSED: 1 = (122)*293 + (-15)*2383 
PASSED: 1 = (-1528)*7422 + (1301)*8717 
PASSED: 1 = (615)*9719 + (-604)*9896 
PASSED: 1 = (-304)*5448 + (959)*1727 
PASSED: 1 = (-139)*4772 + (431)*1539 
PASSED: 1 = (493)*1870 + (-93)*9913 
PASSED: 4 = (628)*5668 + (-565)*6300 
PASSED: 1 = (3461)*7036 + (-2461)*9895 
PASSED: 4 = (-60)*8704 + (137)*3812 
PASSED: 1 = (77)*1323 + (-305)*334 
PASSED: 3 = (662)*7674 + (-1089)*4665 
PASSED: 2 = (3)*5142 + (-2)*7712 
PASSED: 1 = (1840)*8254 + (-2211)*6869 
PASSED: 1 = (2027)*5548 + (-1471)*7645 
PASSED: 1 = (929)*2663 + (-897)*2758 
PASSED: 2 = (-301)*38 + (4)*2860 
PASSED: 2 = (-823)*8724 + (737)*9742 
PASSED: 1 = (-3)*7530 + (29)*779 
PASSED: 1 = (-1499)*2317 + (1144)*3036 
PASSED: 1 = (662)*2191 + (-787)*1843 
PASSED: 1 = (10)*289 + (-27)*107 
PASSED: 1 = (4289)*9041 + (-4336)*8943

No comments: