1use num::integer;
2use rand::prelude::*;
3
4pub fn is_prime(n: u32) -> bool {
5 if n < 2 {
6 return false;
7 }
8
9 if n == 2 {
10 return true;
11 }
12
13 if n.is_multiple_of(2) {
14 return false;
15 }
16
17 for i in 3..(integer::sqrt(n) + 1) {
18 if n.is_multiple_of(i) {
19 return false;
20 }
21 }
22
23 true
24}
25
26pub fn coprime(n: u32) -> u32 {
27 let mut rng = rand::rng();
28
29 loop {
30 let a = rng.random_range(2..n - 1);
31 if gcd(n, a) == 1 {
32 return a;
33 }
34 }
35}
36
37pub fn base_exp(n: u32) -> Option<(u32, u32)> {
38 let s = format!("{:b}", n).chars().count() as u32;
39
40 for i in (2..s).rev() {
41 let a = (n as f64).powf(1.0 / (i as f64)).round() as u32;
42 if a.pow(i) == n {
43 return Some((a, i));
44 }
45 }
46
47 None
48}
49
50pub fn gcd(a: u32, b: u32) -> u32 {
51 integer::gcd(a, b)
52}
53
54pub fn modexp(a: u32, r: u32, n: u32) -> u32 {
55 if a == 0 {
56 return 0;
57 }
58
59 if r == 0 {
60 return 1;
61 }
62
63 let mut p = a;
64 for _ in 1..r {
65 p = (p * a) % n
66 }
67
68 p
69}
70
71pub fn modexp2(a: u32, j: u32, n: u32) -> u32 {
72 if a == 0 {
73 return 0;
74 }
75
76 if j == 0 {
77 return a % n;
78 }
79
80 let mut p = a;
81 for _ in 0..j {
82 p = (p * p) % n
83 }
84
85 p
86}
87
88pub fn continued_fraction(f: f64) -> Vec<u32> {
89 let mut list = vec![];
90 let mut r = f;
91
92 loop {
93 let t = r.trunc();
94 list.push(t as u32);
95
96 let diff = r - t;
97 if diff < 1e-3 {
98 break;
99 }
100
101 r = 1.0 / diff;
102 }
103
104 list
105}
106
107pub fn convergent(cf: &[u32]) -> (u32, u32) {
108 let len = cf.len();
109 if len == 1 {
110 return (cf[0], 1);
111 }
112
113 let mut s = 1;
114 let mut r = cf[len - 1];
115 for i in 2..len {
116 let tmp = s;
117 s = r;
118 r = cf[len - i] * r + tmp;
119 }
120 s += cf[0] * r;
121
122 (s, r)
123}
124
125pub fn parse_float(bin: &[char]) -> f64 {
126 let mut f = 0.0;
127
128 for (i, b) in bin.iter().enumerate() {
129 if *b == '0' {
130 continue;
131 }
132
133 f += 0.5_f64.powf((i + 1) as f64);
134 }
135
136 f
137}
138
139pub fn find_order(a: u32, n: u32, bin: &[char]) -> (u32, u32, bool) {
140 if bin.is_empty() {
141 return (0, 1, false);
142 }
143
144 let fv = parse_float(bin);
145 let cf = continued_fraction(fv);
146
147 for i in 0..cf.len() {
148 let (s, r) = convergent(&cf[0..(i + 1)]);
149
150 if modexp(a, r, n) == 1 {
151 return (s, r, true);
152 }
153 }
154
155 let (s, r) = convergent(&cf);
156 (s, r, false)
157}
158
159pub fn is_trivial(n: u32, factor: &[u32]) -> bool {
160 for p in factor.iter() {
161 if 1 < *p && *p < n && n.is_multiple_of(*p) {
162 return false;
163 }
164 }
165
166 true
167}
168
169pub fn is_odd(n: u32) -> bool {
170 n % 2 == 1
171}
172
173#[test]
174fn test_is_prime() {
175 assert!(is_prime(2));
176 assert!(is_prime(3));
177 assert!(is_prime(5));
178 assert!(is_prime(7));
179 assert!(is_prime(11));
180 assert!(is_prime(13));
181}
182
183#[test]
184fn test_coprime() {
185 assert!(gcd(15, coprime(15)) == 1);
186 assert!(gcd(21, coprime(21)) == 1);
187 assert!(gcd(35, coprime(35)) == 1);
188 assert!(gcd(51, coprime(51)) == 1);
189}
190
191#[test]
192fn test_gcd() {
193 assert_eq!(gcd(15, 2), 1);
194 assert_eq!(gcd(15, 3), 3);
195 assert_eq!(gcd(15, 4), 1);
196 assert_eq!(gcd(15, 5), 5);
197 assert_eq!(gcd(15, 6), 3);
198 assert_eq!(gcd(15, 7), 1);
199 assert_eq!(gcd(15, 8), 1);
200 assert_eq!(gcd(15, 9), 3);
201 assert_eq!(gcd(15, 10), 5);
202 assert_eq!(gcd(15, 11), 1);
203 assert_eq!(gcd(15, 12), 3);
204 assert_eq!(gcd(15, 13), 1);
205 assert_eq!(gcd(15, 14), 1);
206}
207
208#[test]
209fn test_base_exp() {
210 assert_eq!(base_exp(25), Some((5, 2)));
211 assert_eq!(base_exp(27), Some((3, 3)));
212 assert_eq!(base_exp(36), Some((6, 2)));
213 assert_eq!(base_exp(49), Some((7, 2)));
214}
215
216#[test]
217fn test_modexp2() {
218 assert_eq!(modexp2(7, 0, 15), 7);
219 assert_eq!(modexp2(7, 1, 15), 4);
220 assert_eq!(modexp2(7, 2, 15), 1);
221 assert_eq!(modexp2(7, 3, 15), 1);
222 assert_eq!(modexp2(7, 4, 15), 1);
223 assert_eq!(modexp2(7, 5, 15), 1);
224 assert_eq!(modexp2(7, 6, 15), 1);
225 assert_eq!(modexp2(7, 7, 15), 1);
226 assert_eq!(modexp2(7, 8, 15), 1);
227 assert_eq!(modexp2(7, 9, 15), 1);
228 assert_eq!(modexp2(7, 10, 15), 1);
229 assert_eq!(modexp2(7, 11, 15), 1);
230 assert_eq!(modexp2(7, 12, 15), 1);
231 assert_eq!(modexp2(7, 13, 15), 1);
232 assert_eq!(modexp2(7, 14, 15), 1);
233 assert_eq!(modexp2(0, 15, 15), 0);
234}
235
236#[test]
237fn test_continued_fraction() {
238 assert_eq!(continued_fraction(0.0), [0]);
239 assert_eq!(continued_fraction(1.0 / 16.0), [0, 16]);
240 assert_eq!(continued_fraction(4.0 / 16.0), [0, 4]);
241 assert_eq!(continued_fraction(7.0 / 16.0), [0, 2, 3, 1, 1]);
242 assert_eq!(continued_fraction(13.0 / 16.0), [0, 1, 4, 3]);
243 assert_eq!(continued_fraction(0.42857), [0, 2, 2, 1]);
244 assert_eq!(continued_fraction(0.166656494140625), [0, 6]);
245}
246
247#[test]
248fn test_convergent() {
249 assert_eq!(convergent(&continued_fraction(1.0 / 16.0)), (1, 16));
250 assert_eq!(convergent(&continued_fraction(4.0 / 16.0)), (1, 4));
251 assert_eq!(convergent(&continued_fraction(7.0 / 16.0)), (7, 16));
252 assert_eq!(convergent(&continued_fraction(13.0 / 16.0)), (13, 16));
253 assert_eq!(convergent(&continued_fraction(0.42857)), (3, 7));
254 assert_eq!(convergent(&continued_fraction(0.166656494140625)), (1, 6));
255}
256
257#[test]
258fn test_find_order() {
259 assert_eq!(find_order(7, 15, &['0', '1', '0']), (1, 4, true));
260 assert_eq!(find_order(7, 15, &['1', '0', '0']), (1, 2, false));
261 assert_eq!(find_order(7, 15, &['1', '1', '0']), (3, 4, true));
262}