Implement TwicePlus optimization in Fp curves

3 files changed, 439 insertions, 312 deletions

diff --git a/crypto/src/math/ec/ECPoint.cs b/crypto/src/math/ec/ECPoint.cs
index e4ce58d2d..d81558939 100644
--- a/crypto/src/math/ec/ECPoint.cs
+++ b/crypto/src/math/ec/ECPoint.cs
@@ -128,6 +128,16 @@ namespace Org.BouncyCastle.Math.EC
         public abstract ECPoint Twice();
         public abstract ECPoint Multiply(BigInteger b);
 
+        public virtual ECPoint TwicePlus(ECPoint b)
+        {
+            return Twice().Add(b);
+        }
+
+        public virtual ECPoint ThreeTimes()
+        {
+            return TwicePlus(this);
+        }
+
         /**
         * Sets the appropriate <code>ECMultiplier</code>, unless already set. 
         */
@@ -139,7 +149,7 @@ namespace Org.BouncyCastle.Math.EC
                 {
                     if (this.multiplier == null)
                     {
-                        this.multiplier = new FpNafMultiplier();
+                        this.multiplier = new WNafMultiplier();
                     }
                 }
             }
@@ -270,54 +280,184 @@ namespace Org.BouncyCastle.Math.EC
             ECPoint b)
         {
             if (this.IsInfinity)
+            {
                 return b;
-
+            }
             if (b.IsInfinity)
+            {
                 return this;
+            }
+            if (this == b)
+            {
+                return Twice();
+            }
+
+            ECFieldElement X1 = this.x, Y1 = this.y;
+            ECFieldElement X2 = b.x, Y2 = b.y;
 
-            // Check if b = this or b = -this
-            if (this.x.Equals(b.x))
+            ECFieldElement dx = X2.Subtract(X1), dy = Y2.Subtract(Y1);
+
+            if (dx.IsZero)
             {
-                if (this.y.Equals(b.y))
+                if (dy.IsZero)
                 {
-                    // this = b, i.e. this must be doubled
-                    return this.Twice();
+                    // this == b, i.e. this must be doubled
+                    return Twice();
                 }
 
-                Debug.Assert(this.y.Equals(b.y.Negate()));
-
-                // this = -b, i.e. the result is the point at infinity
-                return this.curve.Infinity;
+                // this == -b, i.e. the result is the point at infinity
+                return curve.Infinity;
             }
 
-            ECFieldElement gamma = b.y.Subtract(this.y).Divide(b.x.Subtract(this.x));
-
-            ECFieldElement x3 = gamma.Square().Subtract(this.x).Subtract(b.x);
-            ECFieldElement y3 = gamma.Multiply(this.x.Subtract(x3)).Subtract(this.y);
+            ECFieldElement gamma = dy.Divide(dx);
+            ECFieldElement X3 = gamma.Square().Subtract(X1).Subtract(X2);
+            ECFieldElement Y3 = gamma.Multiply(X1.Subtract(X3)).Subtract(Y1);
 
-            return new FpPoint(curve, x3, y3, withCompression);
+            return new FpPoint(curve, X3, Y3, this.withCompression);
         }
 
         // B.3 pg 62
         public override ECPoint Twice()
         {
-            // Twice identity element (point at infinity) is identity
             if (this.IsInfinity)
+            {
                 return this;
+            }
 
-            // if y1 == 0, then (x1, y1) == (x1, -y1)
-            // and hence this = -this and thus 2(x1, y1) == infinity
-            if (this.y.ToBigInteger().SignValue == 0)
-                return this.curve.Infinity;
+            ECFieldElement Y1 = this.y;
+            if (Y1.IsZero) 
+            {
+                return curve.Infinity;
+            }
 
-            ECFieldElement TWO = this.curve.FromBigInteger(BigInteger.Two);
-            ECFieldElement THREE = this.curve.FromBigInteger(BigInteger.Three);
-            ECFieldElement gamma = this.x.Square().Multiply(THREE).Add(curve.a).Divide(y.Multiply(TWO));
+            ECFieldElement X1 = this.x;
 
-            ECFieldElement x3 = gamma.Square().Subtract(this.x.Multiply(TWO));
-            ECFieldElement y3 = gamma.Multiply(this.x.Subtract(x3)).Subtract(this.y);
+            ECFieldElement X1Squared = X1.Square();
+            ECFieldElement gamma = Three(X1Squared).Add(this.Curve.A).Divide(Two(Y1));
+            ECFieldElement X3 = gamma.Square().Subtract(Two(X1));
+            ECFieldElement Y3 = gamma.Multiply(X1.Subtract(X3)).Subtract(Y1);
 
-            return new FpPoint(curve, x3, y3, this.withCompression);
+            return new FpPoint(curve, X3, Y3, this.withCompression);
+        }
+
+        public override ECPoint TwicePlus(ECPoint b)
+        {
+            if (this == b)
+            {
+                return ThreeTimes();
+            }
+            if (this.IsInfinity)
+            {
+                return b;
+            }
+            if (b.IsInfinity)
+            {
+                return Twice();
+            }
+
+            ECFieldElement Y1 = this.y;
+            if (Y1.IsZero)
+            {
+                return b;
+            }
+
+            ECFieldElement X1 = this.x;
+            ECFieldElement X2 = b.x, Y2 = b.y;
+
+            ECFieldElement dx = X2.Subtract(X1), dy = Y2.Subtract(Y1);
+
+            if (dx.IsZero)
+            {
+                if (dy.IsZero)
+                {
+                    // this == b i.e. the result is 3P
+                    return ThreeTimes();
+                }
+
+                // this == -b, i.e. the result is P
+                return this;
+            }
+
+            /*
+             * Optimized calculation of 2P + Q, as described in "Trading Inversions for
+             * Multiplications in Elliptic Curve Cryptography", by Ciet, Joye, Lauter, Montgomery.
+             */
+
+            ECFieldElement X = dx.Square(), Y = dy.Square();
+            ECFieldElement d = X.Multiply(Two(X1).Add(X2)).Subtract(Y);
+            if (d.IsZero)
+            {
+                return curve.Infinity;
+            }
+
+            ECFieldElement D = d.Multiply(dx);
+            ECFieldElement I = D.Invert();
+            ECFieldElement L1 = d.Multiply(I).Multiply(dy);
+            ECFieldElement L2 = Two(Y1).Multiply(X).Multiply(dx).Multiply(I).Subtract(L1);
+            ECFieldElement X4 = (L2.Subtract(L1)).Multiply(L1.Add(L2)).Add(X2);
+            ECFieldElement Y4 = (X1.Subtract(X4)).Multiply(L2).Subtract(Y1);
+
+            return new FpPoint(curve, X4, Y4, this.withCompression);
+        }
+
+        public override ECPoint ThreeTimes()
+        {
+            if (this.IsInfinity || this.y.IsZero)
+            {
+                return this;
+            }
+
+            ECFieldElement X1 = this.x, Y1 = this.y;
+
+            ECFieldElement _2Y1 = Two(Y1);
+            ECFieldElement X = _2Y1.Square();
+            ECFieldElement Z = Three(X1.Square()).Add(this.Curve.A);
+            ECFieldElement Y = Z.Square();
+
+            ECFieldElement d = Three(X1).Multiply(X).Subtract(Y);
+            if (d.IsZero)
+            {
+                return this.Curve.Infinity;
+            }
+
+            ECFieldElement D = d.Multiply(_2Y1);
+            ECFieldElement I = D.Invert();
+            ECFieldElement L1 = d.Multiply(I).Multiply(Z);
+            ECFieldElement L2 = X.Square().Multiply(I).Subtract(L1);
+
+            ECFieldElement X4 = (L2.Subtract(L1)).Multiply(L1.Add(L2)).Add(X1);
+            ECFieldElement Y4 = (X1.Subtract(X4)).Multiply(L2).Subtract(Y1);
+            return new FpPoint(curve, X4, Y4, this.withCompression);
+        }
+
+        protected virtual ECFieldElement Two(ECFieldElement x)
+        {
+            return x.Add(x);
+        }
+
+        protected virtual ECFieldElement Three(ECFieldElement x)
+        {
+            return Two(x).Add(x);
+        }
+
+        protected virtual ECFieldElement Four(ECFieldElement x)
+        {
+            return Two(Two(x));
+        }
+
+        protected virtual ECFieldElement Eight(ECFieldElement x)
+        {
+            return Four(Two(x));
+        }
+
+        protected virtual ECFieldElement DoubleProductFromSquares(ECFieldElement a, ECFieldElement b,
+            ECFieldElement aSquared, ECFieldElement bSquared)
+        {
+            /*
+             * NOTE: If squaring in the field is faster than multiplication, then this is a quicker
+             * way to calculate 2.A.B, if A^2 and B^2 are already known.
+             */
+            return a.Add(b).Square().Subtract(aSquared).Subtract(bSquared);
         }
 
         // D.3.2 pg 102 (see Note:)
