FFmpeg
rational.c
Go to the documentation of this file.
1 /*
2  * rational numbers
3  * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at>
4  *
5  * This file is part of FFmpeg.
6  *
7  * FFmpeg is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * FFmpeg is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with FFmpeg; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  */
21 
22 #include "libavutil/rational.c"
23 #include "libswscale/rational64.c"
24 #include "libavutil/integer.h"
25 #include "libavutil/intfloat.h"
26 
27 int main(void)
28 {
29  AVRational a,b,r;
30  AVRational64 a64,b64,r64;
31  int i,j,k;
32  static const int64_t numlist[] = {
33  INT64_MIN, INT64_MIN+1, INT64_MAX, INT32_MIN, INT32_MAX, 1,0,-1,
34  123456789, INT32_MAX-1, INT32_MAX+1LL, UINT32_MAX-1, UINT32_MAX, UINT32_MAX+1LL
35  };
36 
37  for (a.num = -2; a.num <= 2; a.num++) {
38  for (a.den = -2; a.den <= 2; a.den++) {
39  for (b.num = -2; b.num <= 2; b.num++) {
40  for (b.den = -2; b.den <= 2; b.den++) {
41  int c = av_cmp_q(a,b);
42  double d = av_q2d(a) == av_q2d(b) ?
43  0 : (av_q2d(a) - av_q2d(b));
44  if (d > 0) d = 1;
45  else if (d < 0) d = -1;
46  else if (d != d) d = INT_MIN;
47  if (c != d)
48  av_log(NULL, AV_LOG_ERROR, "%d/%d %d/%d, %d %f\n", a.num,
49  a.den, b.num, b.den, c,d);
50  r = av_sub_q(av_add_q(b,a), b);
51  if(b.den && (r.num*a.den != a.num*r.den || !r.num != !a.num || !r.den != !a.den))
52  av_log(NULL, AV_LOG_ERROR, "%d/%d ", r.num, r.den);
53  }
54  }
55  }
56  }
57 
58  for (a64.num = -2; a64.num <= 2; a64.num++) {
59  for (a64.den = -2; a64.den <= 2; a64.den++) {
60  for (b64.num = -2; b64.num <= 2; b64.num++) {
61  for (b64.den = -2; b64.den <= 2; b64.den++) {
62  const double adbl = av_q2d_64(a64);
63  const double bdbl = av_q2d_64(b64);
64  const int c = av_cmp_q64(a64,b64);
65  const int d = adbl == bdbl ? 0 :
66  adbl > bdbl ? 1 :
67  adbl < bdbl ? -1 : INT_MIN;
68 
69  if (c != d)
70  av_log(NULL, AV_LOG_ERROR, "%lld/%lld %lld/%lld, %d != %d\n",
71  (long long) a64.num, (long long) a64.den,
72  (long long) b64.num, (long long) b64.den, c,d);
73 
74  // Check arithmetic result
75  if (a64.den && b64.den) {
76  double rdbl;
77 
78  r64 = av_add_q64(a64, b64);
79  rdbl = av_q2d_64(r64);
80  if (rdbl != adbl + bdbl) {
81  av_log(NULL, AV_LOG_ERROR, "%f + %f = %f != %f\n",
82  adbl, bdbl, rdbl, adbl + bdbl);
83  }
84 
85  r64 = av_mul_q64(a64, b64);
86  rdbl = av_q2d_64(r64);
87  if (rdbl != adbl * bdbl) {
88  av_log(NULL, AV_LOG_ERROR, "%f * %f = %f != %f\n",
89  adbl, bdbl, rdbl, adbl * bdbl);
90  }
91  }
92 
93  // Check addition round-trip
94  r64 = av_sub_q64(av_add_q64(a64, b64), b64);
95  if (b64.den && (r64.num*a64.den != a64.num*r64.den ||
96  !r64.num != !a64.num ||
97  !r64.den != !a64.den))
98  {
99  av_log(NULL, AV_LOG_ERROR, "%lld/%lld != %lld/%lld\n",
100  (long long) a64.num, (long long) a64.den,
101  (long long) r64.num, (long long) r64.den);
102  }
103 
104  if (b64.num) {
105  // Check multiplication round-trip
106  r64 = av_div_q64(av_mul_q64(a64, b64), b64);
107  if (b64.den && (r64.num*a64.den != a64.num*r64.den ||
108  !r64.num != !a64.num ||
109  !r64.den != !a64.den))
110  {
111  av_log(NULL, AV_LOG_ERROR, "%lld/%lld != %lld/%lld\n",
112  (long long) a64.num, (long long) a64.den,
113  (long long) r64.num, (long long) r64.den);
114  }
115  }
116  }
117  }
118  }
119  }
120 
121  /* Check overflow behavior and edge cases */
122  static const AVRational unit_mul_q[][3] = {
123  {{INT_MAX, 2}, { 2, 1}, { INT_MAX, 1}},
124  {{INT_MAX, 2}, {-2, 1}, {-INT_MAX, 1}},
125  {{INT_MAX, 2}, { 0, 1}, {0, 1}},
126  {{INT_MIN, 2}, { 2, 1}, {-INT_MAX, 1}}, /* not INT_MIN */
127  {{INT_MIN, 2}, {-2, 1}, { INT_MAX, 1}},
128  {{INT_MIN, 2}, { 0, 1}, {0, 1}},
129  {{INT_MAX >> 8, 1}, {INT_MAX >> 8, 1}, {INT_MAX, 1}},
130  {{1, INT_MAX >> 8}, {1, INT_MAX >> 8}, {0, 1}},
131  {{1, 1}, {0, 0}, {0, 0}},
132  {{0, 1}, {0, 0}, {0, 0}},
133  };
134 
135  for (i = 0; i < FF_ARRAY_ELEMS(unit_mul_q); i++) {
136  for (int c = 0; c < 2; c++) { /* test commutativity */
137  AVRational a = unit_mul_q[i][c ? 1 : 0];
138  AVRational b = unit_mul_q[i][c ? 0 : 1];
139  AVRational c = unit_mul_q[i][2];
140  AVRational r = av_mul_q(a, b);
141  if (r.num != c.num || r.den != c.den) {
142  av_log(NULL, AV_LOG_ERROR, "%d/%d * %d/%d = %d/%d, expected %d/%d\n",
143  a.num, a.den, b.num, b.den, r.num, r.den, c.num, c.den);
144  }
145  }
146  }
147 
148  static const AVRational unit_add_q[][3] = {
149  {{INT_MAX, 1}, { 2, 2}, { INT_MAX, 1}},
150  {{INT_MAX, 1}, {-2, 2}, { INT_MAX - 1, 1}},
151  {{INT_MAX, 1}, { 0, 2}, { INT_MAX, 1}},
152  {{INT_MIN, 1}, { 2, 2}, {-INT_MAX, 1}},
153  {{INT_MIN, 1}, {-2, 2}, {-INT_MAX, 1}},
154  {{INT_MIN, 1}, { 0, 2}, {-INT_MAX, 1}},
155  {{INT_MAX - 10, 1}, {20, 1}, { INT_MAX, 1}},
156  {{2, INT_MAX}, {2, INT_MAX}, {4, INT_MAX}},
157  {{1, 1}, {0, 0}, {0, 0}},
158  {{0, 1}, {0, 0}, {0, 0}},
159  };
160 
161  for (i = 0; i < FF_ARRAY_ELEMS(unit_add_q); i++) {
162  for (int c = 0; c < 2; c++) { /* test commutativity */
163  AVRational a = unit_add_q[i][c ? 1 : 0];
164  AVRational b = unit_add_q[i][c ? 0 : 1];
165  AVRational c = unit_add_q[i][2];
166  AVRational r = av_add_q(a, b);
167  if (r.num != c.num || r.den != c.den) {
168  av_log(NULL, AV_LOG_ERROR, "%d/%d + %d/%d = %d/%d, expected %d/%d\n",
169  a.num, a.den, b.num, b.den, r.num, r.den, c.num, c.den);
170  }
171  }
172  }
173 
174  static const AVRational64 unit_mul_q64[][3] = {
175  {{INT64_MAX, 2}, { 2, 1}, { INT64_MAX, 1}},
176  {{INT64_MAX, 2}, {-2, 1}, {-INT64_MAX, 1}},
177  {{INT64_MAX, 2}, { 0, 1}, {0, 1}},
178  {{INT64_MIN, 2}, { 2, 1}, {-INT64_MAX, 1}}, /* not INT64_MIN */
179  {{INT64_MIN, 2}, {-2, 1}, { INT64_MAX, 1}},
180  {{INT64_MIN, 2}, { 0, 1}, {0, 1}},
181  {{INT64_MAX >> 8, 1}, {INT64_MAX >> 8, 1}, {INT64_MAX, 1}},
182  {{1, INT64_MAX >> 8}, {1, INT64_MAX >> 8}, {0, 1}},
183  {{1, 1}, {0, 0}, {0, 0}},
184  {{0, 1}, {0, 0}, {0, 0}},
185  };
186 
187  for (i = 0; i < FF_ARRAY_ELEMS(unit_mul_q64); i++) {
188  for (int c = 0; c < 2; c++) { /* test commutativity */
189  AVRational64 a = unit_mul_q64[i][c ? 1 : 0];
190  AVRational64 b = unit_mul_q64[i][c ? 0 : 1];
191  AVRational64 c = unit_mul_q64[i][2];
193  if (r.num != c.num || r.den != c.den) {
194  av_log(NULL, AV_LOG_ERROR, "%lld/%lld * %lld/%lld = %lld/%lld, expected %lld/%lld\n",
195  (long long) a.num, (long long) a.den,
196  (long long) b.num, (long long) b.den,
197  (long long) r.num, (long long) r.den,
198  (long long) c.num, (long long) c.den);
199  }
200  }
201  }
202 
203  static const AVRational64 unit_add_q64[][3] = {
204  {{INT64_MAX, 1}, { 2, 2}, { INT64_MAX, 1}},
205  {{INT64_MAX, 1}, {-2, 2}, { INT64_MAX - 1, 1}},
206  {{INT64_MAX, 1}, { 0, 2}, { INT64_MAX, 1}},
207  {{INT64_MIN, 1}, { 2, 2}, {-INT64_MAX, 1}},
208  {{INT64_MIN, 1}, {-2, 2}, {-INT64_MAX, 1}},
209  {{INT64_MIN, 1}, { 0, 2}, {-INT64_MAX, 1}},
210  {{INT64_MAX - 10, 1}, {20, 1}, { INT64_MAX, 1}},
211  {{2, INT64_MAX}, {2, INT64_MAX}, {4, INT64_MAX}},
212  {{1, 1}, {0, 0}, {0, 0}},
213  {{0, 1}, {0, 0}, {0, 0}},
214  };
215 
216  for (i = 0; i < FF_ARRAY_ELEMS(unit_add_q64); i++) {
217  for (int c = 0; c < 2; c++) { /* test commutativity */
218  AVRational64 a = unit_add_q64[i][c ? 1 : 0];
219  AVRational64 b = unit_add_q64[i][c ? 0 : 1];
220  AVRational64 c = unit_add_q64[i][2];
222  if (r.num != c.num || r.den != c.den) {
223  av_log(NULL, AV_LOG_ERROR, "%lld/%lld + %lld/%lld = %lld/%lld, expected %lld/%lld\n",
224  (long long) a.num, (long long) a.den,
225  (long long) b.num, (long long) b.den,
226  (long long) r.num, (long long) r.den,
227  (long long) c.num, (long long) c.den);
228  }
229  }
230  }
231 
232  for (i = 0; i < FF_ARRAY_ELEMS(numlist); i++) {
233  int64_t a = numlist[i];
234 
235  for (j = 0; j < FF_ARRAY_ELEMS(numlist); j++) {
236  int64_t b = numlist[j];
237  if (b<=0)
238  continue;
239  for (k = 0; k < FF_ARRAY_ELEMS(numlist); k++) {
240  int64_t c = numlist[k];
241  int64_t res;
242  AVInteger ai;
243 
244  if (c<=0)
245  continue;
246  res = av_rescale_rnd(a,b,c, AV_ROUND_ZERO);
247 
248  ai = av_mul_i(av_int2i(a), av_int2i(b));
249  ai = av_div_i(ai, av_int2i(c));
250 
251  if (av_cmp_i(ai, av_int2i(INT64_MAX)) > 0 && res == INT64_MIN)
252  continue;
253  if (av_cmp_i(ai, av_int2i(INT64_MIN)) < 0 && res == INT64_MIN)
254  continue;
255  if (av_cmp_i(ai, av_int2i(res)) == 0)
256  continue;
257 
258  // Special exception for INT64_MIN, remove this in case INT64_MIN is handled without off by 1 error
259  if (av_cmp_i(ai, av_int2i(res-1)) == 0 && a == INT64_MIN)
260  continue;
261 
262  av_log(NULL, AV_LOG_ERROR, "%"PRId64" * %"PRId64" / %"PRId64" = %"PRId64" or %"PRId64"\n", a,b,c, res, av_i2int(ai));
263  }
264  }
265  }
266 
267  for (a.num = 1; a.num <= 10; a.num++) {
268  for (a.den = 1; a.den <= 10; a.den++) {
269  if (av_gcd(a.num, a.den) > 1)
270  continue;
271  for (b.num = 1; b.num <= 10; b.num++) {
272  for (b.den = 1; b.den <= 10; b.den++) {
273  int start;
274  if (av_gcd(b.num, b.den) > 1)
275  continue;
276  if (av_cmp_q(b, a) < 0)
277  continue;
278  for (start = 0; start < 10 ; start++) {
279  int acc= start;
280  int i;
281 
282  for (i = 0; i<100; i++) {
283  int exact = start + av_rescale_q(i+1, b, a);
284  acc = av_add_stable(a, acc, b, 1);
285  if (FFABS(acc - exact) > 2) {
286  av_log(NULL, AV_LOG_ERROR, "%d/%d %d/%d, %d %d\n", a.num,
287  a.den, b.num, b.den, acc, exact);
288  return 1;
289  }
290  }
291  }
292  }
293  }
294  }
295  }
296 
297  for (a.den = 1; a.den < 0x100000000U/3; a.den*=3) {
298  for (a.num = -1; a.num < (1<<27); a.num += 1 + a.num/100) {
299  float f = av_int2float(av_q2intfloat(a));
300  float f2 = av_q2d(a);
301  if (fabs(f - f2) > fabs(f)/5000000) {
302  av_log(NULL, AV_LOG_ERROR, "%d/%d %f %f\n", a.num,
303  a.den, f, f2);
304  return 1;
305  }
306 
307  }
308  }
309 
310  return 0;
311 }
r
const char * r
Definition: vf_curves.c:127
av_add_stable
int64_t av_add_stable(AVRational ts_tb, int64_t ts, AVRational inc_tb, int64_t inc)
Add a value to a timestamp.
Definition: mathematics.c:191
int64_t
long long int64_t
Definition: coverity.c:34
main
int main(void)
Definition: rational.c:27
b
#define b
Definition: input.c:43
AV_ROUND_ZERO
@ AV_ROUND_ZERO
Round toward zero.
Definition: mathematics.h:131
av_sub_q
AVRational av_sub_q(AVRational b, AVRational c)
Subtract one rational from another.
Definition: rational.c:101
av_i2int
int64_t av_i2int(AVInteger a)
Convert the given AVInteger to an int64_t.
Definition: integer.c:160
intfloat.h
av_gcd
int64_t av_gcd(int64_t a, int64_t b)
Compute the greatest common divisor of two integer operands.
Definition: mathematics.c:37
av_q2intfloat
uint32_t av_q2intfloat(AVRational q)
Convert an AVRational to a IEEE 32-bit float expressed in fixed-point format.
Definition: rational.c:154
av_int2float
static av_always_inline float av_int2float(uint32_t i)
Reinterpret a 32-bit integer as a float.
Definition: intfloat.h:40
AV_LOG_ERROR
#define AV_LOG_ERROR
Something went wrong and cannot losslessly be recovered.
Definition: log.h:210
FF_ARRAY_ELEMS
#define FF_ARRAY_ELEMS(a)
Definition: sinewin_tablegen.c:29
av_int2i
AVInteger av_int2i(int64_t a)
Convert the given int64_t to an AVInteger.
Definition: integer.c:149
av_cmp_i
int av_cmp_i(AVInteger a, AVInteger b)
Return 0 if a==b, 1 if a>b and -1 if a<b.
Definition: integer.c:87
av_q2d
static double av_q2d(AVRational a)
Convert an AVRational to a double.
Definition: rational.h:104
av_cmp_q64
int av_cmp_q64(AVRational64 a, AVRational64 b)
Compare two 64-bit rationals.
Definition: rational64.c:108
av_rescale_q
int64_t av_rescale_q(int64_t a, AVRational bq, AVRational cq)
Rescale a 64-bit integer by 2 rational numbers.
Definition: mathematics.c:142
av_mul_q64
AVRational64 av_mul_q64(AVRational64 b, AVRational64 c)
Multiply two 64-bit rationals.
Definition: rational64.c:124
FFABS
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition: common.h:74
rational.c
fabs
static __device__ float fabs(float a)
Definition: cuda_runtime.h:182
NULL
#define NULL
Definition: coverity.c:32
AVRational
Rational number (pair of numerator and denominator).
Definition: rational.h:58
c
Undefined Behavior In the C some operations are like signed integer dereferencing freed accessing outside allocated Undefined Behavior must not occur in a C it is not safe even if the output of undefined operations is unused The unsafety may seem nit picking but Optimizing compilers have in fact optimized code on the assumption that no undefined Behavior occurs Optimizing code based on wrong assumptions can and has in some cases lead to effects beyond the output of computations The signed integer overflow problem in speed critical code Code which is highly optimized and works with signed integers sometimes has the problem that often the output of the computation does not c
Definition: undefined.txt:32
av_div_q64
AVRational64 av_div_q64(AVRational64 b, AVRational64 c)
Divide one 64-bit rational by another.
Definition: rational64.c:130
av_rescale_rnd
int64_t av_rescale_rnd(int64_t a, int64_t b, int64_t c, enum AVRounding rnd)
Rescale a 64-bit integer with specified rounding.
Definition: mathematics.c:58
f
f
Definition: af_crystalizer.c:122
av_sub_q64
AVRational64 av_sub_q64(AVRational64 b, AVRational64 c)
Subtract one 64-bit rational from another.
Definition: rational64.c:141
i
#define i(width, name, range_min, range_max)
Definition: cbs_h264.c:63
av_mul_i
AVInteger av_mul_i(AVInteger a, AVInteger b)
Definition: integer.c:66
AVRational64
64-bit Rational number (pair of numerator and denominator).
Definition: rational64.h:52
a
The reader does not expect b to be semantically here and if the code is changed by maybe adding a a division or other the signedness will almost certainly be mistaken To avoid this confusion a new type was SUINT is the C unsigned type but it holds a signed int to use the same example SUINT a
Definition: undefined.txt:41
AVInteger
Definition: integer.h:36
av_cmp_q
static int av_cmp_q(AVRational a, AVRational b)
Compare two rationals.
Definition: rational.h:89
AVRational64::den
int64_t den
Denominator.
Definition: rational64.h:54
av_div_i
AVInteger av_div_i(AVInteger a, AVInteger b)
Return a/b.
Definition: integer.c:143
av_q2d_64
static double av_q2d_64(AVRational64 a)
Convert an AVRational64 to a double.
Definition: rational64.h:90
rational64.c
av_mul_q
AVRational av_mul_q(AVRational b, AVRational c)
Multiply two rationals.
Definition: rational.c:80
av_add_q
AVRational av_add_q(AVRational b, AVRational c)
Add two rationals.
Definition: rational.c:93
AVRational64::num
int64_t num
Numerator.
Definition: rational64.h:53
av_log
#define av_log(a,...)
Definition: tableprint_vlc.h:27
integer.h
av_add_q64
AVRational64 av_add_q64(AVRational64 b, AVRational64 c)
Add two 64-bit rationals.
Definition: rational64.c:135