1/* mpi-mod.c -  Modular reduction
2   Copyright (C) 1998, 1999, 2001, 2002, 2003,
3                 2007  Free Software Foundation, Inc.
4
5   This file is part of Libgcrypt.
6
7   Libgcrypt is free software; you can redistribute it and/or modify
8   it under the terms of the GNU Lesser General Public License as
9   published by the Free Software Foundation; either version 2.1 of
10   the License, or (at your option) any later version.
11
12   Libgcrypt is distributed in the hope that it will be useful,
13   but WITHOUT ANY WARRANTY; without even the implied warranty of
14   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15   GNU Lesser General Public License for more details.
16
17   You should have received a copy of the GNU Lesser General Public
18   License along with this program; if not, write to the Free Software
19   Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301,
20   USA.  */
21
22
23#include <config.h>
24#include <stdio.h>
25#include <stdlib.h>
26
27#include "mpi-internal.h"
28#include "longlong.h"
29#include "g10lib.h"
30
31
32/* Context used with Barrett reduction.  */
33struct barrett_ctx_s
34{
35  gcry_mpi_t m;   /* The modulus - may not be modified. */
36  int m_copied;   /* If true, M needs to be released.  */
37  int k;
38  gcry_mpi_t y;
39  gcry_mpi_t r1;  /* Helper MPI. */
40  gcry_mpi_t r2;  /* Helper MPI. */
41  gcry_mpi_t r3;  /* Helper MPI allocated on demand. */
42};
43
44
45
46void
47_gcry_mpi_mod (gcry_mpi_t rem, gcry_mpi_t dividend, gcry_mpi_t divisor)
48{
49  _gcry_mpi_fdiv_r (rem, dividend, divisor);
50  rem->sign = 0;
51}
52
53
54/* This function returns a new context for Barrett based operations on
55   the modulus M.  This context needs to be released using
56   _gcry_mpi_barrett_free.  If COPY is true M will be transferred to
57   the context and the user may change M.  If COPY is false, M may not
58   be changed until gcry_mpi_barrett_free has been called. */
59mpi_barrett_t
60_gcry_mpi_barrett_init (gcry_mpi_t m, int copy)
61{
62  mpi_barrett_t ctx;
63  gcry_mpi_t tmp;
64
65  mpi_normalize (m);
66  ctx = gcry_xcalloc (1, sizeof *ctx);
67
68  if (copy)
69    {
70      ctx->m = mpi_copy (m);
71      ctx->m_copied = 1;
72    }
73  else
74    ctx->m = m;
75
76  ctx->k = mpi_get_nlimbs (m);
77  tmp = mpi_alloc (ctx->k + 1);
78
79  /* Barrett precalculation: y = floor(b^(2k) / m). */
80  mpi_set_ui (tmp, 1);
81  mpi_lshift_limbs (tmp, 2 * ctx->k);
82  mpi_fdiv_q (tmp, tmp, m);
83
84  ctx->y  = tmp;
85  ctx->r1 = mpi_alloc ( 2 * ctx->k + 1 );
86  ctx->r2 = mpi_alloc ( 2 * ctx->k + 1 );
87
88  return ctx;
89}
90
91void
92_gcry_mpi_barrett_free (mpi_barrett_t ctx)
93{
94  if (ctx)
95    {
96      mpi_free (ctx->y);
97      mpi_free (ctx->r1);
98      mpi_free (ctx->r2);
99      if (ctx->r3)
100        mpi_free (ctx->r3);
101      if (ctx->m_copied)
102        mpi_free (ctx->m);
103      gcry_free (ctx);
104    }
105}
106
107
108/* R = X mod M
109
110   Using Barrett reduction.  Before using this function
111   _gcry_mpi_barrett_init must have been called to do the
112   precalculations.  CTX is the context created by this precalculation
113   and also conveys M.  If the Barret reduction could no be done a
114   starightforward reduction method is used.
115
116   We assume that these conditions are met:
117   Input:  x =(x_2k-1 ...x_0)_b
118 	   m =(m_k-1 ....m_0)_b	  with m_k-1 != 0
119   Output: r = x mod m
120 */
121void
122_gcry_mpi_mod_barrett (gcry_mpi_t r, gcry_mpi_t x, mpi_barrett_t ctx)
123{
124  gcry_mpi_t m = ctx->m;
125  int k = ctx->k;
126  gcry_mpi_t y = ctx->y;
127  gcry_mpi_t r1 = ctx->r1;
128  gcry_mpi_t r2 = ctx->r2;
129
130  mpi_normalize (x);
131  if (mpi_get_nlimbs (x) > 2*k )
132    {
133      mpi_mod (r, x, m);
134      return;
135    }
136
137  /* 1. q1 = floor( x / b^k-1)
138   *    q2 = q1 * y
139   *    q3 = floor( q2 / b^k+1 )
140   * Actually, we don't need qx, we can work direct on r2
141   */
142  mpi_set ( r2, x );
143  mpi_rshift_limbs ( r2, k-1 );
144  mpi_mul ( r2, r2, y );
145  mpi_rshift_limbs ( r2, k+1 );
146
147  /* 2. r1 = x mod b^k+1
148   *	r2 = q3 * m mod b^k+1
149   *	r  = r1 - r2
150   * 3. if r < 0 then  r = r + b^k+1
151   */
152  mpi_set ( r1, x );
153  if ( r1->nlimbs > k+1 ) /* Quick modulo operation.  */
154    r1->nlimbs = k+1;
155  mpi_mul ( r2, r2, m );
156  if ( r2->nlimbs > k+1 ) /* Quick modulo operation. */
157    r2->nlimbs = k+1;
158  mpi_sub ( r, r1, r2 );
159
160  if ( mpi_is_neg( r ) )
161    {
162      if (!ctx->r3)
163        {
164          ctx->r3 = mpi_alloc ( k + 2 );
165          mpi_set_ui (ctx->r3, 1);
166          mpi_lshift_limbs (ctx->r3, k + 1 );
167        }
168      mpi_add ( r, r, ctx->r3 );
169    }
170
171  /* 4. while r >= m do r = r - m */
172  while ( mpi_cmp( r, m ) >= 0 )
173    mpi_sub ( r, r, m );
174
175}
176
177
178void
179_gcry_mpi_mul_barrett (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v,
180                       mpi_barrett_t ctx)
181{
182  gcry_mpi_mul (w, u, v);
183  mpi_mod_barrett (w, w, ctx);
184}
185