Skip to content

Commit 25ab7c5

Browse files
author
Tomáš Kolárik
committed
Implemented simple_continued_fraction(std::vector<Z>)
1 parent d83bd53 commit 25ab7c5

File tree

1 file changed

+39
-18
lines changed

1 file changed

+39
-18
lines changed

include/boost/math/tools/simple_continued_fraction.hpp

+39-18
Original file line numberDiff line numberDiff line change
@@ -37,21 +37,34 @@ class simple_continued_fraction {
3737
public:
3838
typedef Z int_type;
3939

40-
simple_continued_fraction(Real x) {
40+
simple_continued_fraction(std::vector<Z> data) : b_{std::move(data)} {
41+
const size_t size_ = b_.size();
42+
if (size_ == 0) {
43+
throw std::length_error("Array of coefficients is empty.");
44+
}
45+
46+
for (size_t i = 1; i < size_; ++i) {
47+
if (b_[i] <= 0) {
48+
std::ostringstream oss;
49+
oss << "Found a negative partial denominator: b[" << i << "] = " << b_[i] << ".";
50+
throw std::domain_error(oss.str());
51+
}
52+
}
53+
54+
canonicalize();
55+
}
56+
57+
simple_continued_fraction(Real x) : b_{} {
4158
using std::floor;
4259
using std::abs;
4360
using std::sqrt;
4461
using std::isfinite;
45-
const Real orig_x = x;
4662
if (!isfinite(x)) {
4763
throw std::domain_error("Cannot convert non-finites into continued fractions.");
4864
}
4965

50-
constexpr int p = std::numeric_limits<Real>::max_digits10;
51-
if constexpr (p == 2147483647) {
52-
precision_ = orig_x.backend().precision();
53-
} else {
54-
precision_ = p;
66+
if constexpr (std_precision == 2147483647) {
67+
precision_ = x.backend().precision();
5568
}
5669

5770
b_.reserve(50);
@@ -61,6 +74,8 @@ class simple_continued_fraction {
6174
b_.shrink_to_fit();
6275
return;
6376
}
77+
78+
const Real orig_x = x;
6479
x = 1/(x-bj);
6580
Real f = bj;
6681
if (bj == 0) {
@@ -87,14 +102,7 @@ class simple_continued_fraction {
87102
D = 1/D;
88103
f *= (C*D);
89104
}
90-
// Deal with non-uniqueness of continued fractions: [a0; a1, ..., an, 1] = a0; a1, ..., an + 1].
91-
// The shorter representation is considered the canonical representation,
92-
// so if we compute a non-canonical representation, change it to canonical:
93-
if (b_.size() > 2 && b_.back() == 1) {
94-
b_.pop_back();
95-
b_.back() += 1;
96-
}
97-
b_.shrink_to_fit();
105+
canonicalize();
98106

99107
const size_t size_ = b_.size();
100108
for (size_t i = 1; i < size_; ++i) {
@@ -163,9 +171,22 @@ class simple_continued_fraction {
163171
template<typename T, typename Z2>
164172
friend std::ostream& operator<<(std::ostream& out, simple_continued_fraction<T, Z2>& scf);
165173
private:
166-
std::vector<Z> b_{};
174+
static constexpr int std_precision = std::numeric_limits<Real>::max_digits10;
175+
176+
void canonicalize() {
177+
// Deal with non-uniqueness of continued fractions: [a0; a1, ..., an, 1] = a0; a1, ..., an + 1].
178+
// The shorter representation is considered the canonical representation,
179+
// so if we compute a non-canonical representation, change it to canonical:
180+
if (b_.size() > 2 && b_.back() == 1) {
181+
b_.pop_back();
182+
b_.back() += 1;
183+
}
184+
b_.shrink_to_fit();
185+
}
186+
187+
std::vector<Z> b_;
167188

168-
int precision_{};
189+
int precision_{std_precision};
169190
};
170191

171192

@@ -189,7 +210,7 @@ std::ostream& operator<<(std::ostream& out, simple_continued_fraction<Real, Z2>&
189210
template<typename Real, typename Z = std::int64_t>
190211
inline auto simple_continued_fraction_coefficients(Real x)
191212
{
192-
auto temp = simple_continued_fraction(x);
213+
auto temp = simple_continued_fraction<Real, Z>(x);
193214
return temp.get_data();
194215
}
195216

0 commit comments

Comments
 (0)