Math64x64Extended

The Math64x64Extended library provides high-precision mathematical functions using 64.64 fixed-point arithmetic. It extends the ABDK Math library with normal distribution functions essential for Black-Scholes pricing.

Overview

The library provides two categories of functions:

Category	Functions
ABDK Wrappers	`ln`, `exp`, `sqrt`, `add`, `sub`, `mul`, `div`, `neg`, `abs`
Normal Distribution	`normalPdf`, `normalCdf` (Abramowitz-Stegun approximation)

Why 64.64 Fixed-Point?

┌─────────────────────────────────────────────────────────────────┐
│                    64.64 Format (int128)                        │
├─────────────────────────────────────────────────────────────────┤
│  [  64 bits integer part  ] [  64 bits fractional part  ]       │
│                                                                 │
│  • ~18 decimal digits of precision                              │
│  • Range: approximately ±9.2 × 10^18                            │
│  • Deterministic (no floating-point rounding)                   │
│  • Gas-efficient compared to arbitrary precision                │
└─────────────────────────────────────────────────────────────────┘

Advantages over floating-point:

Deterministic: Same input always produces same output across all nodes
No rounding surprises: Fixed-point arithmetic is predictable
Sufficient precision: 18 decimals exceeds financial requirements

Key Concepts

64.64 Representation

A value is represented as: value = int128_representation / 2^64

Decimal Value	64.64 Representation	Calculation
1.0	`0x10000000000000000`	1 × 2^64
0.5	`0x8000000000000000`	0.5 × 2^64
2.0	`0x20000000000000000`	2 × 2^64
0.3989... (1/√2π)	`7359186146747302912`	(1/√2π) × 2^64

Normal Distribution

The cumulative distribution function (CDF) is critical for Black-Scholes:

N(x) = \frac{1}{\sqrt{2\pi}} \int_{-\infty}^{x} e^{-t^2/2} dt

Since this integral has no closed-form solution, we use the Abramowitz-Stegun polynomial approximation.

Constants

Fixed-Point Constants

int128 internal constant ONE = 0x10000000000000000;      // 1.0
int128 internal constant TWO = 0x20000000000000000;      // 2.0
int128 internal constant HALF = 0x8000000000000000;      // 0.5
int128 internal constant NEG_HALF = -0x8000000000000000; // -0.5
int128 internal constant INV_SQRT_2PI_64X64 = 7_359_186_146_747_302_912;  // 1/√2π

Abramowitz-Stegun Coefficients

// Reference: Abramowitz & Stegun, Handbook of Mathematical Functions, Eq. 26.2.17
// N(x) ≈ 1 - n(x)(a1*t + a2*t² + a3*t³ + a4*t⁴ + a5*t⁵) for x ≥ 0
// where t = 1/(1 + p*x), error < 7.5e-8

int128 internal constant AS_P_64X64 = 4_273_038_846_047_820_800;    // p = 0.2316419
int128 internal constant AS_A1_64X64 = 5_891_549_345_779_789_824;   // a1 = 0.319381530
int128 internal constant AS_A2_64X64 = -6_577_440_832_507_964_416;  // a2 = -0.356563782
int128 internal constant AS_A3_64X64 = 32_862_467_576_799_068_160;  // a3 = 1.781477937
int128 internal constant AS_A4_64X64 = -33_596_242_918_879_592_448; // a4 = -1.821255978
int128 internal constant AS_A5_64X64 = 24_539_231_939_563_106_304;  // a5 = 1.330274429

Normal Distribution Functions

`normalPdf`

Standard normal probability density function:

function normalPdf(int128 x) internal pure returns (int128)

n(x) = \frac{1}{\sqrt{2\pi}} e^{-x^2/2}

Use case: Used in Gamma and Theta calculations.

Implementation:

// Calculate -x²/2
int128 xSquared = mul(x, x);
int128 negHalfXSquared = mul(NEG_HALF, xSquared);

// Calculate e^(-x²/2)
int128 expPart = exp(negHalfXSquared);

// Multiply by 1/√2π
return mul(INV_SQRT_2PI_64X64, expPart);

`normalCdf`

Standard normal cumulative distribution function:

function normalCdf(int128 x) internal pure returns (int128)

N(x) = \Phi(x) = \int_{-\infty}^{x} n(t) dt

Approximation: Abramowitz-Stegun equation 26.2.17

Maximum error: < 7.5 × 10^-8
Polynomial approximation for efficiency

Implementation overview:

For $x \geq 0$ :

N(x) \approx 1 - n(x) \cdot (a_1t + a_2t^2 + a_3t^3 + a_4t^4 + a_5t^5)

Where $t = \frac{1}{1 + px}$ and $p = 0.2316419$

For $x < 0$ : Use symmetry $N(-x) = 1 - N(x)$

function normalCdf(int128 x) internal pure returns (int128) {
    // Handle negative x: N(-x) = 1 - N(x)
    bool isNegative = x < 0;
    if (isNegative) {
        x = neg(x);
    }

    // Calculate t = 1 / (1 + p*x)
    int128 t = div(ONE, add(ONE, mul(AS_P_64X64, x)));

    // Calculate polynomial using powers of t
    int128 t2 = mul(t, t);
    int128 t3 = mul(t2, t);
    int128 t4 = mul(t3, t);
    int128 t5 = mul(t4, t);

    int128 poly = mul(AS_A1_64X64, t);
    poly = add(poly, mul(AS_A2_64X64, t2));
    poly = add(poly, mul(AS_A3_64X64, t3));
    poly = add(poly, mul(AS_A4_64X64, t4));
    poly = add(poly, mul(AS_A5_64X64, t5));

    // N(x) = 1 - n(x) * poly
    int128 result = sub(ONE, mul(normalPdf(x), poly));

    // Apply symmetry for negative input
    if (isNegative) {
        result = sub(ONE, result);
    }

    return result;
}

