/[ascend]/trunk/base/generic/system/jacobian.c
ViewVC logotype

Contents of /trunk/base/generic/system/jacobian.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1347 - (show annotations) (download) (as text)
Tue Mar 13 06:40:06 2007 UTC (13 years, 2 months ago) by jpye
File MIME type: text/x-csrc
File size: 4470 byte(s)
Added routine for outputting system as a DOT graph.
Changed an 'assert' to an 'asc_assert' in logrelation.c.
Changed logical relation in sequence.a4c (was causing a crash).
Added Simulation::write(FILE,char*) for outputting *stuff* from a simulation.
1 /* ASCEND modelling environment
2 Copyright (C) 2007 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 */
19 #include "jacobian.h"
20 #include "relman.h"
21 #include "slv_client.h"
22 #include <utilities/ascMalloc.h>
23 #include <utilities/ascPanic.h>
24 #include <general/mathmacros.h>
25
26 /* #define JACOBIAN_DEBUG */
27
28 int system_jacobian(slv_system_t sys
29 , const rel_filter_t *rfilter, const var_filter_t *vfilter, const int safe
30 , struct SystemJacobianStruct *sysjac
31 ){
32 int i,j,n,nsr,nsv,nr,nv;
33 struct var_variable **svars;
34 struct rel_relation **srels;
35 int *vartocol;
36 int res, err=0;
37 double *derivvals;
38 int *derivvars;
39 mtx_coord_t coord;
40 char *relname;
41 #ifdef JACOBIAN_DEBUG
42 char *varname;
43 #endif
44
45 /* first count the rels */
46 nsr = slv_get_num_solvers_rels(sys);
47 nr = slv_count_solvers_rels(sys,rfilter);
48 srels = slv_get_solvers_rel_list(sys);
49
50 /* and vars */
51 nsv = slv_get_num_solvers_vars(sys);
52 nv = slv_count_solvers_vars(sys,vfilter);
53 svars = slv_get_solvers_var_list(sys);
54
55 /* allocate space for lists */
56 sysjac->vars = ASC_NEW_ARRAY(struct var_variable*, nv);
57 sysjac->rels = ASC_NEW_ARRAY(struct rel_relation*, nr);
58 vartocol = ASC_NEW_ARRAY(int,nsv);
59
60 /* now create a the lists of vars and rels, and temp mapping array */
61 n = 0;
62 for(i=0;i<nsr;++i){
63 if(rel_apply_filter(srels[i],rfilter)){
64 sysjac->rels[n++] = srels[i];
65 }
66 }
67 asc_assert(n==nr);
68
69 n = 0;
70 for(i=0;i<nsv;++i){
71 if(var_apply_filter(svars[i],vfilter)){
72 vartocol[i]=n;
73 sysjac->vars[n++] = svars[i];
74 }else{
75 vartocol[i]=-1;
76 }
77 }
78 asc_assert(n==nv);
79
80 CONSOLE_DEBUG("nr = %d",nr);
81 CONSOLE_DEBUG("nv = %d",nv);
82
83 derivvals = ASC_NEW_ARRAY(double,nv);
84 derivvars = ASC_NEW_ARRAY(int,nv);
85
86 /* now create a matrix */
87 sysjac->M = mtx_create();
88 mtx_set_order(sysjac->M, MAX(nv,nr));
89
90 for(i=0;i<nr;++i){
91 #ifdef JACOBIAN_DEBUG
92 relname = rel_make_name(sys,sysjac->rels[i]);
93 CONSOLE_DEBUG("rel '%s'",relname);
94 ASC_FREE(relname);
95 #endif
96
97 res = relman_diff2(sysjac->rels[i], vfilter, derivvals, derivvars, &n, safe);
98 for(j=0;j<n;++j){
99 #ifdef JACOBIAN_DEBUG
100 asc_assert(var_apply_filter(svars[derivvars[j]],vfilter));
101 varname = var_make_name(sys,svars[derivvars[j]]);
102 CONSOLE_DEBUG("var '%s' (var_deriv = %d)",varname,var_deriv(svars[derivvars[j]]));
103 ASC_FREE(varname);
104 #endif
105 mtx_set_value(sysjac->M,mtx_coord(&coord,i,vartocol[derivvars[j]]),derivvals[j]);
106 }
107 if(res){
108 relname = rel_make_name(sys,sysjac->rels[i]);
109 ERROR_REPORTER_HERE(ASC_PROG_ERR,"Error calculating derivatives for relation '%s'",relname);
110 ASC_FREE(relname);
111 err = 1;
112 }
113 }
114
115 sysjac->n_rels = nr;
116 sysjac->n_vars = nv;
117
118 ASC_FREE(derivvals);
119 ASC_FREE(derivvars);
120
121 return err;
122 }
123
124 /* these are the filters that will partitiong our DAE system */
125
126
127 const rel_filter_t system_rfilter_algeb = {
128 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE | REL_DIFFERENTIAL,
129 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE | 0
130 };
131
132 const rel_filter_t system_rfilter_diff = {
133 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE | REL_DIFFERENTIAL,
134 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE | REL_DIFFERENTIAL
135 };
136
137 const rel_filter_t system_rfilter_all = {
138 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE,
139 REL_INCLUDED | REL_EQUALITY | REL_ACTIVE
140 };
141
142
143 const var_filter_t system_vfilter_algeb = {
144 VAR_SVAR | VAR_ACTIVE | VAR_FIXED | VAR_DERIV | VAR_DIFF,
145 VAR_SVAR | VAR_ACTIVE | 0 | 0 | 0
146 };
147
148 const var_filter_t system_vfilter_diff = {
149 VAR_SVAR | VAR_ACTIVE | VAR_FIXED | VAR_DERIV | VAR_DIFF,
150 VAR_SVAR | VAR_ACTIVE | 0 | 0 | VAR_DIFF
151 };
152
153 const var_filter_t system_vfilter_deriv = {
154 VAR_SVAR | VAR_ACTIVE | VAR_FIXED | VAR_DERIV ,
155 VAR_SVAR | VAR_ACTIVE | 0 | VAR_DERIV
156 };
157
158 const var_filter_t system_vfilter_nonderiv = {
159 VAR_SVAR | VAR_ACTIVE | VAR_FIXED | VAR_DERIV ,
160 VAR_SVAR | VAR_ACTIVE | 0 | 0
161 };

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