OR-Tools  8.2
rational_approximation.cc
Go to the documentation of this file.
1// Copyright 2010-2018 Google LLC
2// Licensed under the Apache License, Version 2.0 (the "License");
3// you may not use this file except in compliance with the License.
4// You may obtain a copy of the License at
5//
6// http://www.apache.org/licenses/LICENSE-2.0
7//
8// Unless required by applicable law or agreed to in writing, software
9// distributed under the License is distributed on an "AS IS" BASIS,
10// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11// See the License for the specific language governing permissions and
12// limitations under the License.
13
15
16#include <cmath>
17#include <limits>
18
20
21namespace operations_research {
22
23// Computes a rational approximation numerator/denominator for value x
24// using a continued fraction algorithm. The absolute difference between the
25// output fraction and the input "x" will not exceed "precision".
26Fraction RationalApproximation(const double x, const double precision) {
27 DCHECK_LT(x, std::numeric_limits<double>::infinity());
28 DCHECK_GT(x, -std::numeric_limits<double>::infinity());
29 // All computations are made on long doubles to guarantee the maximum
30 // precision for the approximations. This way, the approximations when
31 // Fractional is float or double are as accurate as possible.
32 long double abs_x = std::abs(x);
33 long double y = abs_x;
34 int64 previous_numerator = 0;
35 int64 previous_denominator = 1;
36 int64 numerator = 1;
37 int64 denominator = 0;
38 while (true) {
39 const int64 term = static_cast<int64>(std::floor(y));
40 const int64 new_numerator = term * numerator + previous_numerator;
41 const int64 new_denominator = term * denominator + previous_denominator;
42 // If there was an overflow, we prefer returning a not-so-good approximation
43 // rather than something that is completely wrong.
44 if (new_numerator < 0 || new_denominator < 0) break;
45 previous_numerator = numerator;
46 previous_denominator = denominator;
47 numerator = new_numerator;
48 denominator = new_denominator;
49 long double numerator_approximation = abs_x * denominator;
50 if (std::abs(numerator_approximation - numerator) <=
51 precision * numerator_approximation) {
52 break;
53 }
54 y = 1 / (y - term);
55 }
56 return Fraction((x < 0) ? -numerator : numerator, denominator);
57}
58} // namespace operations_research
#define DCHECK_GT(val1, val2)
Definition: base/logging.h:890
#define DCHECK_LT(val1, val2)
Definition: base/logging.h:888
int64_t int64
The vehicle routing library lets one model and solve generic vehicle routing problems ranging from th...
std::pair< int64, int64 > Fraction
Fraction RationalApproximation(const double x, const double precision)