[thirdparty/openssl.git] / crypto / bn / bn_gcd.c

/* crypto/bn/bn_gcd.c */
/* Copyright (C) 1995-1998 Eric Young (eay@cryptsoft.com)
 * All rights reserved.
 *
 * This package is an SSL implementation written
 * by Eric Young (eay@cryptsoft.com).
 * The implementation was written so as to conform with Netscapes SSL.
 *
 * This library is free for commercial and non-commercial use as long as
 * the following conditions are aheared to.  The following conditions
 * apply to all code found in this distribution, be it the RC4, RSA,
 * lhash, DES, etc., code; not just the SSL code.  The SSL documentation
 * included with this distribution is covered by the same copyright terms
 * except that the holder is Tim Hudson (tjh@cryptsoft.com).
 *
 * Copyright remains Eric Young's, and as such any Copyright notices in
 * the code are not to be removed.
 * If this package is used in a product, Eric Young should be given attribution
 * as the author of the parts of the library used.
 * This can be in the form of a textual message at program startup or
 * in documentation (online or textual) provided with the package.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 * 1. Redistributions of source code must retain the copyright
 *    notice, this list of conditions and the following disclaimer.
 * 2. Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in the
 *    documentation and/or other materials provided with the distribution.
 * 3. All advertising materials mentioning features or use of this software
 *    must display the following acknowledgement:
 *    "This product includes cryptographic software written by
 *     Eric Young (eay@cryptsoft.com)"
 *    The word 'cryptographic' can be left out if the rouines from the library
 *    being used are not cryptographic related :-).
 * 4. If you include any Windows specific code (or a derivative thereof) from
 *    the apps directory (application code) you must include an acknowledgement:
 *    "This product includes software written by Tim Hudson (tjh@cryptsoft.com)"
 *
 * THIS SOFTWARE IS PROVIDED BY ERIC YOUNG ``AS IS'' AND
 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
 * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
 * SUCH DAMAGE.
 *
 * The licence and distribution terms for any publically available version or
 * derivative of this code cannot be changed.  i.e. this code cannot simply be
 * copied and put under another distribution licence
 * [including the GNU Public Licence.]
 */
/* ====================================================================
 * Copyright (c) 1998-2001 The OpenSSL Project.  All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 *    notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in
 *    the documentation and/or other materials provided with the
 *    distribution.
 *
 * 3. All advertising materials mentioning features or use of this
 *    software must display the following acknowledgment:
 *    "This product includes software developed by the OpenSSL Project
 *    for use in the OpenSSL Toolkit. (http://www.openssl.org/)"
 *
 * 4. The names "OpenSSL Toolkit" and "OpenSSL Project" must not be used to
 *    endorse or promote products derived from this software without
 *    prior written permission. For written permission, please contact
 *    openssl-core@openssl.org.
 *
 * 5. Products derived from this software may not be called "OpenSSL"
 *    nor may "OpenSSL" appear in their names without prior written
 *    permission of the OpenSSL Project.
 *
 * 6. Redistributions of any form whatsoever must retain the following
 *    acknowledgment:
 *    "This product includes software developed by the OpenSSL Project
 *    for use in the OpenSSL Toolkit (http://www.openssl.org/)"
 *
 * THIS SOFTWARE IS PROVIDED BY THE OpenSSL PROJECT ``AS IS'' AND ANY
 * EXPRESSED OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
 * PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE OpenSSL PROJECT OR
 * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
 * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
 * OF THE POSSIBILITY OF SUCH DAMAGE.
 * ====================================================================
 *
 * This product includes cryptographic software written by Eric Young
 * (eay@cryptsoft.com).  This product includes software written by Tim
 * Hudson (tjh@cryptsoft.com).
 *
 */

#include "cryptlib.h"
#include "bn_lcl.h"

static BIGNUM *euclid(BIGNUM *a, BIGNUM *b);

int BN_gcd(BIGNUM *r, const BIGNUM *in_a, const BIGNUM *in_b, BN_CTX *ctx)
{
    BIGNUM *a, *b, *t;
    int ret = 0;

    bn_check_top(in_a);
    bn_check_top(in_b);

    BN_CTX_start(ctx);
    a = BN_CTX_get(ctx);
    b = BN_CTX_get(ctx);
    if (a == NULL || b == NULL)
        goto err;

    if (BN_copy(a, in_a) == NULL)
        goto err;
    if (BN_copy(b, in_b) == NULL)
        goto err;
    a->neg = 0;
    b->neg = 0;

    if (BN_cmp(a, b) < 0) {
        t = a;
        a = b;
        b = t;
    }
    t = euclid(a, b);
    if (t == NULL)
        goto err;

    if (BN_copy(r, t) == NULL)
        goto err;
    ret = 1;
 err:
    BN_CTX_end(ctx);
    bn_check_top(r);
    return (ret);
}

static BIGNUM *euclid(BIGNUM *a, BIGNUM *b)
{
    BIGNUM *t;
    int shifts = 0;

    bn_check_top(a);
    bn_check_top(b);

    /* 0 <= b <= a */
    while (!BN_is_zero(b)) {
        /* 0 < b <= a */

        if (BN_is_odd(a)) {
            if (BN_is_odd(b)) {
                if (!BN_sub(a, a, b))
                    goto err;
                if (!BN_rshift1(a, a))
                    goto err;
                if (BN_cmp(a, b) < 0) {
                    t = a;
                    a = b;
                    b = t;
                }
            } else {            /* a odd - b even */

                if (!BN_rshift1(b, b))
                    goto err;
                if (BN_cmp(a, b) < 0) {
                    t = a;
                    a = b;
                    b = t;
                }
            }
        } else {                /* a is even */

            if (BN_is_odd(b)) {
                if (!BN_rshift1(a, a))
                    goto err;
                if (BN_cmp(a, b) < 0) {
                    t = a;
                    a = b;
                    b = t;
                }
            } else {            /* a even - b even */

                if (!BN_rshift1(a, a))
                    goto err;
                if (!BN_rshift1(b, b))
                    goto err;
                shifts++;
            }
        }
        /* 0 <= b <= a */
    }

    if (shifts) {
        if (!BN_lshift(a, a, shifts))
            goto err;
    }
    bn_check_top(a);
    return (a);
 err:
    return (NULL);
}

