Crypto++
|
00001 // algebra.cpp - written and placed in the public domain by Wei Dai 00002 00003 #include "pch.h" 00004 00005 #ifndef CRYPTOPP_ALGEBRA_CPP // SunCC workaround: compiler could cause this file to be included twice 00006 #define CRYPTOPP_ALGEBRA_CPP 00007 00008 #include "algebra.h" 00009 #include "integer.h" 00010 00011 #include <vector> 00012 00013 NAMESPACE_BEGIN(CryptoPP) 00014 00015 template <class T> const T& AbstractGroup<T>::Double(const Element &a) const 00016 { 00017 return Add(a, a); 00018 } 00019 00020 template <class T> const T& AbstractGroup<T>::Subtract(const Element &a, const Element &b) const 00021 { 00022 // make copy of a in case Inverse() overwrites it 00023 Element a1(a); 00024 return Add(a1, Inverse(b)); 00025 } 00026 00027 template <class T> T& AbstractGroup<T>::Accumulate(Element &a, const Element &b) const 00028 { 00029 return a = Add(a, b); 00030 } 00031 00032 template <class T> T& AbstractGroup<T>::Reduce(Element &a, const Element &b) const 00033 { 00034 return a = Subtract(a, b); 00035 } 00036 00037 template <class T> const T& AbstractRing<T>::Square(const Element &a) const 00038 { 00039 return Multiply(a, a); 00040 } 00041 00042 template <class T> const T& AbstractRing<T>::Divide(const Element &a, const Element &b) const 00043 { 00044 // make copy of a in case MultiplicativeInverse() overwrites it 00045 Element a1(a); 00046 return Multiply(a1, MultiplicativeInverse(b)); 00047 } 00048 00049 template <class T> const T& AbstractEuclideanDomain<T>::Mod(const Element &a, const Element &b) const 00050 { 00051 Element q; 00052 DivisionAlgorithm(result, q, a, b); 00053 return result; 00054 } 00055 00056 template <class T> const T& AbstractEuclideanDomain<T>::Gcd(const Element &a, const Element &b) const 00057 { 00058 Element g[3]={b, a}; 00059 unsigned int i0=0, i1=1, i2=2; 00060 00061 while (!Equal(g[i1], this->Identity())) 00062 { 00063 g[i2] = Mod(g[i0], g[i1]); 00064 unsigned int t = i0; i0 = i1; i1 = i2; i2 = t; 00065 } 00066 00067 return result = g[i0]; 00068 } 00069 00070 template <class T> const typename QuotientRing<T>::Element& QuotientRing<T>::MultiplicativeInverse(const Element &a) const 00071 { 00072 Element g[3]={m_modulus, a}; 00073 Element v[3]={m_domain.Identity(), m_domain.MultiplicativeIdentity()}; 00074 Element y; 00075 unsigned int i0=0, i1=1, i2=2; 00076 00077 while (!Equal(g[i1], Identity())) 00078 { 00079 // y = g[i0] / g[i1]; 00080 // g[i2] = g[i0] % g[i1]; 00081 m_domain.DivisionAlgorithm(g[i2], y, g[i0], g[i1]); 00082 // v[i2] = v[i0] - (v[i1] * y); 00083 v[i2] = m_domain.Subtract(v[i0], m_domain.Multiply(v[i1], y)); 00084 unsigned int t = i0; i0 = i1; i1 = i2; i2 = t; 00085 } 00086 00087 return m_domain.IsUnit(g[i0]) ? m_domain.Divide(v[i0], g[i0]) : m_domain.Identity(); 00088 } 00089 00090 template <class T> T AbstractGroup<T>::ScalarMultiply(const Element &base, const Integer &exponent) const 00091 { 00092 Element result; 00093 SimultaneousMultiply(&result, base, &exponent, 1); 00094 return result; 00095 } 00096 00097 template <class T> T AbstractGroup<T>::CascadeScalarMultiply(const Element &x, const Integer &e1, const Element &y, const Integer &e2) const 00098 { 00099 const unsigned expLen = STDMAX(e1.BitCount(), e2.BitCount()); 00100 if (expLen==0) 00101 return Identity(); 00102 00103 const unsigned w = (expLen <= 46 ? 1 : (expLen <= 260 ? 2 : 3)); 00104 const unsigned tableSize = 1<<w; 00105 std::vector<Element> powerTable(tableSize << w); 00106 00107 powerTable[1] = x; 00108 powerTable[tableSize] = y; 00109 if (w==1) 00110 powerTable[3] = Add(x,y); 00111 else 00112 { 00113 powerTable[2] = Double(x); 00114 powerTable[2*tableSize] = Double(y); 00115 00116 unsigned i, j; 00117 00118 for (i=3; i<tableSize; i+=2) 00119 powerTable[i] = Add(powerTable[i-2], powerTable[2]); 00120 for (i=1; i<tableSize; i+=2) 00121 for (j=i+tableSize; j<(tableSize<<w); j+=tableSize) 00122 powerTable[j] = Add(powerTable[j-tableSize], y); 00123 00124 for (i=3*tableSize; i<(tableSize<<w); i+=2*tableSize) 00125 powerTable[i] = Add(powerTable[i-2*tableSize], powerTable[2*tableSize]); 00126 for (i=tableSize; i<(tableSize<<w); i+=2*tableSize) 00127 for (j=i+2; j<i+tableSize; j+=2) 00128 powerTable[j] = Add(powerTable[j-1], x); 00129 } 00130 00131 Element result; 00132 unsigned power1 = 0, power2 = 0, prevPosition = expLen-1; 00133 bool firstTime = true; 00134 00135 for (int i = expLen-1; i>=0; i--) 00136 { 00137 power1 = 2*power1 + e1.GetBit(i); 00138 power2 = 2*power2 + e2.GetBit(i); 00139 00140 if (i==0 || 2*power1 >= tableSize || 2*power2 >= tableSize) 00141 { 00142 unsigned squaresBefore = prevPosition-i; 00143 unsigned squaresAfter = 0; 00144 prevPosition = i; 00145 while ((power1 || power2) && power1%2 == 0 && power2%2==0) 00146 { 00147 power1 /= 2; 00148 power2 /= 2; 00149 squaresBefore--; 00150 squaresAfter++; 00151 } 00152 if (firstTime) 00153 { 00154 result = powerTable[(power2<<w) + power1]; 00155 firstTime = false; 00156 } 00157 else 00158 { 00159 while (squaresBefore--) 00160 result = Double(result); 00161 if (power1 || power2) 00162 Accumulate(result, powerTable[(power2<<w) + power1]); 00163 } 00164 while (squaresAfter--) 00165 result = Double(result); 00166 power1 = power2 = 0; 00167 } 00168 } 00169 return result; 00170 } 00171 00172 template <class Element, class Iterator> Element GeneralCascadeMultiplication(const AbstractGroup<Element> &group, Iterator begin, Iterator end) 00173 { 00174 if (end-begin == 1) 00175 return group.ScalarMultiply(begin->base, begin->exponent); 00176 else if (end-begin == 2) 00177 return group.CascadeScalarMultiply(begin->base, begin->exponent, (begin+1)->base, (begin+1)->exponent); 00178 else 00179 { 00180 Integer q, t; 00181 Iterator last = end; 00182 --last; 00183 00184 std::make_heap(begin, end); 00185 std::pop_heap(begin, end); 00186 00187 while (!!begin->exponent) 00188 { 00189 // last->exponent is largest exponent, begin->exponent is next largest 00190 t = last->exponent; 00191 Integer::Divide(last->exponent, q, t, begin->exponent); 00192 00193 if (q == Integer::One()) 00194 group.Accumulate(begin->base, last->base); // avoid overhead of ScalarMultiply() 00195 else 00196 group.Accumulate(begin->base, group.ScalarMultiply(last->base, q)); 00197 00198 std::push_heap(begin, end); 00199 std::pop_heap(begin, end); 00200 } 00201 00202 return group.ScalarMultiply(last->base, last->exponent); 00203 } 00204 } 00205 00206 struct WindowSlider 00207 { 00208 WindowSlider(const Integer &expIn, bool fastNegate, unsigned int windowSizeIn=0) 00209 : exp(expIn), windowModulus(Integer::One()), windowSize(windowSizeIn), windowBegin(0), fastNegate(fastNegate), firstTime(true), finished(false) 00210 { 00211 if (windowSize == 0) 00212 { 00213 unsigned int expLen = exp.BitCount(); 00214 windowSize = expLen <= 17 ? 1 : (expLen <= 24 ? 2 : (expLen <= 70 ? 3 : (expLen <= 197 ? 4 : (expLen <= 539 ? 5 : (expLen <= 1434 ? 6 : 7))))); 00215 } 00216 windowModulus <<= windowSize; 00217 } 00218 00219 void FindNextWindow() 00220 { 00221 unsigned int expLen = exp.WordCount() * WORD_BITS; 00222 unsigned int skipCount = firstTime ? 0 : windowSize; 00223 firstTime = false; 00224 while (!exp.GetBit(skipCount)) 00225 { 00226 if (skipCount >= expLen) 00227 { 00228 finished = true; 00229 return; 00230 } 00231 skipCount++; 00232 } 00233 00234 exp >>= skipCount; 00235 windowBegin += skipCount; 00236 expWindow = word32(exp % (word(1) << windowSize)); 00237 00238 if (fastNegate && exp.GetBit(windowSize)) 00239 { 00240 negateNext = true; 00241 expWindow = (word32(1) << windowSize) - expWindow; 00242 exp += windowModulus; 00243 } 00244 else 00245 negateNext = false; 00246 } 00247 00248 Integer exp, windowModulus; 00249 unsigned int windowSize, windowBegin; 00250 word32 expWindow; 00251 bool fastNegate, negateNext, firstTime, finished; 00252 }; 00253 00254 template <class T> 00255 void AbstractGroup<T>::SimultaneousMultiply(T *results, const T &base, const Integer *expBegin, unsigned int expCount) const 00256 { 00257 std::vector<std::vector<Element> > buckets(expCount); 00258 std::vector<WindowSlider> exponents; 00259 exponents.reserve(expCount); 00260 unsigned int i; 00261 00262 for (i=0; i<expCount; i++) 00263 { 00264 assert(expBegin->NotNegative()); 00265 exponents.push_back(WindowSlider(*expBegin++, InversionIsFast(), 0)); 00266 exponents[i].FindNextWindow(); 00267 buckets[i].resize(1<<(exponents[i].windowSize-1), Identity()); 00268 } 00269 00270 unsigned int expBitPosition = 0; 00271 Element g = base; 00272 bool notDone = true; 00273 00274 while (notDone) 00275 { 00276 notDone = false; 00277 for (i=0; i<expCount; i++) 00278 { 00279 if (!exponents[i].finished && expBitPosition == exponents[i].windowBegin) 00280 { 00281 Element &bucket = buckets[i][exponents[i].expWindow/2]; 00282 if (exponents[i].negateNext) 00283 Accumulate(bucket, Inverse(g)); 00284 else 00285 Accumulate(bucket, g); 00286 exponents[i].FindNextWindow(); 00287 } 00288 notDone = notDone || !exponents[i].finished; 00289 } 00290 00291 if (notDone) 00292 { 00293 g = Double(g); 00294 expBitPosition++; 00295 } 00296 } 00297 00298 for (i=0; i<expCount; i++) 00299 { 00300 Element &r = *results++; 00301 r = buckets[i][buckets[i].size()-1]; 00302 if (buckets[i].size() > 1) 00303 { 00304 for (int j = (int)buckets[i].size()-2; j >= 1; j--) 00305 { 00306 Accumulate(buckets[i][j], buckets[i][j+1]); 00307 Accumulate(r, buckets[i][j]); 00308 } 00309 Accumulate(buckets[i][0], buckets[i][1]); 00310 r = Add(Double(r), buckets[i][0]); 00311 } 00312 } 00313 } 00314 00315 template <class T> T AbstractRing<T>::Exponentiate(const Element &base, const Integer &exponent) const 00316 { 00317 Element result; 00318 SimultaneousExponentiate(&result, base, &exponent, 1); 00319 return result; 00320 } 00321 00322 template <class T> T AbstractRing<T>::CascadeExponentiate(const Element &x, const Integer &e1, const Element &y, const Integer &e2) const 00323 { 00324 return MultiplicativeGroup().AbstractGroup<T>::CascadeScalarMultiply(x, e1, y, e2); 00325 } 00326 00327 template <class Element, class Iterator> Element GeneralCascadeExponentiation(const AbstractRing<Element> &ring, Iterator begin, Iterator end) 00328 { 00329 return GeneralCascadeMultiplication<Element>(ring.MultiplicativeGroup(), begin, end); 00330 } 00331 00332 template <class T> 00333 void AbstractRing<T>::SimultaneousExponentiate(T *results, const T &base, const Integer *exponents, unsigned int expCount) const 00334 { 00335 MultiplicativeGroup().AbstractGroup<T>::SimultaneousMultiply(results, base, exponents, expCount); 00336 } 00337 00338 NAMESPACE_END 00339 00340 #endif