Totally remove the supposedly 'faster' variant in

BN_mod_mul_montgomery, which calls bn_sqr_recursive
without much preparation.

bn_sqr_recursive requires the length of its argument to be
a power of 2, which is not always the case here.
There's no reason for not using BN_sqr -- if a simpler
approach to squaring made sense, then why not change
BN_sqr?  (Using BN_sqr should also speed up DH where g is chosen
such that it becomes small [e.g., 2] when converted
to Montgomery representation.)

Case closed :-)

CHANGES

@@ -4,7 +4,7 @@
 
  Changes between 0.9.5a and 0.9.6  [xx XXX 2000]
 
-  *) Disable 'optimized' squaring variant in BN_mod_mul_montgomery,
+  *) Remove 'optimized' squaring variant in BN_mod_mul_montgomery,
      it can return incorrect results.
      (Note: The buggy variant was not enabled in OpenSSL 0.9.5a,
      but it was in 0.9.6-beta[12].)

crypto/bn/bn_mont.c

@@ -85,36 +85,7 @@ int BN_mod_mul_montgomery(BIGNUM *r, BIGNUM *a, BIGNUM *b,
 
 	if (a == b)
 		{
-#if 0 /* buggy -- try squaring  g  (after converting it to Montgomery
-         representation) in the following parameters
-         (but note that squaring 2 or 4 works):
-Diffie-Hellman-Parameters: (1024 bit)
-    prime:
-        00:ff:ff:ff:ff:ff:ff:ff:ff:c9:0f:da:a2:21:68:
-        c2:34:c4:c6:62:8b:80:dc:1c:d1:29:02:4e:08:8a:
-        67:cc:74:02:0b:be:a6:3b:13:9b:22:51:4a:08:79:
-        8e:34:04:dd:ef:95:19:b3:cd:3a:43:1b:30:2b:0a:
-        6d:f2:5f:14:37:4f:e1:35:6d:6d:51:c2:45:e4:85:
-        b5:76:62:5e:7e:c6:f4:4c:42:e9:a6:37:ed:6b:0b:
-        ff:5c:b6:f4:06:b7:ed:ee:38:6b:fb:5a:89:9f:a5:
-        ae:9f:24:11:7c:4b:1f:e6:49:28:66:51:ec:e6:53:
-        81:ff:ff:ff:ff:ff:ff:ff:ff
-    generator: 8 (0x8)
------BEGIN DH PARAMETERS-----
-MIGHAoGBAP//////////yQ/aoiFowjTExmKLgNwc0SkCTgiKZ8x0Agu+pjsTmyJR
-Sgh5jjQE3e+VGbPNOkMbMCsKbfJfFDdP4TVtbVHCReSFtXZiXn7G9ExC6aY37WsL
-/1y29Aa37e44a/taiZ+lrp8kEXxLH+ZJKGZR7OZTgf//////////AgEI
------END DH PARAMETERS-----
-*/
-		bn_wexpand(tmp,a->top*2);
-		bn_wexpand(tmp2,a->top*4);
-		bn_sqr_recursive(tmp->d,a->d,a->top,tmp2->d);
-		tmp->top=a->top*2;
-		while (tmp->top > 0 && tmp->d[tmp->top-1] == 0)
-			tmp->top--;
-#else
 		if (!BN_sqr(tmp,a,ctx)) goto err;
-#endif
 		}
 	else
 		{

crypto/bn/bn_sqr.c

@@ -188,7 +188,7 @@ void bn_sqr_normal(BN_ULONG *r, BN_ULONG *a, int n, BN_ULONG *tmp)
 
 #ifdef BN_RECURSION
 /* r is 2*n words in size,
- * a and b are both n words in size.
+ * a and b are both n words in size.    (There's not actually a 'b' here ...)
  * n must be a power of 2.
  * We multiply and return the result.
  * t must be 2*n words in size