/* solves ax == 1 (mod n) */
static BIGNUM *BN_mod_inverse_no_branch(BIGNUM *in,
                                        const BIGNUM *a, const BIGNUM *n,
                                        BN_CTX *ctx);
BIGNUM *BN_mod_inverse(BIGNUM *in, const BIGNUM *a, const BIGNUM *n,
                       BN_CTX *ctx)
{
    BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
    BIGNUM *ret = NULL;
    int sign;

    if ((BN_get_flags(a, BN_FLG_CONSTTIME) != 0)
        || (BN_get_flags(n, BN_FLG_CONSTTIME) != 0)) {
        return BN_mod_inverse_no_branch(in, a, n, ctx);
    }

    bn_check_top(a);
    bn_check_top(n);

    BN_CTX_start(ctx);
    A = BN_CTX_get(ctx);
    B = BN_CTX_get(ctx);
    X = BN_CTX_get(ctx);
    D = BN_CTX_get(ctx);
    M = BN_CTX_get(ctx);
    Y = BN_CTX_get(ctx);
    T = BN_CTX_get(ctx);
    if (T == NULL)
        goto err;

    if (in == NULL)
        R = BN_new();
    else
        R = in;
    if (R == NULL)
        goto err;

    BN_one(X);
    BN_zero(Y);
    if (BN_copy(B, a) == NULL)
        goto err;
    if (BN_copy(A, n) == NULL)
        goto err;
    A->neg = 0;
    if (B->neg || (BN_ucmp(B, A) >= 0)) {
        if (!BN_nnmod(B, B, A, ctx))
            goto err;
    }
    sign = -1;
        /*-
         * From  B = a mod |n|,  A = |n|  it follows that
         *
         *      0 <= B < A,
         *     -sign*X*a  ==  B   (mod |n|),
         *      sign*Y*a  ==  A   (mod |n|).
         */

    if (BN_is_odd(n) && (BN_num_bits(n) <= (BN_BITS <= 32 ? 450 : 2048))) {
        /*
         * Binary inversion algorithm; requires odd modulus. This is faster
         * than the general algorithm if the modulus is sufficiently small
         * (about 400 .. 500 bits on 32-bit sytems, but much more on 64-bit
         * systems)
         */
        int shift;

        while (!BN_is_zero(B)) {
                        /*-
                         *      0 < B < |n|,
                         *      0 < A <= |n|,
                         * (1) -sign*X*a  ==  B   (mod |n|),
                         * (2)  sign*Y*a  ==  A   (mod |n|)
                         */

            /*
             * Now divide B by the maximum possible power of two in the
             * integers, and divide X by the same value mod |n|. When we're
             * done, (1) still holds.
             */
            shift = 0;
            while (!BN_is_bit_set(B, shift)) { /* note that 0 < B */
                shift++;

                if (BN_is_odd(X)) {
                    if (!BN_uadd(X, X, n))
                        goto err;
                }
                /*
                 * now X is even, so we can easily divide it by two
                 */
                if (!BN_rshift1(X, X))
                    goto err;
            }
            if (shift > 0) {
                if (!BN_rshift(B, B, shift))
                    goto err;
            }

            /*
             * Same for A and Y.  Afterwards, (2) still holds.
             */
            shift = 0;
            while (!BN_is_bit_set(A, shift)) { /* note that 0 < A */
                shift++;

                if (BN_is_odd(Y)) {
                    if (!BN_uadd(Y, Y, n))
                        goto err;
                }
                /* now Y is even */
                if (!BN_rshift1(Y, Y))
                    goto err;
            }
            if (shift > 0) {
                if (!BN_rshift(A, A, shift))
                    goto err;
            }

                        /*-
                         * We still have (1) and (2).
                         * Both  A  and  B  are odd.
                         * The following computations ensure that
                         *
                         *     0 <= B < |n|,
                         *      0 < A < |n|,
                         * (1) -sign*X*a  ==  B   (mod |n|),
                         * (2)  sign*Y*a  ==  A   (mod |n|),
                         *
                         * and that either  A  or  B  is even in the next iteration.
                         */
            if (BN_ucmp(B, A) >= 0) {
                /* -sign*(X + Y)*a == B - A  (mod |n|) */
                if (!BN_uadd(X, X, Y))
                    goto err;
                /*
                 * NB: we could use BN_mod_add_quick(X, X, Y, n), but that
                 * actually makes the algorithm slower
                 */
                if (!BN_usub(B, B, A))
                    goto err;
            } else {
                /*  sign*(X + Y)*a == A - B  (mod |n|) */
                if (!BN_uadd(Y, Y, X))
                    goto err;
                /*
                 * as above, BN_mod_add_quick(Y, Y, X, n) would slow things
                 * down
                 */
                if (!BN_usub(A, A, B))
                    goto err;
            }
        }
    } else {
        /* general inversion algorithm */

        while (!BN_is_zero(B)) {
            BIGNUM *tmp;

                        /*-
                         *      0 < B < A,
                         * (*) -sign*X*a  ==  B   (mod |n|),
                         *      sign*Y*a  ==  A   (mod |n|)
                         */

            /* (D, M) := (A/B, A%B) ... */
            if (BN_num_bits(A) == BN_num_bits(B)) {
                if (!BN_one(D))
                    goto err;
                if (!BN_sub(M, A, B))
                    goto err;
            } else if (BN_num_bits(A) == BN_num_bits(B) + 1) {
                /* A/B is 1, 2, or 3 */
                if (!BN_lshift1(T, B))
                    goto err;
                if (BN_ucmp(A, T) < 0) {
                    /* A < 2*B, so D=1 */
                    if (!BN_one(D))
                        goto err;
                    if (!BN_sub(M, A, B))
                        goto err;
                } else {
                    /* A >= 2*B, so D=2 or D=3 */
                    if (!BN_sub(M, A, T))
                        goto err;
                    if (!BN_add(D, T, B))
                        goto err; /* use D (:= 3*B) as temp */
                    if (BN_ucmp(A, D) < 0) {
                        /* A < 3*B, so D=2 */
                        if (!BN_set_word(D, 2))
                            goto err;
                        /*
                         * M (= A - 2*B) already has the correct value
                         */
                    } else {
                        /* only D=3 remains */
                        if (!BN_set_word(D, 3))
                            goto err;
                        /*
                         * currently M = A - 2*B, but we need M = A - 3*B
                         */
                        if (!BN_sub(M, M, B))
                            goto err;
                    }
                }
            } else {
                if (!BN_div(D, M, A, B, ctx))
                    goto err;
            }

                        /*-
                         * Now
                         *      A = D*B + M;
                         * thus we have
                         * (**)  sign*Y*a  ==  D*B + M   (mod |n|).
                         */

            tmp = A;            /* keep the BIGNUM object, the value does not
                                 * matter */

            /* (A, B) := (B, A mod B) ... */
            A = B;
            B = M;
            /* ... so we have  0 <= B < A  again */

                        /*-
                         * Since the former  M  is now  B  and the former  B  is now  A,
                         * (**) translates into
                         *       sign*Y*a  ==  D*A + B    (mod |n|),
                         * i.e.
                         *       sign*Y*a - D*A  ==  B    (mod |n|).
                         * Similarly, (*) translates into
                         *      -sign*X*a  ==  A          (mod |n|).
                         *
                         * Thus,
                         *   sign*Y*a + D*sign*X*a  ==  B  (mod |n|),
                         * i.e.
                         *        sign*(Y + D*X)*a  ==  B  (mod |n|).
                         *
                         * So if we set  (X, Y, sign) := (Y + D*X, X, -sign),  we arrive back at
                         *      -sign*X*a  ==  B   (mod |n|),
                         *       sign*Y*a  ==  A   (mod |n|).
                         * Note that  X  and  Y  stay non-negative all the time.
                         */

            /*
             * most of the time D is very small, so we can optimize tmp :=
             * D*X+Y
             */
            if (BN_is_one(D)) {
                if (!BN_add(tmp, X, Y))
                    goto err;
            } else {
                if (BN_is_word(D, 2)) {
                    if (!BN_lshift1(tmp, X))
                        goto err;
                } else if (BN_is_word(D, 4)) {
                    if (!BN_lshift(tmp, X, 2))
                        goto err;
                } else if (D->top == 1) {
                    if (!BN_copy(tmp, X))
                        goto err;
                    if (!BN_mul_word(tmp, D->d[0]))
                        goto err;
                } else {
                    if (!BN_mul(tmp, D, X, ctx))
                        goto err;
                }
                if (!BN_add(tmp, tmp, Y))
                    goto err;
            }

            M = Y;              /* keep the BIGNUM object, the value does not
                                 * matter */
            Y = X;
            X = tmp;
            sign = -sign;
        }
    }

        /*-
         * The while loop (Euclid's algorithm) ends when
         *      A == gcd(a,n);
         * we have
         *       sign*Y*a  ==  A  (mod |n|),
         * where  Y  is non-negative.
         */

    if (sign < 0) {
        if (!BN_sub(Y, n, Y))
            goto err;
    }
    /* Now  Y*a  ==  A  (mod |n|).  */

    if (BN_is_one(A)) {
        /* Y*a == 1  (mod |n|) */
        if (!Y->neg && BN_ucmp(Y, n) < 0) {
            if (!BN_copy(R, Y))
                goto err;
        } else {
            if (!BN_nnmod(R, Y, n, ctx))
                goto err;
        }
    } else {
        BNerr(BN_F_BN_MOD_INVERSE, BN_R_NO_INVERSE);
        goto err;
    }
    ret = R;
 err:
    if ((ret == NULL) && (in == NULL))
        BN_free(R);
    BN_CTX_end(ctx);
    bn_check_top(ret);
    return (ret);
}

