annotate doc/interpreter/nonlin.txi @ 7018:fd42779a8428

[project @ 2007-10-13 00:52:12 by jwe]
author jwe
date Sat, 13 Oct 2007 00:52:13 +0000
parents 9398f6a81bdf
children 8fb8d6985395
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
6778
083721ae3dfa [project @ 2007-07-18 17:03:10 by jwe]
jwe
parents: 4167
diff changeset
1 @c Copyright (C) 1996, 1997, 2007 John W. Eaton
7018
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
2 @c
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
3 @c This file is part of Octave.
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
4 @c
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
5 @c Octave is free software; you can redistribute it and/or modify it
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
6 @c under the terms of the GNU General Public License as published by the
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
7 @c Free Software Foundation; either version 3 of the License, or (at
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
8 @c your option) any later version.
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
9 @c
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
10 @c Octave is distributed in the hope that it will be useful, but WITHOUT
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
11 @c ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
12 @c FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
13 @c for more details.
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
14 @c
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
15 @c You should have received a copy of the GNU General Public License
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
16 @c along with Octave; see the file COPYING. If not, see
fd42779a8428 [project @ 2007-10-13 00:52:12 by jwe]
jwe
parents: 6850
diff changeset
17 @c <http://www.gnu.org/licenses/>.
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
18
4167
aae05d51353c [project @ 2002-11-12 02:52:50 by jwe]
jwe
parents: 3368
diff changeset
19 @node Nonlinear Equations
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
20 @chapter Nonlinear Equations
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
21 @cindex nonlinear equations
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
22 @cindex equations, nonlinear
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
23
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
24 Octave can solve sets of nonlinear equations of the form
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
25 @iftex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
26 @tex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
27 $$
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
28 f (x) = 0
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
29 $$
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
30 @end tex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
31 @end iftex
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
32 @ifnottex
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
33
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
34 @example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
35 F (x) = 0
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
36 @end example
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
37 @end ifnottex
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
38
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
39 @noindent
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
40 using the function @code{fsolve}, which is based on the @sc{Minpack}
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
41 subroutine @code{hybrd}. This is an iterative technique so a starting
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
42 point will have to be provided. This also has the consequence that
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
43 convergence is not guarantied even if a solution exists.
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
44
3368
a4cd1e9d9962 [project @ 1999-11-20 17:22:48 by jwe]
jwe
parents: 3294
diff changeset
45 @DOCSTRING(fsolve)
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
46
3368
a4cd1e9d9962 [project @ 1999-11-20 17:22:48 by jwe]
jwe
parents: 3294
diff changeset
47 @DOCSTRING(fsolve_options)
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
48
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
49 Here is a complete example. To solve the set of equations
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
50 @iftex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
51 @tex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
52 $$
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
53 \eqalign{-2x^2 + 3xy + 4\sin(y) - 6 &= 0\cr
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
54 3x^2 - 2xy^2 + 3\cos(x) + 4 &= 0}
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
55 $$
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
56 @end tex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
57 @end iftex
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
58 @ifinfo
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
59
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
60 @example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
61 -2x^2 + 3xy + 4 sin(y) = 6
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
62 3x^2 - 2xy^2 + 3 cos(x) = -4
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
63 @end example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
64 @end ifinfo
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
65
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
66 @noindent
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
67 you first need to write a function to compute the value of the given
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
68 function. For example:
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
69
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
70 @example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
71 function y = f (x)
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
72 y(1) = -2*x(1)^2 + 3*x(1)*x(2) + 4*sin(x(2)) - 6;
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
73 y(2) = 3*x(1)^2 - 2*x(1)*x(2)^2 + 3*cos(x(1)) + 4;
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
74 endfunction
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
75 @end example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
76
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
77 Then, call @code{fsolve} with a specified initial condition to find the
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
78 roots of the system of equations. For example, given the function
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
79 @code{f} defined above,
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
80
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
81 @example
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
82 [x, info] = fsolve (@@f, [1; 2])
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
83 @end example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
84
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
85 @noindent
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
86 results in the solution
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
87
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
88 @example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
89 x =
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
90
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
91 0.57983
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
92 2.54621
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
93
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
94 info = 1
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
95 @end example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
96
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
97 @noindent
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
98 A value of @code{info = 1} indicates that the solution has converged.
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
99
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
100 The function @code{perror} may be used to print English messages
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
101 corresponding to the numeric error codes. For example,
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
102
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
103 @example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
104 @group
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
105 perror ("fsolve", 1)
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
106 @print{} solution converged to requested tolerance
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
107 @end group
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
108 @end example
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
109
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
110 When no Jacobian is supplied (as in the example above) it is approximated
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
111 numerically. This requires more function evaluations, and hence is
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
112 less efficient. In the example above we could compute the Jacobian
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
113 analytically as
3294
bfe1573bd2ae [project @ 1999-10-19 10:06:07 by jwe]
jwe
parents:
diff changeset
114
6850
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
115 @iftex
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
116 @tex
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
117 $$
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
118 \left[\matrix{ {\partial f_1 \over \partial x_1} &
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
119 {\partial f_1 \over \partial x_2} \cr
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
120 {\partial f_2 \over \partial x_1} &
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
121 {\partial f_2 \over \partial x_2} \cr}\right] =
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
122 \left[\matrix{ 3 x_2 - 4 x_1 &
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
123 4 \cos(x_2) + 3 x_1 \cr
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
124 -2 x_2^2 - 3 \sin(x_1) + 6 x_1 &
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
125 -4 x_1 x_2 \cr }\right]
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
126 $$
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
127 @end tex
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
128 which is computed with the following Octave function
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
129 @end iftex
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
130
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
131 @example
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
132 function J = jacobian(x)
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
133 J(1,1) = 3*x(2) - 4*x(1);
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
134 J(1,2) = 4*cos(x(2)) + 3*x(1);
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
135 J(2,1) = -2*x(2)^2 - 3*sin(x(1)) + 6*x(1);
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
136 J(2,2) = -4*x(1)*x(2);
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
137 endfunction
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
138 @end example
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
139
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
140 @noindent
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
141 Using this Jacobian is done with the following code
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
142
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
143 @example
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
144 [x, info] = fsolve (@{@@f, @@jacobian@}, [1; 2]);
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
145 @end example
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
146
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
147 @noindent
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
148 which gives the same solution as before.
9398f6a81bdf [project @ 2007-08-31 17:29:22 by jwe]
jwe
parents: 6778
diff changeset
149