Don't like ads? PRO users don't see any ads ;-)
Public Pastes
Guest

ECM

By: a guest on Nov 3rd, 2012  |  syntax: C++  |  size: 7.79 KB  |  hits: 56  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. #include <iostream>
  2. #include <mpir.h>
  3. #include <vector>
  4.  
  5. using namespace std;
  6.  
  7. void NextValues(mpz_t X, mpz_t Z, mpz_t k, mpz_t n, mpz_t a, mpz_t b);
  8. void x_sub_2i(mpz_t ret, mpz_t r, mpz_t s, mpz_t a, mpz_t b, mpz_t n);
  9. void z_sub_2i(mpz_t ret, mpz_t r, mpz_t s, mpz_t a, mpz_t b, mpz_t n);
  10. void x_sub_2i_plus_1(mpz_t ret, mpz_t r, mpz_t s, mpz_t u, mpz_t v, mpz_t a, mpz_t b, mpz_t n, mpz_t Z);
  11. void z_sub_2i_plus_1(mpz_t ret, mpz_t r, mpz_t s, mpz_t u, mpz_t v, mpz_t n, mpz_t X);
  12. void ECM(mpz_t n, mpz_t X, mpz_t Y, mpz_t a, mpz_t MAX);
  13.  
  14. int main(){
  15.         mpz_t n;
  16.         mpz_t X;
  17.         mpz_t Y;
  18.         mpz_t a;
  19.         mpz_t MAX;
  20.         mpz_t rn;
  21.         mpz_init(n);
  22.         mpz_init(X);
  23.         mpz_init(Y);
  24.         mpz_init(a);
  25.         mpz_init(MAX);
  26.         gmp_randstate_t r_state;
  27.        
  28.         char input[512];
  29.     unsigned long int seed;
  30.    
  31.         cout << "ECM\n";
  32.         cout << "n = ";
  33.         cin >> input;
  34.         mpz_set_str(n, input, 10);
  35.         cout << "Random Seed: ";
  36.         cin >> seed;
  37.         cout << "max = ";
  38.         cin >> input;
  39.         mpz_set_str(MAX, input, 10);
  40.        
  41.         gmp_randinit_default (r_state);
  42.     gmp_randseed_ui(r_state, seed);
  43.         mpz_init(rn);
  44.  
  45.         while(1){
  46.                 mpz_urandomm(rn,r_state,n);
  47.                 mpz_set(X, rn);
  48.                 mpz_urandomm(rn,r_state,n);
  49.                 mpz_set(Y, rn);
  50.                 mpz_urandomm(rn,r_state,n);
  51.                 mpz_set(a, rn);
  52.                 ECM(n, X, Y, a, MAX);
  53.         }
  54.        
  55.         return 0;
  56. }
  57.  
  58. void ECM(mpz_t n, mpz_t X, mpz_t Y, mpz_t a, mpz_t MAX){
  59.         mpz_t b;
  60.         mpz_t Y2;
  61.         mpz_t X3;
  62.         mpz_t aX;
  63.         mpz_t a3;
  64.         mpz_t b2;
  65.         mpz_t g;
  66.         mpz_t Z;
  67.         mpz_t k;
  68.         mpz_t temp1;
  69.         mpz_t temp2;
  70.         mpz_init(b);
  71.         mpz_init(Y2);
  72.         mpz_init(X3);
  73.         mpz_init(aX);
  74.         mpz_init(a3);
  75.         mpz_init(b2);
  76.         mpz_init(g);
  77.         mpz_init(Z);
  78.         mpz_init(k);
  79.         mpz_init(temp1);
  80.         mpz_init(temp2);
  81.         unsigned long long i;
  82.  
  83.         //b = Y^2 - X^3 - aX mod n
  84.         mpz_mul(Y2, Y, Y);
  85.         mpz_pow_ui(X3, X, 3);
  86.         mpz_mul(aX, a, X);
  87.         mpz_sub(temp1, Y2, X3);
  88.         mpz_sub(temp1, temp1, aX);
  89.         mpz_mod(b, temp1, n);
  90.  
  91.         //g = gcd(4a^3 + 27b^2, n)
  92.         mpz_pow_ui(a3, a, 3);
  93.         mpz_mul(b2, b, b);
  94.         mpz_mul_si(temp1, a3, 4);
  95.         mpz_mul_si(temp2, b2, 27);
  96.         mpz_add(temp1, temp1, temp2);
  97.         mpz_gcd(g, temp1, n);
  98.  
  99.         if(mpz_cmp_ui(g, 1) != 0){
  100.                 mpz_out_str(stdout, 10, g);
  101.                 cout << "\n\a";
  102.                 system("pause");
  103.                 return;
  104.         }
  105.  
  106.         mpz_set_si(Z, 1);
  107.         mpz_set_si(k, 2);
  108.  
  109.         while(mpz_cmp(k, MAX) <= 0){
  110.                 for(i = 1; i <= 10; i++){
  111.                         NextValues(X, Z, k, n, a, b);
  112.                         mpz_add_ui(k, k, 1);
  113.                 }
  114.  
  115.                 mpz_gcd(g, Z, n);
  116.                 if(mpz_cmp_ui(g, 1) != 0){
  117.                         mpz_out_str(stdout, 10, g);
  118.                         cout << "\n\a";
  119.                         system("pause");
  120.                         return;
  121.                 }
  122.         }
  123.  
  124.         mpz_clear(b);
  125.         mpz_clear(Y2);
  126.         mpz_clear(X3);
  127.         mpz_clear(aX);
  128.         mpz_clear(a3);
  129.         mpz_clear(b2);
  130.         mpz_clear(g);
  131.         mpz_clear(Z);
  132.         mpz_clear(k);
  133.         mpz_clear(temp1);
  134.         mpz_clear(temp2);
  135. }
  136.  
  137. void NextValues(mpz_t X, mpz_t Z, mpz_t k, mpz_t n, mpz_t a, mpz_t b){
  138.         unsigned int i = 0;
  139.         unsigned int length;
  140.         vector<int> C;
  141.         mpz_t ik;
  142.         mpz_init(ik);
  143.         mpz_set(ik, k); //ik is used to preserve k
  144.         mpz_t temp1;
  145.         mpz_init(temp1);
  146.         mpz_t X1;
  147.         mpz_t Z1;
  148.         mpz_t X2;
  149.         mpz_t Z2;
  150.         mpz_init(X1);
  151.         mpz_init(Z1);
  152.         mpz_init(X2);
  153.         mpz_init(Z2);
  154.         mpz_t U1;
  155.         mpz_t U2;
  156.         mpz_init(U1);
  157.         mpz_init(U2);
  158.         mpz_t temp;
  159.         mpz_init(temp);
  160.  
  161.         //binary expansion of ik
  162.         while(mpz_cmp_ui(ik, 0) > 0){
  163.                 mpz_mod_ui(temp1, ik, 2);
  164.                 if(mpz_cmp_ui(temp1, 0) == 0){
  165.                         C.push_back(0);
  166.                         mpz_div_ui(ik, ik, 2);
  167.                 }
  168.                 else
  169.                 {
  170.                         C.push_back(1);
  171.                         mpz_sub_ui(ik, ik, 1);
  172.                         mpz_div_ui(ik, ik, 2);
  173.                 }
  174.         }
  175.  
  176.         length = C.size();
  177.  
  178.         mpz_set(X1, X);
  179.         mpz_set(Z1, Z);
  180.        
  181.         x_sub_2i(X2, X, Z, a, b, n);
  182.         z_sub_2i(Z2, X, Z, a, b, n);
  183.  
  184.         for(i = length - 1; i >= 1; i--){
  185.                 x_sub_2i_plus_1(U1, X1, Z1, X2, Z2, a, b, n, Z);
  186.                 z_sub_2i_plus_1(U2, X1, Z1, X2, Z2, n, X);
  187.  
  188.                 if(C[i] == 0){
  189.                         x_sub_2i(temp, X1, Z1, a, b, n);
  190.                         z_sub_2i(Z1, X1, Z1, a, b, n);
  191.                         mpz_set(X1, temp);
  192.                         mpz_set(X2, U1);
  193.                         mpz_set(Z2, U2);
  194.                 }
  195.                 else
  196.                 {
  197.                         x_sub_2i(temp, X2, Z2, a, b, n);
  198.                         z_sub_2i(Z2, X2, Z2, a, b, n);
  199.                         mpz_set(X2, temp);
  200.                         mpz_set(X1, U1);
  201.                         mpz_set(Z1, U2);
  202.                 }
  203.         }
  204.        
  205.         C.erase(C.begin(), C.end());
  206.         mpz_set(X, X1);
  207.         mpz_set(Z, Z1);
  208.  
  209.         mpz_clear(ik);
  210.         mpz_clear(temp1);
  211.         mpz_clear(X1);
  212.         mpz_clear(Z1);
  213.         mpz_clear(X2);
  214.         mpz_clear(Z2);
  215.         mpz_clear(U1);
  216.         mpz_clear(U2);
  217.         mpz_clear(temp);
  218.  
  219.         return;
  220. }
  221.  
  222. void x_sub_2i(mpz_t ret, mpz_t r, mpz_t s, mpz_t a, mpz_t b, mpz_t n){
  223.         mpz_t term;
  224.         mpz_t r2;
  225.         mpz_t s2;
  226.         mpz_t value;
  227.         mpz_t term2;
  228.         mpz_t s3;
  229.         mpz_t temp;
  230.         mpz_init(term);
  231.         mpz_init(r2);
  232.         mpz_init(s2);
  233.         mpz_init(value);
  234.         mpz_init(term2);
  235.         mpz_init(s3);
  236.         mpz_init(temp);
  237.  
  238.         //term = r^2 - a * s^2 MOD n
  239.         mpz_mul(r2, r, r);
  240.         mpz_mul(s2, s, s);
  241.         mpz_mul(temp, a, s2);
  242.         mpz_sub(temp, r2, temp);
  243.         mpz_mod(term, temp, n);
  244.  
  245.         //value = term^2 - 8 * b * r * s^3 MOD n
  246.         mpz_mul(term2, term, term);
  247.         mpz_pow_ui(s3, s, 3);
  248.         mpz_mul_si(temp, b, 8);
  249.         mpz_mul(temp, temp, r);
  250.         mpz_mul(temp, temp, s3);
  251.         mpz_sub(temp, term2, temp);
  252.         mpz_mod(value, temp, n);
  253.  
  254.         mpz_set(ret, value);
  255.  
  256.         mpz_clear(term);
  257.         mpz_clear(r2);
  258.         mpz_clear(s2);
  259.         mpz_clear(value);
  260.         mpz_clear(term2);
  261.         mpz_clear(s3);
  262.         mpz_clear(temp);
  263.         return;
  264. }
  265.  
  266. void z_sub_2i(mpz_t ret, mpz_t r, mpz_t s, mpz_t a, mpz_t b, mpz_t n){
  267.         mpz_t term;
  268.         mpz_t r3;
  269.         mpz_t s2;
  270.         mpz_t s3;
  271.         mpz_t value;
  272.         mpz_t temp;
  273.         mpz_init(term);
  274.         mpz_init(r3);
  275.         mpz_init(s2);
  276.         mpz_init(s3);
  277.         mpz_init(value);
  278.         mpz_init(temp);
  279.         mpz_t temp2;
  280.         mpz_init(temp2);
  281.  
  282.         //term = r^3 + a * r * s^2 + b * s^3 MOD n
  283.         mpz_pow_ui(r3, r, 3);
  284.         mpz_pow_ui(s2, s, 2);
  285.         mpz_pow_ui(s3, s, 3);
  286.         mpz_mul(temp, a, r);
  287.         mpz_mul(temp, temp, s2);
  288.         mpz_mul(temp2, b, s3);
  289.         mpz_add(temp, r3, temp);
  290.         mpz_add(temp, temp, temp2);
  291.         mpz_mod(term, temp, n);
  292.  
  293.         //value = 4 * s * term MOD n
  294.         mpz_mul_ui(temp, s, 4);
  295.         mpz_mul(temp, temp, term);
  296.         mpz_mod(value, temp, n);
  297.         mpz_set(ret, value);
  298.  
  299.         mpz_clear(term);
  300.         mpz_clear(r3);
  301.         mpz_clear(s2);
  302.         mpz_clear(s3);
  303.         mpz_clear(value);
  304.         mpz_clear(temp);
  305.         mpz_clear(temp2);
  306.         return;
  307. }
  308.  
  309. void x_sub_2i_plus_1(mpz_t ret, mpz_t r, mpz_t s, mpz_t u, mpz_t v, mpz_t a, mpz_t b, mpz_t n, mpz_t Z){
  310.         mpz_t term1;
  311.         mpz_t term2;
  312.         mpz_t term12;
  313.         mpz_t value;
  314.         mpz_t ru;
  315.         mpz_t asv;
  316.         mpz_t bsv;
  317.         mpz_t rv;
  318.         mpz_t su;
  319.         mpz_t temp1;
  320.         mpz_t temp2;
  321.         mpz_t temp3;
  322.         mpz_t temp4;
  323.         mpz_init(term1);
  324.         mpz_init(term2);
  325.         mpz_init(term12);
  326.         mpz_init(value);
  327.         mpz_init(ru);
  328.         mpz_init(asv);
  329.         mpz_init(bsv);
  330.         mpz_init(rv);
  331.         mpz_init(su);
  332.         mpz_init(temp1);
  333.         mpz_init(temp2);
  334.         mpz_init(temp3);
  335.         mpz_init(temp4);
  336.  
  337.         //term1 = r * u - a * s * v MOD n
  338.         mpz_mul(ru, r, u);
  339.         mpz_mul(asv, a, s);
  340.         mpz_mul(asv, asv, v);
  341.         mpz_sub(temp1, ru, asv);
  342.         mpz_mod(term1, temp1, n);
  343.  
  344.         //term2 = b * s * v * (r * v + s * u) MOD n
  345.         mpz_mul(bsv, b, s);
  346.         mpz_mul(bsv, bsv, v);
  347.         mpz_mul(rv, r, v);
  348.         mpz_mul(su, s, u);
  349.         mpz_add(temp1, rv, su);
  350.         mpz_mul(temp1, bsv, temp1);
  351.         mpz_mod(term2, temp1, n);
  352.  
  353.         //value = Z * ( term1^2 - 4 * term2 ) MOD n
  354.         mpz_pow_ui(term12, term1, 2);
  355.         mpz_mul_ui(temp1, term2, 4);
  356.         mpz_sub(temp1, term12, temp1);
  357.         mpz_mul(temp1, Z, temp1);
  358.         mpz_mod(value, temp1, n);
  359.         mpz_set(ret, value);
  360.  
  361.         mpz_clear(term1);
  362.         mpz_clear(term2);
  363.         mpz_clear(term12);
  364.         mpz_clear(value);
  365.         mpz_clear(ru);
  366.         mpz_clear(asv);
  367.         mpz_clear(bsv);
  368.         mpz_clear(rv);
  369.         mpz_clear(su);
  370.         mpz_clear(temp1);
  371.         mpz_clear(temp2);
  372.         mpz_clear(temp3);
  373.         mpz_clear(temp4);
  374.         return;
  375. }
  376.  
  377. void z_sub_2i_plus_1(mpz_t ret, mpz_t r, mpz_t s, mpz_t u, mpz_t v, mpz_t n, mpz_t X){
  378.         mpz_t term;
  379.         mpz_t us;
  380.         mpz_t rv;
  381.         mpz_t term2;
  382.         mpz_t value;
  383.         mpz_t temp;
  384.         mpz_init(term);
  385.         mpz_init(us);
  386.         mpz_init(rv);
  387.         mpz_init(term2);
  388.         mpz_init(value);
  389.         mpz_init(temp);
  390.        
  391.         //term = u * s - r * v MOD n
  392.         mpz_mul(us, u, s);
  393.         mpz_mul(rv, r, v);
  394.         mpz_sub(temp, us, rv);
  395.         mpz_mod(term, temp, n);
  396.  
  397.         //value = X * term^2 MOD n
  398.         mpz_pow_ui(term2, term, 2);
  399.         mpz_mul(temp, term2, X);
  400.         mpz_mod(value, temp, n);
  401.         mpz_set(ret, value);
  402.  
  403.         mpz_clear(term);
  404.         mpz_clear(us);
  405.         mpz_clear(rv);
  406.         mpz_clear(term2);
  407.         mpz_clear(value);
  408.         mpz_clear(temp);
  409.         return;
  410. }
create a new version of this paste RAW Paste Data