JavaScript must be enabled in order for you to use JSXGraph and JSXGraph reference. However, it seems JavaScript is either disabled or not supported by your browser.

Class Index | File Index

Elements
Classes

Namespace JXG.Math.Numerics


      ↳ JXG.Math.Numerics



Defined in: numerics.js.

Namespace Summary
Constructor Attributes Constructor Name and Description
 
The JXG.Math.Numerics namespace holds numerical algorithms, constants, and variables.
Field Summary
Field Attributes Field Name and Description
<static>  
JXG.Math.Numerics.maxIterationsMinimize
Maximum number of iterations in JXG.Math.Numerics.fminbr
<static>  
JXG.Math.Numerics.maxIterationsRoot
Maximum number of iterations in JXG.Math.Numerics.fzero and JXG.Math.Numerics.chandrupatla
Method Summary
Method Attributes Method Name and Description
<private> <static>  
JXG.Math.Numerics._gaussKronrod(interval, f, n, xgk, wg, wgk, resultObj)
Generic Gauss-Kronrod quadrature algorithm.
<private> <static>  
JXG.Math.Numerics._initCubicPoly(x1, x2, t1, t2)
Determine the coefficients of a cardinal spline polynom, See https://stackoverflow.com/questions/9489736/catmull-rom-curve-with-no-cusps-and-no-self-intersections
<private> <static>  
JXG.Math.Numerics._rescale_error(err, result_abs, result_asc)
Scale error in Gauss Kronrod quadrature.
<private> <static>  
JXG.Math.Numerics._riemannValue(x, f, type, delta)
Evaluate the function term for {@see #riemann}.
<private> <static>  
JXG.Math.Numerics._workspace(interval, n)
Generate workspace object for JXG.Math.Numerics.Qag.
<static>  
JXG.Math.Numerics.backwardSolve(R, b, canModify)
Solves a system of linear equations given by the right triangular matrix R and vector b.
<static>  
JXG.Math.Numerics.bezier(points)
Computes the cubic Bezier curve through a given set of points.
<static>  
JXG.Math.Numerics.bspline(points, order)
Computes the B-spline curve of order k (order = degree+1) through a given set of points.
<static>  
JXG.Math.Numerics.CardinalSpline(points, tau, type)
Computes the cubic cardinal spline curve through a given set of points.
<static>  
JXG.Math.Numerics.CatmullRomSpline(points, type)
Computes the cubic Catmull-Rom spline curve through a given set of points.
<static>  
JXG.Math.Numerics.chandrupatla(f, x0, object)
Find zero of an univariate function f.
<static>  
JXG.Math.Numerics.D(f, obj)
Numerical (symmetric) approximation of derivative.
<static>  
JXG.Math.Numerics.det(mat)
Computes the determinant of a square nxn matrix with the Gauss-Bareiss algorithm.
<static>  
JXG.Math.Numerics.findBracket(f, x0, object)
Given a value x_0, this function tries to find a second value x_1 such that the function f has opposite signs at x_0 and x_1.
<static>  
JXG.Math.Numerics.fminbr(f, x0, context)
Find minimum of an univariate function f.
<static>  
JXG.Math.Numerics.fzero(f, x0, object)
Find zero of an univariate function f.
<static>  
JXG.Math.Numerics.Gauss(A, b)
Solves a system of linear equations given by A and b using the Gauss-Jordan-elimination.
<private> <static>  
JXG.Math.Numerics.gaussBareiss(mat)
Gauss-Bareiss algorithm to compute the determinant of matrix without fractions.
<static>  
JXG.Math.Numerics.GaussKronrod15(interval, f, resultObj)
15 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
<static>  
JXG.Math.Numerics.GaussKronrod21(interval, f, resultObj)
21 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
<static>  
JXG.Math.Numerics.GaussKronrod31(interval, f, resultObj)
31 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
<static>  
JXG.Math.Numerics.GaussLegendre(interval, f, config)
Calculates the integral of function f over interval using Gauss-Legendre quadrature.
<static>  
JXG.Math.Numerics.generalizedNewton(c1, c2, t1ini, t2ini)
Compute an intersection of the curves c1 and c2 with a generalized Newton method.
<static>  
JXG.Math.Numerics.generatePolynomialTerm(coeffs, deg, varname, prec)
Generate a string containing the function term of a polynomial.
<static>  
JXG.Math.Numerics.glomin(f, x0)
Purpose: GLOMIN seeks a global minimum of a function F(X) in an interval [A,B].
<static>  
JXG.Math.Numerics.I(interval, f)
Integral of function f over interval.
<static>  
JXG.Math.Numerics.Jacobi(Ain)
Compute the Eigenvalues and Eigenvectors of a symmetric 3x3 matrix with the Jacobi method Adaption of a FORTRAN program by Ed Wilson, Dec.
<static>  
JXG.Math.Numerics.lagrangePolynomial(p)
Computes the polynomial through a given set of coordinates in Lagrange form.
<static>  
JXG.Math.Numerics.lagrangePolynomialCoefficients(points)
Determine the Lagrange polynomial through an array of points and return the coefficients of the polynomial as array.
<static>  
JXG.Math.Numerics.lagrangePolynomialTerm(points, digits, param, dot)
Determine the Lagrange polynomial through an array of points and return the term of the polynomial as string.
<static>  
JXG.Math.Numerics.Neville(p)
Returns the Lagrange polynomials for curves with equidistant nodes, see Jean-Paul Berrut, Lloyd N.
<static>  
JXG.Math.Numerics.Newton(f, x, context)
Newton's method to find roots of a funtion in one variable.
<static>  
JXG.Math.Numerics.NewtonCotes(interval, f, config)
Calculates the integral of function f over interval using Newton-Cotes-algorithm.
<static>  
JXG.Math.Numerics.polzeros(a, deg, tol, max_it, initial_values)
Determine all roots of a polynomial with real or complex coefficients by using the iterative method attributed to Weierstrass, Durand, Kerner, Aberth, and Ehrlich.
<static>  
JXG.Math.Numerics.Qag(interval, f, config)
Quadrature algorithm qag from QUADPACK.
<static>  
JXG.Math.Numerics.RamerDouglasPeucker(pts, eps)
Implements the Ramer-Douglas-Peucker algorithm.
<static> <deprecated>  
JXG.Math.Numerics.RamerDouglasPeuker(pts, eps)
Old name for the implementation of the Ramer-Douglas-Peucker algorithm.
<static>  
JXG.Math.Numerics.regressionPolynomial(degree, dataX, dataY)
Computes the regression polynomial of a given degree through a given set of coordinates.
<static>  
JXG.Math.Numerics.riemann(f, n, type, start, end)
Helper function to create curve which displays Riemann sums.
<static>  
JXG.Math.Numerics.riemannsum(f, n, type, start, end)
Approximate the integral by Riemann sums.
<static>  
JXG.Math.Numerics.Romberg(interval, f, config)
Calculates the integral of function f over interval using Romberg iteration.
<static>  
JXG.Math.Numerics.root(f, x, context)
Abstract method to find roots of univariate functions, which - for the time being - is an alias for JXG.Math.Numerics.chandrupatla.
<static>  
JXG.Math.Numerics.rungeKutta(butcher, x0, I, N, f)
Solve initial value problems numerically using explicit Runge-Kutta methods.
<static>  
JXG.Math.Numerics.splineDef(x, y)
Calculates second derivatives at the given knots.
<static>  
JXG.Math.Numerics.splineEval(x0, x, y, F)
Evaluate points on spline.
<static>  
JXG.Math.Numerics.Visvalingam(pts, numPoints)
Implements the Visvalingam-Whyatt algorithm.
Namespace Detail
JXG.Math.Numerics
The JXG.Math.Numerics namespace holds numerical algorithms, constants, and variables.
Field Detail
<static> {Number} JXG.Math.Numerics.maxIterationsMinimize
Maximum number of iterations in JXG.Math.Numerics.fminbr
Default Value:
500

<static> {Number} JXG.Math.Numerics.maxIterationsRoot
Maximum number of iterations in JXG.Math.Numerics.fzero and JXG.Math.Numerics.chandrupatla
Default Value:
80
Method Detail
<private> <static> {Number} JXG.Math.Numerics._gaussKronrod(interval, f, n, xgk, wg, wgk, resultObj)
Generic Gauss-Kronrod quadrature algorithm. Internal method used in JXG.Math.Numerics.GaussKronrod15, JXG.Math.Numerics.GaussKronrod21, JXG.Math.Numerics.GaussKronrod31. Taken from QUADPACK.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Number} n
order
{Array} xgk
Kronrod quadrature abscissae
{Array} wg
Weights of the Gauss rule
{Array} wgk
Weights of the Kronrod rule
{Object} resultObj
Object returning resultObj.abserr, resultObj.resabs, resultObj.resasc. See the library QUADPACK for an explanation.
Returns:
{Number} Integral value of f over interval

