Skip to main content

nautilus_model/defi/tick_map/
full_math.rs

1// -------------------------------------------------------------------------------------------------
2//  Copyright (C) 2015-2026 Nautech Systems Pty Ltd. All rights reserved.
3//  https://nautechsystems.io
4//
5//  Licensed under the GNU Lesser General Public License Version 3.0 (the "License");
6//  You may not use this file except in compliance with the License.
7//  You may obtain a copy of the License at https://www.gnu.org/licenses/lgpl-3.0.en.html
8//
9//  Unless required by applicable law or agreed to in writing, software
10//  distributed under the License is distributed on an "AS IS" BASIS,
11//  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12//  See the License for the specific language governing permissions and
13//  limitations under the License.
14// -------------------------------------------------------------------------------------------------
15
16use alloy_primitives::{I256, U160, U256};
17
18pub const Q128: U256 = U256::from_limbs([0, 0, 1, 0]);
19pub const Q96_U160: U160 = U160::from_limbs([0, 1 << 32, 0]);
20pub(crate) const DECIMAL_EXPONENT_MAX: u8 = 77;
21
22/// Contains 512-bit math functions for Uniswap V3 style calculations
23/// Handles "phantom overflow" - allows multiplication and division where
24/// intermediate values overflow 256 bits
25#[derive(Debug)]
26pub struct FullMath;
27
28impl FullMath {
29    /// Calculates floor(a×b÷denominator) with full precision
30    ///
31    /// Follows the Solidity implementation from Uniswap V3's `FullMath` library:
32    /// <https://github.com/Uniswap/v3-core/blob/main/contracts/libraries/FullMath.sol>
33    ///
34    /// # Errors
35    ///
36    /// Returns error if `denominator` is zero or the result would overflow 256 bits.
37    pub fn mul_div(a: U256, b: U256, mut denominator: U256) -> anyhow::Result<U256> {
38        // 512-bit multiply [prod1 prod0] = a * b
39        // Compute the product mod 2**256 and mod 2**256 - 1
40        // then use the Chinese Remainder Theorem to reconstruct
41        // the 512 bit result. The result is stored in two 256
42        // variables such that product = prod1 * 2**256 + prod0
43        let mm = a.mul_mod(b, U256::MAX);
44
45        // Least significant 256 bits of the product
46        let mut prod_0 = a * b;
47        let mut prod_1 = mm - prod_0 - U256::from_limbs([u64::from(mm < prod_0), 0, 0, 0]);
48
49        // Make sure the result is less than 2**256.
50        // Also prevents denominator == 0
51        if denominator <= prod_1 {
52            anyhow::bail!("Result would overflow 256 bits");
53        }
54
55        ///////////////////////////////////////////////
56        // 512 by 256 division.
57        ///////////////////////////////////////////////
58
59        // Make division exact by subtracting the remainder from [prod1 prod0]
60        // Compute remainder using mul_mod
61        let remainder = a.mul_mod(b, denominator);
62
63        // Subtract 256 bit number from 512 bit number
64        prod_1 -= U256::from_limbs([u64::from(remainder > prod_0), 0, 0, 0]);
65        prod_0 -= remainder;
66
67        // Factor powers of two out of denominator
68        // Compute largest power of two divisor of denominator.
69        // Always >= 1.
70        let mut twos = (-denominator) & denominator;
71
72        // Divide denominator by power of two
73        denominator /= twos;
74
75        // Divide [prod1 prod0] by the factors of two
76        prod_0 /= twos;
77
78        // Shift in bits from prod1 into prod0. For this we need
79        // to flip `twos` such that it is 2**256 / twos.
80        // If twos is zero, then it becomes one
81        twos = (-twos) / twos + U256::from(1);
82
83        prod_0 |= prod_1 * twos;
84
85        // Invert denominator mod 2**256
86        // Now that denominator is an odd number, it has an inverse
87        // modulo 2**256 such that denominator * inv = 1 mod 2**256.
88        // Compute the inverse by starting with a seed that is correct
89        // correct for four bits. That is, denominator * inv = 1 mod 2**4
90        let mut inv = (U256::from(3) * denominator) ^ U256::from(2);
91
92        // Now use Newton-Raphson iteration to improve the precision.
93        // Thanks to Hensel's lifting lemma, this also works in modular
94        // arithmetic, doubling the correct bits in each step.
95        inv *= U256::from(2) - denominator * inv; // inverse mod 2**8
96
97        inv *= U256::from(2) - denominator * inv; // inverse mod 2**16
98
99        inv *= U256::from(2) - denominator * inv; // inverse mod 2**32
100
101        inv *= U256::from(2) - denominator * inv; // inverse mod 2**64
102
103        inv *= U256::from(2) - denominator * inv; // inverse mod 2**128
104
105        inv *= U256::from(2) - denominator * inv; // inverse mod 2**256
106
107        // Because the division is now exact we can divide by multiplying
108        // with the modular inverse of denominator. This will give us the
109        // correct result modulo 2**256. Since the preconditions guarantee
110        // that the outcome is less than 2**256, this is the final result.
111        // We don't need to compute the high bits of the result and prod1
112        // is no longer required.
113        let result = prod_0 * inv;
114
115        Ok(result)
116    }
117
118    pub(crate) fn mul_div_scaled(
119        a: U256,
120        b: U256,
121        denominator: U256,
122        scales: &[U256],
123    ) -> anyhow::Result<U256> {
124        let mut quotient = Self::mul_div(a, b, denominator)?;
125        let mut remainder = a.mul_mod(b, denominator);
126
127        for &scale in scales {
128            let scaled_quotient = quotient
129                .checked_mul(scale)
130                .ok_or_else(|| anyhow::anyhow!("Scaled result exceeds 256-bit range"))?;
131            let scaled_remainder = Self::mul_div(remainder, scale, denominator)?;
132            quotient = scaled_quotient
133                .checked_add(scaled_remainder)
134                .ok_or_else(|| anyhow::anyhow!("Scaled result exceeds 256-bit range"))?;
135            remainder = remainder.mul_mod(scale, denominator);
136        }
137
138        Ok(quotient)
139    }
140
141    pub(crate) fn check_decimal_exponent(exponent: u8) -> anyhow::Result<()> {
142        anyhow::ensure!(
143            exponent <= DECIMAL_EXPONENT_MAX,
144            "Decimal exponent {exponent} exceeds supported maximum {DECIMAL_EXPONENT_MAX}"
145        );
146        Ok(())
147    }
148
149    pub(crate) fn pow10(exponent: u8) -> anyhow::Result<U256> {
150        Self::check_decimal_exponent(exponent)?;
151        U256::from(10)
152            .checked_pow(U256::from(exponent))
153            .ok_or_else(|| anyhow::anyhow!("Decimal exponent {exponent} exceeds U256 range"))
154    }
155
156    /// Calculates ceil(a×b÷denominator) with full precision
157    /// Returns `Ok` with the rounded result or an error when rounding cannot be performed safely.
158    ///
159    /// # Errors
160    ///
161    /// Returns error if `denominator` is zero or the rounded result would overflow `U256`.
162    pub fn mul_div_rounding_up(a: U256, b: U256, denominator: U256) -> anyhow::Result<U256> {
163        let result = Self::mul_div(a, b, denominator)?;
164
165        // Check if there's a remainder
166        if a.mul_mod(b, denominator).is_zero() {
167            Ok(result)
168        } else if result == U256::MAX {
169            anyhow::bail!("Result would overflow 256 bits")
170        } else {
171            Ok(result + U256::from(1))
172        }
173    }
174
175    /// Calculates ceil(a÷b) with proper rounding up
176    /// Equivalent to Solidity's divRoundingUp function
177    ///
178    /// # Errors
179    ///
180    /// Returns error if `b` is zero or if the rounded quotient would overflow `U256`.
181    pub fn div_rounding_up(a: U256, b: U256) -> anyhow::Result<U256> {
182        if b.is_zero() {
183            anyhow::bail!("Cannot divide by zero");
184        }
185
186        let quotient = a / b;
187        let remainder = a % b;
188
189        // Add 1 if there's a remainder (equivalent to gt(mod(x, y), 0) in assembly)
190        if remainder > U256::ZERO {
191            // Check for overflow before incrementing
192            if quotient == U256::MAX {
193                anyhow::bail!("Result would overflow 256 bits");
194            }
195            Ok(quotient + U256::from(1))
196        } else {
197            Ok(quotient)
198        }
199    }
200
201    /// Computes the integer square root of a 256-bit unsigned integer using the Babylonian method
202    #[must_use]
203    pub fn sqrt(x: U256) -> U256 {
204        if x.is_zero() {
205            return U256::ZERO;
206        }
207
208        if x == U256::from(1u128) {
209            return U256::from(1u128);
210        }
211
212        let mut z = x;
213        let mut y = (x + U256::from(1u128)) >> 1;
214
215        while y < z {
216            z = y;
217            y = (x / z + z) >> 1;
218        }
219
220        z
221    }
222
223    /// Truncates a U256 value to u128 by extracting the lower 128 bits.
224    ///
225    /// This matches Solidity's `uint128(value)` cast behavior, which discards
226    /// the upper 128 bits. If the value is larger than `u128::MAX`, the upper
227    /// bits are lost.
228    #[must_use]
229    pub fn truncate_to_u128(value: U256) -> u128 {
230        (value & U256::from(u128::MAX)).to::<u128>()
231    }
232
233    /// Converts an I256 signed integer to U256, mimicking Solidity's `uint256(int256)` cast.
234    ///
235    /// This performs a reinterpret cast, preserving the bit pattern:
236    /// - Positive values: returns the value as-is
237    /// - Negative values: returns the two's complement representation as unsigned
238    #[must_use]
239    pub fn truncate_to_u256(value: I256) -> U256 {
240        value.into_raw()
241    }
242
243    /// Converts a U256 unsigned integer to I256, mimicking Solidity's `int256(uint256)` cast.
244    ///
245    /// This performs a reinterpret cast, preserving the bit pattern.
246    /// Solidity's `SafeCast.toInt256()` checks the value fits in `I256::MAX`, then reinterprets.
247    ///
248    /// # Panics
249    /// Panics if the value exceeds `I256::MAX` (matching Solidity's require check)
250    #[must_use]
251    pub fn truncate_to_i256(value: U256) -> I256 {
252        I256::from_raw(value)
253    }
254}
255
256#[cfg(test)]
257mod tests {
258    use rstest::*;
259
260    use super::*;
261
262    #[rstest]
263    fn test_mul_div_reverts_denominator_zero() {
264        // Test that denominator 0 causes error
265        assert!(FullMath::mul_div(Q128, U256::from(5), U256::ZERO).is_err());
266
267        // Test with numerator overflow and denominator 0
268        assert!(FullMath::mul_div(Q128, Q128, U256::ZERO).is_err());
269    }
270
271    #[rstest]
272    fn test_mul_div_reverts_output_overflow() {
273        // Test output overflow: Q128 * Q128 / 1 would overflow
274        assert!(FullMath::mul_div(Q128, Q128, U256::from(1)).is_err());
275
276        // Test overflow with inputs that would cause prod1 >= denominator
277        // MAX * MAX / 1 would definitely overflow
278        assert!(FullMath::mul_div(U256::MAX, U256::MAX, U256::from(1)).is_err());
279
280        // Test with a smaller denominator that should still cause overflow
281        assert!(FullMath::mul_div(U256::MAX, U256::MAX, U256::from(2)).is_err());
282
283        // Test overflow with all max inputs and denominator = MAX - 1
284        assert!(FullMath::mul_div(U256::MAX, U256::MAX, U256::MAX - U256::from(1)).is_err());
285    }
286
287    #[rstest]
288    fn test_mul_div_all_max_inputs() {
289        // MAX * MAX / MAX = MAX
290        let result = FullMath::mul_div(U256::MAX, U256::MAX, U256::MAX).unwrap();
291        assert_eq!(result, U256::MAX);
292    }
293
294    #[rstest]
295    fn test_mul_div_accurate_without_phantom_overflow() {
296        // Calculate Q128 * 0.5 / 1.5 = Q128 / 3
297        let numerator_b = Q128 * U256::from(50) / U256::from(100); // 0.5
298        let denominator = Q128 * U256::from(150) / U256::from(100); // 1.5
299        let expected = Q128 / U256::from(3);
300
301        let result = FullMath::mul_div(Q128, numerator_b, denominator).unwrap();
302        assert_eq!(result, expected);
303    }
304
305    #[rstest]
306    fn test_mul_div_accurate_with_phantom_overflow() {
307        // Calculate Q128 * 35 * Q128 / (8 * Q128) = 35/8 * Q128 = 4.375 * Q128
308        let numerator_b = U256::from(35) * Q128;
309        let denominator = U256::from(8) * Q128;
310        let expected = U256::from(4375) * Q128 / U256::from(1000);
311
312        let result = FullMath::mul_div(Q128, numerator_b, denominator).unwrap();
313        assert_eq!(result, expected);
314    }
315
316    #[rstest]
317    fn test_mul_div_accurate_with_phantom_overflow_repeating_decimal() {
318        // Calculate Q128 * 1000 * Q128 / (3000 * Q128) = 1/3 * Q128
319        let numerator_b = U256::from(1000) * Q128;
320        let denominator = U256::from(3000) * Q128;
321        let expected = Q128 / U256::from(3);
322
323        let result = FullMath::mul_div(Q128, numerator_b, denominator).unwrap();
324        assert_eq!(result, expected);
325    }
326
327    #[rstest]
328    fn test_mul_div_basic_cases() {
329        // Simple case: 100 * 200 / 50 = 400
330        assert_eq!(
331            FullMath::mul_div(U256::from(100), U256::from(200), U256::from(50)).unwrap(),
332            U256::from(400)
333        );
334
335        // Test with 1: a * 1 / b = a / b
336        assert_eq!(
337            FullMath::mul_div(U256::from(1000), U256::from(1), U256::from(4)).unwrap(),
338            U256::from(250)
339        );
340
341        // Test division that results in 0 due to floor
342        assert_eq!(
343            FullMath::mul_div(U256::from(1), U256::from(1), U256::from(3)).unwrap(),
344            U256::ZERO
345        );
346    }
347
348    #[rstest]
349    fn test_mul_div_scaled_preserves_fractional_precision() {
350        let result = FullMath::mul_div_scaled(
351            U256::from(1),
352            U256::from(1),
353            U256::from(3),
354            &[U256::from(10), U256::from(10)],
355        )
356        .unwrap();
357
358        assert_eq!(result, U256::from(33));
359    }
360
361    #[rstest]
362    fn test_pow10_accepts_largest_supported_exponent() {
363        let result = FullMath::pow10(DECIMAL_EXPONENT_MAX).unwrap();
364        let expected = U256::from_str_radix(
365            "100000000000000000000000000000000000000000000000000000000000000000000000000000",
366            10,
367        )
368        .unwrap();
369
370        assert_eq!(result, expected);
371    }
372
373    #[rstest]
374    fn test_pow10_rejects_first_unsupported_exponent() {
375        let exponent = DECIMAL_EXPONENT_MAX + 1;
376        let error = FullMath::pow10(exponent).unwrap_err();
377
378        assert_eq!(
379            error.to_string(),
380            "Decimal exponent 78 exceeds supported maximum 77"
381        );
382    }
383
384    // mul_div_rounding_up tests
385    #[rstest]
386    fn test_mul_div_rounding_up_reverts_denominator_zero() {
387        // Test that denominator 0 causes error
388        assert!(FullMath::mul_div_rounding_up(Q128, U256::from(5), U256::ZERO).is_err());
389
390        // Test with numerator overflow and denominator 0
391        assert!(FullMath::mul_div_rounding_up(Q128, Q128, U256::ZERO).is_err());
392    }
393
394    #[rstest]
395    fn test_mul_div_rounding_up_reverts_output_overflow() {
396        // Test output overflow: Q128 * Q128 / 1 would overflow
397        assert!(FullMath::mul_div_rounding_up(Q128, Q128, U256::from(1)).is_err());
398
399        // Test overflow with all max inputs minus 1 - this should pass since MAX/MAX-1 = ~1
400        // but since there's a remainder, rounding up would still fit in U256
401        // Let's test a case that actually overflows after rounding
402        assert!(FullMath::mul_div_rounding_up(U256::MAX, U256::MAX, U256::from(2)).is_err());
403
404        // Test overflow with all max inputs and denominator = MAX - 1
405        assert!(
406            FullMath::mul_div_rounding_up(U256::MAX, U256::MAX, U256::MAX - U256::from(1)).is_err()
407        );
408    }
409
410    #[rstest]
411    fn test_mul_div_rounding_up_reverts_overflow_after_rounding_case_1() {
412        // Edge case discovered through fuzzing: mul_div succeeds but result is MAX with remainder
413        // so rounding up would overflow
414        let a = U256::from_str_radix("535006138814359", 10).unwrap();
415        let b = U256::from_str_radix(
416            "432862656469423142931042426214547535783388063929571229938474969",
417            10,
418        )
419        .unwrap();
420        let denominator = U256::from(2);
421
422        assert!(FullMath::mul_div_rounding_up(a, b, denominator).is_err());
423    }
424
425    #[rstest]
426    fn test_mul_div_rounding_up_reverts_overflow_after_rounding_case_2() {
427        // Another edge case discovered through fuzzing: tests boundary condition where
428        // mul_div returns MAX-1 but with remainder, so rounding up would cause overflow
429        let a = U256::from_str_radix(
430            "115792089237316195423570985008687907853269984659341747863450311749907997002549",
431            10,
432        )
433        .unwrap();
434        let b = U256::from_str_radix(
435            "115792089237316195423570985008687907853269984659341747863450311749907997002550",
436            10,
437        )
438        .unwrap();
439        let denominator = U256::from_str_radix(
440            "115792089237316195423570985008687907853269984653042931687443039491902864365164",
441            10,
442        )
443        .unwrap();
444
445        assert!(FullMath::mul_div_rounding_up(a, b, denominator).is_err());
446    }
447
448    #[rstest]
449    fn test_mul_div_rounding_up_all_max_inputs() {
450        // MAX * MAX / MAX = MAX (no rounding needed)
451        let result = FullMath::mul_div_rounding_up(U256::MAX, U256::MAX, U256::MAX).unwrap();
452        assert_eq!(result, U256::MAX);
453    }
454
455    #[rstest]
456    fn test_mul_div_rounding_up_accurate_without_phantom_overflow() {
457        // Calculate Q128 * 0.5 / 1.5 = Q128 / 3, but with rounding up
458        let numerator_b = Q128 * U256::from(50) / U256::from(100); // 0.5
459        let denominator = Q128 * U256::from(150) / U256::from(100); // 1.5
460        let expected = Q128 / U256::from(3) + U256::from(1); // Rounded up
461
462        let result = FullMath::mul_div_rounding_up(Q128, numerator_b, denominator).unwrap();
463        assert_eq!(result, expected);
464    }
465
466    #[rstest]
467    fn test_mul_div_rounding_up_accurate_with_phantom_overflow() {
468        // Calculate Q128 * 35 * Q128 / (8 * Q128) = 35/8 * Q128 = 4.375 * Q128
469        // This should be exact (no remainder), so no rounding up needed
470        let numerator_b = U256::from(35) * Q128;
471        let denominator = U256::from(8) * Q128;
472        let expected = U256::from(4375) * Q128 / U256::from(1000);
473
474        let result = FullMath::mul_div_rounding_up(Q128, numerator_b, denominator).unwrap();
475        assert_eq!(result, expected);
476    }
477
478    #[rstest]
479    fn test_mul_div_rounding_up_accurate_with_phantom_overflow_repeating_decimal() {
480        // Calculate Q128 * 1000 * Q128 / (3000 * Q128) = 1/3 * Q128, with rounding up
481        let numerator_b = U256::from(1000) * Q128;
482        let denominator = U256::from(3000) * Q128;
483        let expected = Q128 / U256::from(3) + U256::from(1); // Rounded up due to remainder
484
485        let result = FullMath::mul_div_rounding_up(Q128, numerator_b, denominator).unwrap();
486        assert_eq!(result, expected);
487    }
488
489    #[rstest]
490    fn test_mul_div_rounding_up_basic_cases() {
491        // Test exact division (no rounding needed)
492        assert_eq!(
493            FullMath::mul_div_rounding_up(U256::from(100), U256::from(200), U256::from(50))
494                .unwrap(),
495            U256::from(400)
496        );
497
498        // Test division with remainder (rounding up needed)
499        assert_eq!(
500            FullMath::mul_div_rounding_up(U256::from(1), U256::from(1), U256::from(3)).unwrap(),
501            U256::from(1) // 0 rounded up to 1
502        );
503
504        // Test another rounding case: 7 * 3 / 4 = 21 / 4 = 5.25 -> 6
505        assert_eq!(
506            FullMath::mul_div_rounding_up(U256::from(7), U256::from(3), U256::from(4)).unwrap(),
507            U256::from(6)
508        );
509
510        // Test case with zero result and zero remainder
511        assert_eq!(
512            FullMath::mul_div_rounding_up(U256::ZERO, U256::from(100), U256::from(3)).unwrap(),
513            U256::ZERO
514        );
515    }
516
517    #[rstest]
518    fn test_mul_div_rounding_up_overflow_at_max() {
519        // Test that rounding up when result is already MAX causes overflow
520        // We need a case where mul_div returns MAX but there's a remainder
521        // This is tricky to construct, so we test the boundary condition
522        assert!(FullMath::mul_div_rounding_up(U256::MAX, U256::from(2), U256::from(2)).is_ok());
523
524        // This should succeed: MAX * 1 / 1 = MAX (no remainder)
525        assert_eq!(
526            FullMath::mul_div_rounding_up(U256::MAX, U256::from(1), U256::from(1)).unwrap(),
527            U256::MAX
528        );
529    }
530
531    #[rstest]
532    fn test_truncate_to_u128_preserves_small_values() {
533        // Small value (fits in u128) should be preserved exactly
534        let value = U256::from(12345u128);
535        assert_eq!(FullMath::truncate_to_u128(value), 12345u128);
536
537        // u128::MAX should be preserved
538        let max_value = U256::from(u128::MAX);
539        assert_eq!(FullMath::truncate_to_u128(max_value), u128::MAX);
540    }
541
542    #[rstest]
543    fn test_truncate_to_u128_discards_upper_bits() {
544        // Value = u128::MAX + 1 (sets bit 128)
545        // Lower 128 bits = 0, so result should be 0
546        let value = U256::from(u128::MAX) + U256::from(1);
547        assert_eq!(FullMath::truncate_to_u128(value), 0);
548
549        // Value with both high and low bits set:
550        // High 128 bits = 0xFFFF...FFFF, Low 128 bits = 0x1234
551        let value = (U256::from(u128::MAX) << 128) | U256::from(0x1234u128);
552        assert_eq!(FullMath::truncate_to_u128(value), 0x1234u128);
553    }
554}