/*
 * BN_mod_inverse_no_branch is a special version of BN_mod_inverse. It does
 * not contain branches that may leak sensitive information.
 */
static BIGNUM *BN_mod_inverse_no_branch(BIGNUM *in,
                                        const BIGNUM *a, const BIGNUM *n,
                                        BN_CTX *ctx)
{
    BIGNUM *A, *B, *X, *Y, *M, *D, *T, *R = NULL;
    BIGNUM local_A, local_B;
    BIGNUM *pA, *pB;
    BIGNUM *ret = NULL;
    int sign;

    bn_check_top(a);
    bn_check_top(n);

    BN_CTX_start(ctx);
    A = BN_CTX_get(ctx);
    B = BN_CTX_get(ctx);
    X = BN_CTX_get(ctx);
    D = BN_CTX_get(ctx);
    M = BN_CTX_get(ctx);
    Y = BN_CTX_get(ctx);
    T = BN_CTX_get(ctx);
    if (T == NULL)
        goto err;

    if (in == NULL)
        R = BN_new();
    else
        R = in;
    if (R == NULL)
        goto err;

    BN_one(X);
    BN_zero(Y);
    if (BN_copy(B, a) == NULL)
        goto err;
    if (BN_copy(A, n) == NULL)
        goto err;
    A->neg = 0;

    if (B->neg || (BN_ucmp(B, A) >= 0)) {
        /*
         * Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
         * BN_div_no_branch will be called eventually.
         */
        pB = &local_B;
        BN_with_flags(pB, B, BN_FLG_CONSTTIME);
        if (!BN_nnmod(B, pB, A, ctx))
            goto err;
    }
    sign = -1;
        /*-
         * From  B = a mod |n|,  A = |n|  it follows that
         *
         *      0 <= B < A,
         *     -sign*X*a  ==  B   (mod |n|),
         *      sign*Y*a  ==  A   (mod |n|).
         */

    while (!BN_is_zero(B)) {
        BIGNUM *tmp;

                /*-
                 *      0 < B < A,
                 * (*) -sign*X*a  ==  B   (mod |n|),
                 *      sign*Y*a  ==  A   (mod |n|)
                 */

        /*
         * Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
         * BN_div_no_branch will be called eventually.
         */
        pA = &local_A;
        BN_with_flags(pA, A, BN_FLG_CONSTTIME);

        /* (D, M) := (A/B, A%B) ... */
        if (!BN_div(D, M, pA, B, ctx))
            goto err;

                /*-
                 * Now
                 *      A = D*B + M;
                 * thus we have
                 * (**)  sign*Y*a  ==  D*B + M   (mod |n|).
                 */

        tmp = A;                /* keep the BIGNUM object, the value does not
                                 * matter */

        /* (A, B) := (B, A mod B) ... */
        A = B;
        B = M;
        /* ... so we have  0 <= B < A  again */

                /*-
                 * Since the former  M  is now  B  and the former  B  is now  A,
                 * (**) translates into
                 *       sign*Y*a  ==  D*A + B    (mod |n|),
                 * i.e.
                 *       sign*Y*a - D*A  ==  B    (mod |n|).
                 * Similarly, (*) translates into
                 *      -sign*X*a  ==  A          (mod |n|).
                 *
                 * Thus,
                 *   sign*Y*a + D*sign*X*a  ==  B  (mod |n|),
                 * i.e.
                 *        sign*(Y + D*X)*a  ==  B  (mod |n|).
                 *
                 * So if we set  (X, Y, sign) := (Y + D*X, X, -sign),  we arrive back at
                 *      -sign*X*a  ==  B   (mod |n|),
                 *       sign*Y*a  ==  A   (mod |n|).
                 * Note that  X  and  Y  stay non-negative all the time.
                 */

        if (!BN_mul(tmp, D, X, ctx))
            goto err;
        if (!BN_add(tmp, tmp, Y))
            goto err;

        M = Y;                  /* keep the BIGNUM object, the value does not
                                 * matter */
        Y = X;
        X = tmp;
        sign = -sign;
    }

        /*-
         * The while loop (Euclid's algorithm) ends when
         *      A == gcd(a,n);
         * we have
         *       sign*Y*a  ==  A  (mod |n|),
         * where  Y  is non-negative.
         */

    if (sign < 0) {
        if (!BN_sub(Y, n, Y))
            goto err;
    }
    /* Now  Y*a  ==  A  (mod |n|).  */

    if (BN_is_one(A)) {
        /* Y*a == 1  (mod |n|) */
        if (!Y->neg && BN_ucmp(Y, n) < 0) {
            if (!BN_copy(R, Y))
                goto err;
        } else {
            if (!BN_nnmod(R, Y, n, ctx))
                goto err;
        }
    } else {
        BNerr(BN_F_BN_MOD_INVERSE_NO_BRANCH, BN_R_NO_INVERSE);
        goto err;
    }
    ret = R;
 err:
    if ((ret == NULL) && (in == NULL))
        BN_free(R);
    BN_CTX_end(ctx);
    bn_check_top(ret);
    return (ret);
}
Commit	Line	Data
d02b48c6	1	/* crypto/bn/bn_gcd.c */
58964a49	2	/* Copyright (C) 1995-1998 Eric Young (eay@cryptsoft.com)
d02b48c6 RE	3	* All rights reserved.
	4	*
	5	* This package is an SSL implementation written
	6	* by Eric Young (eay@cryptsoft.com).
	7	* The implementation was written so as to conform with Netscapes SSL.
40720ce3	8	*
d02b48c6 RE	9	* This library is free for commercial and non-commercial use as long as
	10	* the following conditions are aheared to. The following conditions
	11	* apply to all code found in this distribution, be it the RC4, RSA,
	12	* lhash, DES, etc., code; not just the SSL code. The SSL documentation
	13	* included with this distribution is covered by the same copyright terms
	14	* except that the holder is Tim Hudson (tjh@cryptsoft.com).
40720ce3	15	*
d02b48c6 RE	16	* Copyright remains Eric Young's, and as such any Copyright notices in
	17	* the code are not to be removed.
	18	* If this package is used in a product, Eric Young should be given attribution
	19	* as the author of the parts of the library used.
	20	* This can be in the form of a textual message at program startup or
	21	* in documentation (online or textual) provided with the package.
40720ce3	22	*
d02b48c6 RE	23	* Redistribution and use in source and binary forms, with or without
	24	* modification, are permitted provided that the following conditions
	25	* are met:
	26	* 1. Redistributions of source code must retain the copyright
	27	* notice, this list of conditions and the following disclaimer.
	28	* 2. Redistributions in binary form must reproduce the above copyright
	29	* notice, this list of conditions and the following disclaimer in the
	30	* documentation and/or other materials provided with the distribution.
	31	* 3. All advertising materials mentioning features or use of this software
	32	* must display the following acknowledgement:
	33	* "This product includes cryptographic software written by
	34	* Eric Young (eay@cryptsoft.com)"
	35	* The word 'cryptographic' can be left out if the rouines from the library
	36	* being used are not cryptographic related :-).
40720ce3	37	* 4. If you include any Windows specific code (or a derivative thereof) from
d02b48c6 RE	38	* the apps directory (application code) you must include an acknowledgement:
d02b48c6 RE	39	* "This product includes software written by Tim Hudson (tjh@cryptsoft.com)"
40720ce3	40	*
d02b48c6 RE	41	* THIS SOFTWARE IS PROVIDED BY ERIC YOUNG ``AS IS'' AND
	42	* ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
	43	* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
	44	* ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
	45	* FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
	46	* DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
	47	* OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
	48	* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
	49	* LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
	50	* OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
	51	* SUCH DAMAGE.
40720ce3	52	*
d02b48c6 RE	53	* The licence and distribution terms for any publically available version or
	54	* derivative of this code cannot be changed. i.e. this code cannot simply be
	55	* copied and put under another distribution licence
	56	* [including the GNU Public Licence.]
	57	*/