<private> <static> {Array} JXG.Math.Numerics._initCubicPoly(x1, x2, t1, t2)
Determine the coefficients of a cardinal spline polynom, See https://stackoverflow.com/questions/9489736/catmull-rom-curve-with-no-cusps-and-no-self-intersections
Parameters:
{Number} x1
point 1
{Number} x2
point 2
{Number} t1
tangent slope 1
{Number} t2
tangent slope 2
Returns:
{Array} coefficents array c for the polynomial t maps to c[0] + c[1]*t + c[2]*t*t + c[3]*t*t*t

<private> <static> JXG.Math.Numerics._rescale_error(err, result_abs, result_asc)
Scale error in Gauss Kronrod quadrature. Internal method used in JXG.Math.Numerics._gaussKronrod.
Parameters:
err
result_abs
result_asc

<private> <static> {Number} JXG.Math.Numerics._riemannValue(x, f, type, delta)
Evaluate the function term for {@see #riemann}.
Parameters:
{Number} x
function argument
{function} f
JavaScript function returning a number
{String} type
Name of the Riemann sum type, e.g. 'lower', see {@see #riemann}.
{Number} delta
Width of the bars in user coordinates
Returns:
{Number} Upper (delta > 0) or lower (delta < 0) value of the bar containing x of the Riemann sum.

<private> <static> {Object} JXG.Math.Numerics._workspace(interval, n)
Generate workspace object for JXG.Math.Numerics.Qag.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{Number} n
Max. limit
Returns:
{Object} Workspace object

<static> {Array} JXG.Math.Numerics.backwardSolve(R, b, canModify)
Solves a system of linear equations given by the right triangular matrix R and vector b.
Parameters:
{Array} R
Right triangular matrix represented by an array of rows. All entries a_(i,j) with i < j are ignored.
{Array} b
Right hand side of the linear equation system.
{Boolean} canModify Optional, Default: false
If true, the right hand side vector is allowed to be changed by this method.
Returns:
{Array} An array representing a vector that solves the system of linear equations.

<static> {Array} JXG.Math.Numerics.bezier(points)
Computes the cubic Bezier curve through a given set of points.
Parameters:
{Array} points
Array consisting of 3*k+1 JXG.Points. The points at position k with k mod 3 = 0 are the data points, points at position k with k mod 3 = 1 or 2 are the control points.
Returns:
{Array} An array consisting of two functions of one parameter t which return the x resp. y coordinates of the Bezier curve in t, one zero value, and a third function accepting no parameters and returning one third of the length of the points.

<static> {Array} JXG.Math.Numerics.bspline(points, order)
Computes the B-spline curve of order k (order = degree+1) through a given set of points.
Parameters:
{Array} points
Array consisting of JXG.Points.
{Number} order
Order of the B-spline curve.
Returns:
{Array} An Array consisting of four components: Two functions each of one parameter t which return the x resp. y coordinates of the B-spline curve in t, a zero value, and a function simply returning the length of the points array minus one.

<static> {Array} JXG.Math.Numerics.CardinalSpline(points, tau, type)
Computes the cubic cardinal spline curve through a given set of points. The curve is uniformly parametrized. Two artificial control points at the beginning and the end are added. The implementation (especially the centripetal parametrization) is from https://stackoverflow.com/questions/9489736/catmull-rom-curve-with-no-cusps-and-no-self-intersections .
Parameters:
{Array} points
Array consisting of JXG.Points.
{Number|Function} tau
The tension parameter, either a constant number or a function returning a number. This number is between 0 and 1. tau=1/2 give Catmull-Rom splines.
{String} type
(Optional) parameter which allows to choose between "uniform" (default) and "centripetal" parameterization. Thus the two possible values are "uniform" or "centripetal".
Returns:
{Array} An Array consisting of four components: Two functions each of one parameter t which return the x resp. y coordinates of the Catmull-Rom-spline curve in t, a zero value, and a function simply returning the length of the points array minus three.

<static> {Array} JXG.Math.Numerics.CatmullRomSpline(points, type)
Computes the cubic Catmull-Rom spline curve through a given set of points. The curve is uniformly parametrized. The curve is the cardinal spline curve for tau=0.5. Two artificial control points at the beginning and the end are added.
Parameters:
{Array} points
Array consisting of JXG.Points.
{String} type
(Optional) parameter which allows to choose between "uniform" (default) and "centripetal" parameterization. Thus the two possible values are "uniform" or "centripetal".
Returns:
{Array} An Array consisting of four components: Two functions each of one parameter t which return the x resp. y coordinates of the Catmull-Rom-spline curve in t, a zero value, and a function simply returning the length of the points array minus three.

<static> {Number} JXG.Math.Numerics.chandrupatla(f, x0, object)
Find zero of an univariate function f.
Parameters:
{function} f
Function, whose root is to be found
{Array|Number} x0
Start value or start interval enclosing the root
{Object} object
Parent object in case f is method of it
Returns:
{Number} the approximation of the root Algorithm: Chandrupatla's method, see Tirupathi R. Chandrupatla, "A new hybrid quadratic/bisection algorithm for finding the zero of a nonlinear function without using derivatives", Advances in Engineering Software, Volume 28, Issue 3, April 1997, Pages 145-149. If x0 is an array containing lower and upper bound for the zero algorithm 748 is applied. Otherwise, if x0 is a number, the algorithm tries to bracket a zero of f starting from x0. If this fails, we fall back to Newton's method.
See:
JXG.Math.Numerics.root
JXG.Math.Numerics.fzero

<static> {function} JXG.Math.Numerics.D(f, obj)
Numerical (symmetric) approximation of derivative. suspendUpdate is piped through, see JXG.Curve#updateCurve and JXG.Curve#hasPoint.
Parameters:
{function} f
Function in one variable to be differentiated.
{object} obj Optional
Optional object that is treated as "this" in the function body. This is useful, if the function is a method of an object and contains a reference to its parent object via "this".
Returns:
{function} Derivative function of a given function f.

<static> {Number} JXG.Math.Numerics.det(mat)
Computes the determinant of a square nxn matrix with the Gauss-Bareiss algorithm.
Parameters:
{Array} mat
Matrix.
Returns:
{Number} The determinant pf the matrix mat. The empty matrix returns 0.

<static> {Array} JXG.Math.Numerics.findBracket(f, x0, object)
Given a value x_0, this function tries to find a second value x_1 such that the function f has opposite signs at x_0 and x_1. The return values have to be tested if the method succeeded.
Parameters:
{Function} f
Function, whose root is to be found
{Number} x0
Start value
{Object} object
Parent object in case f is method of it
Returns:
{Array} [x_0, f(x_0), x_1, f(x_1)] in case that x_0 <= x_1 or [x_1, f(x_1), x_0, f(x_0)] in case that x_1 < x_0.
See:
JXG.Math.Numerics.fzero
JXG.Math.Numerics.chandrupatla

<static> {Number} JXG.Math.Numerics.fminbr(f, x0, context)
Find minimum of an univariate function f.

Algorithm: G.Forsythe, M.Malcolm, C.Moler, Computer methods for mathematical computations. M., Mir, 1980, p.180 of the Russian edition

Parameters:
{function} f
Function, whose minimum is to be found
{Array} x0
Start interval enclosing the minimum
{Object} context
Parent object in case f is method of it
Returns:
{Number} the approximation of the minimum value position

<static> {Number} JXG.Math.Numerics.fzero(f, x0, object)
Find zero of an univariate function f.
Parameters:
{function} f
Function, whose root is to be found
{Array|Number} x0
Start value or start interval enclosing the root
{Object} object
Parent object in case f is method of it
Returns:
{Number} the approximation of the root Algorithm: Brent's root finder from G.Forsythe, M.Malcolm, C.Moler, Computer methods for mathematical computations. M., Mir, 1980, p.180 of the Russian edition https://www.netlib.org/c/brent.shar If x0 is an array containing lower and upper bound for the zero algorithm 748 is applied. Otherwise, if x0 is a number, the algorithm tries to bracket a zero of f starting from x0. If this fails, we fall back to Newton's method.
See:
JXG.Math.Numerics.chandrupatla
JXG.Math.Numerics.root

<static> {Array} JXG.Math.Numerics.Gauss(A, b)
Solves a system of linear equations given by A and b using the Gauss-Jordan-elimination. The algorithm runs in-place. I.e. the entries of A and b are changed.
Parameters:
{Array} A
Square matrix represented by an array of rows, containing the coefficients of the lineare equation system.
{Array} b
A vector containing the linear equation system's right hand side.
Throws:
{Error}
If a non-square-matrix is given or if b has not the right length or A's rank is not full.
Returns:
{Array} A vector that solves the linear equation system.

<private> <static> JXG.Math.Numerics.gaussBareiss(mat)
Gauss-Bareiss algorithm to compute the determinant of matrix without fractions. See Henri Cohen, "A Course in Computational Algebraic Number Theory (Graduate texts in mathematics; 138)", Springer-Verlag, ISBN 3-540-55640-0 / 0-387-55640-0 Third, Corrected Printing 1996 "Algorithm 2.2.6", pg. 52-53
Parameters:
{Array} mat
Matrix
Returns:
Number

<static> {Number} JXG.Math.Numerics.GaussKronrod15(interval, f, resultObj)
15 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} resultObj
Object returning resultObj.abserr, resultObj.resabs, resultObj.resasc. See the library QUADPACK for an explanation.
Returns:
{Number} Integral value of f over interval

<static> {Number} JXG.Math.Numerics.GaussKronrod21(interval, f, resultObj)
21 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} resultObj
Object returning resultObj.abserr, resultObj.resabs, resultObj.resasc. See the library QUADPACK for an explanation.
Returns:
{Number} Integral value of f over interval

