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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 571 - (hide 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 johnpye 571 /*
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