dcbd0d74	58	/* ====================================================================
7d0d0996	59	* Copyright (c) 1998-2001 The OpenSSL Project. All rights reserved.
dcbd0d74 BM	60	*
	61	* Redistribution and use in source and binary forms, with or without
	62	* modification, are permitted provided that the following conditions
	63	* are met:
	64	*
	65	* 1. Redistributions of source code must retain the above copyright
40720ce3	66	* notice, this list of conditions and the following disclaimer.
dcbd0d74 BM	67	*
	68	* 2. Redistributions in binary form must reproduce the above copyright
	69	* notice, this list of conditions and the following disclaimer in
	70	* the documentation and/or other materials provided with the
	71	* distribution.
	72	*
	73	* 3. All advertising materials mentioning features or use of this
	74	* software must display the following acknowledgment:
	75	* "This product includes software developed by the OpenSSL Project
	76	* for use in the OpenSSL Toolkit. (http://www.openssl.org/)"
	77	*
	78	* 4. The names "OpenSSL Toolkit" and "OpenSSL Project" must not be used to
	79	* endorse or promote products derived from this software without
	80	* prior written permission. For written permission, please contact
	81	* openssl-core@openssl.org.
	82	*
	83	* 5. Products derived from this software may not be called "OpenSSL"
	84	* nor may "OpenSSL" appear in their names without prior written
	85	* permission of the OpenSSL Project.
	86	*
	87	* 6. Redistributions of any form whatsoever must retain the following
	88	* acknowledgment:
	89	* "This product includes software developed by the OpenSSL Project
	90	* for use in the OpenSSL Toolkit (http://www.openssl.org/)"
	91	*
	92	* THIS SOFTWARE IS PROVIDED BY THE OpenSSL PROJECT ``AS IS'' AND ANY
	93	* EXPRESSED OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
	94	* IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
	95	* PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE OpenSSL PROJECT OR
	96	* ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
	97	* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
	98	* NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
	99	* LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
	100	* HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
	101	* STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
	102	* ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
	103	* OF THE POSSIBILITY OF SUCH DAMAGE.
	104	* ====================================================================
	105	*
	106	* This product includes cryptographic software written by Eric Young
	107	* (eay@cryptsoft.com). This product includes software written by Tim
	108	* Hudson (tjh@cryptsoft.com).
	109	*
	110	*/
