Add an option to make the success probability estimation nonlinear

Our "what is the success probability of paying over a channel with
the given liquidity bounds" calculation currently assumes the
probability of where the liquidity lies in a channel is constant
across the entire capacity of a channel. This is obviously a
somewhat dubious assumption given most nodes don't materially
rebalance and flows within the network often push liquidity
"towards the edges".

Here we add an option to consider this when scoring channels during
routefinding. Specifically, if a new `linear_success_probability`
flag is set on `ProbabilisticScoringFeeParameters`, rather than
assuming a PDF of `1` (across the channel's capacity scaled from 0
to 1), we use `(x - 0.5)^6`.

This assumes liquidity is very likely to be near the edges, which
matches experimental results. Further, calculating the CDF (i.e.
integral) between arbitrary points on the PDF is trivial, which we
do as our main scoring function.

While this (finally) introduces floats in our scoring, its not
practical to exponentiate using fixed-precision, and benchmarks
show this is a performance regression, but not a huge one, more
than made up for by the increase in payment success rates.

Author: TheBlueMatt

bench/benches/bench.rs

@@ -13,6 +13,9 @@ criterion_group!(benches,
 	lightning::routing::router::benches::generate_routes_with_probabilistic_scorer,
 	lightning::routing::router::benches::generate_mpp_routes_with_probabilistic_scorer,
 	lightning::routing::router::benches::generate_large_mpp_routes_with_probabilistic_scorer,
+	lightning::routing::router::benches::generate_routes_with_nonlinear_probabilistic_scorer,
+	lightning::routing::router::benches::generate_mpp_routes_with_nonlinear_probabilistic_scorer,
+	lightning::routing::router::benches::generate_large_mpp_routes_with_nonlinear_probabilistic_scorer,
 	lightning::sign::benches::bench_get_secure_random_bytes,
 	lightning::ln::channelmanager::bench::bench_sends,
 	lightning_persister::bench::bench_sends,

lightning/src/routing/router.rs

@@ -7114,6 +7114,38 @@ pub mod benches {
 			"generate_large_mpp_routes_with_probabilistic_scorer");
 	}
 
+	pub fn generate_routes_with_nonlinear_probabilistic_scorer(bench: &mut Criterion) {
+		let logger = TestLogger::new();
+		let network_graph = bench_utils::read_network_graph(&logger).unwrap();
+		let mut params = ProbabilisticScoringParameters::default();
+		params.linear_success_probability = false;
+		let scorer = ProbabilisticScorer::new(params, &network_graph, &logger);
+		generate_routes(bench, &network_graph, scorer, InvoiceFeatures::empty(), 0,
+			"generate_routes_with_nonlinear_probabilistic_scorer");
+	}
+
+	pub fn generate_mpp_routes_with_nonlinear_probabilistic_scorer(bench: &mut Criterion) {
+		let logger = TestLogger::new();
+		let network_graph = bench_utils::read_network_graph(&logger).unwrap();
+		let mut params = ProbabilisticScoringParameters::default();
+		params.linear_success_probability = false;
+		let scorer = ProbabilisticScorer::new(params, &network_graph, &logger);
+		generate_routes(bench, &network_graph, scorer,
+			channelmanager::provided_invoice_features(&UserConfig::default()), 0,
+			"generate_mpp_routes_with_nonlinear_probabilistic_scorer");
+	}
+
+	pub fn generate_large_mpp_routes_with_nonlinear_probabilistic_scorer(bench: &mut Criterion) {
+		let logger = TestLogger::new();
+		let network_graph = bench_utils::read_network_graph(&logger).unwrap();
+		let mut params = ProbabilisticScoringParameters::default();
+		params.linear_success_probability = false;
+		let scorer = ProbabilisticScorer::new(params, &network_graph, &logger);
+		generate_routes(bench, &network_graph, scorer,
+			channelmanager::provided_invoice_features(&UserConfig::default()), 100_000_000,
+			"generate_large_mpp_routes_with_nonlinear_probabilistic_scorer");
+	}
+
 	fn generate_routes<S: Score>(
 		bench: &mut Criterion, graph: &NetworkGraph<&TestLogger>, mut scorer: S,
 		score_params: &S::ScoreParams, features: Bolt11InvoiceFeatures, starting_amount: u64,

lightning/src/routing/scoring.rs

@@ -509,6 +509,28 @@ pub struct ProbabilisticScoringFeeParameters {
 	/// [`base_penalty_msat`]: Self::base_penalty_msat
 	/// [`anti_probing_penalty_msat`]: Self::anti_probing_penalty_msat
 	pub considered_impossible_penalty_msat: u64,
+
+	/// In order to calculate most of the scores above, we must first convert a lower and upper
+	/// bound on the available liquidity in a channel into the probability that we think a payment
+	/// will succeed. That probability is derived from a Probability Density Function for where we
+	/// think the liquidity in a channel likely lies, given such bounds.
+	///
+	/// If this flag is set, that PDF is simply a constant - we assume that the actual available
+	/// liquidity in a channel is just as likely to be at any point between our lower and upper
+	/// bounds.
+	///
+	/// If this flag is *not* set, that PDF is `(x - 0.5*capacity) ^ 6`. That is, we use an
+	/// exponential curve which expects the liquidity of a channel to lie "at the edges". This
+	/// matches experimental results - most routing nodes do not aggressively rebalance their
+	/// channels and flows in the network are often unbalanced, leaving liquidity usually
+	/// unavailable.
+	///
+	/// Thus, for the "best" routes, leave this flag `false`. However, the flag does imply a number
+	/// of floating-point multiplications in the hottest routing code, which may lead to routing
+	/// performance degradation on some machines.
+	///
+	/// Default value: false
+	pub linear_success_probability: bool,
 }
 
 impl Default for ProbabilisticScoringFeeParameters {
@@ -523,6 +545,7 @@ impl Default for ProbabilisticScoringFeeParameters {
 			considered_impossible_penalty_msat: 1_0000_0000_000,
 			historical_liquidity_penalty_multiplier_msat: 10_000,
 			historical_liquidity_penalty_amount_multiplier_msat: 64,
+			linear_success_probability: false,
 		}
 	}
 }
