pay: Use log probability based bias in channel selection

Changelog-Changed: pay: The route selection will now use the log-propability-based channel selection to increase success rate and reduce time to completion

plugins/libplugin-pay.c

@@ -8,6 +8,7 @@
 #include <common/random_select.h>
 #include <common/type_to_string.h>
 #include <errno.h>
+#include <math.h>
 #include <plugins/libplugin-pay.h>
 #include <sys/types.h>
 #include <wire/peer_wire.h>
@@ -707,22 +708,15 @@ static u64 capacity_bias(const struct gossmap *map,
 			 int dir,
 			 struct amount_msat amount)
 {
-	struct amount_msat fee;
 	struct amount_sat capacity;
-
-	/* Median fees are 1000 base, 10 ppm, so scale capacity bias to that */
-	/* Overflow is pretty-much impossible, so ignore. */
-	if (!amount_msat_fee(&fee, amount, 1000, 10))
-		return 0;
+	u64 capmsat, amtmsat = amount.millisatoshis; /* Raw: lengthy math */
 
 	/* Can fail in theory if gossmap changed underneath. */
 	if (!gossmap_chan_get_capacity(map, c, &capacity))
 		return 0;
 
-	/* bias = fee * (amt / (c + 1)) */
-	return fee.millisatoshis    /* Raw: complex math & laziness */
-		* amount.millisatoshis  /* Raw: complex math & laziness */
-		/ (capacity.satoshis*1000 + 1); /* Raw: complex math & laziness */
+	capmsat = capacity.satoshis * 1000; /* Raw: lengthy math */
+	return -log((capmsat + 1 - amtmsat) / (capmsat + 1));
 }
 
 /* Prioritize costs over distance, but bias to larger channels. */
@@ -732,10 +726,13 @@ static u64 route_score(u32 distance,
 		       int dir,
 		       const struct gossmap_chan *c)
 {
-	u64 costs = cost.millisatoshis + risk.millisatoshis /* Raw: score */
-		/* We use global_gossmap (can't still be NULL)
-		 * *without* get_gossmap() which might change topology. */
-		+ capacity_bias(global_gossmap, c, dir, cost);
+	u64 cmsat = cost.millisatoshis; /* Raw: lengthy math */
+	u64 rmsat = risk.millisatoshis; /* Raw: lengthy math */
+	u64 bias = capacity_bias(global_gossmap, c, dir, cost);
+
+	/* Smoothed harmonic mean to avoid division by 0 */
+	u64 costs = (cmsat * rmsat * bias) / (cmsat + rmsat + bias + 1);
+
 	if (costs > 0xFFFFFFFF)
 		costs = 0xFFFFFFFF;
 	return costs;