{{ message }}

Instantly share code, notes, and snippets.

# mikofski/polyDer2D.m

Last active Mar 13, 2020
horners method in 2-D
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 function [fx,fy] = polyDer2D(p,x,y,n,m) %POLYDER2D Evaluate derivatives of 2-D polynomial using Horner's method. % F = POLYDER2D(P,X,Y,N,M) evaluates the derivatives of 2-D polynomial P at % the points specified by X and Y, which must have the same dimensions. The % outputs FX and FY will have the same dimensions as X and Y. N and M specify % the order of X and Y respectively. Polynomial coefficients are in the % following order. % % F(X,Y) = P_1 * X^N * Y^M + P_2 * X^{N-1} * Y^M + ... + P_{N+1} * Y^M + ... % P_{N+2} * X^N * Y^{M-1} + P_{N+3} * X^{N-1} * Y^{M-1} + ... + P_{2*(N+1)} * Y^{M-1} + ... % ... % P_{M*(N+1)+1} * X^N + P_{M*(N+1)+2} * X^{N-1} + ... + P_{(N+1)*(M+1)} % % See also: POLYFITN by John D'Errico on MathWorks MATLAB Central FEX % http://www.mathworks.com/matlabcentral/fileexchange/34765-polyfitn % % Mark Mikofski, http://poquitopicante.blogspot.com % Version 1-0, 2013-04-17 %% LaTex % % $$f\left(x,y\right)=p_1 x^n y^m+p_2 x^{\left(n-1\right)} y^m+\ldots+p_{n+1} y^m+\ldots$$ % % $$p_{n+2}x^ny^{\left(m-1\right)}+p_{n+3}x^{\left(n-1\right)}y^{\left(m-1\right)}+\ldots+p_{2\left(n+1\right)}y^{\left(m-1\right)}+\ldots$$ % % $$\ldots$$ % % $$p_{m\left(n+1\right)+1}*x^n+p_{m\left(n+1\right)+2}*x^{\left(n-1\right)}+\ldots+p_{\left(n+1\right)\left(m+1\right)}$$ %% check input args validateattributes(p,{'numeric'},{'2d','nonempty','real','finite'}, ... 'polyDer2D','p',1) validateattributes(x,{'numeric'},{'nonempty','real','finite'}, ... 'polyDer2D','x',2) validateattributes(y,{'numeric'},{'nonempty','real','finite'}, ... 'polyDer2D','y',3) assert(all(size(x)==size(y)),'polyDer2D:sizeMismatch', ... 'X and Y must be the same size.') % use size of p to set n & m pdims = size(p); if nargin<4 && all(pdims>1) n = pdims(1)-1; m = pdims(2)-1; end validateattributes(n,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyDer2D','n',4) validateattributes(m,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyDer2D','m',5) if all(pdims>1) assert(pdims(1)==n+1,'polyDer2D:xOrderMismatch', ... 'The number of x coefficients doesn''t match the order n.') assert(pdims(2)==m+1,'polyDer2D:yOrderMismatch', ... 'The number of y coefficients doesn''t match the order m.') end p = p(:); Np = numel(p); assert(Np==(n+1)*(m+1),'OrderMismatch', ... ['The number of x & y coefficients doesn''t match the orders ', ... 'n & m.']) %% fx = df/dx fx = n*p(1); for ni = 1:n-1 fx = fx.*x+(n-ni)*p(1+ni); end for mi = 1:m mj = (n+1)*mi+1; gx = n*p(mj); for ni = 1:n-1 gx = gx.*x+(n-ni)*p(mj+ni); end fx = fx.*y+gx; end %% fy = df/dy fy = p(1); for ni = 1:n fy = fy.*x+p(1+ni); end fy = m*fy; for mi = 1:m-1 mj = (n+1)*mi+1; gy = p(mj); for ni = 1:n gy = gy.*x+p(mj+ni); end fy = fy.*y+(m-mi)*gy; end
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 import numpy as np def polyDer2D(p, x, y, n, m): """ polyDer2D(p, x, y, n, m) Evaluate derivatives of a 2-D polynomial using Horner's method. Evaluates the derivatives of 2-D polynomial p at the points specified by x and y, which must be the same dimensions. The outputs (fx, fy) will have the same dimensions as x and y. The order of x and y are specified by n and m, respectively. Parameters ---------- p : array_like Polynomial coefficients in order specified by polyVal2D.html. x : array_like Values of 1st independent variable. y : array_like Values of 2nd independent variable. n : int Order of the 1st independent variable, x. m : int Order of the 2nd independent variable, y. Returns ------- fx : ndarray Derivative with respect to x. fy : ndarray Derivative with respect to y. See Also -------- numpy.polynomial.polynomial.polyval2d : Evaluate a 2-D polynomial at points (x, y). Example -------- >>> print polyVal2D([1,2,3,4,5,6],2,3,2,1) >>> (39, 11) >>> 1*2*2*3 + 2*3 + 4*2*2 + 5 39 >>> 1*(2**2) + 2*2 + 3 11 >>> print polyVal2D([1,2,3,4,5,6,7,8,9],2,3,2,2) >>> (153, 98) >>> 1*2*2*(3**2) + 2*(3**2) + 4*2*2*3 + 5*3 + 7*2*2 + 8 153 >>> 1*2*(2**2)*3 + 2*2*2*3 + 3*2*3 + 4*(2**2) + 5*2 + 6 98 """ # TODO: check input args p = np.array(p) x = np.array(x) y = np.array(y) n = np.array(n) m = np.array(m) # fx = df/dx fx = n * p[0] for ni in np.arange(n - 1): fx = fx * x + (n - ni - 1) * p[1 + ni] for mi in np.arange(m): mj = (n + 1) * (mi + 1) gx = n * p[mj] for ni in np.arange(n - 1): gx = gx * x + (n - ni - 1) * p[mj + 1 + ni] fx = fx * y + gx # fy = df/dy fy = p[0] for ni in np.arange(n): fy = fy * x + p[1 + ni] fy = m * fy for mi in np.arange(m - 1): mj = (n + 1) * (mi + 1) gy = p[mj] for ni in np.arange(n): gy = gy * x + p[mj + 1 + ni] fy = fy * y + (m - mi - 1) * gy return fx, fy
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 function p = polyFit2D(f,x,y,n,m) %POLYFIT2D Fit data with a 2-D polynomial. % P = POLYFIT2D(F,X,Y,N,M) determines the 2-D polynomial P that best fit the % data F(X,Y) in the least squares sense. Polynomial coefficients are in the % following order. % % F(X,Y) = P_1 * X^N * Y^M + P_2 * X^{N-1} * Y^M + ... + P_{N+1} * Y^M + ... % P_{N+2} * X^N * Y^{M-1} + P_{N+3} * X^{N-1} * Y^{M-1} + ... + P_{2*(N+1)} * Y^{M-1} + ... % ... % P_{M*(N+1)+1} * X^N + P_{M*(N+1)+2} * X^{N-1} + ... + P_{(N+1)*(M+1)} % % See also: POLYFITN by John D'Errico on MathWorks MATLAB Central FEX % http://www.mathworks.com/matlabcentral/fileexchange/34765-polyfitn % % Mark Mikofski, http://poquitopicante.blogspot.com % Version 1-0, 2013-03-13 %% LaTex % % $$f\left(x,y\right)=p_1 x^n y^m+p_2 x^{\left(n-1\right)} y^m+\ldots+p_{n+1} y^m+\ldots$$ % % $$p_{n+2}x^ny^{\left(m-1\right)}+p_{n+3}x^{\left(n-1\right)}y^{\left(m-1\right)}+\ldots+p_{2\left(n+1\right)}y^{\left(m-1\right)}+\ldots$$ % % $$\ldots$$ % % $$p_{m\left(n+1\right)+1}*x^n+p_{m\left(n+1\right)+2}*x^{\left(n-1\right)}+\ldots+p_{\left(n+1\right)\left(m+1\right)}$$ %% check input args validateattributes(f,{'numeric'},{'nonempty','real','finite'}, ... 'polyFit2D','f',1) validateattributes(x,{'numeric'},{'nonempty','real','finite'}, ... 'polyFit2D','x',2) validateattributes(y,{'numeric'},{'nonempty','real','finite'}, ... 'polyFit2D','y',3) assert(all(size(x)==size(y)),'polyFit2D:sizeMismatch', ... 'X and Y must be the same size.') assert(all(size(x)==size(f)),'polyFit2D:sizeMismatch', ... 'X and F must be the same size.') validateattributes(n,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyFit2D','n',4) validateattributes(m,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyFit2D','m',5) Npoints = numel(f); assert(all((n+1)*(m+1)
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 \begin{aligned} f\left(x,y\right)=&p_1 x^n y^m+p_2 x^{\left(n-1\right)} y^m+\ldots+p_{n+1} y^m+\ldots \\ &p_{n+2}x^ny^{\left(m-1\right)}+p_{n+3}x^{\left(n-1\right)}y^{\left(m-1\right)}+\ldots+p_{2\left(n+1\right)}y^{\left(m-1\right)}+\ldots \\ &\ldots \\ &p_{m\left(n+1\right)+1}*x^n+p_{m\left(n+1\right)+2}*x^{\left(n-1\right)}+\ldots+p_{\left(n+1\right)\left(m+1\right)} \end{aligned}
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 function f = polyVal2D(p,x,y,n,m) %POLYVAL2D Evaluate a 2-D polynomial using Horner's method. % F = POLYVAL2D(P,X,Y) evaluates the 2-D polynomial P at the points % specified by X and Y, which must have the same dimensions. The output F % will have the same dimensions as X and Y. N and M specify the order of X % and Y respectively. Polynomial coefficients are in the following order. % % F(X,Y) = P_1 * X^N * Y^M + P_2 * X^{N-1} * Y^M + ... + P_{N+1} * Y^M + ... % P_{N+2} * X^N * Y^{M-1} + P_{N+3} * X^{N-1} * Y^{M-1} + ... + P_{2*(N+1)} * Y^{M-1} + ... % ... % P_{M*(N+1)+1} * X^N + P_{M*(N+1)+2} * X^{N-1} + ... + P_{(N+1)*(M+1)} % % % See also: POLYFITN by John D'Errico on MathWorks MATLAB Central FEX % http://www.mathworks.com/matlabcentral/fileexchange/34765-polyfitn % % Mark Mikofski, http://poquitopicante.blogspot.com % Version 1-0, 2013-03-13 %% LaTex % % $$f\left(x,y\right)=p_1 x^n y^m+p_2 x^{\left(n-1\right)} y^m+\ldots+p_{n+1} y^m+\ldots$$ % % $$p_{n+2}x^ny^{\left(m-1\right)}+p_{n+3}x^{\left(n-1\right)}y^{\left(m-1\right)}+\ldots+p_{2\left(n+1\right)}y^{\left(m-1\right)}+\ldots$$ % % $$\ldots$$ % % $$p_{m\left(n+1\right)+1}*x^n+p_{m\left(n+1\right)+2}*x^{\left(n-1\right)}+\ldots+p_{\left(n+1\right)\left(m+1\right)}$$ %% check input args validateattributes(p,{'numeric'},{'2d','nonempty','real','finite'}, ... 'polyVal2D','p',1) validateattributes(x,{'numeric'},{'nonempty','real','finite'}, ... 'polyVal2D','x',2) validateattributes(y,{'numeric'},{'nonempty','real','finite'}, ... 'polyVal2D','y',3) assert(all(size(x)==size(y)),'polyVal2D:sizeMismatch', ... 'X and Y must be the same size.') % use size of p to set n & m pdims = size(p); if nargin<4 && all(pdims>1) n = pdims(1)-1; m = pdims(2)-1; end validateattributes(n,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyVal2D','n',4) validateattributes(m,{'numeric'},{'scalar','integer','positive','<',10}, ... 'polyVal2D','m',5) if all(pdims>1) assert(pdims(1)==n+1,'polyVal2D:xOrderMismatch', ... 'The number of x coefficients doesn''t match the order n.') assert(pdims(2)==m+1,'polyVal2D:yOrderMismatch', ... 'The number of y coefficients doesn''t match the order m.') end %% evaluate polynomial P f = p(1); for ni = 1:n f = f.*x+p(1+ni); end for mi = 1:m mj = (n+1)*mi+1; g = p(mj); for ni = 1:n g = g.*x+p(mj+ni); end f = f.*y+g; end
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 import numpy as np def polyVal2D(p, x, y, n, m): """ polyVal2D(p, x, y, n, m) Evaluate a 2-D polynomial using Horner's method. Evaluates the 2-D polynomial p at the points specified by x and y, which must be the same dimensions. The output f will have the same dimensions as x and y. The order of x and y are specified by n and m, respectively. Parameters ---------- p : array_like Polynomial coefficients in order specified by polyVal2D.html. x : array_like Values of 1st independent variable. y : array_like Values of 2nd independent variable. n : int Order of the 1st independent variable, x. m : int Order of the 2nd independent variable, y. Returns ------- f : ndarray Values of the evaluated 2-D polynomial. See Also -------- numpy.polynomial.polynomial.polyval2d : Evaluate a 2-D polynomial at points (x, y). Example -------- >>> (5**2 * 6**2 + 2 * 5 * 6**2 + 3 * 6**2 + 4 * 5**2 * 6 + ... 5 * 5 * 6 + 6 * 6 + 7 * 5**2 + 8 * 5 + 9) 2378 >>> polyVal2D([1,2,3,4,5,6,7,8,9],5,6,2,2) 2378 >>> 5 * 6 + 2 * 6 + 3 * 5 + 4 61 >>> polyVal2D([1,2,3,4],5,6,1,1) 61 """ # TODO: check input args p = np.array(p) x = np.array(x) y = np.array(y) n = np.array(n) m = np.array(m) # evaluate f(x,y) f = p[0] for ni in np.arange(n): f = f * x + p[1 + ni] for mi in np.arange(m): mj = (n + 1) * (mi + 1) g = p[mj] for ni in np.arange(n): g = g * x + p[mj + 1 + ni] f = f * y + g return f