d02b48c6	111
d02b48c6 RE	112	#include "cryptlib.h"
	113	#include "bn_lcl.h"
	114
d02b48c6	115	static BIGNUM euclid(BIGNUM a, BIGNUM *b);
9b141126	116
cbd48ba6	117	int BN_gcd(BIGNUM r, const BIGNUM in_a, const BIGNUM in_b, BN_CTX ctx)
40720ce3 MC	118	{
	119	BIGNUM a, b, *t;
	120	int ret = 0;
	121
	122	bn_check_top(in_a);
	123	bn_check_top(in_b);
	124
	125	BN_CTX_start(ctx);
	126	a = BN_CTX_get(ctx);
	127	b = BN_CTX_get(ctx);
	128	if (a == NULL \|\| b == NULL)
	129	goto err;
	130
	131	if (BN_copy(a, in_a) == NULL)
	132	goto err;
	133	if (BN_copy(b, in_b) == NULL)
	134	goto err;
	135	a->neg = 0;
	136	b->neg = 0;
	137
	138	if (BN_cmp(a, b) < 0) {
	139	t = a;
	140	a = b;
	141	b = t;
	142	}
	143	t = euclid(a, b);
	144	if (t == NULL)
	145	goto err;
	146
	147	if (BN_copy(r, t) == NULL)
	148	goto err;
	149	ret = 1;
	150	err:
	151	BN_CTX_end(ctx);
	152	bn_check_top(r);
	153	return (ret);
	154	}
