ch24aok.c

稀疏矩阵、链表、图、队列、二叉树、多叉树、排序、遗传算法等的实现

    int iStatus = EXIT_SUCCESS;

    /*
    * This routine divides the single INT pointed at by the first    
    * parameter by the second INT pointed at by the second parameter.
    * The quotient is given in the first of the pair of INTs pointed 
    * at by the third parameter, and the remainder in the second of  
    * that pair.                                                     
    */

    iAnswerOD[0] = *iNumerator / *iDenominator;
    iAnswerOD[1] = *iNumerator % *iDenominator;

    return iStatus;
}

/* Divide a pair of INTs by one INT, getting quotient and remainder   */
int pairDivide ( INT * iNumerator, INT * iDenominator, INT * iAnswerPD )
{
    int iStatus = EXIT_SUCCESS;
    int j = 0;
    INT r = 0;
    INT w = 0;

    /*
    * This routine divides the pair of INTs pointed at by the first  
    * parameter by the single INT pointed at by the second parameter.
    * The answer is given in the first two of the triplet of INTs    
    * pointed at by the third parameter, and the remainder in the    
    * third of that triplet.                                         
    */

    /* The algorithm used is that of Knuth, Vol 3. 4.3.1, Ex 16       */

    r = 0;
    for (j=0;(j<=1);j++)
    {
         w = (r * BASE + *(iNumerator + j));
         *(iAnswerPD+j) = w / *(iDenominator);
         r = w % *(iDenominator);
    }
    *(iAnswerPD+2) = r;

    return iStatus;
}

/* Divide an INT array by one INT, getting just the quotient          */
int cascDivide ( INT * aNumerator, INT * iDenominator, INT ** aAnswer )
{
    /*
    * This routine allocates space for the answer, which the caller  
    * must, at some time, release                                    
    */

    /*
    * This routine divides the long number, in the INT array pointed 
    * to by the first argument, by the single INT pointed to by the  
    * second argument. Space for the result is allocated. The final  
    * output (pointed to by the pointer the third argument points to)
    * is normalized.                                                 
    * The sign of the single INT is assumed to be positive, and the
    * sign of the answer is that of the incoming long array.
    */

    int iStatus = EXIT_SUCCESS;
    int j = 0;
    int k = 0;
    INT r = 0;
    INT w = 0;
    int iSign = 1;
    INT * pInt;
    INT ** ppInt = &pInt;

    if ((aNumerator==NULL) || (iDenominator==NULL) || ((aAnswer)==NULL))
    {
        /*
        * If any of the incoming pointers is NULL, then return with  
        * error indication:                                          
        */
        iStatus = EXIT_FAILURE;
        return iStatus;
    }

    /* The algorithm used is again that of Knuth, Vol 3. 4.3.1, Ex 16 */

    /* Determine the length of the numerator, and remember its sign   */
    r = 0;
    k = abs(*aNumerator);
    aAllocate(k, ppInt);
    if (*aNumerator<0)
        iSign = -1;
    else
        iSign = 1;

    for (j=1; (j<=k); j++)
    {
         w = (r * BASE + *(aNumerator + j));
         *(pInt + j) = w / *(iDenominator);
         r = w % *(iDenominator);
    }

    /*
    * Note that at the very end we have an outgoing, residual        
    * carry. This is simply lost.                                    
    */

    /* Normalize the intermediate answer into the final answer:   */
    /* Note that this also frees the intermediate work area       */
    j = aNormalize ( ppInt, aAnswer );

    /* Apply the incoming sign to the answer:                         */
    if (iSign<0)
        **aAnswer = - (**aAnswer);

    return iStatus;
}

