FFmpeg
Loading...
Searching...
No Matches
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"
24#include "libavutil/integer.h"
25#include "libavutil/intfloat.h"
26
27int main(void)
28{
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;
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}
#define i(width, name, range_min, range_max)
Definition cbs_h264.c:63
#define f(width, name)
Definition cbs_vp8.c:236
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition common.h:74
#define NULL
Definition coverity.c:32
long long int64_t
Definition coverity.c:34
static __device__ float fabs(float a)
AVRational64 av_div_q64(AVRational64 b, AVRational64 c)
Divide one 64-bit rational by another.
Definition rational64.c:130
AVRational64 av_mul_q64(AVRational64 b, AVRational64 c)
Multiply two 64-bit rationals.
Definition rational64.c:124
int av_cmp_q64(AVRational64 a, AVRational64 b)
Compare two 64-bit rationals.
Definition rational64.c:108
AVRational64 av_add_q64(AVRational64 b, AVRational64 c)
Add two 64-bit rationals.
Definition rational64.c:135
static double av_q2d_64(AVRational64 a)
Convert an AVRational64 to a double.
Definition rational64.h:90
AVRational64 av_sub_q64(AVRational64 b, AVRational64 c)
Subtract one 64-bit rational from another.
Definition rational64.c:141
#define AV_LOG_ERROR
Something went wrong and cannot losslessly be recovered.
Definition log.h:210
AVRational av_add_q(AVRational b, AVRational c)
Add two rationals.
Definition rational.c:93
AVRational av_mul_q(AVRational b, AVRational c)
Multiply two rationals.
Definition rational.c:80
static double av_q2d(AVRational a)
Convert an AVRational to a double.
Definition rational.h:104
static int av_cmp_q(AVRational a, AVRational b)
Compare two rationals.
Definition rational.h:89
AVRational av_sub_q(AVRational b, AVRational c)
Subtract one rational from another.
Definition rational.c:101
uint32_t av_q2intfloat(AVRational q)
Convert an AVRational to a IEEE 32-bit float expressed in fixed-point format.
Definition rational.c:154
int64_t av_add_stable(AVRational ts_tb, int64_t ts, AVRational inc_tb, int64_t inc)
Add a value to a timestamp.
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
int64_t av_gcd(int64_t a, int64_t b)
Compute the greatest common divisor of two integer operands.
Definition mathematics.c:37
int64_t av_rescale_q(int64_t a, AVRational bq, AVRational cq)
Rescale a 64-bit integer by 2 rational numbers.
@ AV_ROUND_ZERO
Round toward zero.
int a
#define r
Definition input.c:42
#define b
Definition input.c:43
AVInteger av_div_i(AVInteger a, AVInteger b)
Return a/b.
Definition integer.c:143
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
AVInteger av_int2i(int64_t a)
Convert the given int64_t to an AVInteger.
Definition integer.c:149
AVInteger av_mul_i(AVInteger a, AVInteger b)
Definition integer.c:66
int64_t av_i2int(AVInteger a)
Convert the given AVInteger to an int64_t.
Definition integer.c:160
arbitrary precision integers
static av_always_inline float av_int2float(uint32_t i)
Reinterpret a 32-bit integer as a float.
Definition intfloat.h:40
64-bit rational numbers
rational numbers
#define FF_ARRAY_ELEMS(a)
64-bit Rational number (pair of numerator and denominator).
Definition rational64.h:52
int64_t num
Numerator.
Definition rational64.h:53
int64_t den
Denominator.
Definition rational64.h:54
Rational number (pair of numerator and denominator).
Definition rational.h:58
#define av_log(a,...)
int main(void)
Definition rational.c:27
static double c[64]