diff --git a/docs/algorithm/add.md b/docs/algorithm/add.md new file mode 100644 index 0000000..4089a2a --- /dev/null +++ b/docs/algorithm/add.md @@ -0,0 +1,11 @@ +# Addition +The addition algorithm operates in Base $2^{64}$, using 64-bit unsigned integers (`uint64_t`) to represent each limb of the arbitrary-precision integer. Overflows beyond the 64-bit boundary are detected using standard unsigned arithmetic comparisons. + +1. **Iteration:** The algorithm traverses from the least significant limb (LSB, `index 0`) to the most significant limb (MSB, `index N`). +2. **Carry Propagation:** For each limb, represented by `i`, it calculates: + `sum = aVal + bVal` + An overflow is detected if the intermediate `sum` is smaller than `aVal`. +3. **Carry Handling:** + * The carry for the next iteration (`nextCarry`) is initially flagged as `1` if an overflow occurred in the initial sum. + * The carry from the previous iteration is then added to the `sum`. If this addition causes a secondary overflow (detected by checking `sum < carry`), `nextCarry` is incremented. +4. **Finalization:** If a carry remains after the final limb is processed, an additional limb containing the value `1` is appended to the vector. Finally, the vector is normalized to prune trailing zeros. \ No newline at end of file diff --git a/docs/algorithm/barret-reducer.md b/docs/algorithm/barret-reducer.md new file mode 100644 index 0000000..1d883a3 --- /dev/null +++ b/docs/algorithm/barret-reducer.md @@ -0,0 +1,62 @@ +# Barrett Reduction + +This module implements **Barrett reduction**, a mathematical technique used to speed up modulo operations on large integers. + +--- + +## 1. What is Barrett Reduction? + +In computer science, when working with very large integers (arbitrary-precision arithmetic), traditional division is one of the slowest operations. Calculating a simple remainder $x \pmod m$ normally requires a slow long-division step. + +Barrett reduction solves this by replacing the slow division with a combination of fast **multiplication** and **digit shifting** (limb-shifting). + +### The Core Concept: +To find the remainder of $x$ divided by $m$, we want to find the quotient $q = \lfloor x / m \rfloor$. Instead of dividing $x$ by $m$, we can multiply $x$ by the "inverse" of $m$ (i.e., $1/m$). + +Since we work with integers, we scale this inverse up to make it an integer. Let $b = 2^{64}$ represent our base (since each "limb" or array element in our code is a 64-bit integer), and let $k$ be the number of limbs in the modulus $m$. + +1. **Precomputation:** + We precalculate a scaling constant, $\mu$ (mu): + $$\mu = \lfloor \frac{b^{2k}}{m} \rfloor$$ + This is the only actual division we need to perform, and it is done once when initializing the reducer. + +2. **Estimating the Quotient:** + To calculate the modulo, we estimate the quotient $q$. Instead of full multiplication, we use a slightly adjusted formula to keep the numbers from growing too large: + $$q \approx \lfloor \frac{\lfloor x / b^{k-1} \rfloor \cdot \mu}{b^{k+1}} \rfloor$$ + In code, dividing by $b^{k-1}$ and $b^{k+1}$ is done using simple right-shifts of the array elements (limbs). + +3. **Calculating the Remainder and Correcting:** + Once we have our estimated quotient $q$, we calculate the remainder: + $$r = x - (q \cdot m)$$ + Because we used an approximation, our estimated quotient might be slightly smaller than the actual quotient (by at most 2). Consequently, the temporary remainder $r$ might still be larger than our modulus $m$. We correct this by subtracting $m$ from $r$ until $r < m$ (which is guaranteed to happen in at most two steps). + +--- + +## 2. Estimation and Final Correction + +To understand how the mathematical approximation interacts with the final adjustment loop, we can look at how the estimate behaves. + +### Visualizing the Process: +The following diagram illustrates how the algorithm behaves: + +![Quotient Estimation and Remainder Progress](img/barrett_reduction.png) + +* **Top Graph (Quotient Estimation):** Shows that the estimated quotient (due to rounding down during limb-shifts) occasionally falls slightly below the actual quotient. +* **Bottom Graph (Remainder Progress):** Shows the preliminary remainder $r$. Whenever the estimated quotient is slightly too small, the temporary remainder $r$ crosses the red line (the modulus $m$). In these instances, the correction loop subtracts $m$ until the remainder is brought within the correct bounds. + +--- + +## 3. Algorithm Flowchart + +The step-by-step control flow is represented in the diagram below: + +![Barrett Reduction Flowchart](img/barrett_flowchart.png) + +### Summary of Steps in Code: +1. **Shortcut check:** If the input $x$ is already smaller than the modulus $m$, return $x$ immediately. +2. **First shift:** Shift $x$ right by $k-1$ limbs to get $q_1$. +3. **Multiply:** Multiply $q_1$ by our precomputed value $\mu$ to get $q_2$. +4. **Second shift:** Shift $q_2$ right by $k+1$ limbs to get our final estimated quotient, $q_3$. +5. **Find remainder:** Calculate $r = x - (q_3 \cdot m)$. +6. **Adjust:** While $r \geq m$, subtract $m$ from $r$ (runs at most twice). +7. **Return:** Return the final reduced value $r$. \ No newline at end of file diff --git a/docs/algorithm/div.md b/docs/algorithm/div.md new file mode 100644 index 0000000..aa66016 --- /dev/null +++ b/docs/algorithm/div.md @@ -0,0 +1,166 @@ +# Division +This documentation explains how the **Bitwise Long Division** (Binary Long Division) is implemented. + +## Definitions +* **Dividend (`data`):** The number being divided. +* **Divisor (`divisor`):** The number you are dividing by. +* **Quotient (`quotient`):** This is the result number. How many times the divisor fits entirely into the dividend. +* **Remainder / Modulo (`remaining`):** The leftover amount that is strictly less than the divisor. +* **Dividend Mask (`dividendMask`):** The active "working remainder" used dynamically within the long division loop to evaluate the current subset of bits. + +## The Basic Concept +Binary Long Division works exactly like the long division taught in elementary school, but it is simpler because the quotient at any step can only ever be `1` or `0` (it either fits, or it doesn't). +1. Isolate the highest piece of the dividend. +2. Check if the divisor fits into it. +3. If it fits, write `1` to the quotient, subtract the divisor, and "bring down" the next bit. +4. If it doesn't fit, write `0` to the quotient, don't subtract, and "bring down" the next bit. + +**Quick Binary Example `7 / 2`:** +* `7` (`111`) +* `2` (`10`) + +**The Setup:** +```text + [Quotient] + _____________ +Divisor: 10 | 1 1 1 <-- Dividend + ^ ^ + MSB LSB + (Highest) (Lowest) +``` + +**Step 1: Isolate the highest bit (MSB)** +```text + 0 <-- 10 does not fit into 1. Write 0 to quotient. + _______ + 10 | 1 1 1 + - 0 <-- No subtraction occurs. + --- + 1 <-- Working remainder is 1. +``` + +**Step 2: Bring down the next bit** +```text + 0 1 <-- 10 fits into 11! Write 1 to quotient. + _______ + 10 | 1 1 1 + . | + 1 1 <-- Bring down the 1. Working remainder is now 11. + - 1 0 <-- Subtract divisor (10). + ----- + 1 <-- New working remainder is 1. +``` + +**Step 3: Bring down the last bit (LSB)** +```text + 0 1 1 <-- 10 fits into 11! Write 1 to quotient. + _______ + 10 | 1 1 1 + . | + 1 1 <-- Bring down the final 1. Working remainder is 11. + - 1 0 <-- Subtract divisor (10). + ----- + 1 <-- Final Remainder! +``` + +**Result:** +The quotient is `001` (`3` in decimal). The leftover value at the bottom is our remainder: **`1`**. +*(7 / 2 = 3 R 1)*. + + +## Implementation +The division algorithm is implemented within 5 steps. + +### 1. Initialization +1. **Division by Zero:** If the divisor is `0` or the dividend is `0` (yielding an invalid MSB index), we abort and default the quotient and remainder to `0`. +2. **Finding the MSB:** Using the `getStartBitIndex` function, we get the absolute highest active bit index across the 64-bit limb array and assign it to `initialDividendIndex`. This ensures that we skip processing leading zero-limbs or "ghost" zeros. +3. **Preallocation:** Because the exact bit-length of the quotient correlates directly to the `initialDividendIndex`, the `quotient` vector is pre-allocated to the maximum required size: `(initialDividendIndex / 64) + 1` 64-bit limbs. + +### 2. Working Remainder (`dividendMask`) +`dividendMask` is initialized by taking a value of `0` and using the `addBitFromNumber` helper to copy down the absolute highest bit from the **dividend** at `initialDividendIndex`. + +### 3. Bitwise Evaluation Loop +In a loop, bit by bit is processed from `initialDividendIndex` down to `-1`. +For every bit position, the mathematical power index (`currentQBitIndex`) is evaluated. + +At each step, we check if the `dividendMask` (working remainder) is $\ge$ the `divisor`. + +**Case A: The Mask is $\ge$ the Divisor** +* The divisor "fits" inside the working remainder. +* A `1` is written to the exact corresponding bit inside the pre-allocated `quotient` vector using bitwise operations: `quotient[currentQBitIndex / 64] |= (1ULL << (currentQBitIndex % 64))`. +* If we have reached the end of the dividend (`dividendIndex < 0`), the final modulo subtraction occurs, and the loop terminates. +* Otherwise, the `divisor` is subtracted from the `dividendMask` in-place, the next bit from the dividend is appended to the mask using `addBitFromNumberInPlace`, and `dividendIndex` is decremented. + +**Case B: The Mask is $<$ the Divisor** +* The divisor does not fit. The quotient bit remains `0` (its default pre-allocated state). +* No subtraction occurs. +* If we have reached the end of the dividend (`dividendIndex < 0`), the `dividendMask` is preserved as-is as the remainder, and the loop terminates. +* Otherwise, the next bit from the dividend is appended to the mask using `addBitFromNumberInPlace` to increase its value for the next loop iteration, and `dividendIndex` is decremented. + +### 4. Modulo (The Remainder) +Because this is integer division, there is often a fractional remainder. The `div` function accepts an optional pointer to a `remaining` vector (`ByteArray* remaining`). +When the loop finishes processing the final bit (`dividendIndex < 0`), whatever mathematical value is left inside the `dividendMask` represents the modulo. If the pointer is provided, the mask (or subtracted mask) is copied into it. + +*(Note: Because of this architecture, evaluating `A % B` requires the same computational effort as `A / B`. Therefore, if both the quotient and remainder are needed, they are extracted simultaneously).* + +### 5. Final Normalization +Even though pre-allocation is tightly bound to the `initialDividendIndex`, the final quotient might have leading zeros depending on the magnitude of the divisor. The `div` function concludes by stripping any trailing zero-limbs from the little-endian vector using `normalizeVector` to maintain strict Base $2^{64}$ normalization guarantees. + +--- + +## Visualization: `13 / 3` +* **Dividend:** `13` (`1101`) +* **Divisor:** `3` (`0011`) + +**Initialization:** +* `getStartBitIndex(1101)` returns `3` (the 0-based index of the highest `1` bit). +* `dividendMask` (the working remainder) is initialized to `0000`. +* `quotient` is pre-allocated and initialized to `0000`. + +### Step 1: Processing Bit Index 3 +* **Action:** Shift `dividendMask` left by 1, and bring down Bit 3 of the dividend (`1`). +* **Mask Status:** `dividendMask` becomes `0001` (Decimal: 1). +* **Comparison:** Is `0001` $\ge$ `0011` (Divisor)? $\rightarrow$ **FALSE** +* **Result:** + * Divisor does not fit. No subtraction. + * Bit 3 of `quotient` remains `0`. + * **Current Quotient:** `0000` + * **Current Mask:** `0001` + +### Step 2: Processing Bit Index 2 +* **Action:** Shift `dividendMask` left by 1 (`0001` $\rightarrow$ `0010`), and bring down Bit 2 of the dividend (`1`). +* **Mask Status:** `dividendMask` becomes `0011` (Decimal: 3). +* **Comparison:** Is `0011` $\ge$ `0011` (Divisor)? $\rightarrow$ **TRUE** +* **Result:** + * Divisor fits! + * Set Bit 2 of `quotient` to `1` using bitwise OR. + * Subtract divisor from mask: `0011` - `0011` = `0000`. + * **Current Quotient:** `0100` + * **Current Mask:** `0000` + +### Step 3: Processing Bit Index 1 +* **Action:** Shift `dividendMask` left by 1 (`0000` $\rightarrow$ `0000`), and bring down Bit 1 of the dividend (`0`). +* **Mask Status:** `dividendMask` becomes `0000` (Decimal: 0). +* **Comparison:** Is `0000` $\ge$ `0011` (Divisor)? $\rightarrow$ **FALSE** +* **Result:** + * Divisor does not fit. No subtraction. + * Bit 1 of `quotient` remains `0`. + * **Current Quotient:** `0100` + * **Current Mask:** `0000` + +### Step 4: Processing Bit Index 0 (LSB) +* **Action:** Shift `dividendMask` left by 1 (`0000` $\rightarrow$ `0000`), and bring down Bit 0 of the dividend (`1`). +* **Mask Status:** `dividendMask` becomes `0001` (Decimal: 1). +* **Comparison:** Is `0001` $\ge$ `0011` (Divisor)? $\rightarrow$ **FALSE** +* **Result:** + * Divisor does not fit. No subtraction. + * Bit 0 of `quotient` remains `0`. + * **Current Quotient:** `0100` + * **Current Mask:** `0001` + +### Final Output Evaluation: +The loop terminates because we have processed bit index `0`. +1. **The Quotient:** The pre-allocated quotient vector holds `0100`, which mathematically evaluates to **`4`**. +2. **The Remainder:** The `dividendMask` is left holding `0001`, which mathematically evaluates to **`1`**. If a pointer for the modulo was provided, this value is safely copied over. + +**Conclusion:** `1101 / 0011` = `0100` with a remainder of `0001` ($13 / 3 = 4\text{ R }1$). \ No newline at end of file diff --git a/docs/algorithm/img/barrett_flowchart.png b/docs/algorithm/img/barrett_flowchart.png new file mode 100644 index 0000000..e2b58c5 Binary files /dev/null and b/docs/algorithm/img/barrett_flowchart.png differ diff --git a/docs/algorithm/img/barrett_reduction.png b/docs/algorithm/img/barrett_reduction.png new file mode 100644 index 0000000..5f4c2e3 Binary files /dev/null and b/docs/algorithm/img/barrett_reduction.png differ diff --git a/docs/algorithm/img/cpp_registermodell_sequentiell.png b/docs/algorithm/img/cpp_registermodell_sequentiell.png new file mode 100644 index 0000000..e5c3bbb Binary files /dev/null and b/docs/algorithm/img/cpp_registermodell_sequentiell.png differ diff --git a/docs/algorithm/img/cpp_ueberfuehrung_bytearray_sequentiell.png b/docs/algorithm/img/cpp_ueberfuehrung_bytearray_sequentiell.png new file mode 100644 index 0000000..b140ea5 Binary files /dev/null and b/docs/algorithm/img/cpp_ueberfuehrung_bytearray_sequentiell.png differ diff --git a/docs/algorithm/img/dezimales_modell.png b/docs/algorithm/img/dezimales_modell.png new file mode 100644 index 0000000..66bc001 Binary files /dev/null and b/docs/algorithm/img/dezimales_modell.png differ diff --git a/docs/algorithm/img/teilprodukte_alignment_128bit.png b/docs/algorithm/img/teilprodukte_alignment_128bit.png new file mode 100644 index 0000000..f16a0e0 Binary files /dev/null and b/docs/algorithm/img/teilprodukte_alignment_128bit.png differ diff --git a/docs/algorithm/mul.md b/docs/algorithm/mul.md new file mode 100644 index 0000000..0a0c0c5 --- /dev/null +++ b/docs/algorithm/mul.md @@ -0,0 +1,31 @@ +# Multiplication Algorithm + +This document outlines the sequential multiplication logic implemented in the `BigInt` class. + +--- + +### 1. Decimal Analogy +The multiplication process starts with a standard positional approach. Using a decimal analogy (e.g., $23 \times 45$), the operands are broken down into their positional parts (tens and ones), cross-multiplied, scaled by their place value, and summed to produce the final product. + +![Dezimales Multiplikationsmodell](img/dezimales_modell.png) + +--- + +### 2. 128-Bit Alignment Layout +In the binary domain, multiplying two 64-bit inputs (each split into 32-bit halves) results in four 64-bit partial products ($r1$, $r2$, $r3$, $r4$). This layout diagram illustrates how these partial products are shifted and placed within the target 128-bit result space. + +![128-Bit-Alignment der Teilprodukte](img/teilprodukte_alignment_128bit.png) + +--- + +### 3. C++ Register Model +The sequential C++ register model processes these partial products to form a 128-bit product. It splits the calculation into a lower 64-bit register (`fl`) and a higher 64-bit register (`fh`), sequentially accumulating the terms and tracking any intermediate overflows using a running carry. + +![C++ Registermodell](img/cpp_registermodell_sequentiell.png) + +--- + +### 4. Result Vector Integration +Finally, the calculated 128-bit product (`fl` and `fh`) is integrated back into the target array elements (`result`). This step handles addition with existing cell values, tracks two potential overflows, and generates a carry-over value to propagate into the subsequent step. + +![Sequentielle Überführung ins ByteArray](img/cpp_ueberfuehrung_bytearray_sequentiell.png) \ No newline at end of file diff --git a/docs/algorithm/sub.md b/docs/algorithm/sub.md new file mode 100644 index 0000000..6fcb257 --- /dev/null +++ b/docs/algorithm/sub.md @@ -0,0 +1,7 @@ +# Subtraction +Subtraction implements a "borrow" propagation system to ensure mathematical correctness within an unsigned 64-bit integer space: + +1. **Underflow Protection:** Before processing, `isBigger` is checked. If the subtrahend `b` is larger than the minuend `a`, the operation clamps and returns `0` (or assigns `{0}` in-place) because negative integers are not supported. +2. **Borrow Propagation:** During limb-by-limb iteration, the algorithm calculates `diff = aVal - bVal`. An underflow is flagged (`nextBorrow = true`) if `aVal < bVal`. +3. **Regrouping:** If an underflow occurs, the algorithm effectively borrows $2^{64}$ (the base value) from the next higher limb (`index i + 1`). When subtracting the previous borrow from the current difference, a secondary underflow check (`diff < borrowVal`) ensures the borrow state is correctly propagated to the next iteration. +4. **Normalization:** Post-calculation, `normalizeVector` is called to prune trailing zero-limbs, keeping the internal vector size proportional to the actual magnitude of the result. \ No newline at end of file