ABDK Wrapper Functions

Transcendental Functions

`ln`

Natural logarithm:

function ln(int128 x) internal pure returns (int128)

Requirements: x > 0

Example:

int128 ratio = div(spot, strike);  // S/K
int128 logRatio = ln(ratio);       // ln(S/K)

`exp`

Exponential function:

function exp(int128 x) internal pure returns (int128)

Example:

int128 negRt = neg(mul(rate, timeToExpiry));  // -rT
int128 discount = exp(negRt);                  // e^(-rT)

`sqrt`

Square root:

function sqrt(int128 x) internal pure returns (int128)

Requirements: x >= 0

Example:

int128 sqrtT = sqrt(timeToExpiry);            // √T
int128 volSqrtT = mul(volatility, sqrtT);     // σ√T

Arithmetic Functions

`add` / `sub`

Addition and subtraction:

function add(int128 x, int128 y) internal pure returns (int128)
function sub(int128 x, int128 y) internal pure returns (int128)

`mul` / `div`

Multiplication and division:

function mul(int128 x, int128 y) internal pure returns (int128)
function div(int128 x, int128 y) internal pure returns (int128)

Note: These handle the 64.64 scaling correctly—no manual adjustment needed.

`neg` / `abs`

Negation and absolute value:

function neg(int128 x) internal pure returns (int128)
function abs(int128 x) internal pure returns (int128)

Conversion Functions

`fromInt` / `fromUInt`

Convert integers to 64.64:

function fromInt(int256 x) internal pure returns (int128)
function fromUInt(uint256 x) internal pure returns (int128)

Example:

int128 two = fromInt(2);      // 2.0 in 64.64
int128 hundred = fromUInt(100); // 100.0 in 64.64

`toInt` / `toUInt`

Convert 64.64 to integers (rounds down):

function toInt(int128 x) internal pure returns (int64)
function toUInt(int128 x) internal pure returns (uint64)

Calculation Example

d1 Calculation in Black-Scholes

d_1 = \frac{\ln(S/K) + (r + \sigma^2/2) \cdot T}{\sigma \cdot \sqrt{T}}

// Given: spot, strike, rate, volatility, timeToExpiry (all 64.64)

// Step 1: σ√T
int128 sqrtT = sqrt(timeToExpiry);
int128 volSqrtT = mul(volatility, sqrtT);

// Step 2: ln(S/K)
int128 ratio = div(spot, strike);
int128 logRatio = ln(ratio);

// Step 3: σ²/2
int128 volSquared = mul(volatility, volatility);
int128 volSquaredHalf = mul(volSquared, HALF);

// Step 4: (r + σ²/2) × T
int128 rateAdjusted = add(rate, volSquaredHalf);
int128 rateAdjustedT = mul(rateAdjusted, timeToExpiry);

// Step 5: d1 = [ln(S/K) + (r + σ²/2)T] / (σ√T)
int128 numerator = add(logRatio, rateAdjustedT);
int128 d1 = div(numerator, volSqrtT);

Gas Costs

Function	Approximate Gas
`add` / `sub`	~100
`mul` / `div`	~200
`sqrt`	~3,000
`ln`	~5,000
`exp`	~5,000
`normalPdf`	~6,000
`normalCdf`	~20,000

Note: normalCdf is the most expensive operation, called 2-4 times per option pricing.

Security Considerations

Overflow/Underflow

The 64.64 format has limited range:

Integer part: 64 bits signed (~±9.2 × 10^18)
Practical limit: Values beyond 10^15 may cause issues

Safe input ranges:

Input	Safe Range
Spot/Strike	0.001 to 10^12
Volatility	0.01 to 5.0 (1% to 500%)
Rate	-0.1 to 0.5 (-10% to 50%)
Time	0.001 to 10 years

Approximation Error

The Abramowitz-Stegun approximation has maximum error < 7.5 × 10^-8, which is:

Sufficient for financial calculations (better than 0.00001% error)
Deterministic across all nodes

Edge Cases

Input	Behavior
`normalCdf(0)`	Returns 0.5 (HALF)
`normalCdf(very large)`	Approaches 1.0
`normalCdf(very negative)`	Approaches 0.0
`ln(0)`	Reverts
`sqrt(negative)`	Reverts

Integration Points

Used By

Library/Contract	Purpose
BlackScholesLib	d1/d2 calculation, Greeks
MarginEngine	Stress scenario pricing

Depends On

Library	Purpose
ABDKMath64x64	Core 64.64 arithmetic
Constants	64.64 constant values

Code Reference

Source: packages/contracts/src/utils/Math64x64Extended.sol

Dependency: ABDK Libraries Solidity

int128 constant ONE_64X64 = 0x10000000000000000;
int128 constant TWO_64X64 = 0x20000000000000000;
int128 constant HALF_64X64 = 0x8000000000000000;
int128 constant NEG_HALF_64X64 = -0x8000000000000000;

BlackScholesLib — Uses these functions for pricing
Black-Scholes (Math) — Mathematical background
Normal Distribution — Why N(x) matters

Math64x64Extended

On this page