1use ndarray::arr2;
8use ndarray::Array2;
9use num_bigint::BigUint;
10use num::FromPrimitive;
11
12pub fn fib_with_mod(n: u64, modulo: u64) -> u64 {
26 if n == 0 {
27 return 0;
28 }
29 if n == 1 {
30 return 1;
31 }
32
33 let f = vec![0, 1];
34 let t = arr2(&[
35 [0, 1],
36 [1, 1]
37 ]);
38 let power_t = matrix_power_with_mod(&t, n, modulo);
39 let mut answer = 0;
40 for i in 0..2 {
41 answer = (answer + (power_t[[0, i]] * f[i])) % modulo;
42 }
43 answer
44}
45
46
47pub fn bigfib_with_mod(n: &BigUint, modulo: &BigUint) -> BigUint {
83 let ZERO: BigUint = FromPrimitive::from_u64(0).unwrap();
84 let ONE: BigUint = FromPrimitive::from_u64(1).unwrap();
85 if n == &ZERO || n == &ONE {
86 return n.clone();
87 }
88
89 let f: Vec<BigUint> = vec![ZERO.clone(), ONE.clone()];
90 let t: Array2<BigUint> = arr2(&[
91 [ZERO.clone(), ONE.clone()],
92 [ONE.clone(), ONE.clone()]
93 ]);
94 let power_t = bigfib_matrix_power(&t, n, modulo);
95 let mut answer: BigUint = ZERO.clone();
96 for i in 0..2 {
97 answer = (answer + (&power_t[[0, i]] * &f[i])) % modulo;
98 }
99 return answer;
100}
101
102
103fn bigfib_matrix_power(mat: &Array2<BigUint>, pow: &BigUint, modulo: &BigUint) -> Array2<BigUint> {
104 let ONE: BigUint = FromPrimitive::from_u64(1).unwrap();
105 let TWO: BigUint = FromPrimitive::from_u64(2).unwrap();
106 if pow == &ONE {
107 return mat.clone();
108 }
109 if pow % &TWO == ONE {
110 return bigfib_multiply(
111 &mat,
112 &bigfib_matrix_power(mat, &(pow - ONE), modulo),
113 modulo
114 );
115 }
116 let x = bigfib_matrix_power(mat, &(pow / TWO), modulo);
117 bigfib_multiply(&x, &x, modulo)
118}
119
120
121fn bigfib_multiply(a: &Array2<BigUint>, b: &Array2<BigUint>, modulo: &BigUint) -> Array2<BigUint> {
122 let ZERO: BigUint = FromPrimitive::from_u64(0).unwrap();
123 let mut return_mat: Array2<BigUint> = arr2(&[
124 [ZERO.clone(), ZERO.clone()],
125 [ZERO.clone(), ZERO.clone()]
126 ]);
127
128 for i in 0..2 {
129 for j in 0..2 {
130 for k in 0..2 {
131 let big_val: BigUint = &return_mat[[i, j]] + (&a[[i, k]] * &b[[k, j]]);
132 return_mat[[i, j]] = big_val % modulo;
133 }
134 }
135 }
136 return_mat
137}
138
139
140fn multiply_with_mod(a: &Array2<u64>, b: &Array2<u64>, modulo: u64) -> Array2<u64> {
141 let mut return_mat: Array2<u64> = Array2::zeros((2, 2));
142
143 let big_mod: BigUint = FromPrimitive::from_u64(modulo).unwrap();
144 for i in 0..2 {
145 for j in 0..2 {
146 for k in 0..2 {
147 let mat_ij: BigUint = FromPrimitive::from_u64(return_mat[[i, j]]).unwrap();
148 let a_ik: BigUint = FromPrimitive::from_u64(a[[i, k]]).unwrap();
149 let b_kj: BigUint = FromPrimitive::from_u64(b[[k, j]]).unwrap();
150
151 let big_val: BigUint = (mat_ij + (
152 a_ik * b_kj
153 )) % &big_mod;
154
155 return_mat[[i, j]] = small_big_int_to_u64(&big_val);
156 }
157 }
158 }
159 return_mat
160}
161
162
163fn matrix_power_with_mod(mat: &Array2<u64>, pow: u64, modulo: u64) -> Array2<u64> {
164 if pow == 1 {
165 return mat.clone();
166 }
167 if pow % 2 == 1 {
168 return multiply_with_mod(
169 &mat,
170 &matrix_power_with_mod(
171 mat,
172 pow - 1,
173 modulo
174 ),
175 modulo
176 );
177 }
178 let x = matrix_power_with_mod(mat, pow / 2, modulo);
179 multiply_with_mod(&x, &x, modulo)
180}
181
182
183fn small_big_int_to_u64(big_int: &BigUint) -> u64 {
184 let mut result: u64 = 0;
185
186 let digits = big_int.to_radix_be(10);
187 for i in 0..digits.len() - 1 {
188 result = result + digits[i] as u64;
189 result = result * 10;
190 }
191 result + digits[digits.len() - 1] as u64
192}
193
194#[cfg(test)]
195mod tests {
196 use crate::*;
197
198 #[test]
199 fn test_first_few() {
200 assert_eq!(fib_with_mod(0, 10), 0);
201 assert_eq!(fib_with_mod(1, 10), 1);
202 assert_eq!(fib_with_mod(2, 10), 1);
203 assert_eq!(fib_with_mod(3, 10), 2);
204 assert_eq!(fib_with_mod(4, 10), 3);
205 assert_eq!(fib_with_mod(5, 10), 5);
206 }
207
208 #[test]
209 fn test_modulo() {
210 assert_eq!(fib_with_mod(100, 1_000_000_000), 261_915_075);
211 assert_eq!(fib_with_mod(100, 1_000_000), 915_075);
212 assert_eq!(fib_with_mod(100, 1_000), 75);
213 assert_eq!(fib_with_mod(100, 10), 5);
214 }
215
216 #[test]
217 fn test_big() {
218 assert_eq!(fib_with_mod(1_000_000_000_000_000, 1_000_000), 546_875);
219 assert_eq!(fib_with_mod(1_955_995_342_096_516, u64::MAX), 2_886_946_313_980_141_317);
220 }
221
222 #[test]
223 fn test_bigfib() {
224 let ns = vec![0, 1, 2, 3, 10, 1_000_000_000_000_000, 1_955_995_342_096_516];
225 let modulos = vec![10, 10, 20, 30, 100, 1_000_000, u64::MAX];
226 let expected_results = vec![0, 1, 1, 2, 55, 546_875, 2_886_946_313_980_141_317];
227
228 for i in 0..ns.len() {
229 assert_eq!(
230 bigfib_with_mod(
231 &FromPrimitive::from_u64(ns[i]).unwrap(),
232 &FromPrimitive::from_u64(modulos[i]).unwrap()
233 ),
234 FromPrimitive::from_u64(expected_results[i]).unwrap()
235 );
236 }
237 }
238
239 #[test]
240 fn test_large_bigfib() {
241 let n: BigUint = BigUint::from_slice(&[100u32, 100, 100, 100, 15129, 12319]);
242 let modulo: BigUint = BigUint::from_slice(&[14u32, 12, 1923876123, 12]);
243 let expected_result: BigUint = BigUint::from_slice(&[2743227343u32, 920986447, 1158660944, 5]);
244 assert_eq!(
245 bigfib_with_mod(&n, &modulo),
246 expected_result
247 );
248 }
249}