// ************************************************************************** // // // // eses eses // // eses eses // // eses eseses esesese eses Embedded Systems Group // // ese ese ese ese ese // // ese eseseses eseseses ese Department of Computer Science // // eses eses ese eses // // eses eseses eseseses eses University of Kaiserslautern // // eses eses // // // // ************************************************************************** // // This file contains functions for radix-B arithmetic. A radix-B number is // // thereby represented as a []nat array where B=2^n, e.g., B=2^16 on Abacus // // with 16-bit words. // // ************************************************************************** // bool correct; // ----------------------------------------------------------------------------- // Radix-B Addition y[0..len] = x1[0..len-1] + x2[0..len-1] // ----------------------------------------------------------------------------- // For the correctness of the computation, note that we should compute c' and // y[i] so that c'*B+y[i] = x1[i] + x2[i] + c holds. It is clear by the code // that g*B + s = x1[i] + x2[i] and p*B + y[i] = s + c hold, so that we also // have (g+p)*B + y[i] = x1[i] + x2[i] + c. Thus, we could determine c' simply // as c' = g+p, but we can alternatively also define it as bitwise disjunction // or bitwise xor, since only one of p,g can be 1 as proved below. // We prove that g,p,c∈{0,1} and that g=1 implies p=0: // Since 0 <= x1[i]+x2[i] <= B+(B-2) holds, and B+(B-2) is the machine double // word |0..01|1..101|, it follows that g is one of the machine words |0..0| or // |0..01|, and if g=|0..01|, then s!=|1..1|. If c is also one of |0..0| or // |0..01|, it follows that 0 <= s+c <= B, and the latter is the machine double // word |0..01|0..0|. Hence, p is also either one of the machine words |0..0| or // |0..01|, and the same follows now for c which is obtained as bitwise // disjunction of words that are either |0..0| or |0..01|. // Finally, if g=|0..01|, then s!=|1..1| (as noted above), and therefore when // adding c which is either |0..0| or |0..01|, the sum will not require // additional bits so that p=|0..0| holds in that case. So, g and p cannot be // both |1..1|, and thus all assignments c' = g+p, c' = g|p or c' = !(q<->p) // are equivalent to each other. // ----------------------------------------------------------------------------- procedure NatAddSameLen([]nat x1,x2,y,nat len) { nat c,g,p,s,i; c = 0; // initial carry digit for(i=0..len-1) { // ripple from least to highest digit g,s = x1[i] + x2[i]; // g:generate carry; s:preliminary sum p,y[i] = s + c; // p:propagate carry; add previous carry c = g+p; // define next carry digit 0 or 1 } y[len] = c; return; } // ----------------------------------------------------------------------------- // Radix-B Subtraction y[0..len] = x1[0..len-1] - x2[0..len-1] // ----------------------------------------------------------------------------- // For the correctness of the computation, note that we should compute c' and // y[i] so that -c'*B+y[i] = x1[i] - x2[i] + c holds with c',c∈{0,1}. It is // clear by the code that g*B + s = x1[i] - x2[i] and p*B + y[i] = s - c hold, // so that we also have (g+p)*B + y[i] = x1[i] - x2[i] - c. Thus, we could // determine c' simply as c' = -(g+p), but we have to consider that all this // will be correct in case of overflows as well. // To this end, we will see that g,p are either |0..0| or |1..1| to be read as // two-complement numbers 0 and -1, respectively. To compute x1[i] - x2[i] with // an unsigned subtraction, both numbers x1[i] and x2[i] are extended by |0..0| // to perform the subtraction in full precision. // Since -(B-1) <= x1[i]-x2[i] <= +(B-1) holds, and -(B-1) and +(B-1) are the // machine double words |1..1|0..01| and |0..0|1..1|, respectively, it follows // that g is one of the machine words |1..1| or |0..0|, and if g=|1..1|, then // s!=|0..0|. Note that we should read the words here as 2-complement numbers, // so that g is either -1 or 0, which means that we have to subtract B or not, // respectively, in the next step via c. // If c∈{0,1}, it follows that -1 <= s-c <= B-1, and the latter is the machine // double word |0..0|1..1| while -1 is |1..1|1..1|. Hence, p is also either one // of the machine words |0..0| or |1..1|. Finally, note that since g=|1..1|, // then s!=|0..0|, and therefore, 0<=s-c<=B-1, so that p=|0..0| follows. Hence, // g and p cannot be both |1..1|, so that g+p, g|p and !(g<->p) are all // equivalent and will either produce |0..0| or |1..1|, meaning the integers 0 // and -1, respectively. Changing their signs by subtracting from 0 yields // finally c'∈{0,1}. // ----------------------------------------------------------------------------- procedure NatSubSameLen([]nat x1,x2,y,nat len) { nat c,g,p,s,i; c = 0; // initial carry digit ∈{0,1} for(i=0..len-1) { // ripple from least to highest digit g,s = x1[i] - x2[i]; // g:generate carry; p∈{|0..0|,|1..1|} p,y[i] = s - c; // p:propagate carry; p∈{|0..0|,|1..1|} c = 0-(g+p); // g+p∈{|0..0|,|1..1|}, thus c∈{0,1} } y[len] = g+p; // non-negative iff y[len]=|0..0| return; // negative iff y[len]=|1..1| } // ----------------------------------------------------------------------------- // Radix-B Multiplication y[0..len1+len2-1] = x1[0..len1-1] * x2[0..len2-1] // ----------------------------------------------------------------------------- // For the correctness, note that the maximum of pr[i+j] + x1[i] * x2[j] + c is // (B-1)*B + (B-1) so that the maximum of the sum of the carry digits computed // in the algorithm is B-1 and therefore their addition does not overflow. // ----------------------------------------------------------------------------- procedure NatMul([]nat x1,x2,pr,nat len1,len2) { nat c1,c,p,i,j,k; for(j=0..len2-1) { // compute pr[j+len1..j] = pr[j+len1-1..j] + x1[len1-1..0] * x2[j] c = 0; for(i=0..len1-1) { c1,p = x1[i] * x2[j]; c,p = p + c; c = c + c1; k = i+j; c1,pr[k] = p + pr[k]; c = c + c1; } pr[len1+j] = c; } return; } // ----------------------------------------------------------------------------- // Radix-B Equality checking y = x1[0..len1-1] = x2[0..len2-1] // ----------------------------------------------------------------------------- function NatEquSameLen([]nat x1,x2,nat len) : bool { nat i; bool eq; i = len-1; eq = true; while(i>0 & eq) { eq = eq & x1[i]==x2[i]; i = i-1; } return eq; } // ----------------------------------------------------------------------------- // Radix-B Comparison y = x1[0..len1-1] ? x2[0..len2-1] yields -1 if x1<x2, // 0 if x1==x2 and +1 if x1>x2 // ----------------------------------------------------------------------------- function NatCmpSameLen([]nat x1,x2,nat len) : int { nat i; bool eq; i = len; eq = true; while(i>0 & eq) { i = i-1; eq = eq & x1[i]==x2[i]; } if(x1[i]==x2[i]) return +0; else if(x1[i]<x2[i]) return -1; else return +1; } // ----------------------------------------------------------------------------- // compute sum digits y[0..len] = x1[0..len-1] + x2[0..len-1] // ----------------------------------------------------------------------------- function NatAddCLA ([]nat x1,x2,y,c, []bool g,p, nat len): []nat { nat i,j,s,s1,len1; // loop variables len1 = len-1; // compute this constant once and forall // first row of generate and propagate signals for(i=0..len1) { c[i],y[i] = x1[i] + x2[i]; // preliminary sum and carry digits g[i] = c[i]>0; // generate signal at level 0 p[i] = y[i]==(nat)(-1); // propagate signal at level 0 } // compute generate/propagate signals by forward traversal on lookahead tree s = 1; while (s<len) { i = s-1; // middle index of considered subtree s = s+s; // width of current subtrees s=2^l j = s-1; // leftmost index of considered subtree if(j>len1) j = len1; while(i<len1) { g[j] = g[j] | g[i] & p[j]; p[j] = p[j] & p[i]; i = i+s; j = j+s; if(j>len1) j = len1; } } // compute carry signals c[i] by backward traversal of the lookahead tree c[0] = 0; // rightmost carry signal c[len] = (nat) (g[len1] | p[len1] & (bool)c[0]); // leftmost carry signal while(s>1) { s1 = s/2; i = s1; j = 0; while(i<len) { c[i] = (nat) (g[i-1] | p[i-1] & (bool)c[j]); i = i+s; j = j+s; } s = s1; } // compute the final sum digits with the carry and preliminary sum digits // note that we should have c[i] as either 1 or 0, but due to the casts from // boolean values to type nat, it becomes MaxNat and subtracting this is // the same as subtracting -1 (if we ignore overflows) for(i=1..len) { y[i] = y[i] - c[i]; } return y; } // ---------------------------------------------------------------------------- // NatAdd implements carry lookahead addition of radix-B numbers x1 and x2 which // both have length len. The sum is stored in array y whose length is len+1. The // carry digits are computed with a prefix tree with a depth of log(len) so that // the procedure could be computed in time log(len) using O(len) processing // elements. // ---------------------------------------------------------------------------- function NatAddCLA ([]nat x1,x2,[]nat y,c, []bool g,p, nat len): []nat { nat i,j,s,s1,len1; // loop variables len1 = len-1; // compute this constant once and forall // first row of generate and propagate signals for(i=0..len1) { c[i],y[i] = x1[i] + x2[i]; // preliminary sum and carry digits g[i] = c[i]>0; // generate signal at level 0 p[i] = y[i]==(nat)(-1); // propagate signal at level 0 } // compute generate/propagate signals by forward traversal on lookahead tree s = 1; while (s<len) { i = s-1; // middle index of considered subtree s = s+s; // width of current subtrees s=2^l j = s-1; // leftmost index of considered subtree if(j>len1) j = len1; while(i<len1) { g[j] = g[j] | g[i] & p[j]; p[j] = p[j] & p[i]; i = i+s; j = j+s; if(j>len1) j = len1; } } // compute carry signals c[i] by backward traversal of the lookahead tree c[0] = 0; // rightmost carry signal c[len] = (nat) (g[len1] | p[len1] & (bool)c[0]); // leftmost carry signal while(s>1) { s1 = s/2; i = s1; j = 0; while(i<len) { c[i] = (nat) (g[i-1] | p[i-1] & (bool)c[j]); i = i+s; j = j+s; } s = s1; } // compute the final sum digits with the carry and preliminary sum digits // note that we should have c[i] as either 1 or 0, but due to the casts from // boolean values to type nat, it becomes MaxNat and subtracting this is // the same as subtracting -1 (if we ignore overflows) for(i=1..len) { y[i] = y[i] - c[i]; } return y; } // ----------------------------------------------------------------------------- // compute sum digits y[0..len] = x1[0..len1-1] + x2[0..len2-1] // ----------------------------------------------------------------------------- procedure NatAdd([]nat x1,x2,sum,nat len1,len2) { nat c,g,p,s,i,lmin; // determine set of common digits if(len1<len2) lmin = len1; else lmin = len2; // carry ripple addition for common digits c = 0; for(i=0..lmin-1) { g,s = x1[i] + x2[i]; p,sum[i] = s + c; c = (nat) ((bool) g | (bool) p); } // propagate c through the final digits of the longer number if(len1<len2) { for(i=lmin..len2-1) c,sum[i] = x2[i] + c; sum[len2] = c; } else { for(i=lmin..len1-1) c,sum[i] = x1[i] + c; sum[len1] = c; } return; } thread test { [8]nat a; [6]nat b; [9]nat sm; nat i; for(i=0..5) a[i] = (nat) -1; b[0] = 1; NatAdd(a,b,sm,8,6); }