/[ascend]/branches/jacob/ascend/integrator/aww.c
ViewVC logotype

Contents of /branches/jacob/ascend/integrator/aww.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2933 - (show annotations) (download) (as text)
Mon Jun 1 20:38:42 2015 UTC (9 years, 1 month ago) by jacob
File MIME type: text/x-csrc
File size: 4754 byte(s)
Create a new branch from trunk r2932, with name 'jacob'
To be used by Jacob Shealy during GSOC2015 for development of ideal 
solution capabilities in FPROPS

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, see <http://www.gnu.org/licenses/>.
16 *//** @file
17 C-code implementation of Art Westerberg's DAE integrator
18 *//*
19 by John Pye, Dec 2006.
20 */
21
22 #include <ascend/utilities/config.h>
23 #include <ascend/general/platform.h>
24 #include <ascend/general/panic.h>
25 #include "aww.h"
26
27 const IntegratorInternals integrator_aww_internals = {
28 integrator_aww_create
29 ,integrator_aww_params_default
30 ,integrator_aww_analyse
31 ,integrator_aww_solve
32 ,NULL /* writematrixfn */
33 ,NULL /* debugfn */
34 ,integrator_aww_free
35 ,INTEG_AWW
36 ,"AWW"
37 };
38
39 typedef struct{
40 unsigned nothing;
41 } IntegratorAWWData;
42
43 /*-------------------------------------------------------------
44 PARAMETERS FOR AWW INTEGRATOR
45 */
46
47 enum aww_parameters{
48 AWW_PARAM_METHOD
49 ,AWW_PARAM_RTOL
50 ,AWW_PARAMS_SIZE
51 };
52
53 enum aww_methods{
54 AWW_BDF
55 ,AWW_AM
56 };
57
58 /**
59 Here the full set of parameters is defined, along with upper/lower bounds,
60 etc. The values are stuck into the blsys->params structure.
61
62 @return 0 on success
63 */
64 int integrator_aww_params_default(IntegratorSystem *blsys){
65 asc_assert(blsys!=NULL);
66 asc_assert(blsys->engine==INTEG_AWW);
67 slv_parameters_t *p;
68 p = &(blsys->params);
69
70 slv_destroy_parms(p);
71
72 if(p->parms==NULL){
73 CONSOLE_DEBUG("params NULL");
74 p->parms = ASC_NEW_ARRAY(struct slv_parameter, AWW_PARAMS_SIZE);
75 if(p->parms==NULL)return -1;
76 p->dynamic_parms = 1;
77 }else{
78 CONSOLE_DEBUG("params not NULL");
79 }
80
81 /* reset the number of parameters to zero so that we can check it at the end */
82 p->num_parms = 0;
83
84 slv_param_char(p,AWW_PARAM_METHOD
85 ,(SlvParameterInitChar){{"method"
86 ,"Stepping method",1
87 ,"See Art's notes, sec. 15.2"
88 }, "BDF"}, (char *[]){"BDF","AM",NULL}
89 );
90
91 slv_param_real(p,AWW_PARAM_RTOL
92 ,(SlvParameterInitReal){{"rtol"
93 ,"Scalar relative error tolerance",1
94 ,"Value of the scalar relative error tolerance."
95 }, 1e-4, 0, DBL_MAX }
96 );
97
98 asc_assert(p->num_parms == AWW_PARAMS_SIZE);
99
100 CONSOLE_DEBUG("Created %d params", p->num_parms);
101
102 return 0;
103 }
104
105 void integrator_aww_create(IntegratorSystem *blsys){
106 CONSOLE_DEBUG("ALLOCATING AWW ENGINE DATA");
107 IntegratorAWWData *enginedata;
108 enginedata = ASC_NEW(IntegratorAWWData);
109 enginedata->nothing = 0;
110 blsys->enginedata = (void *)enginedata;
111 integrator_aww_params_default(blsys);
112 }
113
114 void integrator_aww_free(void *enginedata){
115 CONSOLE_DEBUG("DELETING AWW ENGINE DATA");
116 IntegratorAWWData *d = (IntegratorAWWData *)enginedata;
117 ASC_FREE(d);
118 }
119
120 /** return 1 on success */
121 int integrator_aww_analyse(IntegratorSystem *blsys){
122 return 0;
123 }
124
125 /** return 0 on success */
126 int integrator_aww_solve(IntegratorSystem *blsys
127 , unsigned long start_index, unsigned long finish_index
128 ){
129 double initstep, maxstep, rtol;
130 unsigned maxsubsteps;
131 const char *methodname;
132
133 // solve initial point
134
135 // run valuesForInit
136 // run specifyForInit
137 // solve with QRSlv
138 // if not converged throw error
139
140 // 'DELETE SYSTEM'
141 // run valuesForStep
142 // run specifyForStep
143
144 initstep = integrator_get_stepzero(blsys);
145 maxstep = integrator_get_maxstep(blsys);
146 maxsubsteps = integrator_get_maxsubsteps(blsys);
147 rtol = SLV_PARAM_REAL(&(blsys->params),AWW_PARAM_RTOL);
148
149 methodname = SLV_PARAM_CHAR(&(blsys->params),AWW_PARAM_METHOD);
150 CONSOLE_DEBUG("ASSIGNING STEPPING METHOD '%s'",methodname);
151
152 enum aww_methods method;
153 if(strcmp(methodname,"BDF")==0)method = AWW_BDF;
154 else if(strcmp(methodname,"AM")==0)method = AWW_AM;
155 else{
156 ERROR_REPORTER_HERE(ASC_USER_ERROR,"Invalid method name '%s'",methodname);
157 }
158
159 /* initialise the integrator reporter */
160 integrator_output_init(blsys);
161
162 // while not too many solver calls {
163 // capture plotting points
164 // do something with polynomials
165 // for(i=0; i<= poly_oder; ++i){
166 // if(!in the final step){
167 // if(!have hit stopping configtion) UN 'stepX'
168 // SOLVE with QRSlv
169 // if(!converged)abort 'failed to convert
170 // }
171 // if(in the final step)abort 'STOP reached'
172 // RUN setStopConditions
173 // retrieve status of stop confitions
174 // retureve status of 'last step'
175 // }
176 // compute maximum steps for each variable
177 // }
178
179 ERROR_REPORTER_HERE(ASC_PROG_ERR,"Not implemented");
180 return 1;
181 }
182

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