path: root/crypto/bn/rsaz_exp_x2.c

/*
 * Copyright 2020-2021 The OpenSSL Project Authors. All Rights Reserved.
 * Copyright (c) 2020, Intel Corporation. All Rights Reserved.
 *
 * Licensed under the Apache License 2.0 (the "License").  You may not use
 * this file except in compliance with the License.  You can obtain a copy
 * in the file LICENSE in the source distribution or at
 * https://www.openssl.org/source/license.html
 *
 *
 * Originally written by Ilya Albrekht, Sergey Kirillov and Andrey Matyukov
 * Intel Corporation
 *
 */

#include <openssl/opensslconf.h>
#include "rsaz_exp.h"

#ifndef RSAZ_ENABLED
NON_EMPTY_TRANSLATION_UNIT
#else
# include <assert.h>
# include <string.h>

# if defined(__GNUC__)
#  define ALIGN64 __attribute__((aligned(64)))
# elif defined(_MSC_VER)
#  define ALIGN64 __declspec(align(64))
# else
#  define ALIGN64
# endif

# define ALIGN_OF(ptr, boundary) \
    ((unsigned char *)(ptr) + (boundary - (((size_t)(ptr)) & (boundary - 1))))

/* Internal radix */
# define DIGIT_SIZE (52)
/* 52-bit mask */
# define DIGIT_MASK ((uint64_t)0xFFFFFFFFFFFFF)

# define BITS2WORD8_SIZE(x)  (((x) + 7) >> 3)
# define BITS2WORD64_SIZE(x) (((x) + 63) >> 6)

static ossl_inline uint64_t get_digit52(const uint8_t *in, int in_len);
static ossl_inline void put_digit52(uint8_t *out, int out_len, uint64_t digit);
static void to_words52(BN_ULONG *out, int out_len, const BN_ULONG *in,
                       int in_bitsize);
static void from_words52(BN_ULONG *bn_out, int out_bitsize, const BN_ULONG *in);
static ossl_inline void set_bit(BN_ULONG *a, int idx);

/* Number of |digit_size|-bit digits in |bitsize|-bit value */
static ossl_inline int number_of_digits(int bitsize, int digit_size)
{
    return (bitsize + digit_size - 1) / digit_size;
}

typedef void (*AMM52)(BN_ULONG *res, const BN_ULONG *base,
                      const BN_ULONG *exp, const BN_ULONG *m, BN_ULONG k0);
typedef void (*EXP52_x2)(BN_ULONG *res, const BN_ULONG *base,
                         const BN_ULONG *exp[2], const BN_ULONG *m,
                         const BN_ULONG *rr, const BN_ULONG k0[2]);

/*
 * For details of the methods declared below please refer to
 *    crypto/bn/asm/rsaz-avx512.pl
 *
 * Naming notes:
 *  amm = Almost Montgomery Multiplication
 *  ams = Almost Montgomery Squaring
 *  52x20 - data represented as array of 20 digits in 52-bit radix
 *  _x1_/_x2_ - 1 or 2 independent inputs/outputs
 *  _256 suffix - uses 256-bit (AVX512VL) registers
 */

/*AMM = Almost Montgomery Multiplication. */
void RSAZ_amm52x20_x1_256(BN_ULONG *res, const BN_ULONG *base,
                          const BN_ULONG *exp, const BN_ULONG *m,
                          BN_ULONG k0);
void RSAZ_exp52x20_x2_256(BN_ULONG *res, const BN_ULONG *base,
                      const BN_ULONG *exp[2], const BN_ULONG *m,
                      const BN_ULONG *rr, const BN_ULONG k0[2]);
void RSAZ_amm52x20_x2_256(BN_ULONG *out, const BN_ULONG *a,
                          const BN_ULONG *b, const BN_ULONG *m,
                          const BN_ULONG k0[2]);
void ossl_extract_multiplier_2x20_win5(BN_ULONG *red_Y,
                                       const BN_ULONG *red_table,
                                       int red_table_idx, int tbl_idx);

/*
 * Dual Montgomery modular exponentiation using prime moduli of the
 * same bit size, optimized with AVX512 ISA.
 *
 * Input and output parameters for each exponentiation are independent and
 * denoted here by index |i|, i = 1..2.
 *
 * Input and output are all in regular 2^64 radix.
 *
 * Each moduli shall be |factor_size| bit size.
 *
 * NOTE: currently only 2x1024 case is supported.
 *
 *  [out] res|i|      - result of modular exponentiation: array of qword values
 *                      in regular (2^64) radix. Size of array shall be enough
 *                      to hold |factor_size| bits.
 *  [in]  base|i|     - base
 *  [in]  exp|i|      - exponent
 *  [in]  m|i|        - moduli
 *  [in]  rr|i|       - Montgomery parameter RR = R^2 mod m|i|
 *  [in]  k0_|i|      - Montgomery parameter k0 = -1/m|i| mod 2^64
 *  [in]  factor_size - moduli bit size
 *
 * \return 0 in case of failure,
 *         1 in case of success.
 */
