The template parameter is the type of the base and the results. AbstractGroup is, as it's name suggests, abstract. This class will later be instantiated and inherited by a child-class that will implement group operators (such as addition). The usuall use of it will be through a MultiplicativeGroup, and so it will indeed exponentiate - rather then multiply as suggested by it's name.
view code
//The names of the functions used might be confusing. The usuall use of it will be through a // MultiplicativeGroup, so the Add or Accumulat functions usually multiply, and Double multiplies // the base by itself. template <class T> void AbstractGroup<T>::SimultaneousMultiply(T *results, const T &base, const Integer *expBegin, unsigned int expCount) const { std::vector<std::vector<Element> > buckets(expCount); std::vector<WindowSlider> exponents; exponents.reserve(expCount); unsigned int i; for (i=0; i<expCount; i++) { assert(expBegin->NotNegative()); exponents.push_back(WindowSlider(*expBegin++, InversionIsFast(), 0)); exponents[i].FindNextWindow(); buckets[i].resize(1<<(exponents[i].windowSize-1), Identity()); } unsigned int expBitPosition = 0; Element g = base; bool notDone = true; while (notDone) { notDone = false; for (i=0; i<expCount; i++) { if (!exponents[i].finished && expBitPosition == exponents[i].windowBegin) { Element &bucket = buckets[i][exponents[i].expWindow/2]; if (exponents[i].negateNext) Accumulate(bucket, Inverse(g)); else Accumulate(bucket, g); exponents[i].FindNextWindow(); } notDone = notDone || !exponents[i].finished; } if (notDone) { g = Double(g); expBitPosition++; } } for (i=0; i<expCount; i++) { Element &r = *results++; r = buckets[i][buckets[i].size()-1]; if (buckets[i].size() > 1) { for (int j = buckets[i].size()-2; j >= 1; j--) { Accumulate(buckets[i][j], buckets[i][j+1]); Accumulate(r, buckets[i][j]); } Accumulate(buckets[i][0], buckets[i][1]); r = Add(Double(r), buckets[i][0]); } } }hide code
Threading Intentions
Threading seems to have a bound potential in this case. Most of the work in the function, as seen above, is caused by multiplying B with itself n times. These multiplications cannot be neither threaded nor skipped - there is no way to calculate B4 that will be simpler then calculating B2 first. So threading only tries to speed-up the other n/(k+1)+2k-1 multiplications - of filling the buckets and retrieving the results out of them. Also, retrival of the result from the buckets is a serial process.For that reason, beyond the attempts made to create a speedup for dual-core processors, some of the threading attempts were initially directed towards Multi-Processor environments, allowing extra work to take place while ballancing the threads in hope that the overall result would be faster. Trends in todays CPU development seems to be going in the direction of Dual-Processor systems.
back to top
Threading the function - version 1
The idea in this version was to separate the function into two different threads with very little dependency between them, saving work for the main thread but creating alot of extra work altogether. Each thread calculates the exponentiation of the base by a part of the exponent. When all the threads' work is done the main thread multiplies the results from both threads to get the final result. The exponents for each thread are chosen such that: exp1 + exp2 = exponent. So baseexp1*baseexp2 = baseexponent.One thread will handle thebeginning of the exponent, the other will handle the end. The devision point is set as such that the main thread will have the least possible amount of work, but will still have more work then the helper one. This is to make sure that the main thread won't usually have to wait for the helper thread, thus keeping the main idea of minimum dependency between threads.
Say we devide the exponent, of size n, into two parts - the first holds the first j bits of the exponent - and the secound holds j zeros in it's lower-bits, and the following n-j bits extracted from the original exponent. Say that the first exponent is calculated using window length of q (optimal for j random bits), and the secound is calculated using window length of l (optimal for n-j bits).
The first exponent will calculate using: j + j/(q+1) + 2q -1 multiplications.
The secound exponent will calculate using: n + (n-j)/(l+1) + 2l -1 multiplications.
So we're looking for such j that will give us:
j + j/(q+1) + 2q -1 ~ n + (n-j)/(l+1) + 2l -1
with the right expression a little larger then the left one. This expression can be simplified to:
j/(q+1) + 2q = n-j + (n-j)/(l+1) + 2l
We may assume that 2l is relatively neglegable (since j will surely be far closer to 1.0 then to 0.5). To make calculation easier, and because we want the right expression to be a little larger rather then identical - we shall calculate j as such:
j/(q+1) = (n-j)*A - 2q
j = (n*A - 2q) / (A + 1/(k+1))
A is a constant, A~1 and it's value will determine the margin between the work done by the two threads. Calculating k may seem like a problem, because j is unknown and k is dependant on it. However - n > j > 0.5n. So if we mark k to be the optimal window size for calculating an exponent of length n (n is known, and therefore so is k), q can be only either k or k-1.
J is therefore calculated in two iterations - In the first k is used, then q is derived from the resulted j - and then j is calculated from the new value of q (the final windowsize used will be later calculated from the final value of j).
Summery:
The following are average workloads on threads in terms of group multiplication. To give an idea of the numbers in the summery - average values for a random exponent of size n = 1024bit are shown in brackets.original work: n + n/(k + 1) + 2k - 1 {1233}
main thread: n + (n-j)/(l+1) + 2l -1 {1070}
secound thread:j + j/(q+1) + 2q -1 {1044}
back to top
Threading the function - version 2
The idea in this version was to devide the work into producer and consumer threads. No extra work is created except for unavoidable threading overhead. The produced/consumed values are the values of B mensioned earlier. One thread constantly multiplies B by itself and sends it's value to the other thread whenever a beginning of a window is encountered. The other thread multiplies those values into the appropriate buckets. When the work is done, the main thread calculates the result from the buckets' values.Variant - version 2.1
As we've seen before - the more buckets that exponantiation uses - the less work it will have filling up the buckets, but the more it will have multiplying these buckets into a result. The longer the exponent is - the more buckets the program uses. However, when filling the buckets is done in another thread we can use smaller windows to cause the helper thread to work harder - and the main thread (which calculates the result) to work less. This can decrease overall time because the helper thread has far less work then the main one (it does one multiplication for any (k+1) the main thread does). However, this does create more overall-work so this variant is anticipated to work better for Multi-processor environment in which the threads don't battle for CPU time then on HT machine.I've fixed the window-size on 3, that way the main thred still works about 4 times as much as the main thread - and post-computation only takes 7 multiplications. This value was also chosen via experimentation on different values.
Summery
The following are average workloads on threads in terms of group multiplication. To give an idea of the numbers in the summery - average values for a random exponent of size n = 1024bit are shown in brackets.original work: n + n / (k + 1) + 2k - 1 {1233}
main thread - v2: n + 2k -1 {1087}
secound thread - v2:n / (k + 1) {146}
main thread - v2.1: n + 7 {1031}
secound thread - v2.1:n / 4 {256}
back to top
Implementing Threading
Threading was inserted through another layer of inheritation. The class ThreadedMulGroup inherits from AbstractGroup and has it's own implementation for the function SimultaniousExponentiate. Inheriting from ThreadedMulGroup rather then AbstractGroup makes the inherited group use the threaded version. The threads used for exponentiation are implemented inside the class MulThread (another child of HelperThread).Note: An AbstractRing carrys out all it's exponensiations through it's member MultiplicativeGroup. MultiplicativeGroup is inherited from ThreadedeMulGroup, and therefore an AbstractRing uses threading when exponentiating even if it's not inherited from ThreadedMulGroup. An AbstructRing also contains functions and constructors handeling the HelperThread attached to it's member MultiplicativeGroup. This is also true for classes inheriting from AbstractRing.
Using Modular-Root Thread along with Multiplicative threads
Since I only have two CPUs, threading the multiplications done while calculating the two exponentiations will only result in unnecessary threading overhead. These multiplications are therefore done without threading. This was done by setting the MulThread of the ModularArithmetic objects created for these exponentiations to 0. Multiplications done outside the ModularRoot function will still be threaded. back to topResults
Four versions of the threading mechanisms were tested using a new benchmark made for the cause. The versions are: Non-threaded, Element-threaded (threading version 2), DP variant for element threading (threading version 2.1) and exponent-threaded (threading version 1). The benchmark creates modulus multiplications, with a fixed modulus but random bases and exponents. The modulus is the same size as the bases.Each version was run four times. Each four conclusive tests run were of the four different versions, to prevent as much as possible other CPU interference from dammaging the results. Results shown here are the average millisecounds required per operation. In brackets are the stdev values between the four runs.
Following tests examined threading efficiency over Elliptic-Curve multiplications over GF(p). These multiplications use the same SimultaniousMultiply function. All These Test use the same ECP(155) field used in crypto++'s benchmark.
Integer Exponentiation results
| Test | HT results | DP results | |||||||
|---|---|---|---|---|---|---|---|---|---|
| base size |
exponent size |
Non | elements | elements (DP variant) |
exponent | Non | elements | elements (DP variant) |
exponent |
| 512bit | 256bit | 0.69 (0.01) | 0.78 (0) | 0.81 (0.01) | 1.14 (0) | 0.85 (0) | 0.91 (0) | 0.93 (0.01) | 1.16 (0.07) |
| 512bit | 512bit | 1.33 (0.04) | 1.46 (0) | 1.56 (0.01) | 1.97 (0.01) | 1.59 (0) | 1.7 (0.01) | 1.76 (0.01) | 1.92 (0.07) |
| 1024bit | 128bit | 2.19 (0.04) | 2.15 (0.05) | 2.13 (0.01) | 3.69 (0.02) | 2.72 (0.01) | 2.59 (0.19) | 2.45 (0.03) | 3.77 (0.19) |
| 1024bit | 256bit | 4.3 (0) | 4.14 (0.09) | 4.16 (0.02) | 6.61 (0.03) | 5.4 (0.05) | 4.85 (0.31) | 4.65 (0.02) | 6.53 (0.42) |
| 1024bit | 512bit | 8.39 (0.01) | 8.29 (0.25) | 8.43 (0.02) | 12.4 (0.01) | 10.49 (0.06) | 9.99 (0.83) | 9.39 (0.05) | 11.77 (0.68) |
| 1024bit | 1024bit | 16.45 (0.07) | 16.05 (0.4) | 16.6 (0.04) | 23.67 (0.04) | 20.54 (0.06) | 19.48 (1.53) | 18.51 (0.05) | 22.54 (1.59) |
| 2048bit | 2048bit | 50.94 (0.18) | 49.23 (0.34) | 52.5 (1.71) | 72.35 (1.85) | 63.96 (0.17) | 56.65 (0.38) | 56.65 (0.62) | 65.65 (3.24) |
| 3072bit | 3072bit | 229.87 (0.38) | 222.05 (0.37) | 234.09 (1.73) | 322.74 (4.9) | 289.35 (0.99) | 263.95 (1.45) | 264.06 (5.15) | 299.26 (4.67) |
Elliptic curve multiplication by constant results
| HT results | DP results | |||||||
| constant size |
Non | elements | elements (DP variant) |
exponent | Non | elements | elements (DP variant) |
exponent | 8bit | 0.1 (0.02) | 0.12 (0.02) | 0.11 (0.03) | 0.09 (0.05) | 0.13 (0.03) | 0.13 (0.04) | 0.12 (0.03) | 0.14 (0.02) |
|---|---|---|---|---|---|---|---|---|
| 16bit | 0.47 (0) | 0.47 (0.01) | 0.48 (0.01) | 0.4 (0.19) | 0.59 (0.01) | 0.6 (0.01) | 0.59 (0) | 0.6 (0.01) |
| 62bit | 1.79 (0.02) | 1.86 (0.02) | 1.84 (0) | 2.77 (1.36) | 2.27 (0) | 2.32 (0.02) | 2.31 (0.01) | 3.71 (0.13) |
| 124bit | 3.53 (0.03) | 3.65 (0.04) | 3.58 (0.03) | 4.87 (2.4) | 4.47 (0.03) | 4.53 (0) | 4.51 (0.04) | 6.25 (0.31) |
| 248bit | 6.94 (0.03) | 7.07 (0.04) | 7.08 (0.04) | 8.75 (4.31) | 8.72 (0.06) | 8.8 (0.04) | 8.82 (0.03) | 11.06 (0.86) |
Conclusions
Exponent-threading creates too much overhead work, and produces no speedup over the non-threaded version. Element-threading does produce a measurable speedup over integer modulus exponentiation, and the variant meant for a Multiprocessor system does work better for a DP machine then the original version.Due to the results, I have cancelled Elliptic Curve Engine's use of the threaded multiplication.
back to top
final program
Threading profiles
These results were retrieved from the Intel Thread Profiler (using Windows API), on a Hyper-Threaded pentium4 machine. These are the profiles of the critical path of the program when calculating an RSA 2048 decryption (on the right) and a 1024bit base & exponent exponentiation (on the left). It is clear that RSA calculus has much more ballanced threads, and indeed we saw it creates a larger speedup.
back to top
Threading checks
These results were retrieved from the Intel Thread Checker. This tool is used to check, amongst others, that different threads will not access the same data without blocking each other. Using the HelperThread mechanism helps determining this from a software point of view. All the acess to "shared" data is done through this class, that when used properly provides the necessary thread safety. The check was done while executing a benchmark of the final program. Helper Threads are sometimes not needed for a long period of time, which creates the wornings of long wait times of threads. Since the main thread is the only thread in charge of pther Helper Threads, it dies holding related sync objects. This is also not a prblem.
back to top
Speedup results
These results were retrieved when running the original against the final version, on a P4 and a DP computers respectively. Each test was run for five seconds, each benchmark was run 10 times. The results shown and speedup calculated is the average time (in milliseconds) per action. In braces is the stdev calculated between the 10 retrieved results. Only tests that use the SimultaniousMultiply function are shown.| HT machine | DP machine | |||||
|---|---|---|---|---|---|---|
| Test | Original | final | speedup | Original | final (DP variant) |
speedup |
| RSA 1024 Encryption | 0.12(0) | 0.12(0) | 0% | 0.15(0) | 0.15(0) | 0% |
| RSA 1024 Decryption | 3.49(0.03) | 2.79(0.01) | 20% | 5.68(1.16) | 3.45(0.05) | 39% |
| Rabin 1024 Decryption | 4.16(0.19) | 4.1(0.22) | 1% | 5.17(0.27) | 5.14(0.3) | 1% |
| DLIES 1024 Encryption | 2.82(0.12) | 2.88(0.1) | -2% | 3.55(0.22) | 3.26(0.12) | 8% |
| DLIES 1024 Encryption with precomputation | 3.01(0.14) | 2.99(0.07) | 0% | 3.78(0.17) | 3.74(0.07) | 1% |
| DLIES 1024 Decryption | 9.1(0.13) | 9.1(0.13) | 0% | 11.42(0.21) | 10.66(0.18) | 7% |
| RSA 2048 Encryption | 0.29(0.01) | 0.3(0.01) | -1% | 0.37(0.02) | 0.37(0) | 0% |
| RSA 2048 Decryption | 18.28(0.12) | 14.84(0.09) | 19% | 22.84(0.14) | 17.53(1.88) | 23% |
| Rabin 2048 Decryption | 20.6(1.25) | 20.34(1.24) | 1% | 25.83(1.4) | 24.09(2.87) | 7% |
| DLIES 2048 Encryption | 12.45(0.21) | 11.54(0.27) | 7% | 15.62(0.44) | 13.03(0.17) | 17% |
| DLIES 2048 Encryption with precomputation | 11.79(0.27) | 11.69(0.36) | 1% | 14.75(0.44) | 14.57(0.17) | 1% |
| DLIES 2048 Decryption | 54.83(6.16) | 51.22(0.53) | 7% | 71.57(14.11) | 61.66(10.23) | 14% |
| RSA 1024 Signature | 3.61(0.06) | 2.79(0.03) | 23% | 5.53(1.49) | 3.44(0.08) | 38% |
| RSA 1024 Verification | 0.12(0) | 0.12(0) | 0% | 0.15(0) | 0.15(0) | 0% |
| Rabin 1024 Signature | 4.1(0.06) | 4.05(0.06) | 1% | 5.09(0.05) | 5.07(0.06) | 0% |
| Rabin 1024 Verification | 1.11(0.14) | 1.11(0.14) | 0% | 1.42(0.17) | 1.34(0.11) | 6% |
| RW 1024 Signature | 4.63(0.13) | 3.43(0.05) | 26% | 4.32(0.08) | 4.32(0.03) | 0% |
| RW 1024 Verification | 0.06(0) | 0.06(0) | 0% | 0.08(0.01) | 0.07(0.01) | 4% |
| NR 1024 Signature | 1.44(0.01) | 1.46(0.03) | -2% | 1.81(0.09) | 1.65(0.03) | 9% |
| NR 1024 Signature with precomputation | 0.73(0.03) | 0.71(0.03) | 2% | 0.93(0.06) | 0.89(0.03) | 4% |
| NR 1024 Verification | 1.64(0.07) | 1.68(0.16) | -2% | 2.06(0.08) | 2.08(0.08) | -1% |
| NR 1024 Verification with precomputation | 1.17(0.08) | 1.13(0.12) | 4% | 1.47(0.17) | 1.41(0.17) | 4% |
| DSA 1024 Signature | 1.41(0.04) | 1.44(0.01) | -2% | 1.76(0.04) | 1.65(0.02) | 6% |
| DSA 1024 Signature with precomputation | 0.7(0.02) | 0.71(0.02) | 0% | 0.88(0.01) | 0.88(0) | 0% |
| DSA 1024 Verification | 1.63(0.11) | 1.64(0.12) | -1% | 2.03(0.07) | 1.97(0.09) | 3% |
| DSA 1024 Verification with precomputation | 1.17(0.12) | 1.12(0.13) | 4% | 1.47(0.19) | 1.37(0.14) | 7% |
| ESIGN 1023 Signature | 0.33(0.01) | 0.33(0) | 0% | 0.41(0.01) | 0.41(0.01) | 0% |
| ESIGN 1023 Verification | 0.11(0) | 0.11(0.01) | -1% | 0.14(0) | 0.14(0.01) | -1% |
| ESIGN 1536 Signature | 0.65(0.02) | 0.65(0.02) | 0% | 0.82(0.04) | 0.81(0.02) | 0% |
| ESIGN 1536 Verification | 0.27(0.01) | 0.26(0) | 2% | 0.33(0) | 0.33(0.01) | 0% |
| RSA 2048 Signature | 18.21(0.33) | 14.79(0.09) | 19% | 22.73(0.21) | 17.25(1.14) | 24% |
| RSA 2048 Verification | 0.3(0.03) | 0.29(0) | 2% | 0.36(0.02) | 0.36(0.01) | 0% |
| Rabin 2048 Signature | 20.74(0.6) | 20.47(0.17) | 1% | 25.86(0.32) | 24.33(2.42) | 6% |
| Rabin 2048 Verification | 2.89(0.16) | 2.83(0.15) | 2% | 3.59(0.33) | 3.58(0.29) | 0% |
| RW 2048 Signature | 18.43(0.21) | 18.21(0.12) | 1% | 23(0.21) | 21.34(1.52) | 7% |
| RW 2048 Verification | 0.13(0.01) | 0.13(0.01) | 1% | 0.17(0.01) | 0.17(0.02) | 1% |
| NR 2048 Signature | 6.3(0.16) | 5.85(0.1) | 7% | 7.86(0.09) | 6.63(0.08) | 16% |
| NR 2048 Signature with precomputation | 2.1(0.1) | 2.06(0.03) | 2% | 2.64(0.1) | 2.61(0.08) | 1% |
| NR 2048 Verification | 7.18(0.56) | 7.14(0.19) | 1% | 8.94(0.29) | 8.98(0.33) | 0% |
| NR 2048 Verification with precomputation | 3.4(0.22) | 3.27(0.23) | 4% | 4.27(0.28) | 4.19(0.32) | 2% |
| ESIGN 2046 Signature | 0.77(0.02) | 0.77(0.02) | 0% | 0.96(0.02) | 0.97(0.02) | 0% |
| ESIGN 2046 Verification | 0.28(0) | 0.28(0.01) | 0% | 0.35(0.02) | 0.35(0) | 1% |
| DH 1024 Key-Pair Generation | 1.42(0) | 1.44(0.01) | -1% | 1.78(0.02) | 1.65(0.18) | 7% |
| DH 1024 Key-Pair Generation with precomputation | 1.45(0.01) | 1.51(0.01) | -4% | 1.83(0.02) | 1.91(0.12) | -4% |
| DH 1024 Key Agreement | 2.53(0.11) | 2.53(0.09) | 0% | 3.19(0.1) | 3.01(0.04) | 6% |
| DH 2048 Key-Pair Generation | 6.28(0.21) | 5.89(0.67) | 6% | 7.8(0.06) | 6.49(0.09) | 17% |
| DH 2048 Key-Pair Generation with precomputation | 5.96(0.23) | 5.95(0.67) | 0% | 7.34(0.06) | 7.32(0.07) | 0% |
| DH 2048 Key Agreement | 9.06(0.3) | 8.54(0.12) | 6% | 11.3(0.22) | 10.09(0.89) | 11% |
| MQV 1024 Key-Pair Generation | 1.4(0.02) | 1.41(0.01) | -1% | 1.75(0.03) | 1.6(0.14) | 8% |
| MQV 1024 Key-Pair Generation with precomputation | 0.7(0.03) | 0.69(0.02) | 1% | 0.86(0.06) | 0.87(0.02) | -1% |
| MQV 1024 Key Agreement | 2.63(0.04) | 2.62(0.09) | 0% | 3.31(0.04) | 2.91(0.05) | 12% |
| MQV 2048 Key-Pair Generation | 6.21(0.22) | 5.75(0.13) | 7% | 7.75(0.21) | 6.53(0.22) | 16% |
| MQV 2048 Key-Pair Generation with precomputation | 2.04(0.12) | 2.02(0.07) | 1% | 2.54(0.1) | 2.52(0.04) | 1% |
| MQV 2048 Key Agreement | 11.3(0.12) | 10.46(2.18) | 7% | 14.23(0.41) | 10.82(0.11) | 24% |
Source and executables
final codeexecutable for an HT machine
executable for a DP machine
original executable
dat files
test vectors
Thread activation/parametrization is done in the file helperthread.h. To create the DP variant of element threading (version 2.1) set threadWindowSize() to return 3 (or try a diffrent fixed windowsize). To thread without creating extra work (version 2) make this function return 0.
The .dat files need to be in the same folder as the executable to run the benchmark (use "
back to top
Appendix - optimal window size
As seen before, calculating an exponentiation where the exponent is of size n using window size of k takes: n + n / (k + 1) + 2k - 1 multiplications average (1.5n if k=1), or n + n / k + 2k - 1 in the worst case (2n if k=1).We can clearly see that if for a certain n, k is optimal then for anything less then n, only k or less can be optimal, and k or more for anything over k. So, to find out which k is optimal for any n - I searched for the transition points - meaning I searched for a value of n for which using windowsize k will require the same amout of multiplications as using windowsize k+1. I based my calculations on the worst-case scenario. This means that n can be found from this equation:
n + n / k + 2k - 1 = n + n / (k+1) + 2k+1 - 1
n / k = n / (k+1) + 2k
n(k+1) = nk + k(k+1)2k
n = k(k+1)2k
Transision between k=1 and k=2 is calculated seperately:
n + n/2 +22 - 1 = 2n
n=6
And so the function to retrieve the optimal windowSize for a given exponent length is:
//current program - using worst-case scenario: return (expLen <= 6 ? 1 : (expLen <= 24 ? 2 : (expLen <= 96 ? 3 : (expLen <= 320 ? 4 : (expLen <= 960 ? 5 : (expLen <= 2688 ? 6 : 7)))))); //another possibility could be using average-case scenario: return (expLen <= 17 ? 1 : (expLen <= 47 ? 2 : (expLen <= 159 ? 3 : (expLen <= 479 ? 4 : (expLen <= 1343 ? 5 : (expLen <= 3583 ? 6 : 7))))));back to top