1 /* 2 Copyright 2008-2023 3 Matthias Ehmann, 4 Michael Gerhaeuser, 5 Carsten Miller, 6 Bianca Valentin, 7 Alfred Wassermann, 8 Peter Wilfahrt 9 10 This file is part of JSXGraph. 11 12 JSXGraph is free software dual licensed under the GNU LGPL or MIT License. 13 14 You can redistribute it and/or modify it under the terms of the 15 16 * GNU Lesser General Public License as published by 17 the Free Software Foundation, either version 3 of the License, or 18 (at your option) any later version 19 OR 20 * MIT License: https://github.com/jsxgraph/jsxgraph/blob/master/LICENSE.MIT 21 22 JSXGraph is distributed in the hope that it will be useful, 23 but WITHOUT ANY WARRANTY; without even the implied warranty of 24 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 25 GNU Lesser General Public License for more details. 26 27 You should have received a copy of the GNU Lesser General Public License and 28 the MIT License along with JSXGraph. If not, see <https://www.gnu.org/licenses/> 29 and <https://opensource.org/licenses/MIT/>. 30 */ 31 32 /*global JXG: true, define: true, Float32Array: true */ 33 /*jslint nomen: true, plusplus: true, bitwise: true*/ 34 35 /** 36 * @fileoverview In this file the namespace JXG.Math is defined, which is the base namespace 37 * for namespaces like JXG.Math.Numerics, JXG.Math.Plot, JXG.Math.Statistics, JXG.Math.Clip etc. 38 */ 39 import JXG from "../jxg"; 40 import Type from "../utils/type"; 41 42 var undef, 43 /* 44 * Dynamic programming approach for recursive functions. 45 * From "Speed up your JavaScript, Part 3" by Nicholas C. Zakas. 46 * @see JXG.Math.factorial 47 * @see JXG.Math.binomial 48 * http://blog.thejit.org/2008/09/05/memoization-in-javascript/ 49 * 50 * This method is hidden, because it is only used in JXG.Math. If someone wants 51 * to use it in JSXGraph outside of JXG.Math, it should be moved to jsxgraph.js 52 */ 53 memoizer = function (f) { 54 var cache, join; 55 56 if (f.memo) { 57 return f.memo; 58 } 59 60 cache = {}; 61 join = Array.prototype.join; 62 63 f.memo = function () { 64 var key = join.call(arguments); 65 66 // Seems to be a bit faster than "if (a in b)" 67 return cache[key] !== undef ? cache[key] : (cache[key] = f.apply(this, arguments)); 68 }; 69 70 return f.memo; 71 }; 72 73 /** 74 * Math namespace. Contains mathematics related methods which are 75 * specific to JSXGraph or which extend the JavaScript Math class. 76 * @namespace 77 */ 78 JXG.Math = { 79 /** 80 * eps defines the closeness to zero. If the absolute value of a given number is smaller 81 * than eps, it is considered to be equal to zero. 82 * @type Number 83 */ 84 eps: 0.000001, 85 86 /** 87 * Determine the relative difference between two numbers. 88 * @param {Number} a First number 89 * @param {Number} b Second number 90 * @returns {Number} Relative difference between a and b: |a-b| / max(|a|, |b|) 91 */ 92 relDif: function (a, b) { 93 var c = Math.abs(a), 94 d = Math.abs(b); 95 96 d = Math.max(c, d); 97 98 return d === 0.0 ? 0.0 : Math.abs(a - b) / d; 99 }, 100 101 /** 102 * The JavaScript implementation of the % operator returns the symmetric modulo. 103 * mod and "%" are both identical if a >= 0 and m >= 0 but the results differ if a or m < 0. 104 * @param {Number} a 105 * @param {Number} m 106 * @returns {Number} Mathematical modulo <tt>a mod m</tt> 107 */ 108 mod: function (a, m) { 109 return a - Math.floor(a / m) * m; 110 }, 111 112 /** 113 * Initializes a vector as an array with the coefficients set to the given value resp. zero. 114 * @param {Number} n Length of the vector 115 * @param {Number} [init=0] Initial value for each coefficient 116 * @returns {Array} A <tt>n</tt> times <tt>m</tt>-matrix represented by a 117 * two-dimensional array. The inner arrays hold the columns, the outer array holds the rows. 118 */ 119 vector: function (n, init) { 120 var r, i; 121 122 init = init || 0; 123 r = []; 124 125 for (i = 0; i < n; i++) { 126 r[i] = init; 127 } 128 129 return r; 130 }, 131 132 /** 133 * Initializes a matrix as an array of rows with the given value. 134 * @param {Number} n Number of rows 135 * @param {Number} [m=n] Number of columns 136 * @param {Number} [init=0] Initial value for each coefficient 137 * @returns {Array} A <tt>n</tt> times <tt>m</tt>-matrix represented by a 138 * two-dimensional array. The inner arrays hold the columns, the outer array holds the rows. 139 */ 140 matrix: function (n, m, init) { 141 var r, i, j; 142 143 init = init || 0; 144 m = m || n; 145 r = []; 146 147 for (i = 0; i < n; i++) { 148 r[i] = []; 149 150 for (j = 0; j < m; j++) { 151 r[i][j] = init; 152 } 153 } 154 155 return r; 156 }, 157 158 /** 159 * Generates an identity matrix. If n is a number and m is undefined or not a number, a square matrix is generated, 160 * if n and m are both numbers, an nxm matrix is generated. 161 * @param {Number} n Number of rows 162 * @param {Number} [m=n] Number of columns 163 * @returns {Array} A square matrix of length <tt>n</tt> with all coefficients equal to 0 except a_(i,i), i out of (1, ..., n), if <tt>m</tt> is undefined or not a number 164 * or a <tt>n</tt> times <tt>m</tt>-matrix with a_(i,j) = 0 and a_(i,i) = 1 if m is a number. 165 */ 166 identity: function (n, m) { 167 var r, i; 168 169 if (m === undef && typeof m !== "number") { 170 m = n; 171 } 172 173 r = this.matrix(n, m); 174 175 for (i = 0; i < Math.min(n, m); i++) { 176 r[i][i] = 1; 177 } 178 179 return r; 180 }, 181 182 /** 183 * Generates a 4x4 matrix for 3D to 2D projections. 184 * @param {Number} l Left 185 * @param {Number} r Right 186 * @param {Number} t Top 187 * @param {Number} b Bottom 188 * @param {Number} n Near 189 * @param {Number} f Far 190 * @returns {Array} 4x4 Matrix 191 */ 192 frustum: function (l, r, b, t, n, f) { 193 var ret = this.matrix(4, 4); 194 195 ret[0][0] = (n * 2) / (r - l); 196 ret[0][1] = 0; 197 ret[0][2] = (r + l) / (r - l); 198 ret[0][3] = 0; 199 200 ret[1][0] = 0; 201 ret[1][1] = (n * 2) / (t - b); 202 ret[1][2] = (t + b) / (t - b); 203 ret[1][3] = 0; 204 205 ret[2][0] = 0; 206 ret[2][1] = 0; 207 ret[2][2] = -(f + n) / (f - n); 208 ret[2][3] = -(f * n * 2) / (f - n); 209 210 ret[3][0] = 0; 211 ret[3][1] = 0; 212 ret[3][2] = -1; 213 ret[3][3] = 0; 214 215 return ret; 216 }, 217 218 /** 219 * Generates a 4x4 matrix for 3D to 2D projections. 220 * @param {Number} fov Field of view in vertical direction, given in rad. 221 * @param {Number} ratio Aspect ratio of the projection plane. 222 * @param {Number} n Near 223 * @param {Number} f Far 224 * @returns {Array} 4x4 Projection Matrix 225 */ 226 projection: function (fov, ratio, n, f) { 227 var t = n * Math.tan(fov / 2), 228 r = t * ratio; 229 230 return this.frustum(-r, r, -t, t, n, f); 231 }, 232 233 /** 234 * Multiplies a vector vec to a matrix mat: mat * vec. The matrix is interpreted by this function as an array of rows. Please note: This 235 * function does not check if the dimensions match. 236 * @param {Array} mat Two dimensional array of numbers. The inner arrays describe the columns, the outer ones the matrix' rows. 237 * @param {Array} vec Array of numbers 238 * @returns {Array} Array of numbers containing the result 239 * @example 240 * var A = [[2, 1], 241 * [1, 3]], 242 * b = [4, 5], 243 * c; 244 * c = JXG.Math.matVecMult(A, b) 245 * // c === [13, 19]; 246 */ 247 matVecMult: function (mat, vec) { 248 var i, 249 s, 250 k, 251 m = mat.length, 252 n = vec.length, 253 res = []; 254 255 if (n === 3) { 256 for (i = 0; i < m; i++) { 257 res[i] = mat[i][0] * vec[0] + mat[i][1] * vec[1] + mat[i][2] * vec[2]; 258 } 259 } else { 260 for (i = 0; i < m; i++) { 261 s = 0; 262 for (k = 0; k < n; k++) { 263 s += mat[i][k] * vec[k]; 264 } 265 res[i] = s; 266 } 267 } 268 return res; 269 }, 270 271 /** 272 * Computes the product of the two matrices mat1*mat2. 273 * @param {Array} mat1 Two dimensional array of numbers 274 * @param {Array} mat2 Two dimensional array of numbers 275 * @returns {Array} Two dimensional Array of numbers containing result 276 */ 277 matMatMult: function (mat1, mat2) { 278 var i, 279 j, 280 s, 281 k, 282 m = mat1.length, 283 n = m > 0 ? mat2[0].length : 0, 284 m2 = mat2.length, 285 res = this.matrix(m, n); 286 287 for (i = 0; i < m; i++) { 288 for (j = 0; j < n; j++) { 289 s = 0; 290 for (k = 0; k < m2; k++) { 291 s += mat1[i][k] * mat2[k][j]; 292 } 293 res[i][j] = s; 294 } 295 } 296 return res; 297 }, 298 299 /** 300 * Transposes a matrix given as a two dimensional array. 301 * @param {Array} M The matrix to be transposed 302 * @returns {Array} The transpose of M 303 */ 304 transpose: function (M) { 305 var MT, i, j, m, n; 306 307 // number of rows of M 308 m = M.length; 309 // number of columns of M 310 n = M.length > 0 ? M[0].length : 0; 311 MT = this.matrix(n, m); 312 313 for (i = 0; i < n; i++) { 314 for (j = 0; j < m; j++) { 315 MT[i][j] = M[j][i]; 316 } 317 } 318 319 return MT; 320 }, 321 322 /** 323 * Compute the inverse of an nxn matrix with Gauss elimination. 324 * @param {Array} Ain 325 * @returns {Array} Inverse matrix of Ain 326 */ 327 inverse: function (Ain) { 328 var i, 329 j, 330 k, 331 s, 332 ma, 333 r, 334 swp, 335 n = Ain.length, 336 A = [], 337 p = [], 338 hv = []; 339 340 for (i = 0; i < n; i++) { 341 A[i] = []; 342 for (j = 0; j < n; j++) { 343 A[i][j] = Ain[i][j]; 344 } 345 p[i] = i; 346 } 347 348 for (j = 0; j < n; j++) { 349 // pivot search: 350 ma = Math.abs(A[j][j]); 351 r = j; 352 353 for (i = j + 1; i < n; i++) { 354 if (Math.abs(A[i][j]) > ma) { 355 ma = Math.abs(A[i][j]); 356 r = i; 357 } 358 } 359 360 // Singular matrix 361 if (ma <= this.eps) { 362 return []; 363 } 364 365 // swap rows: 366 if (r > j) { 367 for (k = 0; k < n; k++) { 368 swp = A[j][k]; 369 A[j][k] = A[r][k]; 370 A[r][k] = swp; 371 } 372 373 swp = p[j]; 374 p[j] = p[r]; 375 p[r] = swp; 376 } 377 378 // transformation: 379 s = 1.0 / A[j][j]; 380 for (i = 0; i < n; i++) { 381 A[i][j] *= s; 382 } 383 A[j][j] = s; 384 385 for (k = 0; k < n; k++) { 386 if (k !== j) { 387 for (i = 0; i < n; i++) { 388 if (i !== j) { 389 A[i][k] -= A[i][j] * A[j][k]; 390 } 391 } 392 A[j][k] = -s * A[j][k]; 393 } 394 } 395 } 396 397 // swap columns: 398 for (i = 0; i < n; i++) { 399 for (k = 0; k < n; k++) { 400 hv[p[k]] = A[i][k]; 401 } 402 for (k = 0; k < n; k++) { 403 A[i][k] = hv[k]; 404 } 405 } 406 407 return A; 408 }, 409 410 /** 411 * Inner product of two vectors a and b. n is the length of the vectors. 412 * @param {Array} a Vector 413 * @param {Array} b Vector 414 * @param {Number} [n] Length of the Vectors. If not given the length of the first vector is taken. 415 * @returns {Number} The inner product of a and b. 416 */ 417 innerProduct: function (a, b, n) { 418 var i, 419 s = 0; 420 421 if (n === undef || !Type.isNumber(n)) { 422 n = a.length; 423 } 424 425 for (i = 0; i < n; i++) { 426 s += a[i] * b[i]; 427 } 428 429 return s; 430 }, 431 432 /** 433 * Calculates the cross product of two vectors both of length three. 434 * In case of homogeneous coordinates this is either 435 * <ul> 436 * <li>the intersection of two lines</li> 437 * <li>the line through two points</li> 438 * </ul> 439 * @param {Array} c1 Homogeneous coordinates of line or point 1 440 * @param {Array} c2 Homogeneous coordinates of line or point 2 441 * @returns {Array} vector of length 3: homogeneous coordinates of the resulting point / line. 442 */ 443 crossProduct: function (c1, c2) { 444 return [ 445 c1[1] * c2[2] - c1[2] * c2[1], 446 c1[2] * c2[0] - c1[0] * c2[2], 447 c1[0] * c2[1] - c1[1] * c2[0] 448 ]; 449 }, 450 451 /** 452 * Euclidean norm of a vector. 453 * 454 * @param {Array} a Array containing a vector. 455 * @param {Number} n (Optional) length of the array. 456 * @returns {Number} Euclidean norm of the vector. 457 */ 458 norm: function (a, n) { 459 var i, 460 sum = 0.0; 461 462 if (n === undef || !Type.isNumber(n)) { 463 n = a.length; 464 } 465 466 for (i = 0; i < n; i++) { 467 sum += a[i] * a[i]; 468 } 469 470 return Math.sqrt(sum); 471 }, 472 473 /** 474 * Compute a * x + y for a scalar a and vectors x and y. 475 * 476 * @param {Number} a 477 * @param {Array} x 478 * @param {Array} y 479 * @returns 480 */ 481 axpy: function (a, x, y) { 482 var i, 483 le = x.length, 484 p = []; 485 for (i = 0; i < le; i++) { 486 p[i] = a * x[i] + y[i]; 487 } 488 return p; 489 }, 490 491 /** 492 * Compute the factorial of a positive integer. If a non-integer value 493 * is given, the fraction will be ignored. 494 * @function 495 * @param {Number} n 496 * @returns {Number} n! = n*(n-1)*...*2*1 497 */ 498 factorial: memoizer(function (n) { 499 if (n < 0) { 500 return NaN; 501 } 502 503 n = Math.floor(n); 504 505 if (n === 0 || n === 1) { 506 return 1; 507 } 508 509 return n * this.factorial(n - 1); 510 }), 511 512 /** 513 * Computes the binomial coefficient n over k. 514 * @function 515 * @param {Number} n Fraction will be ignored 516 * @param {Number} k Fraction will be ignored 517 * @returns {Number} The binomial coefficient n over k 518 */ 519 binomial: memoizer(function (n, k) { 520 var b, i; 521 522 if (k > n || k < 0) { 523 return NaN; 524 } 525 526 k = Math.round(k); 527 n = Math.round(n); 528 529 if (k === 0 || k === n) { 530 return 1; 531 } 532 533 b = 1; 534 535 for (i = 0; i < k; i++) { 536 b *= n - i; 537 b /= i + 1; 538 } 539 540 return b; 541 }), 542 543 /** 544 * Calculates the cosine hyperbolicus of x. 545 * @function 546 * @param {Number} x The number the cosine hyperbolicus will be calculated of. 547 * @returns {Number} Cosine hyperbolicus of the given value. 548 */ 549 cosh: 550 Math.cosh || 551 function (x) { 552 return (Math.exp(x) + Math.exp(-x)) * 0.5; 553 }, 554 555 /** 556 * Sine hyperbolicus of x. 557 * @function 558 * @param {Number} x The number the sine hyperbolicus will be calculated of. 559 * @returns {Number} Sine hyperbolicus of the given value. 560 */ 561 sinh: 562 Math.sinh || 563 function (x) { 564 return (Math.exp(x) - Math.exp(-x)) * 0.5; 565 }, 566 567 /** 568 * Hyperbolic arc-cosine of a number. 569 * 570 * @param {Number} x 571 * @returns {Number} 572 */ 573 acosh: 574 Math.acosh || 575 function (x) { 576 return Math.log(x + Math.sqrt(x * x - 1)); 577 }, 578 579 /** 580 * Hyperbolic arcsine of a number 581 * @param {Number} x 582 * @returns {Number} 583 */ 584 asinh: 585 Math.asinh || 586 function (x) { 587 if (x === -Infinity) { 588 return x; 589 } 590 return Math.log(x + Math.sqrt(x * x + 1)); 591 }, 592 593 /** 594 * Computes the cotangent of x. 595 * @function 596 * @param {Number} x The number the cotangent will be calculated of. 597 * @returns {Number} Cotangent of the given value. 598 */ 599 cot: function (x) { 600 return 1 / Math.tan(x); 601 }, 602 603 /** 604 * Computes the inverse cotangent of x. 605 * @param {Number} x The number the inverse cotangent will be calculated of. 606 * @returns {Number} Inverse cotangent of the given value. 607 */ 608 acot: function (x) { 609 return (x >= 0 ? 0.5 : -0.5) * Math.PI - Math.atan(x); 610 }, 611 612 /** 613 * Compute n-th real root of a real number. n must be strictly positive integer. 614 * If n is odd, the real n-th root exists and is negative. 615 * For n even, for negative valuees of x NaN is returned 616 * @param {Number} x radicand. Must be non-negative, if n even. 617 * @param {Number} n index of the root. must be strictly positive integer. 618 * @returns {Number} returns real root or NaN 619 * 620 * @example 621 * nthroot(16, 4): 2 622 * nthroot(-27, 3): -3 623 * nthroot(-4, 2): NaN 624 */ 625 nthroot: function (x, n) { 626 var inv = 1 / n; 627 628 if (n <= 0 || Math.floor(n) !== n) { 629 return NaN; 630 } 631 632 if (x === 0.0) { 633 return 0.0; 634 } 635 636 if (x > 0) { 637 return Math.exp(inv * Math.log(x)); 638 } 639 640 // From here on, x is negative 641 if (n % 2 === 1) { 642 return -Math.exp(inv * Math.log(-x)); 643 } 644 645 // x negative, even root 646 return NaN; 647 }, 648 649 /** 650 * Computes cube root of real number 651 * Polyfill for Math.cbrt(). 652 * 653 * @function 654 * @param {Number} x Radicand 655 * @returns {Number} Cube root of x. 656 */ 657 cbrt: 658 Math.cbrt || 659 function (x) { 660 return this.nthroot(x, 3); 661 }, 662 663 /** 664 * Compute base to the power of exponent. 665 * @param {Number} base 666 * @param {Number} exponent 667 * @returns {Number} base to the power of exponent. 668 */ 669 pow: function (base, exponent) { 670 if (base === 0) { 671 if (exponent === 0) { 672 return 1; 673 } 674 return 0; 675 } 676 677 // exponent is an integer 678 if (Math.floor(exponent) === exponent) { 679 return Math.pow(base, exponent); 680 } 681 682 // exponent is not an integer 683 if (base > 0) { 684 return Math.exp(exponent * Math.log(base)); 685 } 686 687 return NaN; 688 }, 689 690 /** 691 * Compute base to the power of the rational exponent m / n. 692 * This function first reduces the fraction m/n and then computes 693 * JXG.Math.pow(base, m/n). 694 * 695 * This function is necessary to have the same results for e.g. 696 * (-8)^(1/3) = (-8)^(2/6) = -2 697 * @param {Number} base 698 * @param {Number} m numerator of exponent 699 * @param {Number} n denominator of exponent 700 * @returns {Number} base to the power of exponent. 701 */ 702 ratpow: function (base, m, n) { 703 var g; 704 if (m === 0) { 705 return 1; 706 } 707 if (n === 0) { 708 return NaN; 709 } 710 711 g = this.gcd(m, n); 712 return this.nthroot(this.pow(base, m / g), n / g); 713 }, 714 715 /** 716 * Logarithm to base 10. 717 * @param {Number} x 718 * @returns {Number} log10(x) Logarithm of x to base 10. 719 */ 720 log10: function (x) { 721 return Math.log(x) / Math.log(10.0); 722 }, 723 724 /** 725 * Logarithm to base 2. 726 * @param {Number} x 727 * @returns {Number} log2(x) Logarithm of x to base 2. 728 */ 729 log2: function (x) { 730 return Math.log(x) / Math.log(2.0); 731 }, 732 733 /** 734 * Logarithm to arbitrary base b. If b is not given, natural log is taken, i.e. b = e. 735 * @param {Number} x 736 * @param {Number} b base 737 * @returns {Number} log(x, b) Logarithm of x to base b, that is log(x)/log(b). 738 */ 739 log: function (x, b) { 740 if (b !== undefined && Type.isNumber(b)) { 741 return Math.log(x) / Math.log(b); 742 } 743 744 return Math.log(x); 745 }, 746 747 /** 748 * The sign() function returns the sign of a number, indicating whether the number is positive, negative or zero. 749 * 750 * @function 751 * @param {Number} x A Number 752 * @returns {Number} This function has 5 kinds of return values, 753 * 1, -1, 0, -0, NaN, which represent "positive number", "negative number", "positive zero", "negative zero" 754 * and NaN respectively. 755 */ 756 sign: 757 Math.sign || 758 function (x) { 759 x = +x; // convert to a number 760 if (x === 0 || isNaN(x)) { 761 return x; 762 } 763 return x > 0 ? 1 : -1; 764 }, 765 766 /** 767 * A square & multiply algorithm to compute base to the power of exponent. 768 * Implementated by Wolfgang Riedl. 769 * 770 * @param {Number} base 771 * @param {Number} exponent 772 * @returns {Number} Base to the power of exponent 773 */ 774 squampow: function (base, exponent) { 775 var result; 776 777 if (Math.floor(exponent) === exponent) { 778 // exponent is integer (could be zero) 779 result = 1; 780 781 if (exponent < 0) { 782 // invert: base 783 base = 1.0 / base; 784 exponent *= -1; 785 } 786 787 while (exponent !== 0) { 788 if (exponent & 1) { 789 result *= base; 790 } 791 792 exponent >>= 1; 793 base *= base; 794 } 795 return result; 796 } 797 798 return this.pow(base, exponent); 799 }, 800 801 /** 802 * Greatest common divisor (gcd) of two numbers. 803 * @see <a href="https://rosettacode.org/wiki/Greatest_common_divisor#JavaScript">rosettacode.org</a> 804 * 805 * @param {Number} a First number 806 * @param {Number} b Second number 807 * @returns {Number} gcd(a, b) if a and b are numbers, NaN else. 808 */ 809 gcd: function (a, b) { 810 var tmp, 811 endless = true; 812 813 a = Math.abs(a); 814 b = Math.abs(b); 815 816 if (!(Type.isNumber(a) && Type.isNumber(b))) { 817 return NaN; 818 } 819 if (b > a) { 820 tmp = a; 821 a = b; 822 b = tmp; 823 } 824 825 while (endless) { 826 a %= b; 827 if (a === 0) { 828 return b; 829 } 830 b %= a; 831 if (b === 0) { 832 return a; 833 } 834 } 835 }, 836 837 /** 838 * Least common multiple (lcm) of two numbers. 839 * 840 * @param {Number} a First number 841 * @param {Number} b Second number 842 * @returns {Number} lcm(a, b) if a and b are numbers, NaN else. 843 */ 844 lcm: function (a, b) { 845 var ret; 846 847 if (!(Type.isNumber(a) && Type.isNumber(b))) { 848 return NaN; 849 } 850 851 ret = a * b; 852 if (ret !== 0) { 853 return ret / this.gcd(a, b); 854 } 855 856 return 0; 857 }, 858 859 /** 860 * Error function, see {@link https://en.wikipedia.org/wiki/Error_function}. 861 * 862 * @see JXG.Math.PropFunc.erf 863 * @param {Number} x 864 * @returns {Number} 865 */ 866 erf: function (x) { 867 return this.ProbFuncs.erf(x); 868 }, 869 870 /** 871 * Complementary error function, i.e. 1 - erf(x). 872 * 873 * @see JXG.Math.erf 874 * @see JXG.Math.PropFunc.erfc 875 * @param {Number} x 876 * @returns {Number} 877 */ 878 erfc: function (x) { 879 return this.ProbFuncs.erfc(x); 880 }, 881 882 /** 883 * Inverse of error function 884 * 885 * @see JXG.Math.erf 886 * @see JXG.Math.PropFunc.erfi 887 * @param {Number} x 888 * @returns {Number} 889 */ 890 erfi: function (x) { 891 return this.ProbFuncs.erfi(x); 892 }, 893 894 /** 895 * Normal distribution function 896 * 897 * @see JXG.Math.PropFunc.ndtr 898 * @param {Number} x 899 * @returns {Number} 900 */ 901 ndtr: function (x) { 902 return this.ProbFuncs.ndtr(x); 903 }, 904 905 /** 906 * Inverse of normal distribution function 907 * 908 * @see JXG.Math.ndtr 909 * @see JXG.Math.PropFunc.ndtri 910 * @param {Number} x 911 * @returns {Number} 912 */ 913 ndtri: function (x) { 914 return this.ProbFuncs.ndtri(x); 915 }, 916 917 /** 918 * Returns sqrt(a * a + b * b) for a variable number of arguments. 919 * This is a naive implementation which might be faster than Math.hypot. 920 * The latter is numerically more stable. 921 * 922 * @param {Number} a Variable number of arguments. 923 * @returns Number 924 */ 925 hypot: function() { 926 var i, le, a, sum; 927 928 le = arguments.length; 929 for (i = 0, sum = 0.0; i < le; i++) { 930 a = arguments[i]; 931 sum += a * a; 932 } 933 return Math.sqrt(sum); 934 }, 935 936 /** 937 * Heaviside unit step function. Returns 0 for x <, 1 for x > 0, and 0.5 for x == 0. 938 * 939 * @param {Number} x 940 * @returns Number 941 */ 942 hstep: function(x) { 943 return (x > 0.0) ? 1 : 944 ((x < 0.0) ? 0.0 : 0.5); 945 }, 946 947 /* ******************** Comparisons and logical operators ************** */ 948 949 /** 950 * Logical test: a < b? 951 * 952 * @param {Number} a 953 * @param {Number} b 954 * @returns {Boolean} 955 */ 956 lt: function (a, b) { 957 return a < b; 958 }, 959 960 /** 961 * Logical test: a <= b? 962 * 963 * @param {Number} a 964 * @param {Number} b 965 * @returns {Boolean} 966 */ 967 leq: function (a, b) { 968 return a <= b; 969 }, 970 971 /** 972 * Logical test: a > b? 973 * 974 * @param {Number} a 975 * @param {Number} b 976 * @returns {Boolean} 977 */ 978 gt: function (a, b) { 979 return a > b; 980 }, 981 982 /** 983 * Logical test: a >= b? 984 * 985 * @param {Number} a 986 * @param {Number} b 987 * @returns {Boolean} 988 */ 989 geq: function (a, b) { 990 return a >= b; 991 }, 992 993 /** 994 * Logical test: a === b? 995 * 996 * @param {Number} a 997 * @param {Number} b 998 * @returns {Boolean} 999 */ 1000 eq: function (a, b) { 1001 return a === b; 1002 }, 1003 1004 /** 1005 * Logical test: a !== b? 1006 * 1007 * @param {Number} a 1008 * @param {Number} b 1009 * @returns {Boolean} 1010 */ 1011 neq: function (a, b) { 1012 return a !== b; 1013 }, 1014 1015 /** 1016 * Logical operator: a && b? 1017 * 1018 * @param {Boolean} a 1019 * @param {Boolean} b 1020 * @returns {Boolean} 1021 */ 1022 and: function (a, b) { 1023 return a && b; 1024 }, 1025 1026 /** 1027 * Logical operator: !a? 1028 * 1029 * @param {Boolean} a 1030 * @returns {Boolean} 1031 */ 1032 not: function (a) { 1033 return !a; 1034 }, 1035 1036 /** 1037 * Logical operator: a || b? 1038 * 1039 * @param {Boolean} a 1040 * @param {Boolean} b 1041 * @returns {Boolean} 1042 */ 1043 or: function (a, b) { 1044 return a || b; 1045 }, 1046 1047 /** 1048 * Logical operator: either a or b? 1049 * 1050 * @param {Boolean} a 1051 * @param {Boolean} b 1052 * @returns {Boolean} 1053 */ 1054 xor: function (a, b) { 1055 return (a || b) && !(a && b); 1056 }, 1057 1058 /* *************************** Normalize *************************** */ 1059 1060 /** 1061 * Normalize the standard form [c, b0, b1, a, k, r, q0, q1]. 1062 * @private 1063 * @param {Array} stdform The standard form to be normalized. 1064 * @returns {Array} The normalized standard form. 1065 */ 1066 normalize: function (stdform) { 1067 var n, 1068 signr, 1069 a2 = 2 * stdform[3], 1070 r = stdform[4] / a2; 1071 1072 stdform[5] = r; 1073 stdform[6] = -stdform[1] / a2; 1074 stdform[7] = -stdform[2] / a2; 1075 1076 if (!isFinite(r)) { 1077 n = this.hypot(stdform[1], stdform[2]); 1078 1079 stdform[0] /= n; 1080 stdform[1] /= n; 1081 stdform[2] /= n; 1082 stdform[3] = 0; 1083 stdform[4] = 1; 1084 } else if (Math.abs(r) >= 1) { 1085 stdform[0] = (stdform[6] * stdform[6] + stdform[7] * stdform[7] - r * r) / (2 * r); 1086 stdform[1] = -stdform[6] / r; 1087 stdform[2] = -stdform[7] / r; 1088 stdform[3] = 1 / (2 * r); 1089 stdform[4] = 1; 1090 } else { 1091 signr = r <= 0 ? -1 : 1; 1092 stdform[0] = 1093 signr * (stdform[6] * stdform[6] + stdform[7] * stdform[7] - r * r) * 0.5; 1094 stdform[1] = -signr * stdform[6]; 1095 stdform[2] = -signr * stdform[7]; 1096 stdform[3] = signr / 2; 1097 stdform[4] = signr * r; 1098 } 1099 1100 return stdform; 1101 }, 1102 1103 /** 1104 * Converts a two dimensional array to a one dimensional Float32Array that can be processed by WebGL. 1105 * @param {Array} m A matrix in a two dimensional array. 1106 * @returns {Float32Array} A one dimensional array containing the matrix in column wise notation. Provides a fall 1107 * back to the default JavaScript Array if Float32Array is not available. 1108 */ 1109 toGL: function (m) { 1110 var v, i, j; 1111 1112 if (typeof Float32Array === "function") { 1113 v = new Float32Array(16); 1114 } else { 1115 v = new Array(16); 1116 } 1117 1118 if (m.length !== 4 && m[0].length !== 4) { 1119 return v; 1120 } 1121 1122 for (i = 0; i < 4; i++) { 1123 for (j = 0; j < 4; j++) { 1124 v[i + 4 * j] = m[i][j]; 1125 } 1126 } 1127 1128 return v; 1129 }, 1130 1131 /** 1132 * Theorem of Vieta: Given a set of simple zeroes x_0, ..., x_n 1133 * of a polynomial f, compute the coefficients s_k, (k=0,...,n-1) 1134 * of the polynomial of the form. See {@link https://de.wikipedia.org/wiki/Elementarsymmetrisches_Polynom}. 1135 * <p> 1136 * f(x) = (x-x_0)*...*(x-x_n) = 1137 * x^n + sum_{k=1}^{n} (-1)^(k) s_{k-1} x^(n-k) 1138 * </p> 1139 * @param {Array} x Simple zeroes of the polynomial. 1140 * @returns {Array} Coefficients of the polynomial. 1141 * 1142 */ 1143 Vieta: function (x) { 1144 var n = x.length, 1145 s = [], 1146 m, 1147 k, 1148 y; 1149 1150 s = x.slice(); 1151 for (m = 1; m < n; ++m) { 1152 y = s[m]; 1153 s[m] *= s[m - 1]; 1154 for (k = m - 1; k >= 1; --k) { 1155 s[k] += s[k - 1] * y; 1156 } 1157 s[0] += y; 1158 } 1159 return s; 1160 } 1161 }; 1162 1163 export default JXG.Math; 1164