See ChangeLog: Mon Jan 24 22:24:38 CET 2000 Werner Koch
[gnupg.git] / mpi / mpi-pow.c
1 /* mpi-pow.c  -  MPI functions
2  *      Copyright (C) 1998 Free Software Foundation, Inc.
3  *      Copyright (C) 1994, 1996 Free Software Foundation, Inc.
4  *
5  * This file is part of GnuPG.
6  *
7  * GnuPG is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 2 of the License, or
10  * (at your option) any later version.
11  *
12  * GnuPG 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 General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software
19  * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20  *
21  * Note: This code is heavily based on the GNU MP Library.
22  *       Actually it's the same code with only minor changes in the
23  *       way the data is stored; this is to support the abstraction
24  *       of an optional secure memory allocation which may be used
25  *       to avoid revealing of sensitive data due to paging etc.
26  *       The GNU MP Library itself is published under the LGPL;
27  *       however I decided to publish this code under the plain GPL.
28  */
29
30 #include <config.h>
31 #include <stdio.h>
32 #include <stdlib.h>
33 #include <assert.h>
34 #include "mpi-internal.h"
35 #include "longlong.h"
36
37
38 /****************
39  * RES = BASE ^ EXP mod MOD
40  */
41 void
42 gcry_mpi_powm( MPI res, MPI base, MPI exp, MPI mod)
43 {
44     mpi_ptr_t  rp, ep, mp, bp;
45     mpi_size_t esize, msize, bsize, rsize;
46     int        esign, msign, bsign, rsign;
47     int        esec,  msec,  bsec,  rsec;
48     mpi_size_t size;
49     int mod_shift_cnt;
50     int negative_result;
51     mpi_ptr_t mp_marker=NULL, bp_marker=NULL, ep_marker=NULL;
52     mpi_ptr_t xp_marker=NULL;
53     int assign_rp=0;
54     mpi_ptr_t tspace = NULL;
55     mpi_size_t tsize=0;   /* to avoid compiler warning */
56                           /* fixme: we should check that the warning is void*/
57
58     esize = exp->nlimbs;
59     msize = mod->nlimbs;
60     size = 2 * msize;
61     esign = exp->sign;
62     msign = mod->sign;
63
64     esec = mpi_is_secure(exp);
65     msec = mpi_is_secure(mod);
66     bsec = mpi_is_secure(base);
67     rsec = mpi_is_secure(res);
68
69     rp = res->d;
70     ep = exp->d;
71
72     if( !msize )
73         msize = 1 / msize;          /* provoke a signal */
74
75     if( !esize ) {
76         /* Exponent is zero, result is 1 mod MOD, i.e., 1 or 0
77          * depending on if MOD equals 1.  */
78         rp[0] = 1;
79         res->nlimbs = (msize == 1 && mod->d[0] == 1) ? 0 : 1;
80         res->sign = 0;
81         goto leave;
82     }
83
84     /* Normalize MOD (i.e. make its most significant bit set) as required by
85      * mpn_divrem.  This will make the intermediate values in the calculation
86      * slightly larger, but the correct result is obtained after a final
87      * reduction using the original MOD value.  */
88     mp = mp_marker = mpi_alloc_limb_space(msize, msec);
89     count_leading_zeros( mod_shift_cnt, mod->d[msize-1] );
90     if( mod_shift_cnt )
91         mpihelp_lshift( mp, mod->d, msize, mod_shift_cnt );
92     else
93         MPN_COPY( mp, mod->d, msize );
94
95     bsize = base->nlimbs;
96     bsign = base->sign;
97     if( bsize > msize ) { /* The base is larger than the module. Reduce it. */
98         /* Allocate (BSIZE + 1) with space for remainder and quotient.
99          * (The quotient is (bsize - msize + 1) limbs.)  */
100         bp = bp_marker = mpi_alloc_limb_space( bsize + 1, bsec );
101         MPN_COPY( bp, base->d, bsize );
102         /* We don't care about the quotient, store it above the remainder,
103          * at BP + MSIZE.  */
104         mpihelp_divrem( bp + msize, 0, bp, bsize, mp, msize );
105         bsize = msize;
106         /* Canonicalize the base, since we are going to multiply with it
107          * quite a few times.  */
108         MPN_NORMALIZE( bp, bsize );
109     }
110     else
111         bp = base->d;
112
113     if( !bsize ) {
114         res->nlimbs = 0;
115         res->sign = 0;
116         goto leave;
117     }
118
119     if( res->alloced < size ) {
120         /* We have to allocate more space for RES.  If any of the input
121          * parameters are identical to RES, defer deallocation of the old
122          * space.  */
123         if( rp == ep || rp == mp || rp == bp ) {
124             rp = mpi_alloc_limb_space( size, rsec );
125             assign_rp = 1;
126         }
127         else {
128             mpi_resize( res, size );
129             rp = res->d;
130         }
131     }
132     else { /* Make BASE, EXP and MOD not overlap with RES.  */
133         if( rp == bp ) {
134             /* RES and BASE are identical.  Allocate temp. space for BASE.  */
135             assert( !bp_marker );
136             bp = bp_marker = mpi_alloc_limb_space( bsize, bsec );
137             MPN_COPY(bp, rp, bsize);
138         }
139         if( rp == ep ) {
140             /* RES and EXP are identical.  Allocate temp. space for EXP.  */
141             ep = ep_marker = mpi_alloc_limb_space( esize, esec );
142             MPN_COPY(ep, rp, esize);
143         }
144         if( rp == mp ) {
145             /* RES and MOD are identical.  Allocate temporary space for MOD.*/
146             assert( !mp_marker );
147             mp = mp_marker = mpi_alloc_limb_space( msize, msec );
148             MPN_COPY(mp, rp, msize);
149         }
150     }
151
152     MPN_COPY( rp, bp, bsize );
153     rsize = bsize;
154     rsign = bsign;
155
156     {
157         mpi_size_t i;
158         mpi_ptr_t xp = xp_marker = mpi_alloc_limb_space( 2 * (msize + 1), msec );
159         int c;
160         mpi_limb_t e;
161         mpi_limb_t carry_limb;
162
163         negative_result = (ep[0] & 1) && base->sign;
164
165         i = esize - 1;
166         e = ep[i];
167         count_leading_zeros (c, e);
168         e = (e << c) << 1;     /* shift the exp bits to the left, lose msb */
169         c = BITS_PER_MPI_LIMB - 1 - c;
170
171         /* Main loop.
172          *
173          * Make the result be pointed to alternately by XP and RP.  This
174          * helps us avoid block copying, which would otherwise be necessary
175          * with the overlap restrictions of mpihelp_divmod. With 50% probability
176          * the result after this loop will be in the area originally pointed
177          * by RP (==RES->d), and with 50% probability in the area originally
178          * pointed to by XP.
179          */
180         for(;;) {
181             while( c ) {
182                 mpi_ptr_t tp;
183                 mpi_size_t xsize;
184
185                 /*mpihelp_mul_n(xp, rp, rp, rsize);*/
186                 if( rsize < KARATSUBA_THRESHOLD )
187                     mpih_sqr_n_basecase( xp, rp, rsize );
188                 else {
189                     if( !tspace ) {
190                         tsize = 2 * rsize;
191                         tspace = mpi_alloc_limb_space( tsize, 0 );
192                     }
193                     else if( tsize < (2*rsize) ) {
194                         mpi_free_limb_space( tspace );
195                         tsize = 2 * rsize;
196                         tspace = mpi_alloc_limb_space( tsize, 0 );
197
198                     }
199                     mpih_sqr_n( xp, rp, rsize, tspace );
200                 }
201
202                 xsize = 2 * rsize;
203                 if( xsize > msize ) {
204                     mpihelp_divrem(xp + msize, 0, xp, xsize, mp, msize);
205                     xsize = msize;
206                 }
207
208                 tp = rp; rp = xp; xp = tp;
209                 rsize = xsize;
210
211                 if( (mpi_limb_signed_t)e < 0 ) {
212                     mpihelp_mul( xp, rp, rsize, bp, bsize );
213                     xsize = rsize + bsize;
214                     if( xsize > msize ) {
215                         mpihelp_divrem(xp + msize, 0, xp, xsize, mp, msize);
216                         xsize = msize;
217                     }
218
219                     tp = rp; rp = xp; xp = tp;
220                     rsize = xsize;
221                 }
222                 e <<= 1;
223                 c--;
224             }
225
226             i--;
227             if( i < 0 )
228                 break;
229             e = ep[i];
230             c = BITS_PER_MPI_LIMB;
231         }
232
233         /* We shifted MOD, the modulo reduction argument, left MOD_SHIFT_CNT
234          * steps.  Adjust the result by reducing it with the original MOD.
235          *
236          * Also make sure the result is put in RES->d (where it already
237          * might be, see above).
238          */
239         if( mod_shift_cnt ) {
240             carry_limb = mpihelp_lshift( res->d, rp, rsize, mod_shift_cnt);
241             rp = res->d;
242             if( carry_limb ) {
243                 rp[rsize] = carry_limb;
244                 rsize++;
245             }
246         }
247         else {
248             MPN_COPY( res->d, rp, rsize);
249             rp = res->d;
250         }
251
252         if( rsize >= msize ) {
253             mpihelp_divrem(rp + msize, 0, rp, rsize, mp, msize);
254             rsize = msize;
255         }
256
257         /* Remove any leading zero words from the result.  */
258         if( mod_shift_cnt )
259             mpihelp_rshift( rp, rp, rsize, mod_shift_cnt);
260         MPN_NORMALIZE (rp, rsize);
261     }
262
263     if( negative_result && rsize ) {
264         if( mod_shift_cnt )
265             mpihelp_rshift( mp, mp, msize, mod_shift_cnt);
266         mpihelp_sub( rp, mp, msize, rp, rsize);
267         rsize = msize;
268         rsign = msign;
269         MPN_NORMALIZE(rp, rsize);
270     }
271     res->nlimbs = rsize;
272     res->sign = rsign;
273
274   leave:
275     if( assign_rp ) mpi_assign_limb_space( res, rp, size );
276     if( mp_marker ) mpi_free_limb_space( mp_marker );
277     if( bp_marker ) mpi_free_limb_space( bp_marker );
278     if( ep_marker ) mpi_free_limb_space( ep_marker );
279     if( xp_marker ) mpi_free_limb_space( xp_marker );
280     if( tspace )    mpi_free_limb_space( tspace );
281 }
282