@@ -576,6 +599,7 @@ impl ProbabilisticScoringFeeParameters {
 			manual_node_penalties: HashMap::new(),
 			anti_probing_penalty_msat: 0,
 			considered_impossible_penalty_msat: 0,
+			linear_success_probability: true,
 		}
 	}
 }
@@ -947,14 +971,42 @@ const BASE_AMOUNT_PENALTY_DIVISOR: u64 = 1 << 30;
 /// seen an HTLC successfully complete over this channel.
 fn success_probability(
 	min_liquidity_msat: u64, amount_msat: u64, max_liquidity_msat: u64, capacity_msat: u64,
-	_params: &ProbabilisticScoringFeeParameters, min_zero_implies_no_successes: bool,
+	params: &ProbabilisticScoringFeeParameters, min_zero_implies_no_successes: bool,
 ) -> (u64, u64) {
 	debug_assert!(min_liquidity_msat <= amount_msat);
 	debug_assert!(amount_msat < max_liquidity_msat);
 	debug_assert!(max_liquidity_msat <= capacity_msat);
 
-	let numerator = max_liquidity_msat - amount_msat;
-	let mut denominator = (max_liquidity_msat - min_liquidity_msat).saturating_add(1);
+	let (numerator, mut denominator) =
+		if params.linear_success_probability {
+			(max_liquidity_msat - amount_msat,
+				(max_liquidity_msat - min_liquidity_msat).saturating_add(1))
+		} else {
+			let capacity = capacity_msat as f64;
+			let min = (min_liquidity_msat as f64) / capacity;
+			let max = (max_liquidity_msat as f64) / capacity;
+			let amount = (amount_msat as f64) / capacity;
+
+			// Assume the channel has a probability density function of (x - 0.5)^6 for values from 0 to 1
+			// (where 1 is the channel's full capacity). The success probability given some liquidity
+			// bounds is thus the integral under the curve from the minimum liquidity to the amount,
+			// divided by the same integral from the minimum to the maximum liquidity bounds.
+			//
+			// For (x - 0.5)^7, this means simply subtracting the two bounds mius 0.5 to the 7th power.
+			let max_pow = (max - 0.5).powi(7);
+			let num = max_pow - (amount - 0.5).powi(7);
+			let den = max_pow - (min - 0.5).powi(7);
+
+			// Because our numerator and denominator max out at 2^-6 we need to multiply them by
+			// quite a large factor to get something useful (ideally in the 2^30 range).
+			const ALMOST_TRILLION: f64 = 1024.0 * 1024.0 * 1024.0 * 64.0;
+			let numerator = (num * ALMOST_TRILLION) as u64 + 1;
+			let denominator = (den * ALMOST_TRILLION) as u64 + 1;
+			debug_assert!(numerator <= 1 << 31, "Got large numerator ({}) from float {}.", numerator, num);
+			debug_assert!(denominator <= 1 << 31, "Got large denominator ({}) from float {}.", denominator, den);
+			(numerator, denominator)
+		};
+
 	if min_zero_implies_no_successes && min_liquidity_msat == 0 &&
 		numerator < u64::max_value() / 3 {
 		// If we have no knowledge of the channel, scale probability down by 75%