pcli/dex_utils/replicate/
xyk.rs

1use crate::dex_utils::replicate::math_utils;
2use anyhow::Context;
3use ndarray::Array;
4use penumbra_sdk_asset::Value;
5use penumbra_sdk_dex::{
6    lp::{position::Position, Reserves},
7    DirectedUnitPair,
8};
9use penumbra_sdk_num::{fixpoint::U128x128, Amount};
10use rand_core::OsRng;
11use tracing::field;
12
13/// The number of positions that is used to replicate the xyk CFMM.
14pub(crate) const NUM_POOLS_PRECISION: usize = 30;
15
16/// Maximum number of iteration that we allow GS to perform.
17const GAUS_SEIDEL_MAX_ITERATION: usize = 10_000;
18
19/// Sample a range of points around a given price
20pub fn sample_prices(current_price: f64, num_points: usize) -> Vec<f64> {
21    crate::dex_utils::replicate::math_utils::sample_to_upper(3.0 * current_price, num_points)
22}
23
24#[tracing::instrument(name = "replicate_xyk")]
25pub fn replicate(
26    pair: &DirectedUnitPair,
27    raw_r1: &Value,
28    current_price: U128x128,
29    fee_bps: u32,
30) -> anyhow::Result<Vec<Position>> {
31    // First, we find the pool invariant using human display units. This means that we
32    // only need to care about scaling into proper denom units right before posting the
33    // positions. On the other hand, we have to unscale the values that we are given.
34    let fp_raw_r1 = U128x128::from(raw_r1.amount.value());
35    let r1_scaling_factor = U128x128::from(pair.start.unit_amount());
36
37    let fp_r1 = (fp_raw_r1 / r1_scaling_factor).context("scaling factor can't be 0")?;
38    let fp_r2 = (current_price * fp_r1).context("should not overflow when multiplying by price")?;
39
40    tracing::debug!(
41        %fp_r1,
42        %fp_r2,
43        "computed respective quantities"
44    );
45
46    let xyk_invariant = (fp_r1 * fp_r2).expect("no overflow when computing curve invariant!");
47
48    let xyk_invariant: f64 = xyk_invariant.try_into()?;
49    tracing::debug!(?xyk_invariant, "computed the total invariant for the PVF");
50
51    let f64_current_price: f64 = current_price.try_into()?;
52
53    let alphas = sample_prices(f64_current_price, NUM_POOLS_PRECISION);
54
55    alphas
56        .iter()
57        .enumerate()
58        .for_each(|(i, alpha)| tracing::debug!(i, alpha, "sampled tick"));
59
60    // TODO(erwan): unused for now, but next refactor will rip out `solve` internals to
61    // take this vector of solutions as an argument so that we can more easily recover from
62    // working with non-singular matrices etc.
63    let _b: Vec<f64> = alphas
64        .iter()
65        .map(|price: &f64| portfolio_value_function(xyk_invariant, *price))
66        .collect();
67
68    let position_ks = solve(&alphas, xyk_invariant, NUM_POOLS_PRECISION)?.to_vec();
69
70    position_ks
71        .iter()
72        .enumerate()
73        .for_each(|(i, pool_invariant)| tracing::debug!(i, pool_invariant, "found solution"));
74
75    let unit_start = pair.start.unit_amount();
76    let unit_end: U128x128 = pair.end.unit_amount().into();
77
78    let positions: Vec<Position> = position_ks
79        .iter()
80        .enumerate()
81        .zip(alphas)
82        .map(|((i, k_i), alpha_i)| {
83            tracing::debug!(i, f64_current_price, k_i, alpha_i, "constructing pool");
84
85            // Case 1: \alpha_i < current_price
86            // Populating ticks that are below the current price, the intuition
87            // is that the positions accumulates the less valuable asset so as
88            // the price trends to \alpha_i, we must provision inventories of
89            // `asset_2`.
90            // \phi(R) = alpha_i * (R_1 = 0)  + 1 * (R_2 = k_i * alpha_i) = k_i * alpha_i
91            // Case 2: \alpha_i >= current_price
92            // Tick is above the current price, therefore we want
93            // to create a one-sided position with price `alpha_i`
94            // that provisions `asset_1`.
95            // \phi(R) = alpha_i * (R_1 = k_i) + 1 * (R_2 = 0) = alpha_i * k_i
96            let approx_p: U128x128 = alpha_i
97                .try_into()
98                .expect("able to convert alpha_i to U128x128");
99            let scaled_p = (approx_p * unit_end).expect("no overflow when scaling p");
100            let p: Amount = scaled_p
101                .round_down()
102                .try_into()
103                .expect("integral after truncating");
104
105            let unscaled_q = Amount::from(1u64);
106            let q = unscaled_q * unit_start;
107
108            if alpha_i < f64_current_price {
109                let r1: Amount = Amount::from(0u64);
110                let approx_r2: U128x128 = (*k_i * pair.end.unit_amount().value() as f64 * alpha_i)
111                    .try_into()
112                    .expect("able to convert k_i * alpha_i to U128x128");
113                let r2: Amount = approx_r2
114                    .round_down()
115                    .try_into()
116                    .expect("integral after truncating");
117
118                tracing::debug!(
119                    i,
120                    k_i,
121                    alpha_i,
122                    f64_current_price,
123                    directed_pair = pair.to_string(),
124                    r1 = field::display(r1),
125                    r2 = field::display(r2),
126                    ?p,
127                    ?q,
128                    "creating position with a tick below the current price"
129                );
130
131                Position::new(
132                    OsRng,
133                    pair.into_directed_trading_pair(),
134                    fee_bps,
135                    p,
136                    q,
137                    Reserves { r1, r2 },
138                )
139            } else {
140                let approx_r1: U128x128 = (*k_i * pair.start.unit_amount().value() as f64)
141                    .try_into()
142                    .expect("able to convert k_i * alpha_i to U128x128");
143                let r1: Amount = approx_r1
144                    .round_down()
145                    .try_into()
146                    .expect("integral after truncating");
147                let r2: Amount = Amount::from(0u64);
148
149                tracing::debug!(
150                    i,
151                    k_i,
152                    alpha_i,
153                    f64_current_price,
154                    directed_pair = pair.to_string(),
155                    %r1,
156                    %r2,
157                    ?p,
158                    ?q,
159                    "creating position with a tick above the current price"
160                );
161
162                Position::new(
163                    OsRng,
164                    pair.into_directed_trading_pair(),
165                    fee_bps,
166                    p,
167                    q,
168                    Reserves { r1, r2 },
169                )
170            }
171        })
172        .collect();
173    Ok(positions)
174}
175
176pub fn solve(
177    alpha: &[f64],
178    k: f64,
179    n: usize,
180) -> anyhow::Result<Array<f64, ndarray::Dim<[usize; 1]>>> {
181    let mut A = Array::zeros((n, n));
182    let mut b = Array::zeros(n);
183
184    for j in 0..n {
185        b[j] = portfolio_value_function(k, alpha[j]);
186
187        for i in 0..j {
188            A[[j, i]] = alpha[i];
189        }
190        for i in j..n {
191            A[[j, i]] = alpha[j];
192        }
193    }
194
195    math_utils::gauss_seidel(
196        A,
197        b,
198        GAUS_SEIDEL_MAX_ITERATION,
199        super::APPROXIMATION_TOLERANCE,
200    )
201}
202
203pub fn portfolio_value_function(invariant_k: f64, price: f64) -> f64 {
204    2.0 * f64::sqrt(invariant_k * price)
205}