simplify2(bitvec VA1,VA2,VX1,VX2[0..n-1]){
  bitvec OVV[0..2n-1][0..2n-1],OV[0..2n-1];
  for each i=0..n-1 do in parallel { 
     if(VA1[i] == 1)
          OVV[i][] = shiftl(VX1,i);
     else OVV[i][]=0;
     if(VA2[i] == 1)
          OVV[n+i][] = shiftl(VX2,i);
     else OVV[n+i][]=0;
  }
  /* parallel sum */
  OV[]=0;
  for i=0..2n-1 do OV += OVV[i];
  return(OV)
}

simplify(bitvec A[0..n-1],X[0..n-1],OVplus[0..2n-1],OVminus[0..2n-1]){
   bitvec VAplus,VAminus,VXplus,VXminus [0..2n-1];
   bitvec V1plus,V1minus,V2plus,V2minus [0..2n-1];

   /* eliminate sequences SQ1s in A and X*/
   VAplus = bitwise_and(bitwise_not(A),shiftl(A,1));  
   VAminus = bitwise_and(A,bitwise_not(shiftl(A,1)));  
   VXplus = bitwise_and(bitwise_not(X),shiftl(X,1));  
   VXminus = bitwise_and(X,bitwise_not(shiftl(X,1)));  

   /* compute the cartesian sum of powers */
   OVplus = simplify2(VAplus,VXplus,VAminus,VXminus);
   OVminus = simplify2(VAplus,VXminus,VAminus,VXplus);
   /* eliminate powers with opposite signs */
   for each i=0..2n-1 do in parallel 
      if(OVplus[i]==1 and OVminus[i]==1)
         OVplus[i] = OVminus[i] = 0;
   
   /* eliminate sequences SQ1s in OVplus and OVminus*/
   V1plus = bitwise_and(bitwise_not(OVplus),shiftl(OVplus,1));  
   V1minus = bitwise_and(bitwise_not(OVminus),shiftl(OVminus,1));  
   V2minus = bitwise_and(OVplus,bitwise_not(shiftl(OVplus,1)));  
   V2plus = bitwise_and(OVminus,bitwise_not(shiftl(OVminus,1)));  

   /* for each sign unify the two vectors */
   OVplus[0] = xor(V1plus[0],V2plus[0]);
   for each i=1..2n-1 do in parallel 
      OVplus[i] = or(
                     xor(V1plus[i],V2plus[i]),
                     and(V1plus[i-1],V2plus[i-1])
		  );
   OVminus[0] = xor(V1minus[0],V2minus[0]);
   for each i=1..2n-1 do in parallel 
      OVminus[i] = or(
                     xor(V1minus[i],V2minus[i]),
                     and(V1minus[i-1],V2minus[i-1])
		  );
   return;
}

mult3(bitvec A[0..n-1],X[0..n-1],Y[0..n-1]){
   bitvec OVplus,OVminus[0..2n-1];
   bitvec OVV[0..n-1][0..2n-1],OV[0..2n-1];   
   
   simplify(A,X,OVplus,OVminus);

   for each i=0..n-1 do in parallel { 
     if(OVplus[i] == 1)
          OVV[i][] = shiftl(Y,i);
     else if(OVminus[i] == 1)
          OVV[i][] = shiftl(-Y,i);  /*-Y==2's complement of Y*/
     else OVV[i][]=0;
  }
  /* parallel sum */
  OV[]=0;
  for i=0..n-1 do OV += OVV[i];
  return(OV)