1 |
#include "incidencematrix.h" |
2 |
|
3 |
#include <stdexcept> |
4 |
#include <iostream> |
5 |
using namespace std; |
6 |
|
7 |
#include "variable.h" |
8 |
#include "relation.h" |
9 |
|
10 |
extern "C"{ |
11 |
#include <utilities/ascConfig.h> |
12 |
#include <compiler/instance_enum.h> |
13 |
|
14 |
#include <solver/var.h> |
15 |
#include <solver/rel.h> |
16 |
#include <solver/discrete.h> |
17 |
#include <solver/conditional.h> |
18 |
#include <solver/logrel.h> |
19 |
#include <solver/bnd.h> |
20 |
#include <solver/mtx.h> |
21 |
#include <solver/linsol.h> |
22 |
#include <solver/linsolqr.h> |
23 |
#include <solver/slv_common.h> |
24 |
#include <solver/slv_types.h> |
25 |
#include <solver/slv_client.h> |
26 |
} |
27 |
|
28 |
IncidencePoint::IncidencePoint(const int&row, const int &col, const IncidencePointType &type) : row(row), col(col), type(type){ |
29 |
// constructor, IncidencePoint |
30 |
} |
31 |
|
32 |
IncidencePoint::IncidencePoint(const IncidencePoint &old) : row(old.row), col(old.col), type(old.type){ |
33 |
// copy ctor |
34 |
} |
35 |
|
36 |
IncidencePoint::IncidencePoint() : row(-1), col(-1), type(IM_NULL){ |
37 |
// default ctor... don't use. need this to keep swig happy for some strange reason. |
38 |
} |
39 |
|
40 |
IncidenceMatrix::IncidenceMatrix(Simulation &sim) : sim(sim){ |
41 |
// constructor |
42 |
is_built = FALSE; |
43 |
} |
44 |
|
45 |
IncidenceMatrix::~IncidenceMatrix(){ |
46 |
if(is_built){ |
47 |
free_incidence_data(&i); |
48 |
} |
49 |
} |
50 |
|
51 |
void |
52 |
IncidenceMatrix::buildPlotData(){ |
53 |
int c=-1; |
54 |
|
55 |
cerr << "BUILDPLOTDATA" << endl; |
56 |
|
57 |
slv_system_t sys = sim.getSystem(); |
58 |
|
59 |
cerr << "GOT SYSTEM DATA" << endl; |
60 |
|
61 |
if(build_incidence_data(sys,&i)) { |
62 |
cerr << "FAILED TO BUILD INCIDENCE DATA" << endl; |
63 |
free_incidence_data(&i); |
64 |
throw runtime_error("IncidenceMatrix::buildPlotData error calculating grid"); |
65 |
return; |
66 |
} |
67 |
|
68 |
for (int r=0; r < i.nprow; r++) { |
69 |
struct rel_relation *rel = i.rlist[i.pr2e[r]]; |
70 |
const struct var_variable **vp = rel_incidence_list(rel); |
71 |
|
72 |
if(rel_active(rel)){ |
73 |
int nvars = rel_n_incidences(rel); |
74 |
if(rel_included(rel)){ |
75 |
for(int v=0; v < nvars; v++ ) { |
76 |
if(var_flags(vp[v]) & VAR_SVAR) { |
77 |
int vndx = var_sindex(vp[v]); |
78 |
c = i.v2pc[vndx]; |
79 |
if (i.vfixed[vndx]) { |
80 |
data.push_back(IncidencePoint(r,c,IM_ACTIVE_FIXED)); |
81 |
}else{ |
82 |
data.push_back(IncidencePoint(r,c,IM_ACTIVE_FREE)); |
83 |
} |
84 |
} |
85 |
} |
86 |
}else{ /* hollow squares */ |
87 |
for(int v=0; v < nvars; v++ ) { |
88 |
if (var_flags(vp[v]) & VAR_SVAR) { |
89 |
int vndx = var_sindex(vp[v]); |
90 |
c = i.v2pc[vndx]; |
91 |
if (i.vfixed[vndx]) { |
92 |
data.push_back(IncidencePoint(r,c,IM_DORMANT_FIXED)); |
93 |
} else { |
94 |
data.push_back(IncidencePoint(r,c,IM_DORMANT_FREE)); |
95 |
} |
96 |
} |
97 |
} |
98 |
} |
99 |
} |
100 |
} |
101 |
|
102 |
is_built = TRUE; |
103 |
} |
104 |
|
105 |
const int & |
106 |
IncidenceMatrix::getNumRows() const{ |
107 |
return i.nprow; |
108 |
} |
109 |
|
110 |
const int & |
111 |
IncidenceMatrix::getNumCols() const{ |
112 |
return i.npcol; |
113 |
} |
114 |
|
115 |
const vector<IncidencePoint> & |
116 |
IncidenceMatrix::getIncidenceData(){ |
117 |
cerr << "GET INCIDENCE DATA" << endl; |
118 |
if(!is_built){ |
119 |
buildPlotData(); |
120 |
} |
121 |
return data; |
122 |
} |
123 |
|
124 |
const Variable |
125 |
IncidenceMatrix::getVariable(const int &col) const{ |
126 |
if(!is_built)throw runtime_error("Not built"); |
127 |
if(col < 0 || col >= getNumCols())throw range_error("Column out of range"); |
128 |
int vindex = i.pc2v[col]; |
129 |
struct var_variable *var = i.vlist[vindex]; |
130 |
|
131 |
return Variable(&sim, var); |
132 |
} |
133 |
|
134 |
const Relation |
135 |
IncidenceMatrix::getRelation(const int &row) const{ |
136 |
if(!is_built)throw runtime_error("Not built"); |
137 |
if(row < 0 || row >= getNumRows())throw range_error("Row out of range"); |
138 |
int rindex = i.pr2e[row]; |
139 |
struct rel_relation *rel = i.rlist[rindex]; |
140 |
return Relation(&sim, rel); |
141 |
} |
142 |
|
143 |
const int |
144 |
IncidenceMatrix::getBlockRow(const int &row) const{ |
145 |
if(!is_built)throw runtime_error("Not built"); |
146 |
if(row < 0 || row >= getNumRows())throw range_error("Row out of range"); |
147 |
const mtx_block_t *bb = slv_get_solvers_blocks(sim.getSystem()); |
148 |
for(int i=0; i < bb->nblocks; ++i){ |
149 |
if(row >= bb->block[i].row.low && row <= bb->block[i].row.high){ |
150 |
return i; |
151 |
} |
152 |
} |
153 |
return -1; |
154 |
} |
155 |
|
156 |
/** |
157 |
Returns location of specified block |
158 |
@param block the block number |
159 |
@return vector(ve row-low, col-low, row-high, col-high) |
160 |
*/ |
161 |
const vector<int> |
162 |
IncidenceMatrix::getBlockLocation(const int &block) const{ |
163 |
if(!is_built)throw runtime_error("Not built"); |
164 |
const mtx_block_t *bb = slv_get_solvers_blocks(sim.getSystem()); |
165 |
if(block < 0 || block >= bb->nblocks){ |
166 |
throw range_error("Invalid block number"); |
167 |
} |
168 |
vector<int> v; |
169 |
v.push_back(bb->block[block].row.low); |
170 |
v.push_back(bb->block[block].col.low); |
171 |
v.push_back(bb->block[block].row.high); |
172 |
v.push_back(bb->block[block].col.high); |
173 |
return v; |
174 |
} |
175 |
|
176 |
const vector<Variable> |
177 |
IncidenceMatrix::getBlockVars(const int &block){ |
178 |
if(!is_built){ |
179 |
buildPlotData(); |
180 |
} |
181 |
vector<Variable> v; |
182 |
const mtx_block_t *bb = slv_get_solvers_blocks(sim.getSystem()); |
183 |
if(block < 0 || block >= bb->nblocks){ |
184 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Block out of range"); |
185 |
return v; |
186 |
} |
187 |
int low = bb->block[block].col.low; |
188 |
int high = bb->block[block].col.high; |
189 |
for(int j=low; j<=high; ++j){ |
190 |
v.push_back(getVariable(j)); |
191 |
} |
192 |
return v; |
193 |
} |
194 |
|
195 |
const vector<Relation> |
196 |
IncidenceMatrix::getBlockRels(const int &block){ |
197 |
CONSOLE_DEBUG("..."); |
198 |
if(!is_built){ |
199 |
buildPlotData(); |
200 |
} |
201 |
vector<Relation> v; |
202 |
const mtx_block_t *bb = slv_get_solvers_blocks(sim.getSystem()); |
203 |
if(block < 0 || block >= bb->nblocks){ |
204 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Block out of range"); |
205 |
return v; |
206 |
} |
207 |
int low = bb->block[block].row.low; |
208 |
int high = bb->block[block].row.high; |
209 |
for(int j=low; j<=high; ++j){ |
210 |
v.push_back(getRelation(j)); |
211 |
} |
212 |
CONSOLE_DEBUG("..."); |
213 |
return v; |
214 |
} |
215 |
|
216 |
const int |
217 |
IncidenceMatrix::getNumBlocks(){ |
218 |
if(!is_built){ |
219 |
buildPlotData(); |
220 |
} |
221 |
const mtx_block_t *bb = slv_get_solvers_blocks(sim.getSystem()); |
222 |
return bb->nblocks; |
223 |
} |
224 |
|
225 |
#ifdef ASC_WITH_MFGRAPH |
226 |
/** |
227 |
Return an MFGRAPH object containing the incidence information for the |
228 |
local block. It is intended that this should be passed back to the python |
229 |
layer as a SWIGed object that Python can decode. It might prove easier |
230 |
just to write a .dot from this C++ layer however. |
231 |
*/ |
232 |
mfg::DrawGraph |
233 |
IncidenceMatrix::getBlockGraph(const int &block) const{ |
234 |
using namespace mfg; |
235 |
|
236 |
DrawGraph g; |
237 |
Node *a = g.CreateNode(); |
238 |
Node *x = g.CreateNode(); |
239 |
Edge *ax = g.CreateEdge(a, x); |
240 |
ax->Attribs()["color"] = "red"; |
241 |
Subgraph *s = g.CreateSubgraph(); |
242 |
Node *b = g.CreateNode(s); |
243 |
b->Attribs()["style"] = "filled"; |
244 |
b->Attribs()["fillcolor"] = "green"; |
245 |
Edge *ab = g.CreateEdge(a, b); |
246 |
return g; |
247 |
} |
248 |
#endif |
249 |
|