www.gusucode.com > symbolic工具箱matlab源码程序 > symbolic/@sym/bernstein.m
function b = bernstein(f, varargin) % BERNSTEIN(f,n,x) with a nonnegative integer n and a % function handle f is the n-th order Bernstein % polynomial approximation % symsum(nchoosek(n,k)*x^k*(1-x)^(n-k)*f(k/n),k,0,n), % of f, evaluated at the point x. % % BERNSTEIN(f,n,x) with a nonnegative integer n and a % symbolic expression or symfun f is the n-th order % Bernstein polynomial approximation of f at the % point x. The expression f is regarded as a univariate % function of the unknown determined via symvar(f,1). % % BERNSTEIN(f,var,n,x) with a nonnegative integer n and % a symbolic expression or symfun f is the n-th order % Bernstein polynomial approximation of f at the % point x. The expression f is regarded as a univariate % function of the unknown var. % % Functions f must accept one scalar input argument and % return a scalar value. % % The symbolic polynomials returned for symbolic x are % numerically stable when substituting numerical values % between 0 and 1 for x. % % Examples: % The second order approximation of the exponential % function by Bernstein polynomials in the symbolic % variable x is given by: % % syms x; % bernstein(exp(x), 2, x) % % This returns: % % (x - 1)^2 + x^2*exp(1) - 2*x*exp(1/2)*(x - 1) % % The following is the 5th order Bernstein approximation % of a function f representing a linear ramp: % % syms x; % f(x) = triangularPulse(1/4, 3/4, Inf, x); % p = bernstein(f, 5, x) % % This returns: % % 7*x^3*(x-1)^2 - 3*x^2*(x-1)^3 - 5*x^4*(x-1) + x^5. % % The call simplify(p) simplifies this to: % % -x^2*(2*x - 3). % % Note that simplification of a symbolic Bernstein polynomial % may produce a representation that cannot be evaluated in a % numerically stable way. Compare the following plot of a % rectangular pulse function, its numerically stable Bernstein % representation f1, and its simplified version f2 that is not % numerically stable: % % f = @(x)rectangularPulse(1/4,3/4,x); % b1 = bernstein(f,100,sym('x')) % b2 = simplify(b1) % f1 = matlabFunction(b1); % f2 = matlabFunction(b2); % x = 0:0.001:1; % plot(x,f(x),x,f1(x),x,f2(x)) % % See also BERNSTEIN, BERNSTEINMATRIX, SYM/NCHOOSEK, % SYM/SYMVAR, SYM/SYMSUM. % Copyright 2013 The MathWorks, Inc. if nargin == 3 n = varargin{1}; x = varargin{2}; if isa(f, 'sym') var = symvar(f, 1); elseif isa(f, 'function_handle') var = sym([]); % we only need a dummy else var = symvar(sym(f),1); end elseif nargin == 4 var = varargin{1}; n = varargin{2}; x = varargin{3}; if isa(f, 'function_handle') error(message('symbolic:sym:bernstein:ExpectingThreeArguments')); end else error(message('symbolic:sym:bernstein:ExpectingThreeOrFourArguments')); end % check argument var if ~isa(var,'sym'), var = sym(var); end if builtin('numel',var) ~= 1, var = normalizesym(var); end if nargin == 4 && ~isequal(var, symvar(var, 1)) error(message('symbolic:sym:bernstein:ExpectingUnknown2')); end % check argument n n = formula(sym(n)); if builtin('numel',n) ~= 1, n = normalizesym(n); end n = feval(symengine, 'specfunc::makeInteger', n); if ~(isscalar(n) && isreal(n) && n == round(n) && n >= 0) || ~isfinite(n) if nargin == 3 error(message('symbolic:sym:bernstein:ExpectingNonNegativeInteger2')); else error(message('symbolic:sym:bernstein:ExpectingNonNegativeInteger3')); end end % check argument x if ~isa(x,'sym'), x = sym(x); end if builtin('numel',x) ~= 1, x = normalizesym(x); end x = formula(x); if ~isscalar(x) if nargin == 3 error(message('symbolic:sym:bernstein:ExpectingScalar3')); else error(message('symbolic:sym:bernstein:ExpectingScalar4')); end end % check argument f if isa(f, 'function_handle') % --------------------------------------------------------- % Create matrix input [f(0/n), f(1/n), ...] from function f % --------------------------------------------------------- nf = nargin(f); if ~(nf == 1 || nf == -1 || nf == -2) error(message('symbolic:sym:bernstein:ExpectingUnivariateFunctionHandle')); end if n == 0 vector = 0; else vector = (0:1/n:1)'; end f = arrayfun(f, sym(vector), 'UniformOutput', false); f = sym([f{:}]); if numel(size(f)) ~=2 || any(size(f) ~= [1 n+1]) error(message('symbolic:sym:bernstein:FunctionMustProduceScalars')); end if any(~isfinite(f)) error(message('symbolic:sym:bernstein:FunctionMustProduceFiniteValues')); end elseif ~isscalar(f) error(message('symbolic:sym:bernstein:ExpectingScalar1')); end if ~isa(f,'sym') f = sym(f); end if builtin('numel',f) ~= 1 f = normalizesym(f); end %---------------- % replace symfuns %---------------- issymfun = isa(f, 'symfun'); if issymfun fargs = argnames(f); % Eliminate var from the arguments of the symfun fargs = privsubsref(fargs, logical(fargs ~= var)); f = formula(f); if ~isscalar(f) error(message('symbolic:sym:bernstein:ExpectingScalar1')); end end if isempty(var) b = mupadmex('bernstein', f.s, n.s, x.s); else b = mupadmex('bernstein', f.s, var.s, n.s, x.s); end % ------------------------- % re-vitalize symfuns % ------------------------- if issymfun && ~isempty(fargs) b = symfun(b, fargs); end