int RSAZ_mod_exp_avx512_x2(BN_ULONG *res1,
                           const BN_ULONG *base1,
                           const BN_ULONG *exp1,
                           const BN_ULONG *m1,
                           const BN_ULONG *rr1,
                           BN_ULONG k0_1,
                           BN_ULONG *res2,
                           const BN_ULONG *base2,
                           const BN_ULONG *exp2,
                           const BN_ULONG *m2,
                           const BN_ULONG *rr2,
                           BN_ULONG k0_2,
                           int factor_size)
{
    int ret = 0;

    /*
     * Number of word-size (BN_ULONG) digits to store exponent in redundant
     * representation.
     */
    int exp_digits = number_of_digits(factor_size + 2, DIGIT_SIZE);
    int coeff_pow = 4 * (DIGIT_SIZE * exp_digits - factor_size);
    BN_ULONG *base1_red, *m1_red, *rr1_red;
    BN_ULONG *base2_red, *m2_red, *rr2_red;
    BN_ULONG *coeff_red;
    BN_ULONG *storage = NULL;
    BN_ULONG *storage_aligned = NULL;
    BN_ULONG storage_len_bytes = 7 * exp_digits * sizeof(BN_ULONG);

    /* AMM = Almost Montgomery Multiplication */
    AMM52 amm = NULL;
    /* Dual (2-exps in parallel) exponentiation */
    EXP52_x2 exp_x2 = NULL;

    const BN_ULONG *exp[2] = {0};
    BN_ULONG k0[2] = {0};

    /* Only 1024-bit factor size is supported now */
    switch (factor_size) {
    case 1024:
        amm = RSAZ_amm52x20_x1_256;
        exp_x2 = RSAZ_exp52x20_x2_256;
        break;
    default:
        goto err;
    }

    storage = (BN_ULONG *)OPENSSL_malloc(storage_len_bytes + 64);
    if (storage == NULL)
        goto err;
    storage_aligned = (BN_ULONG *)ALIGN_OF(storage, 64);

    /* Memory layout for red(undant) representations */
    base1_red = storage_aligned;
    base2_red = storage_aligned + 1 * exp_digits;
    m1_red    = storage_aligned + 2 * exp_digits;
    m2_red    = storage_aligned + 3 * exp_digits;
    rr1_red   = storage_aligned + 4 * exp_digits;
    rr2_red   = storage_aligned + 5 * exp_digits;
    coeff_red = storage_aligned + 6 * exp_digits;

    /* Convert base_i, m_i, rr_i, from regular to 52-bit radix */
    to_words52(base1_red, exp_digits, base1, factor_size);
    to_words52(base2_red, exp_digits, base2, factor_size);
    to_words52(m1_red, exp_digits, m1, factor_size);
    to_words52(m2_red, exp_digits, m2, factor_size);
    to_words52(rr1_red, exp_digits, rr1, factor_size);
    to_words52(rr2_red, exp_digits, rr2, factor_size);

    /*
     * Compute target domain Montgomery converters RR' for each modulus
     * based on precomputed original domain's RR.
     *
     * RR -> RR' transformation steps:
     *  (1) coeff = 2^k
     *  (2) t = AMM(RR,RR) = RR^2 / R' mod m
     *  (3) RR' = AMM(t, coeff) = RR^2 * 2^k / R'^2 mod m
     * where
     *  k = 4 * (52 * digits52 - modlen)
     *  R  = 2^(64 * ceil(modlen/64)) mod m
     *  RR = R^2 mod M
     *  R' = 2^(52 * ceil(modlen/52)) mod m
     *
     *  modlen = 1024: k = 64, RR = 2^2048 mod m, RR' = 2^2080 mod m
     */
    memset(coeff_red, 0, exp_digits * sizeof(BN_ULONG));
    /* (1) in reduced domain representation */
    set_bit(coeff_red, 64 * (int)(coeff_pow / 52) + coeff_pow % 52);

    amm(rr1_red, rr1_red, rr1_red, m1_red, k0_1);     /* (2) for m1 */
    amm(rr1_red, rr1_red, coeff_red, m1_red, k0_1);   /* (3) for m1 */

    amm(rr2_red, rr2_red, rr2_red, m2_red, k0_2);     /* (2) for m2 */
    amm(rr2_red, rr2_red, coeff_red, m2_red, k0_2);   /* (3) for m2 */

    exp[0] = exp1;
    exp[1] = exp2;

    k0[0] = k0_1;
    k0[1] = k0_2;

    exp_x2(rr1_red, base1_red, exp, m1_red, rr1_red, k0)