Fraction vs PyFraction — gcd Refactor & Head-to-Head

---

# 🧮 Fraction vs PyFraction

## Exact rationals in C++ — the gcd refactor & a fair head-to-head

### `fractions.hpp` 🧰 · `pyfractions.hpp` 🐍

@luk036 👨‍💻 · 2026 📅

---

### 📋 Agenda

**Part 1 — The Contenders** 🧩
- Two classes, one name
- Semantics matrix
- Why exactness matters

**Part 2 — How We Measure** 📏
- `-O3` / `/O2` flags
- The gcd tax

**Part 3 — Killing Template Recursion** 💥
- Why the optimizer gave up
- Concrete overloads

]

**Part 4 — Six Optimizations** ⚡
- `from_raw`, zero-gcd ints
- 128-bit compare, `_umul128`

**Part 5 — Head-to-Head** 📊
- g++ · clang-cl · MSVC
- Who wins what

**Part 6 — Size & Lessons** 🎓
- Lines of code
- Verification pipeline

]

---

## 🧩 Part 1 — The Contenders

---

### 🗂️ Two Classes, One Name

<pre>
graph TB
    NS["namespace fractions"]
    S["Fraction[T = int64_t]\nfractions.hpp\nfull API · zero-denominator OK"]
    P["Fraction[T = int64_t]\npyfractions.hpp\nPython-compatible"]
    E["ExtFraction[T]\nextfractions.hpp\nsibling surface"]
    NS --> S
    NS --> P
    NS --> E
    S -.->|"same name!"| P
    style NS fill:#eceff4,stroke:#4c566a,stroke-width:3px
    style S fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style P fill:#d1c4e9,stroke:#6a1b9a,stroke-width:3px
    style E fill:#c8e6c9,stroke:#2e7d32,stroke-width:3px
</pre>

]

> ⚠️ Both now default to `T = int64_t`; both define `fractions::Fraction` — they cannot share a translation unit.

---

### ⚖️ Same Name ≠ Same Contract

> 🎯 Choose by **semantics** first, speed second.

---

### 📐 Why Exact Fractions?

Geometric predicates return a **sign**, not a magnitude:

