ARTS 2.5.10 (git: 2f1c442c)
rational.h
Go to the documentation of this file.
3
4This program is free software; you can redistribute it and/or modify it
6Free Software Foundation; either version 2, or (at your option) any
7later version.
8
9This program is distributed in the hope that it will be useful,
10but WITHOUT ANY WARRANTY; without even the implied warranty of
11MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12GNU General Public License for more details.
13
14You should have received a copy of the GNU General Public License
15along with this program; if not, write to the Free Software
16Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
17USA. */
18
27#ifndef rational_h
28#define rational_h
29
30#include "array.h"
31#include "bifstream.h"
32#include "bofstream.h"
33#include "math_funcs.h"
34#include "matpack.h"
35
36#include <numeric>
37#include <ostream>
38
39
40using std::gcd;
41
43struct Rational {
46
52 constexpr Rational(const Index n = 0, const Index d = 1) noexcept
53 : numer(d ? n / gcd(n, d) : 0), denom(d / gcd(n, d)) {}
54
67 explicit Rational(const String& s);
68
70 void simplify_in_place() noexcept;
71
77 [[nodiscard]] constexpr bool isUndefined() const noexcept { return (denom == 0); }
78
84 [[nodiscard]] constexpr bool isDefined() const noexcept { return not isUndefined(); }
85
92 [[nodiscard]] constexpr bool isIndex(int n=1) const noexcept {
93 return isDefined() and not bool((n*numer) % denom);
94 }
95
103 [[nodiscard]] constexpr Index toIndex(int n=1) const noexcept {
104 return (n*numer) / denom;;
105 }
106
111 [[nodiscard]] constexpr Numeric toNumeric() const noexcept {
112 return Numeric(numer) / Numeric(denom);
113 }
114
122 [[nodiscard]] constexpr int toInt(int n=1) const noexcept { return int(toIndex(n)); }
123
129 constexpr Rational& operator+=(const Rational& a) noexcept {
130 numer = numer * a.denom + a.numer * denom;
131 denom *= a.denom;
132 return *this;
133 }
134
140 constexpr Rational& operator+=(const Index& a) noexcept {
141 numer += denom * a;
142 return *this;
143 }
144
150 constexpr Rational& operator+=(const int& a) noexcept {
151 numer += denom * a;
152 return *this;
153 }
154
160 constexpr Rational& operator-=(const Rational& a) noexcept {
161 numer = numer * a.denom - a.numer * denom;
162 denom *= a.denom;
163 return *this;
164 }
165
171 constexpr Rational& operator-=(const Index& a) noexcept {
172 numer -= denom * a;
173 return *this;
174 }
175
181 constexpr Rational& operator-=(const int& a) noexcept {
182 numer -= denom * a;
183 return *this;
184 }
185
191 constexpr Rational& operator/=(const Rational& a) noexcept {
192 numer *= a.denom;
193 denom *= a.numer;
194 return *this;
195 }
196
202 constexpr Rational& operator/=(const Index& a) noexcept {
203 denom *= a;
204 return *this;
205 }
206
212 constexpr Rational& operator/=(const int& a) noexcept {
213 denom *= a;
214 return *this;
215 }
216
222 constexpr Rational& operator*=(const Rational& a) noexcept {
223 numer *= a.numer;
224 denom *= a.denom;
225 return *this;
226 }
227
233 constexpr Rational& operator*=(const Index& a) noexcept {
234 numer *= a;
235 return *this;
236 }
237
243 constexpr Rational& operator*=(const int& a) noexcept {
244 numer *= a;
245 return *this;
246 }
247
249 constexpr Rational operator++(int) const noexcept{
250 return {numer + denom, denom};
251 }
252
254 constexpr Rational operator--(int) const noexcept {
255 return {numer - denom, denom};
256 }
257
259 constexpr Rational& operator++() noexcept {
260 numer += denom;
261 return *this;
262 }
263
265 constexpr Rational& operator--() noexcept {
266 numer -= denom;
267 return *this;
268 }
269
271 explicit constexpr operator bool() const noexcept {
272 return isDefined() and bool(numer);
273 }
274
277
280
283
286 bif >> numer >> denom;
287 return bif;
288 }
289
291 bofstream& write(bofstream& bof) const {
292 bof << numer << denom;
293 return bof;
294 }
295
297 constexpr Rational& fixSign() noexcept {
298 if (denom < 0) {
299 numer = -numer;
300 denom = -denom;
301 }
302 return *this;
303 }
304
305 friend std::ostream& operator<<(std::ostream& os, const Rational& a);
306
307 friend std::istream& operator>>(std::istream& is, Rational& a);
308}; // Rational;
309
315constexpr Rational reduce_by_gcd(Rational a) noexcept {
316 const Index div = gcd(a.numer, a.denom);
317 if (div)
318 return {a.numer / div, a.denom / div};
319 return a;
320}
321
329constexpr Rational numeric2rational(Numeric x, size_t maxdec=4) noexcept {
330 Index nom=0, denom=1;
331
333 const bool signchange = x < 0;
334 x = signchange ? -x : x;
335
336 // Add numbers by keeping the floor
337 size_t i=0;
338 do {
339 const auto xi=Index(x);
340 nom += xi;
341 x = 10 * (x - Numeric(xi));
342 nom *= 10;
343 denom *= 10;
344 i++;
345 } while (i<=maxdec);
346
347 // Fix possible rounding error
348 if (x >= 5)
349 nom += 10;
350
351 // Change sign or not
352 return signchange ? Rational{-nom, denom} : Rational{nom, denom};
353}
354
355
356// An undefined rational to be used everywhere for a undefined rationals
357#define RATIONAL_UNDEFINED Rational(0, 0)
358
364constexpr Rational operator-(const Rational a) noexcept {
365 return {-a.numer, a.denom};
366}
367
373constexpr Rational operator+(const Rational a) noexcept { return a; }
374
381constexpr Rational operator+(const Rational a, const Rational b) noexcept {
382 return (a.denom == b.denom)
383 ? Rational(a.numer + b.numer, a.denom)
384 : Rational(a.numer * b.denom + b.numer * a.denom,
385 a.denom * b.denom);
386}
387
394constexpr Rational operator+(const Rational a, Index b) noexcept {
395 return {a.numer + b * a.denom, a.denom};
396}
397
404constexpr Rational operator+(const Rational a, int b) noexcept {
405 return {a.numer + b * a.denom, a.denom};
406}
407
414constexpr Rational operator+(Index b, const Rational a) noexcept { return operator+(a, b); }
415
422constexpr Rational operator+(int b, const Rational a) noexcept { return operator+(a, b); }
423
430constexpr Rational operator-(const Rational a, const Rational b) noexcept {
431 return (a.denom == b.denom)
432 ? Rational(a.numer - b.numer, a.denom)
433 : Rational(a.numer * b.denom - b.numer * a.denom,
434 a.denom * b.denom);
435}
436
443constexpr Rational operator-(const Rational a, Index b) noexcept {
444 return {a.numer - b * a.denom, a.denom};
445}
446
453constexpr Rational operator-(const Rational a, int b) noexcept {
454 return {a.numer - b * a.denom, a.denom};
455}
456
463constexpr Rational operator-(Index b, const Rational a) noexcept {
464 return {-a.numer + b * a.denom, a.denom};
465}
466
473constexpr Rational operator-(int b, const Rational a) noexcept {
474 return {-a.numer + b * a.denom, a.denom};
475}
476
483constexpr Rational operator/(const Rational a, const Rational b) noexcept {
484 return {a.numer * b.denom, a.denom * b.numer};
485}
486
493constexpr Rational operator/(const Rational a, Index b) noexcept {
494 return {a.numer, a.denom * b};
495}
496
503constexpr Rational operator/(const Rational a, int b) noexcept {
504 return {a.numer, a.denom * b};
505}
506
513constexpr Rational operator/(Index b, const Rational a) noexcept {
514 return {a.denom * b, a.numer};
515}
516
523constexpr Rational operator/(int b, const Rational a) noexcept {
524 return {a.denom * b, a.numer};
525}
526
533constexpr Rational operator*(const Rational a, const Rational b) noexcept {
534 return {a.numer * b.numer, a.denom * b.denom};
535}
536
543constexpr Rational operator*(const Rational a, Index b) noexcept {
544 return {a.numer * b, a.denom};
545}
546
553constexpr Rational operator*(const Rational a, int b) noexcept {
554 return {a.numer * b, a.denom};
555}
556
563constexpr Rational operator*(Index b, const Rational a) noexcept { return operator*(a, b); }
564
571constexpr Rational operator*(int b, const Rational a) noexcept { return operator*(a, b); }
572
579constexpr Rational operator%(const Rational a, const Rational b) noexcept {
580 return (a.denom == b.denom)
581 ? Rational(a.numer % b.numer, a.denom)
582 : Rational((a.numer * b.denom) % (a.denom * b.numer),
583 a.denom * b.denom);
584}
585
592constexpr Rational operator%(const Rational a, Index b) noexcept {
593 return {a.numer % (a.denom * b), a.denom};
594}
595
602constexpr Rational operator%(const Rational a, int b) noexcept {
603 return {a.numer % (a.denom * b), a.denom};
604}
605
612constexpr Rational operator%(Index b, const Rational a) noexcept {
613 return {(b * a.denom) % a.numer, a.denom};
614}
615
622constexpr Rational operator%(int b, const Rational a) noexcept {
623 return {(b * a.denom) % a.numer, a.denom};
624}
625
633constexpr bool operator==(const Rational a, const Rational b) noexcept {
634 return a.isDefined() and b.isDefined() and
635 a.numer * b.denom == a.denom * b.numer;
636}
637
645constexpr bool operator!=(const Rational a, const Rational b) noexcept {
646 return not (a.isDefined() and b.isDefined() and operator==(a, b));
647}
648
656constexpr bool operator<(const Rational a, const Rational b) noexcept {
657 return a.isDefined() and b.isDefined() and
658 a.numer * b.denom < a.denom * b.numer;
659}
660
668constexpr bool operator>(const Rational a, const Rational b) noexcept { return operator<(b, a) and a.isDefined() and b.isDefined(); }
669
677constexpr bool operator<=(const Rational a, const Rational b) noexcept {
678 return not operator>(a, b) and a.isDefined() and b.isDefined();
679}
680
688constexpr bool operator>=(const Rational a, const Rational b) noexcept {
689 return not operator<(a, b) and a.isDefined() and b.isDefined();
690}
691
698constexpr bool operator!(const Rational a) noexcept { return a.numer and a.isDefined(); }
699
705inline Numeric sqrt(const Rational r) { return std::sqrt(r.toNumeric()); }
706
713inline Numeric pow(const Rational base, Numeric exp) {
714 return std::pow(base.toNumeric(), exp);
715}
716
723inline Numeric pow(Numeric base, const Rational exp) {
724 return std::pow(base, exp.toNumeric());
725}
726
733inline Numeric pow(const Rational base, const Rational exp) {
734 return pow(base, exp.toNumeric());
735}
736
743constexpr bool operator<(const Index a, const Rational b) noexcept {
744 return Rational(a, 1) < b and b.isDefined();
745}
746
753constexpr bool operator<(const int a, const Rational b) noexcept {
754 return Rational(a, 1) < b and b.isDefined();
755}
756
763constexpr bool operator<(const Rational a, const Index b) noexcept {
764 return a < Rational(b, 1) and a.isDefined();
765}
766
773constexpr bool operator<(const Rational a, const int b) noexcept {
774 return a < Rational(b, 1) and a.isDefined();
775}
776
783constexpr bool operator>(const Index a, const Rational b) noexcept {
784 return Rational(a, 1) > b and b.isDefined();
785}
786
793constexpr bool operator>(const int a, const Rational b) noexcept {
794 return Rational(a, 1) > b and b.isDefined();
795}
796
803constexpr bool operator>(const Rational a, const Index b) noexcept {
804 return a > Rational(b, 1) and a.isDefined();
805}
806
813constexpr bool operator>(const Rational a, const int b) noexcept {
814 return a > Rational(b, 1) and a.isDefined();
815}
816
823constexpr bool operator==(const Rational a, const Index b) noexcept {
824 return a == Rational(b, 1) and a.isDefined();
825}
826
833constexpr bool operator==(const Rational a, const int b) noexcept {
834 return a == Rational(b, 1) and a.isDefined();
835}
836
843constexpr bool operator!=(const Rational a, const Index b) noexcept {
844 return not(a == b) and a.isDefined();
845}
846
853constexpr bool operator!=(const Rational a, const int b) noexcept {
854 return not(a == b) and a.isDefined();
855}
856
862constexpr Rational abs(const Rational a) noexcept {
863 return a < 0 ? -a : a;
864}
865
872constexpr Rational max(const Rational a, const Rational b) noexcept {
873 return a < b ? b : a;
874} // Let other operators find out if this is allowed instead
875
882constexpr Rational min(const Rational a, const Rational b) noexcept {
883 return a < b ? a : b;
884} // Let other operators find out if this is allowed instead
885
887
893constexpr Rational operator ""_2(unsigned long long int n) noexcept {
894 return Rational(n, 2);
895};
896
902constexpr bool iseven(const Rational r) noexcept {
903 return 0 == (r % 2);
904}
905
907constexpr Numeric operator*(Rational y, Numeric x) noexcept {return y.toNumeric() * x;}
908constexpr Numeric operator*(Numeric x, Rational y) noexcept {return x * y.toNumeric();}
909
911constexpr Numeric operator/(Rational y, Numeric x) noexcept {return y.toNumeric() / x;}
912constexpr Numeric operator/(Numeric x, Rational y) noexcept {return x / y.toNumeric();}
913
915constexpr Numeric operator+(Rational y, Numeric x) noexcept {return y.toNumeric() + x;}
916constexpr Numeric operator+(Numeric x, Rational y) noexcept {return x + y.toNumeric();}
917
919constexpr Numeric operator-(Rational y, Numeric x) noexcept {return y.toNumeric() - x;}
920constexpr Numeric operator-(Numeric x, Rational y) noexcept {return x - y.toNumeric();}
921
922#endif // rational_h
923
This file contains the definition of Array.
base max(const Array< base > &x)
Max function.
Definition: array.h:145
base min(const Array< base > &x)
Min function.
Definition: array.h:161
This file contains the class declaration of bifstream.
This file contains the class declaration of bofstream.
This can be used to make arrays out of anything.
Definition: array.h:48
Binary output file stream class.
Definition: bifstream.h:43
Binary output file stream class.
Definition: bofstream.h:42
void abs(Sparse &A, const Sparse &B)
Absolute value of sparse matrix elements.
Definition: matpackII.cc:394
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
constexpr bool operator>(const Rational a, const Rational b) noexcept
More than.
Definition: rational.h:668
constexpr Rational operator+(const Rational a) noexcept
Positive.
Definition: rational.h:373
constexpr bool operator>=(const Rational a, const Rational b) noexcept
More than or equal to.
Definition: rational.h:688
Numeric pow(const Rational base, Numeric exp)
Power of.
Definition: rational.h:713
constexpr bool operator!(const Rational a) noexcept
Not.
Definition: rational.h:698
constexpr Rational operator-(const Rational a) noexcept
Negative.
Definition: rational.h:364
constexpr bool operator!=(const Rational a, const Rational b) noexcept
Inequality.
Definition: rational.h:645
constexpr bool operator<(const Rational a, const Rational b) noexcept
Less than.
Definition: rational.h:656
Numeric sqrt(const Rational r)
Square root.
Definition: rational.h:705
constexpr bool iseven(const Rational r) noexcept
Returns true if even integer.
Definition: rational.h:902
constexpr bool operator==(const Rational a, const Rational b) noexcept
Equality.
Definition: rational.h:633
constexpr Rational numeric2rational(Numeric x, size_t maxdec=4) noexcept
Rational from Numeric.
Definition: rational.h:329
constexpr Rational operator%(const Rational a, const Rational b) noexcept
Remainder.
Definition: rational.h:579
constexpr bool operator<=(const Rational a, const Rational b) noexcept
Less than or equal to.
Definition: rational.h:677
constexpr Rational reduce_by_gcd(Rational a) noexcept
Returns the rational reduced by the greates.
Definition: rational.h:315
constexpr Rational operator*(const Rational a, const Rational b) noexcept
Multiplication.
Definition: rational.h:533
constexpr Rational operator/(const Rational a, const Rational b) noexcept
Division.
Definition: rational.h:483
Implements rational numbers to work with other ARTS types.
Definition: rational.h:43
constexpr Rational & operator++() noexcept
Definition: rational.h:259
constexpr bool isUndefined() const noexcept
Is the object not defined.
Definition: rational.h:77
constexpr Rational & operator-=(const int &a) noexcept
Remove from this.
Definition: rational.h:181
constexpr Rational(const Index n=0, const Index d=1) noexcept
Initialization call.
Definition: rational.h:52
constexpr bool isDefined() const noexcept
Is the object defined.
Definition: rational.h:84
constexpr Rational & operator*=(const int &a) noexcept
Multiply by this.
Definition: rational.h:243
constexpr Numeric toNumeric() const noexcept
Converts this to a Numeric.
Definition: rational.h:111
constexpr Rational & operator--() noexcept
Remove one if possible.
Definition: rational.h:265
constexpr Rational & operator-=(const Index &a) noexcept
Remove from this.
Definition: rational.h:171
constexpr Rational & operator+=(const Rational &a) noexcept
Definition: rational.h:129
constexpr Rational & operator/=(const Index &a) noexcept
Divide by this.
Definition: rational.h:202
constexpr Index toIndex(int n=1) const noexcept
Converts the value to index by n-scaled division.
Definition: rational.h:103
constexpr Rational & operator+=(const Index &a) noexcept
Definition: rational.h:140
Definition: rational.h:285
void simplify_in_place() noexcept
Simplify by reducing the values locally.
Definition: rational.cc:86
constexpr Rational & operator*=(const Rational &a) noexcept
Multiply by this.
Definition: rational.h:222
friend std::ostream & operator<<(std::ostream &os, const Rational &a)
Definition: rational.cc:34
friend std::istream & operator>>(std::istream &is, Rational &a)
Definition: rational.cc:45
constexpr Rational & operator/=(const Rational &a) noexcept
Divide by this.
Definition: rational.h:191
constexpr Rational & operator+=(const int &a) noexcept
Definition: rational.h:150
bofstream & write(bofstream &bof) const
Binary write for Rational.
Definition: rational.h:291
constexpr Rational operator--(int) const noexcept
Remove one if possible.
Definition: rational.h:254
Index numer
Definition: rational.h:44
Index denom
Definition: rational.h:45
constexpr Rational & fixSign() noexcept
Makes the sign of denom positive.
Definition: rational.h:297
constexpr Rational operator++(int) const noexcept
Definition: rational.h:249
constexpr bool isIndex(int n=1) const noexcept
Is the object a n-scaled Index.
Definition: rational.h:92
constexpr Rational & operator/=(const int &a) noexcept
Divide by this.
Definition: rational.h:212
constexpr int toInt(int n=1) const noexcept
Converts the value to int by n-scaled division in Index form.
Definition: rational.h:122
constexpr Rational & operator-=(const Rational &a) noexcept
Remove from this.
Definition: rational.h:160
constexpr Rational & operator*=(const Index &a) noexcept
Multiply by this.
Definition: rational.h:233
#define d
#define a
#define b