fork download
  1. #include <iostream>
  2. #include <stdlib.h>
  3. #include <math.h>
  4.  
  5. // Нахождение GCD
  6. long long gcd(long long a, long long b){
  7. if(b == 0)
  8. return a;
  9. return gcd(b, a % b);
  10. }
  11.  
  12. // Двоичное умножение по модулю для защиты от переполнений
  13. long long mul(long long a, long long b, long long m){
  14. if(b == 1)
  15. return a;
  16. if(b % 2 == 0){
  17. long long t = mul(a, b/2, m);
  18. return (2 * t) % m;
  19. }
  20. return (mul(a, b - 1, m) + a) % m;
  21. }
  22.  
  23. // Возведение в степень по модулю
  24. long long modpow (long long a, long long b, long long m){
  25. if(b==0)
  26. return 1;
  27. if(b % 2 == 0){
  28. long long t = modpow (a, b/2, m);
  29. return mul(t, t, m) % m;
  30. }
  31. return (mul(modpow (a, b - 1, m), a, m)) % m;
  32. }
  33.  
  34. int jacobi (int nominator, int denominator)
  35. {
  36. int s = 0;
  37. int remainder, power_of_two, t;
  38. do{
  39. remainder = nominator % denominator;
  40. power_of_two = t = 0;
  41. while(remainder % 2 == 0)
  42. {
  43. power_of_two++;
  44. remainder >>= 1;
  45. }
  46. t = remainder;
  47. s = (s + power_of_two * (denominator * denominator - 1) / 8 + (t - 1) * (denominator - 1) / 4) % 2;
  48. if(t == 1)
  49. return (s) ? -1 : 1;
  50. nominator = denominator;
  51. denominator = t;
  52. }while(t >= 3);
  53. }
  54.  
  55. int main(int argc, char* argv[])
  56. {
  57. long n; /*n > 2, нечетное натуральное число, которое необходимо проверить на простоту*/
  58. long k; /* параметр, определяющий точность теста */
  59. long a;
  60. double t;
  61. int flag = 0;
  62.  
  63. std::cin >> n;
  64. std::cin >> k;
  65. t = 1 - 1 / pow(2.0, k);
  66.  
  67. for (int i=1; i <= k; i++)
  68. {
  69. a = rand() % (n - 2) + 2;
  70. long long p = modpow(a, (n - 1) / 2, n);
  71. int j = jacobi (a, n);
  72.  
  73. if (gcd(a, n) != 1)
  74. {
  75. flag = 1;
  76. break;
  77. }
  78. else if ((p - n) != j)
  79. {
  80. std::cout << "[DEBUG]: modpow = " << p << std::endl;
  81. std::cout << "[DEBUG]: jacobi = " << j << std::endl;
  82. flag = 2;
  83. break;
  84. }
  85. }
  86.  
  87.  
  88. if (flag !=0)
  89. {
  90. std::cout <<"composite \n";
  91. }
  92. else
  93. std::cout<<"prime with possibility " << t << "\n";
  94. return 0;
  95. }
Success #stdin #stdout 0s 3468KB
stdin
311 10
stdout
[DEBUG]: modpow = 1
[DEBUG]: jacobi = 1
composite