$$ \operatorname{orient}(p,q,r) \;=\; \begin{vmatrix}
{\color{#1565c0}p_x} & {\color{#1565c0}p_y} & 1\\[2pt]
{\color{#2e7d32}q_x} & {\color{#2e7d32}q_y} & 1\\[2pt]
{\color{#bf616a}r_x} & {\color{#bf616a}r_y} & 1
\end{vmatrix} $$

- One wrong sign flips topology. 💥
- Rational coordinates are common in projective constructions. 📐
- Exact arithmetic gives the sign with **zero tolerance**. ✅

> In projgeom, correctness is binary — and fractions keep it binary. 🧷

---

## 📏 Part 2 — How We Measure

---

### 🔧 Compiler Flags — a Fair Fight

| Toolchain | Version | Optimization | Notes |
|:--|:--|:--|:--|
| 🐃 g++ | 15.2.0 | `-O3` | full `-std=c++17 -DNDEBUG` |
| 🧑‍💻 clang-cl | LLVM 23.1.3 | `/O2` | forwards to cc1 **`-O3`** (rejects `/O3`) |
| 🪟 MSVC | cl 14.51 (VS 18) | `/O2` | MSVC has no `-O3`; `/O2` is max speed |

- Both classes instantiated as `Fraction<int64_t>` 🧮
- One process, namespace-isolated; nanobench 4.3.11 📏
- `warmup(300)` · `minEpochIterations(50000)` · **median of 3 passes**
- Inputs seeded from `argc` → no constant folding 🎲

<pre>
graph LR
    SEED["runtime seed\n(argc - 1)"] --> IN["inputs\n1/2, 1/3, large pairs"]
    IN --> OP["operation\n+ - * / compare int-op"]
    OP --> SINK["doNotOptimizeAway"]
    SINK --> NB["nanobench\n50k iterations"]
    NB --> MED["median ns/op"]
    style SEED fill:#fff9c4,stroke:#f57f17,stroke-width:3px
    style OP fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style NB fill:#c8e6c9,stroke:#2e7d32,stroke-width:3px
    style MED fill:#d1c4e9,stroke:#6a1b9a,stroke-width:3px
</pre>

]

---

### 🧮 The gcd Tax

$$ \gcd({\color{#1565c0}m},{\color{#2e7d32}n}) \;=\; \gcd({\color{#2e7d32}n},\, {\color{#1565c0}m} \bmod {\color{#2e7d32}n}), \qquad \gcd(m,0)=|m| $$

Every normalization, reduction and comparison eventually spends time in integer modulo.

<pre>
graph LR
    OP["a/b op c/d"] --> N["normalize()"]
    N --> POS["keep_denom_positive()"]
    POS --> RED["reduce()"]
    RED --> G["gcd(numer, denom)"]
    G --> MOD["integer modulo\n~20-40 cycles"]
    style OP fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style G fill:#bf616a,stroke:#8b2c3a,stroke-width:4px,color:#ffffff
    style MOD fill:#f8bbd0,stroke:#c2185b,stroke-width:3px
</pre>

]

Construction alone is **≈ 36 ns** (g++) · **≈ 21.6 ns** (clang-cl) · **≈ 28.8 ns** (MSVC) — almost all of it gcd. 🧱

---

## 💥 Part 3 — Killing Template Recursion

---

### 🧨 The Problem: Template Recursion Can't Unroll

`fractions.hpp` used a **function template** for gcd:

```cpp
template <typename Mn>
CONSTEXPR14 auto gcd_recur(const Mn& m, const Mn& n) -> Mn {
    if (n == 0) return abs(m);
    return gcd_recur(n, m % n);   // recursive call on Mn
}
```

- The recursive call is **dependently typed** — the optimizer sees a call, not a self-loop. 🚫
- `pyfractions.hpp` had the same shape **inside the class template**. 🐍
- Result: the inner `%` chain never gets inlined / partially unrolled. 😴

$$ \gcd(m,n)=\gcd(n,\,m\bmod n) \quad\text{should compile to}\quad
{\color{#2e7d32}\text{idiv} \to \text{test} \to \text{idiv} \to \cdots} $$

---

### 🔨 Fix #1 — Concrete Overloads in `fractions.hpp`

Replace the template with **non-template overloads**, one per integer type:

```cpp
CONSTEXPR14 auto gcd_recur(long long _m, long long _n) -> long long {
    if (_n == 0) return abs(_m);
    return gcd_recur(_n, _m % _n);   // concrete self-call: unrollable
}

CONSTEXPR14 auto gcd(long long _m, long long _n) -> long long {
    if (_m == 0) return abs(_n);
    return gcd_recur(_m, _n);
}
```

Overload set: `unsigned long`, `long`, `unsigned int`, `int`, **`unsigned long long`**, **`long long`** 🧩

- `long long` / `unsigned long long` added so `int64_t` works on Windows/LP64. 🪟
- `template <typename T = int64_t> struct Fraction` — default now `int64_t`. ✅

---

### 🐍 Fix #2 — Move PyFred's gcd Out of the Class Template

The gcd became **free, non-template functions** in `fractions::detail`:

```cpp
namespace fractions {
    namespace detail {
        CONSTEXPR14 auto gcd_recur(long long _m, long long _n) -> long long {
            if (_n == 0) return abs(_m);
            return gcd_recur(_n, _m % _n);   // plain overload, not a template
        }
        CONSTEXPR14 auto gcd(long long _m, long long _n) -> long long {
            if (_m == 0) return abs(_n);
            return gcd_recur(_m, _n);
        }
```

- The class body now calls `detail::gcd(_numerator, _denominator)`. 📞
- No `T` in any gcd signature — nothing left to prevent unrolling. 🎯
- Kept in `detail` so it never collides with `fractions::gcd`. 🧷

<pre>
graph LR
    C["Fraction<int64_t>\nnormalize()"] --> D["detail::gcd(long long, long long)"]
    D --> R["detail::gcd_recur(long long, long long)"]
    R -->|"m % n"| R
    style C fill:#d1c4e9,stroke:#6a1b9a,stroke-width:3px
    style D fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style R fill:#c8e6c9,stroke:#2e7d32,stroke-width:3px
</pre>

]

---

### 🔬 Evidence — g++ `-O3` Assembly of the Constructor

<rect x="416" y="24" width="304" height="150" rx="8" fill="#e8f5e9" stroke="#2e7d32" stroke-width="3"/>
  <text x="432" y="52" font-family="sans-serif" font-size="15" font-weight="bold" fill="#2e7d32">after — concrete overloads</text>
  <g font-family="monospace" font-size="13" fill="#2e3440">
    <text x="432" y="82">gcd_recur(n, m % n)</text>
    <text x="432" y="106">├─ idivq  %rcx</text>
    <text x="432" y="130">└─ 0 call instructions</text>
    <text x="432" y="158" fill="#2e7d32" font-weight="bold">inlined · 3 × idivq</text>
  </g>
</svg>

</p>

> 🧪 `xmake test`: **137 cases / 958 assertions** pass; `clang-format --dry-run -Werror` clean.

---

## ⚡ Part 4 — Six Optimizations

---

### 1️⃣ `operator+`: Skip the Constructor gcd

```cpp
auto denom = da * right;
auto numer = right * na + left * nb;
if (common == 1) {
    return Fraction::from_raw(numer, denom);   // already reduced
}
const auto g = gcd(numer, common);              // only factors of common
return Fraction::from_raw(numer / g, denom / g);
```

$$ \frac{{\color{#1565c0}a}}{{\color{#2e7d32}b}}+\frac{{\color{#b48ead}c}}{{\color{#d08770}d}}
   =\frac{{\color{#1565c0}a}{\color{#d08770}d}+{\color{#2e7d32}b}{\color{#b48ead}c}}
          {{\color{#2e7d32}b}{\color{#d08770}d}} $$

If $\gcd(b,d)=1$, the sum is **already in lowest terms**. ✅

> 🧾 The old `operator+` ran a **second** gcd in the constructor.

---

### 2️⃣ Integer ± : Zero gcds

Because $a/b$ is canonical:

$$ \gcd({\color{#1565c0}a} \pm {\color{#bf616a}n}\,{\color{#2e7d32}b},\; {\color{#2e7d32}b})
   \;=\; \gcd({\color{#1565c0}a},{\color{#2e7d32}b}) \;=\; 1 $$

```cpp
// before: swap + reduce + swap + reduce
// after:
this->_numer += rhs * this->_denom;   // no gcd at all
```

- Overflow behaviour unchanged: the old code already formed $n \cdot b$. 🧮
- PyFraction still routes `+ int` through a full fraction add → the big outlier. 🐍

---

### 3️⃣ Scalar Multiply: One gcd, Not Two

```cpp
std::swap(this->_numer, rhs);
this->reduce();                 // gcd(rhs, denom)
this->_numer *= rhs;
this->keep_denom_positive();    // was normalize(): another gcd
```

After cross-reduction, the second gcd is **always 1**. ✅

The same idea lets PyFraction build results with its private **skip-normalize** constructor:

```cpp
Fraction(T numerator, T denominator, bool skip_normalize);
```

> 📦 `from_raw` (Simple) and `skip_normalize` (Py) are the same trick, two spellings.

---

### 4️⃣ Comparison: Two Multiplies, Zero gcds

$$ \frac{\color{#1565c0}a}{\color{#2e7d32}b} \;<\; \frac{\color{#b48ead}c}{\color{#d08770}d}
   \iff {\color{#1565c0}a}{\color{#d08770}d} \;<\; {\color{#b48ead}c}{\color{#2e7d32}b},
   \qquad {\color{#2e7d32}b},{\color{#d08770}d}>0 $$

<pre>
graph LR
    C["a/b vs c/d"] --> D{"both denominators > 0?"}
    D -->|yes| W{"__SIZEOF_INT128__?"}
    W -->|yes| I["__int128 cross product\nGCC · Clang · clang-cl"]
    W -->|no, MSVC x64| U["_umul128 mirror\nhi/lo comparison"]
    D -->|no: 0 or NaN| O["legacy reduce path"]
    style C fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style I fill:#c8e6c9,stroke:#2e7d32,stroke-width:3px
    style U fill:#c8e6c9,stroke:#2e7d32,stroke-width:3px
    style O fill:#fff9c4,stroke:#f57f17,stroke-width:3px
</pre>

]

---

### 5️⃣ MSVC: Replace Division-Based Overflow Checks

`_mul_overflows` used `max / |b|` — one integer division **per check**, twice per multiply.

```cpp
const unsigned long long lo = _umul128(_abs_u64(a), _abs_u64(b), &hi);
const unsigned long long limit = negative ? positive_max + 1 : positive_max;
return hi != 0 || lo > limit;
```

- `_umul128` is exact and branch-light on MSVC x64. ⚙️
- Keeps the existing `double` fallback for genuine overflow. 🛟

---

### 6️⃣ Portability: `__int128` Under `-Wpedantic -Werror`

```cpp
#if defined(__SIZEOF_INT128__)
    __extension__ typedef __int128 wide;   // GCC/Clang extension, silenced
    return static_cast<wide>(an) * bd < static_cast<wide>(bn) * ad;
#elif defined(_MSC_VER) && defined(_M_X64)
    // _umul128 sign-magnitude comparison
#else
    // legacy reduce-based path
#endif
```

- `using wide = __int128;` **fails** `-Wpedantic -Werror`. ❌
- `__extension__ typedef` is clean. ✅
- `<intrin.h>` only when `_MSC_VER && _M_X64`. 🪟

---

### 🙃 What Did **Not** Work

<text x="128" y="100" text-anchor="end">mul small — with</text>
    <rect x="138" y="86" width="354" height="20" fill="#bf616a" rx="3"/>
    <text x="500" y="101" fill="#bf616a" font-weight="bold">35.4 · +25%</text>

<text x="128" y="140" text-anchor="end">mul large — without</text>
    <rect x="138" y="126" width="310" height="20" fill="#a3be8c" rx="3"/>
    <text x="456" y="141" fill="#2e7d32" font-weight="bold">31.0</text>

<text x="128" y="170" text-anchor="end">mul large — with</text>
    <rect x="138" y="156" width="371" height="20" fill="#bf616a" rx="3"/>
    <text x="517" y="171" fill="#bf616a" font-weight="bold">37.1 · +20%</text>
  </g>
  <text x="16" y="214" font-family="sans-serif" font-size="13" fill="#4c566a">Two extra branches × two reduce() calls per multiply = measurable regression.</text>
  <text x="16" y="236" font-family="sans-serif" font-size="13" fill="#bf616a" font-weight="bold">Hypothesis rejected. We kept the original reduce().</text>
</svg>

</p>

---

### 📦 Applied Patch

> 🔒 Public API preserved (a defaulted `T` is additive). 🧪 Unit tests untouched.

---

## 📊 Part 5 — Head-to-Head

---

### 🏁 Fraction vs PyFraction — median ns/op

.font-sm[
| Operation | 🐃 g++ -O3 | 🧑‍💻 clang-cl /O2 | 🪟 MSVC /O2 |
|:--|--:|--:|--:|
| construction | 35.95 / 35.93 | 21.62 / 21.61 | 28.76 / 28.79 |
| `a + b` large | 35.97 / **21.63** | 22.33 / **13.03** | 35.92 / **21.63** |
| `a * b` large | 35.95 / 36.09 | 21.61 / 21.60 | 35.94 / 36.16 |
| `a / b` small | 28.78 / **21.60** | 17.30 / **12.97** | 29.38 / **21.61** |
| `a < b` large | **2.64** / 3.18 | 5.04 / **4.73** | 7.21 / **6.48** |
| `a == b` | 1.46 / 1.46 | 2.17 / 2.16 | 2.16 / 2.16 |
| `a + 7` | **1.20** / 1.88 | **2.16** / 8.64 | **1.57** / 13.69 |
| `a * 7` | 14.42 / 14.41 | 8.66 / 8.65 | **15.12** / 21.70 |
| complex expr | 174.19 / 173.63 | **108.14** / 123.88 | **159.63** / 179.39 |
]

> 🐍 PyFraction wins add/sub/div. 🧰 Fraction wins integer ± and most comparisons. Construction is a tie.

---

### 📊 PyFraction / Fraction — ratio (lower = Py faster)

<text x="40" y="52" font-family="sans-serif" font-size="14" font-weight="bold" fill="#1565c0">g++ -O3</text>
  <line x1="120" y1="60" x2="120" y2="214" stroke="#4c566a" stroke-width="1.5"/>
  <g font-family="sans-serif" font-size="12" fill="#2e3440">
    <text x="112" y="80" text-anchor="end">add</text>
    <rect x="60" y="70" width="42" height="15" fill="#2e7d32" rx="2"/>
    <text x="108" y="81" fill="#2e7d32" font-weight="bold">0.60</text>
    <text x="112" y="116" text-anchor="end">div</text>
    <rect x="60" y="106" width="52" height="15" fill="#2e7d32" rx="2"/>
    <text x="118" y="117" fill="#2e7d32" font-weight="bold">0.75</text>
    <text x="112" y="152" text-anchor="end">lt</text>
    <rect x="60" y="142" width="83" height="15" fill="#bf616a" rx="2"/>
    <text x="149" y="153" fill="#bf616a" font-weight="bold">1.19</text>
    <text x="112" y="188" text-anchor="end">complex</text>
    <rect x="60" y="178" width="70" height="15" fill="#4c566a" rx="2"/>
    <text x="136" y="189" fill="#4c566a" font-weight="bold">1.00</text>
  </g>

<text x="300" y="52" font-family="sans-serif" font-size="14" font-weight="bold" fill="#6a1b9a">clang-cl</text>
  <line x1="400" y1="60" x2="400" y2="214" stroke="#4c566a" stroke-width="1.5"/>
  <g font-family="sans-serif" font-size="12" fill="#2e3440">
    <text x="392" y="80" text-anchor="end">add</text>
    <rect x="340" y="70" width="41" height="15" fill="#2e7d32" rx="2"/>
    <text x="387" y="81" fill="#2e7d32" font-weight="bold">0.58</text>
    <text x="392" y="116" text-anchor="end">div</text>
    <rect x="340" y="106" width="52" height="15" fill="#2e7d32" rx="2"/>
    <text x="398" y="117" fill="#2e7d32" font-weight="bold">0.75</text>
    <text x="392" y="152" text-anchor="end">lt</text>
    <rect x="340" y="142" width="65" height="15" fill="#2e7d32" rx="2"/>
    <text x="411" y="153" fill="#2e7d32" font-weight="bold">0.93</text>
    <text x="392" y="188" text-anchor="end">complex</text>
    <rect x="340" y="178" width="80" height="15" fill="#bf616a" rx="2"/>
    <text x="426" y="189" fill="#bf616a" font-weight="bold">1.15</text>
  </g>

<text x="570" y="52" font-family="sans-serif" font-size="14" font-weight="bold" fill="#b48ead">MSVC</text>
  <line x1="670" y1="60" x2="670" y2="214" stroke="#4c566a" stroke-width="1.5"/>
  <g font-family="sans-serif" font-size="12" fill="#2e3440">
    <text x="662" y="80" text-anchor="end">add</text>
    <rect x="610" y="70" width="42" height="15" fill="#2e7d32" rx="2"/>
    <text x="658" y="81" fill="#2e7d32" font-weight="bold">0.60</text>
    <text x="662" y="116" text-anchor="end">div</text>
    <rect x="610" y="106" width="52" height="15" fill="#2e7d32" rx="2"/>
    <text x="668" y="117" fill="#2e7d32" font-weight="bold">0.74</text>
    <text x="662" y="152" text-anchor="end">lt</text>
    <rect x="610" y="142" width="63" height="15" fill="#2e7d32" rx="2"/>
    <text x="679" y="153" fill="#2e7d32" font-weight="bold">0.90</text>
    <text x="662" y="188" text-anchor="end">complex</text>
    <rect x="610" y="178" width="78" height="15" fill="#bf616a" rx="2"/>
    <text x="694" y="189" fill="#bf616a" font-weight="bold">1.12</text>
  </g>

<text x="16" y="238" font-family="sans-serif" font-size="12" fill="#4c566a">Green = Py faster · Red = Py slower · off-chart: `+ 7` = 1.57× (g++) · 4.0× (clang-cl) · 8.7× (MSVC).</text>
</svg>

</p>

---

### 🧭 What Wins Where

<pre>
graph TD
    Q["Which operation?"] --> F["fraction op: + - * /"]
    Q --> I["fraction op with int"]
    Q --> C["compare"]
    F --> FW["PyFraction<br/>add/sub/div ~0.58-0.75×"]
    I --> IW["Fraction<br/>± int 1.3-8.7× faster"]
    C --> CW["Fraction on g++<br/>PyFraction on clang-cl/MSVC"]
    style Q fill:#eceff4,stroke:#4c566a,stroke-width:3px
    style FW fill:#d1c4e9,stroke:#6a1b9a,stroke-width:3px
    style IW fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style CW fill:#fff9c4,stroke:#f57f17,stroke-width:3px
</pre>

]

- 🐍 **PyFraction**: fraction `+`, `-`, `/` win on every compiler; also `lt` on clang-cl/MSVC.
- 🧰 **Fraction**: integer `±` dominates everywhere; `lt` on g++; `complex` on clang-cl/MSVC.
- ⚖️ **Ties**: construction, `*`, `==`.
- 🧠 Compiler-dependent — the refactor made construction a tie; arithmetic paths still differ.

---

## 🎓 Part 6 — Size, Verification, Lessons

---

### 📏 Size: Lines of Code

| Header | total | code | comment | blank | size |
|--------|------:|-----:|--------:|------:|-----:|
| `fractions.hpp` 🧰 | 1545 | 497 | **937** | 111 | 54 KB |
| `extfractions.hpp` 🪶 | 1225 | 397 | 748 | 80 | 45 KB |
| `pyfractions.hpp` 🐍 | 1127 | **647** | 334 | 146 | 41 KB |

- 🐍 **PyFraction is the shortest file but the densest code** — 1.3× more executable lines than Fraction.
- 🧰 **Fraction is the longest file** — Doxygen math diagrams (Graphviz + LaTeX) dominate; 497 lines are code.
- 🔒 No public API removed; `Fraction<T = int64_t>` is additive.

> 📐 Fewer total lines ≠ less code. Compare *code* lines, not file length. 🧠

---

### 🧪 Verification Pipeline

<pre>
graph TD
    T["xmake test\n137 cases · 958 assertions"] --> OK["✅"]
    A["g++ -O3 assembly\n0 calls to gcd"] --> OK
    D["bench: 3 passes\nmedian ns/op"] --> OK
    G["GCC 15.2\n-Wpedantic -Wconversion -Werror"] --> OK
    F["clang-format 22\n--dry-run --Werror"] --> OK
    style OK fill:#c8e6c9,stroke:#2e7d32,stroke-width:4px
    style T fill:#e3f2fd,stroke:#1565c0,stroke-width:3px
    style A fill:#d1c4e9,stroke:#6a1b9a,stroke-width:3px
    style D fill:#fff9c4,stroke:#f57f17,stroke-width:3px
</pre>

]

- Special cases: `0/0`, `±1/0`, negative denominators, both-negative comparisons. 🧷
- Tests untouched; public API unchanged. 🔒
- Benchmarks in one process — no cross-binary noise. 📏

---

### 🎓 Lessons

1. 🧮 **gcd is the denominator of every cost.**
2. 💥 **A recursive function template is not unrollable — overloads are.**
3. 🐍 **Moving gcd out of the class template is what made Py construction a tie.**
4. 📏 **Match optimization flags** — `-O3` vs `/O2`; clang-cl's `/O2` is cc1 `-O3`.
5. 🧪 **Fast paths must preserve NaN/∞/throw semantics.**
6. 🧷 **Don't "fix" a bug a test asserts without permission.**
7. 🧩 **Public API stability is a feature.**

$$ \text{speedup} \;=\; \frac{\text{Fraction}}{\text{PyFraction}},
   \qquad {\color{#2e7d32}\text{add large (MSVC)}}:\;
   \frac{35.92}{21.63} \approx {\color{#2e7d32}1.66\times} $$

> The best optimization is the work you can prove is redundant. ✅

---

## Q&A 🎤

### GitHub: `luk036/fractions-cpp`

Slides: `projgeom/py-vs-ext-remark.html` 🧮

---

# 🏁 Thank You

### Questions, benchmarks, or a PR? 🚀