<static> {Number} JXG.Math.Numerics.GaussKronrod31(interval, f, resultObj)
31 point Gauss-Kronrod quadrature algorithm, see the library QUADPACK
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} resultObj
Object returning resultObj.abserr, resultObj.resabs, resultObj.resasc. See the library QUADPACK for an explanation.
Returns:
{Number} Integral value of f over interval

<static> {Number} JXG.Math.Numerics.GaussLegendre(interval, f, config)
Calculates the integral of function f over interval using Gauss-Legendre quadrature.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} config Optional
The algorithm setup. Accepted property is the order n of type number. n is allowed to take values between 2 and 18, default value is 12.
{Number} config.n Optional, Default: 16
Returns:
{Number} Integral value of f over interval
Examples:
function f(x) {
  return x*x;
}

// calculates integral of f from 0 to 2.
var area1 = JXG.Math.Numerics.GaussLegendre([0, 2], f);

// the same with an anonymous function
var area2 = JXG.Math.Numerics.GaussLegendre([0, 2], function (x) { return x*x; });

// use 16 point Gauss-Legendre rule.
var area3 = JXG.Math.Numerics.GaussLegendre([0, 2], f,
                                  {n: 16});

<static> {JXG.Coords} JXG.Math.Numerics.generalizedNewton(c1, c2, t1ini, t2ini)
Compute an intersection of the curves c1 and c2 with a generalized Newton method. We want to find values t1, t2 such that c1(t1) = c2(t2), i.e. (c1_x(t1)-c2_x(t2),c1_y(t1)-c2_y(t2)) = (0,0). We set (e,f) := (c1_x(t1)-c2_x(t2),c1_y(t1)-c2_y(t2)) The Jacobian J is defined by J = (a, b) (c, d) where a = c1_x'(t1) b = -c2_x'(t2) c = c1_y'(t1) d = -c2_y'(t2) The inverse J^(-1) of J is equal to (d, -b)/ (-c, a) / (ad-bc) Then, (t1new, t2new) := (t1,t2) - J^(-1)*(e,f). If the function meetCurveCurve possesses the properties t1memo and t2memo then these are taken as start values for the Newton algorithm. After stopping of the Newton algorithm the values of t1 and t2 are stored in t1memo and t2memo.
Parameters:
{JXG.Curve} c1
Curve, Line or Circle
{JXG.Curve} c2
Curve, Line or Circle
{Number} t1ini
start value for t1
{Number} t2ini
start value for t2
Returns:
{JXG.Coords} intersection point