@@ -335,23 +475,6 @@ namespace Org.BouncyCastle.Math.EC
         {
             return new FpPoint(this.curve, this.x, this.y.Negate(), this.withCompression);
         }
-
-        /**
-         * Sets the default <code>ECMultiplier</code>, unless already set. 
-         */
-        internal override void AssertECMultiplier()
-        {
-            if (this.multiplier == null)
-            {
-                lock (this)
-                {
-                    if (this.multiplier == null)
-                    {
-                        this.multiplier = new WNafMultiplier();
-                    }
-                }
-            }
-        }
     }
 
     /**
diff --git a/crypto/src/math/ec/multiplier/FpNafMultiplier.cs b/crypto/src/math/ec/multiplier/FpNafMultiplier.cs
index f5a98501a..3453e5600 100644
--- a/crypto/src/math/ec/multiplier/FpNafMultiplier.cs
+++ b/crypto/src/math/ec/multiplier/FpNafMultiplier.cs
@@ -1,39 +1,41 @@
 namespace Org.BouncyCastle.Math.EC.Multiplier
 {
-	/**
-	* Class implementing the NAF (Non-Adjacent Form) multiplication algorithm.
-	*/
-	internal class FpNafMultiplier
-		: ECMultiplier
-	{
-		/**
-		* D.3.2 pg 101
-		* @see org.bouncycastle.math.ec.multiplier.ECMultiplier#multiply(org.bouncycastle.math.ec.ECPoint, java.math.BigInteger)
-		*/
-		public ECPoint Multiply(ECPoint p, BigInteger k, PreCompInfo preCompInfo)
-		{
-			// TODO Probably should try to add this
-			// BigInteger e = k.Mod(n); // n == order of p
-			BigInteger e = k;
-			BigInteger h = e.Multiply(BigInteger.Three);
+    /**
+    * Class implementing the NAF (Non-Adjacent Form) multiplication algorithm.
+    */
+    internal class FpNafMultiplier
+        : ECMultiplier
+    {
+        /**
+        * D.3.2 pg 101
+        * @see org.bouncycastle.math.ec.multiplier.ECMultiplier#multiply(org.bouncycastle.math.ec.ECPoint, java.math.BigInteger)
+        */
+        public ECPoint Multiply(ECPoint p, BigInteger k, PreCompInfo preCompInfo)
+        {
+            // TODO Probably should try to add this
+            // BigInteger e = k.Mod(n); // n == order of p
+            BigInteger e = k;
+            BigInteger h = e.Multiply(BigInteger.Three);
 
-			ECPoint neg = p.Negate();
-			ECPoint R = p;
+            ECPoint neg = p.Negate();
+            ECPoint R = p;
 
-			for (int i = h.BitLength - 2; i > 0; --i)
-			{             
-				R = R.Twice();
+            for (int i = h.BitLength - 2; i > 0; --i)
+            {             
+                bool hBit = h.TestBit(i);
+                bool eBit = e.TestBit(i);
 
-				bool hBit = h.TestBit(i);
-				bool eBit = e.TestBit(i);
+                if (hBit == eBit)
+                {
+                    R = R.Twice();
+                }
+                else
+                {
+                    R = R.TwicePlus(hBit ? p : neg);
+                }
+            }
 
-				if (hBit != eBit)
-				{
-					R = R.Add(hBit ? p : neg);
-				}
-			}
-
-			return R;
-		}
-	}
+            return R;
+        }
+    }
 }
