Botan
1.11.15
|
00001 /** 00002 * (C) Copyright Projet SECRET, INRIA, Rocquencourt 00003 * (C) Bhaskar Biswas and Nicolas Sendrier 00004 * (C) 2014 Jack Lloyd 00005 * 00006 * Botan is released under the Simplified BSD License (see license.txt) 00007 * 00008 */ 00009 00010 #include <botan/mceliece.h> 00011 #include <cmath> 00012 00013 namespace Botan { 00014 00015 namespace { 00016 00017 double binomial(size_t n, size_t k) 00018 { 00019 double x = 1; 00020 00021 for(size_t i = 0; i != k; ++i) 00022 { 00023 x *= n - i; 00024 x /= k -i; 00025 } 00026 00027 return x; 00028 } 00029 00030 double log_binomial(size_t n, size_t k) 00031 { 00032 double x = 0; 00033 00034 for(size_t i = 0; i != k; ++i) 00035 { 00036 x += std::log(n - i); 00037 x -= std::log(k - i); 00038 } 00039 00040 return x / std::log(2); 00041 } 00042 00043 double nb_iter(size_t n, size_t k, size_t w, size_t p, size_t l) 00044 { 00045 double x = 2 * log_binomial(k / 2, p); 00046 x += log_binomial(n - k - l, w - 2 * p); 00047 x = log_binomial(n, w) - x; 00048 return x; 00049 } 00050 00051 double cout_iter(size_t n, size_t k, size_t p, size_t l) 00052 { 00053 double x = binomial(k / 2, p); 00054 const size_t i = static_cast<size_t>(std::log(x) / std::log(2)); 00055 double res = 2 * p * (n - k - l) * ldexp(x * x, -l); 00056 00057 // x <- binomial(k/2,p)*2*(2*l+log[2](binomial(k/2,p))) 00058 x *= 2 * (2 * l + i); 00059 00060 // res <- k*(n-k)/2 + 00061 // binomial(k/2,p)*2*(2*l+log[2](binomial(k/2,p))) + 00062 // 2*p*(n-k-l)*binomial(k/2,p)^2/2^l 00063 res += x + k * ((n - k) / 2.0); 00064 00065 return std::log(res) / std::log(2); // convert to bits 00066 } 00067 00068 double cout_total(size_t n, size_t k, size_t w, size_t p, size_t l) 00069 { 00070 return nb_iter(n, k, w, p, l) + cout_iter(n, k, p, l); 00071 } 00072 00073 double best_wf(size_t n, size_t k, size_t w, size_t p) 00074 { 00075 if(p >= k / 2) 00076 return -1; 00077 00078 double min = cout_total(n, k, w, p, 0); 00079 00080 for(size_t l = 1; l < n - k; ++l) 00081 { 00082 const double lwf = cout_total(n, k, w, p, l); 00083 if(lwf < min) 00084 min = lwf; 00085 else 00086 break; 00087 } 00088 00089 return min; 00090 } 00091 00092 } 00093 00094 size_t mceliece_work_factor(size_t n, size_t k, size_t t) 00095 { 00096 double min = cout_total(n, k, t, 0, 0); // correspond a p=1 00097 for(size_t p = 0; p != t / 2; ++p) 00098 { 00099 double lwf = best_wf(n, k + 1, t, p); 00100 if(lwf < 0) 00101 break; 00102 00103 min = std::min(min, lwf); 00104 } 00105 00106 return min; 00107 } 00108 00109 }