File: ca.sh
   1 #!/bin/sh
   2 
   3 # The MIT License (MIT)
   4 #
   5 # Copyright © 2024 pacman64
   6 #
   7 # Permission is hereby granted, free of charge, to any person obtaining a copy
   8 # of this software and associated documentation files (the “Software”), to deal
   9 # in the Software without restriction, including without limitation the rights
  10 # to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
  11 # copies of the Software, and to permit persons to whom the Software is
  12 # furnished to do so, subject to the following conditions:
  13 #
  14 # The above copyright notice and this permission notice shall be included in
  15 # all copies or substantial portions of the Software.
  16 #
  17 # THE SOFTWARE IS PROVIDED “AS IS”, WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  18 # IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  19 # FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
  20 # AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  21 # LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
  22 # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
  23 # SOFTWARE.
  24 
  25 
  26 # ca [expression]
  27 #
  28 #
  29 # CAlculator is an easier-to-use version of bc (basic calculator) where
  30 #
  31 #   - you give the expression as an argument, while bc reads it from stdin
  32 #   - you can use either ** or ^ to raise powers
  33 #   - you can use [ and ] or ( and ) interchangeably
  34 #   - the number of decimals is 20 by default
  35 #   - automatically includes the extended bc math library via option -l
  36 #   - there are several extra predefined values and functions (see below)
  37 #
  38 #
  39 # Predefined functions and values, besides the ones `bc` already comes with
  40 #
  41 #
  42 # nano   a billionth
  43 # micro  a millionth
  44 # milli  a thousandth
  45 #
  46 # kilo   1000
  47 # mega   a million
  48 # giga   a billion
  49 # tera   a thousand billion / trillion
  50 # peta   a million billion / quadrillion
  51 # exa    a quintillion
  52 # zetta  a sixtillion
  53 #
  54 # binkilo   1024
  55 # binmega   1024^2
  56 # bingiga   1024^3
  57 # bintera   1024^4
  58 # binpeta   1024^5
  59 # binexa    1024^6
  60 # binzetta  1024^7
  61 #
  62 # kb  1024
  63 # mb  1024^2
  64 # gb  1024^3
  65 # tb  1024^4
  66 # pb  1024^5
  67 #
  68 # kib  1024
  69 # mib  1024^2
  70 # gib  1024^3
  71 # tib  1024^4
  72 # pib  1024^5
  73 #
  74 # mol      1 mole, the exact avogadro number
  75 # mole     the exact avogadro number
  76 #
  77 # cup2l    liters in 1 cup
  78 # floz2ml  milliliters in 1 fluid ounce
  79 # gal2l    liters in 1 gallon
  80 # in2cm    centimeters in 1 inch
  81 # lb2kg    kilograms in 1 pound
  82 # mi2km    kilometers in 1 mile
  83 # mpg2kpl  kilometers per liter in 1 mile per gallon
  84 # nmi2km   kilometers in 1 nautical mile
  85 # oz2g     grams in 1 weight ounce
  86 # psi2pa   pascals in 1 pound per square-inch
  87 # ton2kg   kilograms in 1 american ton
  88 # yd2m     meters in 1 yard
  89 #
  90 # ftin(feet, inches)  turn a feet/inches pair into meters
  91 # lboz(feet, inches)  turn a pound/ounces pair into kilograms
  92 #
  93 # eu()          approximate Euler's number
  94 # euler()       approximate Euler's number
  95 # pi()          approximate any circle's circumference/diameter ratio
  96 # tau()         2*pi; approximate any circle's circumference/radius ratio
  97 # abs(x)        the absolute value function
  98 # exp(x)        the exponential function; alias for bc's e
  99 # ln(x)         natural logarithm; an alias for bc's l
 100 # log(x)        natural logarithm; an alias for bc's l
 101 # log2(x)       base-2 logarithm
 102 # log10(x)      base-10 logarithm
 103 # sin(x)        the sine function
 104 # cos(x)        the cosine function
 105 # atan(x)       the arc-tangent function
 106 # hypot(x, y)   length of the hypothenuse
 107 # mod(x, y)     modulus function, which is different from the % operator
 108 # r(x, d)       round to the given number of decimal digits; can be negative
 109 # round(x, d)   round to the given number of decimal digits; can be negative
 110 # r0(x)         round to integer
 111 # r1(x)         round to 1 decimal place
 112 # ...
 113 # r20(x)        round to 20 decimal places
 114 # round0(x)     round to integer
 115 # round1(x)     round to 1 decimal place
 116 # ...
 117 # round20(x)    round to 20 decimal places
 118 # f(x)          the factorial function
 119 # fac(x)        the factorial function
 120 # fact(x)       the factorial function
 121 # factorial(x)  the factorial function
 122 # com(n, k)     number of combinations
 123 # comb(n, k)    number of combinations
 124 # p(n, k)       number of permutations
 125 # per(n, k)     number of permutations
 126 # perm(n, k)    number of permutations
 127 # lcm(x, y)     least-common multiple
 128 # gcd(x, y)     greatest-common divisor
 129 # min(x, y)     the minimum of 2 numbers
 130 # max(x, y)     the maximum of 2 numbers
 131 # sgn(x)        the sign (-1, 0, or +1) of a number
 132 # sign(x)       the sign (-1, 0, or +1) of a number
 133 # isprime(x)    check if number is prime
 134 # lgamma(x)     lanczos approximation of the log-gamma function
 135 # lbeta(x, y)   lanczos approximation of the log-beta function
 136 # gamma(x)      approximation of the gamma function
 137 # beta(x)       approximation of the beta function
 138 # logistic(x)   the logistic function
 139 # mix(x, y, k)  interpolate 2 numbers; extrapolates when k is outside [0, 1]
 140 # clamp(x, min, max)   trap a number in the range given
 141 # wrap(x, min, max)    normalize a number in the range given
 142 # unwrap(x, min, max)  expand a normalized number into the range given
 143 
 144 
 145 # handle help options
 146 case "$1" in
 147     -h|--h|-help|--help)
 148         awk '/^# +ca/, /^$/ { gsub(/^# ?/, ""); print }' "$0"
 149         exit 0
 150     ;;
 151 esac
 152 
 153 # handle optional scale-override parameter
 154 scale=20
 155 if [ $# -eq 2 ]; then
 156     scale="$1"
 157     [ $# -gt 0 ] && shift
 158 fi
 159 
 160 # ensure output is all on 1 line, define several funcs and values, then
 161 # inject the expression given as this script's arguments, transformed
 162 # according to the rules described above
 163 {
 164 BC_LINE_LENGTH=0 bc -l << ENDOFSCRIPT
 165 scale = ${scale};
 166 
 167 femto = 0.000000000000001;
 168 pico = 0.000000000001;
 169 nano = 0.000000001;
 170 micro = 0.000001;
 171 milli = 0.001;
 172 
 173 kilo = 1000;
 174 mega = 1000 * kilo;
 175 giga = 1000 * mega;
 176 tera = 1000 * giga;
 177 peta = 1000 * tera;
 178 exa =  1000 * peta;
 179 zetta =  1000 * exa;
 180 
 181 binkilo = 1024;
 182 binmega = 1024 * binkilo;
 183 bingiga = 1024 * binmega;
 184 bintera = 1024 * bingiga;
 185 binpeta = 1024 * bintera;
 186 binexa = 1024 * binpeta;
 187 binzetta = 1024 * binexa;
 188 
 189 kb = 1024;
 190 mb = 1024 * kb;
 191 gb = 1024 * mb;
 192 tb = 1024 * gb;
 193 pb = 1024 * tb;
 194 eb = 1024 * pb;
 195 zb = 1024 * eb;
 196 
 197 kib = 1024;
 198 mib = 1024 * kib;
 199 gib = 1024 * mib;
 200 tib = 1024 * gib;
 201 pib = 1024 * tib;
 202 zib = 1024 * pib;
 203 
 204 mol = 602214076000000000000000;
 205 mole = 602214076000000000000000;
 206 
 207 cup2l = 0.23658824;
 208 floz2l = 0.0295735295625;
 209 floz2ml = 29.5735295625;
 210 ft2m = 0.3048;
 211 gal2l = 3.785411784;
 212 in2cm = 2.54;
 213 lb2kg = 0.45359237;
 214 mi2km = 1.609344;
 215 mpg2kpl = 0.425143707;
 216 nmi2km = 1.852;
 217 oz2g = 28.349523125
 218 psi2pa = 6894.757293168;
 219 ton2kg = 907.18474;
 220 yd2m = 0.9144;
 221 
 222 ga2l = gal2l;
 223 nm2km = nmi2km;
 224 tn2kg = ton2kg;
 225 
 226 hour = 3600;
 227 day = 24 * hour;
 228 week = 7 * day;
 229 
 230 hr = hour;
 231 wk = week;
 232 
 233 define ftin(f, i) {
 234     return (0.3048 * f + 0.0254 * i);
 235 }
 236 
 237 define lboz(l, o) {
 238     return (0.45359237 * l + 0.028349523 * o);
 239 }
 240 
 241 define eu() {
 242     return (e(1));
 243 }
 244 
 245 define euler() {
 246     return (e(1));
 247 }
 248 
 249 define pi() {
 250     return (4*a(1));
 251 }
 252 
 253 define tau() {
 254     return (8*a(1));
 255 }
 256 
 257 define deg(x) {
 258     return (180 * x / pi());
 259 }
 260 
 261 define rad(x) {
 262     return (pi() * x / 180);
 263 }
 264 
 265 define abs(x) {
 266     if (x >= 0) return (x);
 267     return (-x);
 268 }
 269 
 270 define exp(x) {
 271     return (e(x));
 272 }
 273 
 274 define ln(x) {
 275     return (l(x));
 276 }
 277 
 278 define log(x) {
 279     return (l(x));
 280 }
 281 
 282 define log2(x) {
 283     return (l(x) / l(2));
 284 }
 285 
 286 define log10(x) {
 287     return (l(x) / l(10));
 288 }
 289 
 290 define sin(x) {
 291     return (s(x));
 292 }
 293 
 294 define cos(x) {
 295     return (c(x));
 296 }
 297 
 298 define tan(x) {
 299     return (s(x) / c(x));
 300 }
 301 
 302 define cot(x) {
 303     return (c(x) / s(x));
 304 }
 305 
 306 define atan(x) {
 307     return (a(x));
 308 }
 309 
 310 define sinh(x) {
 311     return ((e(x) - e(-x)) / 2);
 312 }
 313 
 314 define cosh(x) {
 315     return ((e(x) + e(-x)) / 2);
 316 }
 317 
 318 define tanh(x) {
 319     return ((e(x) - e(-x)) / (e(x) + e(-x)));
 320 }
 321 
 322 define coth(x) {
 323     return ((e(x) + e(-x)) / (e(x) - e(-x)));
 324 }
 325 
 326 define hypot(x, y) {
 327     return (sqrt(x*x + y*y));
 328 }
 329 
 330 define sinc(x) {
 331     if (x == 0) return (1);
 332     return (s(x) / x);
 333 }
 334 
 335 define min(x, y) {
 336     if (x <= y) return (x);
 337     return (y);
 338 }
 339 
 340 define max(x, y) {
 341     if (x >= y) return (x);
 342     return (y);
 343 }
 344 
 345 define mod(x, y) {
 346     auto s, m;
 347     s = scale;
 348     scale = 0;
 349     m = x % y;
 350     scale = s;
 351     return (m);
 352 }
 353 
 354 define mod1(x) {
 355     return (mod(x, 1));
 356 }
 357 
 358 define modf(x) {
 359     return (mod(x, 1));
 360 }
 361 
 362 define round0(x) {
 363     auto i;
 364     i = x - mod(x, 1);
 365     if (x - i >= 0.5) {
 366         return (i + 1);
 367     }
 368     return (i);
 369 }
 370 
 371 define round(x, d) {
 372     auto k;
 373     k = 10^d;
 374     return (round0(x * k) / k);
 375 }
 376 
 377 define r(x, d) {
 378     return (round(x, d));
 379 }
 380 
 381 define r0(x) {
 382     return (round0(x));
 383 }
 384 
 385 define r1(x) {
 386     return (round(x, 1));
 387 }
 388 
 389 define r2(x) {
 390     return (round(x, 2));
 391 }
 392 
 393 define r3(x) {
 394     return (round(x, 3));
 395 }
 396 
 397 define r4(x) {
 398     return (round(x, 4));
 399 }
 400 
 401 define r5(x) {
 402     return (round(x, 5));
 403 }
 404 
 405 define r6(x) {
 406     return (round(x, 6));
 407 }
 408 
 409 define r7(x) {
 410     return (round(x, 7));
 411 }
 412 
 413 define r8(x) {
 414     return (round(x, 8));
 415 }
 416 
 417 define r9(x) {
 418     return (round(x, 9));
 419 }
 420 
 421 define r10(x) {
 422     return (round(x, 10));
 423 }
 424 
 425 define r11(x) {
 426     return (round(x, 11));
 427 }
 428 
 429 define r12(x) {
 430     return (round(x, 12));
 431 }
 432 
 433 define r13(x) {
 434     return (round(x, 13));
 435 }
 436 
 437 define r14(x) {
 438     return (round(x, 14));
 439 }
 440 
 441 define r15(x) {
 442     return (round(x, 15));
 443 }
 444 
 445 define r16(x) {
 446     return (round(x, 16));
 447 }
 448 
 449 define r17(x) {
 450     return (round(x, 17));
 451 }
 452 
 453 define r18(x) {
 454     return (round(x, 18));
 455 }
 456 
 457 define r19(x) {
 458     return (round(x, 19));
 459 }
 460 
 461 define r20(x) {
 462     return (round(x, 20));
 463 }
 464 
 465 define round1(x) {
 466     return (round(x, 1));
 467 }
 468 
 469 define round2(x) {
 470     return (round(x, 2));
 471 }
 472 
 473 define round3(x) {
 474     return (round(x, 3));
 475 }
 476 
 477 define round4(x) {
 478     return (round(x, 4));
 479 }
 480 
 481 define round5(x) {
 482     return (round(x, 5));
 483 }
 484 
 485 define round6(x) {
 486     return (round(x, 6));
 487 }
 488 
 489 define round7(x) {
 490     return (round(x, 7));
 491 }
 492 
 493 define round8(x) {
 494     return (round(x, 8));
 495 }
 496 
 497 define round9(x) {
 498     return (round(x, 9));
 499 }
 500 
 501 define round10(x) {
 502     return (round(x, 10));
 503 }
 504 
 505 define round11(x) {
 506     return (round(x, 11));
 507 }
 508 
 509 define round12(x) {
 510     return (round(x, 12));
 511 }
 512 
 513 define round13(x) {
 514     return (round(x, 13));
 515 }
 516 
 517 define round14(x) {
 518     return (round(x, 14));
 519 }
 520 
 521 define round15(x) {
 522     return (round(x, 15));
 523 }
 524 
 525 define round16(x) {
 526     return (round(x, 16));
 527 }
 528 
 529 define round17(x) {
 530     return (round(x, 17));
 531 }
 532 
 533 define round18(x) {
 534     return (round(x, 18));
 535 }
 536 
 537 define round19(x) {
 538     return (round(x, 19));
 539 }
 540 
 541 define round20(x) {
 542     return (round(x, 20));
 543 }
 544 
 545 define fac(x) {
 546     auto f, i;
 547     if (x < 0) return (0);
 548     f = 1;
 549     for (i = x; i >= 2; i--) { f *= i; }
 550     return (f);
 551 }
 552 
 553 define f(x) {
 554     return (fac(x));
 555 }
 556 
 557 define fact(x) {
 558     return (fac(x));
 559 }
 560 
 561 define factorial(x) {
 562     return (fac(x));
 563 }
 564 
 565 define per(n, k) {
 566     auto p, i;
 567     if (n < k) return (0);
 568     p = 1;
 569     for (i = n; i >= n - k + 1; i--) { p *= i; }
 570     return (p);
 571 }
 572 
 573 define p(n, k) {
 574     return (per(n, k));
 575 }
 576 
 577 define perm(n, k) {
 578     return (per(n, k));
 579 }
 580 
 581 define com(n, k) {
 582     if (n < k) return (0);
 583     return (per(n, k) / fac(k));
 584 }
 585 
 586 define comb(n, k) {
 587     return (com(n, k));
 588 }
 589 
 590 define gcd(x, y) {
 591     return (x * y / lcm(x, y));
 592 }
 593 
 594 define lcm(x, y) {
 595     auto a, b, z;
 596 
 597     /* the LCM is defined only for positive integers */
 598     # if (mod(x, 1) != 0 || x < 1 || mod(y, 1) != 0 || y < 1) { return 0; }
 599     if (mod(x, 1) != 0) return (0);
 600     if (x < 1) return (0);
 601     if (mod(y, 1) != 0) return (0);
 602     if (y < 1) return (0);
 603 
 604     a = min(x, y);
 605     b = max(x, y);
 606 
 607     z = b;
 608     while (mod(z, a) != 0) { z += b; }
 609     return (z);
 610 }
 611 
 612 define sgn(x) {
 613     if (x > 0) return (1);
 614     if (x < 0) return (-1);
 615     return (0);
 616 }
 617 
 618 define sign(x) {
 619     if (x > 0) return (1);
 620     if (x < 0) return (-1);
 621     return (0);
 622 }
 623 
 624 define isprime(n) {
 625     auto div;
 626 
 627     # if (mod(n, 1) != 0 || n < 2) { return 0; }
 628     if (mod(n, 1) != 0) return (0);
 629     if (n < 2) return (0);
 630 
 631     /* 2 is the only even prime number */
 632     if (mod(n, 2) == 0) {
 633         if (n == 2) return (1);
 634         return (0);
 635     }
 636 
 637     /*
 638         even numbers have already been handled, so only odd numbers
 639         make it here: only check up to the square-root, so the loop
 640         has O(n**1.5) time-complexity
 641     */
 642     for (div = 3; div <= sqrt(n); div += 2) {
 643         if (mod(n, div) == 0) return (0);
 644     }
 645     return (1);
 646 }
 647 
 648 /*
 649 I adapted an implementation of Lanczos' log-gamma approximation I found at
 650 http://introcs.cs.princeton.edu/java/91float/Gamma.java.html
 651 */
 652 define lgamma(x) {
 653     auto tmp, a, b, c, d;
 654     tmp = (x - 0.5) * l(x + 4.5) - (x + 4.5);
 655     a = 76.18009173 / (x + 0) - 86.50532033 / (x + 1);
 656     b = 24.01409822 / (x + 2) - 1.231739516 / (x + 3);
 657     c = 0.00120858003 / (x + 4) - 0.00000536382 / (x + 5);
 658     d = sqrt(2 * 3.14159265358979323844);
 659     return (tmp + l((1.0 + a + b + c) * d));
 660 }
 661 
 662 define lbeta(x, y) {
 663     return (lgamma(x) + lgamma(y) - lgamma(x + y));
 664 }
 665 
 666 define gamma(x) {
 667     return (e(lgamma(x)));
 668 }
 669 
 670 define beta(x, y) {
 671     return (e(lbeta(x, y)));
 672 }
 673 
 674 define logistic(x) {
 675     return (1 / (1 + e(-x)));
 676 }
 677 
 678 define mix(x, y, k) {
 679     # return (x + (1 - k) * (y - x));
 680     return (x + k * y);
 681 }
 682 
 683 define clamp(x, min, max) {
 684     return (max(min(max, x), min));
 685 }
 686 
 687 define wrap(x, min, max) {
 688     return ((x - min) / (max - min));
 689 }
 690 
 691 define unwrap(x, min, max) {
 692     return ((max - min) * x + min);
 693 }
 694 
 695 $(echo "$@" | sed 's-^+--g; s-_--g; s-\*\*-^-g; s-\[-(-g; s-\]-)-g')
 696 ENDOFSCRIPT
 697 } |
 698 # ensure the result shows at least a zero before the decimal dot, then
 699 # rid the result of trailing zero decimals and/or trailing decimal dots
 700 sed -E 's-^\.-0.-; s/^-\./-0./; s-(\.[0-9]+[1-9]+)0+$-\1-; s-\.0*$--'