%%% bpolynomial.mp %%% Copyright 2007 Stephan Hennig % % This work may be distributed and/or modified under the conditions of % the LaTeX Project Public License, either version 1.3 of this license % or (at your option) any later version. The latest version of this % license is in http://www.latex-project.org/lppl.txt % %%% Identify yourself. if known bpolynomial_fileversion: endinput fi; string bpolynomial_fileversion; bpolynomial_fileversion := "v0.5 (2007/12/12)"; message "Loading bpolynomial " & bpolynomial_fileversion; %%% Main user macro for defining polynomials. %%% Arguments are a suffix and the coefficients %%% of the function a*x^3 + b*x^2 + c*x + d. vardef newBPolynomial@#(expr a, b, c, d)= bpolynomial__defineBPolynomial.@#(a, b, c, d); bpolynomial__defineBPolynomial.@#'(0, 3a, 2b, c); bpolynomial__defineBPolynomial.@#''(0, 0, 6a, 2b); bpolynomial__defineBPolynomial.@#'''(0, 0, 0, 6a); enddef; %%% This macro returns the path of a Bezier curve that matches %%% a function a*x^3 + b*x^2 + c*x + d between two points A and D. %%% This macro is the heart of this package and is used by %%% several other macros. %%% Arguments are the coefficients of the polynomial and the %%% start and end point of the graph/path. vardef bpolynomial__getBezierFromPolynomial(expr a, b, c, d, A, D)= save xA,xB,xC,xD,yA,yB,yC,yD; save xl,yl,xr,yr,dx; numeric xA,xB,xC,xD,yA,yB,yC,yD; numeric xl,yl,xr,yr,dx; xl := xpart A; yl := ypart A; xr := xpart D; yr := ypart D; dx := xpart D - xpart A; %%% Original equation system for x values. % xA = xl; % 3(xB - xA) = dx; % 3(xC - 2xB + xA) = 0; % xD - 3xC + 3xB - xA = 0; %%% Modified equation system. xA := xl; xB := xl + dx/3; xC := xr - dx/3; xD := xr; %%% Original equation system for y values. % yA = ((a*xl + b)*xl + c)*xl + d; % 3(yB - yA) = dx*((3a*xl + 2b)*xl + c); % 3(yC - 2yB + yA) = dx*dx*(3a*xl + b); % yD - 3yC + 3yB - yA = a*dx*dx*dx; %%% Modified equation system. yA := yl; 3(yB - yA) = dx*((3a*xl + 2b)*xl + c); 3(yC - 2yB + yA) = dx*dx*(3a*xl + b); yD := yr; %%% Return path A..controls B and C..D. (xA,yA)..controls (xB,yB) and (xC,yC)..(xD,yD) enddef; %%% This macro returns the path of a Bezier curve that matches %%% a function a*x^3 + b*x^2 + c*x + d in the range [xl, xr]. %%% Arguments are the coefficients of the polynomial and the %%% range boundaries of the graph/path. vardef getBezierFromPolynomial(expr a, b, c, d, xl, xr)= bpolynomial__getBezierFromPolynomial(a, b, c, d, (xl, ((a*xl+b)*xl+c)*xl+d), (xr, ((a*xr+b)*xr+c)*xr+d)) enddef; %%% This macro returns the path of a Bezier curve that matches %%% a function a*x^3 + b*x^2 + c*x + d in the range [xl, xr]. %%% Arguments are the coefficients of the polynomial and the %%% range boundaries of the graph/path. vardef getBezierFromSqrRoot(expr u, v, w, xl, xr)= save yl, yr; numeric yl,yr; if (xl >= -v): yl := xl; else: message "Package bpolynomial warning: Replacing lower range boundary " & decimal xl & " by " & decimal -v & "!"; yl := -v; fi if (xr >= -v): yr := xr; else: message "Package bpolynomial warning: Replacing upper range boundary " & decimal xr & " by " & decimal -v & "!"; yr := -v; fi bpolynomial__getBezierFromPolynomial(0, 1/u/u, -2*w/u/u, (w/u)*(w/u)-v, (u*sqrt(yl+v)+w, yl), (u*sqrt(yr+v)+w, yr)) reflectedabout ((0,0),(1,1)) enddef; %%% This macro returns the path of a Bezier curve that matches %%% a function a*x^3 + b*x^2 + c*x + d in the range [xl, xr]. %%% Arguments are the coefficients of the polynomial and the %%% range boundaries of the graph/path. vardef getBezierFromCubRoot(expr u, v, w, xl, xl)= save yl, yr; numeric yl,yr; if (xl >= -v): yl := xl; else: message "Package bpolynomial warning: Replacing lower range boundary " & decimal xl & " by " & decimal -v & "!"; yl := -v; fi if (xr >= -v): yr := xr; else: message "Package bpolynomial warning: Replacing upper range boundary " & decimal xr & " by " & decimal -v & "!"; yr := -v & "!"; fi bpolynomial__getBezierFromPolynomial(1/u/u/u, -3w/u/u/u, 3(w/u)*(w/u)/u, (w/u)*(w/u)*(w/u)-v, (u*((yl+v)**(1/3))+w, yl), (u*((yr+v)**(1/3))+w, yr)) reflectedabout ((0,0),(1,1)) enddef; %%% This internal macro defines a new polynomial. %%% Arguments are a suffix macro and the coefficients %%% of the polynomial a*x^3 + b*x^2 + c*x + d. vardef bpolynomial__defineBPolynomial@#(expr ca,cb,cc,cd)= numeric @#.a, @#.b, @#.c, @#.d; %%% Save coefficients for later access. %%% Variable @#.a refers to coefficient a of polynomial @#. @#.a := ca; @#.b := cb; @#.c := cc; @#.d := cd; %%% This macro returns values of polynomial @#. %%% Argument is an x value. vardef @#.eval(expr x)= (((@#.a*x + @#.b)*x + @#.c)*x + @#.d) enddef; %%% This macro returns the path corresponding to polynomial @# %%% on the intervall [xl, xr]. vardef @#.getPath(expr xl,xr)= bpolynomial__getBezierFromPolynomial(@#.a, @#.b, @#.c, @#.d, (xl, @#.eval(xl)), (xr, @#.eval(xr))) enddef; %%% This macro returns a path tangent to @# at point (x, f(x)) %%% covering the interval [x+xm, x+xp]. vardef @#.getTangent(expr x, xm, xp)= save m, y; numeric m, y; m := (3@#.a*x + 2@#.b)*x + @#.c; y := @#.eval(x); (x+xm, y + m*xm) -- (x+xp, y + m*xp) enddef; enddef; %%% This macro defines a new square root. %%% Arguments are a suffix macro and the parameters %%% of the function u*(x + v)^(1/2) + w. vardef newBSqrRoot@#(expr cu,cv,cw)= numeric @#.a, @#.b, @#.c, @#.d; numeric @#.u, @#.v, @#.w; %%% Save parameters for later access. %%% Variable @#.v refers to parameters of square root @#. %%% Variables @#.a to @#.d store the coefficients of the %%% corresponding polynomial. @#.u := cu; @#.v := cv; @#.w := cw; @#.a := 0; @#.b := 1/cu/cu; @#.c := -2*cw/cu/cu; @#.d := (cw/cu)*(cw/cu)-cv; %%% This macro returns values of polynomial @#. %%% Argument is an x value. vardef @#.eval(expr x)= if (x >= -@#.v): @#.u*sqrt(x + @#.v) + @#.w else: message "Package bpolynomial warning: Cannot evaluate function at x = " & decimal x & "!"; @#.w fi enddef; %%% This macro returns the path corresponding to square root @# %%% on the intervall [yl, yr]. The path of the corresponing %%% polynomial is computed and then transformed. vardef @#.getPath(expr xl,xr)= save yl, yr; numeric yl, yr; if (xl >= -@#.v): yl := xl; else: message "Package bpolynomial warning: Replacing lower range boundary " & decimal xl & " by " & decimal -@#.v & "!"; yl := -@#.v; fi if (xr >= -@#.v): yr := xr; else: message "Package bpolynomial warning: Replacing upper range boundary " & decimal xr & " by " & decimal -@#.v & "!"; yr := -@#.v; fi bpolynomial__getBezierFromPolynomial(@#.a, @#.b, @#.c, @#.d, (@#.eval(yl), yl), (@#.eval(yr), yr)) reflectedabout ((0,0),(1,1)) enddef; %%% This macro returns a path tangent to square root @# %%% at point (x, f(x)) covering the interval [x+xm, x+xp]. vardef @#.getTangent(expr x, epsl, epsr)= save m, y; numeric m, y; if (x >= -@#.v): m := @#.u/(2sqrt(x + @#.v)); y := @#.eval(x); (x+epsl, y + m*epsl) -- (x+epsr, y + m*epsr) else: message "Package bpolynomial warning: Cannot draw tangent at x = " & decimal x & "!"; (-@#.v, @#.w)--(-@#.v, @#.w+1) fi enddef; enddef; %%% This macro defines a new cubic root. %%% Arguments are a suffix macro and the parameters %%% of the function u*(x + v)^(1/3) + w. vardef newBCubRoot@#(expr cu,cv,cw)= numeric @#.a, @#.b, @#.c, @#.d; numeric @#.u, @#.v, @#.w; %%% Save parameters for later access. %%% Variable @#.v refers to parameters of cubic root @#. %%% Variables @#.a to @#.d store the coefficients of the %%% corresponding polynomial. @#.u := cu; @#.v := cv; @#.w := cw; @#.a := 1/cu/cu/cu; @#.b := -3cw/cu/cu/cu; @#.c := 3(cw/cu)*(cw/cu)/cu; @#.d := (cw/cu)*(cw/cu)*(cw/cu)-cv; %%% This macro returns values of polynomial @#. %%% Argument is an x value. vardef @#.eval(expr x)= if (x >= -@#.v): @#.u*((x+@#.v)**(1/3)) + @#.w else: message "Package bpolynomial warning: Cannot evaluate function at x = " & decimal x & "!"; @#.w fi enddef; %%% This macro returns the path corresponding to cubic root @# %%% on the intervall [yl, yr]. The path of the corresponing %%% polynomial is computed and then transformed. vardef @#.getPath(expr xl,xr)= save yl, yr; numeric yl, yr; if (xl >= -@#.v): yl := xl; else: message "Package bpolynomial warning: Replacing lower range boundary " & decimal xl & " by " & decimal -@#.v & "!"; yl := -@#.v; fi if (xr >= -@#.v): yr := xr; else: message "Package bpolynomial warning: Replacing upper range boundary " & decimal xr & " by " & decimal -@#.v & "!"; yr := -@#.v; fi bpolynomial__getBezierFromPolynomial(@#.a, @#.b, @#.c, @#.d, (@#.eval(yl), yl), (@#.eval(yr), yr)) reflectedabout ((0,0),(1,1)) enddef; %%% This macro returns a path tangent to cubic root @# %%% at point (x, f(x)) covering the interval [x+xm, x+xp]. vardef @#.getTangent(expr x, epsl, epsr)= save m, y; numeric m, y; if (x >= -@#.v): m := @#.u/3/((x + @#.v)**(2/3)); y := @#.eval(x); (x+epsl, y + m*epsl) -- (x+epsr, y + m*epsr) else: message "Package bpolynomial warning: Cannot draw tangent at x = " & decimal x & "!"; (-@#.v, @#.w)--(-@#.v, @#.w+1) fi enddef; enddef;