/[ascend]/trunk/models/johnpye/tron/tron.h
ViewVC logotype

Contents of /trunk/models/johnpye/tron/tron.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1500 - (show annotations) (download) (as text)
Thu Jun 21 14:38:59 2007 UTC (17 years, 5 months ago) by jpye
File MIME type: text/x-chdr
File size: 8183 byte(s)
Preliminary work for TRON tie-in.
1 /**
2 This is a header file to be used when calling TRON from C.
3
4 TRON is written by Chih-Jen Lin and Jorge J. More', but is not
5 included in the ASCEND distribution due to licensing restrictions.
6 */
7 #ifndef TRON_H
8 #define TRON_H
9
10 #ifdef linux
11 # define DTRON dtron_
12 # define FNAME_LCASE_DECOR
13 #else
14 # error "System-dependent information required"
15 #endif
16
17 int DTRON(const int *n, double *x, const double *xl,
18 const double *xu, double *f, double *g, double *a,
19 double *adiag, int *acol_ptr, int *arow_ind,
20 double *frtol, double *fatol, double *fmin, double *cgtol,
21 int *itermax, double *delta, char *task, double *b,
22 double *bdiag, int *bcol_ptr, int *brow_ind,
23 double *l, double *ldiag, int *lcol_ptr, int *lrow_ind,
24 double *xc, double *s, int *indfree, int *isave,
25 double *dsave, double *wa, int *iwa
26 );
27 /**<
28 This is the TRON trust region Newton method solver for the
29 solution of large bound-constrained optimization problems
30
31 min { f(x) : xl <= x <= xu }
32
33 where the Hessian matrix is sparse. The user must evaluate the
34 function, gradient, and the Hessian matrix.
35
36 @param n number of variables (in)
37 @param x double precision array of dimension n (in/out)
38 On entry x specifies the vector x.
39 On exit x is the final minimizer.
40
41 @param xl vector of length n, containing lower bound for each var (in)
42 @param xu vector of length n, containing upper bound for each var (in)
43
44 @param f double precision variable.
45 On entry f must contain the function at x.
46 On exit f is unchanged.
47
48 @param g double precision array of dimension n.
49 On entry g must contain the gradient at x.
50 On exit g is unchanged.
51
52 @param a double precision array of dimension nnz.
53 On entry a must contain the strict lower triangular part
54 of A in compressed column storage.
55 On exit a is unchanged.
56
57 @param adiag double precision array of dimension n.
58 On entry adiag must contain the diagonal elements of A.
59 On exit adiag is unchanged.
60
61 @param acol_ptr integer array of dimension n + 1.
62 On entry acol_ptr must contain pointers to the columns of A.
63 The nonzeros in column j of A must be in positions
64 acol_ptr(j), ... , acol_ptr(j+1) - 1.
65 On exit acol_ptr is unchanged.
66
67 @param arow_ind integer array of dimension nnz.
68 On entry arow_ind must contain row indices for the strict
69 lower triangular part of A in compressed column storage.
70 On exit arow_ind is unchanged.
71
72 @param frtol double precision variable.
73 On entry frtol specifies the relative error desired in the
74 function. Convergence occurs if the estimate of the
75 relative error between f(x) and f(xsol), where xsol
76 local minimizer, is less than frtol.
77 On exit frtol is unchanged.
78
79 @param fatol double precision variable.
80 On entry fatol specifies the absolute error desired in the
81 function. Convergence occurs if the estimate of the
82 absolute error between f(x) and f(xsol), where xsol
83 is a local minimizer, is less than fatol.
84 On exit fatol is unchanged.
85
86 @param fmin double precision variable.
87 On entry fmin specifies a lower bound for the function.
88 The subroutine exits with a warning if f < fmin.
89 On exit fmin is unchanged.
90
91 @param cgtol double precision variable.
92 On entry cgtol specifies the convergence criteria for
93 the conjugate gradient method.
94 On exit cgtol is unchanged.
95
96 @param itermax max number of conjugate gradient iterations (input)
97
98 @param delta trust region bound (input)
99
100 @param task character variable of length at least 60.
101 On initial entry task must be set to 'START'.
102 On exit task indicates the required action:
103
104 If task(1:1) = 'F' then evaluate the function at x.
105
106 If task(1:2) = 'GH' then evaluate the gradient and the
107 Hessian matrix at x.
108
109 If task(1:4) = 'CONV' then the search is successful.
110
111 If task(1:4) = 'WARN' then the subroutine is not able
112 to satisfy the convergence conditions. The exit value
113 of x contains the best approximation found.
114
115 ---
116 @param b double precision array of dimension nnz + n*p.
117 On entry b need not be specified.
118 On exit b contains the strict lower triangular part
119 of B in compressed column storage.
120
121 @param bdiag diagonal elements of B (returned into user-allocated space of size n)
122
123 @param bcol_ptr integer array of dimension n + 1.
124 On entry bcol_ptr need not be specified
125 On exit bcol_ptr contains pointers to the columns of B.
126 The nonzeros in column j of B are in the
127 bcol_ptr(j), ... , bcol_ptr(j+1) - 1 positions of b.
128
129 @param brow_ind integer array of dimension nnz.
130 On entry brow_ind need not be specified.
131 On exit brow_ind contains row indices for the strict lower
132 triangular part of B in compressed column storage.
133 ---
134 @param l double precision array of dimension nnz + n*p.
135 On entry l need not be specified.
136 On exit l contains the strict lower triangular part
137 of L in compressed column storage.
138
139 @param ldiag double precision array of dimension n.
140 On entry ldiag need not be specified.
141 On exit ldiag contains the diagonal elements of L.
142
143 @param lcol_ptr integer array of dimension n + 1.
144 On entry lcol_ptr need not be specified.
145 On exit lcol_ptr contains pointers to the columns of L.
146 The nonzeros in column j of L are in the
147 lcol_ptr(j), ... , lcol_ptr(j+1) - 1 positions of l.
148
149 @param lrow_ind integer array of dimension nnz + n*p.
150 On entry lrow_ind need not be specified.
151 On exit lrow_ind contains row indices for the strict lower
152 triangular part of L in compressed column storage.
153 ---
154
155 @param xc double precision working array of dimension n.
156
157 @param s double precision working array of dimension n.
158
159 @param indfree integer working array of dimension n.
160
161 @param isave integer working array of dimension 3.
162
163 @param dsave is double precision working array of dimension 3.
164
165 @param wa double precision work array of dimension 7*n.
166
167 @param iwa integer work array of dimension 3*n.
168
169
170 This subroutine uses reverse communication.
171 The user must choose an initial approximation x to the minimizer,
172 and make an initial call with task set to 'START'.
173 On exit task indicates the required action.
174
175 A typical invocation has the following outline:
176
177 Compute a starting vector x.
178 Compute the sparsity pattern of the Hessian matrix and
179 store in compressed column storage in (acol_ptr,arow_ind).
180
181 char task[60] = "START";
182 while(search){
183
184 if(strcmp(task,'F')==0 || strcmp(task,"START")==0){
185
186 // Evaluate the function at x and store in f.
187
188 }
189 if(strcmp(task,'GH')==0 || strcmp(task,"START")==0){
190
191 // Evaluate the gradient at x and store in g.
192
193 // Evaluate the Hessian at x and store in compressed
194 column storage in (a,adiag,acol_ptr,arow_ind)
195 }
196
197 dtron(n,x,xl,xu,f,g,a,adiag,acol_ptr,arow_ind,
198 frtol,fatol,fmin,cgtol,itermax,delta,task,
199 b,bdiag,bcol_ptr,brow_ind,
200 l,ldiag,lcol_ptr,lrow_ind,
201 xc,s,indfree,
202 isave,dsave,wa,iwa
203 );
204
205 if(strncmp(task,"CONV",4)==0)search = FALSE;
206
207 }
208
209 NOTE: The user must not alter work arrays between calls.
210
211 TRON calls the following routines from MINPACK-2:
212 dcauchy, dspcg, dssyax
213
214 TRON also calls the 'dcopy' routine from Level 1 BLAS:
215 */
216
217 /* the following macros help with dlopening etc */
218
219 #define TRON_DTRON_ARGS (const int *n, double *x, const double *xl, \
220 const double *xu, double *f, double *g, double *a,\
221 double *adiag, int *acol_ptr, int *arow_ind,\
222 double *frtol, double *fatol, double *fmin, double *cgtol,\
223 int *itermax, double *delta, char *task, double *b,\
224 double *bdiag, int *bcol_ptr, int *brow_ind,\
225 double *l, double *ldiag, int *lcol_ptr, int *lrow_ind,\
226 double *xc, double *s, int *indfree, int *isave,\
227 double *dsave, double *wa, int *iwa)
228
229 #define TRON_DTRON_VALS (n,x,xl,xu,f,g,a,adiag,acol_ptr,arow_ind,\
230 frtol,fatol,fmin,cgtol,itermax,delta,task,\
231 b,bdiag,bcol_ptr,brow_ind,\
232 l,ldiag,lcol_ptr,lrow_ind,\
233 xc,s,indfree,\
234 isave,dsave,wa,iwa)
235
236 typedef int dtron_fn_t TRON_DTRON_ARGS;
237
238 #endif /* TRON_H */
239

john.pye@anu.edu.au
ViewVC Help
Powered by ViewVC 1.1.22