1 /*
  2     Copyright 2008-2026
  3         Matthias Ehmann,
  4         Michael Gerhaeuser,
  5         Carsten Miller,
  6         Bianca Valentin,
  7         Andreas Walter,
  8         Alfred Wassermann,
  9         Peter Wilfahrt
 10 
 11     This file is part of JSXGraph.
 12 
 13     JSXGraph is free software dual licensed under the GNU LGPL or MIT License.
 14 
 15     You can redistribute it and/or modify it under the terms of the
 16 
 17       * GNU Lesser General Public License as published by
 18         the Free Software Foundation, either version 3 of the License, or
 19         (at your option) any later version
 20       OR
 21       * MIT License: https://github.com/jsxgraph/jsxgraph/blob/master/LICENSE.MIT
 22 
 23     JSXGraph is distributed in the hope that it will be useful,
 24     but WITHOUT ANY WARRANTY; without even the implied warranty of
 25     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 26     GNU Lesser General Public License for more details.
 27 
 28     You should have received a copy of the GNU Lesser General Public License and
 29     the MIT License along with JSXGraph. If not, see <https://www.gnu.org/licenses/>
 30     and <https://opensource.org/licenses/MIT/>.
 31  */
 32 
 33 /*global JXG: true, define: true*/
 34 /*jslint nomen: true, plusplus: true*/
 35 
 36 /**
 37  * @fileoverview This file contains the Math.Geometry namespace for calculating algebraic/geometric
 38  * stuff like intersection points, angles, midpoint, and so on.
 39  */
 40 
 41 import JXG from "../jxg.js";
 42 import Const from "../base/constants.js";
 43 import Coords from "../base/coords.js";
 44 import Mat from "./math.js";
 45 import Stat from "../math/statistics.js";
 46 import Numerics from "./numerics.js";
 47 import Type from "../utils/type.js";
 48 import Expect from "../utils/expect.js";
 49 
 50 /**
 51  * Math.Geometry namespace definition. This namespace holds geometrical algorithms,
 52  * especially intersection algorithms.
 53  * @name JXG.Math.Geometry
 54  * @exports Mat.Geometry as JXG.Math.Geometry
 55  * @namespace
 56  */
 57 Mat.Geometry = {};
 58 
 59 // the splitting is necessary due to the shortcut for the circumcircleMidpoint method to circumcenter.
 60 
 61 JXG.extend(
 62     Mat.Geometry,
 63     /** @lends JXG.Math.Geometry */ {
 64         /* ***************************************/
 65         /* *** GENERAL GEOMETRIC CALCULATIONS ****/
 66         /* ***************************************/
 67 
 68         /**
 69          * Calculates the angle defined by the points A, B, C.
 70          * @param {JXG.Point|Array} A A point  or [x,y] array.
 71          * @param {JXG.Point|Array} B Another point or [x,y] array.
 72          * @param {JXG.Point|Array} C A circle - no, of course the third point or [x,y] array.
 73          * @deprecated Use {@link JXG.Math.Geometry.rad} instead.
 74          * @see JXG.Math.Geometry.rad
 75          * @see JXG.Math.Geometry.trueAngle
 76          * @returns {Number} The angle in radian measure.
 77          */
 78         angle: function (A, B, C) {
 79             var u,
 80                 v,
 81                 s,
 82                 t,
 83                 a = [],
 84                 b = [],
 85                 c = [];
 86 
 87             JXG.deprecated("Geometry.angle()", "Geometry.rad()");
 88             if (A.coords) {
 89                 a[0] = A.coords.usrCoords[1];
 90                 a[1] = A.coords.usrCoords[2];
 91             } else {
 92                 a[0] = A[0];
 93                 a[1] = A[1];
 94             }
 95 
 96             if (B.coords) {
 97                 b[0] = B.coords.usrCoords[1];
 98                 b[1] = B.coords.usrCoords[2];
 99             } else {
100                 b[0] = B[0];
101                 b[1] = B[1];
102             }
103 
104             if (C.coords) {
105                 c[0] = C.coords.usrCoords[1];
106                 c[1] = C.coords.usrCoords[2];
107             } else {
108                 c[0] = C[0];
109                 c[1] = C[1];
110             }
111 
112             u = a[0] - b[0];
113             v = a[1] - b[1];
114             s = c[0] - b[0];
115             t = c[1] - b[1];
116 
117             return Math.atan2(u * t - v * s, u * s + v * t);
118         },
119 
120         /**
121          * Calculates the angle defined by the three points A, B, C if you're going from A to C around B counterclockwise.
122          * @param {JXG.Point|Array} A Point or [x,y] array
123          * @param {JXG.Point|Array} B Point or [x,y] array
124          * @param {JXG.Point|Array} C Point or [x,y] array
125          * @see JXG.Math.Geometry.rad
126          * @returns {Number} The angle in degrees.
127          */
128         trueAngle: function (A, B, C) {
129             return this.rad(A, B, C) * 57.295779513082323; // *180.0/Math.PI;
130         },
131 
132         /**
133          * Calculates the internal angle defined by the three points A, B, C if you're going from A to C around B counterclockwise.
134          * @param {JXG.Point|Array} A Point or [x,y] array
135          * @param {JXG.Point|Array} B Point or [x,y] array
136          * @param {JXG.Point|Array} C Point or [x,y] array
137          * @see JXG.Math.Geometry.trueAngle
138          * @returns {Number} Angle in radians.
139          */
140         rad: function (A, B, C) {
141             var ax, ay, bx, by, cx, cy, phi;
142 
143             if (A.coords) {
144                 ax = A.coords.usrCoords[1];
145                 ay = A.coords.usrCoords[2];
146             } else {
147                 ax = A[0];
148                 ay = A[1];
149             }
150 
151             if (B.coords) {
152                 bx = B.coords.usrCoords[1];
153                 by = B.coords.usrCoords[2];
154             } else {
155                 bx = B[0];
156                 by = B[1];
157             }
158 
159             if (C.coords) {
160                 cx = C.coords.usrCoords[1];
161                 cy = C.coords.usrCoords[2];
162             } else {
163                 cx = C[0];
164                 cy = C[1];
165             }
166 
167             phi = Math.atan2(cy - by, cx - bx) - Math.atan2(ay - by, ax - bx);
168 
169             if (phi < 0) {
170                 phi += 6.2831853071795862;
171             }
172 
173             return phi;
174         },
175 
176         /**
177          * Calculates a point on the bisection line between the three points A, B, C.
178          * As a result, the bisection line is defined by two points:
179          * Parameter B and the point with the coordinates calculated in this function.
180          * Does not work for ideal points.
181          * @param {JXG.Point} A Point
182          * @param {JXG.Point} B Point
183          * @param {JXG.Point} C Point
184          * @param [board=A.board] Reference to the board
185          * @returns {JXG.Coords} Coordinates of the second point defining the bisection.
186          */
187         angleBisector: function (A, B, C, board) {
188             var phiA,
189                 phiC,
190                 phi,
191                 Ac = A.coords.usrCoords,
192                 Bc = B.coords.usrCoords,
193                 Cc = C.coords.usrCoords,
194                 x,
195                 y;
196 
197             if (!Type.exists(board)) {
198                 board = A.board;
199             }
200 
201             // Parallel lines
202             if (Bc[0] === 0) {
203                 return new Coords(
204                     Const.COORDS_BY_USER,
205                     [1, (Ac[1] + Cc[1]) * 0.5, (Ac[2] + Cc[2]) * 0.5],
206                     board
207                 );
208             }
209 
210             // Non-parallel lines
211             x = Ac[1] - Bc[1];
212             y = Ac[2] - Bc[2];
213             phiA = Math.atan2(y, x);
214 
215             x = Cc[1] - Bc[1];
216             y = Cc[2] - Bc[2];
217             phiC = Math.atan2(y, x);
218 
219             phi = (phiA + phiC) * 0.5;
220 
221             if (phiA > phiC) {
222                 phi += Math.PI;
223             }
224 
225             x = Math.cos(phi) + Bc[1];
226             y = Math.sin(phi) + Bc[2];
227 
228             return new Coords(Const.COORDS_BY_USER, [1, x, y], board);
229         },
230 
231         // /**
232         //  * Calculates a point on the m-section line between the three points A, B, C.
233         //  * As a result, the m-section line is defined by two points:
234         //  * Parameter B and the point with the coordinates calculated in this function.
235         //  * The m-section generalizes the bisector to any real number.
236         //  * For example, the trisectors of an angle are simply the 1/3-sector and the 2/3-sector.
237         //  * Does not work for ideal points.
238         //  * @param {JXG.Point} A Point
239         //  * @param {JXG.Point} B Point
240         //  * @param {JXG.Point} C Point
241         //  * @param {Number} m Number
242         //  * @param [board=A.board] Reference to the board
243         //  * @returns {JXG.Coords} Coordinates of the second point defining the bisection.
244         //  */
245         // angleMsector: function (A, B, C, m, board) {
246         //     var phiA, phiC, phi,
247         //         Ac = A.coords.usrCoords,
248         //         Bc = B.coords.usrCoords,
249         //         Cc = C.coords.usrCoords,
250         //         x, y;
251 
252         //     if (!Type.exists(board)) {
253         //         board = A.board;
254         //     }
255 
256         //     // Parallel lines
257         //     if (Bc[0] === 0) {
258         //         return new Coords(Const.COORDS_BY_USER,
259         //             [1, (Ac[1] + Cc[1]) * m, (Ac[2] + Cc[2]) * m], board);
260         //     }
261 
262         //     // Non-parallel lines
263         //     x = Ac[1] - Bc[1];
264         //     y = Ac[2] - Bc[2];
265         //     phiA =  Math.atan2(y, x);
266 
267         //     x = Cc[1] - Bc[1];
268         //     y = Cc[2] - Bc[2];
269         //     phiC =  Math.atan2(y, x);
270 
271         //     phi = phiA + ((phiC - phiA) * m);
272 
273         //     if (phiA - phiC > Math.PI) {
274         //         phi += 2*m*Math.PI;
275         //     }
276 
277         //     x = Math.cos(phi) + Bc[1];
278         //     y = Math.sin(phi) + Bc[2];
279 
280         //     return new Coords(Const.COORDS_BY_USER, [1, x, y], board);
281         // },
282 
283         /**
284          * Reflects the point along the line.
285          * @param {JXG.Line} line Axis of reflection.
286          * @param {JXG.Point} point Point to reflect.
287          * @param [board=point.board] Reference to the board
288          * @returns {JXG.Coords} Coordinates of the reflected point.
289          */
290         reflection: function (line, point, board) {
291             // (v,w) defines the slope of the line
292             var x0,
293                 y0,
294                 x1,
295                 y1,
296                 v,
297                 w,
298                 mu,
299                 pc = point.coords.usrCoords,
300                 p1c = line.point1.coords.usrCoords,
301                 p2c = line.point2.coords.usrCoords;
302 
303             if (!Type.exists(board)) {
304                 board = point.board;
305             }
306 
307             v = p2c[1] - p1c[1];
308             w = p2c[2] - p1c[2];
309 
310             x0 = pc[1] - p1c[1];
311             y0 = pc[2] - p1c[2];
312 
313             mu = (v * y0 - w * x0) / (v * v + w * w);
314 
315             // point + mu*(-y,x) is the perpendicular foot
316             x1 = pc[1] + 2 * mu * w;
317             y1 = pc[2] - 2 * mu * v;
318 
319             return new Coords(Const.COORDS_BY_USER, [x1, y1], board);
320         },
321 
322         /**
323          * Computes the new position of a point which is rotated
324          * around a second point (called rotpoint) by the angle phi.
325          * @param {JXG.Point} rotpoint Center of the rotation
326          * @param {JXG.Point} point point to be rotated
327          * @param {Number} phi rotation angle in arc length
328          * @param {JXG.Board} [board=point.board] Reference to the board
329          * @returns {JXG.Coords} Coordinates of the new position.
330          */
331         rotation: function (rotpoint, point, phi, board) {
332             var x0,
333                 y0,
334                 c,
335                 s,
336                 x1,
337                 y1,
338                 pc = point.coords.usrCoords,
339                 rotpc = rotpoint.coords.usrCoords;
340 
341             if (!Type.exists(board)) {
342                 board = point.board;
343             }
344 
345             x0 = pc[1] - rotpc[1];
346             y0 = pc[2] - rotpc[2];
347 
348             c = Math.cos(phi);
349             s = Math.sin(phi);
350 
351             x1 = x0 * c - y0 * s + rotpc[1];
352             y1 = x0 * s + y0 * c + rotpc[2];
353 
354             return new Coords(Const.COORDS_BY_USER, [x1, y1], board);
355         },
356 
357         /**
358          * Calculates the coordinates of a point on the perpendicular to the given line through
359          * the given point.
360          * @param {JXG.Line} line A line.
361          * @param {JXG.Point} point Point which is projected to the line.
362          * @param {JXG.Board} [board=point.board] Reference to the board
363          * @returns {Array} Array of length two containing coordinates of a point on the perpendicular to the given line
364          *                  through the given point and boolean flag "change".
365          */
366         perpendicular: function (line, point, board) {
367             var x,
368                 y,
369                 change,
370                 c,
371                 z,
372                 A = line.point1.coords.usrCoords,
373                 B = line.point2.coords.usrCoords,
374                 C = point.coords.usrCoords;
375 
376             if (!Type.exists(board)) {
377                 board = point.board;
378             }
379 
380             // special case: point is the first point of the line
381             if (point === line.point1) {
382                 x = A[1] + B[2] - A[2];
383                 y = A[2] - B[1] + A[1];
384                 z = A[0] * B[0];
385 
386                 if (Math.abs(z) < Mat.eps) {
387                     x = B[2];
388                     y = -B[1];
389                 }
390                 c = [z, x, y];
391                 change = true;
392 
393                 // special case: point is the second point of the line
394             } else if (point === line.point2) {
395                 x = B[1] + A[2] - B[2];
396                 y = B[2] - A[1] + B[1];
397                 z = A[0] * B[0];
398 
399                 if (Math.abs(z) < Mat.eps) {
400                     x = A[2];
401                     y = -A[1];
402                 }
403                 c = [z, x, y];
404                 change = false;
405 
406                 // special case: point lies somewhere else on the line
407             } else if (Math.abs(Mat.innerProduct(C, line.stdform, 3)) < Mat.eps) {
408                 x = C[1] + B[2] - C[2];
409                 y = C[2] - B[1] + C[1];
410                 z = B[0];
411 
412                 if (Math.abs(z) < Mat.eps) {
413                     x = B[2];
414                     y = -B[1];
415                 }
416 
417                 change = true;
418                 if (
419                     Math.abs(z) > Mat.eps &&
420                     Math.abs(x - C[1]) < Mat.eps &&
421                     Math.abs(y - C[2]) < Mat.eps
422                 ) {
423                     x = C[1] + A[2] - C[2];
424                     y = C[2] - A[1] + C[1];
425                     change = false;
426                 }
427                 c = [z, x, y];
428 
429                 // general case: point does not lie on the line
430                 // -> calculate the foot of the dropped perpendicular
431             } else {
432                 c = [0, line.stdform[1], line.stdform[2]];
433                 c = Mat.crossProduct(c, C); // perpendicuar to line
434                 c = Mat.crossProduct(c, line.stdform); // intersection of line and perpendicular
435                 change = true;
436             }
437 
438             return [new Coords(Const.COORDS_BY_USER, c, board), change];
439         },
440 
441         /**
442          * @deprecated Please use {@link JXG.Math.Geometry.circumcenter} instead.
443          */
444         circumcenterMidpoint: function () {
445             JXG.deprecated("Geometry.circumcenterMidpoint()", "Geometry.circumcenter()");
446             this.circumcenter.apply(this, arguments);
447         },
448 
449         /**
450          * Calculates the center of the circumcircle of the three given points.
451          * @param {JXG.Point} point1 Point
452          * @param {JXG.Point} point2 Point
453          * @param {JXG.Point} point3 Point
454          * @param {JXG.Board} [board=point1.board] Reference to the board
455          * @returns {JXG.Coords} Coordinates of the center of the circumcircle of the given points.
456          */
457         circumcenter: function (point1, point2, point3, board) {
458             var u,
459                 v,
460                 m1,
461                 m2,
462                 A = point1.coords.usrCoords,
463                 B = point2.coords.usrCoords,
464                 C = point3.coords.usrCoords;
465 
466             if (!Type.exists(board)) {
467                 board = point1.board;
468             }
469 
470             u = [B[0] - A[0], -B[2] + A[2], B[1] - A[1]];
471             v = [(A[0] + B[0]) * 0.5, (A[1] + B[1]) * 0.5, (A[2] + B[2]) * 0.5];
472             m1 = Mat.crossProduct(u, v);
473 
474             u = [C[0] - B[0], -C[2] + B[2], C[1] - B[1]];
475             v = [(B[0] + C[0]) * 0.5, (B[1] + C[1]) * 0.5, (B[2] + C[2]) * 0.5];
476             m2 = Mat.crossProduct(u, v);
477 
478             return new Coords(Const.COORDS_BY_USER, Mat.crossProduct(m1, m2), board);
479         },
480 
481         /**
482          * Calculates the Euclidean distance for two given arrays of the same length.
483          * @param {Array} array1 Array of Number
484          * @param {Array} array2 Array of Number
485          * @param {Number} [n] Length of the arrays. Default is the minimum length of the given arrays.
486          * @returns {Number} Euclidean distance of the given vectors.
487          */
488         distance: function (array1, array2, n) {
489             var i,
490                 sum = 0;
491 
492             if (!n) {
493                 n = Math.min(array1.length, array2.length);
494             }
495 
496             for (i = 0; i < n; i++) {
497                 sum += (array1[i] - array2[i]) * (array1[i] - array2[i]);
498             }
499 
500             return Math.sqrt(sum);
501         },
502 
503         /**
504          * Calculates Euclidean distance for two given arrays of the same length.
505          * If one of the arrays contains a zero in the first coordinate, and the Euclidean distance
506          * is different from zero it is a point at infinity and we return Infinity.
507          * @param {Array} array1 Array containing elements of type number.
508          * @param {Array} array2 Array containing elements of type number.
509          * @param {Number} [n] Length of the arrays. Default is the minimum length of the given arrays.
510          * @returns {Number} Euclidean (affine) distance of the given vectors.
511          */
512         affineDistance: function (array1, array2, n) {
513             var d;
514 
515             d = this.distance(array1, array2, n);
516 
517             if (
518                 d > Mat.eps &&
519                 (Math.abs(array1[0]) < Mat.eps || Math.abs(array2[0]) < Mat.eps)
520             ) {
521                 return Infinity;
522             }
523 
524             return d;
525         },
526 
527         /**
528          * Affine ratio of three collinear points a, b, c: (c - a) / (b - a).
529          * If r > 1 or r < 0 then c is outside of the segment ab.
530          *
531          * @param {Array|JXG.Coords} a
532          * @param {Array|JXG.Coords} b
533          * @param {Array|JXG.Coords} c
534          * @returns {Number} affine ratio (c - a) / (b - a)
535          */
536         affineRatio: function (a, b, c) {
537             var r = 0.0,
538                 dx;
539 
540             if (Type.exists(a.usrCoords)) {
541                 a = a.usrCoords;
542             }
543             if (Type.exists(b.usrCoords)) {
544                 b = b.usrCoords;
545             }
546             if (Type.exists(c.usrCoords)) {
547                 c = c.usrCoords;
548             }
549 
550             dx = b[1] - a[1];
551 
552             if (Math.abs(dx) > Mat.eps) {
553                 r = (c[1] - a[1]) / dx;
554             } else {
555                 r = (c[2] - a[2]) / (b[2] - a[2]);
556             }
557             return r;
558         },
559 
560         /**
561          * Sort vertices counter clockwise starting with the first point.
562          * Used in Polygon.sutherlandHodgman, Geometry.signedPolygon.
563          *
564          * @param {Array} p An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
565          *
566          * @returns {Array}
567          */
568         sortVertices: function (p) {
569             var ll,
570                 ps = Expect.each(p, Expect.coordsArray),
571                 N = ps.length,
572                 lastPoint = null;
573 
574             // If the last point equals the first point, we take the last point out of the array.
575             // It may be that the several points at the end of the array are equal to the first point.
576             // The polygonal chain is been closed by JSXGraph, but this may also have been done by the user.
577             // Therefore, we use a while loop to pop the last points.
578             while (
579                 ps[0][0] === ps[N - 1][0] &&
580                 ps[0][1] === ps[N - 1][1] &&
581                 ps[0][2] === ps[N - 1][2]
582             ) {
583                 lastPoint = ps.pop();
584                 N--;
585             }
586 
587             ll = ps[0];
588             // Sort ps in increasing order of the angle between a point and the first point ll.
589             // If a point is equal to the first point ll, the angle is defined to be -Infinity.
590             // Otherwise, atan2 would return zero, which is a value which also attained by points
591             // on the same horizontal line.
592             ps.sort(function (a, b) {
593                 var rad1 =
594                         (a[2] === ll[2] && a[1] === ll[1])
595                             ? -Infinity
596                             : Math.atan2(a[2] - ll[2], a[1] - ll[1]),
597                     rad2 =
598                         (b[2] === ll[2] && b[1] === ll[1])
599                             ? -Infinity
600                             : Math.atan2(b[2] - ll[2], b[1] - ll[1]);
601                 return rad1 - rad2;
602             });
603 
604             // If the last point has been taken out of the array, we put it in again.
605             if (lastPoint !== null) {
606                 ps.push(lastPoint);
607             }
608 
609             return ps;
610         },
611 
612         /**
613          * Signed triangle area of the three points given. It can also be used
614          * to test the orientation of the triangle.
615          * <ul>
616          * <li> If the return value is < 0, then the point p2 is left of the line [p1, p3] (i.e p3 is right from [p1, p2]).
617          * <li> If the return value is > 0, then the point p2 is right of the line [p1, p3] (i.e p3 is left from [p1, p2]).
618          * <li> If the return value is = 0, then the points p1, p2, p3 are collinear.
619          * </ul>
620          *
621          * @param {JXG.Point|JXG.Coords|Array} p1
622          * @param {JXG.Point|JXG.Coords|Array} p2
623          * @param {JXG.Point|JXG.Coords|Array} p3
624          *
625          * @returns {Number}
626          */
627         signedTriangle: function (p1, p2, p3) {
628             var A = Expect.coordsArray(p1),
629                 B = Expect.coordsArray(p2),
630                 C = Expect.coordsArray(p3);
631             return 0.5 * ((B[1] - A[1]) * (C[2] - A[2]) - (B[2] - A[2]) * (C[1] - A[1]));
632         },
633 
634         /**
635          * Determine the signed area of a non-self-intersecting polygon.
636          * Surveyor's Formula
637          *
638          * @param {Array} p An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
639          * @param {Boolean} [sort=true]
640          *
641          * @returns {Number}
642          */
643         signedPolygon: function (p, sort) {
644             var i,
645                 N,
646                 A = 0,
647                 ps = Expect.each(p, Expect.coordsArray);
648 
649             if (sort === undefined) {
650                 sort = true;
651             }
652 
653             if (!sort) {
654                 ps = this.sortVertices(ps);
655             } else {
656                 // Make sure the polygon is closed. If it is already closed this won't change the sum because the last
657                 // summand will be 0.
658                 ps.unshift(ps[ps.length - 1]);
659             }
660 
661             N = ps.length;
662 
663             for (i = 1; i < N; i++) {
664                 A += ps[i - 1][1] * ps[i][2] - ps[i][1] * ps[i - 1][2];
665             }
666 
667             return 0.5 * A;
668         },
669 
670         /**
671          * Calculate the complex hull of a point cloud by the Graham scan algorithm.
672          *
673          * @param {Array} points An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
674          *
675          * @returns {Array} List of objects <pre>{i: index, c: coords}</pre> containing the convex hull points
676          *  in form of the index in the original input array and a coords array.
677          *
678          * @example
679          *     // Static example
680          *
681          *     var i, hull,
682          *       p = [],
683          *       q = [];
684          *
685          *     p.push( board.create('point', [4, 0], {withLabel:false }) );
686          *     p.push( board.create('point', [0, 4], {withLabel:false }) );
687          *     p.push( board.create('point', [0, 0], {withLabel:false }) );
688          *     p.push([-1, 0]);
689          *     p.push([-3, -3]);
690          *
691          *     hull = JXG.Math.Geometry.GrahamScan(p);
692          *     for (i = 0; i < hull.length; i++) {
693          *       console.log(hull[i]);
694          *       q.push(hull[i].c);
695          *     }
696          *     board.create('polygon', q);
697          *     // Output:
698          *     // { i: 4, c: [1, -3, 3]}
699          *     // { i: 0, c: [1, 4, 0]}
700          *     // { i: 1, c: [1, 0, 4]}
701          *
702          * </pre><div id="JXGb310b874-595e-4020-b0c2-566482797836" class="jxgbox" style="width: 300px; height: 300px;"></div>
703          * <script type="text/javascript">
704          *     (function() {
705          *         var board = JXG.JSXGraph.initBoard('JXGb310b874-595e-4020-b0c2-566482797836',
706          *             {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false});
707          *         var i, hull,
708          *           p = [],
709          *           q = [];
710          *
711          *         p.push( board.create('point', [4, 0], {withLabel:false }) );
712          *         p.push( board.create('point', [0, 4], {withLabel:false }) );
713          *         p.push( board.create('point', [0, 0], {withLabel:false }) );
714          *         p.push([-1, 0]);
715          *         p.push([-3, -3]);
716          *
717          *         hull = JXG.Math.Geometry.GrahamScan(p);
718          *         for (i = 0; i < hull.length; i++) {
719          *           console.log(hull[i]);
720          *           q.push(hull[i].c);
721          *         }
722          *         board.create('polygon', q);
723          *
724          *     })();
725          *
726          * </script><pre>
727          *
728          */
729         GrahamScan: function (points) {
730             var i, M, o,
731                 mi_idx,
732                 mi_x, mi_y, ma_x, ma_y,
733                 mi_xpy, mi_xmy, ma_xpy, ma_xmy,
734                 mi_x_i, ma_x_i, mi_y_i, ma_y_i,
735                 mi_xpy_i, mi_xmy_i, ma_xpy_i, ma_xmy_i,
736                 v, c,
737                 eps = Mat.eps * Mat.eps,
738                 that = this,
739                 ps_idx = [],
740                 stack = [],
741                 ps = Expect.each(points, Expect.coordsArray), // New array object, i.e. a copy of the input array.
742                 N,
743                 AklToussaint = 1024;  // This is a rough threshold where the heuristic pays off.
744 
745             N = ps.length;
746             if (N === 0) {
747                 return [];
748             }
749 
750             if (N > AklToussaint) {
751                 //
752                 // Akl-Toussaint heuristic
753                 // Determine an irregular convex octagon whose inside can be discarded.
754                 //
755                 mi_x = ps[0][1];
756                 ma_x = mi_x;
757                 mi_y = ps[0][2];
758                 ma_y = mi_y;
759 
760                 mi_xmy = ps[0][1] - ps[0][2];
761                 ma_xmy = mi_xmy;
762                 mi_xpy = ps[0][1] + ps[0][2];
763                 ma_xpy = mi_xpy;
764 
765                 mi_x_i = 0;
766                 ma_x_i = 0;
767                 mi_y_i = 0;
768                 ma_y_i = 0;
769 
770                 mi_xmy_i = 0;
771                 ma_xmy_i = 0;
772                 mi_xpy_i = 0;
773                 ma_xpy_i = 0;
774                 for (i = 1; i < N; i++) {
775                     v = ps[i][1];
776                     if (v < mi_x) {
777                         mi_x = v;
778                         mi_x_i = i;
779                     } else if (v > ma_x) {
780                         ma_x = v;
781                         ma_x_i = i;
782                     }
783 
784                     v = ps[i][2];
785                     if (v < mi_y) {
786                         mi_y = v;
787                         mi_y_i = i;
788                     } else if (v > ma_y) {
789                         ma_y = v;
790                         ma_y_i = i;
791                     }
792 
793                     v = ps[i][1] - ps[i][2];
794                     if (v < mi_xmy) {
795                         mi_xmy = v;
796                         mi_xmy_i = i;
797                     } else if (v > ma_xmy) {
798                         ma_xmy = v;
799                         ma_xmy_i = i;
800                     }
801 
802                     v = ps[i][1] + ps[i][2];
803                     if (v < mi_xpy) {
804                         mi_xpy = v;
805                         mi_xpy_i = i;
806                     } else if (v > ma_xpy) {
807                         ma_xpy = v;
808                         ma_xpy_i = i;
809                     }
810                 }
811             }
812 
813             // Keep track of the indices of the input points.
814             for (i = 0; i < N; i++) {
815                 c = ps[i];
816                 if (N <= AklToussaint ||
817                     // Discard inside of the octagon according to the Akl-Toussaint heuristic
818                     i in [mi_x_i, ma_x_i, mi_y_i, ma_y_i, mi_xpy_i, mi_xmy_i, ma_xpy_i, ma_xmy_i] ||
819                     (mi_x_i !== mi_xmy_i && this.signedTriangle(ps[mi_x_i], ps[mi_xmy_i], c) >= -eps) ||
820                     (mi_xmy_i !== ma_y_i && this.signedTriangle(ps[mi_xmy_i], ps[ma_y_i], c) >= -eps) ||
821                     (ma_y_i !== ma_xpy_i && this.signedTriangle(ps[ma_y_i], ps[ma_xpy_i], c) >= -eps) ||
822                     (ma_xpy_i !== ma_x_i && this.signedTriangle(ps[ma_xpy_i], ps[ma_x_i], c) >= -eps) ||
823                     (ma_x_i !== ma_xmy_i && this.signedTriangle(ps[ma_x_i], ps[ma_xmy_i], c) >= -eps) ||
824                     (ma_xmy_i !== mi_y_i && this.signedTriangle(ps[ma_xmy_i], ps[mi_y_i], c) >= -eps) ||
825                     (mi_y_i !== mi_xpy_i && this.signedTriangle(ps[mi_y_i], ps[mi_xpy_i], c) >= -eps) ||
826                     (mi_xpy_i !== mi_x_i && this.signedTriangle(ps[mi_xpy_i], ps[mi_x_i], c) >= -eps)
827                 ) {
828                     ps_idx.push({
829                         i: i,
830                         c: c
831                     });
832                 }
833             }
834             N = ps_idx.length;
835 
836             // Find the point with the lowest y value
837             mi_idx = 0;
838             mi_x = ps_idx[0].c[1];
839             mi_y = ps_idx[0].c[2];
840             for (i = 1; i < N; i++) {
841                 if ((ps_idx[i].c[2] < mi_y) || (ps_idx[i].c[2] === mi_y && ps_idx[i].c[1] < mi_x)) {
842                     mi_x = ps_idx[i].c[1];
843                     mi_y = ps_idx[i].c[2];
844                     mi_idx = i;
845                 }
846             }
847             ps_idx = Type.swap(ps_idx, mi_idx, 0);
848 
849             // Our origin o, i.e. the first point.
850             o = ps_idx[0].c;
851 
852             // Sort according to the angle around o.
853             ps_idx.sort(function(a_obj, b_obj) {
854                 var a = a_obj.c,
855                     b = b_obj.c,
856                     v = that.signedTriangle(o, a, b);
857 
858                 if (v === 0) {
859                     // if o, a, b are collinear, the point which is further away
860                     // from o is considered greater.
861                     return Mat.hypot(a[1] - o[1], a[2] - o[2]) - Mat.hypot(b[1] - o[1], b[2] - o[2]);
862                 }
863 
864                 // if v < 0, a is to the left of [o, b], i.e. angle(a) > angle(b)
865                 return -v;
866             });
867 
868             // Do the Graham scan.
869             M = 0;
870             for (i = 0; i < N; i++) {
871                 while (M > 1 && this.signedTriangle(stack[M - 2].c, stack[M - 1].c, ps_idx[i].c) <= 0) {
872                     // stack[M - 1] is to the left of stack[M - 1], ps[i]: discard it
873                     stack.pop();
874                     M--;
875                 }
876                 stack.push(ps_idx[i]);
877                 M++;
878             }
879 
880             return stack;
881         },
882 
883         // Original method
884         // GrahamScan: function (points, indices) {
885         //     var i,
886         //         M = 1,
887         //         ps = Expect.each(points, Expect.coordsArray),
888         //         N = ps.length;
889         //     ps = this.sortVertices(ps);
890         //     N = ps.length;
891         //     for (i = 2; i < N; i++) {
892         //         while (this.signedTriangle(ps[M - 1], ps[M], ps[i]) <= 0) {
893         //             if (M > 1) {
894         //                 M -= 1;
895         //             } else if (i === N - 1) {
896         //                 break;
897         //             }
898         //             i += 1;
899         //         }
900         //         M += 1;
901         //         ps = Type.swap(ps, M, i);
902         //         indices = Type.swap(indices, M, i);
903         //     }
904         //     return ps.slice(0, M);
905         // },
906 
907         /**
908          * Calculate the complex hull of a point cloud by the Graham scan algorithm.
909          *
910          * @param {Array} points An array containing {@link JXG.Point}, {@link JXG.Coords}, and/or arrays.
911          * @param {Boolean} [returnCoords=false] If true, return an array of coords. Otherwise return a list of pointers
912          * to the input list elements. That is, if the input is a list of {@link JXG.Point} elements, the returned list
913          * will contain the points that form the convex hull.
914          * @returns {Array} List containing the convex hull. Format depends on returnCoords.
915          * @see JXG.Math.Geometry.GrahamScan
916          *
917          * @example
918          *     // Static example
919          *     var i, hull,
920          *         p = [];
921          *
922          *     p.push( board.create('point', [4, 0], {withLabel:false }) );
923          *     p.push( board.create('point', [0, 4], {withLabel:false }) );
924          *     p.push( board.create('point', [0, 0], {withLabel:false }) );
925          *     p.push( board.create('point', [1, 1], {withLabel:false }) );
926          *     hull = JXG.Math.Geometry.convexHull(p);
927          *     for (i = 0; i < hull.length; i++) {
928          *       hull[i].setAttribute({color: 'blue'});
929          *     }
930          *
931          * </pre><div id="JXGdfc76123-81b8-4250-96f9-419253bd95dd" class="jxgbox" style="width: 300px; height: 300px;"></div>
932          * <script type="text/javascript">
933          *     (function() {
934          *         var board = JXG.JSXGraph.initBoard('JXGdfc76123-81b8-4250-96f9-419253bd95dd',
935          *             {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false});
936          *         var i, hull,
937          *             p = [];
938          *
939          *         p.push( board.create('point', [4, 0], {withLabel:false }) );
940          *         p.push( board.create('point', [0, 4], {withLabel:false }) );
941          *         p.push( board.create('point', [0, 0], {withLabel:false }) );
942          *         p.push( board.create('point', [1, 1], {withLabel:false }) );
943          *         hull = JXG.Math.Geometry.convexHull(p);
944          *         for (i = 0; i < hull.length; i++) {
945          *           hull[i].setAttribute({color: 'blue'});
946          *         }
947          *
948          *     })();
949          *
950          * </script><pre>
951          *
952          * @example
953          *     // Dynamic version using returnCoords==true: drag the points
954          *     var p = [];
955          *
956          *     p.push( board.create('point', [4, 0], {withLabel:false }) );
957          *     p.push( board.create('point', [0, 4], {withLabel:false }) );
958          *     p.push( board.create('point', [0, 0], {withLabel:false }) );
959          *     p.push( board.create('point', [1, 1], {withLabel:false }) );
960          *
961          *     var c = board.create('curve', [[], []], {fillColor: 'yellow', fillOpacity: 0.3});
962          *     c.updateDataArray = function() {
963          *       var i,
964          *         hull = JXG.Math.Geometry.convexHull(p, true);
965          *
966          *       this.dataX = [];
967          *       this.dataY = [];
968          *
969          *       for (i = 0; i < hull.length; i ++) {
970          *         this.dataX.push(hull[i][1]);
971          *         this.dataY.push(hull[i][2]);
972          *       }
973          *       this.dataX.push(hull[0][1]);
974          *       this.dataY.push(hull[0][2]);
975          *     };
976          *     board.update();
977          *
978          * </pre><div id="JXG61e51909-da0b-432f-9aa7-9fb0c8bb01c9" class="jxgbox" style="width: 300px; height: 300px;"></div>
979          * <script type="text/javascript">
980          *     (function() {
981          *         var board = JXG.JSXGraph.initBoard('JXG61e51909-da0b-432f-9aa7-9fb0c8bb01c9',
982          *             {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false});
983          *         var p = [];
984          *
985          *         p.push( board.create('point', [4, 0], {withLabel:false }) );
986          *         p.push( board.create('point', [0, 4], {withLabel:false }) );
987          *         p.push( board.create('point', [0, 0], {withLabel:false }) );
988          *         p.push( board.create('point', [1, 1], {withLabel:false }) );
989          *
990          *         var c = board.create('curve', [[], []], {fillColor: 'yellow', fillOpacity: 0.3});
991          *         c.updateDataArray = function() {
992          *           var i,
993          *             hull = JXG.Math.Geometry.convexHull(p, true);
994          *
995          *           this.dataX = [];
996          *           this.dataY = [];
997          *
998          *           for (i = 0; i < hull.length; i ++) {
999          *             this.dataX.push(hull[i][1]);
1000          *             this.dataY.push(hull[i][2]);
1001          *           }
1002          *           this.dataX.push(hull[0][1]);
1003          *           this.dataY.push(hull[0][2]);
1004          *         };
1005          *         board.update();
1006          *
1007          *
1008          *     })();
1009          *
1010          * </script><pre>
1011          *
1012          */
1013         convexHull: function(points, returnCoords) {
1014             var i, hull,
1015                 res = [];
1016 
1017             hull = this.GrahamScan(points);
1018             for (i = 0; i < hull.length; i++) {
1019                 if (returnCoords) {
1020                     res.push(hull[i].c);
1021                 } else {
1022                     res.push(points[hull[i].i]);
1023                 }
1024             }
1025             return res;
1026         },
1027 
1028         // /**
1029         //  * Determine if a polygon or a path element is convex, non-convex or complex which are defined like this:
1030         //  * <ul>
1031         //  * <li> A polygon is convex if for every pair of points, the line segment connecting them does not intersect
1032         //  * an edge of the polygon in one point.
1033         //  * A single line segment or a a single point is considered as convex. A necessary condition for a polygon
1034         //  * to be convex that the angle sum of its interior angles equals ± 2 π.
1035         //  * <li> A polygon is non-convex, if it does not self-intersect, but is not convex.
1036         //  * <li> A polygon is complex if its the angle sum is not equal to ± 2 π.
1037         //  * That is, there must be self-intersection (contiguous coincident points in the path are not treated as self-intersection).
1038         //  * </ul>
1039         //  * A path  element might be specified as an array of coordinate arrays or {@link JXG.Coords}.
1040         //  *
1041         //  * @param {Array|Polygon|PolygonalChain} points Polygon or list of coordinates
1042         //  * @returns {Number} -1: if complex, 0: if non-convex, 1: if convex
1043         //  */
1044         /**
1045          * Determine if a polygon or a path element is convex:
1046          * <p>
1047          * A polygon is convex if for every pair of points, the line segment connecting them does not intersect
1048          * an edge of the polygon in one point.
1049          * A single line segment, a single point, or the empty set is considered as convex. A necessary condition for a polygon
1050          * to be convex that the angle sum of its interior angles equals ± 2 π.
1051          * <p>
1052          * A path  element might be specified as an array of coordinate arrays or {@link JXG.Coords}.
1053          * See the discussion at <a href="https://stackoverflow.com/questions/471962/how-do-i-efficiently-determine-if-a-polygon-is-convex-non-convex-or-complex">stackoverflow</a>.
1054          *
1055          * @param {Array|Polygon|PolygonalChain} points Polygon or list of coordinates
1056          * @returns {Boolean} true if convex
1057          *
1058          * @example
1059          * var pol = board.create('polygon', [
1060          *     [-1, -1],
1061          *     [3, -1],
1062          *     [4, 2],
1063          *     [3, 3],
1064          *     [0, 4],
1065          *     [-3, 1]
1066          * ], {
1067          *     vertices: {
1068          *         color: 'blue',
1069          *         snapToGrid: true
1070          *     }
1071          * });
1072          *
1073          * console.log(JXG.Math.Geometry.isConvex(pol));
1074          * // > true
1075          *
1076          *
1077          *
1078          * </pre><div id="JXG9b43cc53-15b4-49be-92cc-2a1dfc06665b" class="jxgbox" style="width: 300px; height: 300px;"></div>
1079          * <script type="text/javascript">
1080          *     (function() {
1081          *         var board = JXG.JSXGraph.initBoard('JXG9b43cc53-15b4-49be-92cc-2a1dfc06665b',
1082          *             {boundingbox: [-8, 8, 8,-8], axis: true, showcopyright: false, shownavigation: false});
1083          *     var pol = board.create('polygon', [
1084          *         [-1, -1],
1085          *         [3, -1],
1086          *         [4, 2],
1087          *         [3, 3],
1088          *         [0, 4],
1089          *         [-3, 1]
1090          *     ], {
1091          *         vertices: {
1092          *             color: 'blue',
1093          *             snapToGrid: true
1094          *         }
1095          *     });
1096          *
1097          *     console.log(JXG.Math.Geometry.isConvex(pol));
1098          *
1099          *
1100          *
1101          *     })();
1102          *
1103          * </script><pre>
1104          *
1105          */
1106         isConvex: function(points) {
1107             var ps, le, i,
1108                 eps = Mat.eps * Mat.eps,
1109                 old_x, old_y, old_dir,
1110                 new_x, new_y, new_dir,
1111                 angle,
1112                 orient,
1113                 angle_sum = 0.0;
1114 
1115             if (Type.isArray(points)) {
1116                 ps = Expect.each(points, Expect.coordsArray);
1117             } else if (Type.exists(points.type) && points.type === Const.OBJECT_TYPE_POLYGON) {
1118                 ps = Expect.each(points.vertices, Expect.coordsArray);
1119             }
1120             le = ps.length;
1121             if (le === 0) {
1122                 // Empty set is convex
1123                 return true;
1124             }
1125             if (le < 3) {
1126                 // Segments and points are convex
1127                 return true;
1128             }
1129 
1130             orient = null;
1131             old_x = ps[le - 2][1];
1132             old_y = ps[le - 2][2];
1133             new_x = ps[le - 1][1];
1134             new_y = ps[le - 1][2];
1135             new_dir = Math.atan2(new_y - old_y, new_x - old_x);
1136             for (i = 0; i < le; i++) {
1137                 old_x = new_x;
1138                 old_y = new_y;
1139                 old_dir = new_dir;
1140                 new_x = ps[i][1];
1141                 new_y = ps[i][2];
1142                 if (old_x === new_x && old_y === new_y) {
1143                     // Repeated consecutive points are ignored
1144                     continue;
1145                 }
1146                 new_dir = Math.atan2(new_y - old_y, new_x - old_x);
1147                 angle = new_dir - old_dir;
1148                 if (angle <= -Math.PI) {
1149                     angle += 2 * Math.PI;
1150                 } else if (angle > Math.PI) {
1151                     angle -= 2 * Math.PI;
1152                 }
1153                 if (orient === null) {
1154                     if (angle === 0.0) {
1155                         continue;
1156                     }
1157                     orient = (angle > 0) ? 1 : -1;
1158                 } else {
1159                     if (orient * angle < -eps) {
1160                         return false;
1161                     }
1162                 }
1163                 angle_sum += angle;
1164             }
1165 
1166             if ((Math.abs(angle_sum / (2 * Math.PI)) - 1) < eps) {
1167                 return true;
1168             }
1169             return false;
1170         },
1171 
1172         /**
1173          * A line can be a segment, a straight, or a ray. So it is not always delimited by point1 and point2
1174          * calcStraight determines the visual start point and end point of the line. A segment is only drawn
1175          * from start to end point, a straight line is drawn until it meets the boards boundaries.
1176          * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point.
1177          * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and
1178          * set by this method.
1179          * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set
1180          * by this method.
1181          * @param {Number} margin Optional margin, to avoid the display of the small sides of lines.
1182          * @returns null
1183          * @see Line
1184          * @see JXG.Line
1185          */
1186         calcStraight: function (el, point1, point2, margin) {
1187             var takePoint1,
1188                 takePoint2,
1189                 intersection,
1190                 intersect1,
1191                 intersect2,
1192                 straightFirst,
1193                 straightLast,
1194                 c, p1, p2;
1195 
1196             if (!Type.exists(margin)) {
1197                 // Enlarge the drawable region slightly. This hides the small sides
1198                 // of thick lines in most cases.
1199                 margin = 10;
1200             }
1201 
1202             straightFirst = el.evalVisProp('straightfirst');
1203             straightLast = el.evalVisProp('straightlast');
1204 
1205             // If one of the point is an ideal point in homogeneous coordinates
1206             // drawing of line segments or rays are not possible.
1207             if (Math.abs(point1.scrCoords[0]) < Mat.eps) {
1208                 straightFirst = true;
1209             }
1210             if (Math.abs(point2.scrCoords[0]) < Mat.eps) {
1211                 straightLast = true;
1212             }
1213 
1214             // Do nothing in case of line segments (inside or outside of the board)
1215             if (!straightFirst && !straightLast) {
1216                 return;
1217             }
1218 
1219             // Compute the stdform of the line in screen coordinates.
1220             c = [];
1221             c[0] =
1222                 el.stdform[0] -
1223                 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX +
1224                 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY;
1225             c[1] = el.stdform[1] / el.board.unitX;
1226             c[2] = -el.stdform[2] / el.board.unitY;
1227 
1228             // If p1=p2
1229             if (isNaN(c[0] + c[1] + c[2])) {
1230                 return;
1231             }
1232 
1233             takePoint1 = false;
1234             takePoint2 = false;
1235 
1236             // Line starts at point1 and point1 is inside the board
1237             takePoint1 =
1238                 !straightFirst &&
1239                 Math.abs(point1.usrCoords[0]) >= Mat.eps &&
1240                 point1.scrCoords[1] >= 0.0 &&
1241                 point1.scrCoords[1] <= el.board.canvasWidth &&
1242                 point1.scrCoords[2] >= 0.0 &&
1243                 point1.scrCoords[2] <= el.board.canvasHeight;
1244 
1245             // Line ends at point2 and point2 is inside the board
1246             takePoint2 =
1247                 !straightLast &&
1248                 Math.abs(point2.usrCoords[0]) >= Mat.eps &&
1249                 point2.scrCoords[1] >= 0.0 &&
1250                 point2.scrCoords[1] <= el.board.canvasWidth &&
1251                 point2.scrCoords[2] >= 0.0 &&
1252                 point2.scrCoords[2] <= el.board.canvasHeight;
1253 
1254             // Intersect the line with the four borders of the board.
1255             intersection = this.meetLineBoard(c, el.board, margin);
1256             intersect1 = intersection[0];
1257             intersect2 = intersection[1];
1258 
1259             /**
1260              * At this point we have four points:
1261              * point1 and point2 are the first and the second defining point on the line,
1262              * intersect1, intersect2 are the intersections of the line with border around the board.
1263              */
1264 
1265             /*
1266              * Here we handle rays where both defining points are outside of the board.
1267              */
1268             // If both points are outside and the complete ray is outside we do nothing
1269             if (!takePoint1 && !takePoint2) {
1270                 // Ray starting at point 1
1271                 if (
1272                     !straightFirst &&
1273                     straightLast &&
1274                     !this.isSameDirection(point1, point2, intersect1) &&
1275                     !this.isSameDirection(point1, point2, intersect2)
1276                 ) {
1277                     return;
1278                 }
1279 
1280                 // Ray starting at point 2
1281                 if (
1282                     straightFirst &&
1283                     !straightLast &&
1284                     !this.isSameDirection(point2, point1, intersect1) &&
1285                     !this.isSameDirection(point2, point1, intersect2)
1286                 ) {
1287                     return;
1288                 }
1289             }
1290 
1291             /*
1292              * If at least one of the defining points is outside of the board
1293              * we take intersect1 or intersect2 as one of the end points
1294              * The order is also important for arrows of axes
1295              */
1296             if (!takePoint1) {
1297                 if (!takePoint2) {
1298                     // Two border intersection points are used
1299                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1300                         p1 = intersect1;
1301                         p2 = intersect2;
1302                     } else {
1303                         p2 = intersect1;
1304                         p1 = intersect2;
1305                     }
1306                 } else {
1307                     // One border intersection points is used
1308                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1309                         p1 = intersect1;
1310                     } else {
1311                         p1 = intersect2;
1312                     }
1313                 }
1314             } else {
1315                 if (!takePoint2) {
1316                     // One border intersection points is used
1317                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1318                         p2 = intersect2;
1319                     } else {
1320                         p2 = intersect1;
1321                     }
1322                 }
1323             }
1324 
1325             if (p1) {
1326                 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1));
1327                 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords);
1328             }
1329 
1330             if (p2) {
1331                 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1));
1332                 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords);
1333             }
1334         },
1335 
1336         /**
1337          * A line can be a segment, a straight, or a ray. so it is not always delimited by point1 and point2.
1338          *
1339          * This method adjusts the line's delimiting points taking into account its nature, the viewport defined
1340          * by the board.
1341          *
1342          * A segment is delimited by start and end point, a straight line or ray is delimited until it meets the
1343          * boards boundaries. However, if the line has infinite ticks, it will be delimited by the projection of
1344          * the boards vertices onto itself.
1345          *
1346          * @param {JXG.Line} el Reference to a line object, that needs calculation of start and end point.
1347          * @param {JXG.Coords} point1 Coordinates of the point where line drawing begins. This value is calculated and
1348          * set by this method.
1349          * @param {JXG.Coords} point2 Coordinates of the point where line drawing ends. This value is calculated and set
1350          * by this method.
1351          * @see Line
1352          * @see JXG.Line
1353          */
1354         calcLineDelimitingPoints: function (el, point1, point2) {
1355             var distP1P2,
1356                 boundingBox,
1357                 lineSlope,
1358                 intersect1,
1359                 intersect2,
1360                 straightFirst,
1361                 straightLast,
1362                 c,
1363                 p1,
1364                 p2,
1365                 takePoint1 = false,
1366                 takePoint2 = false;
1367 
1368             straightFirst = el.evalVisProp('straightfirst');
1369             straightLast = el.evalVisProp('straightlast');
1370 
1371             // If one of the point is an ideal point in homogeneous coordinates
1372             // drawing of line segments or rays are not possible.
1373             if (Math.abs(point1.scrCoords[0]) < Mat.eps) {
1374                 straightFirst = true;
1375             }
1376             if (Math.abs(point2.scrCoords[0]) < Mat.eps) {
1377                 straightLast = true;
1378             }
1379 
1380             // Compute the stdform of the line in screen coordinates.
1381             c = [];
1382             c[0] =
1383                 el.stdform[0] -
1384                 (el.stdform[1] * el.board.origin.scrCoords[1]) / el.board.unitX +
1385                 (el.stdform[2] * el.board.origin.scrCoords[2]) / el.board.unitY;
1386             c[1] = el.stdform[1] / el.board.unitX;
1387             c[2] = -el.stdform[2] / el.board.unitY;
1388 
1389             // p1=p2
1390             if (isNaN(c[0] + c[1] + c[2])) {
1391                 return;
1392             }
1393 
1394             takePoint1 = !straightFirst;
1395             takePoint2 = !straightLast;
1396             // Intersect the board vertices on the line to establish the available visual space for the infinite ticks
1397             // Based on the slope of the line we can optimise and only project the two outer vertices
1398 
1399             // boundingBox = [x1, y1, x2, y2] upper left, lower right vertices
1400             boundingBox = el.board.getBoundingBox();
1401             lineSlope = el.getSlope();
1402             if (lineSlope >= 0) {
1403                 // project vertices (x2,y1) (x1, y2)
1404                 intersect1 = this.projectPointToLine(
1405                     { coords: { usrCoords: [1, boundingBox[2], boundingBox[1]] } },
1406                     el,
1407                     el.board
1408                 );
1409                 intersect2 = this.projectPointToLine(
1410                     { coords: { usrCoords: [1, boundingBox[0], boundingBox[3]] } },
1411                     el,
1412                     el.board
1413                 );
1414             } else {
1415                 // project vertices (x1, y1) (x2, y2)
1416                 intersect1 = this.projectPointToLine(
1417                     { coords: { usrCoords: [1, boundingBox[0], boundingBox[1]] } },
1418                     el,
1419                     el.board
1420                 );
1421                 intersect2 = this.projectPointToLine(
1422                     { coords: { usrCoords: [1, boundingBox[2], boundingBox[3]] } },
1423                     el,
1424                     el.board
1425                 );
1426             }
1427 
1428             /**
1429              * we have four points:
1430              * point1 and point2 are the first and the second defining point on the line,
1431              * intersect1, intersect2 are the intersections of the line with border around the board.
1432              */
1433 
1434             /*
1435              * Here we handle rays/segments where both defining points are outside of the board.
1436              */
1437             if (!takePoint1 && !takePoint2) {
1438                 // Segment, if segment does not cross the board, do nothing
1439                 if (!straightFirst && !straightLast) {
1440                     distP1P2 = point1.distance(Const.COORDS_BY_USER, point2);
1441                     // if  intersect1 not between point1 and point2
1442                     if (
1443                         Math.abs(
1444                             point1.distance(Const.COORDS_BY_USER, intersect1) +
1445                             intersect1.distance(Const.COORDS_BY_USER, point2) -
1446                             distP1P2
1447                         ) > Mat.eps
1448                     ) {
1449                         return;
1450                     }
1451                     // if insersect2 not between point1 and point2
1452                     if (
1453                         Math.abs(
1454                             point1.distance(Const.COORDS_BY_USER, intersect2) +
1455                             intersect2.distance(Const.COORDS_BY_USER, point2) -
1456                             distP1P2
1457                         ) > Mat.eps
1458                     ) {
1459                         return;
1460                     }
1461                 }
1462 
1463                 // If both points are outside and the complete ray is outside we do nothing
1464                 // Ray starting at point 1
1465                 if (
1466                     !straightFirst &&
1467                     straightLast &&
1468                     !this.isSameDirection(point1, point2, intersect1) &&
1469                     !this.isSameDirection(point1, point2, intersect2)
1470                 ) {
1471                     return;
1472                 }
1473 
1474                 // Ray starting at point 2
1475                 if (
1476                     straightFirst &&
1477                     !straightLast &&
1478                     !this.isSameDirection(point2, point1, intersect1) &&
1479                     !this.isSameDirection(point2, point1, intersect2)
1480                 ) {
1481                     return;
1482                 }
1483             }
1484 
1485             /*
1486              * If at least one of the defining points is outside of the board
1487              * we take intersect1 or intersect2 as one of the end points
1488              * The order is also important for arrows of axes
1489              */
1490             if (!takePoint1) {
1491                 if (!takePoint2) {
1492                     // Two border intersection points are used
1493                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1494                         p1 = intersect1;
1495                         p2 = intersect2;
1496                     } else {
1497                         p2 = intersect1;
1498                         p1 = intersect2;
1499                     }
1500                 } else {
1501                     // One border intersection points is used
1502                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1503                         p1 = intersect1;
1504                     } else {
1505                         p1 = intersect2;
1506                     }
1507                 }
1508             } else {
1509                 if (!takePoint2) {
1510                     // One border intersection points is used
1511                     if (this.isSameDir(point1, point2, intersect1, intersect2)) {
1512                         p2 = intersect2;
1513                     } else {
1514                         p2 = intersect1;
1515                     }
1516                 }
1517             }
1518 
1519             if (p1) {
1520                 //point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords.slice(1));
1521                 point1.setCoordinates(Const.COORDS_BY_USER, p1.usrCoords);
1522             }
1523 
1524             if (p2) {
1525                 //point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords.slice(1));
1526                 point2.setCoordinates(Const.COORDS_BY_USER, p2.usrCoords);
1527             }
1528         },
1529 
1530         /**
1531          * Calculates the visProp.position corresponding to a given angle.
1532          * @param {number} angle angle in radians. Must be in range (-2pi,2pi).
1533          */
1534         calcLabelQuadrant: function (angle) {
1535             var q;
1536             if (angle < 0) {
1537                 angle += 2 * Math.PI;
1538             }
1539             q = Math.floor((angle + Math.PI / 8) / (Math.PI / 4)) % 8;
1540             return ["rt", "urt", "top", "ulft", "lft", "llft", "lrt"][q];
1541         },
1542 
1543         /**
1544          * The vectors <tt>p2-p1</tt> and <tt>i2-i1</tt> are supposed to be collinear. If their cosine is positive
1545          * they point into the same direction otherwise they point in opposite direction.
1546          * @param {JXG.Coords} p1
1547          * @param {JXG.Coords} p2
1548          * @param {JXG.Coords} i1
1549          * @param {JXG.Coords} i2
1550          * @returns {Boolean} True, if <tt>p2-p1</tt> and <tt>i2-i1</tt> point into the same direction
1551          */
1552         isSameDir: function (p1, p2, i1, i2) {
1553             var dpx = p2.usrCoords[1] - p1.usrCoords[1],
1554                 dpy = p2.usrCoords[2] - p1.usrCoords[2],
1555                 dix = i2.usrCoords[1] - i1.usrCoords[1],
1556                 diy = i2.usrCoords[2] - i1.usrCoords[2];
1557 
1558             if (Math.abs(p2.usrCoords[0]) < Mat.eps) {
1559                 dpx = p2.usrCoords[1];
1560                 dpy = p2.usrCoords[2];
1561             }
1562 
1563             if (Math.abs(p1.usrCoords[0]) < Mat.eps) {
1564                 dpx = -p1.usrCoords[1];
1565                 dpy = -p1.usrCoords[2];
1566             }
1567 
1568             return dpx * dix + dpy * diy >= 0;
1569         },
1570 
1571         /**
1572          * If you're looking from point "start" towards point "s" and you can see the point "p", return true.
1573          * Otherwise return false.
1574          * @param {JXG.Coords} start The point you're standing on.
1575          * @param {JXG.Coords} p The point in which direction you're looking.
1576          * @param {JXG.Coords} s The point that should be visible.
1577          * @returns {Boolean} True, if from start the point p is in the same direction as s is, that means s-start = k*(p-start) with k>=0.
1578          */
1579         isSameDirection: function (start, p, s) {
1580             var dx,
1581                 dy,
1582                 sx,
1583                 sy,
1584                 r = false;
1585 
1586             dx = p.usrCoords[1] - start.usrCoords[1];
1587             dy = p.usrCoords[2] - start.usrCoords[2];
1588 
1589             sx = s.usrCoords[1] - start.usrCoords[1];
1590             sy = s.usrCoords[2] - start.usrCoords[2];
1591 
1592             if (Math.abs(dx) < Mat.eps) {
1593                 dx = 0;
1594             }
1595 
1596             if (Math.abs(dy) < Mat.eps) {
1597                 dy = 0;
1598             }
1599 
1600             if (Math.abs(sx) < Mat.eps) {
1601                 sx = 0;
1602             }
1603 
1604             if (Math.abs(sy) < Mat.eps) {
1605                 sy = 0;
1606             }
1607 
1608             if (dx >= 0 && sx >= 0) {
1609                 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0);
1610             } else if (dx <= 0 && sx <= 0) {
1611                 r = (dy >= 0 && sy >= 0) || (dy <= 0 && sy <= 0);
1612             }
1613 
1614             return r;
1615         },
1616 
1617         /**
1618          * Determinant of three points in the Euclidean plane.
1619          * Zero, if the points are collinear. Used to determine of a point q is left or
1620          * right to a segment defined by points p1 and p2.
1621          * <p>
1622          * Non-homogeneous version.
1623          *
1624          * @param  {Array|JXG.Point} p1 First point or its coordinates of the segment. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1625          * @param  {Array|JXG.Point} p2 Second point or its coordinates of the segment. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1626          * @param  {Array|JXG.Point} q Point or its coordinates. Point object or array of length 3. First (homogeneous) coordinate is equal to 1.
1627          * @return {Number} Signed area of the triangle formed by these three points.
1628          *
1629          * @see JXG.Math.Geometry.windingNumber
1630          */
1631         det3p: function (p1, p2, q) {
1632             var pp1, pp2, qq;
1633 
1634             if (Type.isPoint(p1)) {
1635                 pp1 = p1.Coords(true);
1636             } else {
1637                 pp1 = p1;
1638             }
1639             if (Type.isPoint(p2)) {
1640                 pp2 = p2.Coords(true);
1641             } else {
1642                 pp2 = p2;
1643             }
1644             if (Type.isPoint(q)) {
1645                 qq = q.Coords(true);
1646             } else {
1647                 qq = q;
1648             }
1649 
1650             return (pp1[1] - qq[1]) * (pp2[2] - qq[2]) - (pp2[1] - qq[1]) * (pp1[2] - qq[2]);
1651         },
1652 
1653         /**
1654          * Winding number of a point in respect to a polygon path.
1655          *
1656          * The point is regarded outside if the winding number is zero,
1657          * inside otherwise. The algorithm tries to find degenerate cases, i.e.
1658          * if the point is on the path. This is regarded as "outside".
1659          * If the point is a vertex of the path, it is regarded as "inside".
1660          *
1661          * Implementation of algorithm 7 from "The point in polygon problem for
1662          * arbitrary polygons" by Kai Hormann and Alexander Agathos, Computational Geometry,
1663          * Volume 20, Issue 3, November 2001, Pages 131-144.
1664          *
1665          * @param  {Array} usrCoords Homogenous coordinates of the point
1666          * @param  {Array} path      Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements
1667          * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords.
1668          * @param  {Boolean} [doNotClosePath=false] If true the last point of the path is not connected to the first point.
1669          * This is necessary if the path consists of two or more closed subpaths, e.g. if the figure has a hole.
1670          *
1671          * @return {Number}          Winding number of the point. The point is
1672          *                           regarded outside if the winding number is zero,
1673          *                           inside otherwise.
1674          */
1675         windingNumber: function (usrCoords, path, doNotClosePath) {
1676             var wn = 0,
1677                 le = path.length,
1678                 x = usrCoords[1],
1679                 y = usrCoords[2],
1680                 p0,
1681                 p1,
1682                 p2,
1683                 d,
1684                 sign,
1685                 i,
1686                 off = 0;
1687 
1688             if (le === 0) {
1689                 return 0;
1690             }
1691 
1692             doNotClosePath = doNotClosePath || false;
1693             if (doNotClosePath) {
1694                 off = 1;
1695             }
1696 
1697             // Infinite points are declared outside
1698             if (isNaN(x) || isNaN(y)) {
1699                 return 1;
1700             }
1701 
1702             if (Type.exists(path[0].coords)) {
1703                 p0 = path[0].coords;
1704                 p1 = path[le - 1].coords;
1705             } else {
1706                 p0 = path[0];
1707                 p1 = path[le - 1];
1708             }
1709             // Handle the case if the point is the first vertex of the path, i.e. inside.
1710             if (p0.usrCoords[1] === x && p0.usrCoords[2] === y) {
1711                 return 1;
1712             }
1713 
1714             for (i = 0; i < le - off; i++) {
1715                 // Consider the edge from p1 = path[i] to p2 = path[i+1]isClosedPath
1716                 if (Type.exists(path[i].coords)) {
1717                     p1 = path[i].coords.usrCoords;
1718                     p2 = path[(i + 1) % le].coords.usrCoords;
1719                 } else {
1720                     p1 = path[i].usrCoords;
1721                     p2 = path[(i + 1) % le].usrCoords;
1722                 }
1723 
1724                 // If one of the two points p1, p2 is undefined or infinite,
1725                 // move on.
1726                 if (
1727                     p1[0] === 0 ||
1728                     p2[0] === 0 ||
1729                     isNaN(p1[1]) ||
1730                     isNaN(p2[1]) ||
1731                     isNaN(p1[2]) ||
1732                     isNaN(p2[2])
1733                 ) {
1734                     continue;
1735                 }
1736 
1737                 if (p2[2] === y) {
1738                     if (p2[1] === x) {
1739                         return 1;
1740                     }
1741                     if (p1[2] === y && p2[1] > x === p1[1] < x) {
1742                         return 0;
1743                     }
1744                 }
1745 
1746                 if (p1[2] < y !== p2[2] < y) {
1747                     // Crossing
1748                     sign = 2 * (p2[2] > p1[2] ? 1 : 0) - 1;
1749                     if (p1[1] >= x) {
1750                         if (p2[1] > x) {
1751                             wn += sign;
1752                         } else {
1753                             d = this.det3p(p1, p2, usrCoords);
1754                             if (d === 0) {
1755                                 // Point is on line, i.e. outside
1756                                 return 0;
1757                             }
1758                             if (d > 0 + Mat.eps === p2[2] > p1[2]) {
1759                                 // Right crossing
1760                                 wn += sign;
1761                             }
1762                         }
1763                     } else {
1764                         if (p2[1] > x) {
1765                             d = this.det3p(p1, p2, usrCoords);
1766                             if (d > 0 + Mat.eps === p2[2] > p1[2]) {
1767                                 // Right crossing
1768                                 wn += sign;
1769                             }
1770                         }
1771                     }
1772                 }
1773             }
1774 
1775             return wn;
1776         },
1777 
1778         /**
1779          * Decides if a point (x,y) is inside of a path / polygon.
1780          * Does not work correct if the path has hole. In this case, windingNumber is the preferred method.
1781          * Implements W. Randolf Franklin's pnpoly method.
1782          *
1783          * See <a href="https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html">https://wrf.ecse.rpi.edu/Research/Short_Notes/pnpoly.html</a>.
1784          *
1785          * @param {Number} x_in x-coordinate (screen or user coordinates)
1786          * @param {Number} y_in y-coordinate (screen or user coordinates)
1787          * @param  {Array} path  Array of points / coords determining a path, i.e. the vertices of the polygon / path. The array elements
1788          * do not have to be full points, but have to have a subobject "coords" or should be of type JXG.Coords.
1789          * @param {Number} [coord_type=JXG.COORDS_BY_SCREEN] Type of coordinates used here.
1790          *   Possible values are <b>JXG.COORDS_BY_USER</b> and <b>JXG.COORDS_BY_SCREEN</b>.
1791          *   Default value is JXG.COORDS_BY_SCREEN.
1792          * @param {JXG.Board} board Board object
1793          *
1794          * @returns {Boolean} if (x_in, y_in) is inside of the polygon.
1795          * @see JXG.Polygon#hasPoint
1796          * @see JXG.Polygon#pnpoly
1797          * @see JXG.Math.Geometry.windingNumber
1798          *
1799          * @example
1800          * var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]);
1801          * var p = board.create('point', [4, 3]);
1802          * var txt = board.create('text', [-1, 0.5, function() {
1803          *   return 'Point A is inside of the polygon = ' +
1804          *     JXG.Math.Geometry.pnpoly(p.X(), p.Y(), pol.vertices, JXG.COORDS_BY_USER, board);
1805          * }]);
1806          *
1807          * </pre><div id="JXG4656ed42-f965-4e35-bb66-c334a4529683" class="jxgbox" style="width: 300px; height: 300px;"></div>
1808          * <script type="text/javascript">
1809          *     (function() {
1810          *         var board = JXG.JSXGraph.initBoard('JXG4656ed42-f965-4e35-bb66-c334a4529683',
1811          *             {boundingbox: [-2, 5, 5,-2], axis: true, showcopyright: false, shownavigation: false});
1812          *     var pol = board.create('polygon', [[-1,2], [2,2], [-1,4]]);
1813          *     var p = board.create('point', [4, 3]);
1814          *     var txt = board.create('text', [-1, 0.5, function() {
1815          *     		return 'Point A is inside of the polygon = ' + JXG.Math.Geometry.pnpoly(p.X(), p.Y(), pol.vertices, JXG.COORDS_BY_USER, board);
1816          *     }]);
1817          *
1818          *     })();
1819          *
1820          * </script><pre>
1821          *
1822          */
1823         pnpoly: function (x_in, y_in, path, coord_type, board) {
1824             var i, j, vi, vj, len,
1825                 x, y, crds,
1826                 v = path,
1827                 isIn = false;
1828 
1829             if (coord_type === Const.COORDS_BY_USER) {
1830                 crds = new Coords(Const.COORDS_BY_USER, [x_in, y_in], board);
1831                 x = crds.scrCoords[1];
1832                 y = crds.scrCoords[2];
1833             } else {
1834                 x = x_in;
1835                 y = y_in;
1836             }
1837 
1838             len = path.length;
1839             for (i = 0, j = len - 2; i < len - 1; j = i++) {
1840                 vi = Type.exists(v[i].coords) ? v[i].coords : v[i];
1841                 vj = Type.exists(v[j].coords) ? v[j].coords : v[j];
1842 
1843                 if (
1844                     vi.scrCoords[2] > y !== vj.scrCoords[2] > y &&
1845                     x <
1846                     ((vj.scrCoords[1] - vi.scrCoords[1]) * (y - vi.scrCoords[2])) /
1847                     (vj.scrCoords[2] - vi.scrCoords[2]) +
1848                     vi.scrCoords[1]
1849                 ) {
1850                     isIn = !isIn;
1851                 }
1852             }
1853 
1854             return isIn;
1855         },
1856 
1857         /****************************************/
1858         /****          INTERSECTIONS         ****/
1859         /****************************************/
1860 
1861         /**
1862          * Generate the function which computes the coordinates of the intersection point.
1863          * Primarily used in {@link JXG.Point.createIntersectionPoint}.
1864          * @param {JXG.Board} board object
1865          * @param {JXG.Line,JXG.Circle_JXG.Line,JXG.Circle_Number|Function} el1,el2,i The result will be a intersection point on el1 and el2.
1866          * i determines the intersection point if two points are available: <ul>
1867          *   <li>i==0: use the positive square root,</li>
1868          *   <li>i==1: use the negative square root.</li></ul>
1869          * @param {Boolean} alwaysintersect. Flag that determines if segments and arc can have an outer intersection point
1870          * on their defining line or circle.
1871          * @returns {Function} Function returning a {@link JXG.Coords} object that determines
1872          * the intersection point.
1873          *
1874          * @see JXG.Point.createIntersectionPoint
1875          */
1876         intersectionFunction: function (board, el1, el2, i, j, alwaysintersect) {
1877             var func,
1878                 that = this,
1879                 el1_isArcType = false,
1880                 el2_isArcType = false;
1881 
1882             el1_isArcType =
1883                 el1.elementClass === Const.OBJECT_CLASS_CURVE &&
1884                     (el1.type === Const.OBJECT_TYPE_ARC || el1.type === Const.OBJECT_TYPE_SECTOR)
1885                     ? true
1886                     : false;
1887             el2_isArcType =
1888                 el2.elementClass === Const.OBJECT_CLASS_CURVE &&
1889                     (el2.type === Const.OBJECT_TYPE_ARC || el2.type === Const.OBJECT_TYPE_SECTOR)
1890                     ? true
1891                     : false;
1892 
1893             if (
1894                 (el1.elementClass === Const.OBJECT_CLASS_CURVE ||
1895                     el2.elementClass === Const.OBJECT_CLASS_CURVE) &&
1896                 (el1.elementClass === Const.OBJECT_CLASS_CURVE ||
1897                     el1.elementClass === Const.OBJECT_CLASS_CIRCLE) &&
1898                 (el2.elementClass === Const.OBJECT_CLASS_CURVE ||
1899                     el2.elementClass === Const.OBJECT_CLASS_CIRCLE) /*&&
1900                 !(el1_isArcType && el2_isArcType)*/
1901             ) {
1902                 // curve - curve
1903                 // with the exception that both elements are arc types
1904                 /** @ignore */
1905                 func = function () {
1906                     return that.meetCurveCurve(el1, el2, i, j, el1.board);
1907                 };
1908             } else if (
1909                 (el1.elementClass === Const.OBJECT_CLASS_CURVE &&
1910                     !el1_isArcType &&
1911                     el2.elementClass === Const.OBJECT_CLASS_LINE) ||
1912                 (el2.elementClass === Const.OBJECT_CLASS_CURVE &&
1913                     !el2_isArcType &&
1914                     el1.elementClass === Const.OBJECT_CLASS_LINE)
1915             ) {
1916                 // curve - line (this includes intersections between conic sections and lines)
1917                 // with the exception that the curve is of arc type
1918                 /** @ignore */
1919                 func = function () {
1920                     return that.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect));
1921                 };
1922             } else if (
1923                 el1.type === Const.OBJECT_TYPE_POLYGON ||
1924                 el2.type === Const.OBJECT_TYPE_POLYGON
1925             ) {
1926                 // polygon - other
1927                 // Uses the Greiner-Hormann clipping algorithm
1928                 // Not implemented: polygon - point
1929 
1930                 if (el1.elementClass === Const.OBJECT_CLASS_LINE) {
1931                     // line - path
1932                     /** @ignore */
1933                     func = function () {
1934                         var first1 = el1.evalVisProp('straightfirst'),
1935                             last1 = el1.evalVisProp('straightlast'),
1936                             first2 = el2.evalVisProp('straightfirst'),
1937                             last2 = el2.evalVisProp('straightlast'),
1938                             a_not;
1939 
1940                         a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2));
1941                         return that.meetPolygonLine(el2, el1, i, el1.board, a_not);
1942                     };
1943                 } else if (el2.elementClass === Const.OBJECT_CLASS_LINE) {
1944                     // path - line
1945                     /** @ignore */
1946                     func = function () {
1947                         var first1 = el1.evalVisProp('straightfirst'),
1948                             last1 = el1.evalVisProp('straightlast'),
1949                             first2 = el2.evalVisProp('straightfirst'),
1950                             last2 = el2.evalVisProp('straightlast'),
1951                             a_not;
1952 
1953                         a_not = (!Type.evaluate(alwaysintersect) && (!first1 || !last1 || !first2 || !last2));
1954                         return that.meetPolygonLine(el1, el2, i, el1.board, a_not);
1955                     };
1956                 } else {
1957                     // path - path
1958                     /** @ignore */
1959                     func = function () {
1960                         return that.meetPathPath(el1, el2, i, el1.board);
1961                     };
1962                 }
1963             } else if (
1964                 el1.elementClass === Const.OBJECT_CLASS_LINE &&
1965                 el2.elementClass === Const.OBJECT_CLASS_LINE
1966             ) {
1967                 // line - line, lines may also be segments.
1968                 /** @ignore */
1969                 func = function () {
1970                     var res,
1971                         c,
1972                         first1 = el1.evalVisProp('straightfirst'),
1973                         last1 = el1.evalVisProp('straightlast'),
1974                         first2 = el2.evalVisProp('straightfirst'),
1975                         last2 = el2.evalVisProp('straightlast');
1976 
1977                     /**
1978                      * If one of the lines is a segment or ray and
1979                      * the intersection point should disappear if outside
1980                      * of the segment or ray we call
1981                      * meetSegmentSegment
1982                      */
1983                     if (
1984                         !Type.evaluate(alwaysintersect) &&
1985                         (!first1 || !last1 || !first2 || !last2)
1986                     ) {
1987                         res = that.meetSegmentSegment(
1988                             el1.point1.coords.usrCoords,
1989                             el1.point2.coords.usrCoords,
1990                             el2.point1.coords.usrCoords,
1991                             el2.point2.coords.usrCoords
1992                         );
1993 
1994                         if (
1995                             (!first1 && res[1] < 0) ||
1996                             (!last1 && res[1] > 1) ||
1997                             (!first2 && res[2] < 0) ||
1998                             (!last2 && res[2] > 1)
1999                         ) {
2000                             // Non-existent
2001                             c = [0, NaN, NaN];
2002                         } else {
2003                             c = res[0];
2004                         }
2005 
2006                         return new Coords(Const.COORDS_BY_USER, c, el1.board);
2007                     }
2008 
2009                     return that.meet(el1.stdform, el2.stdform, i, el1.board);
2010                 };
2011             } else {
2012                 // All other combinations of circles and lines,
2013                 // Arc types are treated as circles.
2014                 /** @ignore */
2015                 func = function () {
2016                     var res = that.meet(el1.stdform, el2.stdform, i, el1.board),
2017                         has = true,
2018                         first,
2019                         last,
2020                         r;
2021 
2022                     if (Type.evaluate(alwaysintersect)) {
2023                         return res;
2024                     }
2025                     if (el1.elementClass === Const.OBJECT_CLASS_LINE) {
2026                         first = el1.evalVisProp('straightfirst');
2027                         last = el1.evalVisProp('straightlast');
2028                         if (!first || !last) {
2029                             r = that.affineRatio(el1.point1.coords, el1.point2.coords, res);
2030                             if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) {
2031                                 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
2032                             }
2033                         }
2034                     }
2035                     if (el2.elementClass === Const.OBJECT_CLASS_LINE) {
2036                         first = el2.evalVisProp('straightfirst');
2037                         last = el2.evalVisProp('straightlast');
2038                         if (!first || !last) {
2039                             r = that.affineRatio(el2.point1.coords, el2.point2.coords, res);
2040                             if ((!last && r > 1 + Mat.eps) || (!first && r < 0 - Mat.eps)) {
2041                                 return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
2042                             }
2043                         }
2044                     }
2045                     if (el1_isArcType) {
2046                         has = that.coordsOnArc(el1, res);
2047                         if (has && el2_isArcType) {
2048                             has = that.coordsOnArc(el2, res);
2049                         }
2050                         if (!has) {
2051                             return new Coords(JXG.COORDS_BY_USER, [0, NaN, NaN], el1.board);
2052                         }
2053                     }
2054                     return res;
2055                 };
2056             }
2057 
2058             return func;
2059         },
2060 
2061         otherIntersectionFunction: function (input, others, alwaysintersect, precision) {
2062             var func, board,
2063                 el1, el2,
2064                 that = this;
2065 
2066             el1 = input[0];
2067             el2 = input[1];
2068             board = el1.board;
2069             /** @ignore */
2070             func = function () {
2071                 var i, k, c, d,
2072                     isClose,
2073                     le = others.length,
2074                     eps = Type.evaluate(precision);
2075 
2076                 for (i = le; i >= 0; i--) {
2077                     if (el1.elementClass === Const.OBJECT_CLASS_CIRCLE &&
2078                         [Const.OBJECT_CLASS_CIRCLE, Const.OBJECT_CLASS_LINE].indexOf(el2.elementClass) >= 0) {
2079                         // circle, circle|line
2080                         c = that.meet(el1.stdform, el2.stdform, i, board);
2081                     } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE &&
2082                         [Const.OBJECT_CLASS_CURVE, Const.OBJECT_CLASS_CIRCLE].indexOf(el2.elementClass) >= 0) {
2083                         // curve, circle|curve
2084                         c = that.meetCurveCurve(el1, el2, i, 0, board);
2085                     } else if (el1.elementClass === Const.OBJECT_CLASS_CURVE && el2.elementClass === Const.OBJECT_CLASS_LINE) {
2086                         // curve, line
2087                         if (Type.exists(el1.dataX)) {
2088                             c = JXG.Math.Geometry.meetCurveLine(el1, el2, i, el1.board, Type.evaluate(alwaysintersect));
2089                         } else {
2090                             c = JXG.Math.Geometry.meetCurveLineContinuous(el1, el2, i, el1.board);
2091                         }
2092                     }
2093 
2094                     if (c === undefined) {
2095                         // Intersection point does not exist
2096                         continue;
2097                     }
2098 
2099                     // If the intersection is close to one of the points in other
2100                     // we have to search for another intersection point.
2101                     isClose = false;
2102                     for (k = 0; !isClose && k < le; k++) {
2103                         if (Type.exists(c) && Type.exists(c.distance)) {
2104                             d = c.distance(JXG.COORDS_BY_USER, others[k].coords);
2105                             if (d < eps) {
2106                                 isClose = true;
2107                             }
2108                         }
2109                     }
2110                     if (!isClose) {
2111                         // We are done, the intersection is away from any other
2112                         // intersection point.
2113                         return c;
2114                     }
2115                 }
2116                 // Otherwise we return the last intersection point
2117                 return c;
2118             };
2119             return func;
2120         },
2121 
2122         /**
2123          * Returns true if the coordinates are on the arc element,
2124          * false otherwise. Usually, coords is an intersection
2125          * on the circle line. Now it is decided if coords are on the
2126          * circle restricted to the arc line.
2127          * @param  {Arc} arc arc or sector element
2128          * @param  {JXG.Coords} coords Coords object of an intersection
2129          * @returns {Boolean}
2130          * @private
2131          */
2132         coordsOnArc: function (arc, coords) {
2133             var angle = this.rad(arc.radiuspoint, arc.center, coords.usrCoords.slice(1)),
2134                 alpha = 0.0,
2135                 beta = this.rad(arc.radiuspoint, arc.center, arc.anglepoint),
2136                 ev_s = arc.evalVisProp('selection');
2137 
2138             if (arc.evalVisProp('orientation') === 'clockwise') {
2139                 angle = 2 * Math.PI - angle;
2140                 beta = 2 * Math.PI - beta;
2141             }
2142 
2143             if ((ev_s === "minor" && beta > Math.PI) || (ev_s === "major" && beta < Math.PI)) {
2144                 alpha = beta;
2145                 beta = 2 * Math.PI;
2146             }
2147             if (angle < alpha || angle > beta) {
2148                 return false;
2149             }
2150             return true;
2151         },
2152 
2153         /**
2154          * Computes the intersection of a pair of lines, circles or both.
2155          * It uses the internal data array stdform of these elements.
2156          * @param {Array} el1 stdform of the first element (line or circle)
2157          * @param {Array} el2 stdform of the second element (line or circle)
2158          * @param {Number|Function} i Index of the intersection point that should be returned.
2159          * @param board Reference to the board.
2160          * @returns {JXG.Coords} Coordinates of one of the possible two or more intersection points.
2161          * Which point will be returned is determined by i.
2162          */
2163         meet: function (el1, el2, i, board) {
2164             var result,
2165                 eps = Mat.eps;
2166 
2167             if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) < eps) {
2168                 // line line
2169                 result = this.meetLineLine(el1, el2, i, board);
2170             } else if (Math.abs(el1[3]) >= eps && Math.abs(el2[3]) < eps) {
2171                 // circle line
2172                 result = this.meetLineCircle(el2, el1, i, board);
2173             } else if (Math.abs(el1[3]) < eps && Math.abs(el2[3]) >= eps) {
2174                 // line circle
2175                 result = this.meetLineCircle(el1, el2, i, board);
2176             } else {
2177                 // circle circle
2178                 result = this.meetCircleCircle(el1, el2, i, board);
2179             }
2180 
2181             return result;
2182         },
2183 
2184         /**
2185          * Intersection of the line with the board
2186          * @param  {Array}     line   stdform of the line in screen coordinates
2187          * @param  {JXG.Board} board  reference to a board.
2188          * @param  {Number}    margin optional margin, to avoid the display of the small sides of lines.
2189          * @returns {Array}            [intersection coords 1, intersection coords 2]
2190          */
2191         meetLineBoard: function (line, board, margin) {
2192             // Intersect the line with the four borders of the board.
2193             var s = [],
2194                 intersect1,
2195                 intersect2,
2196                 i, j;
2197 
2198             if (!Type.exists(margin)) {
2199                 margin = 0;
2200             }
2201 
2202             // top
2203             s[0] = Mat.crossProduct(line, [margin, 0, 1]);
2204             // left
2205             s[1] = Mat.crossProduct(line, [margin, 1, 0]);
2206             // bottom
2207             s[2] = Mat.crossProduct(line, [-margin - board.canvasHeight, 0, 1]);
2208             // right
2209             s[3] = Mat.crossProduct(line, [-margin - board.canvasWidth, 1, 0]);
2210 
2211             // Normalize the intersections
2212             for (i = 0; i < 4; i++) {
2213                 if (Math.abs(s[i][0]) > Mat.eps) {
2214                     for (j = 2; j > 0; j--) {
2215                         s[i][j] /= s[i][0];
2216                     }
2217                     s[i][0] = 1.0;
2218                 }
2219             }
2220 
2221             // line is parallel to "left", take "top" and "bottom"
2222             if (Math.abs(s[1][0]) < Mat.eps) {
2223                 intersect1 = s[0]; // top
2224                 intersect2 = s[2]; // bottom
2225                 // line is parallel to "top", take "left" and "right"
2226             } else if (Math.abs(s[0][0]) < Mat.eps) {
2227                 intersect1 = s[1]; // left
2228                 intersect2 = s[3]; // right
2229                 // left intersection out of board (above)
2230             } else if (s[1][2] < 0) {
2231                 intersect1 = s[0]; // top
2232 
2233                 // right intersection out of board (below)
2234                 if (s[3][2] > board.canvasHeight) {
2235                     intersect2 = s[2]; // bottom
2236                 } else {
2237                     intersect2 = s[3]; // right
2238                 }
2239                 // left intersection out of board (below)
2240             } else if (s[1][2] > board.canvasHeight) {
2241                 intersect1 = s[2]; // bottom
2242 
2243                 // right intersection out of board (above)
2244                 if (s[3][2] < 0) {
2245                     intersect2 = s[0]; // top
2246                 } else {
2247                     intersect2 = s[3]; // right
2248                 }
2249             } else {
2250                 intersect1 = s[1]; // left
2251 
2252                 // right intersection out of board (above)
2253                 if (s[3][2] < 0) {
2254                     intersect2 = s[0]; // top
2255                     // right intersection out of board (below)
2256                 } else if (s[3][2] > board.canvasHeight) {
2257                     intersect2 = s[2]; // bottom
2258                 } else {
2259                     intersect2 = s[3]; // right
2260                 }
2261             }
2262 
2263             return [
2264                 new Coords(Const.COORDS_BY_SCREEN, intersect1.slice(1), board),
2265                 new Coords(Const.COORDS_BY_SCREEN, intersect2.slice(1), board)
2266             ];
2267         },
2268 
2269         /**
2270          * Intersection of two lines.
2271          * @param {Array} l1 stdform of the first line
2272          * @param {Array} l2 stdform of the second line
2273          * @param {number} i unused
2274          * @param {JXG.Board} board Reference to the board.
2275          * @returns {JXG.Coords} Coordinates of the intersection point.
2276          */
2277         meetLineLine: function (l1, l2, i, board) {
2278             var s = isNaN(l1[5] + l2[5]) ? [0, 0, 0] : Mat.crossProduct(l1, l2);
2279 
2280             // Make intersection of parallel lines more robust:
2281             if (Math.abs(s[0]) < 1.0e-14) {
2282                 s[0] = 0;
2283             }
2284             return new Coords(Const.COORDS_BY_USER, s, board);
2285         },
2286 
2287         /**
2288          * Intersection of line and circle.
2289          * @param {Array} lin stdform of the line
2290          * @param {Array} circ stdform of the circle
2291          * @param {number|function} i number of the returned intersection point.
2292          *   i==0: use the positive square root,
2293          *   i==1: use the negative square root.
2294          * @param {JXG.Board} board Reference to a board.
2295          * @returns {JXG.Coords} Coordinates of the intersection point
2296          */
2297         meetLineCircle: function (lin, circ, i, board) {
2298             var a, b, c, d, n, A, B, C, k, t;
2299 
2300             // Radius is zero, return center of circle
2301             if (circ[4] < Mat.eps) {
2302                 if (Math.abs(Mat.innerProduct([1, circ[6], circ[7]], lin, 3)) < Mat.eps) {
2303                     return new Coords(Const.COORDS_BY_USER, circ.slice(6, 8), board);
2304                 }
2305 
2306                 return new Coords(Const.COORDS_BY_USER, [NaN, NaN], board);
2307             }
2308             c = circ[0];
2309             b = circ.slice(1, 3);
2310             a = circ[3];
2311             d = lin[0];
2312             n = lin.slice(1, 3);
2313 
2314             // Line is assumed to be normalized. Therefore, nn==1 and we can skip some operations:
2315             /*
2316              var nn = n[0]*n[0]+n[1]*n[1];
2317              A = a*nn;
2318              B = (b[0]*n[1]-b[1]*n[0])*nn;
2319              C = a*d*d - (b[0]*n[0]+b[1]*n[1])*d + c*nn;
2320              */
2321             A = a;
2322             B = b[0] * n[1] - b[1] * n[0];
2323             C = a * d * d - (b[0] * n[0] + b[1] * n[1]) * d + c;
2324 
2325             k = B * B - 4 * A * C;
2326             if (k > -Mat.eps * Mat.eps) {
2327                 k = Math.sqrt(Math.abs(k));
2328                 t = [(-B + k) / (2 * A), (-B - k) / (2 * A)];
2329 
2330                 return Type.evaluate(i) === 0
2331                     ? new Coords(
2332                         Const.COORDS_BY_USER,
2333                         [-t[0] * -n[1] - d * n[0], -t[0] * n[0] - d * n[1]],
2334                         board
2335                     )
2336                     : new Coords(
2337                         Const.COORDS_BY_USER,
2338                         [-t[1] * -n[1] - d * n[0], -t[1] * n[0] - d * n[1]],
2339                         board
2340                     );
2341             }
2342 
2343             return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
2344         },
2345 
2346         /**
2347          * Intersection of two circles.
2348          * @param {Array} circ1 stdform of the first circle
2349          * @param {Array} circ2 stdform of the second circle
2350          * @param {number|function} i number of the returned intersection point.
2351          *   i==0: use the positive square root,
2352          *   i==1: use the negative square root.
2353          * @param {JXG.Board} board Reference to the board.
2354          * @returns {JXG.Coords} Coordinates of the intersection point
2355          */
2356         meetCircleCircle: function (circ1, circ2, i, board) {
2357             var radicalAxis;
2358 
2359             // Radius is zero, return center of circle, if on other circle
2360             if (circ1[4] < Mat.eps) {
2361                 if (
2362                     Math.abs(this.distance(circ1.slice(6, 2), circ2.slice(6, 8)) - circ2[4]) <
2363                     Mat.eps
2364                 ) {
2365                     return new Coords(Const.COORDS_BY_USER, circ1.slice(6, 8), board);
2366                 }
2367 
2368                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
2369             }
2370 
2371             // Radius is zero, return center of circle, if on other circle
2372             if (circ2[4] < Mat.eps) {
2373                 if (
2374                     Math.abs(this.distance(circ2.slice(6, 2), circ1.slice(6, 8)) - circ1[4]) <
2375                     Mat.eps
2376                 ) {
2377                     return new Coords(Const.COORDS_BY_USER, circ2.slice(6, 8), board);
2378                 }
2379 
2380                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
2381             }
2382 
2383             radicalAxis = [
2384                 circ2[3] * circ1[0] - circ1[3] * circ2[0],
2385                 circ2[3] * circ1[1] - circ1[3] * circ2[1],
2386                 circ2[3] * circ1[2] - circ1[3] * circ2[2],
2387                 0,
2388                 1,
2389                 Infinity,
2390                 Infinity,
2391                 Infinity
2392             ];
2393             radicalAxis = Mat.normalize(radicalAxis);
2394 
2395             return this.meetLineCircle(radicalAxis, circ1, i, board);
2396         },
2397 
2398         /**
2399          * Segment-wise search for the nr-th intersection of two curves.
2400          * testSegment is always assumed to be true.
2401          *
2402          * @param {JXG.Curve} c1 Curve, Line or Circle
2403          * @param {JXG.Curve} c2 Curve, Line or Circle
2404          * @param {Number} nr the nr-th intersection point will be returned
2405          * @param {JXG.Board} [board=c1.board] Reference to a board object
2406          * @returns {JXG.Coords} intersection as Coords object
2407          *
2408          * @private
2409          * @see JXG.Math.Geometry.meetCurveCurve
2410          */
2411         meetCurveCurveDiscrete: function (c1, c2, nr, board) {
2412             var co,
2413                 i = Type.evaluate(nr);
2414 
2415             if (c1.bezierDegree === 3 || c2.bezierDegree === 3) {
2416                 co = this.meetBezierCurveRedBlueSegments(c1, c2, i);
2417             } else {
2418                 co = this.meetCurveRedBlueSegments(c1, c2, i);
2419             }
2420             return new Coords(Const.COORDS_BY_USER, co, board);
2421         },
2422 
2423         /**
2424          * Apply Newton-Raphson to search for an intersection of two curves
2425          * in a given range of the first curve.
2426          *
2427          * @param {JXG.Curve} c1 Curve, Line or Circle
2428          * @param {JXG.Curve} c2 Curve, Line or Circle
2429          * @param {Array} range Domain for the search of an intersection. The start value
2430          * for the search is chosen to be inside of that range.
2431          * @param {Boolean} testSegment If true require that t1 and t2 are inside of the allowed bounds.
2432          * @returns {Array} [[z, x, y], t1, t2, t, ||c1[t1]-c2[t2]||**2]. The last entry is set to
2433          * 10000 if the intersection is outside of the given domain (range) for the first curve.
2434          * @private
2435          * @see JXG.Math.Geometry._meetCurveCurveRecursive
2436          * @see JXG.Math.Geometry._meetCurveCurveIterative
2437          * @see JXG.Math.Numerics.generalizedDampedNewton
2438          * @see JXG.Math.Geometry.meetCurveCurveCobyla
2439          */
2440         meetCurveCurveNewton: function (c1, c2, range1, range2, testSegment) {
2441             var t1, t2,
2442                 co, r,
2443                 inphi = (Math.sqrt(5) - 1) * 0.5,
2444                 damp = 0.85, // (
2445                 eps3 = Mat.eps * Mat.eps * Mat.eps,
2446                 eps2 = Mat.eps * Mat.eps,
2447 
2448                 ma1 = c1.maxX(),
2449                 mi1 = c1.minX(),
2450                 ma2 = c2.maxX(),
2451                 mi2 = c2.minX(),
2452 
2453                 F = function(t, n) {
2454                     var f1 = c1.Ft(t[0]),
2455                         f2 = c2.Ft(t[1]),
2456                         e = f1[1] - f2[1],
2457                         f = f1[2] - f2[2];
2458 
2459                     return [e, f];
2460                 },
2461                 D = function(t, n) {
2462                     var h = Mat.eps,
2463                         h2 = 2 * h,
2464                         f1_1 = c1.Ft(t[0] - h),
2465                         f1_2 = c1.Ft(t[0] + h),
2466                         f2_1 = c2.Ft(t[1] - h),
2467                         f2_2 = c2.Ft(t[1] + h);
2468                     return [
2469                         [ (f1_2[1] - f1_1[1]) / h2,
2470                          -(f2_2[1] - f2_1[1]) / h2],
2471                         [ (f1_2[2] - f1_1[2]) / h2,
2472                          -(f2_2[2] - f2_1[2]) / h2]
2473                     ];
2474                 };
2475 
2476             t1 = range1[0] + (range1[1] - range1[0]) * (1 - inphi);
2477             t2 = range2[0] + (range2[1] - range2[0]) * (1 - inphi);
2478 
2479             // Use damped Newton
2480             // r = Numerics.generalizedDampedNewtonCurves(c1, c2, t1, t2, damp, eps3);
2481             r = Numerics.generalizedDampedNewton(F, D, 2, [t1, t2], damp, eps3, 40);
2482             // r: [t1, t2, F2]
2483 
2484             t1 = r[0][0];
2485             t2 = r[0][1];
2486             co = c1.Ft(t1);
2487 
2488             if (
2489                 t1 < range1[0] - Mat.eps || t1 > range1[1] + Mat.eps ||
2490                 t2 < range2[0] - Mat.eps || t2 > range2[1] + Mat.eps ||
2491                 (testSegment &&
2492                     (t1 < mi1 - eps2 || t1 > ma1 + eps2 ||
2493                      t2 < mi2 - eps2 || t2 > ma2 + eps2)
2494                 )
2495             ) {
2496                 // Damped-Newton found solution outside of range
2497                 return [co, t1, t2, 10000];
2498             }
2499 // console.log(t1, r[3])
2500 
2501             return [co, t1, t2, r[1]];
2502         },
2503 
2504         /**
2505          * Return a list of the (at most) first i intersection points of two curves.
2506          * Computed iteratively.
2507          *
2508          * @param {JXG.Curve} c1 Curve, Line or Circle
2509          * @param {JXG.Curve} c2 Curve, Line or Circle
2510          * @param {Number} low Lower bound of the search domain (between [0, 1])
2511          * @param {Number} up Upper bound of the search domain (between [0, 1])
2512          * @param {Number} i Return a list of the first i intersection points
2513          * @param {Boolean} testSegment If true require that t1 and t2 are inside of the allowed bounds.
2514          * @returns {Array} List of the first i intersection points, given by the parameter t.
2515          * @private
2516          * @see JXG.Math.Geometry.meetCurveCurveNewton
2517          * @see JXG.Math.Geometry.meetCurveCurve
2518          */
2519         _meetCurveCurveIterative: function(c1, c2, range1, range2, i, testSegment) {
2520             var ret,
2521                 t1,// t2,
2522                 low1, low2, up1, up2,
2523                 eps = Mat.eps * 100, // Minimum difference between zeros
2524                 j1, j2,
2525                 steps = 20,
2526                 d1, d2,
2527                 zeros = [];
2528 
2529             low1 = range1[0];
2530             up1 = range1[1];
2531             low2 = range2[0];
2532             up2 = range2[1];
2533             if (up1 < low1 || up2 < low2) {
2534                 return [];
2535             }
2536 
2537             // console.log('DO iterative', [low1, up1], [low2, up2])
2538 
2539             d1 = (up1 - low1) / steps;
2540             d2 = (up2 - low2) / steps;
2541             for (j1 = 0; j1 < steps; j1++) {
2542                 for (j2 = 0; j2 < steps; j2++) {
2543 
2544                     ret = this.meetCurveCurveNewton(c1, c2,
2545                         [low1 + j1 * d1, low1 + (j1 + 1) * d1],
2546                         [low2 + j2 * d2, low2 + (j2 + 1) * d2],
2547                         testSegment);
2548 
2549                     if (ret[3] < Mat.eps) {
2550                         t1 = ret[1];
2551                         // t2 = ret[2];
2552                         // console.log("\tFOUND", t1, t2, c1.Ft(t1)[2])
2553                         zeros = zeros.concat([t1]);
2554                         zeros = Type.toUniqueArrayFloat(zeros, eps);
2555                         // console.log(zeros, i)
2556                         if (zeros.length > i) {
2557                             return zeros;
2558                         }
2559                     }
2560                 }
2561             }
2562 
2563             return zeros;
2564         },
2565 
2566         /**
2567          * Compute an intersection of the curves c1 and c2.
2568          * We want to find values t1, t2 such that
2569          * c1(t1) = c2(t2), i.e. (c1_x(t1) - c2_x(t2), c1_y(t1) - c2_y(t2)) = (0, 0).
2570          *
2571          * Available methods:
2572          * <ul>
2573          *  <li> discrete, segment-wise intersections
2574          *  <li> generalized damped Newton-Raphson
2575          * </ul>
2576          *
2577          * Segment-wise intersection is more stable, but has problems with tangent points.
2578          * Damped Newton-Raphson converges very rapidly but sometimes behaves chaotic.
2579          *
2580          * @param {JXG.Curve} c1 Curve, Line or Circle
2581          * @param {JXG.Curve} c2 Curve, Line or Circle
2582          * @param {Number|Function} nr the nr-th intersection point will be returned. For backwards compatibility:
2583          * if method='newton' and nr is not an integer, {@link JXG.Math.Numerics.generalizedNewton} is called
2584          * directly with nr as start value (not recommended).
2585          * @param {Number} t2ini not longer used. Must be supplied and is ignored.
2586          * @param {JXG.Board} [board=c1.board] Reference to a board object.
2587          * @param {String} [method] Intersection method, possible values are 'newton' and 'segment'.
2588          * If both curves are given by functions (assumed to be continuous), 'newton' is the default, otherwise
2589          * 'segment' is the default.
2590          * @returns {JXG.Coords} intersection point
2591          *
2592          * @see JXG.Math.Geometry.meetCurveCurveDiscrete
2593          * @see JXG.Math.Geometry._meetCurveCurveRecursive
2594          * @see JXG.Math.Geometry.meetCurveCurveIterative
2595          */
2596         meetCurveCurve: function (c1, c2, nr, t2ini, board, method, testSegment) {
2597             var co,
2598                 zeros,
2599                 mi1, ma1, mi2, ma2,
2600                 i = Type.evaluate(nr);
2601 
2602             board = board || c1.board;
2603             if (method === 'segment' || Type.exists(c1.dataX) || Type.exists(c2.dataX)) {
2604                 // Discrete data points, i.e. x-coordinates of c1 or c2 are given in an array)
2605                 return this.meetCurveCurveDiscrete(c1, c2, i, board);
2606             }
2607 
2608             // Outdated:
2609             // Backwards compatibility if nr is not a positive integer then
2610             // generalizedNewton is still used.
2611             if (Type.exists(method) && method === 'newton' && i < 0 || parseInt(i) !== i) {
2612                 co = Numerics.generalizedNewton(c1, c2, i, t2ini);
2613                 return new Coords(Const.COORDS_BY_USER, co, board);
2614             }
2615 
2616             // Method 'newton'
2617             mi1 = c1.minX();
2618             ma1 = c1.maxX();
2619             mi2 = c2.minX();
2620             ma2 = c2.maxX();
2621 
2622             // console.time('curvecurve')
2623             zeros = this._meetCurveCurveIterative(c1, c2, [mi1, ma1], [mi2, ma2], i, testSegment);
2624             // console.timeEnd('curvecurve')
2625 
2626             if (zeros.length > i) {
2627                 co = c1.Ft(zeros[i]);
2628             } else {
2629                 return [0, NaN, NaN];
2630             }
2631 
2632             return new Coords(Const.COORDS_BY_USER, co, board);
2633         },
2634 
2635         /**
2636          * Intersection of curve with line,
2637          * Order of input does not matter for el1 and el2.
2638          * From version 0.99.7 on this method calls
2639          * {@link JXG.Math.Geometry.meetCurveLineDiscrete}.
2640          * If higher precision is needed, {@link JXG.Math.Geometry.meetCurveLineContinuous}
2641          * has to be used.
2642          *
2643          * @param {JXG.Curve|JXG.Line} el1 Curve or Line
2644          * @param {JXG.Curve|JXG.Line} el2 Curve or Line
2645          * @param {Number|Function} nr the nr-th intersection point will be returned.
2646          * @param {JXG.Board} [board=el1.board] Reference to a board object.
2647          * @param {Boolean} alwaysIntersect If false just the segment between the two defining points are tested for intersection
2648          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2649          * the ideal point [0,1,0] is returned.
2650          */
2651         meetCurveLine: function (el1, el2, nr, board, alwaysIntersect) {
2652             var v = [0, NaN, NaN],
2653                 cu,
2654                 li;
2655 
2656             if (!Type.exists(board)) {
2657                 board = el1.board;
2658             }
2659 
2660             if (el1.elementClass === Const.OBJECT_CLASS_CURVE) {
2661                 cu = el1;
2662                 li = el2;
2663             } else {
2664                 cu = el2;
2665                 li = el1;
2666             }
2667 
2668             if (Type.exists(cu.dataX)) {
2669                 // We use the discrete version if
2670                 //   the curve is not a parametric curve, e.g. implicit plots
2671                 v = this.meetCurveLineDiscrete(cu, li, nr, board, !alwaysIntersect);
2672             } else {
2673                 v = this.meetCurveCurve(cu, li, nr, 0, board, 'newton', !alwaysIntersect);
2674             }
2675 
2676             return v;
2677         },
2678 
2679         /**
2680          * Intersection of line and curve, continuous case.
2681          * Finds the nr-th intersection point
2682          * Uses {@link JXG.Math.Geometry.meetCurveLineDiscrete} as a first approximation.
2683          * A more exact solution is then found with {@link JXG.Math.Numerics.root}.
2684          *
2685          * @param {JXG.Curve} cu Curve
2686          * @param {JXG.Line} li Line
2687          * @param {NumberFunction} nr Will return the nr-th intersection point.
2688          * @param {JXG.Board} board
2689          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the
2690          * line defined by the segment
2691          * @returns {JXG.Coords} Coords object containing the intersection.
2692          */
2693         meetCurveLineContinuous: function (cu, li, nr, board, testSegment) {
2694             var func0, func1,
2695                 t, v, x, y, z,
2696                 eps = Mat.eps,
2697                 epsLow = Mat.eps,
2698                 steps,
2699                 delta,
2700                 tnew, tmin, fmin,
2701                 i, ft;
2702 
2703             v = this.meetCurveLineDiscrete(cu, li, nr, board, testSegment);
2704             x = v.usrCoords[1];
2705             y = v.usrCoords[2];
2706 
2707             func0 = function (t) {
2708                 var c1, c2;
2709 
2710                 if (t > cu.maxX() || t < cu.minX()) {
2711                     return Infinity;
2712                 }
2713                 c1 = cu.X(t) - x;
2714                 c2 = cu.Y(t) - y;
2715                 return c1 * c1 + c2 * c2;
2716                 // return c1 * (cu.X(t + h) - cu.X(t - h)) + c2 * (cu.Y(t + h) - cu.Y(t - h)) / h;
2717             };
2718 
2719             func1 = function (t) {
2720                 var v = li.stdform[0] + li.stdform[1] * cu.X(t) + li.stdform[2] * cu.Y(t);
2721                 return v * v;
2722             };
2723 
2724             // Find t
2725             steps = 50;
2726             delta = (cu.maxX() - cu.minX()) / steps;
2727             tnew = cu.minX();
2728             fmin = 0.0001; //eps;
2729             tmin = NaN;
2730             for (i = 0; i < steps; i++) {
2731                 t = Numerics.root(func0, [
2732                     Math.max(tnew, cu.minX()),
2733                     Math.min(tnew + delta, cu.maxX())
2734                 ]);
2735                 ft = Math.abs(func0(t));
2736                 if (ft <= fmin) {
2737                     fmin = ft;
2738                     tmin = t;
2739                     if (fmin < eps) {
2740                         break;
2741                     }
2742                 }
2743 
2744                 tnew += delta;
2745             }
2746             t = tmin;
2747             // Compute "exact" t
2748             t = Numerics.root(func1, [
2749                 Math.max(t - delta, cu.minX()),
2750                 Math.min(t + delta, cu.maxX())
2751             ]);
2752 
2753             ft = func1(t);
2754             // Is the point on the line?
2755             if (isNaN(ft) || Math.abs(ft) > epsLow) {
2756                 z = 0.0; //NaN;
2757             } else {
2758                 z = 1.0;
2759             }
2760 
2761             return new Coords(Const.COORDS_BY_USER, [z, cu.X(t), cu.Y(t)], board);
2762         },
2763 
2764         /**
2765          * Intersection of line and curve, discrete case.
2766          * Segments are treated as lines.
2767          * Finding the nr-th intersection point should work for all nr.
2768          * @param {JXG.Curve} cu
2769          * @param {JXG.Line} li
2770          * @param {Number|Function} nr
2771          * @param {JXG.Board} board
2772          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the
2773          * line defined by the segment
2774          *
2775          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2776          * the ideal point [0,1,0] is returned.
2777          */
2778         meetCurveLineDiscrete: function (cu, li, nr, board, testSegment) {
2779             var i, j,
2780                 n = Type.evaluate(nr),
2781                 p1, p2,
2782                 p, q,
2783                 lip1 = li.point1.coords.usrCoords,
2784                 lip2 = li.point2.coords.usrCoords,
2785                 d, res,
2786                 cnt = 0,
2787                 len = cu.numberPoints,
2788                 ev_sf = li.evalVisProp('straightfirst'),
2789                 ev_sl = li.evalVisProp('straightlast');
2790 
2791             // In case, no intersection will be found we will take this
2792             q = new Coords(Const.COORDS_BY_USER, [0, NaN, NaN], board);
2793 
2794             if (lip1[0] === 0.0) {
2795                 lip1 = [1, lip2[1] + li.stdform[2], lip2[2] - li.stdform[1]];
2796             } else if (lip2[0] === 0.0) {
2797                 lip2 = [1, lip1[1] + li.stdform[2], lip1[2] - li.stdform[1]];
2798             }
2799 
2800             p2 = cu.points[0].usrCoords;
2801             for (i = 1; i < len; i += cu.bezierDegree) {
2802                 p1 = p2.slice(0);
2803                 p2 = cu.points[i].usrCoords;
2804                 d = this.distance(p1, p2);
2805 
2806                 // The defining points are not identical
2807                 if (d > Mat.eps) {
2808                     if (cu.bezierDegree === 3) {
2809                         res = this.meetBeziersegmentBeziersegment(
2810                             [
2811                                 cu.points[i - 1].usrCoords.slice(1),
2812                                 cu.points[i].usrCoords.slice(1),
2813                                 cu.points[i + 1].usrCoords.slice(1),
2814                                 cu.points[i + 2].usrCoords.slice(1)
2815                             ],
2816                             [lip1.slice(1), lip2.slice(1)],
2817                             testSegment
2818                         );
2819                     } else {
2820                         res = [this.meetSegmentSegment(p1, p2, lip1, lip2)];
2821                     }
2822 
2823                     for (j = 0; j < res.length; j++) {
2824                         p = res[j];
2825                         if (0 <= p[1] && p[1] <= 1) {
2826                             if (cnt === n) {
2827                                 /**
2828                                  * If the intersection point is not part of the segment,
2829                                  * this intersection point is set to non-existent.
2830                                  * This prevents jumping behavior of the intersection points.
2831                                  * But it may be discussed if it is the desired behavior.
2832                                  */
2833                                 if (
2834                                     testSegment &&
2835                                     ((!ev_sf && p[2] < 0) || (!ev_sl && p[2] > 1))
2836                                 ) {
2837                                     return q; // break;
2838                                 }
2839 
2840                                 q = new Coords(Const.COORDS_BY_USER, p[0], board);
2841                                 return q; // break;
2842                             }
2843                             cnt += 1;
2844                         }
2845                     }
2846                 }
2847             }
2848 
2849             return q;
2850         },
2851 
2852         /**
2853          * Find the n-th intersection point of two curves named red (first parameter) and blue (second parameter).
2854          * We go through each segment of the red curve and search if there is an intersection with a segment of the blue curve.
2855          * This double loop, i.e. the outer loop runs along the red curve and the inner loop runs along the blue curve, defines
2856          * the n-th intersection point. The segments are either line segments or Bezier curves of degree 3. This depends on
2857          * the property bezierDegree of the curves.
2858          * <p>
2859          * This method works also for transformed curves, since only the already
2860          * transformed points are used.
2861          *
2862          * @param {JXG.Curve} red
2863          * @param {JXG.Curve} blue
2864          * @param {Number|Function} nr
2865          */
2866         meetCurveRedBlueSegments: function (red, blue, nr) {
2867             var i,
2868                 j,
2869                 n = Type.evaluate(nr),
2870                 red1,
2871                 red2,
2872                 blue1,
2873                 blue2,
2874                 m,
2875                 minX,
2876                 maxX,
2877                 iFound = 0,
2878                 lenBlue = blue.numberPoints,
2879                 lenRed = red.numberPoints;
2880 
2881             if (lenBlue <= 1 || lenRed <= 1) {
2882                 return [0, NaN, NaN];
2883             }
2884 
2885             for (i = 1; i < lenRed; i++) {
2886                 red1 = red.points[i - 1].usrCoords;
2887                 red2 = red.points[i].usrCoords;
2888                 minX = Math.min(red1[1], red2[1]);
2889                 maxX = Math.max(red1[1], red2[1]);
2890 
2891                 blue2 = blue.points[0].usrCoords;
2892                 for (j = 1; j < lenBlue; j++) {
2893                     blue1 = blue2;
2894                     blue2 = blue.points[j].usrCoords;
2895                     if (
2896                         Math.min(blue1[1], blue2[1]) < maxX &&
2897                         Math.max(blue1[1], blue2[1]) > minX
2898                     ) {
2899                         m = this.meetSegmentSegment(red1, red2, blue1, blue2);
2900                         if (
2901                             m[1] >= 0.0 && m[2] >= 0.0 &&
2902                             // The two segments meet in the interior or at the start points
2903                             ((m[1] < 1.0 && m[2] < 1.0) ||
2904                               // One of the curve is intersected in the very last point
2905                                 (i === lenRed - 1 && m[1] === 1.0) ||
2906                                 (j === lenBlue - 1 && m[2] === 1.0))
2907                         ) {
2908                             if (iFound === n) {
2909                                 return m[0];
2910                             }
2911 
2912                             iFound++;
2913                         }
2914                     }
2915                 }
2916             }
2917 
2918             return [0, NaN, NaN];
2919         },
2920 
2921         /**
2922          * (Virtual) Intersection of two segments.
2923          * @param {Array} p1 First point of segment 1 using normalized homogeneous coordinates [1,x,y]
2924          * @param {Array} p2 Second point or direction of segment 1 using normalized homogeneous coordinates [1,x,y] or point at infinity [0,x,y], respectively
2925          * @param {Array} q1 First point of segment 2 using normalized homogeneous coordinates [1,x,y]
2926          * @param {Array} q2 Second point or direction of segment 2 using normalized homogeneous coordinates [1,x,y] or point at infinity [0,x,y], respectively
2927          * @returns {Array} [Intersection point, t, u] The first entry contains the homogeneous coordinates
2928          * of the intersection point. The second and third entry give the position of the intersection with respect
2929          * to the definiting parameters. For example, the second entry t is defined by: intersection point = p1 + t * deltaP, where
2930          * deltaP = (p2 - p1) when both parameters are coordinates, and deltaP = p2 if p2 is a point at infinity.
2931          * If the two segments are collinear, [[0,0,0], Infinity, Infinity] is returned.
2932          **/
2933         meetSegmentSegment: function (p1, p2, q1, q2) {
2934             var t,
2935                 u,
2936                 i,
2937                 d,
2938                 li1 = Mat.crossProduct(p1, p2),
2939                 li2 = Mat.crossProduct(q1, q2),
2940                 c = Mat.crossProduct(li1, li2);
2941 
2942             if (Math.abs(c[0]) < Mat.eps) {
2943                 return [c, Infinity, Infinity];
2944             }
2945 
2946             // Normalize the intersection coordinates
2947             c[1] /= c[0];
2948             c[2] /= c[0];
2949             c[0] /= c[0];
2950 
2951             // Now compute in principle:
2952             //    t = dist(c - p1) / dist(p2 - p1) and
2953             //    u = dist(c - q1) / dist(q2 - q1)
2954             // However: the points q1, q2, p1, p2 might be ideal points - or in general - the
2955             // coordinates might be not normalized.
2956             // Note that the z-coordinates of p2 and q2 are used to determine whether it should be interpreted
2957             // as a segment coordinate or a direction.
2958             i = Math.abs(p2[1] - p2[0] * p1[1]) < Mat.eps ? 2 : 1;
2959             d = p1[i] / p1[0];
2960             t = (c[i] - d) / (p2[0] !== 0 ? p2[i] / p2[0] - d : p2[i]);
2961 
2962             i = Math.abs(q2[1] - q2[0] * q1[1]) < Mat.eps ? 2 : 1;
2963             d = q1[i] / q1[0];
2964             u = (c[i] - d) / (q2[0] !== 0 ? q2[i] / q2[0] - d : q2[i]);
2965 
2966             return [c, t, u];
2967         },
2968 
2969         /**
2970          * Find the n-th intersection point of two pathes, usually given by polygons. Uses parts of the
2971          * Greiner-Hormann algorithm in JXG.Math.Clip.
2972          *
2973          * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path1
2974          * @param {JXG.Circle|JXG.Curve|JXG.Polygon} path2
2975          * @param {Number|Function} n
2976          * @param {JXG.Board} board
2977          *
2978          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
2979          * the ideal point [0,0,0] is returned.
2980          *
2981          */
2982         meetPathPath: function (path1, path2, nr, board) {
2983             var S, C, len, intersections,
2984                 n = Type.evaluate(nr);
2985 
2986             S = JXG.Math.Clip._getPath(path1, board);
2987             len = S.length;
2988             if (
2989                 len > 0 &&
2990                 this.distance(S[0].coords.usrCoords, S[len - 1].coords.usrCoords, 3) < Mat.eps
2991             ) {
2992                 S.pop();
2993             }
2994 
2995             C = JXG.Math.Clip._getPath(path2, board);
2996             len = C.length;
2997             if (
2998                 len > 0 &&
2999                 this.distance(C[0].coords.usrCoords, C[len - 1].coords.usrCoords, 3) <
3000                 Mat.eps * Mat.eps
3001             ) {
3002                 C.pop();
3003             }
3004 
3005             // Handle cases where at least one of the paths is empty
3006             if (nr < 0 || JXG.Math.Clip.isEmptyCase(S, C, 'intersection')) {
3007                 return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
3008             }
3009 
3010             JXG.Math.Clip.makeDoublyLinkedList(S);
3011             JXG.Math.Clip.makeDoublyLinkedList(C);
3012 
3013             intersections = JXG.Math.Clip.findIntersections(S, C, board)[0];
3014             if (n < intersections.length) {
3015                 return intersections[n].coords;
3016             }
3017             return new Coords(Const.COORDS_BY_USER, [0, 0, 0], board);
3018         },
3019 
3020         /**
3021          * Find the n-th intersection point between a polygon and a line.
3022          * @param {JXG.Polygon} path
3023          * @param {JXG.Line} line
3024          * @param {Number|Function} nr
3025          * @param {JXG.Board} board
3026          * @param {Boolean} alwaysIntersect If false just the segment between the two defining points of the line are tested for intersection.
3027          *
3028          * @returns {JXG.Coords} Intersection point. In case no intersection point is detected,
3029          * the ideal point [0,0,0] is returned.
3030          */
3031         meetPolygonLine: function (path, line, nr, board, alwaysIntersect) {
3032             var i,
3033                 n = Type.evaluate(nr),
3034                 res,
3035                 border,
3036                 crds = [0, 0, 0],
3037                 len = path.borders.length,
3038                 intersections = [];
3039 
3040             for (i = 0; i < len; i++) {
3041                 border = path.borders[i];
3042                 res = this.meetSegmentSegment(
3043                     border.point1.coords.usrCoords,
3044                     border.point2.coords.usrCoords,
3045                     line.point1.coords.usrCoords,
3046                     line.point2.coords.usrCoords
3047                 );
3048 
3049                 if (
3050                     (!alwaysIntersect || (res[2] >= 0 && res[2] < 1)) &&
3051                     res[1] >= 0 &&
3052                     res[1] < 1
3053                 ) {
3054                     intersections.push(res[0]);
3055                 }
3056             }
3057 
3058             if (n >= 0 && n < intersections.length) {
3059                 crds = intersections[n];
3060             }
3061             return new Coords(Const.COORDS_BY_USER, crds, board);
3062         },
3063 
3064         /****************************************/
3065         /****   BEZIER CURVE ALGORITHMS      ****/
3066         /****************************************/
3067 
3068         /**
3069          * Splits a Bezier curve segment defined by four points into
3070          * two Bezier curve segments. Dissection point is t=1/2.
3071          * @param {Array} curve Array of four coordinate arrays of length 2 defining a
3072          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3073          * @returns {Array} Array consisting of two coordinate arrays for Bezier curves.
3074          */
3075         _bezierSplit: function (curve) {
3076             var p0, p1, p2, p00, p22, p000;
3077 
3078             p0 = [(curve[0][0] + curve[1][0]) * 0.5, (curve[0][1] + curve[1][1]) * 0.5];
3079             p1 = [(curve[1][0] + curve[2][0]) * 0.5, (curve[1][1] + curve[2][1]) * 0.5];
3080             p2 = [(curve[2][0] + curve[3][0]) * 0.5, (curve[2][1] + curve[3][1]) * 0.5];
3081 
3082             p00 = [(p0[0] + p1[0]) * 0.5, (p0[1] + p1[1]) * 0.5];
3083             p22 = [(p1[0] + p2[0]) * 0.5, (p1[1] + p2[1]) * 0.5];
3084 
3085             p000 = [(p00[0] + p22[0]) * 0.5, (p00[1] + p22[1]) * 0.5];
3086 
3087             return [
3088                 [curve[0], p0, p00, p000],
3089                 [p000, p22, p2, curve[3]]
3090             ];
3091         },
3092 
3093         /**
3094          * Computes the bounding box [minX, maxY, maxX, minY] of a Bezier curve segment
3095          * from its control points.
3096          * @param {Array} curve Array of four coordinate arrays of length 2 defining a
3097          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3098          * @returns {Array} Bounding box [minX, maxY, maxX, minY]
3099          */
3100         _bezierBbox: function (curve) {
3101             var bb = [];
3102 
3103             if (curve.length === 4) {
3104                 // bezierDegree == 3
3105                 bb[0] = Math.min(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // minX
3106                 bb[1] = Math.max(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // maxY
3107                 bb[2] = Math.max(curve[0][0], curve[1][0], curve[2][0], curve[3][0]); // maxX
3108                 bb[3] = Math.min(curve[0][1], curve[1][1], curve[2][1], curve[3][1]); // minY
3109             } else {
3110                 // bezierDegree == 1
3111                 bb[0] = Math.min(curve[0][0], curve[1][0]); // minX
3112                 bb[1] = Math.max(curve[0][1], curve[1][1]); // maxY
3113                 bb[2] = Math.max(curve[0][0], curve[1][0]); // maxX
3114                 bb[3] = Math.min(curve[0][1], curve[1][1]); // minY
3115             }
3116 
3117             return bb;
3118         },
3119 
3120         /**
3121          * Decide if two Bezier curve segments overlap by comparing their bounding boxes.
3122          * @param {Array} bb1 Bounding box of the first Bezier curve segment
3123          * @param {Array} bb2 Bounding box of the second Bezier curve segment
3124          * @returns {Boolean} true if the bounding boxes overlap, false otherwise.
3125          */
3126         _bezierOverlap: function (bb1, bb2) {
3127             return bb1[2] >= bb2[0] && bb1[0] <= bb2[2] && bb1[1] >= bb2[3] && bb1[3] <= bb2[1];
3128         },
3129 
3130         /**
3131          * Append list of intersection points to a list.
3132          * @private
3133          */
3134         _bezierListConcat: function (L, Lnew, t1, t2) {
3135             var i,
3136                 t2exists = Type.exists(t2),
3137                 start = 0,
3138                 len = Lnew.length,
3139                 le = L.length;
3140 
3141             if (
3142                 le > 0 &&
3143                 len > 0 &&
3144                 ((L[le - 1][1] === 1 && Lnew[0][1] === 0) ||
3145                     (t2exists && L[le - 1][2] === 1 && Lnew[0][2] === 0))
3146             ) {
3147                 start = 1;
3148             }
3149 
3150             for (i = start; i < len; i++) {
3151                 if (t2exists) {
3152                     Lnew[i][2] *= 0.5;
3153                     Lnew[i][2] += t2;
3154                 }
3155 
3156                 Lnew[i][1] *= 0.5;
3157                 Lnew[i][1] += t1;
3158 
3159                 L.push(Lnew[i]);
3160             }
3161         },
3162 
3163         /**
3164          * Find intersections of two Bezier curve segments by recursive subdivision.
3165          * Below maxlevel determine intersections by intersection line segments.
3166          * @param {Array} red Array of four coordinate arrays of length 2 defining the first
3167          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3168          * @param {Array} blue Array of four coordinate arrays of length 2 defining the second
3169          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3170          * @param {Number} level Recursion level
3171          * @returns {Array} List of intersection points (up to nine). Each intersection point is an
3172          * array of length three (homogeneous coordinates) plus preimages.
3173          */
3174         _bezierMeetSubdivision: function (red, blue, level) {
3175             var bbb,
3176                 bbr,
3177                 ar,
3178                 b0,
3179                 b1,
3180                 r0,
3181                 r1,
3182                 m,
3183                 p0,
3184                 p1,
3185                 q0,
3186                 q1,
3187                 L = [],
3188                 maxLev = 5; // Maximum recursion level
3189 
3190             bbr = this._bezierBbox(blue);
3191             bbb = this._bezierBbox(red);
3192 
3193             if (!this._bezierOverlap(bbr, bbb)) {
3194                 return [];
3195             }
3196 
3197             if (level < maxLev) {
3198                 ar = this._bezierSplit(red);
3199                 r0 = ar[0];
3200                 r1 = ar[1];
3201 
3202                 ar = this._bezierSplit(blue);
3203                 b0 = ar[0];
3204                 b1 = ar[1];
3205 
3206                 this._bezierListConcat(
3207                     L,
3208                     this._bezierMeetSubdivision(r0, b0, level + 1),
3209                     0.0,
3210                     0.0
3211                 );
3212                 this._bezierListConcat(
3213                     L,
3214                     this._bezierMeetSubdivision(r0, b1, level + 1),
3215                     0,
3216                     0.5
3217                 );
3218                 this._bezierListConcat(
3219                     L,
3220                     this._bezierMeetSubdivision(r1, b0, level + 1),
3221                     0.5,
3222                     0.0
3223                 );
3224                 this._bezierListConcat(
3225                     L,
3226                     this._bezierMeetSubdivision(r1, b1, level + 1),
3227                     0.5,
3228                     0.5
3229                 );
3230 
3231                 return L;
3232             }
3233 
3234             // Make homogeneous coordinates
3235             q0 = [1].concat(red[0]);
3236             q1 = [1].concat(red[3]);
3237             p0 = [1].concat(blue[0]);
3238             p1 = [1].concat(blue[3]);
3239 
3240             m = this.meetSegmentSegment(q0, q1, p0, p1);
3241 
3242             if (m[1] >= 0.0 && m[2] >= 0.0 && m[1] <= 1.0 && m[2] <= 1.0) {
3243                 return [m];
3244             }
3245 
3246             return [];
3247         },
3248 
3249         /**
3250          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment
3251          */
3252         _bezierLineMeetSubdivision: function (red, blue, level, testSegment) {
3253             var bbb, bbr, ar,
3254                 r0, r1,
3255                 m,
3256                 p0, p1, q0, q1,
3257                 L = [],
3258                 maxLev = 5; // Maximum recursion level
3259 
3260             bbb = this._bezierBbox(blue);
3261             bbr = this._bezierBbox(red);
3262 
3263             if (testSegment && !this._bezierOverlap(bbr, bbb)) {
3264                 return [];
3265             }
3266 
3267             if (level < maxLev) {
3268                 ar = this._bezierSplit(red);
3269                 r0 = ar[0];
3270                 r1 = ar[1];
3271 
3272                 this._bezierListConcat(
3273                     L,
3274                     this._bezierLineMeetSubdivision(r0, blue, level + 1),
3275                     0.0
3276                 );
3277                 this._bezierListConcat(
3278                     L,
3279                     this._bezierLineMeetSubdivision(r1, blue, level + 1),
3280                     0.5
3281                 );
3282 
3283                 return L;
3284             }
3285 
3286             // Make homogeneous coordinates
3287             q0 = [1].concat(red[0]);
3288             q1 = [1].concat(red[3]);
3289             p0 = [1].concat(blue[0]);
3290             p1 = [1].concat(blue[1]);
3291 
3292             m = this.meetSegmentSegment(q0, q1, p0, p1);
3293 
3294             if (m[1] >= 0.0 && m[1] <= 1.0) {
3295                 if (!testSegment || (m[2] >= 0.0 && m[2] <= 1.0)) {
3296                     return [m];
3297                 }
3298             }
3299 
3300             return [];
3301         },
3302 
3303         /**
3304          * Find the nr-th intersection point of two Bezier curve segments.
3305          * @param {Array} red Array of four coordinate arrays of length 2 defining the first
3306          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3307          * @param {Array} blue Array of four coordinate arrays of length 2 defining the second
3308          * Bezier curve segment, i.e. [[x0,y0], [x1,y1], [x2,y2], [x3,y3]].
3309          * @param {Boolean} testSegment Test if intersection has to be inside of the segment or somewhere on the line defined by the segment
3310          * @returns {Array} Array containing the list of all intersection points as homogeneous coordinate arrays plus
3311          * preimages [x,y], t_1, t_2] of the two Bezier curve segments.
3312          *
3313          */
3314         meetBeziersegmentBeziersegment: function (red, blue, testSegment) {
3315             var L, L2, i;
3316 
3317             if (red.length === 4 && blue.length === 4) {
3318                 L = this._bezierMeetSubdivision(red, blue, 0);
3319             } else {
3320                 L = this._bezierLineMeetSubdivision(red, blue, 0, testSegment);
3321             }
3322 
3323             L.sort(function (a, b) {
3324                 return (a[1] - b[1]) * 10000000.0 + (a[2] - b[2]);
3325             });
3326 
3327             L2 = [];
3328             for (i = 0; i < L.length; i++) {
3329                 // Only push entries different from their predecessor
3330                 if (i === 0 || L[i][1] !== L[i - 1][1] || L[i][2] !== L[i - 1][2]) {
3331                     L2.push(L[i]);
3332                 }
3333             }
3334             return L2;
3335         },
3336 
3337         /**
3338          * Find the nr-th intersection point of two Bezier curves, i.e. curves with bezierDegree == 3.
3339          * @param {JXG.Curve} red Curve with bezierDegree == 3
3340          * @param {JXG.Curve} blue Curve with bezierDegree == 3
3341          * @param {Number|Function} nr The number of the intersection point which should be returned.
3342          * @returns {Array} The homogeneous coordinates of the nr-th intersection point.
3343          */
3344         meetBezierCurveRedBlueSegments: function (red, blue, nr) {
3345             var p, i, j, k,
3346                 n = Type.evaluate(nr),
3347                 po, tmp,
3348                 redArr,
3349                 blueArr,
3350                 bbr,
3351                 bbb,
3352                 intersections,
3353                 startRed = 0,
3354                 startBlue = 0,
3355                 lenBlue, lenRed,
3356                 L = [];
3357 
3358             if (blue.numberPoints < blue.bezierDegree + 1 || red.numberPoints < red.bezierDegree + 1) {
3359                 return [0, NaN, NaN];
3360             }
3361             if (red.bezierDegree === 1 && blue.bezierDegree === 3) {
3362                 tmp = red;
3363                 red = blue;
3364                 blue = tmp;
3365             }
3366 
3367             lenBlue = blue.numberPoints - blue.bezierDegree;
3368             lenRed = red.numberPoints - red.bezierDegree;
3369 
3370             // For sectors, we ignore the "legs"
3371             if (red.type === Const.OBJECT_TYPE_SECTOR) {
3372                 startRed = 3;
3373                 lenRed -= 3;
3374             }
3375             if (blue.type === Const.OBJECT_TYPE_SECTOR) {
3376                 startBlue = 3;
3377                 lenBlue -= 3;
3378             }
3379 
3380             for (i = startRed; i < lenRed; i += red.bezierDegree) {
3381                 p = red.points;
3382                 redArr = [p[i].usrCoords.slice(1), p[i + 1].usrCoords.slice(1)];
3383                 if (red.bezierDegree === 3) {
3384                     redArr[2] = p[i + 2].usrCoords.slice(1);
3385                     redArr[3] = p[i + 3].usrCoords.slice(1);
3386                 }
3387 
3388                 bbr = this._bezierBbox(redArr);
3389 
3390                 for (j = startBlue; j < lenBlue; j += blue.bezierDegree) {
3391                     p = blue.points;
3392                     blueArr = [p[j].usrCoords.slice(1), p[j + 1].usrCoords.slice(1)];
3393                     if (blue.bezierDegree === 3) {
3394                         blueArr[2] = p[j + 2].usrCoords.slice(1);
3395                         blueArr[3] = p[j + 3].usrCoords.slice(1);
3396                     }
3397 
3398                     bbb = this._bezierBbox(blueArr);
3399                     if (this._bezierOverlap(bbr, bbb)) {
3400                         intersections = this.meetBeziersegmentBeziersegment(redArr, blueArr);
3401                         if (intersections.length === 0) {
3402                             continue;
3403                         }
3404                         for (k = 0; k < intersections.length; k++) {
3405                             po = intersections[k];
3406                             if (
3407                                 po[1] < -Mat.eps ||
3408                                 po[1] > 1 + Mat.eps ||
3409                                 po[2] < -Mat.eps ||
3410                                 po[2] > 1 + Mat.eps
3411                             ) {
3412                                 continue;
3413                             }
3414                             L.push(po);
3415                         }
3416                         if (L.length > n) {
3417                             return L[n][0];
3418                         }
3419                     }
3420                 }
3421             }
3422             if (L.length > n) {
3423                 return L[n][0];
3424             }
3425 
3426             return [0, NaN, NaN];
3427         },
3428 
3429         bezierSegmentEval: function (t, curve) {
3430             var f,
3431                 x,
3432                 y,
3433                 t1 = 1.0 - t;
3434 
3435             x = 0;
3436             y = 0;
3437 
3438             f = t1 * t1 * t1;
3439             x += f * curve[0][0];
3440             y += f * curve[0][1];
3441 
3442             f = 3.0 * t * t1 * t1;
3443             x += f * curve[1][0];
3444             y += f * curve[1][1];
3445 
3446             f = 3.0 * t * t * t1;
3447             x += f * curve[2][0];
3448             y += f * curve[2][1];
3449 
3450             f = t * t * t;
3451             x += f * curve[3][0];
3452             y += f * curve[3][1];
3453 
3454             return [1.0, x, y];
3455         },
3456 
3457         /**
3458          * Generate the defining points of a 3rd degree bezier curve that approximates
3459          * a circle sector defined by three coordinate points A, B, C, each defined by an array of length three.
3460          * The coordinate arrays are given in homogeneous coordinates.
3461          * @param {Array} A First point
3462          * @param {Array} B Second point (intersection point)
3463          * @param {Array} C Third point
3464          * @param {Boolean} withLegs Flag. If true the legs to the intersection point are part of the curve.
3465          * @param {Number} sgn Wither 1 or -1. Needed for minor and major arcs. In case of doubt, use 1.
3466          */
3467         bezierArc: function (A, B, C, withLegs, sgn) {
3468             var p1, p2, p3, p4,
3469                 r,
3470                 phi, beta, delta,
3471                 // PI2 = Math.PI * 0.5,
3472                 x = B[1],
3473                 y = B[2],
3474                 z = B[0],
3475                 dataX = [],
3476                 dataY = [],
3477                 co, si,
3478                 ax, ay,
3479                 bx, by,
3480                 k, v, d,
3481                 matrix;
3482 
3483             r = this.distance(B, A);
3484 
3485             // x,y, z is intersection point. Normalize it.
3486             x /= z;
3487             y /= z;
3488 
3489             phi = this.rad(A.slice(1), B.slice(1), C.slice(1));
3490             if (sgn === -1) {
3491                 phi = 2 * Math.PI - phi;
3492             }
3493 
3494             // Always divide the arc into four Bezier arcs.
3495             // Otherwise, the position of gliders on this arc
3496             // will be wrong.
3497             delta = phi / 4;
3498 
3499 
3500             p1 = A;
3501             p1[1] /= p1[0];
3502             p1[2] /= p1[0];
3503             p1[0] /= p1[0];
3504 
3505             p4 = p1.slice(0);
3506 
3507             if (withLegs) {
3508                 dataX = [x, x + 0.333 * (p1[1] - x), x + 0.666 * (p1[1] - x), p1[1]];
3509                 dataY = [y, y + 0.333 * (p1[2] - y), y + 0.666 * (p1[2] - y), p1[2]];
3510             } else {
3511                 dataX = [p1[1]];
3512                 dataY = [p1[2]];
3513             }
3514 
3515             while (phi > Mat.eps) {
3516                 // if (phi > PI2) {
3517                 //     beta = PI2;
3518                 //     phi -= PI2;
3519                 // } else {
3520                 //     beta = phi;
3521                 //     phi = 0;
3522                 // }
3523                 if (phi > delta) {
3524                     beta = delta;
3525                     phi -= delta;
3526                 } else {
3527                     beta = phi;
3528                     phi = 0;
3529                 }
3530 
3531                 co = Math.cos(sgn * beta);
3532                 si = Math.sin(sgn * beta);
3533 
3534                 matrix = [
3535                     [1, 0, 0],
3536                     [x * (1 - co) + y * si, co, -si],
3537                     [y * (1 - co) - x * si, si, co]
3538                 ];
3539                 v = Mat.matVecMult(matrix, p1);
3540                 p4 = [v[0] / v[0], v[1] / v[0], v[2] / v[0]];
3541 
3542                 ax = p1[1] - x;
3543                 ay = p1[2] - y;
3544                 bx = p4[1] - x;
3545                 by = p4[2] - y;
3546                 d = Mat.hypot(ax + bx, ay + by);
3547 
3548                 if (Math.abs(by - ay) > Mat.eps) {
3549                     k = ((((ax + bx) * (r / d - 0.5)) / (by - ay)) * 8) / 3;
3550                 } else {
3551                     k = ((((ay + by) * (r / d - 0.5)) / (ax - bx)) * 8) / 3;
3552                 }
3553 
3554                 p2 = [1, p1[1] - k * ay, p1[2] + k * ax];
3555                 p3 = [1, p4[1] + k * by, p4[2] - k * bx];
3556 
3557                 Type.concat(dataX, [p2[1], p3[1], p4[1]]);
3558                 Type.concat(dataY, [p2[2], p3[2], p4[2]]);
3559                 p1 = p4.slice(0);
3560             }
3561 
3562             if (withLegs) {
3563                 Type.concat(dataX, [
3564                     p4[1] + 0.333 * (x - p4[1]),
3565                     p4[1] + 0.666 * (x - p4[1]),
3566                     x
3567                 ]);
3568                 Type.concat(dataY, [
3569                     p4[2] + 0.333 * (y - p4[2]),
3570                     p4[2] + 0.666 * (y - p4[2]),
3571                     y
3572                 ]);
3573             }
3574 
3575             return [dataX, dataY];
3576         },
3577 
3578         /****************************************/
3579         /****           PROJECTIONS          ****/
3580         /****************************************/
3581 
3582         /**
3583          * Calculates the coordinates of the projection of a given point on a given circle. I.o.w. the
3584          * nearest one of the two intersection points of the line through the given point and the circles
3585          * center.
3586          * @param {JXG.Point|JXG.Coords} point Point to project or coords object to project.
3587          * @param {JXG.Circle} circle Circle on that the point is projected.
3588          * @param {JXG.Board} [board=point.board] Reference to the board
3589          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle.
3590          */
3591         projectPointToCircle: function (point, circle, board) {
3592             var dist,
3593                 P,
3594                 x,
3595                 y,
3596                 factor,
3597                 M = circle.center.coords.usrCoords;
3598 
3599             if (!Type.exists(board)) {
3600                 board = point.board;
3601             }
3602 
3603             // gave us a point
3604             if (Type.isPoint(point)) {
3605                 dist = point.coords.distance(Const.COORDS_BY_USER, circle.center.coords);
3606                 P = point.coords.usrCoords;
3607                 // gave us coords
3608             } else {
3609                 dist = point.distance(Const.COORDS_BY_USER, circle.center.coords);
3610                 P = point.usrCoords;
3611             }
3612 
3613             if (Math.abs(dist) < Mat.eps) {
3614                 dist = Mat.eps;
3615             }
3616 
3617             factor = circle.Radius() / dist;
3618             x = M[1] + factor * (P[1] - M[1]);
3619             y = M[2] + factor * (P[2] - M[2]);
3620 
3621             return new Coords(Const.COORDS_BY_USER, [x, y], board);
3622         },
3623 
3624         /**
3625          * Calculates the coordinates of the orthogonal projection of a given point on a given line. I.o.w. the
3626          * intersection point of the given line and its perpendicular through the given point.
3627          * @param {JXG.Point|JXG.Coords} point Point to project.
3628          * @param {JXG.Line} line Line on that the point is projected.
3629          * @param {JXG.Board} [board=point.board|board=line.board] Reference to a board.
3630          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given line.
3631          */
3632         projectPointToLine: function (point, line, board) {
3633             var v = [0, line.stdform[1], line.stdform[2]],
3634                 coords;
3635 
3636             if (!Type.exists(board)) {
3637                 if (Type.exists(point.coords)) {
3638                     board = point.board;
3639                 } else {
3640                     board = line.board;
3641                 }
3642             }
3643 
3644             if (Type.exists(point.coords)) {
3645                 coords = point.coords.usrCoords;
3646             } else {
3647                 coords = point.usrCoords;
3648             }
3649 
3650             v = Mat.crossProduct(v, coords);
3651             return new Coords(Const.COORDS_BY_USER, Mat.crossProduct(v, line.stdform), board);
3652         },
3653 
3654         /**
3655          * Calculates the coordinates of the orthogonal projection of a given coordinate array on a given line
3656          * segment defined by two coordinate arrays.
3657          * @param {Array} p Point to project.
3658          * @param {Array} q1 Start point of the line segment on that the point is projected.
3659          * @param {Array} q2 End point of the line segment on that the point is projected.
3660          * @returns {Array} The coordinates of the projection of the given point on the given segment
3661          * and the factor that determines the projected point as a convex combination of the
3662          * two endpoints q1 and q2 of the segment.
3663          */
3664         projectCoordsToSegment: function (p, q1, q2) {
3665             var t,
3666                 denom,
3667                 s = [q2[1] - q1[1], q2[2] - q1[2]],
3668                 v = [p[1] - q1[1], p[2] - q1[2]];
3669 
3670             /**
3671              * If the segment has length 0, i.e. is a point,
3672              * the projection is equal to that point.
3673              */
3674             if (Math.abs(s[0]) < Mat.eps && Math.abs(s[1]) < Mat.eps) {
3675                 return [q1, 0];
3676             }
3677 
3678             t = Mat.innerProduct(v, s);
3679             denom = Mat.innerProduct(s, s);
3680             t /= denom;
3681 
3682             return [[1, t * s[0] + q1[1], t * s[1] + q1[2]], t];
3683         },
3684 
3685         /**
3686          * Finds the coordinates of the closest point on a Bezier segment of a
3687          * {@link JXG.Curve} to a given coordinate array.
3688          * @param {Array} pos Point to project in homogeneous coordinates.
3689          * @param {JXG.Curve} curve Curve of type "plot" having Bezier degree 3.
3690          * @param {Number} start Number of the Bezier segment of the curve.
3691          * @returns {Array} The coordinates of the projection of the given point
3692          * on the given Bezier segment and the preimage of the curve which
3693          * determines the closest point.
3694          */
3695         projectCoordsToBeziersegment: function (pos, curve, start) {
3696             var t0,
3697                 /** @ignore */
3698                 minfunc = function (t) {
3699                     var z = [1, curve.X(start + t), curve.Y(start + t)];
3700 
3701                     z[1] -= pos[1];
3702                     z[2] -= pos[2];
3703 
3704                     return z[1] * z[1] + z[2] * z[2];
3705                 };
3706 
3707             t0 = JXG.Math.Numerics.fminbr(minfunc, [0.0, 1.0]);
3708 
3709             return [[1, curve.X(t0 + start), curve.Y(t0 + start)], t0];
3710         },
3711 
3712         /**
3713          * Calculates the coordinates of the projection of a given point on a given curve.
3714          * Uses {@link JXG.Math.Geometry.projectCoordsToCurve}.
3715          *
3716          * @param {JXG.Point} point Point to project.
3717          * @param {JXG.Curve} curve Curve on that the point is projected.
3718          * @param {JXG.Board} [board=point.board] Reference to a board.
3719          * @see JXG.Math.Geometry.projectCoordsToCurve
3720          * @returns {Array} [JXG.Coords, position] The coordinates of the projection of the given
3721          * point on the given graph and the relative position on the curve (real number).
3722          */
3723         projectPointToCurve: function (point, curve, board) {
3724             if (!Type.exists(board)) {
3725                 board = point.board;
3726             }
3727 
3728             var x = point.X(),
3729                 y = point.Y(),
3730                 t = point.position,
3731                 result;
3732 
3733             if (!Type.exists(t)) {
3734                 t = curve.evalVisProp('curvetype') === 'functiongraph' ? x : 0.0;
3735             }
3736             result = this.projectCoordsToCurve(x, y, t, curve, board);
3737             // point.position = result[1];
3738 
3739             return result;
3740         },
3741 
3742         /**
3743          * Calculates the coordinates of the projection of a coordinates pair on a given curve. In case of
3744          * function graphs this is the
3745          * intersection point of the curve and the parallel to y-axis through the given point.
3746          * @param {Number} x coordinate to project.
3747          * @param {Number} y coordinate to project.
3748          * @param {Number} t start value for newtons method
3749          * @param {JXG.Curve} curve Curve on that the point is projected.
3750          * @param {JXG.Board} [board=curve.board] Reference to a board.
3751          * @see JXG.Math.Geometry.projectPointToCurve
3752          * @returns {JXG.Coords} Array containing the coordinates of the projection of the given point on the given curve and
3753          * the position on the curve.
3754          */
3755         projectCoordsToCurve: function (x, y, t, curve, board) {
3756             var newCoords, newCoordsObj,
3757                 i, j, mindist, dist, lbda,
3758                 v, coords, d, p1, p2, res, minfunc,
3759                 t_new, f_new, f_old, dy,
3760                 delta, delta1, delta2, steps,
3761                 minX, maxX, minX_glob, maxX_glob,
3762                 infty = Number.POSITIVE_INFINITY;
3763 
3764             if (!Type.exists(board)) {
3765                 board = curve.board;
3766             }
3767 
3768             if (curve.evalVisProp('curvetype') === 'plot') {
3769                 t = 0;
3770                 mindist = infty;
3771                 if (curve.numberPoints === 0) {
3772                     newCoords = [0, 1, 1];
3773                 } else {
3774                     newCoords = [curve.Z(0), curve.X(0), curve.Y(0)];
3775                 }
3776 
3777                 if (curve.numberPoints > 1) {
3778                     v = [1, x, y];
3779                     if (curve.bezierDegree === 3) {
3780                         j = 0;
3781                     } else {
3782                         p1 = [curve.Z(0), curve.X(0), curve.Y(0)];
3783                     }
3784                     for (i = 0; i < curve.numberPoints - 1; i++) {
3785                         if (curve.bezierDegree === 3) {
3786                             res = this.projectCoordsToBeziersegment(v, curve, j);
3787                         } else {
3788                             p2 = [curve.Z(i + 1), curve.X(i + 1), curve.Y(i + 1)];
3789                             res = this.projectCoordsToSegment(v, p1, p2);
3790                         }
3791                         lbda = res[1];
3792                         coords = res[0];
3793 
3794                         if (0.0 <= lbda && lbda <= 1.0) {
3795                             dist = this.distance(coords, v);
3796                             d = i + lbda;
3797                         } else if (lbda < 0.0) {
3798                             coords = p1;
3799                             dist = this.distance(p1, v);
3800                             d = i;
3801                         } else if (lbda > 1.0 && i === curve.numberPoints - 2) {
3802                             coords = p2;
3803                             dist = this.distance(coords, v);
3804                             d = curve.numberPoints - 1;
3805                         }
3806 
3807                         if (dist < mindist) {
3808                             mindist = dist;
3809                             t = d;
3810                             newCoords = coords;
3811                         }
3812 
3813                         if (curve.bezierDegree === 3) {
3814                             j++;
3815                             i += 2;
3816                         } else {
3817                             p1 = p2;
3818                         }
3819                     }
3820                 }
3821 
3822                 newCoordsObj = new Coords(Const.COORDS_BY_USER, newCoords, board);
3823             } else {
3824                 // 'parameter', 'polar', 'functiongraph'
3825 
3826                 minX_glob = curve.minX();
3827                 maxX_glob = curve.maxX();
3828                 minX = minX_glob;
3829                 maxX = maxX_glob;
3830 
3831                 if (curve.evalVisProp('curvetype') === 'functiongraph') {
3832                     // Restrict the possible position of t
3833                     // to the projection of a circle to the x-axis (= t-axis)
3834                     dy = Math.abs(y - curve.Y(x));
3835                     if (!isNaN(dy)) {
3836                         minX = x - dy;
3837                         maxX = x + dy;
3838                     }
3839                 }
3840 
3841                 /**
3842                  * @ignore
3843                  * Find t such that the Euclidean distance between
3844                  * [x, y] and [curve.X(t), curve.Y(t)]
3845                  * is minimized.
3846                  */
3847                 minfunc = function (t) {
3848                     var dx, dy;
3849 
3850                     if (t < minX_glob || t > curve.maxX_glob) {
3851                         return Infinity;
3852                     }
3853                     dx = x - curve.X(t);
3854                     dy = y - curve.Y(t);
3855                     return dx * dx + dy * dy;
3856                 };
3857 
3858                 // Search t which minimizes minfunc(t)
3859                 // in discrete steps
3860                 f_old = minfunc(t);
3861                 steps = 50;
3862                 delta = (maxX - minX) / steps;
3863                 t_new = minX;
3864                 for (i = 0; i < steps; i++) {
3865                     f_new = minfunc(t_new);
3866 
3867                     if (f_new < f_old || f_old === Infinity || isNaN(f_old)) {
3868                         t = t_new;
3869                         f_old = f_new;
3870                     }
3871 
3872                     t_new += delta;
3873                 }
3874 
3875                 // t = Numerics.root(Numerics.D(minfunc), t);
3876 
3877                 // Ensure that minfunc is defined on the
3878                 // enclosing interval [t-delta1, t+delta2]
3879                 delta1 = delta;
3880                 for (i = 0; i < 20 && isNaN(minfunc(t - delta1)); i++, delta1 *= 0.5);
3881                 if (isNaN(minfunc(t - delta1))) {
3882                     delta1 = 0.0;
3883                 }
3884                 delta2 = delta;
3885                 for (i = 0; i < 20 && isNaN(minfunc(t + delta2)); i++, delta2 *= 0.5);
3886                 if (isNaN(minfunc(t + delta2))) {
3887                     delta2 = 0.0;
3888                 }
3889 
3890                 // Finally, apply mathemetical optimization in the determined interval
3891                 t = Numerics.fminbr(minfunc, [
3892                     Math.max(t - delta1, minX),
3893                     Math.min(t + delta2, maxX)
3894                 ]);
3895 
3896                 // Distinction between closed and open curves is not necessary.
3897                 // If closed, the cyclic projection shift will work anyhow
3898                 // if (Math.abs(curve.X(minX) - curve.X(maxX)) < Mat.eps &&
3899                 //     Math.abs(curve.Y(minX) - curve.Y(maxX)) < Mat.eps) {
3900                 //     // Cyclically
3901                 //     if (t < minX) {console.log(t)
3902                 //         t = maxX + t - minX;
3903                 //     }
3904                 //     if (t > maxX) {
3905                 //         t = minX + t - maxX;
3906                 //     }
3907                 // } else {
3908 
3909                 t = t < minX_glob ? minX_glob : t;
3910                 t = t > maxX_glob ? maxX_glob : t;
3911                 // }
3912 
3913                 newCoordsObj = new Coords(
3914                     Const.COORDS_BY_USER,
3915                     [curve.X(t), curve.Y(t)],
3916                     board
3917                 );
3918             }
3919 
3920             return [curve.updateTransform(newCoordsObj), t];
3921         },
3922 
3923         /**
3924          * Calculates the coordinates of the closest orthogonal projection of a given coordinate array onto the
3925          * border of a polygon.
3926          * @param {Array} p Point to project.
3927          * @param {JXG.Polygon} pol Polygon element
3928          * @returns {Array} The coordinates of the closest projection of the given point to the border of the polygon.
3929          */
3930         projectCoordsToPolygon: function (p, pol) {
3931             var i,
3932                 len = pol.vertices.length,
3933                 d_best = Infinity,
3934                 d,
3935                 projection,
3936                 proj,
3937                 bestprojection;
3938 
3939             for (i = 0; i < len - 1; i++) {
3940                 projection = JXG.Math.Geometry.projectCoordsToSegment(
3941                     p,
3942                     pol.vertices[i].coords.usrCoords,
3943                     pol.vertices[i + 1].coords.usrCoords
3944                 );
3945 
3946                 if (0 <= projection[1] && projection[1] <= 1) {
3947                     d = JXG.Math.Geometry.distance(projection[0], p, 3);
3948                     proj = projection[0];
3949                 } else if (projection[1] < 0) {
3950                     d = JXG.Math.Geometry.distance(pol.vertices[i].coords.usrCoords, p, 3);
3951                     proj = pol.vertices[i].coords.usrCoords;
3952                 } else {
3953                     d = JXG.Math.Geometry.distance(pol.vertices[i + 1].coords.usrCoords, p, 3);
3954                     proj = pol.vertices[i + 1].coords.usrCoords;
3955                 }
3956                 if (d < d_best) {
3957                     bestprojection = proj.slice(0);
3958                     d_best = d;
3959                 }
3960             }
3961             return bestprojection;
3962         },
3963 
3964         /**
3965          * Calculates the coordinates of the projection of a given point on a given turtle. A turtle consists of
3966          * one or more curves of curveType 'plot'. Uses {@link JXG.Math.Geometry.projectPointToCurve}.
3967          * @param {JXG.Point} point Point to project.
3968          * @param {JXG.Turtle} turtle on that the point is projected.
3969          * @param {JXG.Board} [board=point.board] Reference to a board.
3970          * @returns {Array} [JXG.Coords, position] Array containing the coordinates of the projection of the given point on the turtle and
3971          * the position on the turtle.
3972          */
3973         projectPointToTurtle: function (point, turtle, board) {
3974             var newCoords,
3975                 t,
3976                 x,
3977                 y,
3978                 i,
3979                 dist,
3980                 el,
3981                 minEl,
3982                 res,
3983                 newPos,
3984                 np = 0,
3985                 npmin = 0,
3986                 mindist = Number.POSITIVE_INFINITY,
3987                 len = turtle.objects.length;
3988 
3989             if (!Type.exists(board)) {
3990                 board = point.board;
3991             }
3992 
3993             // run through all curves of this turtle
3994             for (i = 0; i < len; i++) {
3995                 el = turtle.objects[i];
3996 
3997                 if (el.elementClass === Const.OBJECT_CLASS_CURVE) {
3998                     res = this.projectPointToCurve(point, el);
3999                     newCoords = res[0];
4000                     newPos = res[1];
4001                     dist = this.distance(newCoords.usrCoords, point.coords.usrCoords);
4002 
4003                     if (dist < mindist) {
4004                         x = newCoords.usrCoords[1];
4005                         y = newCoords.usrCoords[2];
4006                         t = newPos;
4007                         mindist = dist;
4008                         minEl = el;
4009                         npmin = np;
4010                     }
4011                     np += el.numberPoints;
4012                 }
4013             }
4014 
4015             newCoords = new Coords(Const.COORDS_BY_USER, [x, y], board);
4016             // point.position = t + npmin;
4017             // return minEl.updateTransform(newCoords);
4018             return [minEl.updateTransform(newCoords), t + npmin];
4019         },
4020 
4021         /**
4022          * Trivial projection of a point to another point.
4023          * @param {JXG.Point} point Point to project (not used).
4024          * @param {JXG.Point} dest Point on that the point is projected.
4025          * @returns {JXG.Coords} The coordinates of the projection of the given point on the given circle.
4026          */
4027         projectPointToPoint: function (point, dest) {
4028             return dest.coords;
4029         },
4030 
4031         /**
4032          *
4033          * @param {JXG.Point|JXG.Coords} point
4034          * @param {JXG.Board} [board]
4035          */
4036         projectPointToBoard: function (point, board) {
4037             var i,
4038                 l,
4039                 c,
4040                 brd = board || point.board,
4041                 // comparison factor, point coord idx, bbox idx, 1st bbox corner x & y idx, 2nd bbox corner x & y idx
4042                 config = [
4043                     // left
4044                     [1, 1, 0, 0, 3, 0, 1],
4045                     // top
4046                     [-1, 2, 1, 0, 1, 2, 1],
4047                     // right
4048                     [-1, 1, 2, 2, 1, 2, 3],
4049                     // bottom
4050                     [1, 2, 3, 0, 3, 2, 3]
4051                 ],
4052                 coords = point.coords || point,
4053                 bbox = brd.getBoundingBox();
4054 
4055             for (i = 0; i < 4; i++) {
4056                 c = config[i];
4057                 if (c[0] * coords.usrCoords[c[1]] < c[0] * bbox[c[2]]) {
4058                     // define border
4059                     l = Mat.crossProduct(
4060                         [1, bbox[c[3]], bbox[c[4]]],
4061                         [1, bbox[c[5]], bbox[c[6]]]
4062                     );
4063                     l[3] = 0;
4064                     l = Mat.normalize(l);
4065 
4066                     // project point
4067                     coords = this.projectPointToLine({ coords: coords }, { stdform: l }, brd);
4068                 }
4069             }
4070 
4071             return coords;
4072         },
4073 
4074         /**
4075          * Calculates the distance of a point to a line. The point and the line are given by homogeneous
4076          * coordinates. For lines this can be line.stdform.
4077          * @param {Array} point Homogeneous coordinates of a point.
4078          * @param {Array} line Homogeneous coordinates of a line ([C,A,B] where A*x+B*y+C*z=0).
4079          * @returns {Number} Distance of the point to the line.
4080          */
4081         distPointLine: function (point, line) {
4082             var a = line[1],
4083                 b = line[2],
4084                 c = line[0],
4085                 nom;
4086 
4087             if (Math.abs(a) + Math.abs(b) < Mat.eps) {
4088                 return Number.POSITIVE_INFINITY;
4089             }
4090 
4091             nom = a * point[1] + b * point[2] + c;
4092             a *= a;
4093             b *= b;
4094 
4095             return Math.abs(nom) / Math.sqrt(a + b);
4096         },
4097 
4098         /**
4099          * Determine the (Euclidean) distance between a point q and a line segment
4100          * defined by two points p1 and p2. In case p1 equals p2, the distance to this
4101          * point is returned.
4102          *
4103          * @param {Array} q Homogeneous coordinates of q
4104          * @param {Array} p1 Homogeneous coordinates of p1
4105          * @param {Array} p2 Homogeneous coordinates of p2
4106          * @returns {Number} Distance of q to line segment [p1, p2]
4107          */
4108         distPointSegment: function (q, p1, p2) {
4109             var x, y, dx, dy,
4110                 den, lbda,
4111                 eps = Mat.eps * Mat.eps,
4112                 huge = 1000000;
4113 
4114             // Difference q - p1
4115             x = q[1] - p1[1];
4116             y = q[2] - p1[2];
4117             x = (x === Infinity) ? huge : (x === -Infinity) ? -huge : x;
4118             y = (y === Infinity) ? huge : (y === -Infinity) ? -huge : y;
4119 
4120             // Difference p2 - p1
4121             dx = p2[1] - p1[1];
4122             dy = p2[2] - p1[2];
4123             dx = (dx === Infinity) ? huge : (dx === -Infinity) ? -huge : dx;
4124             dy = (dy === Infinity) ? huge : (dy === -Infinity) ? -huge : dy;
4125 
4126             // If den==0 then p1 and p2 are identical
4127             // In this case the distance to p1 is returned
4128             den = dx * dx + dy * dy;
4129             if (den > eps) {
4130                 lbda = (x * dx + y * dy) / den;
4131                 if (lbda < 0.0) {
4132                     lbda = 0.0;
4133                 } else if (lbda > 1.0) {
4134                     lbda = 1.0;
4135                 }
4136                 x -= lbda * dx;
4137                 y -= lbda * dy;
4138             }
4139 
4140             return Mat.hypot(x, y);
4141         },
4142 
4143         /* ***************************************/
4144         /* *** 3D CALCULATIONS ****/
4145         /* ***************************************/
4146 
4147         /**
4148          * Generate the function which computes the data of the intersection between
4149          * <ul>
4150          * <li> plane3d, plane3d,
4151          * <li> plane3d, sphere3d,
4152          * <li> sphere3d, plane3d,
4153          * <li> sphere3d, sphere3d
4154          * </ul>
4155          *
4156          * @param {JXG.GeometryElement3D} el1 Plane or sphere element
4157          * @param {JXG.GeometryElement3D} el2 Plane or sphere element
4158          * @returns {Array} of functions needed as input to create the intersecting line or circle.
4159          *
4160          */
4161         intersectionFunction3D: function (view, el1, el2) {
4162             var func,
4163                 that = this;
4164 
4165             if (el1.type === Const.OBJECT_TYPE_PLANE3D) {
4166                 if (el2.type === Const.OBJECT_TYPE_PLANE3D) {
4167                     // func = () => view.intersectionPlanePlane(el1, el2)[i];
4168                     func = view.intersectionPlanePlane(el1, el2);
4169                 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) {
4170                     func = that.meetPlaneSphere(el1, el2);
4171                 }
4172             } else if (el1.type === Const.OBJECT_TYPE_SPHERE3D) {
4173                 if (el2.type === Const.OBJECT_TYPE_PLANE3D) {
4174                     func = that.meetPlaneSphere(el2, el1);
4175                 } else if (el2.type === Const.OBJECT_TYPE_SPHERE3D) {
4176                     func = that.meetSphereSphere(el1, el2);
4177                 }
4178             }
4179 
4180             return func;
4181         },
4182 
4183         /**
4184          * Intersecting point of three planes in 3D. The planes
4185          * are given in Hesse normal form.
4186          *
4187          * @param {Array} n1 Hesse normal form vector of plane 1
4188          * @param {Number} d1 Hesse normal form right hand side of plane 1
4189          * @param {Array} n2 Hesse normal form vector of plane 2
4190          * @param {Number} d2 Hesse normal form right hand side of plane 2
4191          * @param {Array} n3 Hesse normal form vector of plane 1
4192          * @param {Number} d3 Hesse normal form right hand side of plane 3
4193          * @returns {Array} Coordinates array of length 4 of the intersecting point
4194          */
4195         meet3Planes: function (n1, d1, n2, d2, n3, d3) {
4196             var p = [1, 0, 0, 0],
4197                 n31, n12, n23,
4198                 denom,
4199                 i;
4200 
4201             n31 = Mat.crossProduct(n3.slice(1), n1.slice(1));
4202             n12 = Mat.crossProduct(n1.slice(1), n2.slice(1));
4203             n23 = Mat.crossProduct(n2.slice(1), n3.slice(1));
4204 
4205             denom = Mat.innerProduct(n1.slice(1), n23, 3);
4206             for (i = 0; i < 3; i++) {
4207                 p[i + 1] = (d1 * n23[i] + d2 * n31[i] + d3 * n12[i]) / denom;
4208             }
4209 
4210             return p;
4211         },
4212 
4213         /**
4214          * Direction of intersecting line of two planes in 3D.
4215          *
4216          * @param {Array} v11 First vector spanning plane 1 (homogeneous coordinates)
4217          * @param {Array} v12 Second vector spanning plane 1 (homogeneous coordinates)
4218          * @param {Array} v21 First vector spanning plane 2 (homogeneous coordinates)
4219          * @param {Array} v22 Second vector spanning plane 2 (homogeneous coordinates)
4220          * @returns {Array} Coordinates array of length 4 of the direction  (homogeneous coordinates)
4221          */
4222         meetPlanePlane: function (v11, v12, v21, v22) {
4223             var no1,
4224                 no2,
4225                 v, w;
4226 
4227             v = v11.slice(1);
4228             w = v12.slice(1);
4229             no1 = Mat.crossProduct(v, w);
4230 
4231             v = v21.slice(1);
4232             w = v22.slice(1);
4233             no2 = Mat.crossProduct(v, w);
4234 
4235             w = Mat.crossProduct(no1, no2);
4236             w.unshift(0);
4237             return w;
4238         },
4239 
4240         meetPlaneSphere: function (el1, el2) {
4241             var dis = function () {
4242                     return Mat.innerProduct(el1.normal, el2.center.coords, 4) - el1.d;
4243                 };
4244 
4245             return [
4246                 // Center
4247                 function() {
4248                     return Mat.axpy(-dis(), el1.normal, el2.center.coords);
4249                 },
4250                 // Normal
4251                 el1.normal,
4252                 // Radius
4253                 function () {
4254                     // Radius (returns NaN if spheres don't touch)
4255                     var r = el2.Radius(),
4256                         s = dis();
4257                     return Math.sqrt(r * r - s * s);
4258                 }
4259             ];
4260         },
4261 
4262         meetSphereSphere: function (el1, el2) {
4263             var skew = function () {
4264                     var dist = el1.center.distance(el2.center),
4265                         r1 = el1.Radius(),
4266                         r2 = el2.Radius();
4267                     return (r1 - r2) * (r1 + r2) / (dist * dist);
4268                 };
4269             return [
4270                 // Center
4271                 function () {
4272                     var s = skew();
4273                     return [
4274                         1,
4275                         0.5 * ((1 - s) * el1.center.coords[1] + (1 + s) * el2.center.coords[1]),
4276                         0.5 * ((1 - s) * el1.center.coords[2] + (1 + s) * el2.center.coords[2]),
4277                         0.5 * ((1 - s) * el1.center.coords[3] + (1 + s) * el2.center.coords[3])
4278                     ];
4279                 },
4280                 // Normal
4281                 function() {
4282                     return Stat.subtract(el2.center.coords, el1.center.coords);
4283                 },
4284                 // Radius
4285                 function () {
4286                     // Radius (returns NaN if spheres don't touch)
4287                     var dist = el1.center.distance(el2.center),
4288                         r1 = el1.Radius(),
4289                         r2 = el2.Radius(),
4290                         s = skew(),
4291                         rIxnSq = 0.5 * (r1 * r1 + r2 * r2 - 0.5 * dist * dist * (1 + s * s));
4292                     return Math.sqrt(rIxnSq);
4293                 }
4294             ];
4295         },
4296 
4297         /**
4298          * Test if parameters are inside of allowed ranges
4299          *
4300          * @param {Array} params Array of length 1 or 2
4301          * @param {Array} r_u First range
4302          * @param {Array} [r_v] Second range
4303          * @returns Boolean
4304          * @private
4305          */
4306         _paramsOutOfRange: function(params, r_u, r_v) {
4307             return params[0] < r_u[0] || params[0] > r_u[1] ||
4308                 (params.length > 1 && (params[1] < r_v[0] || params[1] > r_v[1]));
4309         },
4310 
4311         /**
4312          * Given the 2D screen coordinates of a point, finds the nearest point on the given
4313          * parametric curve or surface, and returns its view-space coordinates.
4314          * @param {Array} p 3D coordinates for which the closest point on the curve point is searched.
4315          * @param {JXG.Curve3D|JXG.Surface3D} target Parametric curve or surface to project to.
4316          * @param {Array} params New position of point on the target (i.e. it is a return value),
4317          * modified in place during the search, ending up at the nearest point.
4318          * Usually, point.position is supplied for params.
4319          *
4320          * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface.
4321          */
4322         projectCoordsToParametric: function (p, target, n, params) {
4323             // The variables and parameters for the Cobyla constrained
4324             // minimization algorithm are explained in the Cobyla.js comments
4325             var rhobeg,                // initial size of simplex (Cobyla)
4326                 rhoend,                // finial size of simplex (Cobyla)
4327                 iprint = 0,            // no console output (Cobyla)
4328                 maxfun = 200,          // call objective function at most 200 times (Cobyla)
4329                 _minFunc,              // Objective function for Cobyla
4330                 f = Math.random() * 0.01 + 0.5,
4331                 r_u, r_v,
4332                 m = 2 * n;
4333 
4334             // adapt simplex size to parameter range
4335             if (n === 1) {
4336                 r_u = [Type.evaluate(target.range[0]), Type.evaluate(target.range[1])];
4337 
4338                 rhobeg = 0.1 * (r_u[1] - r_u[0]);
4339             } else if (n === 2) {
4340                 r_u = [Type.evaluate(target.range_u[0]), Type.evaluate(target.range_u[1])];
4341                 r_v = [Type.evaluate(target.range_v[0]), Type.evaluate(target.range_v[1])];
4342 
4343                 rhobeg = 0.1 * Math.min(
4344                     r_u[1] - r_u[0],
4345                     r_v[1] - r_v[0]
4346                 );
4347             }
4348             rhoend = rhobeg / 5e6;
4349 
4350             // Minimize distance of the new position to the original position
4351             _minFunc = function (n, m, w, con) {
4352                 var p_new = [
4353                         target.X.apply(target, w),
4354                         target.Y.apply(target, w),
4355                         target.Z.apply(target, w)
4356                     ],
4357                     xDiff = p[0] - p_new[0],
4358                     yDiff = p[1] - p_new[1],
4359                     zDiff = p[2] - p_new[2];
4360 
4361                 if (m >= 2) {
4362                     con[0] =  w[0] - r_u[0];
4363                     con[1] = -w[0] + r_u[1];
4364                 }
4365                 if (m >= 4) {
4366                     con[2] =  w[1] - r_v[0];
4367                     con[3] = -w[1] + r_v[1];
4368                 }
4369 
4370                 return xDiff * xDiff + yDiff * yDiff + zDiff * zDiff;
4371             };
4372 
4373             // First optimization without range constraints to give a smooth draag experience on
4374             // cyclic structures.
4375 
4376             // Set the start values
4377             if (params.length === 0) {
4378                 // If length > 0: take the previous position as start values for the optimization
4379                 params[0] = f * (r_u[0] + r_u[1]);
4380                 if (n === 2) { params[1] = f * (r_v[0] + r_v[1]); }
4381             }
4382             Mat.Nlp.FindMinimum(_minFunc, n, 0, params, rhobeg, rhoend, iprint, maxfun);
4383             // Update p which is used subsequently in _minFunc
4384             p = [target.X.apply(target, params),
4385                 target.Y.apply(target, params),
4386                 target.Z.apply(target, params)
4387             ];
4388 
4389             // If the optimal params are outside of the rang
4390             // Second optimization to obey the range constraints
4391 
4392             if (this._paramsOutOfRange(params, r_u, r_v)) {
4393                 // Set the start values again
4394                 params[0] = f * (r_u[0] + r_u[1]);
4395                 if (n === 2) { params[1] = f * (r_v[0] + r_v[1]); }
4396 
4397                 Mat.Nlp.FindMinimum(_minFunc, n, m, params, rhobeg, rhoend, iprint, maxfun);
4398             }
4399 
4400             return [1,
4401                 target.X.apply(target, params),
4402                 target.Y.apply(target, params),
4403                 target.Z.apply(target, params)
4404             ];
4405         },
4406 
4407         // /**
4408         //  * Given a the screen coordinates of a point, finds the point on the
4409         //  * given parametric curve or surface which is nearest in screen space,
4410         //  * and returns its view-space coordinates.
4411         //  * @param {Array} pScr Screen coordinates to project.
4412         //  * @param {JXG.Curve3D|JXG.Surface3D} target Parametric curve or surface to project to.
4413         //  * @param {Array} params Parameters of point on the target, initially specifying the starting point of
4414         //  * the search. The parameters are modified in place during the search, ending up at the nearest point.
4415         //  * @returns {Array} Array of length 4 containing the coordinates of the nearest point on the curve or surface.
4416         //  */
4417         // projectScreenCoordsToParametric: function (pScr, target, params) {
4418         //     // The variables and parameters for the Cobyla constrained
4419         //     // minimization algorithm are explained in the Cobyla.js comments
4420         //     var rhobeg, // initial size of simplex (Cobyla)
4421         //         rhoend, // finial size of simplex (Cobyla)
4422         //         iprint = 0, // no console output (Cobyla)
4423         //         maxfun = 200, // call objective function at most 200 times (Cobyla)
4424         //         dim = params.length,
4425         //         _minFunc; // objective function (Cobyla)
4426 
4427         //     // adapt simplex size to parameter range
4428         //     if (dim === 1) {
4429         //         rhobeg = 0.1 * (target.range[1] - target.range[0]);
4430         //     } else if (dim === 2) {
4431         //         rhobeg = 0.1 * Math.min(
4432         //             target.range_u[1] - target.range_u[0],
4433         //             target.range_v[1] - target.range_v[0]
4434         //         );
4435         //     }
4436         //     rhoend = rhobeg / 5e6;
4437 
4438         //     // minimize screen distance to cursor
4439         //     _minFunc = function (n, m, w, con) {
4440         //         var c3d = [
4441         //             1,
4442         //             target.X.apply(target, w),
4443         //             target.Y.apply(target, w),
4444         //             target.Z.apply(target, w)
4445         //         ],
4446         //         c2d = target.view.project3DTo2D(c3d),
4447         //         xDiff = pScr[0] - c2d[1],
4448         //         yDiff = pScr[1] - c2d[2];
4449 
4450         //         if (n === 1) {
4451         //             con[0] = w[0] - target.range[0];
4452         //             con[1] = -w[0] + target.range[1];
4453         //         } else if (n === 2) {
4454         //             con[0] = w[0] - target.range_u[0];
4455         //             con[1] = -w[0] + target.range_u[1];
4456         //             con[2] = w[1] - target.range_v[0];
4457         //             con[3] = -w[1] + target.range_v[1];
4458         //         }
4459 
4460         //         return xDiff * xDiff + yDiff * yDiff;
4461         //     };
4462 
4463         //     Mat.Nlp.FindMinimum(_minFunc, dim, 2 * dim, params, rhobeg, rhoend, iprint, maxfun);
4464 
4465         //     return [1, target.X.apply(target, params), target.Y.apply(target, params), target.Z.apply(target, params)];
4466         // },
4467 
4468         project3DTo3DPlane: function (point, normal, foot) {
4469             // TODO: homogeneous 3D coordinates
4470             var sol = [0, 0, 0],
4471                 le,
4472                 d1,
4473                 d2,
4474                 lbda;
4475 
4476             foot = foot || [0, 0, 0];
4477 
4478             le = Mat.norm(normal);
4479             d1 = Mat.innerProduct(point, normal, 3);
4480             d2 = Mat.innerProduct(foot, normal, 3);
4481             // (point - lbda * normal / le) * normal / le == foot * normal / le
4482             // => (point * normal - foot * normal) ==  lbda * le
4483             lbda = (d1 - d2) / le;
4484             sol = Mat.axpy(-lbda, normal, point);
4485 
4486             return sol;
4487         },
4488 
4489         getPlaneBounds: function (v1, v2, q, s, e) {
4490             var s1, s2, e1, e2, mat, rhs, sol;
4491 
4492             if (v1[2] + v2[0] !== 0) {
4493                 mat = [
4494                     [v1[0], v2[0]],
4495                     [v1[1], v2[1]]
4496                 ];
4497                 rhs = [s - q[0], s - q[1]];
4498 
4499                 sol = Numerics.Gauss(mat, rhs);
4500                 s1 = sol[0];
4501                 s2 = sol[1];
4502 
4503                 rhs = [e - q[0], e - q[1]];
4504                 sol = Numerics.Gauss(mat, rhs);
4505                 e1 = sol[0];
4506                 e2 = sol[1];
4507                 return [s1, e1, s2, e2];
4508             }
4509             return null;
4510         },
4511 
4512         /* ***************************************/
4513         /* *** Various ****/
4514         /* ***************************************/
4515 
4516         /**
4517          * Helper function to create curve which displays a Reuleaux polygons.
4518          * @param {Array} points Array of points which should be the vertices of the Reuleaux polygon. Typically,
4519          * these point list is the array vertices of a regular polygon.
4520          * @param {Number} nr Number of vertices
4521          * @returns {Array} An array containing the two functions defining the Reuleaux polygon and the two values
4522          * for the start and the end of the paramtric curve. array may be used as parent array of a
4523          * {@link JXG.Curve}.
4524          *
4525          * @example
4526          * var A = brd.create('point',[-2,-2]);
4527          * var B = brd.create('point',[0,1]);
4528          * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0});
4529          * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3),
4530          *                          {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'});
4531          *
4532          * </pre><div class="jxgbox" id="JXG2543a843-46a9-4372-abc1-94d9ad2db7ac" style="width: 300px; height: 300px;"></div>
4533          * <script type="text/javascript">
4534          * var brd = JXG.JSXGraph.initBoard('JXG2543a843-46a9-4372-abc1-94d9ad2db7ac', {boundingbox: [-5, 5, 5, -5], axis: true, showcopyright:false, shownavigation: false});
4535          * var A = brd.create('point',[-2,-2]);
4536          * var B = brd.create('point',[0,1]);
4537          * var pol = brd.create('regularpolygon',[A,B,3], {withLines:false, fillColor:'none', highlightFillColor:'none', fillOpacity:0.0});
4538          * var reuleauxTriangle = brd.create('curve', JXG.Math.Geometry.reuleauxPolygon(pol.vertices, 3),
4539          *                          {strokeWidth:6, strokeColor:'#d66d55', fillColor:'#ad5544', highlightFillColor:'#ad5544'});
4540          * </script><pre>
4541          */
4542         reuleauxPolygon: function (points, nr) {
4543             var beta,
4544                 pi2 = Math.PI * 2,
4545                 pi2_n = pi2 / nr,
4546                 diag = (nr - 1) / 2,
4547                 d = 0,
4548                 makeFct = function (which, trig) {
4549                     return function (t, suspendUpdate) {
4550                         var t1 = ((t % pi2) + pi2) % pi2,
4551                             j = Math.floor(t1 / pi2_n) % nr;
4552 
4553                         if (!suspendUpdate) {
4554                             d = points[0].Dist(points[diag]);
4555                             beta = Mat.Geometry.rad(
4556                                 [points[0].X() + 1, points[0].Y()],
4557                                 points[0],
4558                                 points[diag % nr]
4559                             );
4560                         }
4561 
4562                         if (isNaN(j)) {
4563                             return j;
4564                         }
4565 
4566                         t1 = t1 * 0.5 + j * pi2_n * 0.5 + beta;
4567 
4568                         return points[j][which]() + d * Math[trig](t1);
4569                     };
4570                 };
4571 
4572             return [makeFct("X", 'cos'), makeFct("Y", 'sin'), 0, pi2];
4573         }
4574 
4575     }
4576 );
4577 
4578 export default Mat.Geometry;
4579