<static> {String} JXG.Math.Numerics.generatePolynomialTerm(coeffs, deg, varname, prec)
Generate a string containing the function term of a polynomial.
Parameters:
{Array} coeffs
Coefficients of the polynomial. The position i belongs to x^i.
{Number} deg
Degree of the polynomial
{String} varname
Name of the variable (usually 'x')
{Number} prec
Precision
Returns:
{String} A string containg the function term of the polynomial.

<static> JXG.Math.Numerics.glomin(f, x0)
Purpose: GLOMIN seeks a global minimum of a function F(X) in an interval [A,B]. Discussion: This function assumes that F(X) is twice continuously differentiable over [A,B] and that F''(X) <= M for all X in [A,B]. Licensing: This code is distributed under the GNU LGPL license. Modified: 17 April 2008 Author: Original FORTRAN77 version by Richard Brent. C version by John Burkardt. https://people.math.sc.edu/Burkardt/c_src/brent/brent.c Reference: Richard Brent, Algorithms for Minimization Without Derivatives, Dover, 2002, ISBN: 0-486-41998-3, LC: QA402.5.B74. Parameters: Input, double A, B, the endpoints of the interval. It must be the case that A < B. Input, double C, an initial guess for the global minimizer. If no good guess is known, C = A or B is acceptable. Input, double M, the bound on the second derivative. Input, double MACHEP, an estimate for the relative machine precision. Input, double E, a positive tolerance, a bound for the absolute error in the evaluation of F(X) for any X in [A,B]. Input, double T, a positive error tolerance. Input, double F (double x ), a user-supplied function whose global minimum is being sought. Output, double *X, the estimated value of the abscissa for which F attains its global minimum value in [A,B]. Output, double GLOMIN, the value F(X).
Parameters:
f
x0

