1 /* 2 Copyright 2008-2023 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"; 34 import Const from "../base/constants"; 35 import Coords from "../base/coords"; 36 import Mat from "./math"; 37 import Extrapolate from "./extrapolate"; 38 import Numerics from "./numerics"; 39 import Statistics from "./statistics"; 40 import Geometry from "./geometry"; 41 import IntervalArithmetic from "./ia"; 42 import Type from "../utils/type"; 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 if (curve.board.updateQuality === curve.board.BOARD_QUALITY_LOW) { 652 depth = Type.evaluate(curve.visProp.recursiondepthlow) || 13; 653 delta = 2; 654 // this.smoothLevel = 5; //depth - 7; 655 this.smoothLevel = depth - 6; 656 this.jumpLevel = 3; 657 } else { 658 depth = Type.evaluate(curve.visProp.recursiondepthhigh) || 17; 659 delta = 2; 660 // smoothLevel has to be small for graphs in a huge interval. 661 // this.smoothLevel = 3; //depth - 7; // 9 662 this.smoothLevel = depth - 9; // 9 663 this.jumpLevel = 2; 664 } 665 this.nanLevel = depth - 4; 666 667 curve.points = []; 668 669 if (this.xterm === "x") { 670 // For function graphs we can restrict the plot interval 671 // to the visible area + plus margin 672 bbox = curve.board.getBoundingBox(); 673 w2 = (bbox[2] - bbox[0]) * 0.3; 674 // h2 = (bbox[1] - bbox[3]) * 0.3; 675 ta = Math.max(mi, bbox[0] - w2); 676 tb = Math.min(ma, bbox[2] + w2); 677 } else { 678 ta = mi; 679 tb = ma; 680 } 681 pa.setCoordinates( 682 Const.COORDS_BY_USER, 683 [curve.X(ta, suspendUpdate), curve.Y(ta, suspendUpdate)], 684 false 685 ); 686 687 // The first function calls of X() and Y() are done. We can now 688 // switch `suspendUpdate` on. If supported by the functions, this 689 // avoids for the rest of the plotting algorithm, evaluation of any 690 // parent elements. 691 suspendUpdate = true; 692 693 pb.setCoordinates( 694 Const.COORDS_BY_USER, 695 [curve.X(tb, suspendUpdate), curve.Y(tb, suspendUpdate)], 696 false 697 ); 698 699 // Find start and end points of the visible area (plus a certain margin) 700 ret_arr = this._findStartPoint(curve, pa.scrCoords, ta, pb.scrCoords, tb); 701 pa.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 702 ta = ret_arr[1]; 703 ret_arr = this._findStartPoint(curve, pb.scrCoords, tb, pa.scrCoords, ta); 704 pb.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 705 tb = ret_arr[1]; 706 707 // Save the visible area. 708 // This can be used in Curve.hasPoint(). 709 this._visibleArea = [ta, tb]; 710 711 // Start recursive plotting algorithm 712 a = pa.copy("scrCoords"); 713 b = pb.copy("scrCoords"); 714 pa._t = ta; 715 curve.points.push(pa); 716 this._lastCrds = pa.copy("scrCoords"); // Used in _insertPoint 717 this._plotRecursive_v2(curve, a, ta, b, tb, depth, delta); 718 pb._t = tb; 719 curve.points.push(pb); 720 721 curve.numberPoints = curve.points.length; 722 //console.timeEnd("plot"); 723 724 return curve; 725 }, 726 727 //---------------------------------------------------------------------- 728 // Plot algorithm v3 729 //---------------------------------------------------------------------- 730 /** 731 * 732 * @param {JXG.Curve} curve JSXGraph curve element 733 * @param {*} pnt 734 * @param {*} t 735 * @param {*} depth 736 * @param {*} limes 737 * @private 738 */ 739 _insertLimesPoint: function (curve, pnt, t, depth, limes) { 740 var p0, p1, p2; 741 742 // Ignore jump point if it follows limes 743 if ( 744 (Math.abs(this._lastUsrCrds[1]) === Infinity && 745 Math.abs(limes.left_x) === Infinity) || 746 (Math.abs(this._lastUsrCrds[2]) === Infinity && Math.abs(limes.left_y) === Infinity) 747 ) { 748 // console.log("SKIP:", pnt.usrCoords, this._lastUsrCrds, limes); 749 return; 750 } 751 752 // // Ignore jump left from limes 753 // if (Math.abs(limes.left_x) > 100 * Math.abs(this._lastUsrCrds[1])) { 754 // x = Math.sign(limes.left_x) * Infinity; 755 // } else { 756 // x = limes.left_x; 757 // } 758 // if (Math.abs(limes.left_y) > 100 * Math.abs(this._lastUsrCrds[2])) { 759 // y = Math.sign(limes.left_y) * Infinity; 760 // } else { 761 // y = limes.left_y; 762 // } 763 // //pnt.setCoordinates(Const.COORDS_BY_USER, [x, y], false); 764 765 // Add points at a jump. pnt contains [NaN, NaN] 766 //console.log("Add", t, pnt.usrCoords, limes, depth) 767 p0 = new Coords(Const.COORDS_BY_USER, [limes.left_x, limes.left_y], curve.board); 768 p0._t = t; 769 curve.points.push(p0); 770 771 if ( 772 !isNaN(limes.left_x) && 773 !isNaN(limes.left_y) && 774 !isNaN(limes.right_x) && 775 !isNaN(limes.right_y) && 776 (Math.abs(limes.left_x - limes.right_x) > Mat.eps || 777 Math.abs(limes.left_y - limes.right_y) > Mat.eps) 778 ) { 779 p1 = new Coords(Const.COORDS_BY_SCREEN, pnt, curve.board); 780 p1._t = t; 781 curve.points.push(p1); 782 } 783 784 p2 = new Coords(Const.COORDS_BY_USER, [limes.right_x, limes.right_y], curve.board); 785 p2._t = t; 786 curve.points.push(p2); 787 this._lastScrCrds = p2.copy("scrCoords"); 788 this._lastUsrCrds = p2.copy("usrCoords"); 789 }, 790 791 /** 792 * Add a point to the curve plot. If the new point is too close to the previously inserted point, 793 * it is skipped. 794 * Used in {@link JXG.Curve._plotRecursive}. 795 * 796 * @private 797 * @param {JXG.Curve} curve JSXGraph curve element 798 * @param {JXG.Coords} pnt Coords to add to the list of points 799 */ 800 _insertPoint: function (curve, pnt, t, depth, limes) { 801 var last_is_real = !isNaN(this._lastScrCrds[1] + this._lastScrCrds[2]), // The last point was real 802 point_is_real = !isNaN(pnt[1] + pnt[2]), // New point is real point 803 cw = curve.board.canvasWidth, 804 ch = curve.board.canvasHeight, 805 p, 806 near = 0.8, 807 off = 500; 808 809 if (Type.exists(limes)) { 810 this._insertLimesPoint(curve, pnt, t, depth, limes); 811 return; 812 } 813 814 // Check if point has real coordinates and 815 // coordinates are not too far away from canvas. 816 point_is_real = 817 point_is_real && 818 pnt[1] > -off && 819 pnt[2] > -off && 820 pnt[1] < cw + off && 821 pnt[2] < ch + off; 822 823 // Prevent two consecutive NaNs 824 if (!last_is_real && !point_is_real) { 825 return; 826 } 827 828 // Prevent two consecutive points which are too close 829 if ( 830 point_is_real && 831 last_is_real && 832 Math.abs(pnt[1] - this._lastScrCrds[1]) < near && 833 Math.abs(pnt[2] - this._lastScrCrds[2]) < near 834 ) { 835 return; 836 } 837 838 // Prevent two consecutive points at infinity (either direction) 839 if ( 840 (Math.abs(pnt[1]) === Infinity && Math.abs(this._lastUsrCrds[1]) === Infinity) || 841 (Math.abs(pnt[2]) === Infinity && Math.abs(this._lastUsrCrds[2]) === Infinity) 842 ) { 843 return; 844 } 845 846 //console.log("add", t, pnt.usrCoords, depth) 847 // Add regular point 848 p = new Coords(Const.COORDS_BY_SCREEN, pnt, curve.board); 849 p._t = t; 850 curve.points.push(p); 851 this._lastScrCrds = p.copy("scrCoords"); 852 this._lastUsrCrds = p.copy("usrCoords"); 853 }, 854 855 /** 856 * Compute distances in screen coordinates between the points ab, 857 * ac, cb, and cd, where d = (a + b)/2. 858 * cd is used for the smoothness test, ab, ac, cb are used to detect jumps, cusps and poles. 859 * 860 * @private 861 * @param {Array} a Screen coordinates of the left interval bound 862 * @param {Array} b Screen coordinates of the right interval bound 863 * @param {Array} c Screen coordinates of the bisection point at (ta + tb) / 2 864 * @returns {Array} array of distances in screen coordinates between: ab, ac, cb, and cd. 865 */ 866 _triangleDists: function (a, b, c) { 867 var d, d_ab, d_ac, d_cb, d_cd; 868 869 d = [a[0] * b[0], (a[1] + b[1]) * 0.5, (a[2] + b[2]) * 0.5]; 870 871 d_ab = Geometry.distance(a, b, 3); 872 d_ac = Geometry.distance(a, c, 3); 873 d_cb = Geometry.distance(c, b, 3); 874 d_cd = Geometry.distance(c, d, 3); 875 876 return [d_ab, d_ac, d_cb, d_cd]; 877 }, 878 879 /** 880 * Test if the function is undefined on an interval: 881 * If the interval borders a and b are undefined, 20 random values 882 * are tested if they are undefined, too. 883 * Only if all values are undefined, we declare the function to be undefined in this interval. 884 * 885 * @private 886 * @param {JXG.Curve} curve JSXGraph curve element 887 * @param {Array} a Screen coordinates of the left interval bound 888 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 889 * @param {Array} b Screen coordinates of the right interval bound 890 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 891 */ 892 _isUndefined: function (curve, a, ta, b, tb) { 893 var t, i, pnt; 894 895 if (!isNaN(a[1] + a[2]) || !isNaN(b[1] + b[2])) { 896 return false; 897 } 898 899 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 900 901 for (i = 0; i < 20; ++i) { 902 t = ta + Math.random() * (tb - ta); 903 pnt.setCoordinates( 904 Const.COORDS_BY_USER, 905 [curve.X(t, true), curve.Y(t, true)], 906 false 907 ); 908 if (!isNaN(pnt.scrCoords[0] + pnt.scrCoords[1] + pnt.scrCoords[2])) { 909 return false; 910 } 911 } 912 913 return true; 914 }, 915 916 /** 917 * Decide if a path segment is too far from the canvas that we do not need to draw it. 918 * @private 919 * @param {Array} a Screen coordinates of the start point of the segment 920 * @param {Array} ta Curve parameter of a (unused). 921 * @param {Array} b Screen coordinates of the end point of the segment 922 * @param {Array} tb Curve parameter of b (unused). 923 * @param {JXG.Board} board 924 * @returns {Boolean} True if the segment is too far away from the canvas, false otherwise. 925 */ 926 _isOutside: function (a, ta, b, tb, board) { 927 var off = 500, 928 cw = board.canvasWidth, 929 ch = board.canvasHeight; 930 931 return !!( 932 (a[1] < -off && b[1] < -off) || 933 (a[2] < -off && b[2] < -off) || 934 (a[1] > cw + off && b[1] > cw + off) || 935 (a[2] > ch + off && b[2] > ch + off) 936 ); 937 }, 938 939 /** 940 * Decide if a point of a curve is too far from the canvas that we do not need to draw it. 941 * @private 942 * @param {Array} a Screen coordinates of the point 943 * @param {JXG.Board} board 944 * @returns {Boolean} True if the point is too far away from the canvas, false otherwise. 945 */ 946 _isOutsidePoint: function (a, board) { 947 var off = 500, 948 cw = board.canvasWidth, 949 ch = board.canvasHeight; 950 951 return !!(a[1] < -off || a[2] < -off || a[1] > cw + off || a[2] > ch + off); 952 }, 953 954 /** 955 * For a curve c(t) defined on the interval [ta, tb] find the first point 956 * which is in the visible area of the board (plus some outside margin). 957 * <p> 958 * This method is necessary to restrict the recursive plotting algorithm 959 * {@link JXG.Curve._plotRecursive} to the visible area and not waste 960 * recursion to areas far outside of the visible area. 961 * <p> 962 * This method can also be used to find the last visible point 963 * by reversing the input parameters. 964 * 965 * @param {JXG.Curve} curve JSXGraph curve element 966 * @param {Array} ta Curve parameter of a. 967 * @param {Array} b Screen coordinates of the end point of the segment (unused) 968 * @param {Array} tb Curve parameter of b 969 * @return {Array} Array of length two containing the screen ccordinates of 970 * the starting point and the curve parameter at this point. 971 * @private 972 */ 973 _findStartPoint: function (curve, a, ta, b, tb) { 974 // The code below is too unstable. 975 // E.g. [function(t) { return Math.pow(t, 2) * (t + 5) * Math.pow(t - 5, 2); }, -8, 8] 976 // Therefore, we return here. 977 return [a, ta]; 978 979 // var i, 980 // delta, 981 // tc, 982 // td, 983 // z, 984 // isFound, 985 // w2, 986 // h2, 987 // pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 988 // steps = 40, 989 // eps = 0.01, 990 // fnX1, 991 // fnX2, 992 // fnY1, 993 // fnY2, 994 // bbox = curve.board.getBoundingBox(); 995 996 // if (true || !this._isOutsidePoint(a, curve.board)) { 997 // return [a, ta]; 998 // } 999 // w2 = (bbox[2] - bbox[0]) * 0.3; 1000 // h2 = (bbox[1] - bbox[3]) * 0.3; 1001 // bbox[0] -= w2; 1002 // bbox[1] += h2; 1003 // bbox[2] += w2; 1004 // bbox[3] -= h2; 1005 1006 // delta = (tb - ta) / steps; 1007 // tc = ta + delta; 1008 // isFound = false; 1009 1010 // fnX1 = function (t) { 1011 // return curve.X(t, true) - bbox[0]; 1012 // }; 1013 // fnY1 = function (t) { 1014 // return curve.Y(t, true) - bbox[1]; 1015 // }; 1016 // fnX2 = function (t) { 1017 // return curve.X(t, true) - bbox[2]; 1018 // }; 1019 // fnY2 = function (t) { 1020 // return curve.Y(t, true) - bbox[3]; 1021 // }; 1022 // for (i = 0; i < steps; ++i) { 1023 // // Left border 1024 // z = bbox[0]; 1025 // td = Numerics.root(fnX1, [tc - delta, tc], curve); 1026 // // td = Numerics.fzero(fnX1, [tc - delta, tc], this); 1027 // // console.log("A", tc - delta, tc, td, Math.abs(this.X(td, true) - z)); 1028 // if (Math.abs(curve.X(td, true) - z) < eps) { 1029 // //} * Math.abs(z)) { 1030 // isFound = true; 1031 // break; 1032 // } 1033 // // Top border 1034 // z = bbox[1]; 1035 // td = Numerics.root(fnY1, [tc - delta, tc], curve); 1036 // // td = Numerics.fzero(fnY1, [tc - delta, tc], this); 1037 // // console.log("B", tc - delta, tc, td, Math.abs(this.Y(td, true) - z)); 1038 // if (Math.abs(curve.Y(td, true) - z) < eps) { 1039 // // * Math.abs(z)) { 1040 // isFound = true; 1041 // break; 1042 // } 1043 // // Right border 1044 // z = bbox[2]; 1045 // td = Numerics.root(fnX2, [tc - delta, tc], curve); 1046 // // td = Numerics.fzero(fnX2, [tc - delta, tc], this); 1047 // // console.log("C", tc - delta, tc, td, Math.abs(this.X(td, true) - z)); 1048 // if (Math.abs(curve.X(td, true) - z) < eps) { 1049 // // * Math.abs(z)) { 1050 // isFound = true; 1051 // break; 1052 // } 1053 // // Bottom border 1054 // z = bbox[3]; 1055 // td = Numerics.root(fnY2, [tc - delta, tc], curve); 1056 // // td = Numerics.fzero(fnY2, [tc - delta, tc], this); 1057 // // console.log("D", tc - delta, tc, td, Math.abs(this.Y(td, true) - z)); 1058 // if (Math.abs(curve.Y(td, true) - z) < eps) { 1059 // // * Math.abs(z)) { 1060 // isFound = true; 1061 // break; 1062 // } 1063 // tc += delta; 1064 // } 1065 // if (isFound) { 1066 // pnt.setCoordinates( 1067 // Const.COORDS_BY_USER, 1068 // [curve.X(td, true), curve.Y(td, true)], 1069 // false 1070 // ); 1071 // return [pnt.scrCoords, td]; 1072 // } 1073 // console.log("TODO _findStartPoint", curve.Y.toString(), tc); 1074 // pnt.setCoordinates(Const.COORDS_BY_USER, [curve.X(ta, true), curve.Y(ta, true)], false); 1075 // return [pnt.scrCoords, ta]; 1076 }, 1077 1078 /** 1079 * Investigate a function term at the bounds of intervals where 1080 * the function is not defined, e.g. log(x) at x = 0. 1081 * 1082 * c is inbetween a and b 1083 * 1084 * @param {JXG.Curve} curve JSXGraph curve element 1085 * @param {Array} a Screen coordinates of the left interval bound 1086 * @param {Array} b Screen coordinates of the right interval bound 1087 * @param {Array} c Screen coordinates of the bisection point at (ta + tb) / 2 1088 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 1089 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 1090 * @param {Number} tc (ta + tb) / 2 = tc. Parameter which evaluates to b, i.e. [1, X(tc), Y(tc)] = c in screen coordinates 1091 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 1092 * @returns {JXG.Boolean} true if the point is inserted and the recursion should stop, false otherwise. 1093 * 1094 * @private 1095 */ 1096 _getBorderPos: function (curve, ta, a, tc, c, tb, b) { 1097 var t, pnt, p, j, 1098 max_it = 30, 1099 is_undef = false, 1100 t_good, t_bad; 1101 1102 pnt = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false); 1103 j = 0; 1104 // Bisect a, b and c until the point t_real is inside of the definition interval 1105 // and as close as possible at the boundary. 1106 // (t_real2 is/was the second closest point). 1107 // There are four cases: 1108 // a | c | b 1109 // --------------- 1110 // inf | R | R 1111 // R | R | inf 1112 // inf | inf | R 1113 // R | inf | inf 1114 // 1115 if (isNaN(a[1] + a[2]) && !isNaN(c[1] + c[2])) { 1116 t_bad = ta; 1117 t_good = tc; 1118 } else if (isNaN(b[1] + b[2]) && !isNaN(c[1] + c[2])) { 1119 t_bad = tb; 1120 t_good = tc; 1121 } else if (isNaN(c[1] + c[2]) && !isNaN(b[1] + b[2])) { 1122 t_bad = tc; 1123 t_good = tb; 1124 } else if (isNaN(c[1] + c[2]) && !isNaN(a[1] + a[2])) { 1125 t_bad = tc; 1126 t_good = ta; 1127 } else { 1128 return false; 1129 } 1130 do { 1131 t = 0.5 * (t_good + t_bad); 1132 pnt.setCoordinates( 1133 Const.COORDS_BY_USER, 1134 [curve.X(t, true), curve.Y(t, true)], 1135 false 1136 ); 1137 p = pnt.usrCoords; 1138 is_undef = isNaN(p[1] + p[2]); 1139 if (is_undef) { 1140 t_bad = t; 1141 } else { 1142 t_good = t; 1143 } 1144 ++j; 1145 } while (j < max_it && Math.abs(t_good - t_bad) > Mat.eps); 1146 return t; 1147 }, 1148 1149 /** 1150 * 1151 * @param {JXG.Curve} curve JSXGraph curve element 1152 * @param {Number} ta 1153 * @param {Number} tb 1154 */ 1155 _getCuspPos: function (curve, ta, tb) { 1156 var a = [curve.X(ta, true), curve.Y(ta, true)], 1157 b = [curve.X(tb, true), curve.Y(tb, true)], 1158 max_func = function (t) { 1159 var c = [curve.X(t, true), curve.Y(t, true)]; 1160 return -( 1161 Mat.hypot(a[0] - c[0], a[1] - c[1]) + 1162 Mat.hypot(b[0] - c[0], b[1] - c[1]) 1163 ); 1164 }; 1165 1166 return Numerics.fminbr(max_func, [ta, tb], curve); 1167 }, 1168 1169 /** 1170 * 1171 * @param {JXG.Curve} curve JSXGraph curve element 1172 * @param {Number} ta 1173 * @param {Number} tb 1174 */ 1175 _getJumpPos: function (curve, ta, tb) { 1176 var max_func = function (t) { 1177 var e = Mat.eps * Mat.eps, 1178 c1 = [curve.X(t, true), curve.Y(t, true)], 1179 c2 = [curve.X(t + e, true), curve.Y(t + e, true)]; 1180 return -Math.abs((c2[1] - c1[1]) / (c2[0] - c1[0])); 1181 }; 1182 1183 return Numerics.fminbr(max_func, [ta, tb], curve); 1184 }, 1185 1186 /** 1187 * 1188 * @param {JXG.Curve} curve JSXGraph curve element 1189 * @param {Number} t 1190 * @private 1191 */ 1192 _getLimits: function (curve, t) { 1193 var res, 1194 step = 2 / (curve.maxX() - curve.minX()), 1195 x_l, 1196 x_r, 1197 y_l, 1198 y_r; 1199 1200 // From left 1201 res = Extrapolate.limit(t, -step, curve.X); 1202 x_l = res[0]; 1203 if (res[1] === "infinite") { 1204 x_l = Math.sign(x_l) * Infinity; 1205 } 1206 1207 res = Extrapolate.limit(t, -step, curve.Y); 1208 y_l = res[0]; 1209 if (res[1] === "infinite") { 1210 y_l = Math.sign(y_l) * Infinity; 1211 } 1212 1213 // From right 1214 res = Extrapolate.limit(t, step, curve.X); 1215 x_r = res[0]; 1216 if (res[1] === "infinite") { 1217 x_r = Math.sign(x_r) * Infinity; 1218 } 1219 1220 res = Extrapolate.limit(t, step, curve.Y); 1221 y_r = res[0]; 1222 if (res[1] === "infinite") { 1223 y_r = Math.sign(y_r) * Infinity; 1224 } 1225 1226 return { 1227 left_x: x_l, 1228 left_y: y_l, 1229 right_x: x_r, 1230 right_y: y_r, 1231 t: t 1232 }; 1233 }, 1234 1235 /** 1236 * 1237 * @param {JXG.Curve} curve JSXGraph curve element 1238 * @param {Array} a 1239 * @param {Number} tc 1240 * @param {Array} c 1241 * @param {Number} tb 1242 * @param {Array} b 1243 * @param {String} may_be_special 1244 * @param {Number} depth 1245 * @private 1246 */ 1247 _getLimes: function (curve, ta, a, tc, c, tb, b, may_be_special, depth) { 1248 var t; 1249 1250 if (may_be_special === "border") { 1251 t = this._getBorderPos(curve, ta, a, tc, c, tb, b); 1252 } else if (may_be_special === "cusp") { 1253 t = this._getCuspPos(curve, ta, tb); 1254 } else if (may_be_special === "jump") { 1255 t = this._getJumpPos(curve, ta, tb); 1256 } 1257 return this._getLimits(curve, t); 1258 }, 1259 1260 /** 1261 * Recursive interval bisection algorithm for curve plotting. 1262 * Used in {@link JXG.Curve.updateParametricCurve}. 1263 * @private 1264 * @param {JXG.Curve} curve JSXGraph curve element 1265 * @param {Array} a Screen coordinates of the left interval bound 1266 * @param {Number} ta Parameter which evaluates to a, i.e. [1, X(ta), Y(ta)] = a in screen coordinates 1267 * @param {Array} b Screen coordinates of the right interval bound 1268 * @param {Number} tb Parameter which evaluates to b, i.e. [1, X(tb), Y(tb)] = b in screen coordinates 1269 * @param {Number} depth Actual recursion depth. The recursion stops if depth is equal to 0. 1270 * @param {Number} delta If the distance of the bisection point at (ta + tb) / 2 from the point (a + b) / 2 is less then delta, 1271 * the segment [a,b] is regarded as straight line. 1272 * @returns {JXG.Curve} Reference to the curve object. 1273 */ 1274 _plotNonRecursive: function (curve, a, ta, b, tb, d) { 1275 var tc, 1276 c, 1277 ds, 1278 mindepth = 0, 1279 limes = null, 1280 a_nan, 1281 b_nan, 1282 isSmooth = false, 1283 may_be_special = "", 1284 x, 1285 y, 1286 oc, 1287 depth, 1288 ds0, 1289 stack = [], 1290 stack_length = 0, 1291 item; 1292 1293 oc = curve.board.origin.scrCoords; 1294 stack[stack_length++] = [a, ta, b, tb, d, Infinity]; 1295 while (stack_length > 0) { 1296 // item = stack.pop(); 1297 item = stack[--stack_length]; 1298 a = item[0]; 1299 ta = item[1]; 1300 b = item[2]; 1301 tb = item[3]; 1302 depth = item[4]; 1303 ds0 = item[5]; 1304 1305 isSmooth = false; 1306 may_be_special = ""; 1307 limes = null; 1308 //console.log(stack.length, item) 1309 1310 if (curve.points.length > 65536) { 1311 return; 1312 } 1313 1314 if (depth < this.nanLevel) { 1315 // Test if the function is undefined in the whole interval [ta, tb] 1316 if (this._isUndefined(curve, a, ta, b, tb)) { 1317 continue; 1318 } 1319 // Test if the graph is far outside the visible are for the interval [ta, tb] 1320 if (this._isOutside(a, ta, b, tb, curve.board)) { 1321 continue; 1322 } 1323 } 1324 1325 tc = (ta + tb) * 0.5; 1326 1327 // Screen coordinates of point at tc 1328 x = curve.X(tc, true); 1329 y = curve.Y(tc, true); 1330 c = [1, oc[1] + x * curve.board.unitX, oc[2] - y * curve.board.unitY]; 1331 ds = this._triangleDists(a, b, c); // returns [d_ab, d_ac, d_cb, d_cd] 1332 1333 a_nan = isNaN(a[1] + a[2]); 1334 b_nan = isNaN(b[1] + b[2]); 1335 if ((a_nan && !b_nan) || (!a_nan && b_nan)) { 1336 may_be_special = "border"; 1337 } else if ( 1338 ds[0] > 0.66 * ds0 || 1339 ds[0] < this.cusp_threshold * (ds[1] + ds[2]) || 1340 ds[1] > 5 * ds[2] || 1341 ds[2] > 5 * ds[1] 1342 ) { 1343 may_be_special = "cusp"; 1344 } else if ( 1345 ds[2] > this.jump_threshold * ds[0] || 1346 ds[1] > this.jump_threshold * ds[0] || 1347 ds[0] === Infinity || 1348 ds[1] === Infinity || 1349 ds[2] === Infinity 1350 ) { 1351 may_be_special = "jump"; 1352 } 1353 isSmooth = 1354 may_be_special === "" && 1355 depth < this.smoothLevel && 1356 ds[3] < this.smooth_threshold; 1357 1358 if (depth < this.testLevel && !isSmooth) { 1359 if (may_be_special === "") { 1360 isSmooth = true; 1361 } else { 1362 limes = this._getLimes(curve, ta, a, tc, c, tb, b, may_be_special, depth); 1363 } 1364 } 1365 1366 if (limes !== null) { 1367 c = [1, NaN, NaN]; 1368 this._insertPoint(curve, c, tc, depth, limes); 1369 } else if (depth <= mindepth || isSmooth) { 1370 this._insertPoint(curve, c, tc, depth, null); 1371 } else { 1372 stack[stack_length++] = [c, tc, b, tb, depth - 1, ds[0]]; 1373 stack[stack_length++] = [a, ta, c, tc, depth - 1, ds[0]]; 1374 } 1375 } 1376 1377 return this; 1378 }, 1379 1380 /** 1381 * Updates the data points of a parametric curve. This version is used if {@link JXG.Curve#plotVersion} is <tt>3</tt>. 1382 * This is an experimental plot version, <b>not recommended</b> to be used. 1383 * @param {JXG.Curve} curve JSXGraph curve element 1384 * @param {Number} mi Left bound of curve 1385 * @param {Number} ma Right bound of curve 1386 * @returns {JXG.Curve} Reference to the curve object. 1387 */ 1388 updateParametricCurve_v3: function (curve, mi, ma) { 1389 var ta, 1390 tb, 1391 a, 1392 b, 1393 suspendUpdate = false, 1394 pa = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 1395 pb = new Coords(Const.COORDS_BY_USER, [0, 0], curve.board, false), 1396 depth, 1397 w2, // h2, 1398 bbox, 1399 ret_arr; 1400 1401 // console.log("-----------------------------------------------------------"); 1402 // console.time("plot"); 1403 if (curve.board.updateQuality === curve.board.BOARD_QUALITY_LOW) { 1404 depth = Type.evaluate(curve.visProp.recursiondepthlow) || 14; 1405 } else { 1406 depth = Type.evaluate(curve.visProp.recursiondepthhigh) || 17; 1407 } 1408 1409 // smoothLevel has to be small for graphs in a huge interval. 1410 this.smoothLevel = 7; //depth - 10; 1411 this.nanLevel = depth - 4; 1412 this.testLevel = 4; 1413 this.cusp_threshold = 0.5; 1414 this.jump_threshold = 0.99; 1415 this.smooth_threshold = 2; 1416 1417 curve.points = []; 1418 1419 if (curve.xterm === "x") { 1420 // For function graphs we can restrict the plot interval 1421 // to the visible area +plus margin 1422 bbox = curve.board.getBoundingBox(); 1423 w2 = (bbox[2] - bbox[0]) * 0.3; 1424 //h2 = (bbox[1] - bbox[3]) * 0.3; 1425 ta = Math.max(mi, bbox[0] - w2); 1426 tb = Math.min(ma, bbox[2] + w2); 1427 } else { 1428 ta = mi; 1429 tb = ma; 1430 } 1431 pa.setCoordinates( 1432 Const.COORDS_BY_USER, 1433 [curve.X(ta, suspendUpdate), curve.Y(ta, suspendUpdate)], 1434 false 1435 ); 1436 1437 // The first function calls of X() and Y() are done. We can now 1438 // switch `suspendUpdate` on. If supported by the functions, this 1439 // avoids for the rest of the plotting algorithm, evaluation of any 1440 // parent elements. 1441 suspendUpdate = true; 1442 1443 pb.setCoordinates( 1444 Const.COORDS_BY_USER, 1445 [curve.X(tb, suspendUpdate), curve.Y(tb, suspendUpdate)], 1446 false 1447 ); 1448 1449 // Find start and end points of the visible area (plus a certain margin) 1450 ret_arr = this._findStartPoint(curve, pa.scrCoords, ta, pb.scrCoords, tb); 1451 pa.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 1452 ta = ret_arr[1]; 1453 ret_arr = this._findStartPoint(curve, pb.scrCoords, tb, pa.scrCoords, ta); 1454 pb.setCoordinates(Const.COORDS_BY_SCREEN, ret_arr[0], false); 1455 tb = ret_arr[1]; 1456 1457 // Save the visible area. 1458 // This can be used in Curve.hasPoint(). 1459 this._visibleArea = [ta, tb]; 1460 1461 // Start recursive plotting algorithm 1462 a = pa.copy("scrCoords"); 1463 b = pb.copy("scrCoords"); 1464 pa._t = ta; 1465 curve.points.push(pa); 1466 this._lastScrCrds = pa.copy("scrCoords"); // Used in _insertPoint 1467 this._lastUsrCrds = pa.copy("usrCoords"); // Used in _insertPoint 1468 1469 this._plotNonRecursive(curve, a, ta, b, tb, depth); 1470 1471 pb._t = tb; 1472 curve.points.push(pb); 1473 1474 curve.numberPoints = curve.points.length; 1475 // console.timeEnd("plot"); 1476 // console.log("number of points:", this.numberPoints); 1477 1478 return curve; 1479 }, 1480 1481 //---------------------------------------------------------------------- 1482 // Plot algorithm v4 1483 //---------------------------------------------------------------------- 1484 1485 _criticalInterval: function (vec, le, level) { 1486 var i, 1487 j, 1488 le1, 1489 med, 1490 sgn, 1491 sgnChange, 1492 isGroup = false, 1493 abs_vec, 1494 last = -Infinity, 1495 very_small = false, 1496 smooth = false, 1497 group = 0, 1498 groups = [], 1499 types = [], 1500 positions = []; 1501 1502 abs_vec = Statistics.abs(vec); 1503 med = Statistics.median(abs_vec); 1504 1505 if (med < 1.0e-7) { 1506 med = 1.0e-7; 1507 very_small = true; 1508 } else { 1509 med *= this.criticalThreshold; 1510 } 1511 1512 //console.log("Median", med); 1513 for (i = 0; i < le; i++) { 1514 // Start a group if not yet done and 1515 // add position to group 1516 if (abs_vec[i] > med /*&& abs_vec[i] > 0.01*/) { 1517 positions.push({ i: i, v: vec[i], group: group }); 1518 last = i; 1519 if (!isGroup) { 1520 isGroup = true; 1521 } 1522 } else { 1523 if (isGroup && i > last + 4) { 1524 // End the group 1525 if (positions.length > 0) { 1526 groups.push(positions.slice(0)); 1527 } 1528 positions = []; 1529 isGroup = false; 1530 group++; 1531 } 1532 } 1533 } 1534 if (isGroup) { 1535 if (positions.length > 1) { 1536 groups.push(positions.slice(0)); 1537 } 1538 } 1539 1540 if (very_small && groups.length === 0) { 1541 smooth = true; 1542 } 1543 1544 // Decide if there is a singular critical point 1545 // or if a whole interval is problematic. 1546 // The latter is the case if the differences have many sign changes. 1547 for (j = 0; j < groups.length; j++) { 1548 types[j] = "point"; 1549 le1 = groups[j].length; 1550 if (le1 < 64) { 1551 continue; 1552 } 1553 sgnChange = 0; 1554 sgn = Math.sign(groups[j][0].v); 1555 for (i = 1; i < le1; i++) { 1556 if (Math.sign(groups[j][i].v) !== sgn) { 1557 sgnChange++; 1558 sgn = Math.sign(groups[j][i].v); 1559 } 1560 } 1561 if (sgnChange * 6 > le1) { 1562 types[j] = "interval"; 1563 } 1564 } 1565 1566 return { smooth: smooth, groups: groups, types: types }; 1567 }, 1568 1569 Component: function () { 1570 this.left_isNaN = false; 1571 this.right_isNaN = false; 1572 this.left_t = null; 1573 this.right_t = null; 1574 this.t_values = []; 1575 this.x_values = []; 1576 this.y_values = []; 1577 this.len = 0; 1578 }, 1579 1580 findComponents: function (curve, mi, ma, steps) { 1581 var i, t, h, 1582 x, y, 1583 components = [], 1584 comp, 1585 comp_nr = 0, 1586 cnt = 0, 1587 cntNaNs = 0, 1588 comp_started = false, 1589 suspended = false; 1590 1591 h = (ma - mi) / steps; 1592 components[comp_nr] = new this.Component(); 1593 comp = components[comp_nr]; 1594 1595 for (i = 0, t = mi; i <= steps; i++, t += h) { 1596 x = curve.X(t, suspended); 1597 y = curve.Y(t, suspended); 1598 1599 if (isNaN(x) || isNaN(y)) { 1600 cntNaNs++; 1601 // Wait for - at least - two consecutive NaNs 1602 // This avoids starting a new component if 1603 // the function value has infinity as intermediate value. 1604 if (cntNaNs > 1 && comp_started) { 1605 // Finalize a component 1606 comp.right_isNaN = true; 1607 comp.right_t = t - h; 1608 comp.len = cnt; 1609 1610 // Prepare a new component 1611 comp_started = false; 1612 comp_nr++; 1613 components[comp_nr] = new this.Component(); 1614 comp = components[comp_nr]; 1615 cntNaNs = 0; 1616 } 1617 } else { 1618 // Now there is a non-NaN entry. 1619 if (!comp_started) { 1620 // Start the component 1621 comp_started = true; 1622 cnt = 0; 1623 if (cntNaNs > 0) { 1624 comp.left_t = t - h; 1625 comp.left_isNaN = true; 1626 } 1627 } 1628 cntNaNs = 0; 1629 // Add the value to the component 1630 comp.t_values[cnt] = t; 1631 comp.x_values[cnt] = x; 1632 comp.y_values[cnt] = y; 1633 cnt++; 1634 } 1635 if (i === 0) { 1636 suspended = true; 1637 } 1638 } 1639 if (comp_started) { 1640 comp.len = cnt; 1641 } else { 1642 components.pop(); 1643 } 1644 1645 return components; 1646 }, 1647 1648 getPointType: function (curve, pos, t_approx, t_values, x_table, y_table, len) { 1649 var x_values = x_table[0], 1650 y_values = y_table[0], 1651 full_len = t_values.length, 1652 result = { 1653 idx: pos, 1654 t: t_approx, //t_values[pos], 1655 x: x_values[pos], 1656 y: y_values[pos], 1657 type: "other" 1658 }; 1659 1660 if (pos < 5) { 1661 result.type = "borderleft"; 1662 result.idx = 0; 1663 result.t = t_values[0]; 1664 result.x = x_values[0]; 1665 result.y = y_values[0]; 1666 1667 // console.log('Border left', result.t); 1668 return result; 1669 } 1670 if (pos > len - 6) { 1671 result.type = "borderright"; 1672 result.idx = full_len - 1; 1673 result.t = t_values[full_len - 1]; 1674 result.x = x_values[full_len - 1]; 1675 result.y = y_values[full_len - 1]; 1676 1677 // console.log('Border right', result.t, full_len - 1); 1678 return result; 1679 } 1680 1681 return result; 1682 }, 1683 1684 newtonApprox: function (idx, t, h, level, table) { 1685 var i, 1686 s = 0.0; 1687 for (i = level; i > 0; i--) { 1688 s = ((s + table[i][idx]) * (t - (i - 1) * h)) / i; 1689 } 1690 return s + table[0][idx]; 1691 }, 1692 1693 thiele: function (t, recip, t_values, idx, degree) { 1694 var i, 1695 v = 0.0; 1696 for (i = degree; i > 1; i--) { 1697 v = (t - t_values[idx + i]) / (recip[i][idx + 1] - recip[i - 2][idx + 1] + v); 1698 } 1699 return recip[0][idx + 1] + (t - t_values[idx + 1]) / (recip[1][idx + 1] + v); 1700 }, 1701 1702 differenceMethodExperiments: function (component, curve) { 1703 var i, 1704 level, 1705 le, 1706 up, 1707 t_values = component.t_values, 1708 x_values = component.x_values, 1709 y_values = component.y_values, 1710 x_diffs = [], 1711 y_diffs = [], 1712 x_slopes = [], 1713 y_slopes = [], 1714 x_table = [], 1715 y_table = [], 1716 x_recip = [], 1717 y_recip = [], 1718 h, 1719 numerator, 1720 // x_med, y_med, 1721 foundCriticalPoint = 0, 1722 pos, 1723 ma, 1724 j, 1725 v, 1726 groups, 1727 criticalPoints = []; 1728 1729 h = t_values[1] - t_values[0]; 1730 x_table.push([]); 1731 y_table.push([]); 1732 x_recip.push([]); 1733 y_recip.push([]); 1734 le = y_values.length; 1735 for (i = 0; i < le; i++) { 1736 x_table[0][i] = x_values[i]; 1737 y_table[0][i] = y_values[i]; 1738 x_recip[0][i] = x_values[i]; 1739 y_recip[0][i] = y_values[i]; 1740 } 1741 1742 x_table.push([]); 1743 y_table.push([]); 1744 x_recip.push([]); 1745 y_recip.push([]); 1746 numerator = h; 1747 le = y_values.length - 1; 1748 for (i = 0; i < le; i++) { 1749 x_diffs[i] = x_values[i + 1] - x_values[i]; 1750 y_diffs[i] = y_values[i + 1] - y_values[i]; 1751 x_slopes[i] = x_diffs[i]; 1752 y_slopes[i] = y_diffs[i]; 1753 x_table[1][i] = x_diffs[i]; 1754 y_table[1][i] = y_diffs[i]; 1755 x_recip[1][i] = numerator / x_diffs[i]; 1756 y_recip[1][i] = numerator / y_diffs[i]; 1757 } 1758 le--; 1759 1760 up = Math.min(8, y_values.length - 1); 1761 for (level = 1; level < up; level++) { 1762 x_table.push([]); 1763 y_table.push([]); 1764 x_recip.push([]); 1765 y_recip.push([]); 1766 numerator *= h; 1767 for (i = 0; i < le; i++) { 1768 x_diffs[i] = x_diffs[i + 1] - x_diffs[i]; 1769 y_diffs[i] = y_diffs[i + 1] - y_diffs[i]; 1770 x_table[level + 1][i] = x_diffs[i]; 1771 y_table[level + 1][i] = y_diffs[i]; 1772 x_recip[level + 1][i] = 1773 numerator / (x_recip[level][i + 1] - x_recip[level][i]) + 1774 x_recip[level - 1][i + 1]; 1775 y_recip[level + 1][i] = 1776 numerator / (y_recip[level][i + 1] - y_recip[level][i]) + 1777 y_recip[level - 1][i + 1]; 1778 } 1779 1780 // if (level == 1) { 1781 // console.log("bends level=", level, y_diffs.toString()); 1782 // } 1783 1784 // Store point location which may be centered around 1785 // critical points. 1786 // If the lebvel is suitable, step out of the loop. 1787 groups = this._criticalPoints(y_diffs, le, level); 1788 if (groups === false) { 1789 // Its seems, the degree of the polynomial is equal to level 1790 console.log("Polynomial of degree", level); 1791 groups = []; 1792 break; 1793 } 1794 if (groups.length > 0) { 1795 foundCriticalPoint++; 1796 if (foundCriticalPoint > 1 && level % 2 === 0) { 1797 break; 1798 } 1799 } 1800 le--; 1801 } 1802 1803 // console.log("Last diffs", y_diffs, "level", level); 1804 1805 // Analyze the groups which have been found. 1806 for (i = 0; i < groups.length; i++) { 1807 // console.log("Group", i, groups[i]) 1808 // Identify the maximum difference, i.e. the center of the "problem" 1809 ma = -Infinity; 1810 for (j = 0; j < groups[i].length; j++) { 1811 v = Math.abs(groups[i][j].v); 1812 if (v > ma) { 1813 ma = v; 1814 pos = j; 1815 } 1816 } 1817 pos = Math.floor(groups[i][pos].i + level / 2); 1818 // Analyze the critical point 1819 criticalPoints.push( 1820 this.getPointType( 1821 curve, 1822 pos, 1823 t_values, 1824 x_values, 1825 y_values, 1826 x_slopes, 1827 y_slopes, 1828 le + 1 1829 ) 1830 ); 1831 } 1832 1833 return [criticalPoints, x_table, y_table, x_recip, y_recip]; 1834 }, 1835 1836 getCenterOfCriticalInterval: function (group, degree, t_values) { 1837 var ma, 1838 j, 1839 pos, 1840 v, 1841 num = 0.0, 1842 den = 0.0, 1843 h = t_values[1] - t_values[0], 1844 pos_mean, 1845 range = []; 1846 1847 // Identify the maximum difference, i.e. the center of the "problem" 1848 // If there are several equal maxima, store the positions 1849 // in the array range and determine the center of the array. 1850 1851 ma = -Infinity; 1852 range = []; 1853 for (j = 0; j < group.length; j++) { 1854 v = Math.abs(group[j].v); 1855 if (v > ma) { 1856 range = [j]; 1857 ma = v; 1858 pos = j; 1859 } else if (ma === v) { 1860 range.push(j); 1861 } 1862 } 1863 if (range.length > 0) { 1864 pos_mean = 1865 range.reduce(function (total, val) { 1866 return total + val; 1867 }, 0) / range.length; 1868 pos = Math.floor(pos_mean); 1869 pos_mean += group[0].i; 1870 } 1871 1872 if (ma < Infinity) { 1873 for (j = 0; j < group.length; j++) { 1874 num += Math.abs(group[j].v) * group[j].i; 1875 den += Math.abs(group[j].v); 1876 } 1877 pos_mean = num / den; 1878 } 1879 pos_mean += degree / 2; 1880 return [ 1881 group[pos].i + degree / 2, 1882 pos_mean, 1883 t_values[Math.floor(pos_mean)] + h * (pos_mean - Math.floor(pos_mean)) 1884 ]; 1885 }, 1886 1887 differenceMethod: function (component, curve) { 1888 var i, 1889 level, 1890 le, 1891 up, 1892 t_values = component.t_values, 1893 x_values = component.x_values, 1894 y_values = component.y_values, 1895 x_table = [], 1896 y_table = [], 1897 foundCriticalPoint = 0, 1898 degree_x = -1, 1899 degree_y = -1, 1900 pos, 1901 res, 1902 res_x, 1903 res_y, 1904 t_approx, 1905 groups = [], 1906 types, 1907 criticalPoints = []; 1908 1909 le = y_values.length; 1910 // x_table.push([]); 1911 // y_table.push([]); 1912 // for (i = 0; i < le; i++) { 1913 // x_table[0][i] = x_values[i]; 1914 // y_table[0][i] = y_values[i]; 1915 // } 1916 x_table.push(new Float64Array(x_values)); 1917 y_table.push(new Float64Array(y_values)); 1918 1919 le--; 1920 up = Math.min(12, le); 1921 for (level = 0; level < up; level++) { 1922 // Old style method: 1923 // x_table.push([]); 1924 // y_table.push([]); 1925 // for (i = 0; i < le; i++) { 1926 // x_table[level + 1][i] = x_table[level][i + 1] - x_table[level][i]; 1927 // y_table[level + 1][i] = y_table[level][i + 1] - y_table[level][i]; 1928 // } 1929 // New method: 1930 x_table.push(new Float64Array(le)); 1931 y_table.push(new Float64Array(le)); 1932 x_table[level + 1] = x_table[level].map(function (v, idx, arr) { 1933 return arr[idx + 1] - v; 1934 }); 1935 y_table[level + 1] = y_table[level].map(function (v, idx, arr) { 1936 return arr[idx + 1] - v; 1937 }); 1938 1939 // Store point location which may be centered around critical points. 1940 // If the level is suitable, step out of the loop. 1941 res_y = this._criticalInterval(y_table[level + 1], le, level); 1942 if (res_y.smooth === true) { 1943 // Its seems, the degree of the polynomial is equal to level 1944 // If the values in level + 1 are zero, it might be a polynomial of degree level. 1945 // Seems to work numerically stable until degree 6. 1946 degree_y = level; 1947 groups = []; 1948 } 1949 res_x = this._criticalInterval(x_table[level + 1], le, level); 1950 if (degree_x === -1 && res_x.smooth === true) { 1951 // Its seems, the degree of the polynomial is equal to level 1952 // If the values in level + 1 are zero, it might be a polynomial of degree level. 1953 // Seems to work numerically stable until degree 6. 1954 degree_x = level; 1955 } 1956 if (degree_y >= 0) { 1957 break; 1958 } 1959 1960 if (res_y.groups.length > 0) { 1961 foundCriticalPoint++; 1962 if (foundCriticalPoint > 2 && (level + 1) % 2 === 0) { 1963 groups = res_y.groups; 1964 types = res_y.types; 1965 break; 1966 } 1967 } 1968 le--; 1969 } 1970 1971 // console.log("Last diffs", y_table[Math.min(level + 1, up)], "level", level + 1); 1972 // Analyze the groups which have been found. 1973 for (i = 0; i < groups.length; i++) { 1974 if (types[i] === "interval") { 1975 continue; 1976 } 1977 // console.log("Group", i, groups[i], types[i], level + 1) 1978 res = this.getCenterOfCriticalInterval(groups[i], level + 1, t_values); 1979 pos = res_y[0]; 1980 pos = Math.floor(res[1]); 1981 t_approx = res[2]; 1982 // console.log("Critical points:", groups, res, pos) 1983 1984 // Analyze the type of the critical point 1985 // Result is of type 'borderleft', borderright', 'other' 1986 criticalPoints.push( 1987 this.getPointType(curve, pos, t_approx, t_values, x_table, y_table, le + 1) 1988 ); 1989 } 1990 1991 // if (level === up) { 1992 // console.log("No convergence!"); 1993 // } else { 1994 // console.log("Convergence level", level); 1995 // } 1996 return [criticalPoints, x_table, y_table, degree_x, degree_y]; 1997 }, 1998 1999 _insertPoint_v4: function (curve, crds, t, doLog) { 2000 var p, 2001 prev = null, 2002 x, 2003 y, 2004 near = 0.8; 2005 2006 if (curve.points.length > 0) { 2007 prev = curve.points[curve.points.length - 1].scrCoords; 2008 } 2009 2010 // Add regular point 2011 p = new Coords(Const.COORDS_BY_USER, crds, curve.board); 2012 2013 if (prev !== null) { 2014 x = p.scrCoords[1] - prev[1]; 2015 y = p.scrCoords[2] - prev[2]; 2016 if (x * x + y * y < near * near) { 2017 // Math.abs(p.scrCoords[1] - prev[1]) < near && 2018 // Math.abs(p.scrCoords[2] - prev[2]) < near) { 2019 return; 2020 } 2021 } 2022 2023 p._t = t; 2024 curve.points.push(p); 2025 }, 2026 2027 getInterval: function (curve, ta, tb) { 2028 var t_int, x_int, y_int; 2029 2030 //console.log('critical point', ta, tb); 2031 IntervalArithmetic.disable(); 2032 2033 t_int = IntervalArithmetic.Interval(ta, tb); 2034 curve.board.mathLib = IntervalArithmetic; 2035 curve.board.mathLibJXG = IntervalArithmetic; 2036 x_int = curve.X(t_int, true); 2037 y_int = curve.Y(t_int, true); 2038 curve.board.mathLib = Math; 2039 curve.board.mathLibJXG = JXG.Math; 2040 2041 //console.log(x_int, y_int); 2042 return y_int; 2043 }, 2044 2045 sign: function (v) { 2046 if (v < 0) { 2047 return -1; 2048 } 2049 if (v > 0) { 2050 return 1; 2051 } 2052 return 0; 2053 }, 2054 2055 handleBorder: function (curve, comp, group, x_table, y_table) { 2056 var idx = group.idx, 2057 t, 2058 t1, 2059 t2, 2060 size = 32, 2061 y_int, 2062 x, 2063 y, 2064 lo, 2065 hi, 2066 i, 2067 components2, 2068 le, 2069 h; 2070 2071 // console.log("HandleBorder at t =", t_approx); 2072 // console.log("component:", comp) 2073 // console.log("Group:", group); 2074 2075 h = comp.t_values[1] - comp.t_values[0]; 2076 if (group.type === "borderleft") { 2077 t = comp.left_isNaN ? comp.left_t : group.t - h; 2078 t1 = t; 2079 t2 = t1 + h; 2080 } else if (group.type === "borderright") { 2081 t = comp.right_isNaN ? comp.right_t : group.t + h; 2082 t2 = t; 2083 t1 = t2 - h; 2084 } else { 2085 console.log("No bordercase!!!"); 2086 } 2087 2088 components2 = this.findComponents(curve, t1, t2, size); 2089 if (components2.length === 0) { 2090 return; 2091 } 2092 if (group.type === "borderleft") { 2093 t1 = components2[0].left_t; 2094 t2 = components2[0].t_values[0]; 2095 h = components2[0].t_values[1] - components2[0].t_values[0]; 2096 t1 = t1 === null ? t2 - h : t1; 2097 t = t1; 2098 y_int = this.getInterval(curve, t1, t2); 2099 if (Type.isObject(y_int)) { 2100 lo = y_int.lo; 2101 hi = y_int.hi; 2102 2103 x = curve.X(t, true); 2104 y = y_table[1][idx] < 0 ? hi : lo; 2105 this._insertPoint_v4(curve, [1, x, y], t); 2106 } 2107 } 2108 2109 le = components2[0].t_values.length; 2110 for (i = 0; i < le; i++) { 2111 t = components2[0].t_values[i]; 2112 x = components2[0].x_values[i]; 2113 y = components2[0].y_values[i]; 2114 this._insertPoint_v4(curve, [1, x, y], t); 2115 } 2116 2117 if (group.type === "borderright") { 2118 t1 = components2[0].t_values[le - 1]; 2119 t2 = components2[0].right_t; 2120 h = components2[0].t_values[1] - components2[0].t_values[0]; 2121 t2 = t2 === null ? t1 + h : t2; 2122 2123 t = t2; 2124 y_int = this.getInterval(curve, t1, t2); 2125 if (Type.isObject(y_int)) { 2126 lo = y_int.lo; 2127 hi = y_int.hi; 2128 x = curve.X(t, true); 2129 y = y_table[1][idx] > 0 ? hi : lo; 2130 this._insertPoint_v4(curve, [1, x, y], t); 2131 } 2132 } 2133 }, 2134 2135 _seconditeration_v4: function (curve, comp, group, x_table, y_table) { 2136 var i, t1, t2, ret, components2, comp2, idx, groups2, g, x_table2, y_table2, start, le; 2137 2138 // Look at two points, hopefully left and right from the critical point 2139 t1 = comp.t_values[group.idx - 2]; 2140 t2 = comp.t_values[group.idx + 2]; 2141 components2 = this.findComponents(curve, t1, t2, 64); 2142 for (idx = 0; idx < components2.length; idx++) { 2143 comp2 = components2[idx]; 2144 ret = this.differenceMethod(comp2, curve); 2145 groups2 = ret[0]; 2146 x_table2 = ret[1]; 2147 y_table2 = ret[2]; 2148 start = 0; 2149 for (g = 0; g <= groups2.length; g++) { 2150 if (g === groups2.length) { 2151 le = comp2.len; 2152 } else { 2153 le = groups2[g].idx; 2154 } 2155 2156 // Insert all uncritical points until next critical point 2157 for (i = start; i < le; i++) { 2158 if (!isNaN(comp2.x_values[i]) && !isNaN(comp2.y_values[i])) { 2159 this._insertPoint_v4( 2160 curve, 2161 [1, comp2.x_values[i], comp2.y_values[i]], 2162 comp2.t_values[i] 2163 ); 2164 } 2165 } 2166 // Handle next critical point 2167 if (g < groups2.length) { 2168 this.handleSingularity(curve, comp2, groups2[g], x_table2, y_table2); 2169 start = groups2[g].idx + 1; 2170 } 2171 } 2172 le = comp2.len; 2173 if (idx < components2.length - 1) { 2174 this._insertPoint_v4(curve, [1, NaN, NaN], comp2.right_t); 2175 } 2176 } 2177 return this; 2178 }, 2179 2180 _recurse_v4: function (curve, t1, t2, x1, y1, x2, y2, level) { 2181 var tol = 2, 2182 t = (t1 + t2) * 0.5, 2183 x = curve.X(t, true), 2184 y = curve.Y(t, true), 2185 dx, 2186 dy; 2187 2188 //console.log("Level", level) 2189 if (level === 0) { 2190 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2191 return; 2192 } 2193 // console.log("R", t1, t2) 2194 dx = (x - x1) * curve.board.unitX; 2195 dy = (y - y1) * curve.board.unitY; 2196 // console.log("D1", Math.sqrt(dx * dx + dy * dy)) 2197 if (Mat.hypot(dx, dy) > tol) { 2198 this._recurse_v4(curve, t1, t, x1, y1, x, y, level - 1); 2199 } else { 2200 this._insertPoint_v4(curve, [1, x, y], t); 2201 } 2202 dx = (x - x2) * curve.board.unitX; 2203 dy = (y - y2) * curve.board.unitY; 2204 // console.log("D2", Math.sqrt(dx * dx + dy * dy), x-x2, y-y2) 2205 if (Mat.hypot(dx, dy) > tol) { 2206 this._recurse_v4(curve, t, t2, x, y, x2, y2, level - 1); 2207 } else { 2208 this._insertPoint_v4(curve, [1, x, y], t); 2209 } 2210 }, 2211 2212 handleSingularity: function (curve, comp, group, x_table, y_table) { 2213 var idx = group.idx, 2214 t, 2215 t1, 2216 t2, 2217 y_int, 2218 i1, 2219 i2, 2220 x, 2221 y, 2222 lo, 2223 hi, 2224 d_lft, 2225 d_rgt, 2226 d_thresh = 100, 2227 di1 = 5, 2228 di2 = 3, 2229 d1, 2230 d2; 2231 2232 t = group.t; 2233 console.log("HandleSingularity at t =", t); 2234 // console.log(comp.t_values[idx - 1], comp.y_values[idx - 1], comp.t_values[idx + 1], comp.y_values[idx + 1]); 2235 // console.log(group); 2236 2237 // Look at two points, hopefully left and right from the critical point 2238 t1 = comp.t_values[idx - di1]; 2239 t2 = comp.t_values[idx + di1]; 2240 2241 y_int = this.getInterval(curve, t1, t2); 2242 if (Type.isObject(y_int)) { 2243 lo = y_int.lo; 2244 hi = y_int.hi; 2245 } else { 2246 if (y_table[0][idx - 1] < y_table[0][idx + 1]) { 2247 lo = y_table[0][idx - 1]; 2248 hi = y_table[0][idx + 1]; 2249 } else { 2250 lo = y_table[0][idx + 1]; 2251 hi = y_table[0][idx - 1]; 2252 } 2253 } 2254 2255 x = curve.X(t, true); 2256 2257 d_lft = 2258 (y_table[0][idx - di2] - y_table[0][idx - di1]) / 2259 (comp.t_values[idx - di2] - comp.t_values[idx - di1]); 2260 d_rgt = 2261 (y_table[0][idx + di2] - y_table[0][idx + di1]) / 2262 (comp.t_values[idx + di2] - comp.t_values[idx + di1]); 2263 2264 console.log(":::", d_lft, d_rgt); 2265 2266 //this._insertPoint_v4(curve, [1, NaN, NaN], 0); 2267 2268 if (d_lft < -d_thresh) { 2269 // Left branch very steep downwards -> add the minimum 2270 this._insertPoint_v4(curve, [1, x, lo], t, true); 2271 if (d_rgt <= d_thresh) { 2272 // Right branch not very steep upwards -> interrupt the curve 2273 // I.e. it looks like -infty / (finite or infty) and not like -infty / -infty 2274 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2275 } 2276 } else if (d_lft > d_thresh) { 2277 // Left branch very steep upwards -> add the maximum 2278 this._insertPoint_v4(curve, [1, x, hi], t); 2279 if (d_rgt >= -d_thresh) { 2280 // Right branch not very steep downwards -> interrupt the curve 2281 // I.e. it looks like infty / (finite or -infty) and not like infty / infty 2282 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2283 } 2284 } else { 2285 if (lo === -Infinity) { 2286 this._insertPoint_v4(curve, [1, x, lo], t, true); 2287 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2288 } 2289 if (hi === Infinity) { 2290 this._insertPoint_v4(curve, [1, NaN, NaN], t); 2291 this._insertPoint_v4(curve, [1, x, hi], t, true); 2292 } 2293 2294 if (group.t < comp.t_values[idx]) { 2295 i1 = idx - 1; 2296 i2 = idx; 2297 } else { 2298 i1 = idx; 2299 i2 = idx + 1; 2300 } 2301 t1 = comp.t_values[i1]; 2302 t2 = comp.t_values[i2]; 2303 this._recurse_v4( 2304 curve, 2305 t1, 2306 t2, 2307 x_table[0][i1], 2308 y_table[0][i1], 2309 x_table[0][i2], 2310 y_table[0][i2], 2311 10 2312 ); 2313 2314 // x = (x_table[0][idx] - x_table[0][idx - 1]) * curve.board.unitX; 2315 // y = (y_table[0][idx] - y_table[0][idx - 1]) * curve.board.unitY; 2316 // d1 = Math.sqrt(x * x + y * y); 2317 // x = (x_table[0][idx + 1] - x_table[0][idx]) * curve.board.unitX; 2318 // y = (y_table[0][idx + 1] - y_table[0][idx]) * curve.board.unitY; 2319 // d2 = Math.sqrt(x * x + y * y); 2320 2321 // console.log("end", t1, t2, t); 2322 // if (true || (d1 > 2 || d2 > 2)) { 2323 2324 // console.log(d1, d2, y_table[0][idx]) 2325 // // Finite jump 2326 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2327 // } else { 2328 // if (lo !== -Infinity && hi !== Infinity) { 2329 // // Critical point which can be ignored 2330 // this._insertPoint_v4(curve, [1, x_table[0][idx], y_table[0][idx]], comp.t_values[idx]); 2331 // } else { 2332 // if (lo === -Infinity) { 2333 // this._insertPoint_v4(curve, [1, x, lo], t, true); 2334 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2335 // } 2336 // if (hi === Infinity) { 2337 // this._insertPoint_v4(curve, [1, NaN, NaN], t); 2338 // this._insertPoint_v4(curve, [1, x, hi], t, true); 2339 // } 2340 // } 2341 // } 2342 } 2343 if (d_rgt < -d_thresh) { 2344 // Right branch very steep downwards -> add the maximum 2345 this._insertPoint_v4(curve, [1, x, hi], t); 2346 } else if (d_rgt > d_thresh) { 2347 // Right branch very steep upwards -> add the minimum 2348 this._insertPoint_v4(curve, [1, x, lo], t); 2349 } 2350 }, 2351 2352 /** 2353 * Number of equidistant points where the function is evaluated 2354 */ 2355 steps: 1021, //2053, // 1021, 2356 2357 /** 2358 * If the absolute maximum of the set of differences is larger than 2359 * criticalThreshold * median of these values, it is regarded as critical point. 2360 * @see JXG.Math.Plot#_criticalInterval 2361 */ 2362 criticalThreshold: 1000, 2363 2364 plot_v4: function (curve, ta, tb, steps) { 2365 var i, 2366 j, 2367 le, 2368 components, 2369 idx, 2370 comp, 2371 groups, 2372 g, 2373 start, 2374 ret, 2375 x_table, 2376 y_table, 2377 t, 2378 t1, 2379 t2, 2380 good, 2381 bad, 2382 x_int, 2383 y_int, 2384 degree_x, 2385 degree_y, 2386 h = (tb - ta) / steps, 2387 Ypl = function (x) { 2388 return curve.Y(x, true); 2389 }, 2390 Ymi = function (x) { 2391 return -curve.Y(x, true); 2392 }, 2393 h2 = h * 0.5; 2394 2395 components = this.findComponents(curve, ta, tb, steps); 2396 for (idx = 0; idx < components.length; idx++) { 2397 comp = components[idx]; 2398 ret = this.differenceMethod(comp, curve); 2399 groups = ret[0]; 2400 x_table = ret[1]; 2401 y_table = ret[2]; 2402 degree_x = ret[3]; 2403 degree_y = ret[4]; 2404 2405 // if (degree_x >= 0) { 2406 // console.log("x polynomial of degree", degree_x); 2407 // } 2408 // if (degree_y >= 0) { 2409 // console.log("y polynomial of degree", degree_y); 2410 // } 2411 if (groups.length === 0 || groups[0].type !== "borderleft") { 2412 groups.unshift({ 2413 idx: 0, 2414 t: comp.t_values[0], 2415 x: comp.x_values[0], 2416 y: comp.y_values[0], 2417 type: "borderleft" 2418 }); 2419 } 2420 if (groups[groups.length - 1].type !== "borderright") { 2421 le = comp.t_values.length; 2422 groups.push({ 2423 idx: le - 1, 2424 t: comp.t_values[le - 1], 2425 x: comp.x_values[le - 1], 2426 y: comp.y_values[le - 1], 2427 type: "borderright" 2428 }); 2429 } 2430 2431 start = 0; 2432 for (g = 0; g <= groups.length; g++) { 2433 if (g === groups.length) { 2434 le = comp.len; 2435 } else { 2436 le = groups[g].idx - 1; 2437 } 2438 2439 good = 0; 2440 bad = 0; 2441 // Insert all uncritical points until next critical point 2442 for (i = start; i < le - 2; i++) { 2443 this._insertPoint_v4( 2444 curve, 2445 [1, comp.x_values[i], comp.y_values[i]], 2446 comp.t_values[i] 2447 ); 2448 j = Math.max(0, i - 2); 2449 // Add more points in critical intervals 2450 if ( 2451 //degree_y === -1 && // No polynomial 2452 i >= start + 3 && 2453 i < le - 3 && // Do not do this if too close to a critical point 2454 y_table.length > 3 && 2455 Math.abs(y_table[2][i]) > 0.2 * Math.abs(y_table[0][i]) 2456 ) { 2457 t = comp.t_values[i]; 2458 h2 = h * 0.25; 2459 y_int = this.getInterval(curve, t, t + h); 2460 if (Type.isObject(y_int)) { 2461 if (y_table[2][i] > 0) { 2462 this._insertPoint_v4(curve, [1, t + h2, y_int.lo], t + h2); 2463 } else { 2464 this._insertPoint_v4( 2465 curve, 2466 [1, t + h - h2, y_int.hi], 2467 t + h - h2 2468 ); 2469 } 2470 } else { 2471 t1 = Numerics.fminbr(Ypl, [t, t + h]); 2472 t2 = Numerics.fminbr(Ymi, [t, t + h]); 2473 if (t1 < t2) { 2474 this._insertPoint_v4( 2475 curve, 2476 [1, curve.X(t1, true), curve.Y(t1, true)], 2477 t1 2478 ); 2479 this._insertPoint_v4( 2480 curve, 2481 [1, curve.X(t2, true), curve.Y(t2, true)], 2482 t2 2483 ); 2484 } else { 2485 this._insertPoint_v4( 2486 curve, 2487 [1, curve.X(t2, true), curve.Y(t2, true)], 2488 t2 2489 ); 2490 this._insertPoint_v4( 2491 curve, 2492 [1, curve.X(t1, true), curve.Y(t1, true)], 2493 t1 2494 ); 2495 } 2496 } 2497 bad++; 2498 } else { 2499 good++; 2500 } 2501 } 2502 // console.log("GOOD", good, "BAD", bad); 2503 2504 // Handle next critical point 2505 if (g < groups.length) { 2506 //console.log("critical point / interval", groups[g]); 2507 2508 i = groups[g].idx; 2509 if (groups[g].type === "borderleft" || groups[g].type === "borderright") { 2510 this.handleBorder(curve, comp, groups[g], x_table, y_table); 2511 } else { 2512 this._seconditeration_v4(curve, comp, groups[g], x_table, y_table); 2513 } 2514 2515 start = groups[g].idx + 1 + 1; 2516 } 2517 } 2518 2519 le = comp.len; 2520 if (idx < components.length - 1) { 2521 this._insertPoint_v4(curve, [1, NaN, NaN], comp.right_t); 2522 } 2523 } 2524 }, 2525 2526 /** 2527 * Updates the data points of a parametric curve, plotVersion 4. This version is used if {@link JXG.Curve#plotVersion} is <tt>4</tt>. 2528 * @param {JXG.Curve} curve JSXGraph curve element 2529 * @param {Number} mi Left bound of curve 2530 * @param {Number} ma Right bound of curve 2531 * @returns {JXG.Curve} Reference to the curve object. 2532 */ 2533 updateParametricCurve_v4: function (curve, mi, ma) { 2534 var ta, tb, w2, bbox; 2535 2536 if (curve.xterm === "x") { 2537 // For function graphs we can restrict the plot interval 2538 // to the visible area +plus margin 2539 bbox = curve.board.getBoundingBox(); 2540 w2 = (bbox[2] - bbox[0]) * 0.3; 2541 // h2 = (bbox[1] - bbox[3]) * 0.3; 2542 ta = Math.max(mi, bbox[0] - w2); 2543 tb = Math.min(ma, bbox[2] + w2); 2544 } else { 2545 ta = mi; 2546 tb = ma; 2547 } 2548 2549 curve.points = []; 2550 2551 //console.log("--------------------"); 2552 this.plot_v4(curve, ta, tb, this.steps); 2553 2554 curve.numberPoints = curve.points.length; 2555 //console.log(curve.numberPoints); 2556 }, 2557 2558 //---------------------------------------------------------------------- 2559 // Plot algorithm alias 2560 //---------------------------------------------------------------------- 2561 2562 /** 2563 * Updates the data points of a parametric curve, alias for {@link JXG.Curve#updateParametricCurve_v2}. 2564 * This is needed for backwards compatibility, if this method has been 2565 * used directly in an application. 2566 * @param {JXG.Curve} curve JSXGraph curve element 2567 * @param {Number} mi Left bound of curve 2568 * @param {Number} ma Right bound of curve 2569 * @returns {JXG.Curve} Reference to the curve object. 2570 * 2571 * @see JXG.Curve#updateParametricCurve_v2 2572 */ 2573 updateParametricCurve: function (curve, mi, ma) { 2574 return this.updateParametricCurve_v2(curve, mi, ma); 2575 } 2576 }; 2577 2578 export default Mat.Plot; 2579