1 ## Copyright (C) 2003 Willem J. Atsma <watsma@users.sf.net>
3 ## This program is free software; you can redistribute it and/or modify
4 ## it under the terms of the GNU General Public License as published by
5 ## the Free Software Foundation; either version 2 of the License, or
6 ## (at your option) any later version.
8 ## This program is distributed in the hope that it will be useful,
9 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
10 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 ## GNU General Public License for more details.
13 ## You should have received a copy of the GNU General Public License
14 ## along with this program; If not, see <http://www.gnu.org/licenses/>.
17 ## @deftypefn {Function File} {[@var{x}, @var{inf}, @var{msg}] =} symfsolve (@dots{})
18 ## Solve a set of symbolic equations using @command{fsolve}. There are a number of
19 ## ways in which this function can be called.
21 ## This solves for all free variables, initial values set to 0:
25 ## x=sym("x"); y=sym("y");
26 ## f=x^2+3*x-1; g=x*y-y^2+3;
27 ## a = symfsolve(f,g);
30 ## This solves for x and y and sets the initial values to 1 and 5 respectively:
33 ## a = symfsolve(f,g,x,1,y,5);
34 ## a = symfsolve(f,g,@{x==1,y==5@});
35 ## a = symfsolve(f,g,[1 5]);
38 ## In all the previous examples vector a holds the results: x=a(1), y=a(2).
39 ## If initial conditions are specified with variables, the latter determine
43 ## a = symfsolve(f,g,@{y==1,x==2@}); # here y=a(1), x=a(2)
46 ## The system of equations to solve for can be given as separate arguments or
47 ## as a single cell-array:
50 ## a = symfsolve(@{f,g@},@{y==1,x==2@}); # here y=a(1), x=a(2)
53 ## If the variables are not specified explicitly with the initial conditions,
54 ## they are placed in alphabetic order. The system of equations can be comma-
55 ## separated or given in a cell-array. The return-values are those of
56 ## fsolve; @var{x} holds the found roots.
60 function [ x,inf,msg ] = symfsolve (varargin)
62 ## separate variables and equations
66 if iscell(varargin{1})
67 if !strcmp(typeinfo(varargin{1}{1}),"ex")
68 error("First argument must be (a cell-array of) symbolic expressions.")
75 tmp = disp(varargin{i});
76 if( iscell(varargin{i}) || ...
77 all(isalnum(tmp) || tmp=="_" || tmp==",") || ...
78 !strcmp(typeinfo(varargin{i}),"ex") )
81 eqns{end+1} = varargin{i};
82 arg_count = arg_count+1;
87 error("No equations specified.")
90 ## make a list with all variables from equations
95 evars = findsymbols(tmp);
98 ## After the equations may follow initial values. The formats are:
100 ## x,0,y,0.3,z,-3,...
101 ## {x==0, y==0.3, z==-3 ...}
102 ## none - default of al zero initial values
108 elseif (nargin-arg_count)>1
109 if mod(nargin-arg_count,2)
110 error("Initial value symbol-value pairs don't match up.")
112 for i=(arg_count+1):2:nargin
113 tmp = disp(varargin{i});
114 if all(isalnum(tmp) | tmp=="_" | tmp==",")
115 vars{end+1} = varargin{i};
116 X0((i-arg_count+1)/2)=varargin{i+1};
118 error("Error in symbol-value pair arguments.")
121 nvars = length(vars);
123 nvars = length(varargin{arg_count+1});
125 error("The number of initial conditions does not match the number of free variables.")
127 if iscell(varargin{arg_count+1})
128 ## cell-array of relations - this should work for a list of strings ("x==3") too.
130 tmp = disp(varargin{arg_count+1}{i});
131 vars{end+1} = {sym(strtok(tmp,"=="))};
132 X0(i) = str2num(tmp((findstr(tmp,"==")+2):length(tmp)));
135 ## straight numbers, match up with symbols in alphabetic order
137 X0 = varargin{arg_count+1};
141 ## X0 is now a vector, vars a list of variables.
142 ## create temporary function:
143 symfn = sprintf("function Y=symfn(X) ");
145 symfn = [symfn sprintf("%s=X(%d); ",disp(vars{i}),i)];
148 symfn = [symfn sprintf("Y(%d)=%s; ",i,disp(eqns{i}))];
150 symfn = [symfn sprintf("endfunction")];
153 [x,inf,msg] = fsolve("symfn",X0);
160 % f = x ^ 2 + 3 * x - 1;
161 % g = x * y - y ^ 2 + 3;
163 % assert (symfsolve (f, g), [0.30278; -1.58727]', 1e-5);
165 % assert (symfsolve (f, g, x, 1, y, 5), [0.30278; 1.89004]', 1e-5);
167 % assert (symfsolve (f, g, {x==1,y==5}), [0.30278; 1.89004]', 1e-5);
169 % assert (symfsolve (f, g, [1 5]), [0.30278; 1.89004]', 1e-5);
171 % assert (symfsolve ({f, g}, {y==1,x==2}), [1.89004; 0.30278]', 1e-5);