Skip to main content

kestrel_chartkit/valuation/
volatility.rs

1//! Volatility over strike and maturity.
2
3#[cfg(feature = "serde")]
4use serde::{Deserialize, Serialize};
5
6use crate::finance::{year_fraction, Date, DayCountConvention};
7
8use super::ValuationContextError;
9
10/// A volatility surface given as a grid of quoted points.
11///
12/// Maturities run along one axis and strikes along the other, and every cell holds the volatility
13/// quoted there. Between the grid lines the surface is bilinear: linear in maturity, linear in
14/// strike, in that order. Outside them it is held flat at the edge value — the same convention
15/// the yield curve uses, stated here for the same reason: it is an assumption, and an assumption
16/// that is written down can be argued with.
17///
18/// The grid is deliberately the whole model. A parametric fit — a smile function with its own
19/// coefficients — would produce a smooth surface everywhere, including where nothing was quoted,
20/// and its shape would come from the parametrisation rather than from the market. Interpolating
21/// what was actually quoted keeps the two apart.
22///
23/// `validity` bounds how far a caller may ask beyond the quoted range before the surface refuses
24/// to answer. Flat extrapolation one strike beyond the grid is a mild assumption; ten strikes
25/// beyond it is a fabrication, and the difference should not be silent.
26#[derive(Debug, Clone, PartialEq)]
27#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
28pub struct VolatilitySurface {
29    reference: Date,
30    day_count: DayCountConvention,
31    /// Ascending maturities in years from the reference date.
32    maturities: Vec<f64>,
33    /// Ascending strikes.
34    strikes: Vec<f64>,
35    /// `volatilities[maturity index][strike index]`.
36    volatilities: Vec<Vec<f64>>,
37    validity: SurfaceValidity,
38}
39
40/// How far outside the quoted grid the surface still answers.
41#[derive(Debug, Clone, Copy, PartialEq)]
42#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
43pub struct SurfaceValidity {
44    /// Multiples of the quoted strike range that may be exceeded on either side. `0.0` refuses
45    /// any strike outside the grid, `f64::INFINITY` never refuses.
46    pub strike_tolerance: f64,
47    /// Years beyond the last quoted maturity that may still be asked for.
48    pub maturity_tolerance_years: f64,
49}
50
51impl Default for SurfaceValidity {
52    fn default() -> Self {
53        // A tenth of the quoted strike width and a quarter year: close enough to the grid that
54        // holding the edge value flat is a defensible reading of the quotes.
55        Self {
56            strike_tolerance: 0.1,
57            maturity_tolerance_years: 0.25,
58        }
59    }
60}
61
62impl VolatilitySurface {
63    /// Builds a surface from a grid.
64    ///
65    /// `maturities` are years from `reference`, `strikes` are prices, and `volatilities` is one
66    /// row per maturity. Both axes must be finite, strictly ascending and non-empty, every
67    /// volatility finite and non-negative, and the grid rectangular.
68    pub fn new(
69        reference: Date,
70        maturities: Vec<f64>,
71        strikes: Vec<f64>,
72        volatilities: Vec<Vec<f64>>,
73        day_count: DayCountConvention,
74        validity: SurfaceValidity,
75    ) -> Result<Self, ValuationContextError> {
76        if maturities.is_empty() || strikes.is_empty() {
77            return Err(ValuationContextError::InvalidCurve(
78                "a surface needs at least one maturity and one strike",
79            ));
80        }
81        if ascending_violation(&maturities) || ascending_violation(&strikes) {
82            return Err(ValuationContextError::InvalidCurve(
83                "maturities and strikes must be finite and strictly ascending",
84            ));
85        }
86        if maturities.iter().any(|t| *t < 0.0) || strikes.iter().any(|k| *k <= 0.0) {
87            return Err(ValuationContextError::InvalidCurve(
88                "maturities must be non-negative and strikes positive",
89            ));
90        }
91        if volatilities.len() != maturities.len()
92            || volatilities.iter().any(|row| row.len() != strikes.len())
93        {
94            return Err(ValuationContextError::InvalidCurve(
95                "the volatility grid must be rectangular and match both axes",
96            ));
97        }
98        if volatilities
99            .iter()
100            .any(|row| row.iter().any(|v| !v.is_finite() || *v < 0.0))
101        {
102            return Err(ValuationContextError::InvalidCurve(
103                "volatilities must be finite and non-negative",
104            ));
105        }
106
107        Ok(Self {
108            reference,
109            day_count,
110            maturities,
111            strikes,
112            volatilities,
113            validity,
114        })
115    }
116
117    /// A surface with one volatility everywhere — the compatible special case, so a caller who
118    /// only has one number keeps working.
119    pub fn flat(
120        reference: Date,
121        volatility: f64,
122        day_count: DayCountConvention,
123    ) -> Result<Self, ValuationContextError> {
124        Self::new(
125            reference,
126            vec![0.0],
127            vec![1.0],
128            vec![vec![volatility]],
129            day_count,
130            SurfaceValidity {
131                strike_tolerance: f64::INFINITY,
132                maturity_tolerance_years: f64::INFINITY,
133            },
134        )
135    }
136
137    pub fn reference_date(&self) -> Date {
138        self.reference
139    }
140
141    /// Volatility for a maturity in years and a strike.
142    pub fn volatility_at(&self, time: f64, strike: f64) -> Result<f64, ValuationContextError> {
143        if !time.is_finite() || time < 0.0 || !strike.is_finite() || strike <= 0.0 {
144            return Err(ValuationContextError::TimeOutsideCurve);
145        }
146        self.check_validity(time, strike)?;
147
148        let (row_low, row_high, row_weight) = bracket(&self.maturities, time);
149        let (col_low, col_high, col_weight) = bracket(&self.strikes, strike);
150
151        let interpolate = |row: usize| {
152            let low = self.volatilities[row][col_low];
153            let high = self.volatilities[row][col_high];
154            low + col_weight * (high - low)
155        };
156        let low = interpolate(row_low);
157        let high = interpolate(row_high);
158        Ok(low + row_weight * (high - low))
159    }
160
161    /// Volatility for an expiry date and a strike.
162    pub fn volatility(&self, expiry: Date, strike: f64) -> Result<f64, ValuationContextError> {
163        if expiry < self.reference {
164            return Err(ValuationContextError::TimeOutsideCurve);
165        }
166        self.volatility_at(
167            year_fraction(self.reference, expiry, self.day_count),
168            strike,
169        )
170    }
171
172    fn check_validity(&self, time: f64, strike: f64) -> Result<(), ValuationContextError> {
173        let last_maturity = self.maturities[self.maturities.len() - 1];
174        if time > last_maturity + self.validity.maturity_tolerance_years {
175            return Err(ValuationContextError::OutsideSurfaceValidity);
176        }
177
178        let (first_strike, last_strike) = (self.strikes[0], self.strikes[self.strikes.len() - 1]);
179        let width = (last_strike - first_strike).max(f64::MIN_POSITIVE);
180        let allowance = width * self.validity.strike_tolerance;
181        if strike < first_strike - allowance || strike > last_strike + allowance {
182            return Err(ValuationContextError::OutsideSurfaceValidity);
183        }
184        Ok(())
185    }
186}
187
188fn ascending_violation(values: &[f64]) -> bool {
189    values.iter().any(|v| !v.is_finite()) || values.windows(2).any(|w| w[0] >= w[1])
190}
191
192/// Neighbouring indices and the weight between them; outside the axis both indices are the edge,
193/// which is what holds the value flat.
194fn bracket(axis: &[f64], value: f64) -> (usize, usize, f64) {
195    if value <= axis[0] {
196        return (0, 0, 0.0);
197    }
198    let last = axis.len() - 1;
199    if value >= axis[last] {
200        return (last, last, 0.0);
201    }
202    let upper = axis.partition_point(|entry| *entry <= value).max(1);
203    let lower = upper - 1;
204    let weight = (value - axis[lower]) / (axis[upper] - axis[lower]);
205    (lower, upper, weight)
206}