Mercurial > hg > octave-avbm

comparison src/fsolve.cc @ 1:78fd87e624cb

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