d02b48c6	155
6b691a5c	156	static BIGNUM euclid(BIGNUM a, BIGNUM *b)
40720ce3 MC	157	{
	158	BIGNUM *t;
	159	int shifts = 0;
	160
	161	bn_check_top(a);
	162	bn_check_top(b);
	163
	164	/* 0 <= b <= a */
	165	while (!BN_is_zero(b)) {
	166	/* 0 < b <= a */
	167
	168	if (BN_is_odd(a)) {
	169	if (BN_is_odd(b)) {
	170	if (!BN_sub(a, a, b))
	171	goto err;
	172	if (!BN_rshift1(a, a))
	173	goto err;
	174	if (BN_cmp(a, b) < 0) {
	175	t = a;
	176	a = b;
	177	b = t;
	178	}
	179	} else { /* a odd - b even */
	180
	181	if (!BN_rshift1(b, b))
	182	goto err;
	183	if (BN_cmp(a, b) < 0) {
	184	t = a;
	185	a = b;
	186	b = t;
	187	}
	188	}
	189	} else { /* a is even */
	190
	191	if (BN_is_odd(b)) {
	192	if (!BN_rshift1(a, a))
	193	goto err;
	194	if (BN_cmp(a, b) < 0) {
	195	t = a;
	196	a = b;
	197	b = t;
	198	}
	199	} else { /* a even - b even */
	200
	201	if (!BN_rshift1(a, a))
	202	goto err;
	203	if (!BN_rshift1(b, b))
	204	goto err;
	205	shifts++;
	206	}
	207	}
	208	/* 0 <= b <= a */
	209	}
	210
	211	if (shifts) {
	212	if (!BN_lshift(a, a, shifts))
	213	goto err;
	214	}
	215	bn_check_top(a);
	216	return (a);
	217	err:
	218	return (NULL);
	219	}
