Line data Source code
1 : /*
2 : This file is part of TALER
3 : Copyright (C) 2023 Taler Systems SA
4 :
5 : TALER is free software; you can redistribute it and/or modify it under the
6 : terms of the GNU General Public License as published by the Free Software
7 : Foundation; either version 3, or (at your option) any later version.
8 :
9 : TALER is distributed in the hope that it will be useful, but WITHOUT ANY
10 : WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11 : A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12 :
13 : You should have received a copy of the GNU General Public License along with
14 : TALER; see the file COPYING. If not, see
15 : <http://www.gnu.org/licenses/>
16 : */
17 : /**
18 : * @file lib/exchange_api_stefan.c
19 : * @brief calculations on the STEFAN curve
20 : * @author Christian Grothoff
21 : */
22 : #include "taler/taler_json_lib.h"
23 : #include <gnunet/gnunet_curl_lib.h>
24 : #include "exchange_api_handle.h"
25 : #include <math.h>
26 :
27 :
28 : /**
29 : * Determine smallest denomination in @a keys.
30 : *
31 : * @param keys exchange response to evaluate
32 : * @return NULL on error (no denominations)
33 : */
34 : static const struct TALER_Amount *
35 54460 : get_unit (const struct TALER_EXCHANGE_Keys *keys)
36 : {
37 54460 : const struct TALER_Amount *min = NULL;
38 :
39 108920 : for (unsigned int i = 0; i<keys->num_denom_keys; i++)
40 : {
41 54460 : const struct TALER_EXCHANGE_DenomPublicKey *dk
42 54460 : = &keys->denom_keys[i];
43 :
44 54460 : if ( (NULL == min) ||
45 0 : (1 == TALER_amount_cmp (min,
46 : /* > */
47 : &dk->value)) )
48 54460 : min = &dk->value;
49 : }
50 54460 : GNUNET_break (NULL != min);
51 54460 : return min;
52 : }
53 :
54 :
55 : /**
56 : * Convert amount to double for STEFAN curve evaluation.
57 : *
58 : * @param a input amount
59 : * @return (rounded) amount as a double
60 : */
61 : static double
62 180492 : amount_to_double (const struct TALER_Amount *a)
63 : {
64 180492 : double d = (double) a->value;
65 :
66 180492 : d += a->fraction / ((double) TALER_AMOUNT_FRAC_BASE);
67 180492 : return d;
68 : }
69 :
70 :
71 : /**
72 : * Convert double to amount for STEFAN curve evaluation.
73 : *
74 : * @param dv input amount
75 : * @param currency deisred currency
76 : * @param[out] rval (rounded) amount as a double
77 : */
78 : static void
79 33072 : double_to_amount (double dv,
80 : const char *currency,
81 : struct TALER_Amount *rval)
82 : {
83 : double rem;
84 :
85 33072 : GNUNET_assert (GNUNET_OK ==
86 : TALER_amount_set_zero (currency,
87 : rval));
88 33072 : if (! (dv > 0.0))
89 : {
90 : /* Negative or NaN: converting this to the unsigned @e value
91 : would be undefined behaviour; the zero amount is the best
92 : (and safe) answer here. */
93 0 : GNUNET_break (0);
94 0 : return;
95 : }
96 33072 : if (dv >= (double) TALER_AMOUNT_MAX_VALUE)
97 : {
98 : /* Out of range for a Taler amount; saturate instead of
99 : invoking undefined behaviour in the conversion below. */
100 0 : GNUNET_break (0);
101 0 : rval->value = TALER_AMOUNT_MAX_VALUE - 1;
102 0 : rval->fraction = TALER_AMOUNT_FRAC_BASE - 1;
103 0 : return;
104 : }
105 33072 : rval->value = floorl (dv);
106 33072 : rem = dv - ((double) rval->value);
107 33072 : if (rem < 0.0)
108 0 : rem = 0.0;
109 33072 : rem *= TALER_AMOUNT_FRAC_BASE;
110 33072 : rval->fraction = floorl (rem);
111 33072 : if (rval->fraction >= TALER_AMOUNT_FRAC_BASE)
112 : {
113 : /* Strange, multiplication overflowed our range,
114 : round up value instead */
115 0 : rval->fraction = 0;
116 0 : rval->value += 1;
117 : }
118 : }
119 :
120 :
121 : enum GNUNET_GenericReturnValue
122 23762 : TALER_EXCHANGE_keys_stefan_b2n (
123 : const struct TALER_EXCHANGE_Keys *keys,
124 : const struct TALER_Amount *brut,
125 : struct TALER_Amount *net)
126 : {
127 : const struct TALER_Amount *min;
128 23762 : double log_d = amount_to_double (&keys->stefan_log);
129 23762 : double lin_d = keys->stefan_lin;
130 23762 : double abs_d = amount_to_double (&keys->stefan_abs);
131 23762 : double bru_d = amount_to_double (brut);
132 : double min_d;
133 : double fee_d;
134 : double net_d;
135 :
136 23762 : if (TALER_amount_is_zero (brut))
137 : {
138 2376 : *net = *brut;
139 2376 : return GNUNET_NO;
140 : }
141 21386 : min = get_unit (keys);
142 21386 : if (NULL == min)
143 0 : return GNUNET_SYSERR;
144 21386 : if (1.0 <= keys->stefan_lin)
145 : {
146 : /* This cannot work, linear STEFAN fee estimate always
147 : exceed any gross amount. */
148 2 : GNUNET_break_op (0);
149 2 : return GNUNET_SYSERR;
150 : }
151 21384 : min_d = amount_to_double (min);
152 21384 : if (min_d <= 0.0)
153 : {
154 : /* smallest denomination must have a positive value,
155 : otherwise the logarithm below is undefined */
156 0 : GNUNET_break_op (0);
157 0 : return GNUNET_SYSERR;
158 : }
159 21384 : fee_d = abs_d
160 21384 : + log_d * log2 (bru_d / min_d)
161 21384 : + lin_d * bru_d;
162 21384 : if (fee_d < 0.0)
163 : {
164 : /* log2() is negative for gross amounts below the smallest
165 : denomination; a fee is never negative, and without this
166 : clamp the net amount would exceed the gross amount */
167 0 : fee_d = 0.0;
168 : }
169 21384 : if (fee_d > bru_d)
170 : {
171 4848 : GNUNET_assert (GNUNET_OK ==
172 : TALER_amount_set_zero (brut->currency,
173 : net));
174 4848 : return GNUNET_NO;
175 : }
176 16536 : net_d = bru_d - fee_d;
177 16536 : double_to_amount (net_d,
178 16536 : brut->currency,
179 : net);
180 16536 : return GNUNET_OK;
181 : }
182 :
183 :
184 : /**
185 : * Our function
186 : * f(x) := ne + ab + lo * log2(x/mi) + li * x - x
187 : * for #newton().
188 : */
189 : static double
190 50076 : eval_f (double mi,
191 : double ab,
192 : double lo,
193 : double li,
194 : double ne,
195 : double x)
196 : {
197 50076 : return ne + ab + lo * log2 (x / mi) + li * x - x;
198 : }
199 :
200 :
201 : /**
202 : * Our function
203 : * f'(x) := lo / log(2) / x + li - 1
204 : * for #newton().
205 : */
206 : static double
207 50076 : eval_fp (double mi,
208 : double lo,
209 : double li,
210 : double ne,
211 : double x)
212 : {
213 50076 : return lo / log (2) / x + li - 1;
214 : }
215 :
216 :
217 : /**
218 : * Use Newton's method to find x where f(x)=0.
219 : *
220 : * @return x where "eval_f(x)==0".
221 : */
222 : static double
223 16536 : newton (double mi,
224 : double ab,
225 : double lo,
226 : double li,
227 : double ne)
228 : {
229 16536 : const double eps = 0.00000001; /* max error allowed */
230 16536 : double min_ab = ne + ab; /* result cannot be smaller than this! */
231 : /* compute lower bounds by various heuristics */
232 16536 : double min_ab_li = min_ab + min_ab * li;
233 16536 : double min_ab_li_lo = min_ab_li + log2 (min_ab_li / mi) * lo;
234 16536 : double min_ab_lo = min_ab + log2 (min_ab / mi) * lo;
235 16536 : double min_ab_lo_li = min_ab_lo + min_ab_lo * li;
236 : /* take global lower bound */
237 16536 : double x_min = GNUNET_MAX (min_ab_lo_li,
238 : min_ab_li_lo);
239 16536 : double x = x_min; /* use lower bound as starting point */
240 :
241 : /* Objective: invert
242 : ne := br - ab - lo * log2 (br/mi) - li * br
243 : to find 'br'.
244 : Method: use Newton's method to find root of:
245 : f(x) := ne + ab + lo * log2 (x/mi) + li * x - x
246 : using also
247 : f'(x) := lo / log(2) / x + li - 1
248 : */
249 : /* Loop to abort in case of divergence;
250 : 100 is already very high, 2-4 is normal! */
251 50076 : for (unsigned int i = 0; i<100; i++)
252 : {
253 50076 : double fx = eval_f (mi, ab, lo, li, ne, x);
254 50076 : double fxp = eval_fp (mi, lo, li, ne, x);
255 50076 : double x_new = x - fx / fxp;
256 :
257 50076 : if (fabs (x - x_new) <= eps)
258 : {
259 16536 : GNUNET_log (GNUNET_ERROR_TYPE_INFO,
260 : "Needed %u rounds from %f to result BRUT %f => NET: %f\n",
261 : i,
262 : x_min,
263 : x_new,
264 : x_new - ab - li * x_new - lo * log2 (x / mi));
265 16536 : return x_new;
266 : }
267 33540 : if (x_new < x_min)
268 : {
269 0 : GNUNET_log (GNUNET_ERROR_TYPE_WARNING,
270 : "Divergence, obtained very bad estimate %f after %u rounds!\n",
271 : x_new,
272 : i);
273 0 : return x_min;
274 : }
275 33540 : x = x_new;
276 : }
277 0 : GNUNET_log (GNUNET_ERROR_TYPE_WARNING,
278 : "Slow convergence, returning bad estimate %f!\n",
279 : x);
280 0 : return x;
281 : }
282 :
283 :
284 : enum GNUNET_GenericReturnValue
285 23762 : TALER_EXCHANGE_keys_stefan_n2b (
286 : const struct TALER_EXCHANGE_Keys *keys,
287 : const struct TALER_Amount *net,
288 : struct TALER_Amount *brut)
289 : {
290 : const struct TALER_Amount *min;
291 23762 : double lin_d = keys->stefan_lin;
292 23762 : double log_d = amount_to_double (&keys->stefan_log);
293 23762 : double abs_d = amount_to_double (&keys->stefan_abs);
294 23762 : double net_d = amount_to_double (net);
295 : double min_d;
296 : double brut_d;
297 :
298 23762 : if (TALER_amount_is_zero (net))
299 : {
300 7224 : *brut = *net;
301 7224 : return GNUNET_NO;
302 : }
303 16538 : min = get_unit (keys);
304 16538 : if (NULL == min)
305 0 : return GNUNET_SYSERR;
306 16538 : if (1.0 <= keys->stefan_lin)
307 : {
308 : /* This cannot work, linear STEFAN fee estimate always
309 : exceed any gross amount. */
310 2 : GNUNET_break_op (0);
311 2 : return GNUNET_SYSERR;
312 : }
313 16536 : min_d = amount_to_double (min);
314 16536 : brut_d = newton (min_d,
315 : abs_d,
316 : log_d,
317 : lin_d,
318 : net_d);
319 16536 : double_to_amount (brut_d,
320 16536 : net->currency,
321 : brut);
322 16536 : return GNUNET_OK;
323 : }
324 :
325 :
326 : void
327 17724 : TALER_EXCHANGE_keys_stefan_round (
328 : const struct TALER_EXCHANGE_Keys *keys,
329 : struct TALER_Amount *val)
330 : {
331 : const struct TALER_Amount *min;
332 : uint32_t mod;
333 : uint32_t frac;
334 : uint32_t lim;
335 :
336 17724 : if (0 == val->fraction)
337 : {
338 : /* rounding of non-fractions not supported */
339 1188 : return;
340 : }
341 16536 : min = get_unit (keys);
342 16536 : if (NULL == min)
343 0 : return;
344 16536 : if (0 == min->fraction)
345 : {
346 0 : frac = TALER_AMOUNT_FRAC_BASE;
347 : }
348 : else
349 : {
350 16536 : frac = min->fraction;
351 : }
352 16536 : lim = frac / 2;
353 16536 : mod = val->fraction % frac;
354 16536 : if (mod < lim)
355 0 : val->fraction -= mod; /* round down */
356 : else
357 16536 : val->fraction += frac - mod; /* round up */
358 16536 : GNUNET_assert (GNUNET_SYSERR !=
359 : TALER_amount_normalize (val));
360 : }
|