stdlib/private/rational.ajla

   1 {*
   2  * Copyright (C) 2024 Mikulas Patocka
   3  *
   4  * This file is part of Ajla.
   5  *
   6  * Ajla is free software: you can redistribute it and/or modify it under the
   7  * terms of the GNU General Public License as published by the Free Software
   8  * Foundation, either version 3 of the License, or (at your option) any later
   9  * version.
  10  *
  11  * Ajla is distributed in the hope that it will be useful, but WITHOUT ANY
  12  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  13  * A PARTICULAR PURPOSE. See the GNU General Public License for more details.
  14  *
  15  * You should have received a copy of the GNU General Public License along with
  16  * Ajla. If not, see <https://www.gnu.org/licenses/>.
  17  *}
  18
  19 private unit private.rational;
  20
  21 fn rational_zero : rational;
  22 fn rational_one : rational;
  23 fn rational_pi : rational;
  24 fn rational_neg(i1 : rational) : rational;
  25 fn rational_add(i1 i2 : rational) : rational;
  26 fn rational_subtract(i1 i2 : rational) : rational;
  27 fn rational_multiply(i1 i2 : rational) : rational;
  28 fn rational_recip(i1 : rational) : rational;
  29 fn rational_divide(i1 i2 : rational) : rational;
  30 fn rational_modulo(i1 i2 : rational) : rational;
  31 fn rational_power(i1 i2 : rational) : rational;
  32 fn rational_ldexp(i1 i2 : rational) : rational;
  33 fn rational_atan2(i1 i2 : rational) : rational;
  34 fn rational_sqrt(i1 : rational) : rational;
  35 fn rational_cbrt(i1 : rational) : rational;
  36 fn rational_sin(i1 : rational) : rational;
  37 fn rational_cos(i1 : rational) : rational;
  38 fn rational_tan(i1 : rational) : rational;
  39 fn rational_asin(i1 : rational) : rational;
  40 fn rational_acos(i1 : rational) : rational;
  41 fn rational_atan(i1 : rational) : rational;
  42 fn rational_sinh(i1 : rational) : rational;
  43 fn rational_cosh(i1 : rational) : rational;
  44 fn rational_tanh(i1 : rational) : rational;
  45 fn rational_asinh(i1 : rational) : rational;
  46 fn rational_acosh(i1 : rational) : rational;
  47 fn rational_atanh(i1 : rational) : rational;
  48 fn rational_log2(i1 : rational) : rational;
  49 fn rational_log(i1 : rational) : rational;
  50 fn rational_log10(i1 : rational) : rational;
  51 fn rational_exp2(i1 : rational) : rational;
  52 fn rational_exp(i1 : rational) : rational;
  53 fn rational_exp10(i1 : rational) : rational;
  54 fn rational_round(i1 : rational) : rational;
  55 fn rational_ceil(i1 : rational) : rational;
  56 fn rational_floor(i1 : rational) : rational;
  57 fn rational_trunc(i1 : rational) : rational;
  58 fn rational_fract(i1 : rational) : rational;
  59 fn rational_mantissa(i1 : rational) : rational;
  60 fn rational_exponent(i1 : rational) : rational;
  61 fn rational_next_number(i1 : rational) : rational;
  62 fn rational_prev_number(i1 : rational) : rational;
  63 fn rational_is_negative(i1 : rational) : bool;
  64 fn rational_is_infinity(i1 : rational) : bool;
  65 fn rational_equal(i1 i2 : rational) : bool;
  66 fn rational_less(i1 i2 : rational) : bool;
  67 fn rational_to_int(i1 : rational) : int;
  68 fn rational_from_int(i1 : int) : rational;
  69 fn rational_to_bytes(i1 : rational) : bytes;
  70 fn rational_to_bytes_base_precision(i1 : rational, b : int, digits : int) : bytes;
  71 fn rational_from_bytes(a1 : bytes) : rational;
  72 fn rational_from_bytes_base(a1 : bytes, b : int) : rational;
  73
  74 implementation
  75
  76 uses exception;
  77 uses private.show;
  78
  79 fn gcd(a b : int) : int
  80 [
  81         a := abs(a);
  82         b := abs(b);
  83         while b <> 0 do
  84                 a, b := b, a mod b;
  85         return a;
  86 ]
  87
  88 fn rational_normalize(i : rational) : rational
  89 [
  90         if i.den = 0 then [
  91                 if i.num = 0 then
  92                         abort exception_make(rational, ec_sync, error_nan, 0, true);
  93                 if i.num > 0 then
  94                         i.num := 1;
  95                 else
  96                         i.num := -1;
  97                 return i;
  98         ]
  99         var g := gcd(i.num, i.den);
 100         i.num div= g;
 101         i.den div= g;
 102         if i.den < 0 and i.num <> 0 then [
 103                 i.num := -i.num;
 104                 i.den := -i.den;
 105         ]
 106         return i;
 107 ]
 108
 109 fn rational_zero : rational
 110 [
 111         return rational.[ num : 0, den : 1 ];
 112 ]
 113
 114 fn rational_one : rational
 115 [
 116         return rational.[ num : 1, den : 1 ];
 117 ]
 118
 119 fn rational_pi : rational
 120 [
 121         return instance_real_number_real64.pi;
 122 ]
 123
 124 fn rational_neg(i1 : rational) : rational
 125 [
 126         if i1.num = 0 then
 127                 return rational.[
 128                         num : i1.num,
 129                         den : -i1.den,
 130                 ];
 131         return rational.[
 132                 num : -i1.num,
 133                 den : i1.den,
 134         ];
 135 ]
 136
 137 fn rational_add(i1 i2 : rational) : rational
 138 [
 139         var r := rational.[
 140                 num : i1.num * i2.den + i1.den * i2.num,
 141                 den : i1.den * i2.den,
 142         ];
 143         return rational_normalize(r);
 144 ]
 145
 146 fn rational_subtract(i1 i2 : rational) : rational
 147 [
 148         var r := rational.[
 149                 num : i1.num * i2.den - i1.den * i2.num,
 150                 den : i1.den * i2.den,
 151         ];
 152         return rational_normalize(r);
 153 ]
 154
 155 fn rational_multiply(i1 i2 : rational) : rational
 156 [
 157         var r := rational.[
 158                 num : i1.num * i2.num,
 159                 den : i1.den * i2.den,
 160         ];
 161         return rational_normalize(r);
 162 ]
 163
 164 fn rational_recip(i1 : rational) : rational
 165 [
 166         var r := rational.[
 167                 num : i1.den,
 168                 den : i1.num,
 169         ];
 170         return r;
 171 ]
 172
 173 fn rational_divide(i1 i2 : rational) : rational
 174 [
 175         var r := rational.[
 176                 num : i1.num * i2.den,
 177                 den : i1.den * i2.num,
 178         ];
 179         return rational_normalize(r);
 180 ]
 181
 182 fn rational_modulo(i1 i2 : rational) : rational
 183 [
 184         var r := trunc(i1 / i2);
 185         return i1 - i2 * r;
 186 ]
 187
 188 fn rational_power(i1 i2 : rational) : rational
 189 [
 190         if i2.den = 1 then [
 191                 i1.num := ipower(i1.num, i2.num);
 192                 i1.den := ipower(i1.den, i2.num);
 193                 return i1;
 194         ]
 195         var r1 : real64 := i1;
 196         var r2 : real64 := i2;
 197         return power(r1, r2);
 198 ]
 199
 200 fn rational_ldexp(i1 i2 : rational) : rational
 201 [
 202         if i2.den = 1 then [
 203                 if i2.num >= 0 then [
 204                         i1.num shl= i2.num;
 205                         return rational_normalize(i1);
 206                 ] else [
 207                         i1.den shl= -i2.num;
 208                         return rational_normalize(i1);
 209                 ]
 210         ]
 211         var r1 : real64 := i1;
 212         var r2 : real64 := i2;
 213         return ldexp(r1, r2);
 214 ]
 215
 216 fn rational_atan2(i1 i2 : rational) : rational
 217 [
 218         var r1 : real64 := i1;
 219         var r2 : real64 := i2;
 220         return atan2(r1, r2);
 221 ]
 222
 223 fn rational_sqrt(i1 : rational) : rational [ var r : real64 := i1; return sqrt(r); ]
 224 fn rational_cbrt(i1 : rational) : rational [ var r : real64 := i1; return cbrt(r); ]
 225 fn rational_sin(i1 : rational) : rational [ var r : real64 := i1; return sin(r); ]
 226 fn rational_cos(i1 : rational) : rational [ var r : real64 := i1; return cos(r); ]
 227 fn rational_tan(i1 : rational) : rational [ var r : real64 := i1; return tan(r); ]
 228 fn rational_asin(i1 : rational) : rational [ var r : real64 := i1; return asin(r); ]
 229 fn rational_acos(i1 : rational) : rational [ var r : real64 := i1; return acos(r); ]
 230 fn rational_atan(i1 : rational) : rational [ var r : real64 := i1; return atan(r); ]
 231 fn rational_sinh(i1 : rational) : rational [ var r : real64 := i1; return sinh(r); ]
 232 fn rational_cosh(i1 : rational) : rational [ var r : real64 := i1; return cosh(r); ]
 233 fn rational_tanh(i1 : rational) : rational [ var r : real64 := i1; return tanh(r); ]
 234 fn rational_asinh(i1 : rational) : rational [ var r : real64 := i1; return asinh(r); ]
 235 fn rational_acosh(i1 : rational) : rational [ var r : real64 := i1; return acosh(r); ]
 236 fn rational_atanh(i1 : rational) : rational [ var r : real64 := i1; return atanh(r); ]
 237 fn rational_log2(i1 : rational) : rational [ var r : real64 := i1; return log2(r); ]
 238 fn rational_log(i1 : rational) : rational [ var r : real64 := i1; return log(r); ]
 239 fn rational_log10(i1 : rational) : rational [ var r : real64 := i1; return log10(r); ]
 240 fn rational_exp2(i1 : rational) : rational [ var r : real64 := i1; return exp2(r); ]
 241 fn rational_exp(i1 : rational) : rational [ var r : real64 := i1; return exp(r); ]
 242 fn rational_exp10(i1 : rational) : rational [ var r : real64 := i1; return exp10(r); ]
 243
 244 fn rational_round(i1 : rational) : rational
 245 [
 246         if abs(rational_fract(i1)) = 0.5 then [
 247                 var i : int := i1;
 248                 if i bt 0 = is_negative i1 then
 249                         return i - (i and 1);
 250         ]
 251         return rational_floor(i1 + 0.5);
 252 ]
 253
 254 fn rational_ceil(i1 : rational) : rational
 255 [
 256         var i : int := i1;
 257         if i1 > i then
 258                 i += 1;
 259         return i;
 260 ]
 261
 262 fn rational_floor(i1 : rational) : rational
 263 [
 264         var i : int := i1;
 265         if i1 < i then
 266                 i -= 1;
 267         return i;
 268 ]
 269
 270 fn rational_trunc(i1 : rational) : rational
 271 [
 272         var i : int := i1;
 273         return i;
 274 ]
 275
 276 fn rational_fract(i1 : rational) : rational
 277 [
 278         var r := rational_trunc(i1);
 279         return i1 - r;
 280 ]
 281
 282 fn rational_mantissa(i1 : rational) : rational [ var r : real64 := i1; return mantissa(r); ]
 283 fn rational_exponent(i1 : rational) : rational [ var r : real64 := i1; return exponent(r); ]
 284 fn rational_next_number(i1 : rational) : rational [ var r : real64 := i1; return next_number(r); ]
 285 fn rational_prev_number(i1 : rational) : rational [ var r : real64 := i1; return prev_number(r); ]
 286
 287 fn rational_is_negative(i1 : rational) : bool
 288 [
 289         if i1.num = 0 then
 290                 return i1.den < 0;
 291         return i1.num < 0;
 292 ]
 293
 294 fn rational_is_infinity(i1 : rational) : bool
 295 [
 296         return i1.den = 0;
 297 ]
 298
 299 fn rational_equal(i1 i2 : rational) : bool
 300 [
 301         return i1.num * i2.den = i1.den * i2.num;
 302 ]
 303
 304 fn rational_less(i1 i2 : rational) : bool
 305 [
 306         return i1.num * i2.den < i1.den * i2.num;
 307 ]
 308
 309 fn rational_to_int(i1 : rational) : int
 310 [
 311         if i1.den = 0 then
 312                 abort exception_make(int, ec_sync, error_infinity, 0, true);
 313         return i1.num div i1.den;
 314 ]
 315
 316 fn rational_from_int(i1 : int) : rational
 317 [
 318         return rational.[
 319                 num : i1,
 320                 den : 1,
 321         ];
 322 ]
 323
 324 fn rational_to_bytes(i1 : rational) : bytes
 325 [
 326         return real_to_bytes(i1);
 327 ]
 328
 329 fn rational_to_bytes_base_precision(i1 : rational, b : int, digits : int) : bytes
 330 [
 331         return real_to_bytes_base_precision(i1, b, digits);
 332 ]
 333
 334 fn rational_from_bytes(a1 : bytes) : rational
 335 [
 336         return bytes_to_real(rational, instance_real_number_rational, a1);
 337 ]
 338
 339 fn rational_from_bytes_base(a1 : bytes, b : int) : rational
 340 [
 341         return bytes_to_real_base(rational, instance_real_number_rational, a1, b);
 342 ]