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