<static> {Number} JXG.Math.Numerics.I(interval, f)
Integral of function f over interval.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
Returns:
{Number} The value of the integral of f over interval
See:
JXG.Math.Numerics.NewtonCotes
JXG.Math.Numerics.Romberg
JXG.Math.Numerics.Qag

<static> {Array} JXG.Math.Numerics.Jacobi(Ain)
Compute the Eigenvalues and Eigenvectors of a symmetric 3x3 matrix with the Jacobi method Adaption of a FORTRAN program by Ed Wilson, Dec. 25, 1990
Parameters:
{Array} Ain
A symmetric 3x3 matrix.
Returns:
{Array} [A,V] the matrices A and V. The diagonal of A contains the Eigenvalues, V contains the Eigenvectors.

<static> {function} JXG.Math.Numerics.lagrangePolynomial(p)
Computes the polynomial through a given set of coordinates in Lagrange form. Returns the Lagrange polynomials, see Jean-Paul Berrut, Lloyd N. Trefethen: Barycentric Lagrange Interpolation, SIAM Review, Vol 46, No 3, (2004) 501-517.

It possesses the method getTerm() which returns the string containing the function term of the polynomial.

Parameters:
{Array} p
Array of JXG.Points
Returns:
{function} A function of one parameter which returns the value of the polynomial, whose graph runs through the given points.
Examples:
var p = [];
p[0] = board.create('point', [-1,2], {size:4});
p[1] = board.create('point', [0,3], {size:4});
p[2] = board.create('point', [1,1], {size:4});
p[3] = board.create('point', [3,-1], {size:4});
var f = JXG.Math.Numerics.lagrangePolynomial(p);
var graph = board.create('functiongraph', [f,-10, 10], {strokeWidth:3});


					
					
var points = [];
points[0] = board.create('point', [-1,2], {size:4});
points[1] = board.create('point', [0, 0], {size:4});
points[2] = board.create('point', [2, 1], {size:4});

var f = JXG.Math.Numerics.lagrangePolynomial(points);
var graph = board.create('functiongraph', [f,-10, 10], {strokeWidth:3});
var txt = board.create('text', [-3, -4,  () => f.getTerm(2, 't', ' * ')], {fontSize: 16});
var txt2 = board.create('text', [-3, -6,  () => f.getCoefficients()], {fontSize: 12});


					
					

<static> {Function} JXG.Math.Numerics.lagrangePolynomialCoefficients(points)
Determine the Lagrange polynomial through an array of points and return the coefficients of the polynomial as array. The leading coefficient is at position 0.
Parameters:
{Array} points
Array of JXG.Points
Returns:
{Function} returning the coefficients of the Lagrange polynomial through the supplied points.
Examples:
var points = [];
points[0] = board.create('point', [-1,2], {size:4});
points[1] = board.create('point', [0, 0], {size:4});
points[2] = board.create('point', [2, 1], {size:4});

var f = JXG.Math.Numerics.lagrangePolynomial(points);
var graph = board.create('functiongraph', [f,-10, 10], {strokeWidth:3});

var f_arr = JXG.Math.Numerics.lagrangePolynomialCoefficients(points);
var txt = board.create('text', [1, -4, f_arr], {fontSize: 10});


					
					

<static> {Function} JXG.Math.Numerics.lagrangePolynomialTerm(points, digits, param, dot)
Determine the Lagrange polynomial through an array of points and return the term of the polynomial as string.
Parameters:
{Array} points
Array of JXG.Points
{Number} digits
Number of decimal digits of the coefficients
{String} param
Name of the parameter. Default: 'x'.
{String} dot
Multiplication symbol. Default: ' * '.
Returns:
{Function} returning the Lagrange polynomial term through the supplied points as string
Examples:
var points = [];
points[0] = board.create('point', [-1,2], {size:4});
points[1] = board.create('point', [0, 0], {size:4});
points[2] = board.create('point', [2, 1], {size:4});

var f = JXG.Math.Numerics.lagrangePolynomial(points);
var graph = board.create('functiongraph', [f,-10, 10], {strokeWidth:3});

var f_txt = JXG.Math.Numerics.lagrangePolynomialTerm(points, 2, 't', ' * ');
var txt = board.create('text', [-3, -4, f_txt], {fontSize: 16});


					
					

