/[ascend]/trunk/pygtk/simulation.cpp
ViewVC logotype

Annotation of /trunk/pygtk/simulation.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 841 - (hide annotations) (download) (as text)
Sun Sep 10 06:13:44 2006 UTC (13 years, 5 months ago) by johnpye
File MIME type: text/x-c++src
File size: 21177 byte(s)
Fixed crash in simulation.cpp when freeable vars are not returned.
1 johnpye 669 /* 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     */
19 johnpye 132 #include <iostream>
20     #include <iomanip>
21     #include <stdexcept>
22     #include <sstream>
23     using namespace std;
24    
25 johnpye 506 #include "config.h"
26    
27 johnpye 132 extern "C"{
28     #include <utilities/ascConfig.h>
29 johnpye 506 #include <utilities/error.h>
30 johnpye 132 #include <utilities/ascSignal.h>
31     #include <utilities/ascMalloc.h>
32     #include <general/dstring.h>
33 johnpye 196 #include <general/tm_time.h>
34 johnpye 132 #include <compiler/instance_enum.h>
35     #include <compiler/fractions.h>
36     #include <compiler/compiler.h>
37     #include <compiler/dimen.h>
38     #include <compiler/symtab.h>
39     #include <compiler/instance_io.h>
40     #include <compiler/instantiate.h>
41     #include <compiler/bintoken.h>
42     #include <compiler/instance_enum.h>
43     #include <compiler/instquery.h>
44     #include <compiler/check.h>
45     #include <compiler/name.h>
46     #include <compiler/pending.h>
47    
48     #include <utilities/readln.h>
49     #include <solver/mtx.h>
50     #include <solver/slv_types.h>
51     #include <solver/var.h>
52     #include <solver/rel.h>
53     #include <solver/discrete.h>
54     #include <solver/conditional.h>
55     #include <solver/logrel.h>
56     #include <solver/bnd.h>
57     #include <solver/calc.h>
58     #include <solver/relman.h>
59     #include <solver/slv_common.h>
60     #include <solver/linsol.h>
61     #include <solver/linsolqr.h>
62     #include <solver/slv_client.h>
63     #include <solver/system.h>
64     #include <solver/slv_interface.h>
65     #include <solver/slvDOF.h>
66     #include <solver/slv3.h>
67     #include <solver/slv_stdcalls.h>
68 johnpye 328 #include <solver/slv_server.h>
69 johnpye 132 }
70    
71     #include "simulation.h"
72     #include "solver.h"
73 johnpye 208 #include "solverparameters.h"
74 johnpye 132 #include "name.h"
75 johnpye 233 #include "incidencematrix.h"
76 johnpye 237 #include "variable.h"
77 johnpye 307 #include "solverstatus.h"
78 johnpye 310 #include "solverreporter.h"
79 johnpye 132
80     /**
81     Create an instance of a type (call compiler etc)
82    
83     @TODO fix mutex on compile command filenames
84     */
85     Simulation::Simulation(Instance *i, const SymChar &name) : Instanc(i, name), simroot(GetSimulationRoot(i),SymChar("simroot")){
86 johnpye 153 is_built = false;
87 johnpye 132 // Create an Instance object for the 'simulation root' (we'll call
88 johnpye 190 // it the 'simulation model') and it can be fetched using 'getModel()'
89 johnpye 132 // any time later.
90     //simroot = Instanc(GetSimulationRoot(i),name);
91     }
92    
93 johnpye 164 Simulation::Simulation(const Simulation &old) : Instanc(old), simroot(old.simroot){
94     is_built = old.is_built;
95     sys = old.sys;
96     bin_srcname = old.bin_srcname;
97     bin_objname = old.bin_objname;
98     bin_libname = old.bin_libname;
99     bin_cmd = old.bin_cmd;
100 johnpye 190 bin_rm = old.bin_rm;
101 johnpye 775 sing = NULL;
102 johnpye 164 }
103    
104     Simulation::~Simulation(){
105     //CONSOLE_DEBUG("Deleting simulation %s", getName().toString());
106     }
107    
108 johnpye 132 Instanc &
109     Simulation::getModel(){
110     if(!simroot.getInternalType()){
111     throw runtime_error("Simulation::getModel: simroot.getInternalType()is NULL");
112     }
113     return simroot;
114     }
115    
116 johnpye 775
117     slv_system_structure *
118     Simulation::getSystem(){
119     if(!sys)throw runtime_error("Can't getSystem: simulation not yet built");
120     return sys;
121 johnpye 132 }
122    
123 johnpye 295
124 johnpye 775 const string
125     Simulation::getInstanceName(const Instanc &i) const{
126     char *n;
127     n = WriteInstanceNameString(i.getInternalType(),simroot.getInternalType());
128     string s(n);
129     ascfree(n);
130     return s;
131     }
132 johnpye 295
133 johnpye 775 const int
134     Simulation::getNumVars(){
135     return slv_get_num_solvers_vars(getSystem());
136 johnpye 295 }
137    
138    
139 johnpye 775 void
140     Simulation::write(){
141     simroot.write();
142 johnpye 295 }
143    
144 johnpye 775 //------------------------------------------------------------------------------
145     // RUNNING MODEL 'METHODS'
146 johnpye 772
147 johnpye 295 void
148 johnpye 132 Simulation::run(const Method &method){
149 johnpye 774 Instanc &model = getModel();
150     this->run(method,model);
151     }
152    
153     void
154     Simulation::run(const Method &method, Instanc &model){
155    
156 johnpye 132 cerr << "RUNNING PROCEDURE " << method.getName() << endl;
157     Nam name = Nam(method.getSym());
158 johnpye 154 //cerr << "CREATED NAME '" << name.getName() << "'" << endl;
159 johnpye 132 Proc_enum pe;
160     pe = Initialize(
161 johnpye 774 &*(model.getInternalType()) ,name.getInternalType(), "__not_named__"
162 johnpye 132 ,ASCERR
163     ,0, NULL, NULL
164     );
165    
166     if(pe == Proc_all_ok){
167 johnpye 190 ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Method '%s' was run (check above for errors)\n",method.getName());
168 johnpye 132 //cerr << "METHOD " << method.getName() << " COMPLETED OK" << endl;
169     }else{
170     stringstream ss;
171     ss << "Simulation::run: Method '" << method.getName() << "' returned error: ";
172     switch(pe){
173     case Proc_CallOK: ss << "Call OK"; break;
174     case Proc_CallError: ss << "Error occurred in call"; break;
175     case Proc_CallReturn: ss << "Request that caller return (OK)"; break;
176     case Proc_CallBreak: ss << "Break out of enclosing loop"; break;
177     case Proc_CallContinue: ss << "Skip to next iteration"; break;
178    
179     case Proc_break: ss << "Break"; break;
180     case Proc_continue: ss << "Continue"; break;
181     case Proc_fallthru: ss << "Fall-through"; break;
182     case Proc_return: ss << "Return"; break;
183     case Proc_stop: ss << "Stop"; break;
184     case Proc_stack_exceeded: ss << "Stack exceeded"; break;
185     case Proc_stack_exceeded_this_frame: ss << "Stack exceeded this frame"; break;
186     case Proc_case_matched: ss << "Case matched"; break;
187     case Proc_case_unmatched: ss << "Case unmatched"; break;
188    
189     case Proc_case_undefined_value: ss << "Undefined value in case"; break;
190     case Proc_case_boolean_mismatch: ss << "Boolean mismatch in case"; break;
191     case Proc_case_integer_mismatch: ss << "Integer mismatch in case"; break;
192     case Proc_case_symbol_mismatch: ss << "Symbol mismatch in case"; break;
193     case Proc_case_wrong_index: ss << "Wrong index in case"; break;
194     case Proc_case_wrong_value: ss << "Wrong value in case"; break;
195     case Proc_case_extra_values: ss << "Extra values in case"; break;
196     case Proc_bad_statement: ss << "Bad statement"; break;
197     case Proc_bad_name: ss << "Bad name"; break;
198     case Proc_for_duplicate_index: ss << "Duplicate index"; break;
199     case Proc_for_set_err: ss << "For set error"; break;
200     case Proc_for_not_set: ss << "For not set"; break;
201     case Proc_illegal_name_use: ss << "Illegal name use"; break;
202     case Proc_name_not_found: ss << "Name not found"; break;
203     case Proc_instance_not_found: ss << "Instance not found"; break;
204     case Proc_type_not_found: ss << "Type not found"; break;
205     case Proc_illegal_type_use: ss << "Illegal use"; break;
206     case Proc_proc_not_found: ss << "Method not found"; break;
207     case Proc_if_expr_error_typeconflict: ss << "Type conflict in 'if' expression"; break;
208     case Proc_if_expr_error_nameunfound: ss << "Name not found in 'if' expression"; break;
209     case Proc_if_expr_error_incorrectname: ss << "Incorrect name in 'if' expression"; break;
210     case Proc_if_expr_error_undefinedvalue: ss << "Undefined value in 'if' expression"; break;
211     case Proc_if_expr_error_dimensionconflict: ss << "Dimension conflict in 'if' expression"; break;
212     case Proc_if_expr_error_emptychoice: ss << "Empty choice in 'if' expression"; break;
213     case Proc_if_expr_error_emptyintersection: ss << "Empty intersection in 'if' expression"; break;
214     case Proc_if_expr_error_confused: ss << "Confused in 'if' expression"; break;
215     case Proc_if_real_expr: ss << "Real-valued result in 'if' expression"; break;
216     case Proc_if_integer_expr: ss << "Integeter-valued result in 'if' expression"; break;
217     case Proc_if_symbol_expr: ss << "Symbol-valued result in 'if' expression"; break;
218     case Proc_if_set_expr: ss << "Set-valued result in 'if' expression"; break;
219     case Proc_if_not_logical: ss << "If expression is not logical"; break;
220     case Proc_user_interrupt: ss << "User interrupt"; break;
221     case Proc_infinite_loop: ss << "Infinite loop"; break;
222     case Proc_declarative_constant_assignment: ss << "Declarative constant assignment"; break;
223     case Proc_nonsense_assignment: ss << "Nonsense assginment (bogus)"; break;
224     case Proc_nonconsistent_assignment: ss << "Inconsistent assignment"; break;
225     case Proc_nonatom_assignment: ss << "Non-atom assignment"; break;
226     case Proc_nonboolean_assignment: ss << "Non-boolean assignment"; break;
227     case Proc_noninteger_assignment: ss << "Non-integer assignment"; break;
228     case Proc_nonreal_assignment: ss << "Non-real assignment"; break;
229     case Proc_nonsymbol_assignment: ss << "Non-symbol assignment"; break;
230     case Proc_lhs_error: ss << "Left-hand-side error"; break;
231     case Proc_rhs_error: ss << "Right-hand-side error"; break;
232     case Proc_unknown_error: ss << "Unknown error"; break;
233     default:
234     ss << "Invalid error code";
235     }
236    
237     ss << " (" << int(pe) << ")";
238     throw runtime_error(ss.str());
239     }
240     }
241    
242 johnpye 775 //-----------------------------------------------------------------------------
243     // CHECKING METHODS
244    
245 johnpye 772 /**
246 johnpye 775 Check that all the analysis went OK: solver lists are all there, etc...?
247    
248     Can't return anything here because of limitations in the C API
249 johnpye 772 */
250 johnpye 775 void
251     Simulation::checkInstance(){
252     cerr << "CHECKING SIMULATION INSTANCE" << endl;
253 johnpye 132 Instance *i1 = getModel().getInternalType();
254 johnpye 164 CheckInstance(stderr, &*i1);
255 johnpye 775 cerr << "DONE CHECKING INSTANCE" << endl;
256     }
257 johnpye 772
258 johnpye 775 /**
259     @return 1 = underspecified, 2 = square, 3 = structurally singular, 4 = overspecified
260     */
261     enum StructuralStatus
262     Simulation::checkDoF() const{
263     cerr << "CHECKING DOF..." << endl;
264     int dof, status;
265     if(!sys){
266     throw runtime_error("System not yet built");
267     }
268     slvDOF_status(sys, &status, &dof);
269     switch(status){
270     case ASCXX_DOF_UNDERSPECIFIED:
271     case ASCXX_DOF_SQUARE:
272     case ASCXX_DOF_OVERSPECIFIED:
273     case ASCXX_DOF_STRUCT_SINGULAR:
274     return (enum StructuralStatus)status;
275     case 5:
276     throw runtime_error("Unable to resolve degrees of freedom"); break;
277     default:
278     throw runtime_error("Invalid return status from slvDOF_status");
279     }
280 johnpye 190 }
281 johnpye 132
282 johnpye 775 /**
283     Check consistency
284    
285     @TODO what is the difference between this and checkStructuralSingularity?
286    
287     @return list of freeable variables. List will be empty if sys is consistent.
288     */
289     vector<Variable>
290     Simulation::getFreeableVariables(){
291     vector<Variable> v;
292    
293     cerr << "CHECKING CONSISTENCY..." << endl;
294 johnpye 841 int *fixedarrayptr=NULL;
295 johnpye 775
296 johnpye 841 if(!sys){
297     throw runtime_error("System not yet built");
298     }
299    
300 johnpye 775 int res = consistency_analysis(sys, &fixedarrayptr);
301    
302     if(res==1){
303     cerr << "STRUCTURALLY CONSISTENT" << endl;
304     }else{
305 johnpye 841 if(fixedarrayptr ==NULL){
306     ERROR_REPORTER_HERE(ASC_USER_ERROR,"STRUCTURALLY INCONSISTENT");
307     throw runtime_error("Invalid constistency analysis result returned!");
308     }
309    
310     struct var_variable **vp = slv_get_master_var_list(sys);
311 johnpye 775 for(int i=0; fixedarrayptr[i]!=-1; ++i){
312     v.push_back( Variable(this, vp[fixedarrayptr[i]]) );
313     }
314     }
315     return v;
316     }
317    
318     /** Returns TRUE if all is OK (not singular) */
319     bool
320     Simulation::checkStructuralSingularity(){
321     int *vil;
322     int *ril;
323     int *fil;
324    
325     cerr << "RETREIVING slfDOF_structsing INFO" << endl;
326    
327     int res = slvDOF_structsing(sys, mtx_FIRST, &vil, &ril, &fil);
328     struct var_variable **varlist = slv_get_solvers_var_list(sys);
329     struct rel_relation **rellist = slv_get_solvers_rel_list(sys);
330    
331     if(this->sing){
332     cerr << "DELETING OLD SINGULATING INFO" << endl;
333     delete this->sing;
334     this->sing = NULL;
335     }
336    
337     if(res==1){
338     CONSOLE_DEBUG("processing singularity data...");
339     sing = new SingularityInfo();
340    
341     // pull in the lists of vars and rels, and the freeable vars:
342     for(int i=0; ril[i]!=-1; ++i){
343     sing->rels.push_back( Relation(this, rellist[ril[i]]) );
344     }
345    
346     for(int i=0; vil[i]!=-1; ++i){
347     sing->vars.push_back( Variable(this, varlist[vil[i]]) );
348     }
349    
350     for(int i=0; fil[i]!=-1; ++i){
351     sing->freeablevars.push_back( Variable(this, varlist[fil[i]]) );
352     }
353    
354     // we're done with those lists now
355     ASC_FREE(vil);
356     ASC_FREE(ril);
357     ASC_FREE(fil);
358    
359     if(sing->isSingular()){
360     CONSOLE_DEBUG("singularity found");
361     this->sing = sing;
362     return FALSE;
363     }
364     CONSOLE_DEBUG("no singularity");
365     delete sing;
366     return TRUE;
367     }else{
368     if(res==0){
369     throw runtime_error("Unable to determine singularity lists");
370     }else{
371     throw runtime_error("Invalid return from slvDOF_structsing.");
372     }
373     }
374     }
375    
376     /**
377     If the checkStructuralSingularity analysis has been done,
378     this funciton will let you access the SingularityInfo data that was
379     stored.
380     */
381     const SingularityInfo &
382     Simulation::getSingularityInfo() const{
383     if(sing==NULL){
384     throw runtime_error("No singularity info present");
385     }
386     return *sing;
387     }
388    
389     //------------------------------------------
390     // ASSIGNING SOLVER TO SIMULATION
391    
392 johnpye 132 void
393 johnpye 775 Simulation::setSolver(Solver &solver){
394     cerr << "SETTING SOLVER ON SIMULATION TO " << solver.getName() << endl;
395    
396     if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");
397     // Update the solver object because sometimes an alternative solver can be returned, apparently.
398    
399     int selected = slv_select_solver(sys, solver.getIndex());
400     //cerr << "Simulation::setSolver: slv_select_solver returned " << selected << endl;
401    
402     if(selected<0){
403     ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Failed to select solver");
404     throw runtime_error("Failed to select solver");
405     }
406    
407     if(selected!=solver.getIndex()){
408     solver = Solver(slv_solver_name(selected));
409     ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Substitute solver '%s' (index %d) selected.\n", solver.getName().c_str(), selected);
410     }
411    
412     if( slv_eligible_solver(sys) <= 0){
413     ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Inelegible solver '%s'", solver.getName().c_str() );
414     throw runtime_error("Inelegible solver");
415     }
416     }
417    
418     const Solver
419     Simulation::getSolver() const{
420     int index = slv_get_selected_solver(sys);
421     //cerr << "Simulation::getSolver: index = " << index << endl;
422     if(index<0)throw runtime_error("No solver selected");
423    
424     return Solver(slv_solver_name(index));
425     }
426    
427     //------------------------------------------------------------------------------
428     // BUILD THE SYSTEM (SEND IT TO THE SOLVER)
429    
430     /**
431     Build the system (send it to the solver)
432     */
433     void
434 johnpye 132 Simulation::build(){
435     cerr << "BUILDING SIMULATION..." << endl;
436 johnpye 164 Instance *i1 = getModel().getInternalType();
437     sys = system_build(&*i1);
438     if(!sys){
439 johnpye 132 throw runtime_error("Unable to build system");
440     }
441 johnpye 153 is_built = true;
442 johnpye 132 cerr << "...DONE BUILDING" << endl;
443     }
444    
445 johnpye 775
446     //------------------------------------------------------------------------------
447     // SOLVER CONFIGURATION PARAMETERS
448    
449     /**
450     Get solver parameters struct wrapped up as a SolverParameters class.
451     */
452     SolverParameters
453     Simulation::getSolverParameters() const{
454     if(!sys)throw runtime_error("Can't getSolverParameters: Simulation system has not been built yet.");
455    
456     slv_parameters_t p;
457     slv_get_parameters(sys,&p);
458     return SolverParameters(p);
459     }
460    
461     /**
462     Update the solver parameters by passing a new set back
463     */
464     void
465     Simulation::setSolverParameters(SolverParameters &P){
466     if(!sys)throw runtime_error("Can't set solver parameters: simulation has not been built yet.");
467     slv_set_parameters(sys, &(P.getInternalType()));
468     }
469    
470     //------------------------------------------------------------------------------
471     // PRE-SOLVE DIAGNOSTICS
472    
473     /**
474     Get a list of variables to fix to make an underspecified system
475     become square. Also seems to return stuff when you have a structurally
476     singuler system.
477     */
478 johnpye 132 vector<Variable>
479     Simulation::getFixableVariables(){
480     cerr << "GETTING FIXABLE VARIABLES..." << endl;
481     vector<Variable> vars;
482    
483 johnpye 164 if(!sys){
484 johnpye 132 throw runtime_error("Simulation system not yet built");
485     }
486    
487 johnpye 316 int32 *vip; /** TODO ensure 32 bit integers are used */
488 johnpye 132
489     // Get IDs of elegible variables in array at vip...
490 johnpye 316 if(!slvDOF_eligible(sys,&vip)){
491 johnpye 190 ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"No fixable variables found.");
492 johnpye 132 }else{
493     struct var_variable **vp = slv_get_solvers_var_list(sys);
494    
495     if(vp==NULL){
496     throw runtime_error("Simulation variable list is null");
497     }
498 johnpye 190
499 johnpye 132 // iterate through this list until we find a -1:
500     int i=0;
501 johnpye 316 int var_index = vip[i];
502 johnpye 132 while(var_index >= 0){
503     struct var_variable *var = vp[var_index];
504 johnpye 237 vars.push_back( Variable(this, var) );
505 johnpye 132 ++i;
506 johnpye 316 var_index = vip[i];
507 johnpye 132 }
508 johnpye 190 ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"Found %d fixable variables.",i);
509 johnpye 316 ascfree(vip);
510 johnpye 132 }
511    
512     return vars;
513     }
514    
515 johnpye 775 /**
516     Get the list of variables near their bounds. Helps to indentify why
517     you might be having non-convergence problems.
518     */
519 johnpye 328 vector<Variable>
520     Simulation::getVariablesNearBounds(const double &epsilon){
521     cerr << "GETTING VARIABLES NEAR BOUNDS..." << endl;
522     vector<Variable> vars;
523    
524     if(!sys){
525     throw runtime_error("Simulation system not yet built");
526     }
527    
528     int *vip;
529     if(slv_near_bounds(sys,epsilon,&vip)){
530     struct var_variable **vp = slv_get_solvers_var_list(sys);
531     struct var_variable *var;
532     cerr << "VARS FOUND NEAR BOUNDS" << endl;
533     int nlow = vip[0];
534     int nhigh = vip[1];
535     int lim1 = 2 + nlow;
536     for(int i=2; i<lim1; ++i){
537     var = vp[vip[i]];
538     char *var_name = var_make_name(sys,var);
539     cerr << "AT LOWER BOUND: " << var_name << endl;
540     ascfree(var_name);
541     vars.push_back(Variable(this,var));
542     };
543     int lim2 = lim1 + nhigh;
544     for(int i=lim1; i<lim2; ++i){
545     var = vp[vip[i]];
546     char *var_name = var_make_name(sys,var);
547     cerr << "AT UPPER BOUND: " << var_name << endl;
548     ascfree(var_name);
549     vars.push_back(Variable(this,var));
550     }
551     }
552     ascfree(vip);
553     return vars;
554     }
555 johnpye 502
556 johnpye 775
557     bool
558     SingularityInfo::isSingular() const{
559     if(vars.size()||rels.size()){
560     return true;
561     }
562     return false;
563     }
564    
565     //------------------------------------------------------------------------------
566     // SOLVING
567    
568     /**
569     Solve the system through to convergence. This function is hardwired with
570     a maximum of 1000 iterations, but will interrupt itself when the 'stop'
571     condition comes back from the SolverReporter.
572     */
573 johnpye 132 void
574 johnpye 310 Simulation::solve(Solver solver, SolverReporter &reporter){
575 johnpye 153 if(!is_built){
576     throw runtime_error("Simulation::solver: simulation is not yet built, can't start solving.");
577     }
578    
579 johnpye 132 cerr << "SIMULATION::SOLVE STARTING..." << endl;
580     enum inst_t k = getModel().getKind();
581     if(k!=MODEL_INST)throw runtime_error("Can't solve: not an instance of type MODEL_INST");
582 johnpye 190
583 johnpye 164 Instance *i1 = getInternalType();
584     int npend = NumberPendingInstances(&*i1);
585 johnpye 190 if(npend)throw runtime_error("Can't solve: There are still %d pending instances");
586 johnpye 132
587     if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");
588 johnpye 190
589 johnpye 132 cerr << "SIMULATION::SOLVE: SET SOLVER..." << endl;
590     setSolver(solver);
591 johnpye 190
592    
593 johnpye 207 cerr << "PRESOLVING SYSTEM...";
594 johnpye 132 slv_presolve(sys);
595 johnpye 207 cerr << "DONE" << endl;
596 johnpye 190
597 johnpye 132 cerr << "SOLVING SYSTEM..." << endl;
598 johnpye 196 // Add some stuff here for cleverer iteration....
599 johnpye 207 unsigned niter = 1000;
600 johnpye 714 //double updateinterval = 0.02;
601 johnpye 132
602 johnpye 196 double starttime = tm_cpu_time();
603 johnpye 714 //double lastupdate = starttime;
604 johnpye 307 SolverStatus status;
605 johnpye 714 //int solved_vars=0;
606 johnpye 196 bool stop=false;
607    
608 johnpye 307 status.getSimulationStatus(*this);
609 johnpye 314 reporter.report(&status);
610 johnpye 307
611 johnpye 746 for(unsigned iter = 1; iter <= niter && !stop; ++iter){
612 johnpye 285
613 johnpye 307 if(status.isReadyToSolve()){
614 johnpye 196 slv_iterate(sys);
615     }
616 johnpye 307
617     status.getSimulationStatus(*this);
618 johnpye 502
619 johnpye 317 if(reporter.report(&status)){
620     stop = true;
621     }
622 johnpye 196 }
623 johnpye 317
624 johnpye 230 double elapsed = tm_cpu_time() - starttime;
625 johnpye 196
626 johnpye 319
627 johnpye 307 activeblock = status.getCurrentBlockNum();
628 johnpye 285
629 johnpye 319 reporter.finalise(&status);
630    
631 johnpye 322 // Just a little bit of console output:
632    
633 johnpye 307 if(status.isOK()){
634 johnpye 260 cerr << "... SOLVED, STATUS OK" << endl;
635 johnpye 132 }else{
636 johnpye 321 cerr << "... SOLVER FAILED" << endl;
637 johnpye 132 }
638    
639 johnpye 307 cerr << "SOLVER PERFORMED " << status.getIterationNum() << " ITERATIONS IN " << elapsed << "s" << endl;
640 johnpye 132 }
641    
642 johnpye 775 //------------------------------------------------------------------------------
643     // POST-SOLVE DIAGNOSTICS
644 johnpye 132
645 johnpye 775 const int
646     Simulation::getActiveBlock() const{
647     return activeblock;
648 johnpye 132 }
649    
650 johnpye 208 /**
651 johnpye 775 Return an IncidenceMatrix built from the current state of the solver system.
652 johnpye 132
653 johnpye 775 This will actually return something meaningful even before solve.
654 johnpye 225 */
655 johnpye 233 IncidenceMatrix
656     Simulation::getIncidenceMatrix(){
657     return IncidenceMatrix(*this);
658     }
659 johnpye 252
660 johnpye 775 /**
661     This function looks at all the variables in the solve's list and updates
662     the variable status for the corresponding instances.
663 johnpye 317
664 johnpye 775 It does this by using the 'interface pointer' in the Instance, see
665     the C-API function GetInterfacePtr.
666 johnpye 317
667 johnpye 775 This is used to display visually which variables have been solved, which
668     ones have not yet been attempted, and which ones were active when the solver
669     failed (ASCXX_VAR_ACTIVE).
670     */
671 johnpye 255 void
672     Simulation::processVarStatus(){
673    
674     // this is a cheap function call:
675     const mtx_block_t *bb = slv_get_solvers_blocks(getSystem());
676 johnpye 297
677 johnpye 255 var_variable **vlist = slv_get_solvers_var_list(getSystem());
678     int nvars = slv_get_num_solvers_vars(getSystem());
679    
680     slv_status_t status;
681     slv_get_status(getSystem(), &status);
682    
683 johnpye 297 if(status.block.number_of == 0){
684     cerr << "Variable statuses can't be set: block structure not yet determined." << endl;
685     return;
686     }
687    
688 johnpye 255 int activeblock = status.block.current_block;
689     int low = bb->block[activeblock].col.low;
690     int high = bb->block[activeblock].col.high;
691 johnpye 258 bool allsolved = status.converged;
692 johnpye 255 for(int c=0; c < nvars; ++c){
693     var_variable *v = vlist[c];
694     Instanc i((Instance *)var_instance(v));
695     VarStatus s = ASCXX_VAR_STATUS_UNKNOWN;
696 johnpye 258 if(i.isFixed()){
697     s = ASCXX_VAR_FIXED;
698     }else if(var_incident(v) && var_active(v)){
699     if(allsolved || c < low){
700 johnpye 255 s = ASCXX_VAR_SOLVED;
701     }else if(c <= high){
702     s = ASCXX_VAR_ACTIVE;
703     }else{
704     s = ASCXX_VAR_UNSOLVED;
705     }
706     }
707     i.setVarStatus(s);
708     }
709     }
710 johnpye 285

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