用mathematica写的miller rabin算法

mathematica

mathematica 发表于 2008-11-30 18:47
上面的那个数是一个素数。这个例子给出一个强伪素数，它是前46个素数
的强伪素数。可以看出，最后的结果判 ...

观察下面的结果
{2,1,-1,prime}
{3,0,-1,prime}
{4,0,-1,prime}
{5,0,1,prime}
{6,1,-1,prime}
{7,1,-1,prime}
{8,1,-1,prime}
{9,0,1,prime}
{10,1,-1,prime}
{11,1,-1,prime}
{12,0,1,prime}
{13,1,-1,prime}
{14,composite}
提取
{2,1,-1,prime}
{7,1,-1,prime}
这个意味着
\(x^2\equiv-1\pmod{n}\)
\(y^2\equiv-1\pmod{n}\)
那么
\((x*y)^2\equiv1^2\pmod{n}\)
并且x与y互质,
所以可以利用这个同余等式来分解质因数
分别计算
x*y-1与n的最大公约数
x*y+1与n的最大公约数
然后就能得到质因子
4009582166394996054183064520845468530051881660411325087745062047380032\
1707011962427162231915972197335821631650853581669691452338139171692875\
27980445796800452592031836601

代码如下:

1. Clear["Global`*"];(*Clear all variables*)
   
2. n=803837457453639491257079614341942108138837688287558145837488917522\
   
3. 2974273765333652186502336163960045457915042023603208766569966760987284\
   
4. 0439654082329287387918508691668573282677617710293896977394701670823042\
   
5. 8687109997439976544144845341155872450633409279022275296229414984230688\
   
6. 1685404326457534018329786111298960644845216191652872597534901;
   
7. m=n-1;s=0;
   
8. While[Mod[m,2]==0,m=m/2;s=s+1];
   
9. kk=PowerMod[2*7,m,n]
   
10. ga=GCD[kk-1,n]
    

复制代码

mathematica

1. Clear["Global`*"];
   
2. MR[n0_,a0_]:=Module[{n=n0,a=a0,s,m,t1,k},
   
3.     s=0;m=n-1;While[Mod[m,2]==0,m=m/2;s=s+1];
   
4.     t1=PowerMod[a,m,n];
   
5.     If[t1==1,Return[True]];
   
6.     k=0;While[k<s-1&&t1!=n-1,k=k+1;t1=Mod[t1^2,n]];
   
7.     If[t1==n-1,Return[True],Return[False]]
   
8. ]
   

复制代码

这是我重新写的miller rabin算法的代码,又简化了!

mathematica

mathematica 发表于 2019-2-24 09:33
这是我重新写的miller rabin算法的代码,又简化了!

1. Clear["Global`*"];
   
2. n=2887148238050771212671429597130393991977609459279722700926516024197432\
   
3. 3037991527331163289831446392259419778031109293496555784189494417409338\
   
4. 0561511397999942154241693397290542371100275104208013496673175515285922\
   
5. 6962916775325475044445856101949404200039904432116776619949629539250452\
   
6. 6987193290703735640322737012784538991261203092448414947289768854060249\
   
7. 76768122077071687938121709811322297802059565867;
   
8. (*miller rabin子函数*)
   
9. MR[n0_,a0_]:=Module[{n=n0,a=a0,s,m,t1,k},
   
10.     s=0;m=n-1;While[Mod[m,2]==0,m=m/2;s=s+1];
    
11.     t1=PowerMod[a,m,n];
    
12.     If[t1==1,Return[True]];
    
13.     k=0;While[k<s-1&&t1!=n-1,k=k+1;t1=Mod[t1^2,n]];
    
14.     If[t1==n-1,Return[True],Return[False]]
    
15. ]
    
16. Do[Print[{k,MR[n,k]}],{k,1,307}]
    

复制代码

运行结果

1. {1,True}
   
2. {2,True}
   
3. {3,True}
   
4. {4,True}
   
5. {5,True}
   
6. {6,True}
   
7. {7,True}
   
