xref: /openbsd/lib/libcrypto/bn/bn_mul.c (revision 913ec974)
15b37fcf3Sryker /* crypto/bn/bn_mul.c */
25b37fcf3Sryker /* Copyright (C) 1995-1998 Eric Young (eay@cryptsoft.com)
35b37fcf3Sryker  * All rights reserved.
45b37fcf3Sryker  *
55b37fcf3Sryker  * This package is an SSL implementation written
65b37fcf3Sryker  * by Eric Young (eay@cryptsoft.com).
75b37fcf3Sryker  * The implementation was written so as to conform with Netscapes SSL.
85b37fcf3Sryker  *
95b37fcf3Sryker  * This library is free for commercial and non-commercial use as long as
105b37fcf3Sryker  * the following conditions are aheared to.  The following conditions
115b37fcf3Sryker  * apply to all code found in this distribution, be it the RC4, RSA,
125b37fcf3Sryker  * lhash, DES, etc., code; not just the SSL code.  The SSL documentation
135b37fcf3Sryker  * included with this distribution is covered by the same copyright terms
145b37fcf3Sryker  * except that the holder is Tim Hudson (tjh@cryptsoft.com).
155b37fcf3Sryker  *
165b37fcf3Sryker  * Copyright remains Eric Young's, and as such any Copyright notices in
175b37fcf3Sryker  * the code are not to be removed.
185b37fcf3Sryker  * If this package is used in a product, Eric Young should be given attribution
195b37fcf3Sryker  * as the author of the parts of the library used.
205b37fcf3Sryker  * This can be in the form of a textual message at program startup or
215b37fcf3Sryker  * in documentation (online or textual) provided with the package.
225b37fcf3Sryker  *
235b37fcf3Sryker  * Redistribution and use in source and binary forms, with or without
245b37fcf3Sryker  * modification, are permitted provided that the following conditions
255b37fcf3Sryker  * are met:
265b37fcf3Sryker  * 1. Redistributions of source code must retain the copyright
275b37fcf3Sryker  *    notice, this list of conditions and the following disclaimer.
285b37fcf3Sryker  * 2. Redistributions in binary form must reproduce the above copyright
295b37fcf3Sryker  *    notice, this list of conditions and the following disclaimer in the
305b37fcf3Sryker  *    documentation and/or other materials provided with the distribution.
315b37fcf3Sryker  * 3. All advertising materials mentioning features or use of this software
325b37fcf3Sryker  *    must display the following acknowledgement:
335b37fcf3Sryker  *    "This product includes cryptographic software written by
345b37fcf3Sryker  *     Eric Young (eay@cryptsoft.com)"
355b37fcf3Sryker  *    The word 'cryptographic' can be left out if the rouines from the library
365b37fcf3Sryker  *    being used are not cryptographic related :-).
375b37fcf3Sryker  * 4. If you include any Windows specific code (or a derivative thereof) from
385b37fcf3Sryker  *    the apps directory (application code) you must include an acknowledgement:
395b37fcf3Sryker  *    "This product includes software written by Tim Hudson (tjh@cryptsoft.com)"
405b37fcf3Sryker  *
415b37fcf3Sryker  * THIS SOFTWARE IS PROVIDED BY ERIC YOUNG ``AS IS'' AND
425b37fcf3Sryker  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
435b37fcf3Sryker  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
445b37fcf3Sryker  * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
455b37fcf3Sryker  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
465b37fcf3Sryker  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
475b37fcf3Sryker  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
485b37fcf3Sryker  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
495b37fcf3Sryker  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
505b37fcf3Sryker  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
515b37fcf3Sryker  * SUCH DAMAGE.
525b37fcf3Sryker  *
535b37fcf3Sryker  * The licence and distribution terms for any publically available version or
545b37fcf3Sryker  * derivative of this code cannot be changed.  i.e. this code cannot simply be
555b37fcf3Sryker  * copied and put under another distribution licence
565b37fcf3Sryker  * [including the GNU Public Licence.]
575b37fcf3Sryker  */
585b37fcf3Sryker 
595b37fcf3Sryker #include <stdio.h>
605b37fcf3Sryker #include "cryptlib.h"
615b37fcf3Sryker #include "bn_lcl.h"
625b37fcf3Sryker 
63*913ec974Sbeck #ifdef BN_RECURSION
64*913ec974Sbeck /* r is 2*n2 words in size,
65*913ec974Sbeck  * a and b are both n2 words in size.
66*913ec974Sbeck  * n2 must be a power of 2.
67*913ec974Sbeck  * We multiply and return the result.
68*913ec974Sbeck  * t must be 2*n2 words in size
69*913ec974Sbeck  * We calulate
70*913ec974Sbeck  * a[0]*b[0]
71*913ec974Sbeck  * a[0]*b[0]+a[1]*b[1]+(a[0]-a[1])*(b[1]-b[0])
72*913ec974Sbeck  * a[1]*b[1]
73*913ec974Sbeck  */
74*913ec974Sbeck void bn_mul_recursive(BN_ULONG *r, BN_ULONG *a, BN_ULONG *b, int n2,
75*913ec974Sbeck 	     BN_ULONG *t)
765b37fcf3Sryker 	{
77*913ec974Sbeck 	int n=n2/2,c1,c2;
78*913ec974Sbeck 	unsigned int neg,zero;
79*913ec974Sbeck 	BN_ULONG ln,lo,*p;
805b37fcf3Sryker 
81*913ec974Sbeck #ifdef BN_COUNT
82*913ec974Sbeck printf(" bn_mul_recursive %d * %d\n",n2,n2);
83*913ec974Sbeck #endif
84*913ec974Sbeck #ifdef BN_MUL_COMBA
85*913ec974Sbeck /*	if (n2 == 4)
865b37fcf3Sryker 		{
87*913ec974Sbeck 		bn_mul_comba4(r,a,b);
88*913ec974Sbeck 		return;
89*913ec974Sbeck 		}
90*913ec974Sbeck 	else */ if (n2 == 8)
91*913ec974Sbeck 		{
92*913ec974Sbeck 		bn_mul_comba8(r,a,b);
93*913ec974Sbeck 		return;
94*913ec974Sbeck 		}
95*913ec974Sbeck #endif
96*913ec974Sbeck 	if (n2 < BN_MUL_RECURSIVE_SIZE_NORMAL)
97*913ec974Sbeck 		{
98*913ec974Sbeck 		/* This should not happen */
99*913ec974Sbeck 		bn_mul_normal(r,a,n2,b,n2);
100*913ec974Sbeck 		return;
101*913ec974Sbeck 		}
102*913ec974Sbeck 	/* r=(a[0]-a[1])*(b[1]-b[0]) */
103*913ec974Sbeck 	c1=bn_cmp_words(a,&(a[n]),n);
104*913ec974Sbeck 	c2=bn_cmp_words(&(b[n]),b,n);
105*913ec974Sbeck 	zero=neg=0;
106*913ec974Sbeck 	switch (c1*3+c2)
107*913ec974Sbeck 		{
108*913ec974Sbeck 	case -4:
109*913ec974Sbeck 		bn_sub_words(t,      &(a[n]),a,      n); /* - */
110*913ec974Sbeck 		bn_sub_words(&(t[n]),b,      &(b[n]),n); /* - */
111*913ec974Sbeck 		break;
112*913ec974Sbeck 	case -3:
113*913ec974Sbeck 		zero=1;
114*913ec974Sbeck 		break;
115*913ec974Sbeck 	case -2:
116*913ec974Sbeck 		bn_sub_words(t,      &(a[n]),a,      n); /* - */
117*913ec974Sbeck 		bn_sub_words(&(t[n]),&(b[n]),b,      n); /* + */
118*913ec974Sbeck 		neg=1;
119*913ec974Sbeck 		break;
120*913ec974Sbeck 	case -1:
121*913ec974Sbeck 	case 0:
122*913ec974Sbeck 	case 1:
123*913ec974Sbeck 		zero=1;
124*913ec974Sbeck 		break;
125*913ec974Sbeck 	case 2:
126*913ec974Sbeck 		bn_sub_words(t,      a,      &(a[n]),n); /* + */
127*913ec974Sbeck 		bn_sub_words(&(t[n]),b,      &(b[n]),n); /* - */
128*913ec974Sbeck 		neg=1;
129*913ec974Sbeck 		break;
130*913ec974Sbeck 	case 3:
131*913ec974Sbeck 		zero=1;
132*913ec974Sbeck 		break;
133*913ec974Sbeck 	case 4:
134*913ec974Sbeck 		bn_sub_words(t,      a,      &(a[n]),n);
135*913ec974Sbeck 		bn_sub_words(&(t[n]),&(b[n]),b,      n);
136*913ec974Sbeck 		break;
1375b37fcf3Sryker 		}
1385b37fcf3Sryker 
139*913ec974Sbeck #ifdef BN_MUL_COMBA
140*913ec974Sbeck 	if (n == 4)
1415b37fcf3Sryker 		{
142*913ec974Sbeck 		if (!zero)
143*913ec974Sbeck 			bn_mul_comba4(&(t[n2]),t,&(t[n]));
144*913ec974Sbeck 		else
145*913ec974Sbeck 			memset(&(t[n2]),0,8*sizeof(BN_ULONG));
146*913ec974Sbeck 
147*913ec974Sbeck 		bn_mul_comba4(r,a,b);
148*913ec974Sbeck 		bn_mul_comba4(&(r[n2]),&(a[n]),&(b[n]));
1495b37fcf3Sryker 		}
150*913ec974Sbeck 	else if (n == 8)
151*913ec974Sbeck 		{
152*913ec974Sbeck 		if (!zero)
153*913ec974Sbeck 			bn_mul_comba8(&(t[n2]),t,&(t[n]));
154*913ec974Sbeck 		else
155*913ec974Sbeck 			memset(&(t[n2]),0,16*sizeof(BN_ULONG));
156*913ec974Sbeck 
157*913ec974Sbeck 		bn_mul_comba8(r,a,b);
158*913ec974Sbeck 		bn_mul_comba8(&(r[n2]),&(a[n]),&(b[n]));
159*913ec974Sbeck 		}
160*913ec974Sbeck 	else
161*913ec974Sbeck #endif
162*913ec974Sbeck 		{
163*913ec974Sbeck 		p= &(t[n2*2]);
164*913ec974Sbeck 		if (!zero)
165*913ec974Sbeck 			bn_mul_recursive(&(t[n2]),t,&(t[n]),n,p);
166*913ec974Sbeck 		else
167*913ec974Sbeck 			memset(&(t[n2]),0,n2*sizeof(BN_ULONG));
168*913ec974Sbeck 		bn_mul_recursive(r,a,b,n,p);
169*913ec974Sbeck 		bn_mul_recursive(&(r[n2]),&(a[n]),&(b[n]),n,p);
1705b37fcf3Sryker 		}
1715b37fcf3Sryker 
172*913ec974Sbeck 	/* t[32] holds (a[0]-a[1])*(b[1]-b[0]), c1 is the sign
173*913ec974Sbeck 	 * r[10] holds (a[0]*b[0])
174*913ec974Sbeck 	 * r[32] holds (b[1]*b[1])
175*913ec974Sbeck 	 */
1765b37fcf3Sryker 
177*913ec974Sbeck 	c1=(int)(bn_add_words(t,r,&(r[n2]),n2));
1785b37fcf3Sryker 
179*913ec974Sbeck 	if (neg) /* if t[32] is negative */
1805b37fcf3Sryker 		{
181*913ec974Sbeck 		c1-=(int)(bn_sub_words(&(t[n2]),t,&(t[n2]),n2));
1825b37fcf3Sryker 		}
1835b37fcf3Sryker 	else
1845b37fcf3Sryker 		{
185*913ec974Sbeck 		/* Might have a carry */
186*913ec974Sbeck 		c1+=(int)(bn_add_words(&(t[n2]),&(t[n2]),t,n2));
1875b37fcf3Sryker 		}
1885b37fcf3Sryker 
189*913ec974Sbeck 	/* t[32] holds (a[0]-a[1])*(b[1]-b[0])+(a[0]*b[0])+(a[1]*b[1])
190*913ec974Sbeck 	 * r[10] holds (a[0]*b[0])
191*913ec974Sbeck 	 * r[32] holds (b[1]*b[1])
192*913ec974Sbeck 	 * c1 holds the carry bits
193*913ec974Sbeck 	 */
194*913ec974Sbeck 	c1+=(int)(bn_add_words(&(r[n]),&(r[n]),&(t[n2]),n2));
195*913ec974Sbeck 	if (c1)
1965b37fcf3Sryker 		{
197*913ec974Sbeck 		p= &(r[n+n2]);
198*913ec974Sbeck 		lo= *p;
199*913ec974Sbeck 		ln=(lo+c1)&BN_MASK2;
200*913ec974Sbeck 		*p=ln;
201*913ec974Sbeck 
202*913ec974Sbeck 		/* The overflow will stop before we over write
203*913ec974Sbeck 		 * words we should not overwrite */
204*913ec974Sbeck 		if (ln < (BN_ULONG)c1)
205*913ec974Sbeck 			{
206*913ec974Sbeck 			do	{
207*913ec974Sbeck 				p++;
208*913ec974Sbeck 				lo= *p;
209*913ec974Sbeck 				ln=(lo+1)&BN_MASK2;
210*913ec974Sbeck 				*p=ln;
211*913ec974Sbeck 				} while (ln == 0);
212*913ec974Sbeck 			}
213*913ec974Sbeck 		}
2145b37fcf3Sryker 	}
2155b37fcf3Sryker 
216*913ec974Sbeck /* n+tn is the word length
217*913ec974Sbeck  * t needs to be n*4 is size, as does r */
218*913ec974Sbeck void bn_mul_part_recursive(BN_ULONG *r, BN_ULONG *a, BN_ULONG *b, int tn,
219*913ec974Sbeck 	     int n, BN_ULONG *t)
2205b37fcf3Sryker 	{
221*913ec974Sbeck 	int i,j,n2=n*2;
222*913ec974Sbeck 	unsigned int c1;
223*913ec974Sbeck 	BN_ULONG ln,lo,*p;
2245b37fcf3Sryker 
225*913ec974Sbeck #ifdef BN_COUNT
226*913ec974Sbeck printf(" bn_mul_part_recursive %d * %d\n",tn+n,tn+n);
227*913ec974Sbeck #endif
228*913ec974Sbeck 	if (n < 8)
2295b37fcf3Sryker 		{
230*913ec974Sbeck 		i=tn+n;
231*913ec974Sbeck 		bn_mul_normal(r,a,i,b,i);
232*913ec974Sbeck 		return;
2335b37fcf3Sryker 		}
2345b37fcf3Sryker 
235*913ec974Sbeck 	/* r=(a[0]-a[1])*(b[1]-b[0]) */
236*913ec974Sbeck 	bn_sub_words(t,      a,      &(a[n]),n); /* + */
237*913ec974Sbeck 	bn_sub_words(&(t[n]),b,      &(b[n]),n); /* - */
2385b37fcf3Sryker 
239*913ec974Sbeck /*	if (n == 4)
2405b37fcf3Sryker 		{
241*913ec974Sbeck 		bn_mul_comba4(&(t[n2]),t,&(t[n]));
242*913ec974Sbeck 		bn_mul_comba4(r,a,b);
243*913ec974Sbeck 		bn_mul_normal(&(r[n2]),&(a[n]),tn,&(b[n]),tn);
244*913ec974Sbeck 		memset(&(r[n2+tn*2]),0,sizeof(BN_ULONG)*(n2-tn*2));
245*913ec974Sbeck 		}
246*913ec974Sbeck 	else */ if (n == 8)
247*913ec974Sbeck 		{
248*913ec974Sbeck 		bn_mul_comba8(&(t[n2]),t,&(t[n]));
249*913ec974Sbeck 		bn_mul_comba8(r,a,b);
250*913ec974Sbeck 		bn_mul_normal(&(r[n2]),&(a[n]),tn,&(b[n]),tn);
251*913ec974Sbeck 		memset(&(r[n2+tn*2]),0,sizeof(BN_ULONG)*(n2-tn*2));
252*913ec974Sbeck 		}
253*913ec974Sbeck 	else
254*913ec974Sbeck 		{
255*913ec974Sbeck 		p= &(t[n2*2]);
256*913ec974Sbeck 		bn_mul_recursive(&(t[n2]),t,&(t[n]),n,p);
257*913ec974Sbeck 		bn_mul_recursive(r,a,b,n,p);
258*913ec974Sbeck 		i=n/2;
259*913ec974Sbeck 		/* If there is only a bottom half to the number,
260*913ec974Sbeck 		 * just do it */
261*913ec974Sbeck 		j=tn-i;
262*913ec974Sbeck 		if (j == 0)
263*913ec974Sbeck 			{
264*913ec974Sbeck 			bn_mul_recursive(&(r[n2]),&(a[n]),&(b[n]),i,p);
265*913ec974Sbeck 			memset(&(r[n2+i*2]),0,sizeof(BN_ULONG)*(n2-i*2));
266*913ec974Sbeck 			}
267*913ec974Sbeck 		else if (j > 0) /* eg, n == 16, i == 8 and tn == 11 */
268*913ec974Sbeck 				{
269*913ec974Sbeck 				bn_mul_part_recursive(&(r[n2]),&(a[n]),&(b[n]),
270*913ec974Sbeck 					j,i,p);
271*913ec974Sbeck 				memset(&(r[n2+tn*2]),0,
272*913ec974Sbeck 					sizeof(BN_ULONG)*(n2-tn*2));
273*913ec974Sbeck 				}
274*913ec974Sbeck 		else /* (j < 0) eg, n == 16, i == 8 and tn == 5 */
275*913ec974Sbeck 			{
276*913ec974Sbeck 			memset(&(r[n2]),0,sizeof(BN_ULONG)*n2);
277*913ec974Sbeck 			if (tn < BN_MUL_RECURSIVE_SIZE_NORMAL)
278*913ec974Sbeck 				{
279*913ec974Sbeck 				bn_mul_normal(&(r[n2]),&(a[n]),tn,&(b[n]),tn);
280*913ec974Sbeck 				}
281*913ec974Sbeck 			else
282*913ec974Sbeck 				{
283*913ec974Sbeck 				for (;;)
284*913ec974Sbeck 					{
285*913ec974Sbeck 					i/=2;
286*913ec974Sbeck 					if (i < tn)
287*913ec974Sbeck 						{
288*913ec974Sbeck 						bn_mul_part_recursive(&(r[n2]),
289*913ec974Sbeck 							&(a[n]),&(b[n]),
290*913ec974Sbeck 							tn-i,i,p);
291*913ec974Sbeck 						break;
292*913ec974Sbeck 						}
293*913ec974Sbeck 					else if (i == tn)
294*913ec974Sbeck 						{
295*913ec974Sbeck 						bn_mul_recursive(&(r[n2]),
296*913ec974Sbeck 							&(a[n]),&(b[n]),
297*913ec974Sbeck 							i,p);
298*913ec974Sbeck 						break;
299*913ec974Sbeck 						}
300*913ec974Sbeck 					}
301*913ec974Sbeck 				}
302*913ec974Sbeck 			}
3035b37fcf3Sryker 		}
3045b37fcf3Sryker 
305*913ec974Sbeck 	/* t[32] holds (a[0]-a[1])*(b[1]-b[0]), c1 is the sign
306*913ec974Sbeck 	 * r[10] holds (a[0]*b[0])
307*913ec974Sbeck 	 * r[32] holds (b[1]*b[1])
308*913ec974Sbeck 	 */
3095b37fcf3Sryker 
310*913ec974Sbeck 	c1=(int)(bn_add_words(t,r,&(r[n2]),n2));
311*913ec974Sbeck 	c1-=(int)(bn_sub_words(&(t[n2]),t,&(t[n2]),n2));
3125b37fcf3Sryker 
313*913ec974Sbeck 	/* t[32] holds (a[0]-a[1])*(b[1]-b[0])+(a[0]*b[0])+(a[1]*b[1])
314*913ec974Sbeck 	 * r[10] holds (a[0]*b[0])
315*913ec974Sbeck 	 * r[32] holds (b[1]*b[1])
316*913ec974Sbeck 	 * c1 holds the carry bits
317*913ec974Sbeck 	 */
318*913ec974Sbeck 	c1+=(int)(bn_add_words(&(r[n]),&(r[n]),&(t[n2]),n2));
319*913ec974Sbeck 	if (c1)
320*913ec974Sbeck 		{
321*913ec974Sbeck 		p= &(r[n+n2]);
322*913ec974Sbeck 		lo= *p;
323*913ec974Sbeck 		ln=(lo+c1)&BN_MASK2;
324*913ec974Sbeck 		*p=ln;
3255b37fcf3Sryker 
326*913ec974Sbeck 		/* The overflow will stop before we over write
327*913ec974Sbeck 		 * words we should not overwrite */
328*913ec974Sbeck 		if (ln < c1)
329*913ec974Sbeck 			{
330*913ec974Sbeck 			do	{
331*913ec974Sbeck 				p++;
332*913ec974Sbeck 				lo= *p;
333*913ec974Sbeck 				ln=(lo+1)&BN_MASK2;
334*913ec974Sbeck 				*p=ln;
335*913ec974Sbeck 				} while (ln == 0);
336*913ec974Sbeck 			}
337*913ec974Sbeck 		}
338*913ec974Sbeck 	}
3395b37fcf3Sryker 
340*913ec974Sbeck /* a and b must be the same size, which is n2.
341*913ec974Sbeck  * r needs to be n2 words and t needs to be n2*2
342*913ec974Sbeck  */
343*913ec974Sbeck void bn_mul_low_recursive(BN_ULONG *r, BN_ULONG *a, BN_ULONG *b, int n2,
344*913ec974Sbeck 	     BN_ULONG *t)
345*913ec974Sbeck 	{
346*913ec974Sbeck 	int n=n2/2;
3475b37fcf3Sryker 
348*913ec974Sbeck #ifdef BN_COUNT
349*913ec974Sbeck printf(" bn_mul_low_recursive %d * %d\n",n2,n2);
350*913ec974Sbeck #endif
3515b37fcf3Sryker 
352*913ec974Sbeck 	bn_mul_recursive(r,a,b,n,&(t[0]));
353*913ec974Sbeck 	if (n >= BN_MUL_LOW_RECURSIVE_SIZE_NORMAL)
354*913ec974Sbeck 		{
355*913ec974Sbeck 		bn_mul_low_recursive(&(t[0]),&(a[0]),&(b[n]),n,&(t[n2]));
356*913ec974Sbeck 		bn_add_words(&(r[n]),&(r[n]),&(t[0]),n);
357*913ec974Sbeck 		bn_mul_low_recursive(&(t[0]),&(a[n]),&(b[0]),n,&(t[n2]));
358*913ec974Sbeck 		bn_add_words(&(r[n]),&(r[n]),&(t[0]),n);
359*913ec974Sbeck 		}
360*913ec974Sbeck 	else
361*913ec974Sbeck 		{
362*913ec974Sbeck 		bn_mul_low_normal(&(t[0]),&(a[0]),&(b[n]),n);
363*913ec974Sbeck 		bn_mul_low_normal(&(t[n]),&(a[n]),&(b[0]),n);
364*913ec974Sbeck 		bn_add_words(&(r[n]),&(r[n]),&(t[0]),n);
365*913ec974Sbeck 		bn_add_words(&(r[n]),&(r[n]),&(t[n]),n);
366*913ec974Sbeck 		}
367*913ec974Sbeck 	}
3685b37fcf3Sryker 
369*913ec974Sbeck /* a and b must be the same size, which is n2.
370*913ec974Sbeck  * r needs to be n2 words and t needs to be n2*2
371*913ec974Sbeck  * l is the low words of the output.
372*913ec974Sbeck  * t needs to be n2*3
373*913ec974Sbeck  */
374*913ec974Sbeck void bn_mul_high(BN_ULONG *r, BN_ULONG *a, BN_ULONG *b, BN_ULONG *l, int n2,
375*913ec974Sbeck 	     BN_ULONG *t)
376*913ec974Sbeck 	{
377*913ec974Sbeck 	int i,n;
378*913ec974Sbeck 	int c1,c2;
379*913ec974Sbeck 	int neg,oneg,zero;
380*913ec974Sbeck 	BN_ULONG ll,lc,*lp,*mp;
381*913ec974Sbeck 
382*913ec974Sbeck #ifdef BN_COUNT
383*913ec974Sbeck printf(" bn_mul_high %d * %d\n",n2,n2);
384*913ec974Sbeck #endif
385*913ec974Sbeck 	n=n2/2;
386*913ec974Sbeck 
387*913ec974Sbeck 	/* Calculate (al-ah)*(bh-bl) */
388*913ec974Sbeck 	neg=zero=0;
389*913ec974Sbeck 	c1=bn_cmp_words(&(a[0]),&(a[n]),n);
390*913ec974Sbeck 	c2=bn_cmp_words(&(b[n]),&(b[0]),n);
391*913ec974Sbeck 	switch (c1*3+c2)
392*913ec974Sbeck 		{
393*913ec974Sbeck 	case -4:
394*913ec974Sbeck 		bn_sub_words(&(r[0]),&(a[n]),&(a[0]),n);
395*913ec974Sbeck 		bn_sub_words(&(r[n]),&(b[0]),&(b[n]),n);
396*913ec974Sbeck 		break;
397*913ec974Sbeck 	case -3:
398*913ec974Sbeck 		zero=1;
399*913ec974Sbeck 		break;
400*913ec974Sbeck 	case -2:
401*913ec974Sbeck 		bn_sub_words(&(r[0]),&(a[n]),&(a[0]),n);
402*913ec974Sbeck 		bn_sub_words(&(r[n]),&(b[n]),&(b[0]),n);
403*913ec974Sbeck 		neg=1;
404*913ec974Sbeck 		break;
405*913ec974Sbeck 	case -1:
406*913ec974Sbeck 	case 0:
407*913ec974Sbeck 	case 1:
408*913ec974Sbeck 		zero=1;
409*913ec974Sbeck 		break;
410*913ec974Sbeck 	case 2:
411*913ec974Sbeck 		bn_sub_words(&(r[0]),&(a[0]),&(a[n]),n);
412*913ec974Sbeck 		bn_sub_words(&(r[n]),&(b[0]),&(b[n]),n);
413*913ec974Sbeck 		neg=1;
414*913ec974Sbeck 		break;
415*913ec974Sbeck 	case 3:
416*913ec974Sbeck 		zero=1;
417*913ec974Sbeck 		break;
418*913ec974Sbeck 	case 4:
419*913ec974Sbeck 		bn_sub_words(&(r[0]),&(a[0]),&(a[n]),n);
420*913ec974Sbeck 		bn_sub_words(&(r[n]),&(b[n]),&(b[0]),n);
421*913ec974Sbeck 		break;
422*913ec974Sbeck 		}
423*913ec974Sbeck 
424*913ec974Sbeck 	oneg=neg;
425*913ec974Sbeck 	/* t[10] = (a[0]-a[1])*(b[1]-b[0]) */
426*913ec974Sbeck 	/* r[10] = (a[1]*b[1]) */
427*913ec974Sbeck #ifdef BN_MUL_COMBA
428*913ec974Sbeck 	if (n == 8)
429*913ec974Sbeck 		{
430*913ec974Sbeck 		bn_mul_comba8(&(t[0]),&(r[0]),&(r[n]));
431*913ec974Sbeck 		bn_mul_comba8(r,&(a[n]),&(b[n]));
432*913ec974Sbeck 		}
433*913ec974Sbeck 	else
434*913ec974Sbeck #endif
435*913ec974Sbeck 		{
436*913ec974Sbeck 		bn_mul_recursive(&(t[0]),&(r[0]),&(r[n]),n,&(t[n2]));
437*913ec974Sbeck 		bn_mul_recursive(r,&(a[n]),&(b[n]),n,&(t[n2]));
438*913ec974Sbeck 		}
439*913ec974Sbeck 
440*913ec974Sbeck 	/* s0 == low(al*bl)
441*913ec974Sbeck 	 * s1 == low(ah*bh)+low((al-ah)*(bh-bl))+low(al*bl)+high(al*bl)
442*913ec974Sbeck 	 * We know s0 and s1 so the only unknown is high(al*bl)
443*913ec974Sbeck 	 * high(al*bl) == s1 - low(ah*bh+s0+(al-ah)*(bh-bl))
444*913ec974Sbeck 	 * high(al*bl) == s1 - (r[0]+l[0]+t[0])
445*913ec974Sbeck 	 */
446*913ec974Sbeck 	if (l != NULL)
447*913ec974Sbeck 		{
448*913ec974Sbeck 		lp= &(t[n2+n]);
449*913ec974Sbeck 		c1=(int)(bn_add_words(lp,&(r[0]),&(l[0]),n));
450*913ec974Sbeck 		}
451*913ec974Sbeck 	else
452*913ec974Sbeck 		{
453*913ec974Sbeck 		c1=0;
454*913ec974Sbeck 		lp= &(r[0]);
455*913ec974Sbeck 		}
456*913ec974Sbeck 
457*913ec974Sbeck 	if (neg)
458*913ec974Sbeck 		neg=(int)(bn_sub_words(&(t[n2]),lp,&(t[0]),n));
459*913ec974Sbeck 	else
460*913ec974Sbeck 		{
461*913ec974Sbeck 		bn_add_words(&(t[n2]),lp,&(t[0]),n);
462*913ec974Sbeck 		neg=0;
463*913ec974Sbeck 		}
464*913ec974Sbeck 
465*913ec974Sbeck 	if (l != NULL)
466*913ec974Sbeck 		{
467*913ec974Sbeck 		bn_sub_words(&(t[n2+n]),&(l[n]),&(t[n2]),n);
468*913ec974Sbeck 		}
469*913ec974Sbeck 	else
470*913ec974Sbeck 		{
471*913ec974Sbeck 		lp= &(t[n2+n]);
472*913ec974Sbeck 		mp= &(t[n2]);
473*913ec974Sbeck 		for (i=0; i<n; i++)
474*913ec974Sbeck 			lp[i]=((~mp[i])+1)&BN_MASK2;
475*913ec974Sbeck 		}
476*913ec974Sbeck 
477*913ec974Sbeck 	/* s[0] = low(al*bl)
478*913ec974Sbeck 	 * t[3] = high(al*bl)
479*913ec974Sbeck 	 * t[10] = (a[0]-a[1])*(b[1]-b[0]) neg is the sign
480*913ec974Sbeck 	 * r[10] = (a[1]*b[1])
481*913ec974Sbeck 	 */
482*913ec974Sbeck 	/* R[10] = al*bl
483*913ec974Sbeck 	 * R[21] = al*bl + ah*bh + (a[0]-a[1])*(b[1]-b[0])
484*913ec974Sbeck 	 * R[32] = ah*bh
485*913ec974Sbeck 	 */
486*913ec974Sbeck 	/* R[1]=t[3]+l[0]+r[0](+-)t[0] (have carry/borrow)
487*913ec974Sbeck 	 * R[2]=r[0]+t[3]+r[1](+-)t[1] (have carry/borrow)
488*913ec974Sbeck 	 * R[3]=r[1]+(carry/borrow)
489*913ec974Sbeck 	 */
490*913ec974Sbeck 	if (l != NULL)
491*913ec974Sbeck 		{
492*913ec974Sbeck 		lp= &(t[n2]);
493*913ec974Sbeck 		c1= (int)(bn_add_words(lp,&(t[n2+n]),&(l[0]),n));
494*913ec974Sbeck 		}
495*913ec974Sbeck 	else
496*913ec974Sbeck 		{
497*913ec974Sbeck 		lp= &(t[n2+n]);
498*913ec974Sbeck 		c1=0;
499*913ec974Sbeck 		}
500*913ec974Sbeck 	c1+=(int)(bn_add_words(&(t[n2]),lp,  &(r[0]),n));
501*913ec974Sbeck 	if (oneg)
502*913ec974Sbeck 		c1-=(int)(bn_sub_words(&(t[n2]),&(t[n2]),&(t[0]),n));
503*913ec974Sbeck 	else
504*913ec974Sbeck 		c1+=(int)(bn_add_words(&(t[n2]),&(t[n2]),&(t[0]),n));
505*913ec974Sbeck 
506*913ec974Sbeck 	c2 =(int)(bn_add_words(&(r[0]),&(r[0]),&(t[n2+n]),n));
507*913ec974Sbeck 	c2+=(int)(bn_add_words(&(r[0]),&(r[0]),&(r[n]),n));
508*913ec974Sbeck 	if (oneg)
509*913ec974Sbeck 		c2-=(int)(bn_sub_words(&(r[0]),&(r[0]),&(t[n]),n));
510*913ec974Sbeck 	else
511*913ec974Sbeck 		c2+=(int)(bn_add_words(&(r[0]),&(r[0]),&(t[n]),n));
512*913ec974Sbeck 
513*913ec974Sbeck 	if (c1 != 0) /* Add starting at r[0], could be +ve or -ve */
514*913ec974Sbeck 		{
515*913ec974Sbeck 		i=0;
516*913ec974Sbeck 		if (c1 > 0)
517*913ec974Sbeck 			{
518*913ec974Sbeck 			lc=c1;
519*913ec974Sbeck 			do	{
520*913ec974Sbeck 				ll=(r[i]+lc)&BN_MASK2;
521*913ec974Sbeck 				r[i++]=ll;
522*913ec974Sbeck 				lc=(lc > ll);
523*913ec974Sbeck 				} while (lc);
524*913ec974Sbeck 			}
525*913ec974Sbeck 		else
526*913ec974Sbeck 			{
527*913ec974Sbeck 			lc= -c1;
528*913ec974Sbeck 			do	{
529*913ec974Sbeck 				ll=r[i];
530*913ec974Sbeck 				r[i++]=(ll-lc)&BN_MASK2;
531*913ec974Sbeck 				lc=(lc > ll);
532*913ec974Sbeck 				} while (lc);
533*913ec974Sbeck 			}
534*913ec974Sbeck 		}
535*913ec974Sbeck 	if (c2 != 0) /* Add starting at r[1] */
536*913ec974Sbeck 		{
537*913ec974Sbeck 		i=n;
538*913ec974Sbeck 		if (c2 > 0)
539*913ec974Sbeck 			{
540*913ec974Sbeck 			lc=c2;
541*913ec974Sbeck 			do	{
542*913ec974Sbeck 				ll=(r[i]+lc)&BN_MASK2;
543*913ec974Sbeck 				r[i++]=ll;
544*913ec974Sbeck 				lc=(lc > ll);
545*913ec974Sbeck 				} while (lc);
546*913ec974Sbeck 			}
547*913ec974Sbeck 		else
548*913ec974Sbeck 			{
549*913ec974Sbeck 			lc= -c2;
550*913ec974Sbeck 			do	{
551*913ec974Sbeck 				ll=r[i];
552*913ec974Sbeck 				r[i++]=(ll-lc)&BN_MASK2;
553*913ec974Sbeck 				lc=(lc > ll);
554*913ec974Sbeck 				} while (lc);
555*913ec974Sbeck 			}
556*913ec974Sbeck 		}
5575b37fcf3Sryker 	}
5585b37fcf3Sryker #endif
559*913ec974Sbeck 
560*913ec974Sbeck int BN_mul(BIGNUM *r, BIGNUM *a, BIGNUM *b, BN_CTX *ctx)
561*913ec974Sbeck 	{
562*913ec974Sbeck 	int top,al,bl;
563*913ec974Sbeck 	BIGNUM *rr;
564*913ec974Sbeck #ifdef BN_RECURSION
565*913ec974Sbeck 	BIGNUM *t;
566*913ec974Sbeck 	int i,j,k;
567*913ec974Sbeck #endif
568*913ec974Sbeck 
569*913ec974Sbeck #ifdef BN_COUNT
570*913ec974Sbeck printf("BN_mul %d * %d\n",a->top,b->top);
571*913ec974Sbeck #endif
572*913ec974Sbeck 
573*913ec974Sbeck 	bn_check_top(a);
574*913ec974Sbeck 	bn_check_top(b);
575*913ec974Sbeck 	bn_check_top(r);
576*913ec974Sbeck 
577*913ec974Sbeck 	al=a->top;
578*913ec974Sbeck 	bl=b->top;
579*913ec974Sbeck 	r->neg=a->neg^b->neg;
580*913ec974Sbeck 
581*913ec974Sbeck 	if ((al == 0) || (bl == 0))
582*913ec974Sbeck 		{
583*913ec974Sbeck 		BN_zero(r);
584*913ec974Sbeck 		return(1);
585*913ec974Sbeck 		}
586*913ec974Sbeck 	top=al+bl;
587*913ec974Sbeck 
588*913ec974Sbeck 	if ((r == a) || (r == b))
589*913ec974Sbeck 		rr= &(ctx->bn[ctx->tos+1]);
590*913ec974Sbeck 	else
591*913ec974Sbeck 		rr=r;
592*913ec974Sbeck 
593*913ec974Sbeck #if defined(BN_MUL_COMBA) || defined(BN_RECURSION)
594*913ec974Sbeck 	if (al == bl)
595*913ec974Sbeck 		{
596*913ec974Sbeck #  ifdef BN_MUL_COMBA
597*913ec974Sbeck /*		if (al == 4)
598*913ec974Sbeck 			{
599*913ec974Sbeck 			if (bn_wexpand(rr,8) == NULL) return(0);
600*913ec974Sbeck 			rr->top=8;
601*913ec974Sbeck 			bn_mul_comba4(rr->d,a->d,b->d);
602*913ec974Sbeck 			goto end;
603*913ec974Sbeck 			}
604*913ec974Sbeck 		else */ if (al == 8)
605*913ec974Sbeck 			{
606*913ec974Sbeck 			if (bn_wexpand(rr,16) == NULL) return(0);
607*913ec974Sbeck 			rr->top=16;
608*913ec974Sbeck 			bn_mul_comba8(rr->d,a->d,b->d);
609*913ec974Sbeck 			goto end;
610*913ec974Sbeck 			}
611*913ec974Sbeck 		else
612*913ec974Sbeck #  endif
613*913ec974Sbeck #ifdef BN_RECURSION
614*913ec974Sbeck 		if (al < BN_MULL_SIZE_NORMAL)
615*913ec974Sbeck #endif
616*913ec974Sbeck 			{
617*913ec974Sbeck 			if (bn_wexpand(rr,top) == NULL) return(0);
618*913ec974Sbeck 			rr->top=top;
619*913ec974Sbeck 			bn_mul_normal(rr->d,a->d,al,b->d,bl);
620*913ec974Sbeck 			goto end;
621*913ec974Sbeck 			}
622*913ec974Sbeck #  ifdef BN_RECURSION
623*913ec974Sbeck 		goto symetric;
624*913ec974Sbeck #  endif
625*913ec974Sbeck 		}
626*913ec974Sbeck #endif
627*913ec974Sbeck #ifdef BN_RECURSION
628*913ec974Sbeck 	else if ((al < BN_MULL_SIZE_NORMAL) || (bl < BN_MULL_SIZE_NORMAL))
629*913ec974Sbeck 		{
630*913ec974Sbeck 		if (bn_wexpand(rr,top) == NULL) return(0);
631*913ec974Sbeck 		rr->top=top;
632*913ec974Sbeck 		bn_mul_normal(rr->d,a->d,al,b->d,bl);
633*913ec974Sbeck 		goto end;
634*913ec974Sbeck 		}
635*913ec974Sbeck 	else
636*913ec974Sbeck 		{
637*913ec974Sbeck 		i=(al-bl);
638*913ec974Sbeck 		if ((i ==  1) && !BN_get_flags(b,BN_FLG_STATIC_DATA))
639*913ec974Sbeck 			{
640*913ec974Sbeck 			bn_wexpand(b,al);
641*913ec974Sbeck 			b->d[bl]=0;
642*913ec974Sbeck 			bl++;
643*913ec974Sbeck 			goto symetric;
644*913ec974Sbeck 			}
645*913ec974Sbeck 		else if ((i ==  -1) && !BN_get_flags(a,BN_FLG_STATIC_DATA))
646*913ec974Sbeck 			{
647*913ec974Sbeck 			bn_wexpand(a,bl);
648*913ec974Sbeck 			a->d[al]=0;
649*913ec974Sbeck 			al++;
650*913ec974Sbeck 			goto symetric;
651*913ec974Sbeck 			}
652*913ec974Sbeck 		}
653*913ec974Sbeck #endif
654*913ec974Sbeck 
655*913ec974Sbeck 	/* asymetric and >= 4 */
656*913ec974Sbeck 	if (bn_wexpand(rr,top) == NULL) return(0);
657*913ec974Sbeck 	rr->top=top;
658*913ec974Sbeck 	bn_mul_normal(rr->d,a->d,al,b->d,bl);
659*913ec974Sbeck 
660*913ec974Sbeck #ifdef BN_RECURSION
661*913ec974Sbeck 	if (0)
662*913ec974Sbeck 		{
663*913ec974Sbeck symetric:
664*913ec974Sbeck 		/* symetric and > 4 */
665*913ec974Sbeck 		/* 16 or larger */
666*913ec974Sbeck 		j=BN_num_bits_word((BN_ULONG)al);
667*913ec974Sbeck 		j=1<<(j-1);
668*913ec974Sbeck 		k=j+j;
669*913ec974Sbeck 		t= &(ctx->bn[ctx->tos]);
670*913ec974Sbeck 		if (al == j) /* exact multiple */
671*913ec974Sbeck 			{
672*913ec974Sbeck 			bn_wexpand(t,k*2);
673*913ec974Sbeck 			bn_wexpand(rr,k*2);
674*913ec974Sbeck 			bn_mul_recursive(rr->d,a->d,b->d,al,t->d);
675*913ec974Sbeck 			}
676*913ec974Sbeck 		else
677*913ec974Sbeck 			{
678*913ec974Sbeck 			bn_wexpand(a,k);
679*913ec974Sbeck 			bn_wexpand(b,k);
680*913ec974Sbeck 			bn_wexpand(t,k*4);
681*913ec974Sbeck 			bn_wexpand(rr,k*4);
682*913ec974Sbeck 			for (i=a->top; i<k; i++)
683*913ec974Sbeck 				a->d[i]=0;
684*913ec974Sbeck 			for (i=b->top; i<k; i++)
685*913ec974Sbeck 				b->d[i]=0;
686*913ec974Sbeck 			bn_mul_part_recursive(rr->d,a->d,b->d,al-j,j,t->d);
687*913ec974Sbeck 			}
688*913ec974Sbeck 		rr->top=top;
689*913ec974Sbeck 		}
690*913ec974Sbeck #endif
691*913ec974Sbeck #if defined(BN_MUL_COMBA) || defined(BN_RECURSION)
692*913ec974Sbeck end:
693*913ec974Sbeck #endif
694*913ec974Sbeck 	bn_fix_top(rr);
695*913ec974Sbeck 	if (r != rr) BN_copy(r,rr);
696*913ec974Sbeck 	return(1);
697*913ec974Sbeck 	}
698*913ec974Sbeck 
699*913ec974Sbeck void bn_mul_normal(BN_ULONG *r, BN_ULONG *a, int na, BN_ULONG *b, int nb)
700*913ec974Sbeck 	{
701*913ec974Sbeck 	BN_ULONG *rr;
702*913ec974Sbeck 
703*913ec974Sbeck #ifdef BN_COUNT
704*913ec974Sbeck printf(" bn_mul_normal %d * %d\n",na,nb);
705*913ec974Sbeck #endif
706*913ec974Sbeck 
707*913ec974Sbeck 	if (na < nb)
708*913ec974Sbeck 		{
709*913ec974Sbeck 		int itmp;
710*913ec974Sbeck 		BN_ULONG *ltmp;
711*913ec974Sbeck 
712*913ec974Sbeck 		itmp=na; na=nb; nb=itmp;
713*913ec974Sbeck 		ltmp=a;   a=b;   b=ltmp;
714*913ec974Sbeck 
715*913ec974Sbeck 		}
716*913ec974Sbeck 	rr= &(r[na]);
717*913ec974Sbeck 	rr[0]=bn_mul_words(r,a,na,b[0]);
718*913ec974Sbeck 
719*913ec974Sbeck 	for (;;)
720*913ec974Sbeck 		{
721*913ec974Sbeck 		if (--nb <= 0) return;
722*913ec974Sbeck 		rr[1]=bn_mul_add_words(&(r[1]),a,na,b[1]);
723*913ec974Sbeck 		if (--nb <= 0) return;
724*913ec974Sbeck 		rr[2]=bn_mul_add_words(&(r[2]),a,na,b[2]);
725*913ec974Sbeck 		if (--nb <= 0) return;
726*913ec974Sbeck 		rr[3]=bn_mul_add_words(&(r[3]),a,na,b[3]);
727*913ec974Sbeck 		if (--nb <= 0) return;
728*913ec974Sbeck 		rr[4]=bn_mul_add_words(&(r[4]),a,na,b[4]);
729*913ec974Sbeck 		rr+=4;
730*913ec974Sbeck 		r+=4;
731*913ec974Sbeck 		b+=4;
732*913ec974Sbeck 		}
733*913ec974Sbeck 	}
734*913ec974Sbeck 
735*913ec974Sbeck void bn_mul_low_normal(BN_ULONG *r, BN_ULONG *a, BN_ULONG *b, int n)
736*913ec974Sbeck 	{
737*913ec974Sbeck #ifdef BN_COUNT
738*913ec974Sbeck printf(" bn_mul_low_normal %d * %d\n",n,n);
739*913ec974Sbeck #endif
740*913ec974Sbeck 	bn_mul_words(r,a,n,b[0]);
741*913ec974Sbeck 
742*913ec974Sbeck 	for (;;)
743*913ec974Sbeck 		{
744*913ec974Sbeck 		if (--n <= 0) return;
745*913ec974Sbeck 		bn_mul_add_words(&(r[1]),a,n,b[1]);
746*913ec974Sbeck 		if (--n <= 0) return;
747*913ec974Sbeck 		bn_mul_add_words(&(r[2]),a,n,b[2]);
748*913ec974Sbeck 		if (--n <= 0) return;
749*913ec974Sbeck 		bn_mul_add_words(&(r[3]),a,n,b[3]);
750*913ec974Sbeck 		if (--n <= 0) return;
751*913ec974Sbeck 		bn_mul_add_words(&(r[4]),a,n,b[4]);
752*913ec974Sbeck 		r+=4;
753*913ec974Sbeck 		b+=4;
754*913ec974Sbeck 		}
755*913ec974Sbeck 	}
756*913ec974Sbeck 
757