<static> {Array} JXG.Math.Numerics.Neville(p)
Returns the Lagrange polynomials for curves with equidistant nodes, see Jean-Paul Berrut, Lloyd N. Trefethen: Barycentric Lagrange Interpolation, SIAM Review, Vol 46, No 3, (2004) 501-517. The graph of the parametric curve [x(t),y(t)] runs through the given points.
Parameters:
{Array} p
Array of JXG.Points
Returns:
{Array} An array consisting of two functions x(t), y(t) which define a parametric curve f(t) = (x(t), y(t)), a number x1 (which equals 0) and a function x2 defining the curve's domain. That means the curve is defined between x1 and x2(). x2 returns the (length of array p minus one).
Examples:
var p = [];

p[0] = board.create('point', [0, -2], {size:2, name: 'C(a)'});
p[1] = board.create('point', [-1.5, 5], {size:2, name: ''});
p[2] = board.create('point', [1, 4], {size:2, name: ''});
p[3] = board.create('point', [3, 3], {size:2, name: 'C(b)'});

// Curve
var fg = JXG.Math.Numerics.Neville(p);
var graph = board.create('curve', fg, {strokeWidth:3, strokeOpacity:0.5});


					
					

<static> {Number} JXG.Math.Numerics.Newton(f, x, context)
Newton's method to find roots of a funtion in one variable.
Parameters:
{function} f
We search for a solution of f(x)=0.
{Number} x
initial guess for the root, i.e. start value.
{Object} context
optional object that is treated as "this" in the function body. This is useful if the function is a method of an object and contains a reference to its parent object via "this".
Returns:
{Number} A root of the function f.

<static> {Number} JXG.Math.Numerics.NewtonCotes(interval, f, config)
Calculates the integral of function f over interval using Newton-Cotes-algorithm.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} config Optional
The algorithm setup. Accepted properties are number_of_nodes of type number and integration_type with value being either 'trapez', 'simpson', or 'milne'.
{Number} config.number_of_nodes Optional, Default: 28
{String} config.integration_type Optional, Default: 'milne'
Possible values are 'milne', 'simpson', 'trapez'
Throws:
{Error}
If config.number_of_nodes doesn't match config.integration_type an exception is thrown. If you want to use simpson rule respectively milne rule config.number_of_nodes must be dividable by 2 respectively 4.
Returns:
{Number} Integral value of f over interval
Examples:
function f(x) {
  return x*x;
}

// calculates integral of f from 0 to 2.
var area1 = JXG.Math.Numerics.NewtonCotes([0, 2], f);

// the same with an anonymous function
var area2 = JXG.Math.Numerics.NewtonCotes([0, 2], function (x) { return x*x; });

// use trapez rule with 16 nodes
var area3 = JXG.Math.Numerics.NewtonCotes([0, 2], f,
                                  {number_of_nodes: 16, integration_type: 'trapez'});

<static> {Array} JXG.Math.Numerics.polzeros(a, deg, tol, max_it, initial_values)
Determine all roots of a polynomial with real or complex coefficients by using the iterative method attributed to Weierstrass, Durand, Kerner, Aberth, and Ehrlich. In particular, the iteration method with cubic convergence is used that is usually attributed to Ehrlich-Aberth.

The returned roots are sorted with respect to their real values.

This method makes use of the JSXGraph classes JXG.Complex and JXG.C to handle complex numbers.

Parameters:
{Array} a
Array of coefficients of the polynomial a[0] + a[1]*x+ a[2]*x**2... The coefficients are of type Number or JXG.Complex.
{Number} deg Optional
Optional degree of the polynomial. Otherwise all entries are taken, with leading zeros removed.
{Number} tol Optional, Default: Number.EPSILON
Approximation tolerance
{Number} max_it Optional, Default: 30
Maximum number of iterations
{Array} initial_values Optional, Default: null
Array of initial values for the roots. If not given, starting values are determined by the method of Ozawa.
Returns:
{Array} Array of complex numbers (of JXG.Complex) approximating the roots of the polynomial.
See:
JXG.Complex
JXG.C
Examples:
// Polynomial p(z) = -1 + 1z^2
var i, roots,
    p = [-1, 0, 1];

roots = JXG.Math.Numerics.polzeros(p);
for (i = 0; i < roots.length; i++) {
    console.log(i, roots[i].toString());
}
// Output:
  0 -1 + -3.308722450212111e-24i
  1 1 + 0i
// Polynomial p(z) = -1 + 3z - 9z^2 + z^3 - 8z^6 + 9z^7 - 9z^8 + z^9
var i, roots,
    p = [-1, 3, -9, 1, 0, 0, -8, 9, -9, 1];

roots = JXG.Math.Numerics.polzeros(p);
for (i = 0; i < roots.length; i++) {
    console.log(i, roots[i].toString());
}
// Output:
0 -0.7424155888401961 + 0.4950476539211721i
1 -0.7424155888401961 + -0.4950476539211721i
2 0.16674869833354108 + 0.2980502714610669i
3 0.16674869833354108 + -0.29805027146106694i
4 0.21429002063640837 + 1.0682775088132996i
5 0.21429002063640842 + -1.0682775088132999i
6 0.861375497926218 + -0.6259177003583295i
7 0.8613754979262181 + 0.6259177003583295i
8 8.000002743888055 + -1.8367099231598242e-40i