8. {8,True}
   
9. {9,True}
   
10. {10,True}
    
11. {11,True}
    
12. {12,True}
    
13. {13,True}
    
14. {14,True}
    
15. {15,True}
    
16. {16,True}
    
17. {17,True}
    
18. {18,True}
    
19. {19,True}
    
20. {20,True}
    
21. {21,True}
    
22. {22,True}
    
23. {23,True}
    
24. {24,True}
    
25. {25,True}
    
26. {26,True}
    
27. {27,True}
    
28. {28,True}
    
29. {29,True}
    
30. {30,True}
    
31. {31,True}
    
32. {32,True}
    
33. {33,True}
    
34. {34,True}
    
35. {35,True}
    
36. {36,True}
    
37. {37,True}
    
38. {38,True}
    
39. {39,True}
    
40. {40,True}
    
41. {41,True}
    
42. {42,True}
    
43. {43,True}
    
44. {44,True}
    
45. {45,True}
    
46. {46,True}
    
47. {47,True}
    
48. {48,True}
    
49. {49,True}
    
50. {50,True}
    
51. {51,True}
    
52. {52,True}
    
53. {53,True}
    
54. {54,True}
    
55. {55,True}
    
56. {56,True}
    
57. {57,True}
    
58. {58,True}
    
59. {59,True}
    
60. {60,True}
    
61. {61,True}
    
62. {62,True}
    
63. {63,True}
    
64. {64,True}
    
65. {65,True}
    
66. {66,True}
    
67. {67,True}
    
68. {68,True}
    
69. {69,True}
    
70. {70,True}
    
71. {71,True}
    
72. {72,True}
    
73. {73,True}
    
74. {74,True}
    
75. {75,True}
    
76. {76,True}
    
77. {77,True}
    
78. {78,True}
    
79. {79,True}
    
80. {80,True}
    
81. {81,True}
    
82. {82,True}
    
83. {83,True}
    
84. {84,True}
    
85. {85,True}
    
86. {86,True}
    
87. {87,True}
    
88. {88,True}
    
89. {89,True}
    
90. {90,True}
    
91. {91,True}
    
92. {92,True}
    
93. {93,True}
    
94. {94,True}
    
95. {95,True}
    
96. {96,True}
    
97. {97,True}
    
98. {98,True}
    
99. {99,True}
    
100. {100,True}
     
101. {101,True}
     
102. {102,True}
     
103. {103,True}
     
104. {104,True}
     
105. {105,True}
     
106. {106,True}
     
107. {107,True}
     
108. {108,True}
     
109. {109,True}
     
110. {110,True}
     
111. {111,True}
     
112. {112,True}
     
113. {113,True}
     
114. {114,True}
     
115. {115,True}
     
116. {116,True}
     
117. {117,True}
     
118. {118,True}
     
119. {119,True}
     
120. {120,True}
     
121. {121,True}
     
122. {122,True}
     
123. {123,True}
     
124. {124,True}
     
125. {125,True}
     
126. {126,True}
     
127. {127,True}
     
128. {128,True}
     
129. {129,True}
     
130. {130,True}
     
131. {131,True}
     
132. {132,True}
     
133. {133,True}
     
134. {134,True}
     
135. {135,True}
     
136. {136,True}
     
137. {137,True}
     
138. {138,True}
     
139. {139,True}
     
140. {140,True}
     
141. {141,True}
     
142. {142,True}
     
143. {143,True}
     
144. {144,True}
     
145. {145,True}
     
146. {146,True}
     
147. {147,True}
     
148. {148,True}
     
149. {149,True}
     
150. {150,True}
     
151. {151,True}
     
152. {152,True}
     
153. {153,True}
     
154. {154,True}
     
155. {155,True}
     
156. {156,True}
     
157. {157,True}
     
158. {158,True}
     
159. {159,True}
     
160. {160,True}
     
161. {161,True}
     
162. {162,True}
     
163. {163,True}
     
164. {164,True}
     
165. {165,True}
     