/* Divide an INT array by an INT array, getting just the quotient     */
int aDivide ( INT * aNumerator, INT * aDenominator, INT ** aAnswer )
{
    /*
    * This routine allocates space for the answer, which the caller  
    * must, at some time, release                                    
    */

    /*
    * This is the most complex of the routines, and uses the algorithm  
    * described by Knuth in "The Art of Computer Programming".   
    * This is an extension of the simple "pencil and paper" technique
    * that we were taught at school - and many of us have since forgotten.
    */

    int iStatus = EXIT_SUCCESS;
    int i = 0;
    int j = 0;
    int k = 0;
    int l = 0;
    int m = 0;
    int n = 0;
    int borrow = 0;
    int iNmax;
    int iDmax;
    INT carry = 0;
    INT iBase = BASE;
    INT wN[2];
    INT aW[3];
    INT aX[3];
    INT aY[3];
    INT aZ[3];
    INT bigD = 0;
    INT littleD = 0;
    /* Space for a local copy of the numerator. This copy will be     */
    /* long enough to be longer than the denominator.                 */
    INT * plNumerator = NULL;
    INT ** pplNumerator = &plNumerator;
    /* Space for a local copy of the Normalized denominator. This will*/
    /* ensure that the leading INT does not have value zero:          */
    INT *plDenominator = NULL;
    INT ** pplDenominator = & plDenominator;
    /* Space for a local copy of the resultant quotient. This copy is */
    /* prior to the required normalization.                           */
    INT * plQuotient = NULL;
    INT ** pplQuotient = &plQuotient;
    /* Space for three work areas required during the calculation.    */
    INT * plSubtrahend = NULL;
    INT ** pplSubtrahend = &plSubtrahend;
    INT * plMinuend = NULL;
    INT ** pplMinuend = &plMinuend;
    INT * pInter = NULL;
    INT ** ppInter = &pInter;

    /*
    * Normalization for division means (a) ensuring that the numerator
    * is longer than the denominator, and (b) the leading "digit" 
    * of the denominator is greater than floor(BASE/2). We can get   
    * this by multiplying both numerator and denominator by "D",     
    * where "D = floor((b-1)/d" and "d" is the leading "digit" of the
    * denominator. "floor" means "round to the next integer below".
    */

    /* Find the length required for the denominator copy, and allocate*/
    iDmax = abs(aDenominator[0]);
    k = aAllocate ( iDmax, ppInter );

    /* Copy from argument to local copy:                              */
    for (i=1;(i<=iDmax);i++)
        *(pInter+i) = *(aDenominator+i);

    /* And normalize the local copy of the denominator:               */
    k = aNormalize ( ppInter, ppInter );

    /*
    * Get the length of the normalized copy, and compute the scaling 
    * factor by which to multiply both the numerator and denominator: 
    */
    iDmax = *pInter;
    littleD = *(pInter+1);
    bigD = (BASE - 1) / littleD;

    /*
    * Get a scaled (multiplied-up) local copy of the denominator, and
    * release the previous work area:                                
    */
    cascMultiply ( &bigD, pInter, pplDenominator );
    free ( pInter );

    /*
    * Find out how much space to allocate for a local copy of the 
    * numerator, and get the space, multiplying up at the same time:    
    */
    pplNumerator = &plNumerator;
    plNumerator = NULL;
    cascMultiply ( &bigD, aNumerator, ppInter );

    /*
    * Perform "Division Normalize", which inserts a leading zero, and
    * release the intermediate work area:                            
    */
    k = aDivNormalize ( pInter, pplNumerator );
    free ( pInter );

    /*
    * Get the lengths of the normalized numerator, and compute the   
    * length of the quotient:                                        
    * Knuth D1:                                                      
    */
    iNmax = *plNumerator;
    n = iDmax;
    m = iNmax - iDmax;

    /*
    * The length of the quotient was one greater than the difference 
    * in length of the numerator and denominator. Remember that this 
    * routine requires the numerator to be at least as long as the   
    * denominator. The difference in length (and not the "difference 
    * plus one") is required later in the calculation. We have already
    * "lengthened" the numerator, so we need use only the real 
    * difference in length now:                                      
    */

    /* Allocate space for the quotient pre-normalization:             */
    pplQuotient = &plQuotient;
    k = aAllocate ( m, pplQuotient );

    /* Allocate space for holding the intermediate subtraction        */
    pplSubtrahend = &plSubtrahend;
    aAllocate ( n, pplSubtrahend );

    /* Allocate space for holding numerator portion for subtraction:  */
    pplMinuend = &plMinuend;
    aAllocate ( n + 1, pplMinuend );

    /*
    * Division takes place from left to right: we start with the most
    * significant digits - those on the left:                        
    * Knuth D2:                                                      
    */
    for (j=1;(j<=m);j++)
    {
        /* Knuth D3:                                                  */
        wN[0] = *(plNumerator + j); 
        wN[1] = (j==(m+1)) ? 0 : *(plNumerator + j + 1); 
        k = pairDivide ( wN, plDenominator + 1, aW );

        /*
        * Now test for the cases where our putative quotient is too  
        * high. Note that it will sometimes be one too high, and very
        * rarely two too high.                                       
        * If quotient is equal to base, it is too high:              
        */
D3A:    if ((aW[0]==0) && (aW[1]==BASE))
            goto D3B;
        /* If quotient is equal to zero, it cannot be too high:       */
        if ((aW[0]==0) && (aW[1]==0))
            goto D3C;
        /* If there are no more places in the denominator, then there */
        /* is nothing else against which to test, so it's not too high*/
        if (*plDenominator<2)
            goto D3C;

        /* Now we have to compute the inequality of Knuth step D3     */
        pairMultiply (*(plDenominator+2),aW[0],aX);
        pairMultiply (*(plDenominator+2),aW[1],aY);
        aY[0] += aY[1] / BASE; 
        aY[1] %= BASE;
        aX[1] += aY[0] / BASE; 
        aY[0] %= BASE;
        aZ[0] = aW[2];        
        aZ[1] = *(plNumerator + j - m + 3); 
        carry = aZ[1] / BASE; 
        aZ[1] %= BASE;

        if (aX[1]>carry)
            goto D3B;
        if (aY[0]>aZ[0])
            goto D3B;
        if (aY[0]<aZ[0]) 
            goto D3C;
        if (aY[1]>aZ[1]) 
            goto D3B;
        goto D3C;

D3B:    {
            /*
            * The putative quotient is too high. Reduce it, and go   
            * back round to see whether it is still too high:         
            */
            aW[1] -= 1;
            aW[2] += *(plDenominator + 1);
            if (aW[2]<BASE)
                goto D3A;
        }
D3C:    /*
        * This quotient digit is correct.                            
        * Knuth D4:                                                  
        * We now have to subtract from the local copy of the numerator
        * the denominator times the last quotient "digit":        
        */
        if (aW[1]!=0)
        {
            for (i=1; (i<=(n+1)); i++)
            {
                *(plMinuend + i) = *(plNumerator + m - j + i);
                if (i<(n + 1))
                    *(plSubtrahend + i) = *(plDenominator + i);
            }

            /* cascade multiply denominator by quotient digit */
            cascMultiply ( &aW[1], plDenominator, ppInter );

            /* subtract:                                      */
            if ( *pInter > *plMinuend )
                borrow = 1;
            else
            {
                carry = 0;
                k = *pInter;
                l = (*plMinuend) - k;
                for (i=k; (i>=0); i--)
                {
                    if ((i + l)>0)
                    {
                        *(plMinuend + i + l) -= carry; 
                        *(plMinuend + i + l) -= (i>0) ? *(pInter+i) : 0 ;
                        carry = 0;
                        if (*(plMinuend + i + l)<0)
                        {
                            *(plMinuend + i + l) += BASE;
                            carry = 1;
                        }
                        else
                            carry = 0;
                    }   /* end of "if/then" safe to make replacement  */
                }   /* end of "for/i" computing digits to subtract    */

                /*
                * Test the sign of the answer, and if it is negative then
                * set the "borrow" flag. If it is positive or zero then
                * replace a portion of the local copy of the numerator:
                */
                if (*(plMinuend+1)<0)
                    borrow = 1;
                else
                {
                    borrow = 0;
                    /* replace a portion of our local numerator copy  */
                    for (i=n+1; (i>0); i--)
                    {
                        if ((i + l) <= iNmax)
                           *(plNumerator + i + l) = *(plMinuend + i);
                    }   /* end of "for/i" replace portion of numerator    */
                }   /* end of "if/else" numerator portion was replaced    */
            }   /* end of "if/else" lengths allowed numerator subtraction */

            /* Knuth D5:                                                  */
            *(plQuotient + j) = aW[1];

            /* If the remainder was negative then add back:               */
            if (borrow != 0)
            {
                /* Knuth D6:                                              */