<static> {Number} JXG.Math.Numerics.Qag(interval, f, config)
Quadrature algorithm qag from QUADPACK. Internal method used in JXG.Math.Numerics.GaussKronrod15, JXG.Math.Numerics.GaussKronrod21, JXG.Math.Numerics.GaussKronrod31.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} config Optional
The algorithm setup. Accepted propert are max. recursion limit of type number, and epsrel and epsabs, the relative and absolute required precision of type number. Further, q the internal quadrature sub-algorithm of type function.
{Number} config.limit Optional, Default: 15
{Number} config.epsrel Optional, Default: 0.0000001
{Number} config.epsabs Optional, Default: 0.0000001
{Number} config.q Optional, Default: JXG.Math.Numerics.GaussKronrod15
Returns:
{Number} Integral value of f over interval
Examples:
function f(x) {
  return x*x;
}

// calculates integral of f from 0 to 2.
var area1 = JXG.Math.Numerics.Qag([0, 2], f);

// the same with an anonymous function
var area2 = JXG.Math.Numerics.Qag([0, 2], function (x) { return x*x; });

// use JXG.Math.Numerics.GaussKronrod31 rule as sub-algorithm.
var area3 = JXG.Math.Numerics.Quag([0, 2], f,
                                  {q: JXG.Math.Numerics.GaussKronrod31});

<static> {Array} JXG.Math.Numerics.RamerDouglasPeucker(pts, eps)
Implements the Ramer-Douglas-Peucker algorithm. It discards points which are not necessary from the polygonal line defined by the point array pts. The computation is done in screen coordinates. Average runtime is O(nlog(n)), worst case runtime is O(n^2), where n is the number of points.
Parameters:
{Array} pts
Array of JXG.Coords
{Number} eps
If the absolute value of a given number x is smaller than eps it is considered to be equal 0.
Returns:
{Array} An array containing points which represent an apparently identical curve as the points of pts do, but contains fewer points.

<static> JXG.Math.Numerics.RamerDouglasPeuker(pts, eps)
Old name for the implementation of the Ramer-Douglas-Peucker algorithm.
Parameters:
pts
eps
Deprecated:
Use JXG.Math.Numerics.RamerDouglasPeucker

<static> {function} JXG.Math.Numerics.regressionPolynomial(degree, dataX, dataY)
Computes the regression polynomial of a given degree through a given set of coordinates. Returns the regression polynomial function.
Parameters:
{Number|function|Slider} degree
number, function or slider. Either
{Array} dataX
Array containing either the x-coordinates of the data set or both coordinates in an array of JXG.Points or JXG.Coords. In the latter case, the dataY parameter will be ignored.
{Array} dataY
Array containing the y-coordinates of the data set,
Returns:
{function} A function of one parameter which returns the value of the regression polynomial of the given degree. It possesses the method getTerm() which returns the string containing the function term of the polynomial. The function returned will throw an exception, if the data set is malformed.

<static> {Array} JXG.Math.Numerics.riemann(f, n, type, start, end)
Helper function to create curve which displays Riemann sums. Compute coordinates for the rectangles showing the Riemann sum.

In case of type "simpson" and "trapezoidal", the horizontal line approximating the function value is replaced by a parabola or a secant. IN case of "simpson", the parabola is approximated visually by a polygonal chain of fixed step width.

Parameters:
{Function|Array} f
Function or array of two functions. If f is a function the integral of this function is approximated by the Riemann sum. If f is an array consisting of two functions the area between the two functions is filled by the Riemann sum bars.
{Number} n
number of rectangles.
{String} type
Type of approximation. Possible values are: 'left', 'right', 'middle', 'lower', 'upper', 'random', 'simpson', or 'trapezoidal'. "simpson" is Simpson's 1/3 rule.
{Number} start
Left border of the approximation interval
{Number} end
Right border of the approximation interval
Returns:
{Array} An array of two arrays containing the x and y coordinates for the rectangles showing the Riemann sum. This array may be used as parent array of a JXG.Curve. The third parameteris the riemann sum, i.e. the sum of the volumes of all rectangles.

<static> {Number} JXG.Math.Numerics.riemannsum(f, n, type, start, end)
Approximate the integral by Riemann sums. Compute the area described by the riemann sum rectangles. If there is an element of type Riemannsum, then it is more efficient to use the method JXG.Curve.Value() of this element instead.
Parameters:
{Function_Array} f
Function or array of two functions. If f is a function the integral of this function is approximated by the Riemann sum. If f is an array consisting of two functions the area between the two functions is approximated by the Riemann sum.
{Number} n
number of rectangles.
{String} type
Type of approximation. Possible values are: 'left', 'right', 'middle', 'lower', 'upper', 'random', 'simpson' or 'trapezoidal'.
{Number} start
Left border of the approximation interval
{Number} end
Right border of the approximation interval
Returns:
{Number} The sum of the areas of the rectangles.

<static> {Number} JXG.Math.Numerics.Romberg(interval, f, config)
Calculates the integral of function f over interval using Romberg iteration.
Parameters:
{Array} interval
The integration interval, e.g. [0, 3].
{function} f
A function which takes one argument of type number and returns a number.
{Object} config Optional
The algorithm setup. Accepted properties are max_iterations of type number and precision eps.
{Number} config.max_iterations Optional, Default: 20
{Number} config.eps Optional, Default: 0.0000001
Returns:
{Number} Integral value of f over interval
Examples:
function f(x) {
  return x*x;
}

// calculates integral of f from 0 to 2.
var area1 = JXG.Math.Numerics.Romberg([0, 2], f);

// the same with an anonymous function
var area2 = JXG.Math.Numerics.Romberg([0, 2], function (x) { return x*x; });

