solve - 符号方程解

2支棒棒糖 · 发表于 2018-9-20 14:50:06

function varargout = solve(varargin)
%求各种类型方程组的解析解,
%当找不到解析解时,将自动寻求原点附近的一个近似解, 并可达任意精度.
%用法
% solve('方程1','方程2',...,'方程N','变量1','变量2',...,'变量M')
%SOLVE Symbolic solution of algebraic equations.
% SOLVE('eqn1','eqn2',...,'eqnN')
% SOLVE('eqn1','eqn2',...,'eqnN','var1,var2,...,varN')
% SOLVE('eqn1','eqn2',...,'eqnN','var1','var2',...'varN')
%
% The eqns are symbolic expressions or strings specifying equations. The
% vars are symbolic variables or strings specifying the unknown variables.
% SOLVE seeks zeros of the expressions or solutions of the equations.
% If not specified, the unknowns in the system are determined by FINDSYM.
% If no analytical solution is found and the number of equations equals
% the number of dependent variables, a numeric solution is attempted.
%
% Three different types of output are possible. For one equation and one
% output, the resulting solution is returned, with multiple solutions to
% a nonlinear equation in a symbolic vector. For several equations and
% an equal number of outputs, the results are sorted in lexicographic
% order and assigned to the outputs. For several equations and a single
% output, a structure containing the solutions is returned.
%
% Examples:
%
% solve('p*sin(x) = r') chooses 'x' as the unknown and returns
%
% ans =
% asin(r/p)
%
% [x,y] = solve('x^2 + x*y + y = 3','x^2 - 4*x + 3 = 0') returns
%
% x =
% [ 1]
% [ 3]
%
% y =
% [ 1]
% [ -3/2]
%
% S = solve('x^2*y^2 - 2*x - 1 = 0','x^2 - y^2 - 1 = 0') returns
% the solutions in a structure.
%
% S =
% x: [8x1 sym]
% y: [8x1 sym]
%
% [u,v] = solve('a*u^2 + v^2 = 0','u - v = 1') regards 'a' as a
% parameter and solves the two equations for u and v.
%
% S = solve('a*u^2 + v^2','u - v = 1','a,u') regards 'v' as a
% parameter, solves the two equations, and returns S.a and S.u.
%
% [a,u,v] = solve('a*u^2 + v^2','u - v = 1','a^2 - 5*a + 6') solves
% the three equations for a, u and v.
%
% [x,y] = solve('sin(x+y)-exp(x)*y = 0','x^2-y = 2') cannot find
% an analytic solution, so returns a numeric solution.
%
% See also DSOLVE.
% $Revision: 1.13 $ $Date: 1997/11/29 01:06:35 $
% Set the Explicit solution (for equations of order 4 or less)
% option on in the Maple Workspace.
maple('_EnvExplicit := true;');
% Collect input arguments together in either equation or variable lists.
eqns = [];
vars = [];
for k = 1:nargin
v = varargin{k};
if all(isletter(v) | ('0'<=v & v<='9') | v == '_' | v == ',')
vars = [vars ',' v];
else
[t,stat] = maple(v);
if stat
error(['''' v ''' is not a valid expression or equation.'])
end
eqns = [eqns ',' v];
end
end
if isempty(eqns)
warning('List of equations is empty.')
varargout = cell(1,nargout);
varargout{1} = sym([]);
return
else
eqns(1) = '{'; eqns(end+1) = '}';
end
neqns = sum(commas(eqns)) + 1;
% Determine default variables and sort variable list.
if isempty(vars)
vars = ['[' findsym(sym(eqns),neqns) ']'];
else
vars(1) = '['; vars(end+1) = ']';
end
v = vars;
[vars,stat] = maple('sort',v);
if stat
error(['''' v ''' is not a valid variable list.'])
end
vars(1) = '{'; vars(end) = '}';
nvars = sum(commas(vars)) + 1;
if neqns ~= nvars
warning([int2str(neqns) ' equations in ' int2str(nvars) ' variables.'])
end
% Seek analytic solution.
[R,stat] = maple('solve',eqns,vars);
% If no analytic solution, seek numerical solution.
if (isempty(R) & (nvars == neqns))
[R,stat] = maple('fsolve',eqns,vars);
end
if stat
error(R)
end
% If still no solution, give up.
if isempty(R) | findstr(R,'fsolve')
warning('Explicit solution could not be found.');
varargout = cell(1,nargout);
varargout{1} = sym([]);
return
end
% Expand any RootOf.
while ~isempty(findstr(R,'RootOf'))
k = min(findstr(R,'RootOf'));
p = findstr(R,'{'); p = max(p(p<k));
q = findstr(R,'}'); q = min(q(q>k));
s = R(p:q);
t = s(min(findstr(s,'RootOf'))+6:end);
e = min(find(cumsum((t=='(')-(t==')'))==0));
if isempty(find(commas(t(2:e-1))))
% RootOf with one argument, possibly an imbedded RootOf.
[s,stat] = maple('allvalues',s,'dependent');
if stat
error(s)
end
else
% RootOf with a good estimate as the second argument.
s = maple('evalf',s);
end
R = [R(1:p-1) s R(q+1:end)];
end
% Parse the result.
if nvars == 1 & nargout <= 1
% One variable and at most one output.
% Return a single scalar or vector sym.
S = sym([]);
c = find(commas(R) | R == '}');
for p = find(R == '=')
q = min(c(c>p));
t = trim(R(p+1:q-1)); % The solution (xk)
S = [S; sym(t)];
end
varargout{1} = S;
else
% Several variables.
% Create a skeleton structure.
c = [1 find(commas(vars)) length(vars)];
S = [];
for j = 1:nvars
v = trim(vars(c(j)+1:c(j+1)-1));
S = setfield(S,v,[]);
end
% Complete the structure.
c = [1 find(commas(R) | R == '{' | R == '}') length(R)];
for p = find(R == '=')
q = max(c(c<p));
v = trim(R(q+1:p-1)); % The variable (x)
q = min(c(c>p));
t = trim(R(p+1:q-1)); % The solution (xk)
S = setfield(S,v,[getfield(S,v); sym(t)]);
end
if nargout <= 1
% At most one output, return the structure.
varargout{1} = S;
elseif nargout == nvars
% Same number of outputs as variables.
% Match results in lexicographic order to outputs.
v = fieldnames(S);
for j = 1:nvars
varargout{j} = getfield(S,v{j});
end
else
error([int2str(nvars) ' variables does not match ' ...
int2str(nargout) ' outputs.'])
end
end
%-------------------------
function s = trim(s);
%TRIM TRIM(s) deletes any leading or trailing blanks.
while s(1) == ' ', s(1) = []; end
while s(end) == ' ', s(end) = []; end
%-------------------------
function c = commas(s)
%COMMAS COMMAS(s) is true for commas not inside parentheses.
p = cumsum((s == '(') - (s == ')'));
c = (s == ',') & (p == 0);

复制代码

		自动登录	找回密码
密码			立即注册

solve - 符号方程解

站长推荐 /1