166. {166,True}
     
167. {167,True}
     
168. {168,True}
     
169. {169,True}
     
170. {170,True}
     
171. {171,True}
     
172. {172,True}
     
173. {173,True}
     
174. {174,True}
     
175. {175,True}
     
176. {176,True}
     
177. {177,True}
     
178. {178,True}
     
179. {179,True}
     
180. {180,True}
     
181. {181,True}
     
182. {182,True}
     
183. {183,True}
     
184. {184,True}
     
185. {185,True}
     
186. {186,True}
     
187. {187,True}
     
188. {188,True}
     
189. {189,True}
     
190. {190,True}
     
191. {191,True}
     
192. {192,True}
     
193. {193,True}
     
194. {194,True}
     
195. {195,True}
     
196. {196,True}
     
197. {197,True}
     
198. {198,True}
     
199. {199,True}
     
200. {200,True}
     
201. {201,True}
     
202. {202,True}
     
203. {203,True}
     
204. {204,True}
     
205. {205,True}
     
206. {206,True}
     
207. {207,True}
     
208. {208,True}
     
209. {209,True}
     
210. {210,True}
     
211. {211,True}
     
212. {212,True}
     
213. {213,True}
     
214. {214,True}
     
215. {215,True}
     
216. {216,True}
     
217. {217,True}
     
218. {218,True}
     
219. {219,True}
     
220. {220,True}
     
221. {221,True}
     
222. {222,True}
     
223. {223,True}
     
224. {224,True}
     
225. {225,True}
     
226. {226,True}
     
227. {227,True}
     
228. {228,True}
     
229. {229,True}
     
230. {230,True}
     
231. {231,True}
     
232. {232,True}
     
233. {233,True}
     
234. {234,True}
     
235. {235,True}
     
236. {236,True}
     
237. {237,True}
     
238. {238,True}
     
239. {239,True}
     
240. {240,True}
     
241. {241,True}
     
242. {242,True}
     
243. {243,True}
     
244. {244,True}
     
245. {245,True}
     
246. {246,True}
     
247. {247,True}
     
248. {248,True}
     
249. {249,True}
     
250. {250,True}
     
251. {251,True}
     
252. {252,True}
     
253. {253,True}
     
254. {254,True}
     
255. {255,True}
     
256. {256,True}
     
257. {257,True}
     
258. {258,True}
     
259. {259,True}
     
260. {260,True}
     
261. {261,True}
     
262. {262,True}
     
263. {263,True}
     
264. {264,True}
     
265. {265,True}
     
266. {266,True}
     
267. {267,True}
     
268. {268,True}
     
269. {269,True}
     
270. {270,True}
     
271. {271,True}
     
272. {272,True}
     
273. {273,True}
     
274. {274,True}
     
275. {275,True}
     
276. {276,True}
     
277. {277,True}
     
278. {278,True}
     
279. {279,True}
     
280. {280,True}
     
281. {281,True}
     
282. {282,True}
     
283. {283,True}
     
284. {284,True}
     
285. {285,True}
     
286. {286,True}
     
287. {287,True}
     
288. {288,True}
     
289. {289,True}
     
290. {290,True}
     
291. {291,True}
     
292. {292,True}
     
293. {293,True}
     
294. {294,True}
     
295. {295,True}
     
296. {296,True}
     
297. {297,True}
     
298. {298,True}
     
299. {299,True}
     
300. {300,True}
     
301. {301,True}
     
302. {302,True}
     
303. {303,True}
     
304. {304,True}
     
305. {305,True}
     
306. {306,True}
     
307. {307,False}
     

复制代码

mathematica

https://bbs.emath.ac.cn/forum.ph ... 734&fromuid=865
用此处的lucas判定办法,就回得到false的答案!

nyy

好地方 发表于 2008-11-30 20:31
这个伪素数够强，不过越强的伪素数其实越脆弱，略施小计，只需一眨眼的时间就将它分解为下面两个素数的乘积 ...

能说说你是怎么分解的吗？我需要的是方法