// ************************************************************************** //
//                                                                            //
//    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 version of the B-complement division algorithm is an optimization     //
// of IntDivModSeq0 that consists in the removal of the local variable z.     //
// Instead, the partial remainders z[N+i..i] are stored in r[N-1..0]@x[i..0]. //
// ** Note: We are only interested in even values of B **                     // 
// ************************************************************************** //

macro B = 4;    // base of the radix numbers
macro M = 4;    // number of digits used for x
macro N = 3;    // number of digits used for y
macro alpha(x) = (x<B/2 ? +x : +x-B);
macro gamma(y) = (y<0 ? y+B : y);
macro dval(x,i,k) = (i==k-1 ? alpha(x[i]) : +x[i]);
macro natval(x,m) = sum(i=0..m-1) (x[i] * exp(B,i));
macro intval(x,k) = sum(i=0..k-1) (dval(x,i,k) * exp(B,i));
macro z(i,j) = (j==0 ? x[i] : r[j-1]); // should be z[j+i] of IntDivModSeq0

module IntDivModSeq1([M]nat{B} ?x,[N]nat{B} ?y,[M+1]nat{B} q,[N]nat{B} r,event !rdy) {
    [M+N]nat{B} z;    // partial remainders
    [N]nat{B} s;      // auxiliary variable for sums
    [N+1]int{B} c;    // auxiliary variable for carry digits
    int{2} sgn_y;     // sign of y
    event isNegative; // is set when subtraction yields negative result

    // apply digit extension to x and put the extended digits in r
    for(j=0..N-1)
        r[j] = (x[M-1]<B/2 ? 0 : B-1);

    //-------------------------------------------------------------------------
    // determine leftmost digit of the quotient and the first partial remainder
    //-------------------------------------------------------------------------
    sgn_y = (y[N-1]<B/2 ? +1 :  -1);
    q[M]  = (x[M-1]<B/2 ?  0 : B-1);
    // compute r[N-1..0] = r[N-1..0] - alpha(q[M]) * sgn_y * y
    c[0] = 0;
    for(j=0..N-1)
        // note that -(B-1)<=sm<=2*B-1, thus -1<=sm/B<=1
        let(yJ = (j==N-1 ? alpha(y[j]) : y[j]))
        let(sm = r[j] + c[j] - alpha(q[M]) * sgn_y * yJ) {
        next(r[j]) = sm % B;
        c[j+1] = sm / B;
        }
    w0: pause;

    //-------------------------------------------------------------------------
    // computing the remaining quotient digits
    //-------------------------------------------------------------------------
    for(iup=0..M-1)
        let(i=M-1-iup) {// thus i=M-1..0
        // try digits q[i] = B-1,B-2,... until we reach one where 
        // z[N+i..i] - q[i] * sgn_y * y[N-1..0] remains non-negative
        next(q[i]) = B-1;
        do {
            w1: pause;
            // compute s[N..0] = z[N+i..i] - q[i] * sgn_y * y[N-1..0]
            c[0] = 0;
            for(j=0..N-1)
                // note that -B*(B-1)<=sm<=B-1, thus -(B-1)<=sm/B<=0
                let(yJ = (j==N-1 ? alpha(y[j]) : y[j]))
                let(sm = z(i,j) + c[j] - q[i] * sgn_y * yJ) {
                s[j]   = sm % B;
                c[j+1] = sm / B;
                }
            isNegative = (z(i,N) + c[N] < 0);
            if(isNegative) 
                next(q[i]) = q[i] - 1;
        } while(isNegative);

        // Here, q[i] is the largest digit with non-negative s[N..0].
        // We therefore commit the subtraction in that we copy the preliminary
        // result from s to z before computing the next digit q[i].
        w2: pause;
        for(j=0..N-1)
            r[j] = s[j];
    }

    emit(rdy);
    // final assertion
    if(intval(y,N)!=0) {
        assert(intval(x,M) == intval(q,M+1) * sgn_y * intval(y,N) + intval(r,N));
        assert(intval(r,N) < abs(intval(y,N)));
    }
}