/[ascend]/trunk/tcltk/generic/interface/Integrators.h
ViewVC logotype

Contents of /trunk/tcltk/generic/interface/Integrators.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 571 - (show annotations) (download) (as text)
Tue May 9 00:14:59 2006 UTC (14 years, 8 months ago) by johnpye
File MIME type: text/x-chdr
File size: 19639 byte(s)
Renaming 'tcltk98' to 'tcltk', continued...
1 /*
2 * Integrators.h
3 * by Kirk Abbott and Ben Allan
4 * Created: 1/94
5 * Version: $Revision: 1.16 $
6 * Version control file: $RCSfile: Integrators.h,v $
7 * Date last modified: $Date: 2003/08/23 18:43:06 $
8 * Last modified by: $Author: ballan $
9 *
10 * This file is part of the ASCEND Tcl/Tk interface
11 *
12 * Copyright 1997, Carnegie Mellon University
13 *
14 * The ASCEND Tcl/Tk interface is free software; you can redistribute
15 * it and/or modify it under the terms of the GNU General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
18 *
19 * The ASCEND Tcl/Tk interface is distributed in hope that it will be
20 * useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
22 * General Public License for more details.
23 *
24 * You should have received a copy of the GNU General Public License
25 * along with the program; if not, write to the Free Software Foundation,
26 * Inc., 675 Mass Ave, Cambridge, MA 02139 USA. Check the file named
27 * COPYING. COPYING is found in ../compiler.
28 */
29
30 /** @file
31 * Implementation of the Integration Interface.
32 * <pre>
33 * Author: Boyd T. Safrit
34 *
35 * Changes: 8/95 Added support code for blsode, an
36 * experimental user interface
37 * which identifies states, derivatives,
38 * observations, and tolerances by flags.
39 * Interface redesign (ASCEND side) by
40 * Jennifer Stokes and (C side) Ben Allan.
41 * The only required vector left is that
42 * defining the time steps to calculate
43 * results at. All state and observation
44 * data are (optionally) logged to files.
45 * 5/96 Removed old lsode interface routines.
46 * The Piela version of lsode is dead.
47 * 8/96 Removed the var == instance assumption from
48 * the implementation.
49 *
50 * Description: This module defines the general integration
51 * auxillaries for Ascend.
52 *
53 * Requires #include "tcl.h"
54 * #include "utilities/ascConfig.h
55 * #include "general/list.h"
56 * #include "compiler/instance_enum.h"
57 * #include "compiler/dimen.h"
58 * #include "solver/slv_client.h"
59 * #include "compiler/interface.h"
60 * #include <stdio.h> (for current INTEG_DEBUG functions)
61 * </pre>
62 */
63
64 #ifndef ASCTK_INTEGRATORS_H
65 #define ASCTK_INTEGRATORS_H
66
67 /**
68 * Define as TRUE to enable debug message printing.
69 * @todo this needs to go away.
70 */
71 #define INTEG_DEBUG FALSE
72 #define print_debug(msg, value) if (INTEG_DEBUG) FPRINTF(stdout, msg, value)
73 /**< Print a debug message with value if INTEG_DEBUG is TRUE. */
74 #define print_debugstring(msg) if (INTEG_DEBUG) FPRINTF(stdout, msg)
75 /**< Print a debug message string if INTEG_DEBUG is TRUE. */
76
77 /** The integrators that are supported */
78 enum Integrator_type {
79 BLSODE,
80 /* SDASSL, */
81 UNKNOWN
82 };
83
84 struct Integ_DiffData {
85 long n_eqns; /**< dimension of state vector */
86 int *input_indices; /**< vector of state vars indexes */
87 int *output_indices; /**< vector of derivative var indexes */
88 struct var_variable **y_vars; /**< NULL terminated list of states vars */
89 struct var_variable **ydot_vars; /**< NULL terminated list of derivative vars*/
90 struct rel_relation **rlist; /**< NULL terminated list of relevant rels
91 to be differentiated */
92 double **dydx_dx; /**< change in derivatives wrt states
93 I prefer to call this: d(ydot)/dy */
94 };
95
96 extern struct Integ_DiffData g_intg_diff;
97 /**<
98 * global creature for C's FEX and JEX like functions to access data
99 * when called from Fortran.
100 */
101
102 /**
103 * Experimental ivp problem description token. anyone making a copy of
104 * the y, ydot, or obs pointers who plans to free that pointer later
105 * should increment the reference count for that pointer so if blsys
106 * is destroyed later we don't end up double freeing it. Note this
107 * scheme will only work for the first copy and could be eliminated
108 * if we decoupled blsode from lsode as we will once things restabilize.
109 */
110 struct Integ_system_t {
111 /* temporaries for build. these elements will be NULL to clients */
112 struct gl_list_t *indepvars; /**< all apparent independent vars */
113 struct gl_list_t *dynvars; /**< all state and deriv instances plus indices */
114 struct gl_list_t *obslist; /**< observation instance plus indices */
115 struct gl_list_t *states; /**< ordered list of state variables and indices*/
116 struct gl_list_t *derivs; /**< ordered list of derivative (ydot) insts */
117 /* persistent, these elements are valid to C clients. */
118 struct var_variable *x; /**< independent variable */
119 struct var_variable **y; /**< array form of states */
120 struct var_variable **ydot; /**< array form of derivatives */
121 struct var_variable **obs; /**< array form of observed variables */
122 long *y_id; /**< array form of y/ydot user indices */
123 long *obs_id; /**< array form of obs user indices */
124 long n_y;
125 long n_obs;
126 int ycount; /**< number of external references to y */
127 int ydotcount; /**< number of external references to ydot */
128 int obscount; /**< number of external references to obs */
129 int maxsteps; /**< most steps between mesh poins */
130 double stepzero; /**< initial step length, SI units. */
131 double stepmin; /**< shortest step length, SI units. */
132 double stepmax; /**< longest step length, SI units. */
133 };
134
135 extern int Asc_IntegInstIntegrable(struct Instance *inst,
136 enum Integrator_type integrator);
137 /**<
138 * Takes an instance and a type of integrator and determines
139 * if the instance can be integrated using this integrator.
140 */
141
142 extern int Asc_IntegCheckStatus(slv_status_t status);
143 /**<
144 * Takes a status and checks for convergence. If converged then return 1
145 * else return 0 and print error message if appropriate.
146 */
147
148 /* functions affecting the logging of data during integration */
149 extern FILE *Asc_IntegOpenYFile(void);
150 /**<
151 * Returns the pointer to a file to which state
152 * variables should be written at each recorded time interval. Other
153 * info may be written here as well, but that is up to the designer
154 * of the specific integrator interface. The filename used will be
155 * the one most recently set via Asc_IntegSetYFileCmd(). If no
156 * name (or an empty name) has been set via Asc_IntegSetYFileCmd(),
157 * the file returned will be NULL.<br><br>
158 *
159 * This function may return NULL, in which case do not write to
160 * them as there was some problem opening the file.
161 * If opened successfully, the file returned will have DATASET and
162 * the open time stamped at the bottom. This will separate multiple
163 * runs or changes in the same run.<br><br>
164 *
165 * Closing the files is your job. When finished with the file, you
166 * should call Asc_IntegReleaseYFile() before closing it.
167 */
168 extern FILE *Asc_IntegOpenObsFile(void);
169 /**<
170 * Returns the pointer to a file to which observation
171 * variables should be written at each recorded time interval. Other
172 * info may be written here as well, but that is up to the designer
173 * of the specific integrator interface. The filename used will be
174 * the one most recently set via Asc_IntegSetObsFileCmd(). If no
175 * name (or an empty name) has been set via Asc_IntegSetObsFileCmd(),
176 * the file returned will be NULL.<br><br>
177 *
178 * This function may return NULL, in which case do not write to
179 * them as there was some problem opening the file.
180 * If opened successfully, the file returned will have DATASET and
181 * the open time stamped at the bottom. This will separate multiple
182 * runs or changes in the same run.<br><br>
183 *
184 * Closing the files is your job. When finished with the file, you
185 * should call Asc_IntegReleaseObsFile() before closing it.
186 */
187 extern FILE *Asc_IntegGetYFile(void);
188 /**<
189 * Returns the current FILE * for writing state variables,
190 * or NULL if none.
191 */
192 extern FILE *Asc_IntegGetObsFile(void);
193 /**<
194 * Returns the current FILE * for writing observation variables,
195 * or NULL if none.
196 */
197 extern void Asc_IntegReleaseYFile(void);
198 extern void Asc_IntegReleaseObsFile(void);
199 /**<
200 * Releases the internally-stored FILE * for observation variables.
201 * This does not close the file (which you still need to do).
202 * This function should be called just before closing the files,
203 * though it may be called sooner if you want to keep the file pointer
204 * returned all to yourself.
205 */
206
207 extern void Asc_IntegPrintYHeader(FILE *fp, struct Integ_system_t *blsys);
208 /**<
209 * Prints Y header info to the file given.
210 * If FILE is NULL, returns immediately.
211 * If FILE ok, but blsys NULL prints warning to stderr and returns.
212 * At the moment, the written header info is:
213 * A vertical listing of
214 * ode_id/obs_id instance_name_as_seen_in_solver.d UNITS
215 * then a row (which will line up with printed values)
216 * of the ode_id/obs_id.
217 */
218 extern void Asc_IntegPrintObsHeader(FILE *fp, struct Integ_system_t *blsys);
219 /**<
220 * Prints observation header info to the file given.
221 * If FILE is NULL, returns immediately.
222 * If FILE ok, but blsys NULL prints warning to stderr and returns.
223 * At the moment header info is:
224 * A vertical listing of
225 * ode_id/obs_id instance_name_as_seen_in_solver.d UNITS
226 * then a row (which will line up with printed values)
227 * of the ode_id/obs_id.
228 */
229
230 extern void Asc_IntegPrintYLine(FILE *fp, struct Integ_system_t *blsys);
231 /**<
232 * Prints a Y line to the file given.
233 * If FILE is NULL, returns immediately.
234 * If FILE ok, but blsys NULL prints warning to stderr and returns.
235 * Each line is a set of values, identified by the header, following
236 * the value of the independent variable. These lines are of arbitrary
237 * length since it is expected that some other program will be
238 * reading the output, not a human.
239 */
240 extern void Asc_IntegPrintObsLine(FILE *fp, struct Integ_system_t *blsys);
241 /**<
242 * Prints an observation line to the file given.
243 * If FILE is NULL, returns immediately.
244 * If FILE ok, but blsys NULL prints warning to stderr and returns.
245 * Each line is a set of values, identified by the header, following
246 * the value of the independent variable. These lines are of arbitrary
247 * length since it is expected that some other program will be
248 * reading the output, not a human.
249 */
250
251 extern int Asc_IntegSetYFileCmd(ClientData cdata, Tcl_Interp *interp,
252 int argc, CONST84 char *argv[]);
253 /**<
254 * Set the next filename to be used for integration state variable
255 * logging. Does not check the filesystem sanity of the given name.<br><br>
256 *
257 * Registered as: integrate_set_y_file filename
258 */
259
260 extern int Asc_IntegSetObsFileCmd(ClientData cdata, Tcl_Interp *interp,
261 int argc, CONST84 char *argv[]);
262 /**<
263 * Set the next filename to be used for integration observation variable
264 * logging. Does not check the filesystem sanity of the given name.<br><br>
265 *
266 * Registered as: integrate_set_obs_file filename
267 */
268
269 extern int Asc_IntegSetFileUnitsCmd(ClientData cdata, Tcl_Interp *interp,
270 int argc, CONST84 char *argv[]);
271 /**<
272 * Sets output to be in SI (internal) units or in the user set display
273 * units. If display is selected and the variable to be printed cannot
274 * be displayed in the display units because of floating point errors
275 * in conversion, the value will be printed in SI units and a warning
276 * issued to stderr. The user is then expected to deal with the problem
277 * of inconsistency in his output.
278 * Only the first letter of the argument is significant.<br><br>
279 *
280 * Registered as: integrate_logunits <display,si>
281 */
282
283 extern int Asc_IntegSetFileFormatCmd(ClientData cdata, Tcl_Interp *interp,
284 int argc, CONST84 char *argv[]);
285 /**<
286 * Sets output data to be in fixed column width (with extra white) or in
287 * space separated columns whose width is variable from line to line.
288 * Only the first letter of the argument is significant.<br><br>
289 *
290 * Registered as: integrate_logformat <fixed,variable>
291 */
292
293 extern int Asc_IntegSetXSamples(long ns, double *values, dim_type *dimp);
294 /**<
295 * <!-- Asc_IntegSetXSamples(ns, values, dimp); -->
296 * Sets values of time samples to the values given (ns of them) and
297 * keeps both the dim pointer and vector given. The vector and dimp
298 * given may be freed by calling this again, but xsamples owns
299 * them until then. If called with ns<1 or values==NULL, will free any
300 * previously captured values/dimp. If called with dimp==NULL, will
301 * assume WildDimension. Don't call this with a dimp we can't free later.
302 * Return is 1 if for some reason we can't set as expected, 0 otherwise.
303 */
304
305 extern double Asc_IntegGetStepZero(struct Integ_system_t *sys);
306 /**<
307 * Returns the length of the initial step user specified,
308 * or 0.0 if none was set.
309 */
310
311 extern double Asc_IntegGetStepMin(struct Integ_system_t *sys);
312 /**<
313 * Returns the length of the longest allowable step.
314 * or 0.0 if none was set by user.
315 */
316
317 extern double Asc_IntegGetStepMax(struct Integ_system_t *sys);
318 /**<
319 * Returns the length of the shortest allowable step.
320 * or 0.0 if none was set by user.
321 */
322
323 extern int Asc_IntegGetMaxSteps(struct Integ_system_t *sys);
324 /**<
325 * Returns the most internal steps allowed between
326 * two time samples, or 0 if none was set by user.
327 */
328
329 extern long Asc_IntegGetNSamples(void);
330 /**<
331 * Returns the number of values currently stored in xsamples.
332 */
333
334 extern double Asc_IntegGetXSamplei(long i);
335 /**<
336 * Returns the value stored in xsamples[i]. Will whine if
337 * if xsample[i] does not exist. No, there is no way to get
338 * back the pointer to the xsamples vector.
339 */
340 extern void Asc_IntegSetXSamplei(long i, double val);
341 /**<
342 * Sets the value stored in xsamples[i] to val. Will whine if
343 * if xsample[i] does not exist. No, there is no way to get
344 * back the pointer to the xsamples vector.
345 */
346
347 extern dim_type *Asc_IntegGetXSamplesDim(void);
348 /**<
349 * Returns a pointer to the dimensionality of the samples currently
350 * stored, or NULL if none are stored. Do not free this pointer: we
351 * own it.
352 */
353
354 extern int Asc_IntegGetXSamplesCmd(ClientData cdata, Tcl_Interp *interp,
355 int argc, CONST84 char *argv[]);
356 /**<
357 * Returns the current vector of independent variable samples as fed to
358 * smarter integrators such as lsode. Returns 2 element list:
359 * {units of x} {list of doubles in the units indicated in the first
360 * element}. In the event that display is given but the sample values
361 * (which are stored in si) cannot be converted to display units, the
362 * units of x in the first element will change and the SI values will
363 * be returned in the second element. If display is not given, samples
364 * will be returned in SI units.
365 * The units used will be dimensionally consistent with those specified
366 * when integrate_set_samples was called.
367 * Only the first letter of the optional argument display is significant.<br><br>
368 *
369 * Registered as: integrate_get_samples [display]
370 */
371
372 extern int Asc_IntegSetXSamplesCmd(ClientData data, Tcl_Interp *interp,
373 int argc, CONST84 char *argv[]);
374 /**<
375 * Takes the values given and sets the intervals for smarter
376 * integrators accordingly. Values will be converted to SI values
377 * using the units given. ?,* is a legal unit, in which case the
378 * values will be assumed to be SI requiring no conversion.
379 * If no arguments are given, we will free memory allocated for storage
380 * of the list, otherwise at least 2 values are required.<br><br>
381 *
382 * Note: It is the users responsibility to make sure that the dimens/
383 * units this is called with match the independent variable in the ivp.
384 * C clients are given the power to verify this.
385 * Also, the units given must already exist: we do not create them.<br><br>
386 *
387 * Registered as: integrate_set_samples <units> <value [value...] value>
388 */
389
390 extern int Asc_IntegInstIntegrableCmd(ClientData cdata, Tcl_Interp *interp,
391 int argc, CONST84 char *argv[]);
392 /**<
393 * Returns "1" if instance is integrable by the named integrator, "0"
394 * otherwise. Currently valid names are: lsode
395 * The instance examined is "current", "search", or "solver".<br><br>
396 *
397 * Registered as: integrate_able <instance> <NAME>
398 */
399
400 /*
401 * The following functions fetch and set parts with names specific to
402 * the type definitions in ivp.lib, the ASCEND initial value problem
403 * model file. They are for use by any ivp solver interface.
404 * All names are given at the scope of ivp.
405 */
406
407 /*
408 * Problem size parameters.
409 */
410
411 /*
412 * Parts of type definition derivatives refinements.
413 */
414
415 extern double Asc_IntegGetDX(void);
416 /**<
417 * <!-- IntegGet_d_x(value) -->
418 * Gets d.x value.
419 */
420
421 extern void Asc_IntegSetDX(double value);
422 /**<
423 * <!-- IntegSet_d_x(value) -->
424 * Sets d.x to value.
425 */
426
427 extern double *Asc_IntegGetDY(double *vector);
428 /**<
429 * <!-- vector=IntegGetDY(vector) -->
430 * Return the vector d.y.
431 * If vector given is NULL, allocates vector, which the caller then owns.
432 * Vector, if given, should be IntegGet_d_neq()+1 long.
433 */
434
435 extern void Asc_IntegSetDY(double *vector);
436 /**<
437 * <!-- IntegSet_d_y(vector) -->
438 * Sets d.y[] to values in vector.
439 */
440
441 extern double *Asc_IntegGetDDydx(double *vector);
442 /**<
443 * <!-- IntegGet_d_dydx(vector) -->
444 * Returns the vector d.dydx.
445 * If vector given is NULL, allocates vector, which the caller then owns.
446 * Vector, if given, should be IntegGet_d_neq()+1 long.
447 */
448
449 extern void Asc_IntegSetDDydx(double *vector);
450 /**<
451 * <!-- IntegSet_d_dydx(vector) -->
452 * Sets d.dydx[] to values in vector.
453 */
454
455 extern double *Asc_IntegGetDObs(double *vector);
456 /**<
457 * <!-- IntegGet_d_obs(vector) -->
458 * Returns the vector d.obs.
459 * Vector should be of sufficient length (g_intginst_d_n_obs+1).
460 * If NULL is given a vector is allocated and is the caller's
461 * responsibility to deallocate.
462 */
463
464 /*
465 * END of the general ivp type utilities.
466 */
467
468 extern int Asc_IntegSetupCmd(ClientData cdata,Tcl_Interp *interp,
469 int argc, CONST84 char *argv[]);
470 /**<
471 * Set up the integrator.
472 * itype is one of the supported integrator types (currently:
473 * BLSODE)
474 * n1 is the start index to run the integrator from.
475 * n2 is the end index to stop at.
476 *
477 * Registered as: integrateSetup itype n1 n2
478 */
479
480 extern int Asc_IntegCleanupCmd(ClientData cdata, Tcl_Interp *interp,
481 int argc, CONST84 char *argv[]);
482 /**<
483 * Do any C level housekeeping need after calling integration.
484 * Always use integrate_setup/integrate_cleanup in pairs.<br><br>
485 *
486 * Registered as: integrate_cleanup
487 */
488
489 #endif /* ASCTK_INTEGRATORS_H */
490

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