solver/models/
avellaneda_drift.rs

1use super::control::{ControlProblem, StateDerivatives};
2use super::market_making::MarketMakingControl;
3use super::normal_cdf;
4use crate::numeric::finite_difference::discretization::{DimensionKind, Transport};
5use crate::numeric::finite_difference::pde::PdeProblem;
6
7/// Avellaneda-Stoikov Model with Drift
8///
9/// Adds a drift term $\mu$ to the mid-price process: $dS_t = \mu dt + \sigma dW_t$.
10/// The optimal quotes are skewed based on the expected future price movement.
11#[derive(Clone)]
12pub struct AvellanedaDrift {
13    pub gamma: f64,
14    pub sigma: f64,
15    pub kappa: f64,
16    pub a: f64,
17    /// Price drift parameter ($\mu$).
18    pub mu: f64,
19}
20
21impl AvellanedaDrift {
22    /// Computes optimal spreads, adjusted for drift.
23    pub fn get_spreads(&self, derivs: &StateDerivatives<2>) -> (f64, f64) {
24        let base_spread = (1.0 / self.gamma) * (1.0 + self.gamma / self.kappa).ln();
25
26        let dv_dq_buy = derivs.fwd[0];
27        let dv_dq_sell = derivs.bwd[0];
28
29        let mut delta_bid = base_spread - dv_dq_buy;
30        let mut delta_ask = base_spread + dv_dq_sell;
31
32        let min_spread = -5.0;
33        delta_bid = delta_bid.max(min_spread);
34        delta_ask = delta_ask.max(min_spread);
35
36        (delta_bid, delta_ask)
37    }
38
39    /// The base arrival intensity used for intensity-to-spread conversion.
40    pub fn fill_rate_base(&self, _state: &[f64; 2]) -> f64 {
41        self.a
42    }
43
44    /// The fill-rate decay parameter used for intensity-to-spread conversion.
45    pub fn fill_rate_decay(&self) -> f64 {
46        self.kappa
47    }
48}
49
50impl ControlProblem<2> for AvellanedaDrift {
51    type Control = MarketMakingControl;
52
53    fn optimize(&self, _t: f64, _state: &[f64; 2], derivs: &StateDerivatives<2>) -> Self::Control {
54        let base_spread = (1.0 / self.gamma) * (1.0 + self.gamma / self.kappa).ln();
55        let delta_bid = (base_spread - derivs.fwd[0]).max(-5.0);
56        let delta_ask = (base_spread + derivs.bwd[0]).max(-5.0);
57
58        MarketMakingControl::new(
59            self.a * (-self.kappa * delta_bid).exp(),
60            self.a * (-self.kappa * delta_ask).exp(),
61        )
62    }
63
64    fn running_reward(&self, _t: f64, _state: &[f64; 2], control: &Self::Control) -> f64 {
65        // Optimized Hamiltonian H* = (lambda_bid + lambda_ask) / (gamma + kappa).
66        (control.bid_intensity + control.ask_intensity) / (self.gamma + self.kappa)
67    }
68
69    fn bsde_driver(
70        &self,
71        t: f64,
72        state: &[f64; 2],
73        control: &Self::Control,
74        derivs: &StateDerivatives<2>,
75        _dt: f64,
76    ) -> f64 {
77        self.driver(t, state, control, derivs)
78    }
79
80    fn generator(
81        &self,
82        _t: f64,
83        state: &[f64; 2],
84        _control: &Self::Control,
85        _derivs: &StateDerivatives<2>,
86    ) -> f64 {
87        // Risk penalty plus drift effect q * mu.
88        let q = state[0];
89        -0.5 * self.gamma * self.sigma.powi(2) * q.powi(2) + q * self.mu
90    }
91
92    fn terminal(&self, _state: &[f64; 2]) -> f64 {
93        0.0
94    }
95
96    fn discount_rate(&self, _state: &[f64; 2]) -> f64 {
97        0.0
98    }
99
100    fn constant_discount_rate(&self) -> Option<f64> {
101        Some(0.0)
102    }
103
104    fn next_step(&self, _t: f64, state: &[f64; 2], dt: f64, noise: &[f64; 2]) -> [f64; 2] {
105        let mut next = *state;
106
107        let u = normal_cdf(noise[0]);
108        let base_spread = (1.0 / self.gamma) * (1.0 + self.gamma / self.kappa).ln();
109        let delta_bid = base_spread;
110        let delta_ask = base_spread;
111        let lambda_bid = self.a * (-self.kappa * delta_bid).exp();
112        let lambda_ask = self.a * (-self.kappa * delta_ask).exp();
113        let p_bid = (lambda_bid * dt).clamp(0.0, 1.0);
114        let p_ask = (lambda_ask * dt).clamp(0.0, 1.0);
115        if u < p_bid {
116            next[0] += 1.0;
117        } else if u > 1.0 - p_ask {
118            next[0] -= 1.0;
119        }
120
121        next[1] += self.mu * dt + self.sigma * dt.sqrt() * noise[1];
122        next
123    }
124
125    fn is_reduced_value(&self) -> bool {
126        true
127    }
128
129    fn is_diffusion_dimension(&self, dim: usize) -> bool {
130        dim == 1
131    }
132
133    fn gradient_step(&self, _dim: usize) -> f64 {
134        1.0
135    }
136}
137
138impl PdeProblem<2> for AvellanedaDrift {
139    fn dimension_kind(&self, dim: usize) -> DimensionKind {
140        match dim {
141            0 => DimensionKind::DiscreteJump,
142            _ => DimensionKind::Diffusion,
143        }
144    }
145
146    fn transport(
147        &self,
148        _t: f64,
149        state: &[f64; 2],
150        control: &Self::Control,
151        derivs: &StateDerivatives<2>,
152    ) -> Transport<2> {
153        let q = state[0];
154        let hamiltonian =
155            (control.bid_intensity + control.ask_intensity) / (self.gamma + self.kappa);
156        let jump_transport =
157            control.bid_intensity * derivs.fwd[0] - control.ask_intensity * derivs.bwd[0];
158        let local = -0.5 * self.gamma * self.sigma.powi(2) * q.powi(2) + q * self.mu;
159        Transport::new(
160            [control.bid_intensity, 0.0],
161            [control.ask_intensity, 0.0],
162            hamiltonian + local - jump_transport,
163        )
164    }
165}
solver/models/avellaneda_drift.rs

solver/models/
avellaneda_drift.rs