dcbd0d74	220
d02b48c6	221	/* solves ax == 1 (mod n) */
283aedf4	222	static BIGNUM BN_mod_inverse_no_branch(BIGNUM in,
40720ce3 MC	223	const BIGNUM a, const BIGNUM n,
	224	BN_CTX *ctx);
	225	BIGNUM BN_mod_inverse(BIGNUM in, const BIGNUM a, const BIGNUM n,
	226	BN_CTX *ctx)
	227	{
	228	BIGNUM A, B, X, Y, M, D, T, R = NULL;
	229	BIGNUM *ret = NULL;
	230	int sign;
	231
	232	if ((BN_get_flags(a, BN_FLG_CONSTTIME) != 0)
	233	\|\| (BN_get_flags(n, BN_FLG_CONSTTIME) != 0)) {
	234	return BN_mod_inverse_no_branch(in, a, n, ctx);
	235	}
	236
	237	bn_check_top(a);
	238	bn_check_top(n);
	239
	240	BN_CTX_start(ctx);
	241	A = BN_CTX_get(ctx);
	242	B = BN_CTX_get(ctx);
	243	X = BN_CTX_get(ctx);
	244	D = BN_CTX_get(ctx);
	245	M = BN_CTX_get(ctx);
	246	Y = BN_CTX_get(ctx);
	247	T = BN_CTX_get(ctx);
	248	if (T == NULL)
	249	goto err;
	250
	251	if (in == NULL)
	252	R = BN_new();
	253	else
	254	R = in;
	255	if (R == NULL)
	256	goto err;
	257
	258	BN_one(X);
	259	BN_zero(Y);
	260	if (BN_copy(B, a) == NULL)
	261	goto err;
	262	if (BN_copy(A, n) == NULL)
	263	goto err;
	264	A->neg = 0;
	265	if (B->neg \|\| (BN_ucmp(B, A) >= 0)) {
	266	if (!BN_nnmod(B, B, A, ctx))
	267	goto err;
	268	}
	269	sign = -1;
	270	/*-
	271	* From B = a mod \|n\|, A = \|n\| it follows that
	272	*
	273	* 0 <= B < A,
	274	* -signXa == B (mod \|n\|),
	275	* signYa == A (mod \|n\|).
	276	*/
	277
	278	if (BN_is_odd(n) && (BN_num_bits(n) <= (BN_BITS <= 32 ? 450 : 2048))) {
	279	/*
	280	* Binary inversion algorithm; requires odd modulus. This is faster
	281	* than the general algorithm if the modulus is sufficiently small
	282	* (about 400 .. 500 bits on 32-bit sytems, but much more on 64-bit
	283	* systems)
	284	*/
	285	int shift;
	286
287	while (!BN_is_zero(B)) {
288	/*-
289	* 0 < B < \|n\|,
290	* 0 < A <= \|n\|,
291	* (1) -signXa == B (mod \|n\|),
292	* (2) signYa == A (mod \|n\|)
293	*/
294
295	/*
296	* Now divide B by the maximum possible power of two in the
297	* integers, and divide X by the same value mod \|n\|. When we're
298	* done, (1) still holds.
299	*/
300	shift = 0;
301	while (!BN_is_bit_set(B, shift)) { /* note that 0 < B */
302	shift++;
303
304	if (BN_is_odd(X)) {
305	if (!BN_uadd(X, X, n))
306	goto err;
307	}
308	/*
309	* now X is even, so we can easily divide it by two
310	*/
311	if (!BN_rshift1(X, X))
312	goto err;
313	}
314	if (shift > 0) {
315	if (!BN_rshift(B, B, shift))
316	goto err;
317	}
318
319	/*
320	* Same for A and Y. Afterwards, (2) still holds.
321	*/
322	shift = 0;
323	while (!BN_is_bit_set(A, shift)) { /* note that 0 < A */
324	shift++;
325
326	if (BN_is_odd(Y)) {
327	if (!BN_uadd(Y, Y, n))
328	goto err;
329	}
330	/* now Y is even */
331	if (!BN_rshift1(Y, Y))
332	goto err;
333	}
334	if (shift > 0) {
335	if (!BN_rshift(A, A, shift))
336	goto err;
337	}
338
339	/*-
340	* We still have (1) and (2).
341	* Both A and B are odd.
342	* The following computations ensure that
343	*
344	* 0 <= B < \|n\|,
345	* 0 < A < \|n\|,
346	* (1) -signXa == B (mod \|n\|),
347	* (2) signYa == A (mod \|n\|),
348	*
349	* and that either A or B is even in the next iteration.
350	*/
351	if (BN_ucmp(B, A) >= 0) {
352	/* -sign(X + Y)a == B - A (mod \|n\|) */
353	if (!BN_uadd(X, X, Y))
354	goto err;
355	/*
356	* NB: we could use BN_mod_add_quick(X, X, Y, n), but that
357	* actually makes the algorithm slower
358	*/
359	if (!BN_usub(B, B, A))
360	goto err;
361	} else {
362	/* sign(X + Y)a == A - B (mod \|n\|) */
363	if (!BN_uadd(Y, Y, X))
364	goto err;
365	/*
366	* as above, BN_mod_add_quick(Y, Y, X, n) would slow things
367	* down
368	*/
369	if (!BN_usub(A, A, B))
370	goto err;
371	}
372	}
373	} else {
374	/* general inversion algorithm */
375
376	while (!BN_is_zero(B)) {
377	BIGNUM *tmp;
378
379	/*-
380	* 0 < B < A,
381	* () -signX*a == B (mod \|n\|),
382	* signYa == A (mod \|n\|)
383	*/
384
385	/* (D, M) := (A/B, A%B) ... */
386	if (BN_num_bits(A) == BN_num_bits(B)) {
387	if (!BN_one(D))
388	goto err;
389	if (!BN_sub(M, A, B))
390	goto err;
391	} else if (BN_num_bits(A) == BN_num_bits(B) + 1) {
392	/* A/B is 1, 2, or 3 */
393	if (!BN_lshift1(T, B))
394	goto err;
395	if (BN_ucmp(A, T) < 0) {
396	/* A < 2B, so D=1 /
397	if (!BN_one(D))
398	goto err;
399	if (!BN_sub(M, A, B))
400	goto err;
401	} else {
402	/* A >= 2B, so D=2 or D=3 /
403	if (!BN_sub(M, A, T))
404	goto err;
405	if (!BN_add(D, T, B))
406	goto err; /* use D (:= 3B) as temp /
407	if (BN_ucmp(A, D) < 0) {
408	/* A < 3B, so D=2 /
409	if (!BN_set_word(D, 2))
410	goto err;
411	/*
412	* M (= A - 2*B) already has the correct value
413	*/
414	} else {
415	/* only D=3 remains */
416	if (!BN_set_word(D, 3))
417	goto err;
418	/*
419	* currently M = A - 2B, but we need M = A - 3B
420	*/
421	if (!BN_sub(M, M, B))
422	goto err;
423	}
424	}
425	} else {
426	if (!BN_div(D, M, A, B, ctx))
427	goto err;
428	}
429
430	/*-
431	* Now
432	* A = D*B + M;
433	* thus we have
434	* (*) signYa == DB + M (mod \|n\|).
435	*/
436
437	tmp = A; /* keep the BIGNUM object, the value does not
438	* matter */
439
440	/* (A, B) := (B, A mod B) ... */
441	A = B;
442	B = M;
443	/* ... so we have 0 <= B < A again */
444
445	/*-
446	* Since the former M is now B and the former B is now A,
447	* (**) translates into
448	* signYa == D*A + B (mod \|n\|),
449	* i.e.
450	* signYa - D*A == B (mod \|n\|).
451	* Similarly, (*) translates into
452	* -signXa == A (mod \|n\|).
453	*
454	* Thus,
455	* signYa + DsignX*a == B (mod \|n\|),
456	* i.e.
457	* sign(Y + DX)*a == B (mod \|n\|).
458	*
459	* So if we set (X, Y, sign) := (Y + D*X, X, -sign), we arrive back at
460	* -signXa == B (mod \|n\|),
461	* signYa == A (mod \|n\|).
462	* Note that X and Y stay non-negative all the time.
463	*/
464
465	/*
466	* most of the time D is very small, so we can optimize tmp :=
467	* D*X+Y
468	*/
469	if (BN_is_one(D)) {
470	if (!BN_add(tmp, X, Y))
471	goto err;
472	} else {
473	if (BN_is_word(D, 2)) {
474	if (!BN_lshift1(tmp, X))
475	goto err;
476	} else if (BN_is_word(D, 4)) {
477	if (!BN_lshift(tmp, X, 2))
478	goto err;
479	} else if (D->top == 1) {
480	if (!BN_copy(tmp, X))
481	goto err;
482	if (!BN_mul_word(tmp, D->d[0]))
483	goto err;
484	} else {
485	if (!BN_mul(tmp, D, X, ctx))
486	goto err;
487	}
488	if (!BN_add(tmp, tmp, Y))
489	goto err;
490	}
491
492	M = Y; /* keep the BIGNUM object, the value does not
493	* matter */
494	Y = X;
495	X = tmp;
496	sign = -sign;
497	}
498	}
499
500	/*-
501	* The while loop (Euclid's algorithm) ends when
502	* A == gcd(a,n);
503	* we have
504	* signYa == A (mod \|n\|),
505	* where Y is non-negative.
506	*/
507
508	if (sign < 0) {
509	if (!BN_sub(Y, n, Y))
510	goto err;
511	}
512	/* Now Ya == A (mod \|n\|). /
513
514	if (BN_is_one(A)) {
515	/* Ya == 1 (mod \|n\|) /
516	if (!Y->neg && BN_ucmp(Y, n) < 0) {
517	if (!BN_copy(R, Y))
518	goto err;
519	} else {
520	if (!BN_nnmod(R, Y, n, ctx))
521	goto err;
522	}
523	} else {
524	BNerr(BN_F_BN_MOD_INVERSE, BN_R_NO_INVERSE);
525	goto err;
526	}
527	ret = R;
528	err:
529	if ((ret == NULL) && (in == NULL))
530	BN_free(R);
531	BN_CTX_end(ctx);
532	bn_check_top(ret);
533	return (ret);
534	}
535
536	/*
537	* BN_mod_inverse_no_branch is a special version of BN_mod_inverse. It does
538	* not contain branches that may leak sensitive information.
7cdb8158	539	*/
283aedf4	540	static BIGNUM BN_mod_inverse_no_branch(BIGNUM in,
40720ce3 MC	541	const BIGNUM a, const BIGNUM n,
	542	BN_CTX *ctx)
	543	{
	544	BIGNUM A, B, X, Y, M, D, T, R = NULL;
	545	BIGNUM local_A, local_B;
	546	BIGNUM pA, pB;
	547	BIGNUM *ret = NULL;
	548	int sign;
	549
	550	bn_check_top(a);
	551	bn_check_top(n);
	552
	553	BN_CTX_start(ctx);
	554	A = BN_CTX_get(ctx);
	555	B = BN_CTX_get(ctx);
	556	X = BN_CTX_get(ctx);
	557	D = BN_CTX_get(ctx);
	558	M = BN_CTX_get(ctx);
	559	Y = BN_CTX_get(ctx);
	560	T = BN_CTX_get(ctx);
	561	if (T == NULL)
	562	goto err;
	563
	564	if (in == NULL)
	565	R = BN_new();
	566	else
	567	R = in;
	568	if (R == NULL)
	569	goto err;
	570
	571	BN_one(X);
	572	BN_zero(Y);
	573	if (BN_copy(B, a) == NULL)
	574	goto err;
	575	if (BN_copy(A, n) == NULL)
	576	goto err;
	577	A->neg = 0;
	578
	579	if (B->neg \|\| (BN_ucmp(B, A) >= 0)) {
	580	/*
	581	* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
	582	* BN_div_no_branch will be called eventually.
	583	*/
	584	pB = &local_B;
	585	BN_with_flags(pB, B, BN_FLG_CONSTTIME);
	586	if (!BN_nnmod(B, pB, A, ctx))
	587	goto err;
	588	}
	589	sign = -1;
	590	/*-
	591	* From B = a mod \|n\|, A = \|n\| it follows that
	592	*
	593	* 0 <= B < A,
	594	* -signXa == B (mod \|n\|),
	595	* signYa == A (mod \|n\|).
	596	*/
	597
	598	while (!BN_is_zero(B)) {
	599	BIGNUM *tmp;
	600
	601	/*-
	602	* 0 < B < A,
	603	* () -signX*a == B (mod \|n\|),
	604	* signYa == A (mod \|n\|)
605	*/
606
607	/*
608	* Turn BN_FLG_CONSTTIME flag on, so that when BN_div is invoked,
609	* BN_div_no_branch will be called eventually.
610	*/
611	pA = &local_A;
612	BN_with_flags(pA, A, BN_FLG_CONSTTIME);
613
614	/* (D, M) := (A/B, A%B) ... */
615	if (!BN_div(D, M, pA, B, ctx))
616	goto err;
617
618	/*-
619	* Now
620	* A = D*B + M;
621	* thus we have
622	* (*) signYa == DB + M (mod \|n\|).
623	*/
624
625	tmp = A; /* keep the BIGNUM object, the value does not
626	* matter */
627
628	/* (A, B) := (B, A mod B) ... */
629	A = B;
630	B = M;
631	/* ... so we have 0 <= B < A again */
632
633	/*-
634	* Since the former M is now B and the former B is now A,
635	* (**) translates into
636	* signYa == D*A + B (mod \|n\|),
637	* i.e.
638	* signYa - D*A == B (mod \|n\|).
639	* Similarly, (*) translates into
640	* -signXa == A (mod \|n\|).
641	*
642	* Thus,
643	* signYa + DsignX*a == B (mod \|n\|),
644	* i.e.
645	* sign(Y + DX)*a == B (mod \|n\|).
646	*
647	* So if we set (X, Y, sign) := (Y + D*X, X, -sign), we arrive back at
648	* -signXa == B (mod \|n\|),
649	* signYa == A (mod \|n\|).
650	* Note that X and Y stay non-negative all the time.
651	*/
652
653	if (!BN_mul(tmp, D, X, ctx))
654	goto err;
655	if (!BN_add(tmp, tmp, Y))
656	goto err;
657
658	M = Y; /* keep the BIGNUM object, the value does not
659	* matter */
660	Y = X;
661	X = tmp;
662	sign = -sign;
663	}
664
665	/*-
666	* The while loop (Euclid's algorithm) ends when
667	* A == gcd(a,n);
668	* we have
669	* signYa == A (mod \|n\|),
670	* where Y is non-negative.
671	*/
672
673	if (sign < 0) {
674	if (!BN_sub(Y, n, Y))
675	goto err;
676	}
677	/* Now Ya == A (mod \|n\|). /
678
679	if (BN_is_one(A)) {
680	/* Ya == 1 (mod \|n\|) /
681	if (!Y->neg && BN_ucmp(Y, n) < 0) {
682	if (!BN_copy(R, Y))
683	goto err;
684	} else {
685	if (!BN_nnmod(R, Y, n, ctx))
686	goto err;
687	}
688	} else {
689	BNerr(BN_F_BN_MOD_INVERSE_NO_BRANCH, BN_R_NO_INVERSE);
690	goto err;
691	}
692	ret = R;
693	err:
694	if ((ret == NULL) && (in == NULL))
695	BN_free(R);
696	BN_CTX_end(ctx);
697	bn_check_top(ret);
698	return (ret);
699	}