diff --git a/crypto/src/math/ec/multiplier/WNafMultiplier.cs b/crypto/src/math/ec/multiplier/WNafMultiplier.cs
index b5cf34ba8..570ceee95 100644
--- a/crypto/src/math/ec/multiplier/WNafMultiplier.cs
+++ b/crypto/src/math/ec/multiplier/WNafMultiplier.cs
@@ -2,240 +2,242 @@ using System;
 
 namespace Org.BouncyCastle.Math.EC.Multiplier
 {
-	/**
-	* Class implementing the WNAF (Window Non-Adjacent Form) multiplication
-	* algorithm.
-	*/
-	internal class WNafMultiplier
-		: ECMultiplier 
-	{
-		/**
-		* Computes the Window NAF (non-adjacent Form) of an integer.
-		* @param width The width <code>w</code> of the Window NAF. The width is
-		* defined as the minimal number <code>w</code>, such that for any
-		* <code>w</code> consecutive digits in the resulting representation, at
-		* most one is non-zero.
-		* @param k The integer of which the Window NAF is computed.
-		* @return The Window NAF of the given width, such that the following holds:
-		* <code>k = &#8722;<sub>i=0</sub><sup>l-1</sup> k<sub>i</sub>2<sup>i</sup>
-		* </code>, where the <code>k<sub>i</sub></code> denote the elements of the
-		* returned <code>sbyte[]</code>.
-		*/
-		public sbyte[] WindowNaf(sbyte width, BigInteger k)
-		{
-			// The window NAF is at most 1 element longer than the binary
-			// representation of the integer k. sbyte can be used instead of short or
-			// int unless the window width is larger than 8. For larger width use
-			// short or int. However, a width of more than 8 is not efficient for
-			// m = log2(q) smaller than 2305 Bits. Note: Values for m larger than
-			// 1000 Bits are currently not used in practice.
-			sbyte[] wnaf = new sbyte[k.BitLength + 1];
-
-			// 2^width as short and BigInteger
-			short pow2wB = (short)(1 << width);
-			BigInteger pow2wBI = BigInteger.ValueOf(pow2wB);
-
-			int i = 0;
-
-			// The actual length of the WNAF
-			int length = 0;
-
-			// while k >= 1
-			while (k.SignValue > 0)
-			{
-				// if k is odd
-				if (k.TestBit(0))
-				{
-					// k Mod 2^width
-					BigInteger remainder = k.Mod(pow2wBI);
-
-					// if remainder > 2^(width - 1) - 1
-					if (remainder.TestBit(width - 1))
-					{
-						wnaf[i] = (sbyte)(remainder.IntValue - pow2wB);
-					}
-					else
-					{
-						wnaf[i] = (sbyte)remainder.IntValue;
-					}
-					// wnaf[i] is now in [-2^(width-1), 2^(width-1)-1]
-
-					k = k.Subtract(BigInteger.ValueOf(wnaf[i]));
-					length = i;
-				}
-				else
-				{
-					wnaf[i] = 0;
-				}
-
-				// k = k/2
-				k = k.ShiftRight(1);
-				i++;
-			}
-
-			length++;
-
-			// Reduce the WNAF array to its actual length
-			sbyte[] wnafShort = new sbyte[length];
-			Array.Copy(wnaf, 0, wnafShort, 0, length);
-			return wnafShort;
-		}
-
-		/**
-		* Multiplies <code>this</code> by an integer <code>k</code> using the
-		* Window NAF method.
-		* @param k The integer by which <code>this</code> is multiplied.
-		* @return A new <code>ECPoint</code> which equals <code>this</code>
-		* multiplied by <code>k</code>.
-		*/
-		public ECPoint Multiply(ECPoint p, BigInteger k, PreCompInfo preCompInfo)
-		{
-			WNafPreCompInfo wnafPreCompInfo;
-
-			if ((preCompInfo != null) && (preCompInfo is WNafPreCompInfo))
-			{
-				wnafPreCompInfo = (WNafPreCompInfo)preCompInfo;
-			}
-			else
-			{
-				// Ignore empty PreCompInfo or PreCompInfo of incorrect type
-				wnafPreCompInfo = new WNafPreCompInfo();
-			}
-
-			// floor(log2(k))
-			int m = k.BitLength;
-
-			// width of the Window NAF
-			sbyte width;
-
-			// Required length of precomputation array
-			int reqPreCompLen;
-
-			// Determine optimal width and corresponding length of precomputation
-			// array based on literature values
-			if (m < 13)
-			{
-				width = 2;
-				reqPreCompLen = 1;
-			}
-			else
-			{
-				if (m < 41)
-				{
-					width = 3;
-					reqPreCompLen = 2;
-				}
-				else
-				{
-					if (m < 121)
-					{
-						width = 4;
-						reqPreCompLen = 4;
-					}
-					else
-					{
-						if (m < 337)
-						{
-							width = 5;
-							reqPreCompLen = 8;
-						}
-						else
-						{
-							if (m < 897)
-							{
-								width = 6;
-								reqPreCompLen = 16;
-							}
-							else
-							{
-								if (m < 2305)
-								{
-									width = 7;
-									reqPreCompLen = 32;
-								}
-								else 
-								{
-									width = 8;
-									reqPreCompLen = 127;
-								}
-							}
-						}
-					}
-				}
-			}
-
-			// The length of the precomputation array
-			int preCompLen = 1;
-
-			ECPoint[] preComp = wnafPreCompInfo.GetPreComp();
-			ECPoint twiceP = wnafPreCompInfo.GetTwiceP();
-
-			// Check if the precomputed ECPoints already exist
-			if (preComp == null)
-			{
-				// Precomputation must be performed from scratch, create an empty
-				// precomputation array of desired length
-				preComp = new ECPoint[]{ p };
-			}
-			else
-			{
-				// Take the already precomputed ECPoints to start with
-				preCompLen = preComp.Length;
-			}
-
-			if (twiceP == null)
-			{
-				// Compute twice(p)
-				twiceP = p.Twice();
-			}
-
-			if (preCompLen < reqPreCompLen)
-			{
-				// Precomputation array must be made bigger, copy existing preComp
-				// array into the larger new preComp array
-				ECPoint[] oldPreComp = preComp;
-				preComp = new ECPoint[reqPreCompLen];
-				Array.Copy(oldPreComp, 0, preComp, 0, preCompLen);
-
-				for (int i = preCompLen; i < reqPreCompLen; i++)
-				{
-					// Compute the new ECPoints for the precomputation array.
-					// The values 1, 3, 5, ..., 2^(width-1)-1 times p are
-					// computed
-					preComp[i] = twiceP.Add(preComp[i - 1]);
-				}            
-			}
-
-			// Compute the Window NAF of the desired width
-			sbyte[] wnaf = WindowNaf(width, k);
-			int l = wnaf.Length;
-
-			// Apply the Window NAF to p using the precomputed ECPoint values.
-			ECPoint q = p.Curve.Infinity;
-			for (int i = l - 1; i >= 0; i--)
-			{
-				q = q.Twice();
-
-				if (wnaf[i] != 0)
-				{
-					if (wnaf[i] > 0)
-					{
-						q = q.Add(preComp[(wnaf[i] - 1)/2]);
-					}
-					else
-					{
-						// wnaf[i] < 0
-						q = q.Subtract(preComp[(-wnaf[i] - 1)/2]);
-					}
-				}
-			}
-
-			// Set PreCompInfo in ECPoint, such that it is available for next
-			// multiplication.
-			wnafPreCompInfo.SetPreComp(preComp);
-			wnafPreCompInfo.SetTwiceP(twiceP);
-			p.SetPreCompInfo(wnafPreCompInfo);
-			return q;
-		}
-	}
+    /**
+    * Class implementing the WNAF (Window Non-Adjacent Form) multiplication
+    * algorithm.
+    */
+    internal class WNafMultiplier
+        : ECMultiplier 
+    {
+        /**
+        * Computes the Window NAF (non-adjacent Form) of an integer.
+        * @param width The width <code>w</code> of the Window NAF. The width is
+        * defined as the minimal number <code>w</code>, such that for any
+        * <code>w</code> consecutive digits in the resulting representation, at
+        * most one is non-zero.
+        * @param k The integer of which the Window NAF is computed.
+        * @return The Window NAF of the given width, such that the following holds:
+        * <code>k = &#8722;<sub>i=0</sub><sup>l-1</sup> k<sub>i</sub>2<sup>i</sup>
+        * </code>, where the <code>k<sub>i</sub></code> denote the elements of the
+        * returned <code>sbyte[]</code>.
+        */
+        public sbyte[] WindowNaf(sbyte width, BigInteger k)
+        {
+            // The window NAF is at most 1 element longer than the binary
+            // representation of the integer k. sbyte can be used instead of short or
+            // int unless the window width is larger than 8. For larger width use
+            // short or int. However, a width of more than 8 is not efficient for
+            // m = log2(q) smaller than 2305 Bits. Note: Values for m larger than
+            // 1000 Bits are currently not used in practice.
+            sbyte[] wnaf = new sbyte[k.BitLength + 1];
+
+            // 2^width as short and BigInteger
+            short pow2wB = (short)(1 << width);
+            BigInteger pow2wBI = BigInteger.ValueOf(pow2wB);
+
+            int i = 0;
+
+            // The actual length of the WNAF
+            int length = 0;
+
+            // while k >= 1
+            while (k.SignValue > 0)
+            {
+                // if k is odd
+                if (k.TestBit(0))
+                {
+                    // k Mod 2^width
+                    BigInteger remainder = k.Mod(pow2wBI);
+
+                    // if remainder > 2^(width - 1) - 1
+                    if (remainder.TestBit(width - 1))
+                    {
+                        wnaf[i] = (sbyte)(remainder.IntValue - pow2wB);
+                    }
+                    else
+                    {
+                        wnaf[i] = (sbyte)remainder.IntValue;
+                    }
+                    // wnaf[i] is now in [-2^(width-1), 2^(width-1)-1]
+
+                    k = k.Subtract(BigInteger.ValueOf(wnaf[i]));
+                    length = i;
+                }
+                else
+                {
+                    wnaf[i] = 0;
+                }
+
+                // k = k/2
+                k = k.ShiftRight(1);
+                i++;
+            }
+
+            length++;
+
+            // Reduce the WNAF array to its actual length
+            sbyte[] wnafShort = new sbyte[length];
+            Array.Copy(wnaf, 0, wnafShort, 0, length);
+            return wnafShort;
+        }
+
+        /**
+        * Multiplies <code>this</code> by an integer <code>k</code> using the
+        * Window NAF method.
+        * @param k The integer by which <code>this</code> is multiplied.
+        * @return A new <code>ECPoint</code> which equals <code>this</code>
+        * multiplied by <code>k</code>.
+        */
+        public ECPoint Multiply(ECPoint p, BigInteger k, PreCompInfo preCompInfo)
+        {
+            WNafPreCompInfo wnafPreCompInfo;
+
+            if ((preCompInfo != null) && (preCompInfo is WNafPreCompInfo))
+            {
+                wnafPreCompInfo = (WNafPreCompInfo)preCompInfo;
+            }
+            else
+            {
+                // Ignore empty PreCompInfo or PreCompInfo of incorrect type
+                wnafPreCompInfo = new WNafPreCompInfo();
+            }
+
+            // floor(log2(k))
+            int m = k.BitLength;
+
+            // width of the Window NAF
+            sbyte width;
+
+            // Required length of precomputation array
+            int reqPreCompLen;
+
+            // Determine optimal width and corresponding length of precomputation
+            // array based on literature values
+            if (m < 13)
+            {
+                width = 2;
+                reqPreCompLen = 1;
+            }
+            else
+            {
+                if (m < 41)
+                {
+                    width = 3;
+                    reqPreCompLen = 2;
+                }
+                else
+                {
+                    if (m < 121)
+                    {
+                        width = 4;
+                        reqPreCompLen = 4;
+                    }
+                    else
+                    {
+                        if (m < 337)
+                        {
+                            width = 5;
+                            reqPreCompLen = 8;
+                        }
+                        else
+                        {
+                            if (m < 897)
+                            {
+                                width = 6;
+                                reqPreCompLen = 16;
+                            }
+                            else
+                            {
+                                if (m < 2305)
+                                {
+                                    width = 7;
+                                    reqPreCompLen = 32;
+                                }
+                                else 
+                                {
+                                    width = 8;
+                                    reqPreCompLen = 127;
+                                }
+                            }
+                        }
+                    }
+                }
+            }
+
+            // The length of the precomputation array
+            int preCompLen = 1;
+
+            ECPoint[] preComp = wnafPreCompInfo.GetPreComp();
+            ECPoint twiceP = wnafPreCompInfo.GetTwiceP();
+
+            // Check if the precomputed ECPoints already exist
+            if (preComp == null)
+            {
+                // Precomputation must be performed from scratch, create an empty
+                // precomputation array of desired length
+                preComp = new ECPoint[]{ p };
+            }
+            else
+            {
+                // Take the already precomputed ECPoints to start with
+                preCompLen = preComp.Length;
+            }
+
+            if (twiceP == null)
+            {
+                // Compute twice(p)
+                twiceP = p.Twice();
+            }
+
+            if (preCompLen < reqPreCompLen)
+            {
+                // Precomputation array must be made bigger, copy existing preComp
+                // array into the larger new preComp array
+                ECPoint[] oldPreComp = preComp;
+                preComp = new ECPoint[reqPreCompLen];
+                Array.Copy(oldPreComp, 0, preComp, 0, preCompLen);
+
+                for (int i = preCompLen; i < reqPreCompLen; i++)
+                {
+                    // Compute the new ECPoints for the precomputation array.
+                    // The values 1, 3, 5, ..., 2^(width-1)-1 times p are
+                    // computed
+                    preComp[i] = twiceP.Add(preComp[i - 1]);
+                }            
+            }
+
+            // Compute the Window NAF of the desired width
+            sbyte[] wnaf = WindowNaf(width, k);
+            int l = wnaf.Length;
+
+            // Apply the Window NAF to p using the precomputed ECPoint values.
+            ECPoint q = p.Curve.Infinity;
+            for (int i = l - 1; i >= 0; i--)
+            {
+                int wi = wnaf[i];
+                if (wi == 0)
+                {
+                    q = q.Twice();
+                }
+                else
+                {
+                    int n = System.Math.Abs(wi);
+                    ECPoint r = preComp[n >> 1];
+                    if (wi < 0)
+                    {
+                        r = r.Negate();
+                    }
+
+                    q = q.TwicePlus(r);
+                }
+            }
+
+            // Set PreCompInfo in ECPoint, such that it is available for next
+            // multiplication.
+            wnafPreCompInfo.SetPreComp(preComp);
+            wnafPreCompInfo.SetTwiceP(twiceP);
+            p.SetPreCompInfo(wnafPreCompInfo);
+            return q;
+        }
+    }
 }