// use trapez rule with maximum of 16 iterations or stop if the precision 0.0001 has been reached.
var area3 = JXG.Math.Numerics.Romberg([0, 2], f,
                                  {max_iterations: 16, eps: 0.0001});

<static> {Number} JXG.Math.Numerics.root(f, x, context)
Abstract method to find roots of univariate functions, which - for the time being - is an alias for JXG.Math.Numerics.chandrupatla.
Parameters:
{function} f
We search for a solution of f(x)=0.
{Number|Array} x
initial guess for the root, i.e. starting value, or start interval enclosing the root.
{Object} context
optional object that is treated as "this" in the function body. This is useful if the function is a method of an object and contains a reference to its parent object via "this".
Returns:
{Number} A root of the function f.
See:
JXG.Math.Numerics.chandrupatla
JXG.Math.Numerics.fzero

<static> {Array} JXG.Math.Numerics.rungeKutta(butcher, x0, I, N, f)
Solve initial value problems numerically using explicit Runge-Kutta methods. See https://en.wikipedia.org/wiki/Runge-Kutta_methods for more information on the algorithm.
Parameters:
{object|String} butcher
Butcher tableau describing the Runge-Kutta method to use. This can be either a string describing a Runge-Kutta method with a Butcher tableau predefined in JSXGraph like 'euler', 'heun', 'rk4' or an object providing the structure
{
    s: <Number>,
    A: <matrix>,
    b: <Array>,
    c: <Array>
}
which corresponds to the Butcher tableau structure shown here: https://en.wikipedia.org/w/index.php?title=List_of_Runge%E2%80%93Kutta_methods&oldid=357796696 . Default is 'euler'.
{Array} x0
Initial value vector. Even if the problem is one-dimensional, the initial value has to be given in an array.
{Array} I
Interval on which to integrate.
{Number} N
Number of integration intervals, i.e. there are N+1 evaluation points.
{function} f
Function describing the right hand side of the first order ordinary differential equation, i.e. if the ode is given by the equation
dx/dt = f(t, x(t))
. So, f has to take two parameters, a number t and a vector x, and has to return a vector of the same length as x has.
Returns:
{Array} An array of vectors describing the solution of the ode on the given interval I.
Examples:
// A very simple autonomous system dx(t)/dt = x(t);
var f = function(t, x) {
    return [x[0]];
}

// Solve it with initial value x(0) = 1 on the interval [0, 2]
// with 20 evaluation points.
var data = JXG.Math.Numerics.rungeKutta('heun', [1], [0, 2], 20, f);

// Prepare data for plotting the solution of the ode using a curve.
var dataX = [];
var dataY = [];
var h = 0.1;        // (I[1] - I[0])/N  = (2-0)/20
var i;
for(i=0; i<data.length; i++) {
    dataX[i] = i*h;
    dataY[i] = data[i][0];
}
var g = board.create('curve', [dataX, dataY], {strokeWidth:'2px'});

					
					

<static> {Array} JXG.Math.Numerics.splineDef(x, y)
Calculates second derivatives at the given knots.
Parameters:
{Array} x
x values of knots
{Array} y
y values of knots
Returns:
{Array} Second derivatives of the interpolated function at the knots.
See:
#splineEval

<static> {Number|Array} JXG.Math.Numerics.splineEval(x0, x, y, F)
Evaluate points on spline.
Parameters:
{Number|Array} x0
A single float value or an array of values to evaluate
{Array} x
x values of knots
{Array} y
y values of knots
{Array} F
Second derivatives at knots, calculated by JXG.Math.Numerics.splineDef
Returns:
{Number|Array} A single value or an array, depending on what is given as x0.
See:
#splineDef

<static> {Array} JXG.Math.Numerics.Visvalingam(pts, numPoints)
Implements the Visvalingam-Whyatt algorithm. See M. Visvalingam, J. D. Whyatt: "Line generalisation by repeated elimination of the smallest area", C.I.S.R.G Discussion paper 10, July 1992 The algorithm discards points which are not necessary from the polygonal line defined by the point array pts (consisting of type JXG.Coords).
Parameters:
{Array} pts
Array of JXG.Coords
{Number} numPoints
Number of remaining intermediate points. The first and the last point of the original points will be taken in any case.
Returns:
{Array} An array containing points which approximates the curve defined by pts.
Examples:
    var i, p = [];
    for (i = 0; i < 5; ++i) {
        p.push(board.create('point', [Math.random() * 12 - 6, Math.random() * 12 - 6]));
    }

    // Plot a cardinal spline curve
    var splineArr = JXG.Math.Numerics.CardinalSpline(p, 0.5);
    var cu1 = board.create('curve', splineArr, {strokeColor: 'green'});

    var c = board.create('curve', [[0],[0]], {strokeWidth: 2, strokeColor: 'black'});
    c.updateDataArray = function() {
        var i, len, points;

        // Reduce number of intermediate points with Visvakingam-Whyatt to 6
        points = JXG.Math.Numerics.Visvalingam(cu1.points, 6);
        // Plot the remaining points
        len = points.length;
        this.dataX = [];
        this.dataY = [];
        for (i = 0; i < len; i++) {
            this.dataX.push(points[i].usrCoords[1]);
            this.dataY.push(points[i].usrCoords[2]);
        }
    };
    board.update();


					
					

Documentation generated by JsDoc Toolkit 2.4.0 on Fri Mar 08 2024 12:21:01 GMT+0100 (Mitteleuropäische Normalzeit)