annotate src/fsolve.cc @ 114:a20da3b44f08

[project @ 1993-09-18 09:06:35 by jwe]
author jwe
date Sat, 18 Sep 1993 09:06:35 +0000
parents d68036bcad4c
children d1c5e5edbf1e
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
7
d68036bcad4c [project @ 1993-08-08 03:53:43 by jwe]
jwe
parents: 1
diff changeset
1 // f-fsolve.cc -*- C++ -*-
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
2 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
3
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
4 Copyright (C) 1993 John W. Eaton
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
5
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
6 This file is part of Octave.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
7
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
8 Octave is free software; you can redistribute it and/or modify it
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
9 under the terms of the GNU General Public License as published by the
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
10 Free Software Foundation; either version 2, or (at your option) any
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
11 later version.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
12
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
13 Octave is distributed in the hope that it will be useful, but WITHOUT
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
16 for more details.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
17
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
18 You should have received a copy of the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
19 along with Octave; see the file COPYING. If not, write to the Free
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
20 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
21
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
22 */
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
23
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
24 #ifdef __GNUG__
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
25 #pragma implementation
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
26 #endif
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
27
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
28 #include "NLEqn.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
29
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
30 #include "tree-const.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
31 #include "variables.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
32 #include "gripes.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
33 #include "error.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
34 #include "utils.h"
7
d68036bcad4c [project @ 1993-08-08 03:53:43 by jwe]
jwe
parents: 1
diff changeset
35 #include "f-fsolve.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
36
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
37 // Global pointer for user defined function required by hybrd1.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
38 static tree *fsolve_fcn;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
39
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
40 #ifdef WITH_DLD
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
41 tree_constant *
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
42 builtin_fsolve (tree_constant *args, int nargin, int nargout)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
43 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
44 return fsolve (args, nargin, nargout);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
45 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
46 #endif
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
47
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
48 int
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
49 hybrd_info_to_fsolve_info (int info)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
50 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
51 switch (info)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
52 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
53 case 0:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
54 info = -1;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
55 break;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
56 case 1:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
57 break;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
58 case 2:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
59 info = 4;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
60 break;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
61 case 3:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
62 case 4:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
63 case 5:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
64 info = 3;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
65 break;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
66 default:
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
67 panic_impossible ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
68 break;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
69 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
70 return info;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
71 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
72
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
73 ColumnVector
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
74 fsolve_user_function (ColumnVector& x)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
75 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
76 ColumnVector retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
77
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
78 int n = x.capacity ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
79
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
80 // tree_constant name = tree_constant (fsolve_fcn->name ());
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
81 tree_constant *args = new tree_constant [2];
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
82 // args[0] = name;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
83
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
84 if (n > 1)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
85 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
86 Matrix m (n, 1);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
87 for (int i = 0; i < n; i++)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
88 m (i, 0) = x.elem (i);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
89 tree_constant vars (m);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
90 args[1] = vars;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
91 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
92 else
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
93 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
94 double d = x.elem (0);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
95 tree_constant vars (d);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
96 args[1] = vars;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
97 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
98
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
99 if (fsolve_fcn != NULL_TREE)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
100 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
101 tree_constant *tmp = fsolve_fcn->eval (args, 2, 1, 0);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
102 delete [] args;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
103 if (tmp != NULL_TREE_CONST && tmp[0].is_defined ())
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
104 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
105 retval = tmp[0].to_vector ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
106 delete [] tmp;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
107 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
108 else
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
109 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
110 delete [] tmp;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
111 gripe_user_supplied_eval ("fsolve");
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
112 jump_to_top_level ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
113 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
114 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
115
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
116 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
117 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
118
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
119 tree_constant *
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
120 fsolve (tree_constant *args, int nargin, int nargout)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
121 {
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
122 // Assumes that we have been given the correct number of arguments.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
123
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
124 tree_constant *retval = NULL_TREE_CONST;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
125
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
126 fsolve_fcn = is_valid_function (args[1], "fsolve", 1);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
127 if (fsolve_fcn == NULL_TREE
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
128 || takes_correct_nargs (fsolve_fcn, 2, "fsolve", 1) != 1)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
129 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
130
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
131 ColumnVector x = args[2].to_vector ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
132
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
133 if (nargin > 3)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
134 message ("fsolve", "ignoring optional arguments...");
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
135
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
136 if (nargout > 2)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
137 message ("fsolve", "can't compute path output yet...");
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
138
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
139 NLFunc foo_fcn (fsolve_user_function);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
140 NLEqn foo (x, foo_fcn);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
141
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
142 int info;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
143 ColumnVector soln = foo.solve (info);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
144
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
145 info = hybrd_info_to_fsolve_info (info);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
146
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
147 retval = new tree_constant [nargout+1];
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
148 retval[0] = tree_constant (soln, 1);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
149
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
150 if (nargout > 1)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
151 retval[1] = tree_constant ((double) info);
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
152
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
153 if (nargout > 2)
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
154 retval[2] = tree_constant ();
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
155
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
156 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
157 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
158
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
159 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
160 ;;; Local Variables: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
161 ;;; mode: C++ ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
162 ;;; page-delimiter: "^/\\*" ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
163 ;;; End: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
164 */
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
165