--- pgp/src/mpilib.c	2018/04/24 16:39:10	1.1.1.3
+++ pgp/src/mpilib.c	2018/04/24 16:40:00	1.1.1.4
@@ -538,13 +538,18 @@ int mp_mult(register unitptr prod,
 		     C9AB 6095 C62F 3756 3843 E4D0 3950 7AD9
 		We'll have to fix this some day.  In the meantime, we'll
 		just have the compiler skip over this code.
+
+		BUG NOTE: Possibly fixed.  Needs testing.
 	*/
 	int bits;
 	register unit bitmask;
 	unitptr product, mplier, temp;
 	short mprec,mprec2;
 	unit mplicand[MAX_UNIT_PRECISION];
-	mp_init(prod,0);		/* better clear full width--double precision */
+
+	/* better clear full width--double precision */
+	mp_init(prod+tohigher(global_precision),0);
+
 	if (testeq(multiplicand,0))
 		return(0);	/* zero multiplicand means zero product */
 
@@ -1765,107 +1770,112 @@ int mp_modexp(register unitptr expout,re
 	return(0);		/* normal return */
 }	/* mp_modexp */
 
-
-
-
-/*********************************************************************
-	RSA-specific routines follow.  These are the only functions that
-	are specific to the RSA public key cryptosystem.  The other
-	multiprecision integer math functions may be used for non-RSA
-	applications.
-
-	The RSA public key cryptosystem is patented by the Massachusetts
-	Institute of Technology (U.S. patent #4,405,829).  This patent does
-	not apply outside the USA.  Public Key Partners (PKP) holds the
-	exclusive commercial license to sell and sub-license the RSA public
-	key cryptosystem.  The author of this software implementation of the
-	RSA algorithm is providing this software for educational use only.
-	Licensing this algorithm from PKP is the responsibility of you, the
-	user, not Philip Zimmermann, the author of this software.  The author
-	assumes no liability for any breach of patent law resulting from the
-	unlicensed use of this software by the user.
-*********************************************************************/
-
-
-int rsa_decrypt(unitptr M, unitptr C,
-	unitptr d, unitptr p, unitptr q, unitptr u)
-	/*	Uses Quisquater's Chinese Remainder Theorem shortcut
-		for RSA decryption.
-		M is the output plaintext message.
-		C is the input ciphertext message.
-		d is the secret decryption exponent.
-		p and q are the secret prime factors of n.
-		u is the multiplicative inverse of p, mod q.
-		Note that u is precomputed on the assumption that p<q.
-		n, the common modulus, is not used here because of the
-		Chinese Remainder Theorem shortcut.
+int mp_modexp_crt(unitptr expout, unitptr expin,
+	unitptr p, unitptr q, unitptr ep, unitptr eq, unitptr u)
+	/*	This is a faster modexp for moduli with a known
+		factorisation into two relatively prime factors p and q,
+		and an input relatively prime to the modulus,
+		the Chinese Remainder Theorem to do the computation
+		mod p and mod q, and then combine the results.  This
+		relies on a number of precomputed values, but does not
+		actually require the modulus n or the exponent e.
+
+		expout = expin ^ e mod (p*q).
+		We form this by evaluating
+		p2 = (expin ^ e) mod p and
+		q2 = (expin ^ e) mod q
+		and then combining the two by the CRT.
+
+		Two optimisations of this are possible.  First, we can
+		reduce expin modulo p and q before starting.
+
+		Second, since we know the factorisation of p and q
+		(trivially derived from the factorisation of n = p*q),
+		and expin is relatively prime to both p and q,
+		we can use Euler's theorem, expin^phi(m) = 1 (mod m),
+		to throw away multiples of phi(p) or phi(q) in e.
+		Letting ep = e mod phi(p) and
+		        eq = e mod phi(q)
+		then combining these two speedups, we only need to evaluate
+		p2 = ((expin mod p) ^ ep) mod p and
+		q2 = ((expin mod q) ^ eq) mod q.
+
+		Now we need to apply the CRT.  Starting with
+		expout = p2 (mod p) and
+		expout = q2 (mod q)
+		we can say that expout = p2 + p * k, and if we assume that
+		0 <= p2 < p, then 0 <= expout < p*q for some 0 <= k < q.
+		Since we want expout = q2 (mod q), then
+		p*k = q2-p2 (mod q).  Since p and q are relatively
+		prime, p has a multiplicative inverse u mod q.  In other
+		words,
+		       u = 1/p (mod q).
+		Multiplying by u on both sides gives k = u*(q2-p2) (mod q).
+		Since we want 0 <= k < q, we can thus find k as
+		k = (u * (q2-p2)) mod q.
+
+		Once we have k, evaluating p2 + p * k is easy, and
+		that gives us the result.
+
+		In the detailed implementation, there is a temporary, temp,
+		used to hold intermediate results, p2 is held in expout,
+		and q2 is used as a temporary in the final step when it is
+		no longer needed.  With that, you should be able to
+		understand the code below.
 	*/
 {
-	unit p2[MAX_UNIT_PRECISION];
 	unit q2[MAX_UNIT_PRECISION];
-	unit temp1[MAX_UNIT_PRECISION];
-	unit temp2[MAX_UNIT_PRECISION];
+	unit temp[MAX_UNIT_PRECISION];
 	int status;
 
-	mp_init(M,1);	/* initialize result in case of error */
-
-	if (mp_compare(p,q) >= 0)	/* ensure that p<q */
-	{	/* swap the pointers p and q */
-		unitptr t;
-		t = p;  p = q; q = t;
-	}
-
-/*	Rather than decrypting by computing modexp with full mod n
-	precision, compute a shorter modexp with mod p precision... */
+/*	First, compiute p2 (physically held in M) */
 
-/*	p2 = [ (C mod p)**( d mod (p-1) ) ] mod p		*/
-	mp_move(temp1,p);
-	mp_dec(temp1);			/* temp1 = p-1 */
-	mp_mod(temp2,d,temp1);	/* temp2 = d mod (p-1) ) */
-	mp_mod(temp1,C,p);		/* temp1 = C mod p  */
-	status = mp_modexp(p2,temp1,temp2,p);
+/*	p2 = [ (expin mod p) ^ ep ] mod p		*/
+	mp_mod(temp,expin,p);		/* temp = expin mod p  */
+	status = mp_modexp(expout,temp,ep,p);
 	if (status < 0)	/* mp_modexp returned an error. */
+	{	mp_init(expout,1);
 		return(status);	/* error return */
+	}
 
-/*	Now compute a short modexp with mod q precision... */
+/*	And the same thing for q2 */
 
-/*	q2 = [ (C mod q)**( d mod (q-1) ) ] mod q		*/
-	mp_move(temp1,q);
-	mp_dec(temp1);			/* temp1 = q-1 */
-	mp_mod(temp2,d,temp1);	/* temp2 = d mod (q-1) */
-	mp_mod(temp1,C,q);		/* temp1 = C mod q  */
-	status = mp_modexp(q2,temp1,temp2,q);
+/*	q2 = [ (expin mod q) ^ eq ] mod q		*/
+	mp_mod(temp,expin,q);		/* temp = expin mod q  */
+	status = mp_modexp(q2,temp,eq,q);
 	if (status < 0)	/* mp_modexp returned an error. */
+	{	mp_init(expout,1);
 		return(status);	/* error return */
+	}
 
 /*	Now use the multiplicative inverse u to glue together the
-	two halves, saving a lot of time by avoiding a full mod n
-	exponentiation.  At key generation time, u was computed
-	with the secret key as the multiplicative inverse of
-	p, mod q such that:  (p*u) mod q = 1, assuming p<q.
-*/
-	if (mp_compare(p2,q2) == 0)	/* only happens if C<p */
-		mp_move(M,p2);
-	else
-	{	/*	q2 = q2 - p2;  if q2<0 then q2 = q2 + q  */
-		if (mp_sub(q2,p2))	/* if the result went negative... */
-			mp_add(q2,q);		/* add q to q2 */
-
-		/*	M = p2 + ( p * [(q2*u) mod q] )	*/
-		mp_mult(temp1,q2,u);		/* temp1 = q2*u  */
-		mp_mod(temp2,temp1,q);		/* temp2 = temp1 mod q */
-		mp_mult(temp1,p,temp2);	/* temp1 = p * temp2 */
-		mp_add(temp1,p2);		/* temp1 = temp1 + p2 */
-		mp_move(M,temp1);		/* M = temp1 */
+	two halves.
+*/
+#if 0
+/* This optimisation is useful if you commonly get small results,
+   but for random results and large q, the odds of (1/q) of it
+   being useful do not warrant the test.
+*/
+	if (mp_compare(expout,q2) != 0)
+	{
+#endif
+	/*	Find q2-p2 mod q */
+	if (mp_sub(q2,expout))	/* if the result went negative */
+		mp_add(q2,q);	/* add q to q2 */
+
+	/*	expout = p2 + ( p * [(q2*u) mod q] )	*/
+	mp_mult(temp,q2,u);		/* q2*u  */
+	mp_mod(q2,temp,q);		/* (q2*u) mod q */
+	mp_mult(temp,p,q2);		/* p * [(q2*u) mod q] */
+	mp_add(expout,temp);		/* expout = p2 + p * [...] */
+#if 0
 	}
+#endif
 
-	mp_burn(p2);	/* burn the evidence on the stack...*/
-	mp_burn(q2);
-	mp_burn(temp1);
-	mp_burn(temp2);
+	mp_burn(q2);	/* burn the evidence on the stack...*/
+	mp_burn(temp);
 	return(0);	/* normal return */
-}	/* rsa_decrypt */
+}	/* mp_modexp_crt */
 
 
 /****************** end of MPI library ****************************/
-