/[ascend]/trunk/base/generic/solver/integrator.h
ViewVC logotype

Contents of /trunk/base/generic/solver/integrator.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 941 - (show annotations) (download) (as text)
Fri Nov 24 10:46:32 2006 UTC (13 years, 11 months ago) by johnpye
File MIME type: text/x-chdr
File size: 15311 byte(s)
Changed integrator_set_engine to return 0 on success.
Fixed Integrator::setEngine to throw range_error / IndexError on invalid selection.
Test suite contains testIDA that works now (more tests yet to come)
1 /* ASCEND modelling environment
2 Copyright (C) 2006 Carnegie Mellon University
3
4 This program is free software; you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation; either version 2, or (at your option)
7 any later version.
8
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13
14 You should have received a copy of the GNU General Public License
15 along with this program; if not, write to the Free Software
16 Foundation, Inc., 59 Temple Place - Suite 330,
17 Boston, MA 02111-1307, USA.
18 *//** @file
19 Implementation of non-GUI-specific parts of the Integration Interface.
20
21 Here we are integrating a model with independent variable t for a specified
22 set of timesteps, and and observing the values of certain of the
23 independent variables.
24
25 You can specify the values of the time steps that you want via the
26 'samplelist' interface.
27
28 You can have your results reported however you like using the
29 'IntegratorReporter' interface.
30
31 There's nothing here yet to explicitly support DAEs -- that's the next task.
32
33 (old annotation:)
34 The following functions fetch and set parts with names specific to
35 the type definitions in ivp.lib, the ASCEND initial value problem
36 model file. They are for use by any ivp solver interface.
37 All names are given at the scope of ivp.
38
39 *//*
40 by John Pye, May 2006.
41
42 Derived from tcltk/.../Integrators.h (with heavy editing).
43 Removed all global variables and created a 'IntegratorReporter' interface
44 for reporting results of integration to whatever interface you're using.
45 */
46
47 #ifndef ASC_INTEGRATOR_H
48 #define ASC_INTEGRATOR_H
49
50 #include <utilities/ascConfig.h>
51 #include <utilities/ascMalloc.h>
52 #include <general/list.h>
53 #include <general/dstring.h>
54 #include <compiler/compiler.h>
55 #include <compiler/instance_enum.h>
56
57 #include <compiler/fractions.h>
58 #include <compiler/dimen.h>
59 #include <compiler/units.h>
60 #include <compiler/module.h>
61 #include <compiler/library.h>
62 #include <compiler/expr_types.h>
63 #include <compiler/child.h>
64 #include <compiler/type_desc.h>
65 #include <compiler/atomvalue.h>
66 #include <compiler/instance_name.h>
67 #include <compiler/instquery.h>
68 #include <compiler/parentchild.h>
69 #include <compiler/symtab.h>
70 #include <compiler/instance_io.h>
71
72 #include "slv_types.h"
73 #include "mtx.h"
74 #include "var.h"
75 #include "rel.h"
76 #include "discrete.h"
77 #include "conditional.h"
78 #include "logrel.h"
79 #include "bnd.h"
80 #include "slv_common.h"
81 #include "samplelist.h"
82
83 /*---------------------------*/
84 /**
85 Struct containin the list of supported integrators
86 */
87 typedef enum{
88 INTEG_UNKNOWN,
89 INTEG_LSODE,
90 INTEG_IDA
91 } IntegratorEngine;
92
93 /*------------------------
94 abstraction of the integrator output interface
95 */
96 struct IntegratorSystemStruct;
97
98 /**
99 Initialisation. This hook allows initialisation of the GUI or reporting
100 mechanism to be performed when integration begins
101 */
102 typedef int IntegratorOutputInitFn(struct IntegratorSystemStruct *);
103
104 /** Status report. This hook allows raw data to be output as integration
105 proceeds, or for a GUI to perform status updates. An integrator should check
106 the return status on this one, as this is the suggested way to perform
107 GUI interruption of the integrator.
108
109 @return 1 on success, 0 on user interrupt
110 */
111 typedef int IntegratorOutputWriteFn(struct IntegratorSystemStruct *);
112
113 /**
114 Observation reporting. This hook should be implemented to record
115 observations in a way that can be presented to the use, recorded in a
116 file, etc.
117 */
118 typedef int IntegratorOutputWriteObsFn(struct IntegratorSystemStruct *);
119
120 /**
121 Finalisation. This hook can be used to terminate recording of observations,
122 close files, terminate GUI status reporting, etc.
123 */
124 typedef int IntegratorOutputCloseFn(struct IntegratorSystemStruct *);
125
126 /**
127 This struct allows arbitrary functions to be used for the reporting
128 of integrator progress.
129 */
130 typedef struct{
131 IntegratorOutputInitFn *init;
132 IntegratorOutputWriteFn *write;
133 IntegratorOutputWriteObsFn *write_obs;
134 IntegratorOutputCloseFn *close;
135 } IntegratorReporter;
136
137 /*------------------------------------*/
138 /**
139 Initial value problem description struct. Anyone making a copy of
140 the y, ydot, or obs pointers who plans to free that pointer later
141 should increment the reference count for that pointer so if blsys
142 is destroyed later we don't end up double freeing it. Note this
143 scheme will only work for the first copy and could be eliminated
144 if we decoupled blsode from lsode as we will once things restabilize.
145
146 JP: adding to this structure the instance being solved, the engine
147 that we're solving it with, and a struct containing the output function ptrs
148 */
149 struct IntegratorSystemStruct{
150 struct Instance *instance; /**< not sure if this one is really necessary... -- JP */
151 slv_system_t system; /**< the system that we're integrating in ASCEND */
152 IntegratorEngine engine; /**< enum containing the ID of the integrator engine we're using */
153 IntegratorReporter *reporter;/**< functions for reporting integration results */
154 SampleList *samples; /**< pointer to the list of samples. we *don't own* this **/
155 void *enginedata; /**< space where the integrator engine can store stuff */
156 void *clientdata; /**< any stuff that the GUI/CLI needs to associate with this */
157
158 int nstates; /**< was a local global in integrator.c, moved it here. */
159 int nderivs; /**< ditto, as for nstates */
160
161 /* temporaries for build. these elements will be NULL to clients */
162 struct gl_list_t *indepvars; /**< all apparent independent vars */
163 struct gl_list_t *dynvars; /**< all state and deriv instances plus indices */
164 struct gl_list_t *obslist; /**< observation instance plus indices */
165 struct gl_list_t *states; /**< ordered list of state variables and indices*/
166 struct gl_list_t *derivs; /**< ordered list of derivative (ydot) insts */
167
168 /* persistent, these elements are valid to C clients. */
169 struct var_variable *x; /**< independent variable */
170 struct var_variable **y; /**< array form of states */
171 struct var_variable **ydot; /**< array form of derivatives */
172 struct var_variable **obs; /**< array form of observed variables */
173 long *y_id; /**< array form of y/ydot user indices, for DAEs we use negatives here for derivative vars */
174 long *obs_id; /**< array form of obs user indices */
175 long n_y;
176 long n_obs;
177 long currentstep; /**< current step number (also @see integrator_getnsamples) */
178 int ycount; /**< number of external references to y */
179 int ydotcount; /**< number of external references to ydot */
180 int obscount; /**< number of external references to obs */
181 int maxsubsteps; /**< most steps between mesh poins */
182 double stepzero; /**< initial step length, SI units. */
183 double minstep; /**< shortest step length, SI units. */
184 double maxstep; /**< longest step length, SI units. */
185 };
186
187 typedef struct IntegratorSystemStruct IntegratorSystem;
188
189 /*------------------------------------------------------------------------------
190 PUBLIC INTERFACE (for use by the GUI/CLI)]
191 */
192
193 ASC_DLLSPEC(IntegratorSystem *) integrator_new(slv_system_t sys, struct Instance *inst);
194
195 ASC_DLLSPEC(int) integrator_analyse(IntegratorSystem *blsys);
196
197 ASC_DLLSPEC(int) integrator_solve(IntegratorSystem *blsys, long i0, long i1);
198 /**<
199 Takes the type of integrator and sets up the global variables into the
200 current integration instance.
201
202 @return 1 on success, 0 on failure.
203 */
204
205 ASC_DLLSPEC(void) integrator_free(IntegratorSystem *blsys);
206 /**<
207 Deallocates any memory used and sets all integration global points to NULL.
208 */
209
210 ASC_DLLSPEC(int) integrator_set_engine(IntegratorSystem *blsys, IntegratorEngine engine);
211 /**
212 Sets the engine for this integrator. Checks that the integrator can be used
213 on the given system.
214
215 @return 0 on success,
216 1 if unknown engine,
217 2 if invalid engine (not suitable for this present problem)
218 */
219
220 ASC_DLLSPEC(IntegratorEngine) integrator_get_engine(const IntegratorSystem *blsys);
221 /**
222 Returns the engine (ID) selected for use in this IntegratorSystem (may return
223 INTEG_UNKNOWN if none or invalid setting).
224 */
225
226 ASC_DLLSPEC(int) integrator_set_reporter(IntegratorSystem *blsys
227 , IntegratorReporter *r);
228 /**<
229 Use this function from your interface to tell the integrator how it needs
230 to make its output. You need to pass in a struct containing the required
231 function pointers. (We will wrap IntegratorReporter and access it from
232 Python, eventually)
233 */
234
235 ASC_DLLSPEC(int) integrator_find_indep_var(IntegratorSystem *blsys);
236 /**<
237 Attempt to locate the independent variable. It will be stored in the blsys
238 structure if found.
239 @return 1 if found, 0 if not found.
240 */
241
242 extern int integrator_checkstatus(slv_status_t status);
243 /**<
244 * Takes a status and checks for convergence. If converged then return 1
245 * else return 0 and print error message if appropriate.
246 */
247
248 ASC_DLLSPEC(void) integrator_set_samples(IntegratorSystem *blsys, SampleList *samples);
249 /**<
250 Sets values of time samples to the values given (ns of them) and
251 keeps both the dim pointer and vector given. The vector and dimp
252 given may be freed by calling this again, but xsamples owns
253 them until then. If called with ns<1 or values==NULL, will free any
254 previously captured values/dimp. If called with dimp==NULL, will
255 assume WildDimension. Don't call this with a dimp we can't free later.
256 Return is 1 if for some reason we can't set as expected, 0 otherwise.
257 */
258
259 ASC_DLLSPEC(void) integrator_set_stepzero(IntegratorSystem *blsys, double);
260 ASC_DLLSPEC(double) integrator_get_stepzero(IntegratorSystem *blsys);
261 /**<
262 Returns the length of the initial step user specified,
263 or 0.0 if none was set.
264 */
265
266 ASC_DLLSPEC(void) integrator_set_minstep(IntegratorSystem *blsys, double);
267 ASC_DLLSPEC(double) integrator_get_minstep(IntegratorSystem *blsys);
268 /**<
269 Returns the length of the longest allowable step.
270 or 0.0 if none was set by user.
271 */
272
273 ASC_DLLSPEC(void) integrator_set_maxstep(IntegratorSystem *blsys, double);
274 ASC_DLLSPEC(double) integrator_get_maxstep(IntegratorSystem *blsys);
275 /**<
276 Returns the length of the shortest allowable step.
277 or 0.0 if none was set by user.
278 */
279
280 ASC_DLLSPEC(void) integrator_set_maxsubsteps(IntegratorSystem *blsys, int);
281 ASC_DLLSPEC(int) integrator_get_maxsubsteps(IntegratorSystem *blsys);
282 /**<
283 Returns the most internal steps allowed between
284 two time samples, or 0 if none was set by user.
285 */
286
287 ASC_DLLSPEC(long) integrator_getnsamples(IntegratorSystem *blsys);
288 /**<
289 Returns the number of values currently stored in xsamples.
290 */
291
292 ASC_DLLSPEC(long) integrator_getcurrentstep(IntegratorSystem *blsys);
293 /**<
294 Returns the current step number (blsys->currentstep). Should be reset to
295 zero inside your solver's integrator_*_solve method.
296 */
297
298 extern double integrator_getsample(IntegratorSystem *blsys, long i);
299 /**< The following functions fetch and set parts with names specific to
300 the type definitions in ivp.lib, the ASCEND initial value problem
301 model file. They are for use by any ivp solver interface.
302 All names are given at the scope of ivp.
303
304 Returns the value stored in xsamples[i]. Will whine if
305 if xsample[i] does not exist. No, there is no way to get
306 back the pointer to the xsamples vector.
307 */
308
309 extern void integrator_setsample(IntegratorSystem *blsys, long i, double val);
310 /**<
311 Sets the value stored in xsamples[i] to val. Will whine if
312 if xsample[i] does not exist. No, there is no way to get
313 back the pointer to the xsamples vector.
314 */
315
316 extern const dim_type *integrator_getsampledim(IntegratorSystem *blsys);
317 /**<
318 Returns a pointer to the dimensionality of the samples currently
319 stored, or NULL if none are stored. Do not free this pointer: we
320 own it.
321 */
322
323 /*------------------------------------------------------------------------------
324 BACKEND INTERFACE (for use by the integrator engine)
325 */
326
327 /* Problem size parameters. */
328
329 /* Parts of type definition derivatives refinements. */
330
331 ASC_DLLSPEC(double) integrator_get_t(IntegratorSystem *blsys);
332 /**<
333 Gets value of the independent variable value from the ASCEND model
334 (retrieves the value from the ASCEND compiler's instance hierarchy)
335 */
336
337 extern void integrator_set_t(IntegratorSystem *blsys, double value);
338 /**<
339 Sets value of the independent variable in the model (inserts the value in
340 the correct place in the compiler's instance hierarchy via the var_variable
341 API).
342 */
343
344 extern double *integrator_get_y(IntegratorSystem *blsys, double *vector);
345 /**<
346 Gets the value of vector 'y' from the model (what 'y' is depends on your
347 particular integrator).
348
349 If vector given is NULL, allocates vector, which the caller then owns.
350 Vector, if given, should be IntegGet_d_neq()+1 long.
351 */
352
353 extern void integrator_set_y(IntegratorSystem *blsys, double *vector);
354 /**<
355 Sets d.y[] to values in vector.
356 */
357
358 extern double *integrator_get_ydot(IntegratorSystem *blsys, double *vector);
359 /**<
360 Returns the vector ydot (derivatives of the 'states')
361
362 If vector given is NULL, allocates vector, which the caller then owns.
363 Vector, if given, should be IntegGet_d_neq()+1 long.
364 */
365
366 extern void integrator_set_ydot(IntegratorSystem *blsys, double *vector);
367 /**<
368 Sets d.dydx[] to values in vector.
369 */
370
371 ASC_DLLSPEC(double *) integrator_get_observations(IntegratorSystem *blsys, double *vector);
372 /**<
373 Returns the vector d.obs.
374 Vector should be of sufficient length (g_intginst_d_n_obs+1).
375 If NULL is given a vector is allocated and is the caller's
376 responsibility to deallocate.
377 */
378
379 ASC_DLLSPEC(struct var_variable *) integrator_get_observed_var(IntegratorSystem *blsys, const long i);
380 /**<
381 Returns the var_variable contained in the ith position in the observed variable list.
382 */
383
384 ASC_DLLSPEC(struct var_variable *) integrator_get_independent_var(IntegratorSystem *blsys);
385 /**<
386 Return a pointer to the variable identified as the independent variable.
387 */
388
389 /*-------------------------------
390 Stuff to facilitate output to the interface
391 */
392
393 /**
394 This call should be used to get the file streams ready, output column
395 headings etc.
396 */
397 extern int integrator_output_init(IntegratorSystem *blsys);
398
399 /**
400 This call should be used to output the immediate integration results from
401 a single timestep. For an ODE solver, this means just the 'y' values but
402 perhaps not the values of the algebraic varaibles, which must be calculated
403 separately. They'll be stored in the Integ_system_t somehow (not yet
404 known how)
405
406 @return 1 on success, return 0 on failure (integration will be cancelled)
407 */
408 extern int integrator_output_write(IntegratorSystem *blsys);
409
410 /**
411 Write out the 'observed values' for the integration. In the case of ODE
412 integration, we assume that the values of all the algabraic variables are
413 also now calculated.
414 */
415 extern int integrator_output_write_obs(IntegratorSystem *blsys);
416
417 /**
418 This call will close file stream and perhaps perform some kind of
419 user notification or screen update, etc.
420 */
421 extern int integrator_output_close(IntegratorSystem *blsys);
422
423 #endif

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