1 /* 2 Copyright 2008-2026 3 Matthias Ehmann, 4 Michael Gerhaeuser, 5 Carsten Miller, 6 Alfred Wassermann 7 8 This file is part of JSXGraph. 9 10 JSXGraph is free software dual licensed under the GNU LGPL or MIT License. 11 12 You can redistribute it and/or modify it under the terms of the 13 14 * GNU Lesser General Public License as published by 15 the Free Software Foundation, either version 3 of the License, or 16 (at your option) any later version 17 OR 18 * MIT License: https://github.com/jsxgraph/jsxgraph/blob/master/LICENSE.MIT 19 20 JSXGraph is distributed in the hope that it will be useful, 21 but WITHOUT ANY WARRANTY; without even the implied warranty of 22 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 23 GNU Lesser General Public License for more details. 24 25 You should have received a copy of the GNU Lesser General Public License and 26 the MIT License along with JSXGraph. If not, see <https://www.gnu.org/licenses/> 27 and <https://opensource.org/licenses/MIT/>. 28 */ 29 30 /*global JXG: true, define: true*/ 31 /*jslint nomen: true, plusplus: true*/ 32 33 import JXG from "../jxg.js"; 34 import Const from "../base/constants.js"; 35 import Coords from "../base/coords.js"; 36 import Mat from "./math.js"; 37 import Extrapolate from "./extrapolate.js"; 38 import Numerics from "./numerics.js"; 39 import Statistics from "./statistics.js"; 40 import Geometry from "./geometry.js"; 41 import IntervalArithmetic from "./ia.js"; 42 import Type from "../utils/type.js"; 43 44 /** 45 * Functions for plotting of curves. 46 * @name JXG.Math.Plot 47 * @exports Mat.Plot as JXG.Math.Plot 48 * @namespace 49 */ 50 Mat.Plot = { 51 /** 52 * Check if at least one point on the curve is finite and real. 53 **/ 54 checkReal: function (points) { 55 var b = false, 56 i, 57 p, 58 len = points.length; 59 60 for (i = 0; i < len; i++) { 61 p = points[i].usrCoords; 62 if (!isNaN(p[1]) && !isNaN(p[2]) && Math.abs(p[0]) > Mat.eps) { 63 b = true; 64 break; 65 } 66 } 67 return b; 68 }, 69 70 //---------------------------------------------------------------------- 71 // Plot algorithm v0 72 //---------------------------------------------------------------------- 73 /** 74 * Updates the data points of a parametric curve. This version is used if {@link JXG.Curve#doadvancedplot} is <tt>false</tt>. 75 * @param {JXG.Curve} curve JSXGraph curve element 76 * @param {Number} mi Left bound of curve 77 * @param {Number} ma Right bound of curve 78 * @param {Number} len Number of data points 79 * @returns {JXG.Curve} Reference to the curve object. 80 */ 81 updateParametricCurveNaive: function (curve, mi, ma, len) { 82 var i, 83 t, 84 suspendUpdate = false, 85 stepSize = (ma - mi) / len; 86 87 for (i = 0; i < len; i++) { 88 t = mi + i * stepSize; 89 // The last parameter prevents rounding in usr2screen(). 90 curve.points[i].setCoordinates( 91 Const.COORDS_BY_USER, 92 [curve.X(t, suspendUpdate), curve.Y(t, suspendUpdate)], 93 false 94 ); 95 curve.points[i]._t = t; 96 suspendUpdate = true; 97 } 98 return curve; 99 }, 100 101 //---------------------------------------------------------------------- 102 // Plot algorithm v1 103 //---------------------------------------------------------------------- 104 /** 105 * Crude and cheap test if the segment defined by the two points <tt>(x0, y0)</tt> and <tt>(x1, y1)</tt> is 106 * outside the viewport of the board. All parameters have to be given in screen coordinates. 107 * 108 * @private 109 * @deprecated 110 * @param {Number} x0 111 * @param {Number} y0 112 * @param {Number} x1 113 * @param {Number} y1 114 * @param {JXG.Board} board 115 * @returns {Boolean} <tt>true</tt> if the given segment is outside the visible area. 116 */ 117 isSegmentOutside: function (x0, y0, x1, y1, board) { 118 return ( 119 (y0 < 0 && y1 < 0) || 120 (y0 > board.canvasHeight && y1 > board.canvasHeight) || 121 (x0 < 0 && x1 < 0) || 122 (x0 > board.canvasWidth && x1 > board.canvasWidth) 123 ); 124 }, 125 126 /** 127 * Compares the absolute value of <tt>dx</tt> with <tt>MAXX</tt> and the absolute value of <tt>dy</tt> 128 * with <tt>MAXY</tt>. 129 * 130 * @private 131 * @deprecated 132 * @param {Number} dx 133 * @param {Number} dy 134 * @param {Number} MAXX 135 * @param {Number} MAXY 136 * @returns {Boolean} <tt>true</tt>, if <tt>|dx| < MAXX</tt> and <tt>|dy| < MAXY</tt>. 137 */ 138 isDistOK: function (dx, dy, MAXX, MAXY) { 139 return Math.abs(dx) < MAXX && Math.abs(dy) < MAXY && !isNaN(dx + dy); 140 }, 141 142 /** 143 * @private 144 * @deprecated 145 */ 146 isSegmentDefined: function (x0, y0, x1, y1) { 147 return !(isNaN(x0 + y0) && isNaN(x1 + y1)); 148 }, 149 150 /** 151 * Updates the data points of a parametric curve. This version is used if {@link JXG.Curve#doadvancedplot} is <tt>true</tt>. 152 * Since 0.99 this algorithm is deprecated. It still can be used if {@link JXG.Curve#doadvancedplotold} is <tt>true</tt>. 153 * 154 * @deprecated 155 * @param {JXG.Curve} curve JSXGraph curve element 156 * @param {Number} mi Left bound of curve 157 * @param {Number} ma Right bound of curve 158 * @returns {JXG.Curve} Reference to the curve object. 159 */ 160 updateParametricCurveOld: function (curve, mi, ma) { 161 var i, t, d, x, y, 162 x0, y0,// t0, 163 top, 164 depth, 165 MAX_DEPTH, 166 MAX_XDIST, 167 MAX_YDIST, 168 suspendUpdate = false, 169 po = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 170 dyadicStack = [], 171 depthStack = [], 172 pointStack = [], 173 divisors = [], 174 distOK = false, 175 j = 0, 176 distFromLine = function (p1, p2, p0) { 177 var lbda, 178 x0 = p0[1] - p1[1], 179 y0 = p0[2] - p1[2], 180 x1 = p2[0] - p1[1], 181 y1 = p2[1] - p1[2], 182 den = x1 * x1 + y1 * y1; 183 184 if (den >= Mat.eps) { 185 lbda = (x0 * x1 + y0 * y1) / den; 186 if (lbda > 0) { 187 if (lbda <= 1) { 188 x0 -= lbda * x1; 189 y0 -= lbda * y1; 190 // lbda = 1.0; 191 } else { 192 x0 -= x1; 193 y0 -= y1; 194 } 195 } 196 } 197 return Mat.hypot(x0, y0); 198 }; 199 200 JXG.deprecated("Curve.updateParametricCurveOld()"); 201 202 if (curve.board.updateQuality === curve.board.BOARD_QUALITY_LOW) { 203 MAX_DEPTH = 15; 204 MAX_XDIST = 10; // 10 205 MAX_YDIST = 10; // 10 206 } else { 207 MAX_DEPTH = 21; 208 MAX_XDIST = 0.7; // 0.7 209 MAX_YDIST = 0.7; // 0.7 210 } 211 212 divisors[0] = ma - mi; 213 for (i = 1; i < MAX_DEPTH; i++) { 214 divisors[i] = divisors[i - 1] * 0.5; 215 } 216 217 i = 1; 218 dyadicStack[0] = 1; 219 depthStack[0] = 0; 220 221 t = mi; 222 po.setCoordinates( 223 Const.COORDS_BY_USER, 224 [curve.X(t, suspendUpdate), curve.Y(t, suspendUpdate)], 225 false 226 ); 227 228 // Now, there was a first call to the functions defining the curve. 229 // Defining elements like sliders have been evaluated. 230 // Therefore, we can set suspendUpdate to false, so that these defining elements 231 // need not be evaluated anymore for the rest of the plotting. 232 suspendUpdate = true; 233 x0 = po.scrCoords[1]; 234 y0 = po.scrCoords[2]; 235 // t0 = t; 236 237 t = ma; 238 po.setCoordinates( 239 Const.COORDS_BY_USER, 240 [curve.X(t, suspendUpdate), curve.Y(t, suspendUpdate)], 241 false 242 ); 243 x = po.scrCoords[1]; 244 y = po.scrCoords[2]; 245 246 pointStack[0] = [x, y]; 247 248 top = 1; 249 depth = 0; 250 251 curve.points = []; 252 curve.points[j++] = new Coords(Const.COORDS_BY_SCREEN, [x0, y0], curve.board, false); 253 254 do { 255 distOK = 256 this.isDistOK(x - x0, y - y0, MAX_XDIST, MAX_YDIST) || 257 this.isSegmentOutside(x0, y0, x, y, curve.board); 258 while ( 259 depth < MAX_DEPTH && 260 (!distOK || depth < 6) && 261 (depth <= 7 || this.isSegmentDefined(x0, y0, x, y)) 262 ) { 263 // We jump out of the loop if 264 // * depth>=MAX_DEPTH or 265 // * (depth>=6 and distOK) or 266 // * (depth>7 and segment is not defined) 267 268 dyadicStack[top] = i; 269 depthStack[top] = depth; 270 pointStack[top] = [x, y]; 271 top += 1; 272 273 i = 2 * i - 1; 274 // Here, depth is increased and may reach MAX_DEPTH 275 depth++; 276 // In that case, t is undefined and we will see a jump in the curve. 277 t = mi + i * divisors[depth]; 278 279 po.setCoordinates( 280 Const.COORDS_BY_USER, 281 [curve.X(t, suspendUpdate), curve.Y(t, suspendUpdate)], 282 false, 283 true 284 ); 285 x = po.scrCoords[1]; 286 y = po.scrCoords[2]; 287 distOK = 288 this.isDistOK(x - x0, y - y0, MAX_XDIST, MAX_YDIST) || 289 this.isSegmentOutside(x0, y0, x, y, curve.board); 290 } 291 292 if (j > 1) { 293 d = distFromLine( 294 curve.points[j - 2].scrCoords, 295 [x, y], 296 curve.points[j - 1].scrCoords 297 ); 298 if (d < 0.015) { 299 j -= 1; 300 } 301 } 302 303 curve.points[j] = new Coords(Const.COORDS_BY_SCREEN, [x, y], curve.board, false); 304 curve.points[j]._t = t; 305 j += 1; 306 307 x0 = x; 308 y0 = y; 309 // t0 = t; 310 311 top -= 1; 312 x = pointStack[top][0]; 313 y = pointStack[top][1]; 314 depth = depthStack[top] + 1; 315 i = dyadicStack[top] * 2; 316 } while (top > 0 && j < 500000); 317 318 curve.numberPoints = curve.points.length; 319 320 return curve; 321 }, 322 323 //---------------------------------------------------------------------- 324 // Plot algorithm v2 325 //---------------------------------------------------------------------- 326 327 /** 328 * Add a point to the curve plot. If the new point is too close to the previously inserted point, 329 * it is skipped. 330 * Used in {@link JXG.Curve._plotRecursive}. 331 * 332 * @private 333 * @param {JXG.Coords} pnt Coords to add to the list of points 334 */ 335 _insertPoint_v2: function (curve, pnt, t) { 336 var lastReal = !isNaN(this._lastCrds[1] + this._lastCrds[2]), // The last point was real 337 newReal = !isNaN(pnt.scrCoords[1] + pnt.scrCoords[2]), // New point is real point 338 cw = curve.board.canvasWidth, 339 ch = curve.board.canvasHeight, 340 off = 500; 341 342 newReal = 343 newReal && 344 pnt.scrCoords[1] > -off && 345 pnt.scrCoords[2] > -off && 346 pnt.scrCoords[1] < cw + off && 347 pnt.scrCoords[2] < ch + off; 348 349 /* 350 * Prevents two consecutive NaNs or points wich are too close 351 */ 352 if ( 353 (!newReal && lastReal) || 354 (newReal && 355 (!lastReal || 356 Math.abs(pnt.scrCoords[1] - this._lastCrds[1]) > 0.7 || 357 Math.abs(pnt.scrCoords[2] - this._lastCrds[2]) > 0.7)) 358 ) { 359 pnt._t = t; 360 curve.points.push(pnt); 361 this._lastCrds = pnt.copy('scrCoords'); 362 } 363 }, 364 365 /** 366 * Check if there is a single NaN function value at t0. 367 * @param {*} curve 368 * @param {*} t0 369 * @returns {Boolean} true if there is a second NaN point close by, false otherwise 370 */ 371 neighborhood_isNaN_v2: function (curve, t0) { 372 var is_undef, 373 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 374 t, 375 p; 376 377 t = t0 + Mat.eps; 378 pnt.setCoordinates(Const.COORDS_BY_USER, [curve.X(t, true), curve.Y(t, true)], false); 379 p = pnt.usrCoords; 380 is_undef = isNaN(p[1] + p[2]); 381 if (!is_undef) { 382 t = t0 - Mat.eps; 383 pnt.setCoordinates( 384 Const.COORDS_BY_USER, 385 [curve.X(t, true), curve.Y(t, true)], 386 false 387 ); 388 p = pnt.usrCoords; 389 is_undef = isNaN(p[1] + p[2]); 390 if (!is_undef) { 391 return false; 392 } 393 } 394 return true; 395 }, 396 397 /** 398 * Investigate a function term at the bounds of intervals where 399 * the function is not defined, e.g. log(x) at x = 0. 400 * 401 * c is between a and b 402 * @private 403 * @param {JXG.Curve} curve JSXGraph curve element 404 * @param {Array} a Screen coordinates of the left interval bound 405 * @param {Array} b Screen coordinates of the right interval bound 406 * @param {Array} c Screen coordinates of the bisection point at (ta + tb) / 2 407 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 408 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 409 * @param {Number} tc (ta + tb) / 2 = tc. Parameter which evaluates to b, i.e. [1, X(tc), Y(tc)] = c in screen coordinates 410 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 411 * @returns {JXG.Boolean} true if the point is inserted and the recursion should stop, false otherwise. 412 */ 413 _borderCase: function (curve, a, b, c, ta, tb, tc, depth) { 414 var t, pnt, p, 415 p_good = null, 416 j, 417 max_it = 30, 418 is_undef = false, 419 t_nan, t_real;// t_real2; 420 // dx, dy, 421 // vx, vy, vx2, vy2; 422 // asymptote; 423 424 if (depth <= 1) { 425 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 426 // Test if there is a single undefined point. 427 // If yes, we ignore it. 428 if ( 429 isNaN(a[1] + a[2]) && 430 !isNaN(c[1] + c[2]) && 431 !this.neighborhood_isNaN_v2(curve, ta) 432 ) { 433 return false; 434 } 435 if ( 436 isNaN(b[1] + b[2]) && 437 !isNaN(c[1] + c[2]) && 438 !this.neighborhood_isNaN_v2(curve, tb) 439 ) { 440 return false; 441 } 442 if ( 443 isNaN(c[1] + c[2]) && 444 (!isNaN(a[1] + a[2]) || !isNaN(b[1] + b[2])) && 445 !this.neighborhood_isNaN_v2(curve, tc) 446 ) { 447 return false; 448 } 449 450 j = 0; 451 // Bisect a, b and c until the point t_real is inside of the definition interval 452 // and as close as possible at the boundary. 453 // t_real2 is the second closest point. 454 do { 455 // There are four cases: 456 // a | c | b 457 // --------------- 458 // inf | R | R 459 // R | R | inf 460 // inf | inf | R 461 // R | inf | inf 462 // 463 if (isNaN(a[1] + a[2]) && !isNaN(c[1] + c[2])) { 464 t_nan = ta; 465 t_real = tc; 466 // t_real2 = tb; 467 } else if (isNaN(b[1] + b[2]) && !isNaN(c[1] + c[2])) { 468 t_nan = tb; 469 t_real = tc; 470 // t_real2 = ta; 471 } else if (isNaN(c[1] + c[2]) && !isNaN(b[1] + b[2])) { 472 t_nan = tc; 473 t_real = tb; 474 // t_real2 = tb + (tb - tc); 475 } else if (isNaN(c[1] + c[2]) && !isNaN(a[1] + a[2])) { 476 t_nan = tc; 477 t_real = ta; 478 // t_real2 = ta - (tc - ta); 479 } else { 480 return false; 481 } 482 t = 0.5 * (t_nan + t_real); 483 pnt.setCoordinates( 484 Const.COORDS_BY_USER, 485 [curve.X(t, true), curve.Y(t, true)], 486 false 487 ); 488 p = pnt.usrCoords; 489 490 is_undef = isNaN(p[1] + p[2]); 491 if (is_undef) { 492 t_nan = t; 493 } else { 494 // t_real2 = t_real; 495 t_real = t; 496 } 497 ++j; 498 } while (is_undef && j < max_it); 499 500 // If bisection was successful, take this point. 501 // Useful only for general curves, for function graph 502 // the code below overwrite p_good from here. 503 if (j < max_it) { 504 p_good = p.slice(); 505 c = p.slice(); 506 t_real = t; 507 } 508 509 // OK, bisection has been done now. 510 // t_real contains the closest inner point to the border of the interval we could find. 511 // t_real2 is the second nearest point to this boundary. 512 // Now we approximate the derivative by computing the slope of the line through these two points 513 // and test if it is "infinite", i.e larger than 400 in absolute values. 514 // 515 // vx = curve.X(t_real, true); 516 // vx2 = curve.X(t_real2, true); 517 // vy = curve.Y(t_real, true); 518 // vy2 = curve.Y(t_real2, true); 519 // dx = (vx - vx2) / (t_real - t_real2); 520 // dy = (vy - vy2) / (t_real - t_real2); 521 522 if (p_good !== null) { 523 this._insertPoint_v2( 524 curve, 525 new Coords(Const.COORDS_BY_USER, p_good, curve.board, false) 526 ); 527 return true; 528 } 529 } 530 return false; 531 }, 532 533 /** 534 * Recursive interval bisection algorithm for curve plotting. 535 * Used in {@link JXG.Curve.updateParametricCurve}. 536 * @private 537 * @deprecated 538 * @param {JXG.Curve} curve JSXGraph curve element 539 * @param {Array} a Screen coordinates of the left interval bound 540 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 541 * @param {Array} b Screen coordinates of the right interval bound 542 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 543 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 544 * @param {Number} delta If the distance of the bisection point at (ta + tb) / 2 from the point (a + b) / 2 is less then delta, 545 * the segment [a,b] is regarded as straight line. 546 * @returns {JXG.Curve} Reference to the curve object. 547 */ 548 _plotRecursive_v2: function (curve, a, ta, b, tb, depth, delta) { 549 var tc, 550 c, 551 ds, 552 mindepth = 0, 553 isSmooth, 554 isJump, 555 isCusp, 556 cusp_threshold = 0.5, 557 jump_threshold = 0.99, 558 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 559 560 if (curve.numberPoints > 65536) { 561 return; 562 } 563 564 // Test if the function is undefined in an interval 565 if (depth < this.nanLevel && this._isUndefined(curve, a, ta, b, tb)) { 566 return this; 567 } 568 569 if (depth < this.nanLevel && this._isOutside(a, ta, b, tb, curve.board)) { 570 return this; 571 } 572 573 tc = (ta + tb) * 0.5; 574 pnt.setCoordinates(Const.COORDS_BY_USER, [curve.X(tc, true), curve.Y(tc, true)], false); 575 c = pnt.scrCoords; 576 577 if (this._borderCase(curve, a, b, c, ta, tb, tc, depth)) { 578 return this; 579 } 580 581 ds = this._triangleDists(a, b, c); // returns [d_ab, d_ac, d_cb, d_cd] 582 583 isSmooth = depth < this.smoothLevel && ds[3] < delta; 584 585 isJump = 586 ( 587 depth <= this.jumpLevel && (isNaN(ds[0]) || isNaN(ds[1]) || isNaN(ds[2])) 588 ) || ( 589 depth < this.jumpLevel && 590 ( 591 ds[2] > jump_threshold * ds[0] || 592 ds[1] > jump_threshold * ds[0] || 593 ds[0] === Infinity || 594 ds[1] === Infinity || 595 ds[2] === Infinity 596 ) 597 ); 598 599 isCusp = depth < this.smoothLevel + 2 && ds[0] < cusp_threshold * (ds[1] + ds[2]); 600 601 if (isCusp) { 602 mindepth = 0; 603 isSmooth = false; 604 } 605 606 --depth; 607 608 if (isJump) { 609 this._insertPoint_v2( 610 curve, 611 new Coords(Const.COORDS_BY_SCREEN, [NaN, NaN], curve.board, false), 612 tc 613 ); 614 } else if (depth <= mindepth || isSmooth) { 615 this._insertPoint_v2(curve, pnt, tc); 616 //if (this._borderCase(a, b, c, ta, tb, tc, depth)) {} 617 } else { 618 this._plotRecursive_v2(curve, a, ta, c, tc, depth, delta); 619 620 if (!isNaN(pnt.scrCoords[1] + pnt.scrCoords[2])) { 621 this._insertPoint_v2(curve, pnt, tc); 622 } 623 624 this._plotRecursive_v2(curve, c, tc, b, tb, depth, delta); 625 } 626 627 return this; 628 }, 629 630 /** 631 * Updates the data points of a parametric curve. This version is used if {@link JXG.Curve#plotVersion} is <tt>3</tt>. 632 * 633 * @param {JXG.Curve} curve JSXGraph curve element 634 * @param {Number} mi Left bound of curve 635 * @param {Number} ma Right bound of curve 636 * @returns {JXG.Curve} Reference to the curve object. 637 */ 638 updateParametricCurve_v2: function (curve, mi, ma) { 639 var ta, tb, 640 a, b, 641 suspendUpdate = false, 642 pa = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 643 pb = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 644 depth, 645 delta, 646 w2, 647 // h2, 648 bbox, ret_arr; 649 650 //console.time('plot'); 651 // Switching BOARD_QUALITY_LOW/HIGH makes gliders jump 652 // if (curve.board.updateQuality === curve.board.BOARD_QUALITY_LOW) { 653 // depth = curve.evalVisProp('recursiondepthlow') || 13; 654 // delta = 2; 655 // // this.smoothLevel = 5; //depth - 7; 656 // this.smoothLevel = depth - 6; 657 // this.jumpLevel = 3; 658 // } else { 659 depth = curve.evalVisProp('recursiondepthhigh') || 17; 660 delta = 2; 661 // smoothLevel has to be small for graphs in a huge interval. 662 // this.smoothLevel = 3; //depth - 7; // 9 663 this.smoothLevel = depth - 9; // 9 664 this.jumpLevel = 2; 665 // } 666 this.nanLevel = depth - 4; 667 668 curve.points = []; 669 670 if (this.xterm === 'x') { 671 // For function graphs we can restrict the plot interval 672 // to the visible area + plus margin 673 bbox = curve.board.getBoundingBox(); 674 w2 = (bbox[2] - bbox[0]) * 0.3; 675 // h2 = (bbox[1] - bbox[3]) * 0.3; 676 ta = Math.max(mi, bbox[0] - w2); 677 tb = Math.min(ma, bbox[2] + w2); 678 } else { 679 ta = mi; 680 tb = ma; 681 } 682 pa.setCoordinates( 683 Const.COORDS_BY_USER, 684 [curve.X(ta, suspendUpdate), curve.Y(ta, suspendUpdate)], 685 false 686 ); 687 688 // The first function calls of X() and Y() are done. We can now 689 // switch `suspendUpdate` on. If supported by the functions, this 690 // avoids for the rest of the plotting algorithm, evaluation of any 691 // parent elements. 692 suspendUpdate = true; 693 694 pb.setCoordinates( 695 Const.COORDS_BY_USER, 696 [curve.X(tb, suspendUpdate), curve.Y(tb, suspendUpdate)], 697 false 698 ); 699 700 // Find start and end points of the visible area (plus a certain margin) 701 ret_arr = this._findStartPoint(curve, pa.scrCoords, ta, pb.scrCoords, tb); 702 pa.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 703 ta = ret_arr[1]; 704 ret_arr = this._findStartPoint(curve, pb.scrCoords, tb, pa.scrCoords, ta); 705 pb.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 706 tb = ret_arr[1]; 707 708 // Save the visible area. 709 // This can be used in Curve.hasPoint(). 710 this._visibleArea = [ta, tb]; 711 712 // Start recursive plotting algorithm 713 a = pa.copy('scrCoords'); 714 b = pb.copy('scrCoords'); 715 pa._t = ta; 716 curve.points.push(pa); 717 this._lastCrds = pa.copy('scrCoords'); // Used in _insertPoint 718 this._plotRecursive_v2(curve, a, ta, b, tb, depth, delta); 719 pb._t = tb; 720 curve.points.push(pb); 721 722 curve.numberPoints = curve.points.length; 723 //console.timeEnd('plot'); 724 725 return curve; 726 }, 727 728 //---------------------------------------------------------------------- 729 // Plot algorithm v3 730 //---------------------------------------------------------------------- 731 /** 732 * 733 * @param {JXG.Curve} curve JSXGraph curve element 734 * @param {*} pnt 735 * @param {*} t 736 * @param {*} depth 737 * @param {*} limes 738 * @private 739 */ 740 _insertLimesPoint: function (curve, pnt, t, depth, limes) { 741 var p0, p1, p2; 742 743 // Ignore jump point if it follows limes 744 if ( 745 (Math.abs(this._lastUsrCrds[1]) === Infinity && 746 Math.abs(limes.left_x) === Infinity) || 747 (Math.abs(this._lastUsrCrds[2]) === Infinity && Math.abs(limes.left_y) === Infinity) 748 ) { 749 // console.log("SKIP:", pnt.usrCoords, this._lastUsrCrds, limes); 750 return; 751 } 752 753 // // Ignore jump left from limes 754 // if (Math.abs(limes.left_x) > 100 * Math.abs(this._lastUsrCrds[1])) { 755 // x = Math.sign(limes.left_x) * Infinity; 756 // } else { 757 // x = limes.left_x; 758 // } 759 // if (Math.abs(limes.left_y) > 100 * Math.abs(this._lastUsrCrds[2])) { 760 // y = Math.sign(limes.left_y) * Infinity; 761 // } else { 762 // y = limes.left_y; 763 // } 764 // //pnt.setCoordinates(Const.COORDS_BY_USER, [x, y], false); 765 766 // Add points at a jump. pnt contains [NaN, NaN] 767 //console.log("Add", t, pnt.usrCoords, limes, depth) 768 p0 = new Coords(Const.COORDS_BY_USER, [limes.left_x, limes.left_y], curve.board); 769 p0._t = t; 770 curve.points.push(p0); 771 772 if ( 773 !isNaN(limes.left_x) && 774 !isNaN(limes.left_y) && 775 !isNaN(limes.right_x) && 776 !isNaN(limes.right_y) && 777 (Math.abs(limes.left_x - limes.right_x) > Mat.eps || 778 Math.abs(limes.left_y - limes.right_y) > Mat.eps) 779 ) { 780 p1 = new Coords(Const.COORDS_BY_SCREEN, pnt, curve.board); 781 p1._t = t; 782 curve.points.push(p1); 783 } 784 785 p2 = new Coords(Const.COORDS_BY_USER, [limes.right_x, limes.right_y], curve.board); 786 p2._t = t; 787 curve.points.push(p2); 788 this._lastScrCrds = p2.copy('scrCoords'); 789 this._lastUsrCrds = p2.copy('usrCoords'); 790 }, 791 792 /** 793 * Add a point to the curve plot. If the new point is too close to the previously inserted point, 794 * it is skipped. 795 * Used in {@link JXG.Curve._plotRecursive}. 796 * 797 * @private 798 * @param {JXG.Curve} curve JSXGraph curve element 799 * @param {JXG.Coords} pnt Coords to add to the list of points 800 */ 801 _insertPoint: function (curve, pnt, t, depth, limes) { 802 var last_is_real = !isNaN(this._lastScrCrds[1] + this._lastScrCrds[2]), // The last point was real 803 point_is_real = !isNaN(pnt[1] + pnt[2]), // New point is real point 804 cw = curve.board.canvasWidth, 805 ch = curve.board.canvasHeight, 806 p, 807 near = 0.8, 808 off = 500; 809 810 if (Type.exists(limes)) { 811 this._insertLimesPoint(curve, pnt, t, depth, limes); 812 return; 813 } 814 815 // Check if point has real coordinates and 816 // coordinates are not too far away from canvas. 817 point_is_real = 818 point_is_real && 819 pnt[1] > -off && 820 pnt[2] > -off && 821 pnt[1] < cw + off && 822 pnt[2] < ch + off; 823 824 // Prevent two consecutive NaNs 825 if (!last_is_real && !point_is_real) { 826 return; 827 } 828 829 // Prevent two consecutive points which are too close 830 if ( 831 point_is_real && 832 last_is_real && 833 Math.abs(pnt[1] - this._lastScrCrds[1]) < near && 834 Math.abs(pnt[2] - this._lastScrCrds[2]) < near 835 ) { 836 return; 837 } 838 839 // Prevent two consecutive points at infinity (either direction) 840 if ( 841 (Math.abs(pnt[1]) === Infinity && Math.abs(this._lastUsrCrds[1]) === Infinity) || 842 (Math.abs(pnt[2]) === Infinity && Math.abs(this._lastUsrCrds[2]) === Infinity) 843 ) { 844 return; 845 } 846 847 //console.log("add", t, pnt.usrCoords, depth) 848 // Add regular point 849 p = new Coords(Const.COORDS_BY_SCREEN, pnt, curve.board); 850 p._t = t; 851 curve.points.push(p); 852 this._lastScrCrds = p.copy('scrCoords'); 853 this._lastUsrCrds = p.copy('usrCoords'); 854 }, 855 856 /** 857 * Compute distances in screen coordinates between the points ab, 858 * ac, cb, and cd, where d = (a + b)/2. 859 * cd is used for the smoothness test, ab, ac, cb are used to detect jumps, cusps and poles. 860 * 861 * @private 862 * @param {Array} a Screen coordinates of the left interval bound 863 * @param {Array} b Screen coordinates of the right interval bound 864 * @param {Array} c Screen coordinates of the bisection point at (ta + tb) / 2 865 * @returns {Array} array of distances in screen coordinates between: ab, ac, cb, and cd. 866 */ 867 _triangleDists: function (a, b, c) { 868 var d, d_ab, d_ac, d_cb, d_cd; 869 870 d = [a[0] * b[0], (a[1] + b[1]) * 0.5, (a[2] + b[2]) * 0.5]; 871 872 d_ab = Geometry.distance(a, b, 3); 873 d_ac = Geometry.distance(a, c, 3); 874 d_cb = Geometry.distance(c, b, 3); 875 d_cd = Geometry.distance(c, d, 3); 876 877 return [d_ab, d_ac, d_cb, d_cd]; 878 }, 879 880 /** 881 * Test if the function is undefined on an interval: 882 * If the interval borders a and b are undefined, 20 random values 883 * are tested if they are undefined, too. 884 * Only if all values are undefined, we declare the function to be undefined in this interval. 885 * 886 * @private 887 * @param {JXG.Curve} curve JSXGraph curve element 888 * @param {Array} a Screen coordinates of the left interval bound 889 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 890 * @param {Array} b Screen coordinates of the right interval bound 891 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 892 */ 893 _isUndefined: function (curve, a, ta, b, tb) { 894 var t, i, pnt; 895 896 if (!isNaN(a[1] + a[2]) || !isNaN(b[1] + b[2])) { 897 return false; 898 } 899 900 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 901 902 for (i = 0; i < 20; ++i) { 903 t = ta + Math.random() * (tb - ta); 904 pnt.setCoordinates( 905 Const.COORDS_BY_USER, 906 [curve.X(t, true), curve.Y(t, true)], 907 false 908 ); 909 if (!isNaN(pnt.scrCoords[0] + pnt.scrCoords[1] + pnt.scrCoords[2])) { 910 return false; 911 } 912 } 913 914 return true; 915 }, 916 917 /** 918 * Decide if a path segment is too far from the canvas that we do not need to draw it. 919 * @private 920 * @param {Array} a Screen coordinates of the start point of the segment 921 * @param {Array} ta Curve parameter of a (unused). 922 * @param {Array} b Screen coordinates of the end point of the segment 923 * @param {Array} tb Curve parameter of b (unused). 924 * @param {JXG.Board} board 925 * @returns {Boolean} True if the segment is too far away from the canvas, false otherwise. 926 */ 927 _isOutside: function (a, ta, b, tb, board) { 928 var off = 500, 929 cw = board.canvasWidth, 930 ch = board.canvasHeight; 931 932 return !!( 933 (a[1] < -off && b[1] < -off) || 934 (a[2] < -off && b[2] < -off) || 935 (a[1] > cw + off && b[1] > cw + off) || 936 (a[2] > ch + off && b[2] > ch + off) 937 ); 938 }, 939 940 /** 941 * Decide if a point of a curve is too far from the canvas that we do not need to draw it. 942 * @private 943 * @param {Array} a Screen coordinates of the point 944 * @param {JXG.Board} board 945 * @returns {Boolean} True if the point is too far away from the canvas, false otherwise. 946 */ 947 _isOutsidePoint: function (a, board) { 948 var off = 500, 949 cw = board.canvasWidth, 950 ch = board.canvasHeight; 951 952 return !!(a[1] < -off || a[2] < -off || a[1] > cw + off || a[2] > ch + off); 953 }, 954 955 /** 956 * For a curve c(t) defined on the interval [ta, tb] find the first point 957 * which is in the visible area of the board (plus some outside margin). 958 * <p> 959 * This method is necessary to restrict the recursive plotting algorithm 960 * {@link JXG.Curve._plotRecursive} to the visible area and not waste 961 * recursion to areas far outside of the visible area. 962 * <p> 963 * This method can also be used to find the last visible point 964 * by reversing the input parameters. 965 * 966 * @param {JXG.Curve} curve JSXGraph curve element 967 * @param {Array} ta Curve parameter of a. 968 * @param {Array} b Screen coordinates of the end point of the segment (unused) 969 * @param {Array} tb Curve parameter of b 970 * @return {Array} Array of length two containing the screen ccordinates of 971 * the starting point and the curve parameter at this point. 972 * @private 973 */ 974 _findStartPoint: function (curve, a, ta, b, tb) { 975 // The code below is too unstable. 976 // E.g. [function(t) { return Math.pow(t, 2) * (t + 5) * Math.pow(t - 5, 2); }, -8, 8] 977 // Therefore, we return here. 978 return [a, ta]; 979 980 // var i, 981 // delta, 982 // tc, 983 // td, 984 // z, 985 // isFound, 986 // w2, 987 // h2, 988 // pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 989 // steps = 40, 990 // eps = 0.01, 991 // fnX1, 992 // fnX2, 993 // fnY1, 994 // fnY2, 995 // bbox = curve.board.getBoundingBox(); 996 997 // if (true || !this._isOutsidePoint(a, curve.board)) { 998 // return [a, ta]; 999 // } 1000 // w2 = (bbox[2] - bbox[0]) * 0.3; 1001 // h2 = (bbox[1] - bbox[3]) * 0.3; 1002 // bbox[0] -= w2; 1003 // bbox[1] += h2; 1004 // bbox[2] += w2; 1005 // bbox[3] -= h2; 1006 1007 // delta = (tb - ta) / steps; 1008 // tc = ta + delta; 1009 // isFound = false; 1010 1011 // fnX1 = function (t) { 1012 // return curve.X(t, true) - bbox[0]; 1013 // }; 1014 // fnY1 = function (t) { 1015 // return curve.Y(t, true) - bbox[1]; 1016 // }; 1017 // fnX2 = function (t) { 1018 // return curve.X(t, true) - bbox[2]; 1019 // }; 1020 // fnY2 = function (t) { 1021 // return curve.Y(t, true) - bbox[3]; 1022 // }; 1023 // for (i = 0; i < steps; ++i) { 1024 // // Left border 1025 // z = bbox[0]; 1026 // td = Numerics.root(fnX1, [tc - delta, tc], curve); 1027 // // td = Numerics.fzero(fnX1, [tc - delta, tc], this); 1028 // // console.log("A", tc - delta, tc, td, Math.abs(this.X(td, true) - z)); 1029 // if (Math.abs(curve.X(td, true) - z) < eps) { 1030 // //} * Math.abs(z)) { 1031 // isFound = true; 1032 // break; 1033 // } 1034 // // Top border 1035 // z = bbox[1]; 1036 // td = Numerics.root(fnY1, [tc - delta, tc], curve); 1037 // // td = Numerics.fzero(fnY1, [tc - delta, tc], this); 1038 // // console.log("B", tc - delta, tc, td, Math.abs(this.Y(td, true) - z)); 1039 // if (Math.abs(curve.Y(td, true) - z) < eps) { 1040 // // * Math.abs(z)) { 1041 // isFound = true; 1042 // break; 1043 // } 1044 // // Right border 1045 // z = bbox[2]; 1046 // td = Numerics.root(fnX2, [tc - delta, tc], curve); 1047 // // td = Numerics.fzero(fnX2, [tc - delta, tc], this); 1048 // // console.log("C", tc - delta, tc, td, Math.abs(this.X(td, true) - z)); 1049 // if (Math.abs(curve.X(td, true) - z) < eps) { 1050 // // * Math.abs(z)) { 1051 // isFound = true; 1052 // break; 1053 // } 1054 // // Bottom border 1055 // z = bbox[3]; 1056 // td = Numerics.root(fnY2, [tc - delta, tc], curve); 1057 // // td = Numerics.fzero(fnY2, [tc - delta, tc], this); 1058 // // console.log("D", tc - delta, tc, td, Math.abs(this.Y(td, true) - z)); 1059 // if (Math.abs(curve.Y(td, true) - z) < eps) { 1060 // // * Math.abs(z)) { 1061 // isFound = true; 1062 // break; 1063 // } 1064 // tc += delta; 1065 // } 1066 // if (isFound) { 1067 // pnt.setCoordinates( 1068 // Const.COORDS_BY_USER, 1069 // [curve.X(td, true), curve.Y(td, true)], 1070 // false 1071 // ); 1072 // return [pnt.scrCoords, td]; 1073 // } 1074 // console.log("TODO _findStartPoint", curve.Y.toString(), tc); 1075 // pnt.setCoordinates(Const.COORDS_BY_USER, [curve.X(ta, true), curve.Y(ta, true)], false); 1076 // return [pnt.scrCoords, ta]; 1077 }, 1078 1079 /** 1080 * Investigate a function term at the bounds of intervals where 1081 * the function is not defined, e.g. log(x) at x = 0. 1082 * 1083 * c is inbetween a and b 1084 * 1085 * @param {JXG.Curve} curve JSXGraph curve element 1086 * @param {Array} a Screen coordinates of the left interval bound 1087 * @param {Array} b Screen coordinates of the right interval bound 1088 * @param {Array} c Screen coordinates of the bisection point at (ta + tb) / 2 1089 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 1090 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 1091 * @param {Number} tc (ta + tb) / 2 = tc. Parameter which evaluates to b, i.e. [1, X(tc), Y(tc)] = c in screen coordinates 1092 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 1093 * @returns {JXG.Boolean} true if the point is inserted and the recursion should stop, false otherwise. 1094 * 1095 * @private 1096 */ 1097 _getBorderPos: function (curve, ta, a, tc, c, tb, b) { 1098 var t, pnt, p, j, 1099 max_it = 30, 1100 is_undef = false, 1101 t_good, t_bad; 1102 1103 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 1104 j = 0; 1105 // Bisect a, b and c until the point t_real is inside of the definition interval 1106 // and as close as possible at the boundary. 1107 // (t_real2 is/was the second closest point). 1108 // There are four cases: 1109 // a | c | b 1110 // --------------- 1111 // inf | R | R 1112 // R | R | inf 1113 // inf | inf | R 1114 // R | inf | inf 1115 // 1116 if (isNaN(a[1] + a[2]) && !isNaN(c[1] + c[2])) { 1117 t_bad = ta; 1118 t_good = tc; 1119 } else if (isNaN(b[1] + b[2]) && !isNaN(c[1] + c[2])) { 1120 t_bad = tb; 1121 t_good = tc; 1122 } else if (isNaN(c[1] + c[2]) && !isNaN(b[1] + b[2])) { 1123 t_bad = tc; 1124 t_good = tb; 1125 } else if (isNaN(c[1] + c[2]) && !isNaN(a[1] + a[2])) { 1126 t_bad = tc; 1127 t_good = ta; 1128 } else { 1129 return false; 1130 } 1131 do { 1132 t = 0.5 * (t_good + t_bad); 1133 pnt.setCoordinates( 1134 Const.COORDS_BY_USER, 1135 [curve.X(t, true), curve.Y(t, true)], 1136 false 1137 ); 1138 p = pnt.usrCoords; 1139 is_undef = isNaN(p[1] + p[2]); 1140 if (is_undef) { 1141 t_bad = t; 1142 } else { 1143 t_good = t; 1144 } 1145 ++j; 1146 } while (j < max_it && Math.abs(t_good - t_bad) > Mat.eps); 1147 return t; 1148 }, 1149 1150 /** 1151 * 1152 * @param {JXG.Curve} curve JSXGraph curve element 1153 * @param {Number} ta 1154 * @param {Number} tb 1155 */ 1156 _getCuspPos: function (curve, ta, tb) { 1157 var a = [curve.X(ta, true), curve.Y(ta, true)], 1158 b = [curve.X(tb, true), curve.Y(tb, true)], 1159 max_func = function (t) { 1160 var c = [curve.X(t, true), curve.Y(t, true)]; 1161 return -( 1162 Mat.hypot(a[0] - c[0], a[1] - c[1]) + 1163 Mat.hypot(b[0] - c[0], b[1] - c[1]) 1164 ); 1165 }; 1166 1167 return Numerics.fminbr(max_func, [ta, tb], curve); 1168 }, 1169 1170 /** 1171 * 1172 * @param {JXG.Curve} curve JSXGraph curve element 1173 * @param {Number} ta 1174 * @param {Number} tb 1175 */ 1176 _getJumpPos: function (curve, ta, tb) { 1177 var max_func = function (t) { 1178 var e = Mat.eps * Mat.eps, 1179 c1 = [curve.X(t, true), curve.Y(t, true)], 1180 c2 = [curve.X(t + e, true), curve.Y(t + e, true)]; 1181 return -Math.abs((c2[1] - c1[1]) / (c2[0] - c1[0])); 1182 }; 1183 1184 return Numerics.fminbr(max_func, [ta, tb], curve); 1185 }, 1186 1187 /** 1188 * 1189 * @param {JXG.Curve} curve JSXGraph curve element 1190 * @param {Number} t 1191 * @private 1192 */ 1193 _getLimits: function (curve, t) { 1194 var res, 1195 step = 2 / (curve.maxX() - curve.minX()), 1196 x_l, 1197 x_r, 1198 y_l, 1199 y_r; 1200 1201 // From left 1202 res = Extrapolate.limit(t, -step, curve.X); 1203 x_l = res[0]; 1204 if (res[1] === 'infinite') { 1205 x_l = Math.sign(x_l) * Infinity; 1206 } 1207 1208 res = Extrapolate.limit(t, -step, curve.Y); 1209 y_l = res[0]; 1210 if (res[1] === 'infinite') { 1211 y_l = Math.sign(y_l) * Infinity; 1212 } 1213 1214 // From right 1215 res = Extrapolate.limit(t, step, curve.X); 1216 x_r = res[0]; 1217 if (res[1] === 'infinite') { 1218 x_r = Math.sign(x_r) * Infinity; 1219 } 1220 1221 res = Extrapolate.limit(t, step, curve.Y); 1222 y_r = res[0]; 1223 if (res[1] === 'infinite') { 1224 y_r = Math.sign(y_r) * Infinity; 1225 } 1226 1227 return { 1228 left_x: x_l, 1229 left_y: y_l, 1230 right_x: x_r, 1231 right_y: y_r, 1232 t: t 1233 }; 1234 }, 1235 1236 /** 1237 * 1238 * @param {JXG.Curve} curve JSXGraph curve element 1239 * @param {Array} a 1240 * @param {Number} tc 1241 * @param {Array} c 1242 * @param {Number} tb 1243 * @param {Array} b 1244 * @param {String} may_be_special 1245 * @param {Number} depth 1246 * @private 1247 */ 1248 _getLimes: function (curve, ta, a, tc, c, tb, b, may_be_special, depth) { 1249 var t; 1250 1251 if (may_be_special === 'border') { 1252 t = this._getBorderPos(curve, ta, a, tc, c, tb, b); 1253 } else if (may_be_special === 'cusp') { 1254 t = this._getCuspPos(curve, ta, tb); 1255 } else if (may_be_special === 'jump') { 1256 t = this._getJumpPos(curve, ta, tb); 1257 } 1258 return this._getLimits(curve, t); 1259 }, 1260 1261 /** 1262 * Recursive interval bisection algorithm for curve plotting. 1263 * Used in {@link JXG.Curve.updateParametricCurve}. 1264 * @private 1265 * @param {JXG.Curve} curve JSXGraph curve element 1266 * @param {Array} a Screen coordinates of the left interval bound 1267 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 1268 * @param {Array} b Screen coordinates of the right interval bound 1269 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 1270 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 1271 * @param {Number} delta If the distance of the bisection point at (ta + tb) / 2 from the point (a + b) / 2 is less then delta, 1272 * the segment [a,b] is regarded as straight line. 1273 * @returns {JXG.Curve} Reference to the curve object. 1274 */ 1275 _plotNonRecursive: function (curve, a, ta, b, tb, d) { 1276 var tc, 1277 c, 1278 ds, 1279 mindepth = 0, 1280 limes = null, 1281 a_nan, 1282 b_nan, 1283 isSmooth = false, 1284 may_be_special = "", 1285 x, 1286 y, 1287 oc, 1288 depth, 1289 ds0, 1290 stack = [], 1291 stack_length = 0, 1292 item; 1293 1294 oc = curve.board.origin.scrCoords; 1295 stack[stack_length++] = [a, ta, b, tb, d, Infinity]; 1296 while (stack_length > 0) { 1297 // item = stack.pop(); 1298 item = stack[--stack_length]; 1299 a = item[0]; 1300 ta = item[1]; 1301 b = item[2]; 1302 tb = item[3]; 1303 depth = item[4]; 1304 ds0 = item[5]; 1305 1306 isSmooth = false; 1307 may_be_special = ""; 1308 limes = null; 1309 //console.log(stack.length, item) 1310 1311 if (curve.points.length > 65536) { 1312 return; 1313 } 1314 1315 if (depth < this.nanLevel) { 1316 // Test if the function is undefined in the whole interval [ta, tb] 1317 if (this._isUndefined(curve, a, ta, b, tb)) { 1318 continue; 1319 } 1320 // Test if the graph is far outside the visible are for the interval [ta, tb] 1321 if (this._isOutside(a, ta, b, tb, curve.board)) { 1322 continue; 1323 } 1324 } 1325 1326 tc = (ta + tb) * 0.5; 1327 1328 // Screen coordinates of point at tc 1329 x = curve.X(tc, true); 1330 y = curve.Y(tc, true); 1331 c = [1, oc[1] + x * curve.board.unitX, oc[2] - y * curve.board.unitY]; 1332 ds = this._triangleDists(a, b, c); // returns [d_ab, d_ac, d_cb, d_cd] 1333 1334 a_nan = isNaN(a[1] + a[2]); 1335 b_nan = isNaN(b[1] + b[2]); 1336 if ((a_nan && !b_nan) || (!a_nan && b_nan)) { 1337 may_be_special = 'border'; 1338 } else if ( 1339 ds[0] > 0.66 * ds0 || 1340 ds[0] < this.cusp_threshold * (ds[1] + ds[2]) || 1341 ds[1] > 5 * ds[2] || 1342 ds[2] > 5 * ds[1] 1343 ) { 1344 may_be_special = 'cusp'; 1345 } else if ( 1346 ds[2] > this.jump_threshold * ds[0] || 1347 ds[1] > this.jump_threshold * ds[0] || 1348 ds[0] === Infinity || 1349 ds[1] === Infinity || 1350 ds[2] === Infinity 1351 ) { 1352 may_be_special = 'jump'; 1353 } 1354 isSmooth = 1355 may_be_special === "" && 1356 depth < this.smoothLevel && 1357 ds[3] < this.smooth_threshold; 1358 1359 if (depth < this.testLevel && !isSmooth) { 1360 if (may_be_special === "") { 1361 isSmooth = true; 1362 } else { 1363 limes = this._getLimes(curve, ta, a, tc, c, tb, b, may_be_special, depth); 1364 } 1365 } 1366 1367 if (limes !== null) { 1368 c = [1, NaN, NaN]; 1369 this._insertPoint(curve, c, tc, depth, limes); 1370 } else if (depth <= mindepth || isSmooth) { 1371 this._insertPoint(curve, c, tc, depth, null); 1372 } else { 1373 stack[stack_length++] = [c, tc, b, tb, depth - 1, ds[0]]; 1374 stack[stack_length++] = [a, ta, c, tc, depth - 1, ds[0]]; 1375 } 1376 } 1377 1378 return this; 1379 }, 1380 1381 /** 1382 * Updates the data points of a parametric curve. This version is used if {@link JXG.Curve#plotVersion} is <tt>3</tt>. 1383 * This is an experimental plot version, <b>not recommended</b> to be used. 1384 * @param {JXG.Curve} curve JSXGraph curve element 1385 * @param {Number} mi Left bound of curve 1386 * @param {Number} ma Right bound of curve 1387 * @returns {JXG.Curve} Reference to the curve object. 1388 */ 1389 updateParametricCurve_v3: function (curve, mi, ma) { 1390 var ta, 1391 tb, 1392 a, 1393 b, 1394 suspendUpdate = false, 1395 pa = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 1396 pb = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 1397 depth, 1398 w2, // h2, 1399 bbox, 1400 ret_arr; 1401 1402 // console.log("-----------------------------------------------------------"); 1403 // console.time('plot'); 1404 // if (curve.board.updateQuality === curve.board.BOARD_QUALITY_LOW) { 1405 // depth = curve.evalVisProp('recursiondepthlow') || 14; 1406 // } else { 1407 depth = curve.evalVisProp('recursiondepthhigh') || 17; 1408 // } 1409 1410 // smoothLevel has to be small for graphs in a huge interval. 1411 this.smoothLevel = 7; //depth - 10; 1412 this.nanLevel = depth - 4; 1413 this.testLevel = 4; 1414 this.cusp_threshold = 0.5; 1415 this.jump_threshold = 0.99; 1416 this.smooth_threshold = 2; 1417 1418 curve.points = []; 1419 1420 if (curve.xterm === 'x') { 1421 // For function graphs we can restrict the plot interval 1422 // to the visible area +plus margin 1423 bbox = curve.board.getBoundingBox(); 1424 w2 = (bbox[2] - bbox[0]) * 0.3; 1425 //h2 = (bbox[1] - bbox[3]) * 0.3; 1426 ta = Math.max(mi, bbox[0] - w2); 1427 tb = Math.min(ma, bbox[2] + w2); 1428 } else { 1429 ta = mi; 1430 tb = ma; 1431 } 1432 pa.setCoordinates( 1433 Const.COORDS_BY_USER, 1434 [curve.X(ta, suspendUpdate), curve.Y(ta, suspendUpdate)], 1435 false 1436 ); 1437 1438 // The first function calls of X() and Y() are done. We can now 1439 // switch `suspendUpdate` on. If supported by the functions, this 1440 // avoids for the rest of the plotting algorithm, evaluation of any 1441 // parent elements. 1442 suspendUpdate = true; 1443 1444 pb.setCoordinates( 1445 Const.COORDS_BY_USER, 1446 [curve.X(tb, suspendUpdate), curve.Y(tb, suspendUpdate)], 1447 false 1448 ); 1449 1450 // Find start and end points of the visible area (plus a certain margin) 1451 ret_arr = this._findStartPoint(curve, pa.scrCoords, ta, pb.scrCoords, tb); 1452 pa.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 1453 ta = ret_arr[1]; 1454 ret_arr = this._findStartPoint(curve, pb.scrCoords, tb, pa.scrCoords, ta); 1455 pb.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 1456 tb = ret_arr[1]; 1457 1458 // Save the visible area. 1459 // This can be used in Curve.hasPoint(). 1460 this._visibleArea = [ta, tb]; 1461 1462 // Start recursive plotting algorithm 1463 a = pa.copy('scrCoords'); 1464 b = pb.copy('scrCoords'); 1465 pa._t = ta; 1466 curve.points.push(pa); 1467 this._lastScrCrds = pa.copy('scrCoords'); // Used in _insertPoint 1468 this._lastUsrCrds = pa.copy('usrCoords'); // Used in _insertPoint 1469 1470 this._plotNonRecursive(curve, a, ta, b, tb, depth); 1471 1472 pb._t = tb; 1473 curve.points.push(pb); 1474 1475 curve.numberPoints = curve.points.length; 1476 // console.timeEnd('plot'); 1477 // console.log("number of points:", this.numberPoints); 1478 1479 return curve; 1480 }, 1481 1482 //---------------------------------------------------------------------- 1483 // Plot algorithm v4 1484 //---------------------------------------------------------------------- 1485 1486 /** 1487 * TODO 1488 * @param {Array} vec 1489 * @param {Number} le 1490 * @param {Number} level 1491 * @returns Object 1492 * @private 1493 */ 1494 _criticalInterval: function (vec, le, level) { 1495 var i, 1496 j, 1497 le1, 1498 med, 1499 sgn, 1500 sgnChange, 1501 isGroup = false, 1502 abs_vec, 1503 last = -Infinity, 1504 very_small = false, 1505 smooth = false, 1506 group = 0, 1507 groups = [], 1508 types = [], 1509 positions = []; 1510 1511 abs_vec = Statistics.abs(vec); 1512 med = Statistics.median(abs_vec); 1513 1514 if (med < 1.0e-7) { 1515 med = 1.0e-7; 1516 very_small = true; 1517 } else { 1518 med *= this.criticalThreshold; 1519 } 1520 1521 //console.log("Median", med); 1522 for (i = 0; i < le; i++) { 1523 // Start a group if not yet done and 1524 // add position to group 1525 if (abs_vec[i] > med /*&& abs_vec[i] > 0.01*/) { 1526 positions.push({ i: i, v: vec[i], group: group }); 1527 last = i; 1528 if (!isGroup) { 1529 isGroup = true; 1530 } 1531 } else { 1532 if (isGroup && i > last + 4) { 1533 // End the group 1534 if (positions.length > 0) { 1535 groups.push(positions.slice(0)); 1536 } 1537 positions = []; 1538 isGroup = false; 1539 group++; 1540 } 1541 } 1542 } 1543 if (isGroup) { 1544 if (positions.length > 1) { 1545 groups.push(positions.slice(0)); 1546 } 1547 } 1548 1549 if (very_small && groups.length === 0) { 1550 smooth = true; 1551 } 1552 1553 // Decide if there is a singular critical point 1554 // or if a whole interval is problematic. 1555 // The latter is the case if the differences have many sign changes. 1556 for (j = 0; j < groups.length; j++) { 1557 types[j] = 'point'; 1558 le1 = groups[j].length; 1559 if (le1 < 64) { 1560 continue; 1561 } 1562 sgnChange = 0; 1563 sgn = Math.sign(groups[j][0].v); 1564 for (i = 1; i < le1; i++) { 1565 if (Math.sign(groups[j][i].v) !== sgn) { 1566 sgnChange++; 1567 sgn = Math.sign(groups[j][i].v); 1568 } 1569 } 1570 if (sgnChange * 6 > le1) { 1571 types[j] = 'interval'; 1572 } 1573 } 1574 1575 return { smooth: smooth, groups: groups, types: types }; 1576 }, 1577 1578 Component: function () { 1579 this.left_isNaN = false; 1580 this.right_isNaN = false; 1581 this.left_t = null; 1582 this.right_t = null; 1583 this.t_values = []; 1584 this.x_values = []; 1585 this.y_values = []; 1586 this.len = 0; 1587 }, 1588 1589 findComponents: function (curve, mi, ma, steps) { 1590 var i, t, h, 1591 x, y, 1592 components = [], 1593 comp, 1594 comp_nr = 0, 1595 cnt = 0, 1596 cntNaNs = 0, 1597 comp_started = false, 1598 suspended = false; 1599 1600 h = (ma - mi) / steps; 1601 components[comp_nr] = new this.Component(); 1602 comp = components[comp_nr]; 1603 1604 for (i = 0, t = mi; i <= steps; i++, t += h) { 1605 x = curve.X(t, suspended); 1606 y = curve.Y(t, suspended); 1607 1608 if (isNaN(x) || isNaN(y)) { 1609 cntNaNs++; 1610 // Wait for - at least - two consecutive NaNs 1611 // This avoids starting a new component if 1612 // the function value has infinity as intermediate value. 1613 if (cntNaNs > 1 && comp_started) { 1614 // Finalize a component 1615 comp.right_isNaN = true; 1616 comp.right_t = t - h; 1617 comp.len = cnt; 1618 1619 // Prepare a new component 1620 comp_started = false; 1621 comp_nr++; 1622 components[comp_nr] = new this.Component(); 1623 comp = components[comp_nr]; 1624 cntNaNs = 0; 1625 } 1626 } else { 1627 // Now there is a non-NaN entry. 1628 if (!comp_started) { 1629 // Start the component 1630 comp_started = true; 1631 cnt = 0; 1632 if (cntNaNs > 0) { 1633 comp.left_t = t - h; 1634 comp.left_isNaN = true; 1635 } 1636 } 1637 cntNaNs = 0; 1638 // Add the value to the component 1639 comp.t_values[cnt] = t; 1640 comp.x_values[cnt] = x; 1641 comp.y_values[cnt] = y; 1642 cnt++; 1643 } 1644 if (i === 0) { 1645 suspended = true; 1646 } 1647 } 1648 if (comp_started) { 1649 comp.len = cnt; 1650 } else { 1651 components.pop(); 1652 } 1653 1654 return components; 1655 }, 1656 1657 getPointType: function (curve, pos, t_approx, t_values, x_table, y_table, len) { 1658 var x_values = x_table[0], 1659 y_values = y_table[0], 1660 full_len = t_values.length, 1661 result = { 1662 idx: pos, 1663 t: t_approx, //t_values[pos], 1664 x: x_values[pos], 1665 y: y_values[pos], 1666 type: "other" 1667 }; 1668 1669 if (pos < 5) { 1670 result.type = 'borderleft'; 1671 result.idx = 0; 1672 result.t = t_values[0]; 1673 result.x = x_values[0]; 1674 result.y = y_values[0]; 1675 1676 // console.log('Border left', result.t); 1677 return result; 1678 } 1679 if (pos > len - 6) { 1680 result.type = 'borderright'; 1681 result.idx = full_len - 1; 1682 result.t = t_values[full_len - 1]; 1683 result.x = x_values[full_len - 1]; 1684 result.y = y_values[full_len - 1]; 1685 1686 // console.log('Border right', result.t, full_len - 1); 1687 return result; 1688 } 1689 1690 return result; 1691 }, 1692 1693 newtonApprox: function (idx, t, h, level, table) { 1694 var i, 1695 s = 0.0; 1696 for (i = level; i > 0; i--) { 1697 s = ((s + table[i][idx]) * (t - (i - 1) * h)) / i; 1698 } 1699 return s + table[0][idx]; 1700 }, 1701 1702 // Thiele's interpolation formula, 1703 // https://en.wikipedia.org/wiki/Thiele%27s_interpolation_formula 1704 // unused 1705 thiele: function (t, recip, t_values, idx, degree) { 1706 var i, 1707 v = 0.0; 1708 for (i = degree; i > 1; i--) { 1709 v = (t - t_values[idx + i]) / (recip[i][idx + 1] - recip[i - 2][idx + 1] + v); 1710 } 1711 return recip[0][idx + 1] + (t - t_values[idx + 1]) / (recip[1][idx + 1] + v); 1712 }, 1713 1714 differenceMethodExperiments: function (component, curve) { 1715 var i, 1716 level, 1717 le, 1718 up, 1719 t_values = component.t_values, 1720 x_values = component.x_values, 1721 y_values = component.y_values, 1722 x_diffs = [], 1723 y_diffs = [], 1724 x_slopes = [], 1725 y_slopes = [], 1726 x_table = [], 1727 y_table = [], 1728 x_recip = [], 1729 y_recip = [], 1730 h, 1731 numerator, 1732 // x_med, y_med, 1733 foundCriticalPoint = 0, 1734 pos, 1735 ma, 1736 j, 1737 v, 1738 groups, 1739 criticalPoints = []; 1740 1741 h = t_values[1] - t_values[0]; 1742 x_table.push([]); 1743 y_table.push([]); 1744 x_recip.push([]); 1745 y_recip.push([]); 1746 le = y_values.length; 1747 for (i = 0; i < le; i++) { 1748 x_table[0][i] = x_values[i]; 1749 y_table[0][i] = y_values[i]; 1750 x_recip[0][i] = x_values[i]; 1751 y_recip[0][i] = y_values[i]; 1752 } 1753 1754 x_table.push([]); 1755 y_table.push([]); 1756 x_recip.push([]); 1757 y_recip.push([]); 1758 numerator = h; 1759 le = y_values.length - 1; 1760 for (i = 0; i < le; i++) { 1761 x_diffs[i] = x_values[i + 1] - x_values[i]; 1762 y_diffs[i] = y_values[i + 1] - y_values[i]; 1763 x_slopes[i] = x_diffs[i]; 1764 y_slopes[i] = y_diffs[i]; 1765 x_table[1][i] = x_diffs[i]; 1766 y_table[1][i] = y_diffs[i]; 1767 x_recip[1][i] = numerator / x_diffs[i]; 1768 y_recip[1][i] = numerator / y_diffs[i]; 1769 } 1770 le--; 1771 1772 up = Math.min(8, y_values.length - 1); 1773 for (level = 1; level < up; level++) { 1774 x_table.push([]); 1775 y_table.push([]); 1776 x_recip.push([]); 1777 y_recip.push([]); 1778 numerator *= h; 1779 for (i = 0; i < le; i++) { 1780 x_diffs[i] = x_diffs[i + 1] - x_diffs[i]; 1781 y_diffs[i] = y_diffs[i + 1] - y_diffs[i]; 1782 x_table[level + 1][i] = x_diffs[i]; 1783 y_table[level + 1][i] = y_diffs[i]; 1784 x_recip[level + 1][i] = 1785 numerator / (x_recip[level][i + 1] - x_recip[level][i]) + 1786 x_recip[level - 1][i + 1]; 1787 y_recip[level + 1][i] = 1788 numerator / (y_recip[level][i + 1] - y_recip[level][i]) + 1789 y_recip[level - 1][i + 1]; 1790 } 1791 1792 // if (level == 1) { 1793 // console.log("bends level=", level, y_diffs.toString()); 1794 // } 1795 1796 // Store point location which may be centered around 1797 // critical points. 1798 // If the level is suitable, step out of the loop. 1799 groups = this._criticalPoints(y_diffs, le, level); 1800 if (groups === false) { 1801 // Its seems, the degree of the polynomial is equal to level 1802 console.log("Polynomial of degree", level); 1803 groups = []; 1804 break; 1805 } 1806 if (groups.length > 0) { 1807 foundCriticalPoint++; 1808 if (foundCriticalPoint > 1 && level % 2 === 0) { 1809 break; 1810 } 1811 } 1812 le--; 1813 } 1814 1815 // console.log("Last diffs", y_diffs, "level", level); 1816 1817 // Analyze the groups which have been found. 1818 for (i = 0; i < groups.length; i++) { 1819 // console.log("Group", i, groups[i]) 1820 // Identify the maximum difference, i.e. the center of the "problem" 1821 ma = -Infinity; 1822 for (j = 0; j < groups[i].length; j++) { 1823 v = Math.abs(groups[i][j].v); 1824 if (v > ma) { 1825 ma = v; 1826 pos = j; 1827 } 1828 } 1829 pos = Math.floor(groups[i][pos].i + level / 2); 1830 // Analyze the critical point 1831 criticalPoints.push( 1832 this.getPointType( 1833 curve, 1834 pos, 1835 t_values, 1836 x_values, 1837 y_values, 1838 x_slopes, 1839 y_slopes, 1840 le + 1 1841 ) 1842 ); 1843 } 1844 1845 return [criticalPoints, x_table, y_table, x_recip, y_recip]; 1846 }, 1847 1848 getCenterOfCriticalInterval: function (group, degree, t_values) { 1849 var ma, 1850 j, 1851 pos, 1852 v, 1853 num = 0.0, 1854 den = 0.0, 1855 h = t_values[1] - t_values[0], 1856 pos_mean, 1857 range = []; 1858 1859 // Identify the maximum difference, i.e. the center of the "problem" 1860 // If there are several equal maxima, store the positions 1861 // in the array range and determine the center of the array. 1862 1863 ma = -Infinity; 1864 range = []; 1865 for (j = 0; j < group.length; j++) { 1866 v = Math.abs(group[j].v); 1867 if (v > ma) { 1868 range = [j]; 1869 ma = v; 1870 pos = j; 1871 } else if (ma === v) { 1872 range.push(j); 1873 } 1874 } 1875 if (range.length > 0) { 1876 pos_mean = 1877 range.reduce(function (total, val) { 1878 return total + val; 1879 }, 0) / range.length; 1880 pos = Math.floor(pos_mean); 1881 pos_mean += group[0].i; 1882 } 1883 1884 if (ma < Infinity) { 1885 for (j = 0; j < group.length; j++) { 1886 num += Math.abs(group[j].v) * group[j].i; 1887 den += Math.abs(group[j].v); 1888 } 1889 pos_mean = num / den; 1890 } 1891 pos_mean += degree / 2; 1892 return [ 1893 group[pos].i + degree / 2, 1894 pos_mean, 1895 t_values[Math.floor(pos_mean)] + h * (pos_mean - Math.floor(pos_mean)) 1896 ]; 1897 }, 1898 1899 differenceMethod: function (component, curve) { 1900 var i, 1901 level, 1902 le, 1903 up, 1904 t_values = component.t_values, 1905 x_values = component.x_values, 1906 y_values = component.y_values, 1907 x_table = [], 1908 y_table = [], 1909 foundCriticalPoint = 0, 1910 degree_x = -1, 1911 degree_y = -1, 1912 pos, 1913 res, 1914 res_x, 1915 res_y, 1916 t_approx, 1917 groups = [], 1918 types, 1919 criticalPoints = []; 1920 1921 le = y_values.length; 1922 // x_table.push([]); 1923 // y_table.push([]); 1924 // for (i = 0; i < le; i++) { 1925 // x_table[0][i] = x_values[i]; 1926 // y_table[0][i] = y_values[i]; 1927 // } 1928 x_table.push(new Float64Array(x_values)); 1929 y_table.push(new Float64Array(y_values)); 1930 1931 le--; 1932 up = Math.min(12, le); 1933 for (level = 0; level < up; level++) { 1934 // Old style method: 1935 // x_table.push([]); 1936 // y_table.push([]); 1937 // for (i = 0; i < le; i++) { 1938 // x_table[level + 1][i] = x_table[level][i + 1] - x_table[level][i]; 1939 // y_table[level + 1][i] = y_table[level][i + 1] - y_table[level][i]; 1940 // } 1941 // New method: 1942 x_table.push(new Float64Array(le)); 1943 y_table.push(new Float64Array(le)); 1944 x_table[level + 1] = x_table[level].map(function (v, idx, arr) { 1945 return arr[idx + 1] - v; 1946 }); 1947 y_table[level + 1] = y_table[level].map(function (v, idx, arr) { 1948 return arr[idx + 1] - v; 1949 }); 1950 1951 // Store point location which may be centered around critical points. 1952 // If the level is suitable, step out of the loop. 1953 res_y = this._criticalInterval(y_table[level + 1], le, level); 1954 if (res_y.smooth === true) { 1955 // Its seems, the degree of the polynomial is equal to level 1956 // If the values in level + 1 are zero, it might be a polynomial of degree level. 1957 // Seems to work numerically stable until degree 6. 1958 degree_y = level; 1959 groups = []; 1960 } 1961 res_x = this._criticalInterval(x_table[level + 1], le, level); 1962 if (degree_x === -1 && res_x.smooth === true) { 1963 // Its seems, the degree of the polynomial is equal to level 1964 // If the values in level + 1 are zero, it might be a polynomial of degree level. 1965 // Seems to work numerically stable until degree 6. 1966 degree_x = level; 1967 } 1968 if (degree_y >= 0) { 1969 break; 1970 } 1971 1972 if (res_y.groups.length > 0) { 1973 foundCriticalPoint++; 1974 if (foundCriticalPoint > 2 && (level + 1) % 2 === 0) { 1975 groups = res_y.groups; 1976 types = res_y.types; 1977 break; 1978 } 1979 } 1980 le--; 1981 } 1982 1983 // console.log("Last diffs", y_table[Math.min(level + 1, up)], "level", level + 1); 1984 // Analyze the groups which have been found. 1985 for (i = 0; i < groups.length; i++) { 1986 if (types[i] === 'interval') { 1987 continue; 1988 } 1989 // console.log("Group", i, groups[i], types[i], level + 1) 1990 res = this.getCenterOfCriticalInterval(groups[i], level + 1, t_values); 1991 pos = res_y[0]; 1992 pos = Math.floor(res[1]); 1993 t_approx = res[2]; 1994 // console.log("Critical points:", groups, res, pos) 1995 1996 // Analyze the type of the critical point 1997 // Result is of type 'borderleft', borderright', 'other' 1998 criticalPoints.push( 1999 this.getPointType(curve, pos, t_approx, t_values, x_table, y_table, le + 1) 2000 ); 2001 } 2002 2003 // if (level === up) { 2004 // console.log("No convergence!"); 2005 // } else { 2006 // console.log("Convergence level", level); 2007 // } 2008 return [criticalPoints, x_table, y_table, degree_x, degree_y]; 2009 }, 2010 2011 _insertPoint_v4: function (curve, crds, t, doLog) { 2012 var p, 2013 prev = null, 2014 x, 2015 y, 2016 near = 0.8; 2017 2018 if (curve.points.length > 0) { 2019 prev = curve.points[curve.points.length - 1].scrCoords; 2020 } 2021 2022 // Add regular point 2023 p = new Coords(Const.COORDS_BY_USER, crds, curve.board); 2024 2025 if (prev !== null) { 2026 x = p.scrCoords[1] - prev[1]; 2027 y = p.scrCoords[2] - prev[2]; 2028 if (x * x + y * y < near * near) { 2029 // Math.abs(p.scrCoords[1] - prev[1]) < near && 2030 // Math.abs(p.scrCoords[2] - prev[2]) < near) { 2031 return; 2032 } 2033 } 2034 2035 p._t = t; 2036 curve.points.push(p); 2037 }, 2038 2039 getInterval: function (curve, ta, tb) { 2040 var t_int, 2041 // x_int, 2042 y_int; 2043 2044 //console.log('critical point', ta, tb); 2045 IntervalArithmetic.disable(); 2046 2047 t_int = IntervalArithmetic.Interval(ta, tb); 2048 curve.board.mathLib = IntervalArithmetic; 2049 curve.board.mathLibJXG = IntervalArithmetic; 2050 // x_int = curve.X(t_int, true); 2051 y_int = curve.Y(t_int, true); 2052 curve.board.mathLib = Math; 2053 curve.board.mathLibJXG = JXG.Math; 2054 2055 //console.log(x_int, y_int); 2056 return y_int; 2057 }, 2058 2059 sign: function (v) { 2060 if (v < 0) { 2061 return -1; 2062 } 2063 if (v > 0) { 2064 return 1; 2065 } 2066 return 0; 2067 }, 2068 2069 handleBorder: function (curve, comp, group, x_table, y_table) { 2070 var idx = group.idx, 2071 t, 2072 t1, 2073 t2, 2074 size = 32, 2075 y_int, 2076 x, 2077 y, 2078 lo, 2079 hi, 2080 i, 2081 components2, 2082 le, 2083 h; 2084 2085 // console.log("HandleBorder at t =", t_approx); 2086 // console.log("component:", comp) 2087 // console.log("Group:", group); 2088 2089 h = comp.t_values[1] - comp.t_values[0]; 2090 if (group.type === 'borderleft') { 2091 t = comp.left_isNaN ? comp.left_t : group.t - h; 2092 t1 = t; 2093 t2 = t1 + h; 2094 } else if (group.type === 'borderright') { 2095 t = comp.right_isNaN ? comp.right_t : group.t + h; 2096 t2 = t; 2097 t1 = t2 - h; 2098 } else { 2099 console.log("No bordercase!!!"); 2100 } 2101 2102 components2 = this.findComponents(curve, t1, t2, size); 2103 if (components2.length === 0) { 2104 return; 2105 } 2106 if (group.type === 'borderleft') { 2107 t1 = components2[0].left_t; 2108 t2 = components2[0].t_values[0]; 2109 h = components2[0].t_values[1] - components2[0].t_values[0]; 2110 t1 = t1 === null ? t2 - h : t1; 2111 t = t1; 2112 y_int = this.getInterval(curve, t1, t2); 2113 if (Type.isObject(y_int)) { 2114 lo = y_int.lo; 2115 hi = y_int.hi; 2116 2117 x = curve.X(t, true); 2118 y = y_table[1][idx] < 0 ? hi : lo; 2119 this._insertPoint_v4(curve, [1, x, y], t); 2120 } 2121 } 2122 2123 le = components2[0].t_values.length; 2124 for (i = 0; i < le; i++) { 2125 t = components2[0].t_values[i]; 2126 x = components2[0].x_values[i]; 2127 y = components2[0].y_values[i]; 2128 this._insertPoint_v4(curve, [1, x, y], t); 2129 } 2130 2131 if (group.type === 'borderright') { 2132 t1 = components2[0].t_values[le - 1]; 2133 t2 = components2[0].right_t; 2134 h = components2[0].t_values[1] - components2[0].t_values[0]; 2135 t2 = t2 === null ? t1 + h : t2; 2136 2137 t = t2; 2138 y_int = this.getInterval(curve, t1, t2); 2139 if (Type.isObject(y_int)) { 2140 lo = y_int.lo; 2141 hi = y_int.hi; 2142 x = curve.X(t, true); 2143 y = y_table[1][idx] > 0 ? hi : lo; 2144 this._insertPoint_v4(curve, [1, x, y], t); 2145 } 2146 } 2147 }, 2148 2149 _seconditeration_v4: function (curve, comp, group, x_table, y_table) { 2150 var i, t1, t2, ret, components2, comp2, idx, groups2, g, x_table2, y_table2, start, le; 2151 2152 // Look at two points, hopefully left and right from the critical point 2153 t1 = comp.t_values[group.idx - 2]; 2154 t2 = comp.t_values[group.idx + 2]; 2155 components2 = this.findComponents(curve, t1, t2, 64); 2156 for (idx = 0; idx < components2.length; idx++) { 2157 comp2 = components2[idx]; 2158 ret = this.differenceMethod(comp2, curve); 2159 groups2 = ret[0]; 2160 x_table2 = ret[1]; 2161 y_table2 = ret[2]; 2162 start = 0; 2163 for (g = 0; g <= groups2.length; g++) { 2164 if (g === groups2.length) { 2165 le = comp2.len; 2166 } else { 2167 le = groups2[g].idx; 2168 } 2169 2170 // Insert all uncritical points until next critical point 2171 for (i = start; i < le; i++) { 2172 if (!isNaN(comp2.x_values[i]) && !isNaN(comp2.y_values[i])) { 2173 this._insertPoint_v4( 2174 curve, 2175 [1, comp2.x_values[i], comp2.y_values[i]], 2176 comp2.t_values[i] 2177 ); 2178 } 2179 } 2180 // Handle next critical point 2181 if (g < groups2.length) { 2182 this.handleSingularity(curve, comp2, groups2[g], x_table2, y_table2); 2183 start = groups2[g].idx + 1; 2184 } 2185 } 2186 le = comp2.len; 2187 if (idx < components2.length - 1) { 2188 this._insertPoint_v4(curve, [1, NaN, NaN], comp2.right_t); 2189 } 2190 } 2191 return this; 2192 }, 2193 2194 _recurse_v4: function (curve, t1, t2, x1, y1, x2, y2, level) { 2195 var tol = 2, 2196 t = (t1 + t2) * 0.5, 2197 x = curve.X(t, true), 2198 y = curve.Y(t, true), 2199 dx, 2200 dy; 2201 2202 //console.log("Level", level) 2203 if (level === 0) { 2204 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2205 return; 2206 } 2207 // console.log("R", t1, t2) 2208 dx = (x - x1) * curve.board.unitX; 2209 dy = (y - y1) * curve.board.unitY; 2210 // console.log("D1", Math.sqrt(dx * dx + dy * dy)) 2211 if (Mat.hypot(dx, dy) > tol) { 2212 this._recurse_v4(curve, t1, t, x1, y1, x, y, level - 1); 2213 } else { 2214 this._insertPoint_v4(curve, [1, x, y], t); 2215 } 2216 dx = (x - x2) * curve.board.unitX; 2217 dy = (y - y2) * curve.board.unitY; 2218 // console.log("D2", Math.sqrt(dx * dx + dy * dy), x-x2, y-y2) 2219 if (Mat.hypot(dx, dy) > tol) { 2220 this._recurse_v4(curve, t, t2, x, y, x2, y2, level - 1); 2221 } else { 2222 this._insertPoint_v4(curve, [1, x, y], t); 2223 } 2224 }, 2225 2226 handleSingularity: function (curve, comp, group, x_table, y_table) { 2227 var idx = group.idx, 2228 t, 2229 t1, 2230 t2, 2231 y_int, 2232 i1, 2233 i2, 2234 x, 2235 // y, 2236 lo, 2237 hi, 2238 d_lft, 2239 d_rgt, 2240 d_thresh = 100, 2241 // d1, 2242 // d2, 2243 di1 = 5, 2244 di2 = 3; 2245 2246 t = group.t; 2247 console.log("HandleSingularity at t =", t); 2248 // console.log(comp.t_values[idx - 1], comp.y_values[idx - 1], comp.t_values[idx + 1], comp.y_values[idx + 1]); 2249 // console.log(group); 2250 2251 // Look at two points, hopefully left and right from the critical point 2252 t1 = comp.t_values[idx - di1]; 2253 t2 = comp.t_values[idx + di1]; 2254 2255 y_int = this.getInterval(curve, t1, t2); 2256 if (Type.isObject(y_int)) { 2257 lo = y_int.lo; 2258 hi = y_int.hi; 2259 } else { 2260 if (y_table[0][idx - 1] < y_table[0][idx + 1]) { 2261 lo = y_table[0][idx - 1]; 2262 hi = y_table[0][idx + 1]; 2263 } else { 2264 lo = y_table[0][idx + 1]; 2265 hi = y_table[0][idx - 1]; 2266 } 2267 } 2268 2269 x = curve.X(t, true); 2270 2271 d_lft = 2272 (y_table[0][idx - di2] - y_table[0][idx - di1]) / 2273 (comp.t_values[idx - di2] - comp.t_values[idx - di1]); 2274 d_rgt = 2275 (y_table[0][idx + di2] - y_table[0][idx + di1]) / 2276 (comp.t_values[idx + di2] - comp.t_values[idx + di1]); 2277 2278 console.log(":::", d_lft, d_rgt); 2279 2280 //this._insertPoint_v4(curve, [1, NaN, NaN], 0); 2281 2282 if (d_lft < -d_thresh) { 2283 // Left branch very steep downwards -> add the minimum 2284 this._insertPoint_v4(curve, [1, x, lo], t, true); 2285 if (d_rgt <= d_thresh) { 2286 // Right branch not very steep upwards -> interrupt the curve 2287 // I.e. it looks like -infty / (finite or infty) and not like -infty / -infty 2288 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2289 } 2290 } else if (d_lft > d_thresh) { 2291 // Left branch very steep upwards -> add the maximum 2292 this._insertPoint_v4(curve, [1, x, hi], t); 2293 if (d_rgt >= -d_thresh) { 2294 // Right branch not very steep downwards -> interrupt the curve 2295 // I.e. it looks like infty / (finite or -infty) and not like infty / infty 2296 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2297 } 2298 } else { 2299 if (lo === -Infinity) { 2300 this._insertPoint_v4(curve, [1, x, lo], t, true); 2301 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2302 } 2303 if (hi === Infinity) { 2304 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2305 this._insertPoint_v4(curve, [1, x, hi], t, true); 2306 } 2307 2308 if (group.t < comp.t_values[idx]) { 2309 i1 = idx - 1; 2310 i2 = idx; 2311 } else { 2312 i1 = idx; 2313 i2 = idx + 1; 2314 } 2315 t1 = comp.t_values[i1]; 2316 t2 = comp.t_values[i2]; 2317 this._recurse_v4( 2318 curve, 2319 t1, 2320 t2, 2321 x_table[0][i1], 2322 y_table[0][i1], 2323 x_table[0][i2], 2324 y_table[0][i2], 2325 10 2326 ); 2327 2328 // x = (x_table[0][idx] - x_table[0][idx - 1]) * curve.board.unitX; 2329 // y = (y_table[0][idx] - y_table[0][idx - 1]) * curve.board.unitY; 2330 // d1 = Math.sqrt(x * x + y * y); 2331 // x = (x_table[0][idx + 1] - x_table[0][idx]) * curve.board.unitX; 2332 // y = (y_table[0][idx + 1] - y_table[0][idx]) * curve.board.unitY; 2333 // d2 = Math.sqrt(x * x + y * y); 2334 2335 // console.log("end", t1, t2, t); 2336 // if (true || (d1 > 2 || d2 > 2)) { 2337 2338 // console.log(d1, d2, y_table[0][idx]) 2339 // // Finite jump 2340 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2341 // } else { 2342 // if (lo !== -Infinity && hi !== Infinity) { 2343 // // Critical point which can be ignored 2344 // this._insertPoint_v4(curve, [1, x_table[0][idx], y_table[0][idx]], comp.t_values[idx]); 2345 // } else { 2346 // if (lo === -Infinity) { 2347 // this._insertPoint_v4(curve, [1, x, lo], t, true); 2348 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2349 // } 2350 // if (hi === Infinity) { 2351 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2352 // this._insertPoint_v4(curve, [1, x, hi], t, true); 2353 // } 2354 // } 2355 // } 2356 } 2357 if (d_rgt < -d_thresh) { 2358 // Right branch very steep downwards -> add the maximum 2359 this._insertPoint_v4(curve, [1, x, hi], t); 2360 } else if (d_rgt > d_thresh) { 2361 // Right branch very steep upwards -> add the minimum 2362 this._insertPoint_v4(curve, [1, x, lo], t); 2363 } 2364 }, 2365 2366 /** 2367 * Number of equidistant points where the function is evaluated 2368 */ 2369 steps: 1021, //2053, // 1021, 2370 2371 /** 2372 * If the absolute maximum of the set of differences is larger than 2373 * criticalThreshold * median of these values, it is regarded as critical point. 2374 * @see JXG.Math.Plot._criticalInterval 2375 */ 2376 criticalThreshold: 1000, 2377 2378 plot_v4: function (curve, ta, tb, steps) { 2379 var i, 2380 // j, 2381 le, 2382 components, 2383 idx, 2384 comp, 2385 groups, 2386 g, 2387 start, 2388 ret, 2389 x_table, y_table, 2390 t, t1, t2, 2391 // good, 2392 // bad, 2393 // x_int, 2394 y_int, 2395 // degree_x, 2396 // degree_y, 2397 h = (tb - ta) / steps, 2398 Ypl = function (x) { 2399 return curve.Y(x, true); 2400 }, 2401 Ymi = function (x) { 2402 return -curve.Y(x, true); 2403 }, 2404 h2 = h * 0.5; 2405 2406 components = this.findComponents(curve, ta, tb, steps); 2407 for (idx = 0; idx < components.length; idx++) { 2408 comp = components[idx]; 2409 ret = this.differenceMethod(comp, curve); 2410 groups = ret[0]; 2411 x_table = ret[1]; 2412 y_table = ret[2]; 2413 2414 // degree_x = ret[3]; 2415 // degree_y = ret[4]; 2416 // if (degree_x >= 0) { 2417 // console.log("x polynomial of degree", degree_x); 2418 // } 2419 // if (degree_y >= 0) { 2420 // console.log("y polynomial of degree", degree_y); 2421 // } 2422 if (groups.length === 0 || groups[0].type !== 'borderleft') { 2423 groups.unshift({ 2424 idx: 0, 2425 t: comp.t_values[0], 2426 x: comp.x_values[0], 2427 y: comp.y_values[0], 2428 type: "borderleft" 2429 }); 2430 } 2431 if (groups[groups.length - 1].type !== 'borderright') { 2432 le = comp.t_values.length; 2433 groups.push({ 2434 idx: le - 1, 2435 t: comp.t_values[le - 1], 2436 x: comp.x_values[le - 1], 2437 y: comp.y_values[le - 1], 2438 type: "borderright" 2439 }); 2440 } 2441 2442 start = 0; 2443 for (g = 0; g <= groups.length; g++) { 2444 if (g === groups.length) { 2445 le = comp.len; 2446 } else { 2447 le = groups[g].idx - 1; 2448 } 2449 2450 // good = 0; 2451 // bad = 0; 2452 // Insert all uncritical points until next critical point 2453 for (i = start; i < le - 2; i++) { 2454 this._insertPoint_v4( 2455 curve, 2456 [1, comp.x_values[i], comp.y_values[i]], 2457 comp.t_values[i] 2458 ); 2459 // j = Math.max(0, i - 2); 2460 // Add more points in critical intervals 2461 if ( 2462 //degree_y === -1 && // No polynomial 2463 i >= start + 3 && 2464 i < le - 3 && // Do not do this if too close to a critical point 2465 y_table.length > 3 && 2466 Math.abs(y_table[2][i]) > 0.2 * Math.abs(y_table[0][i]) 2467 ) { 2468 t = comp.t_values[i]; 2469 h2 = h * 0.25; 2470 y_int = this.getInterval(curve, t, t + h); 2471 if (Type.isObject(y_int)) { 2472 if (y_table[2][i] > 0) { 2473 this._insertPoint_v4(curve, [1, t + h2, y_int.lo], t + h2); 2474 } else { 2475 this._insertPoint_v4( 2476 curve, 2477 [1, t + h - h2, y_int.hi], 2478 t + h - h2 2479 ); 2480 } 2481 } else { 2482 t1 = Numerics.fminbr(Ypl, [t, t + h]); 2483 t2 = Numerics.fminbr(Ymi, [t, t + h]); 2484 if (t1 < t2) { 2485 this._insertPoint_v4( 2486 curve, 2487 [1, curve.X(t1, true), curve.Y(t1, true)], 2488 t1 2489 ); 2490 this._insertPoint_v4( 2491 curve, 2492 [1, curve.X(t2, true), curve.Y(t2, true)], 2493 t2 2494 ); 2495 } else { 2496 this._insertPoint_v4( 2497 curve, 2498 [1, curve.X(t2, true), curve.Y(t2, true)], 2499 t2 2500 ); 2501 this._insertPoint_v4( 2502 curve, 2503 [1, curve.X(t1, true), curve.Y(t1, true)], 2504 t1 2505 ); 2506 } 2507 } 2508 // bad++; 2509 // } else { 2510 // good++; 2511 } 2512 } 2513 // console.log("GOOD", good, "BAD", bad); 2514 2515 // Handle next critical point 2516 if (g < groups.length) { 2517 //console.log("critical point / interval", groups[g]); 2518 2519 i = groups[g].idx; 2520 if (groups[g].type === "borderleft" || groups[g].type === 'borderright') { 2521 this.handleBorder(curve, comp, groups[g], x_table, y_table); 2522 } else { 2523 this._seconditeration_v4(curve, comp, groups[g], x_table, y_table); 2524 } 2525 2526 start = groups[g].idx + 1 + 1; 2527 } 2528 } 2529 2530 le = comp.len; 2531 if (idx < components.length - 1) { 2532 this._insertPoint_v4(curve, [1, NaN, NaN], comp.right_t); 2533 } 2534 } 2535 }, 2536 2537 /** 2538 * Updates the data points of a parametric curve, plotVersion 4. This version is used if {@link JXG.Curve#plotVersion} is <tt>4</tt>. 2539 * @param {JXG.Curve} curve JSXGraph curve element 2540 * @param {Number} mi Left bound of curve 2541 * @param {Number} ma Right bound of curve 2542 * @returns {JXG.Curve} Reference to the curve object. 2543 */ 2544 updateParametricCurve_v4: function (curve, mi, ma) { 2545 var ta, tb, w2, bbox; 2546 2547 if (curve.xterm === 'x') { 2548 // For function graphs we can restrict the plot interval 2549 // to the visible area +plus margin 2550 bbox = curve.board.getBoundingBox(); 2551 w2 = (bbox[2] - bbox[0]) * 0.3; 2552 // h2 = (bbox[1] - bbox[3]) * 0.3; 2553 ta = Math.max(mi, bbox[0] - w2); 2554 tb = Math.min(ma, bbox[2] + w2); 2555 } else { 2556 ta = mi; 2557 tb = ma; 2558 } 2559 2560 curve.points = []; 2561 2562 //console.log("--------------------"); 2563 this.plot_v4(curve, ta, tb, this.steps); 2564 2565 curve.numberPoints = curve.points.length; 2566 //console.log(curve.numberPoints); 2567 }, 2568 2569 //---------------------------------------------------------------------- 2570 // Plot algorithm alias 2571 //---------------------------------------------------------------------- 2572 2573 /** 2574 * Updates the data points of a parametric curve, alias for {@link JXG.Curve#updateParametricCurve_v2}. 2575 * This is needed for backwards compatibility, if this method has been 2576 * used directly in an application. 2577 * @param {JXG.Curve} curve JSXGraph curve element 2578 * @param {Number} mi Left bound of curve 2579 * @param {Number} ma Right bound of curve 2580 * @returns {JXG.Curve} Reference to the curve object. 2581 * 2582 * @see JXG.Curve#updateParametricCurve_v2 2583 */ 2584 updateParametricCurve: function (curve, mi, ma) { 2585 return this.updateParametricCurve_v2(curve, mi, ma); 2586 } 2587 }; 2588 2589 export default Mat.Plot; 2590