pygtk/interface/simulation.cpp

#include <iostream>
#include <iomanip>
#include <stdexcept>
#include <sstream>
using namespace std;

extern "C"{
#include <utilities/ascConfig.h>
#include <utilities/ascSignal.h>
#include <utilities/ascMalloc.h>
#include <general/dstring.h>
#include <general/tm_time.h>
#include <compiler/instance_enum.h>
#include <compiler/fractions.h>
#include <compiler/compiler.h>
#include <compiler/dimen.h>
#include <compiler/symtab.h>
#include <compiler/instance_io.h>
#include <compiler/instantiate.h>
#include <compiler/bintoken.h>
#include <compiler/instance_enum.h>
#include <compiler/instquery.h>
#include <compiler/check.h>
#include <compiler/name.h>
#include <compiler/pending.h>

#include <utilities/readln.h>
#include <solver/mtx.h>
#include <solver/slv_types.h>
#include <solver/var.h>
#include <solver/rel.h>
#include <solver/discrete.h>
#include <solver/conditional.h>
#include <solver/logrel.h>
#include <solver/bnd.h>
#include <solver/calc.h>
#include <solver/relman.h>
#include <solver/slv_common.h>
#include <solver/linsol.h>
#include <solver/linsolqr.h>
#include <solver/slv_client.h>
#include <solver/system.h>
#include <solver/slv_interface.h>
#include <solver/slvDOF.h>
#include <solver/slv3.h>
#include <solver/slv_stdcalls.h>
}

#include "simulation.h"
#include "solver.h"
#include "solverparameters.h"
#include "name.h"
#include "incidencematrix.h"

/**
    Create an instance of a type (call compiler etc)

    @TODO fix mutex on compile command filenames
*/
Simulation::Simulation(Instance *i, const SymChar &name) : Instanc(i, name), simroot(GetSimulationRoot(i),SymChar("simroot")){
    is_built = false;
    // Create an Instance object for the 'simulation root' (we'll call
    // it the 'simulation model') and it can be fetched using 'getModel()'
    // any time later.
    //simroot = Instanc(GetSimulationRoot(i),name);
}

Simulation::Simulation(const Simulation &old) : Instanc(old), simroot(old.simroot){
    is_built = old.is_built;
    sys = old.sys;
    bin_srcname = old.bin_srcname;
    bin_objname = old.bin_objname;
    bin_libname = old.bin_libname;
    bin_cmd = old.bin_cmd;
    bin_rm = old.bin_rm;
}

Simulation::~Simulation(){
    //CONSOLE_DEBUG("Deleting simulation %s", getName().toString());
}

Instanc &
Simulation::getModel(){
    if(!simroot.getInternalType()){
        throw runtime_error("Simulation::getModel: simroot.getInternalType()is NULL");
    }
    return simroot;
}

void
Simulation::checkDoF() const{
        cerr << "CHECKING DOF..." << endl;
        int dof, status;
        if(!sys){
                throw runtime_error("System not yet built");
        }
        slvDOF_status(sys, &status, &dof);
        switch(status){
                case 1: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Underspecified; %d degrees of freedom",dof); break;
                case 2: ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"Square"); break;
                case 3: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Structurally singular"); break;
                case 4: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Overspecified"); break;
                case 5:
                        throw runtime_error("Unable to resolve degrees of freedom"); break;
                default:
                        throw runtime_error("Invalid return status from slvDOF_status");
        }
}

void
Simulation::run(const Method &method){
    cerr << "RUNNING PROCEDURE " << method.getName() << endl;
    Nam name = Nam(method.getSym());
    //cerr << "CREATED NAME '" << name.getName() << "'" << endl;
    Proc_enum pe;
    pe = Initialize(
        &*(getModel().getInternalType()) ,name.getInternalType(), "__not_named__"
        ,ASCERR
        ,0, NULL, NULL
    );

    if(pe == Proc_all_ok){
        ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Method '%s' was run (check above for errors)\n",method.getName());
        //cerr << "METHOD " << method.getName() << " COMPLETED OK" << endl;
    }else{
        stringstream ss;
        ss << "Simulation::run: Method '" << method.getName() << "' returned error: ";
        switch(pe){
            case Proc_CallOK: ss << "Call OK"; break;
            case Proc_CallError: ss << "Error occurred in call"; break;
            case Proc_CallReturn: ss << "Request that caller return (OK)"; break;
            case Proc_CallBreak: ss << "Break out of enclosing loop"; break;
            case Proc_CallContinue: ss << "Skip to next iteration"; break;

            case Proc_break: ss << "Break"; break;
            case Proc_continue: ss << "Continue"; break;
            case Proc_fallthru: ss << "Fall-through"; break;
            case Proc_return: ss << "Return"; break;
            case Proc_stop: ss << "Stop"; break;
            case Proc_stack_exceeded: ss << "Stack exceeded"; break;
            case Proc_stack_exceeded_this_frame: ss << "Stack exceeded this frame"; break;
            case Proc_case_matched: ss << "Case matched"; break;
            case Proc_case_unmatched: ss << "Case unmatched"; break;

            case Proc_case_undefined_value: ss << "Undefined value in case"; break;
            case Proc_case_boolean_mismatch: ss << "Boolean mismatch in case"; break;
            case Proc_case_integer_mismatch: ss << "Integer mismatch in case"; break;
            case Proc_case_symbol_mismatch: ss << "Symbol mismatch in case"; break;
            case Proc_case_wrong_index: ss << "Wrong index in case"; break;
            case Proc_case_wrong_value: ss << "Wrong value in case"; break;
            case Proc_case_extra_values: ss << "Extra values in case"; break;
            case Proc_bad_statement: ss << "Bad statement"; break;
            case Proc_bad_name: ss << "Bad name"; break;
            case Proc_for_duplicate_index: ss << "Duplicate index"; break;
            case Proc_for_set_err: ss << "For set error"; break;
            case Proc_for_not_set: ss << "For not set"; break;
            case Proc_illegal_name_use: ss << "Illegal name use"; break;
            case Proc_name_not_found: ss << "Name not found"; break;
            case Proc_instance_not_found: ss << "Instance not found"; break;
            case Proc_type_not_found: ss << "Type not found"; break;
            case Proc_illegal_type_use: ss << "Illegal use"; break;
            case Proc_proc_not_found: ss << "Method not found"; break;
            case Proc_if_expr_error_typeconflict: ss << "Type conflict in 'if' expression"; break;
            case Proc_if_expr_error_nameunfound: ss << "Name not found in 'if' expression"; break;
            case Proc_if_expr_error_incorrectname: ss << "Incorrect name in 'if' expression"; break;
            case Proc_if_expr_error_undefinedvalue: ss << "Undefined value in 'if' expression"; break;
            case Proc_if_expr_error_dimensionconflict: ss << "Dimension conflict in 'if' expression"; break;
            case Proc_if_expr_error_emptychoice: ss << "Empty choice in 'if' expression"; break;
            case Proc_if_expr_error_emptyintersection: ss << "Empty intersection in 'if' expression"; break;
            case Proc_if_expr_error_confused: ss << "Confused in 'if' expression"; break;
            case Proc_if_real_expr: ss << "Real-valued result in 'if' expression"; break;
            case Proc_if_integer_expr: ss << "Integeter-valued result in 'if' expression"; break;
            case Proc_if_symbol_expr: ss << "Symbol-valued result in 'if' expression"; break;
            case Proc_if_set_expr: ss << "Set-valued result in 'if' expression"; break;
            case Proc_if_not_logical: ss << "If expression is not logical"; break;
            case Proc_user_interrupt: ss << "User interrupt"; break;
            case Proc_infinite_loop: ss << "Infinite loop"; break;
            case Proc_declarative_constant_assignment: ss << "Declarative constant assignment"; break;
            case Proc_nonsense_assignment: ss << "Nonsense assginment (bogus)"; break;
            case Proc_nonconsistent_assignment: ss << "Inconsistent assignment"; break;
            case Proc_nonatom_assignment: ss << "Non-atom assignment"; break;
            case Proc_nonboolean_assignment: ss << "Non-boolean assignment"; break;
            case Proc_noninteger_assignment: ss << "Non-integer assignment"; break;
            case Proc_nonreal_assignment: ss << "Non-real assignment"; break;
            case Proc_nonsymbol_assignment: ss << "Non-symbol assignment"; break;
            case Proc_lhs_error: ss << "Left-hand-side error"; break;
            case Proc_rhs_error: ss << "Right-hand-side error"; break;
            case Proc_unknown_error: ss << "Unknown error"; break;
            default:
                ss << "Invalid error code";
        }


        ss << " (" << int(pe) << ")";
        throw runtime_error(ss.str());
    }
}

const bool
Simulation::check(){
    cerr << "CHECKING SIMULATION" << endl;
    Instance *i1 = getModel().getInternalType();
    CheckInstance(stderr, &*i1);
    cerr << "...DONE CHECKING" << endl;
}

void
Simulation::build(){
    cerr << "BUILDING SIMULATION..." << endl;
    Instance *i1 = getModel().getInternalType();
    sys = system_build(&*i1);
    if(!sys){
        throw runtime_error("Unable to build system");
    }
    is_built = true;
    cerr << "...DONE BUILDING" << endl;
}

vector<Variable>
Simulation::getFixableVariables(){
    cerr << "GETTING FIXABLE VARIABLES..." << endl;
    vector<Variable> vars;
    vars.reserve(100);

    if(!sys){
        throw runtime_error("Simulation system not yet built");
    }

    int32 **vip; /** TODO ensure 32 bit integers are used */

    // Get IDs of elegible variables in array at vip...
    if(!slvDOF_eligible(sys,vip)){
        ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"No fixable variables found.");
    }else{
        //cerr << "FIXABLE VARS FOUND" << endl;
        struct var_variable **vp = slv_get_solvers_var_list(sys);

        /*struct var_variable *first_var = vp[0];
        char *first_var_name = var_make_name(sys,first_var);
        cerr << "FIRST SYS VAR IS NAMED " << var_make_name(s,first_var) << endl;
        ascfree(first_var_name);*/

        if(vp==NULL){
            throw runtime_error("Simulation variable list is null");
        }

        // iterate through this list until we find a -1:
        int i=0;
        int var_index = (*vip)[i];
        while(var_index >= 0){
            //cerr << "FOUND VARIABLE var_index = " << var_index << endl;
            struct var_variable *var = vp[var_index];
            //cerr << "VARIABLE " << var_index << " IS ELIGIBLE" << endl;
            char *var_name = var_make_name(sys,var);
            //cerr << "ELIGIBLE VAR: " << var_name << endl;
            ascfree(var_name);
            vars.push_back(Variable( sys,var ));
            ++i;
            var_index = (*vip)[i];
        }
        ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"Found %d fixable variables.",i);
        //cerr << "END ELEGIBLE VARS LIST" << endl;
        ascfree(*vip);
        //cerr << "FREED VIP LIST" << endl;
    }

    //cerr << "FINISHED WITH FINDING ELEGIBLE VARIABLES" << endl;
    return vars;
}

    
void
Simulation::solve(Solver solver){
    if(!is_built){
        throw runtime_error("Simulation::solver: simulation is not yet built, can't start solving.");
    }

    cerr << "SIMULATION::SOLVE STARTING..." << endl;
    enum inst_t k = getModel().getKind();
    if(k!=MODEL_INST)throw runtime_error("Can't solve: not an instance of type MODEL_INST");

    Instance *i1 = getInternalType();
    int npend = NumberPendingInstances(&*i1);
    if(npend)throw runtime_error("Can't solve: There are still %d pending instances");

    if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");

    cerr << "SIMULATION::SOLVE: SET SOLVER..." << endl;
    setSolver(solver);


    cerr << "PRESOLVING SYSTEM...";
    slv_presolve(sys);
    cerr << "DONE" << endl;

    cerr << "SOLVING SYSTEM..." << endl;
    // Add some stuff here for cleverer iteration....
    unsigned niter = 1000;
    double updateinterval = 0.02;

    double starttime = tm_cpu_time();
    double lastupdate = starttime;
    slv_status_t status;
    int solved_vars=0;
    bool stop=false;

    for(int iter = 1; iter <= niter && !stop; ++iter){
        slv_get_status(sys,&status);
        if(status.ready_to_solve){
            slv_iterate(sys);
        }
        if(tm_cpu_time() - lastupdate > updateinterval && iter > 0){
            if(solved_vars < status.block.previous_total_size){
                solved_vars = status.block.previous_total_size;
                CONSOLE_DEBUG("Solved %5d vars. Now block %3d/%-3d...", solved_vars, status.block.current_block, status.block.number_of);
            }else{
                CONSOLE_DEBUG("Total %5d iterations (%5d in current block of %d vars)", iter, status.block.iteration, status.block.current_size);
            }
            lastupdate = tm_cpu_time();
        }
    }
    double elapsed = tm_cpu_time() - starttime;

    if(status.ok){
        cerr << "... DONE SOLVING SYSTEM" << endl;
    }else{
        ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Solver failed");
    }

    cerr << "SOLVER PERFORMED " << status.iteration << " ITERATIONS IN " << elapsed << "s" << endl;

    if(status.iteration_limit_exceeded){
        ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Exceeded interation limit");
    }

    if(status.converged){
        ERROR_REPORTER_NOLINE(ASC_USER_SUCCESS,"Solver converged: %d iterations (%.2f s)"
            ,status.iteration,elapsed);
    }else{
        ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Solver not converged after %d iterations (%.2f s).",status.iteration,elapsed);
    }

}

void
Simulation::write(){
    simroot.write();
}

//------------------------------------------
// ASSIGNING SOLVER TO SIMULATION

void
Simulation::setSolver(Solver &solver){
    cerr << "SETTING SOLVER ON SIMULATION TO " << solver.getName() << endl;

    if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");
    // Update the solver object because sometimes an alternative solver can be returned, apparently.

    int selected = slv_select_solver(sys, solver.getIndex());
    //cerr << "Simulation::setSolver: slv_select_solver returned " << selected << endl;

    if(selected<0){
        ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Failed to select solver");
        throw runtime_error("Failed to select solver");
    }

    if(selected!=solver.getIndex()){
        solver = Solver(slv_solver_name(selected));
        ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Substitute solver '%s' (index %d) selected.\n", solver.getName().c_str(), selected);
    }

    if( slv_eligible_solver(sys) <= 0){
        ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Inelegible solver '%s'", solver.getName().c_str() );
        throw runtime_error("Inelegible solver");
    }
}

const Solver
Simulation::getSolver() const{
    int index = slv_get_selected_solver(sys);
    //cerr << "Simulation::getSolver: index = " << index << endl;
    if(index<0)throw runtime_error("No solver selected");

    return Solver(slv_solver_name(index));
}


/**
    Get solver parameters struct wrapped up as a SolverParameters class.
*/
SolverParameters
Simulation::getSolverParameters() const{
    if(!sys)throw runtime_error("Can't getSolverParameters: Simulation system has not been built yet.");

    slv_parameters_t p;
    slv_get_parameters(sys,&p);
    return SolverParameters(p);
}

/**
    Update the solver parameters by passing a new set back
*/
void
Simulation::setSolverParameters(SolverParameters &P){
    if(!sys)throw runtime_error("Can't set solver parameters: simulation has not been built yet.");
    slv_set_parameters(sys, &(P.getInternalType()));
}

slv_system_structure *
Simulation::getSystem(){
    if(!sys)throw runtime_error("Can't getSystem: simulation not yet built");
    return sys;
}

IncidenceMatrix
Simulation::getIncidenceMatrix(){
    return IncidenceMatrix(*this);
}
1	johnpye	132	#include <iostream>
2			#include <iomanip>
3			#include <stdexcept>
4			#include <sstream>
5			using namespace std;
6
7			extern "C"{
8			#include <utilities/ascConfig.h>
9			#include <utilities/ascSignal.h>
10			#include <utilities/ascMalloc.h>
11			#include <general/dstring.h>
12	johnpye	196	#include <general/tm_time.h>
13	johnpye	132	#include <compiler/instance_enum.h>
14			#include <compiler/fractions.h>
15			#include <compiler/compiler.h>
16			#include <compiler/dimen.h>
17			#include <compiler/symtab.h>
18			#include <compiler/instance_io.h>
19			#include <compiler/instantiate.h>
20			#include <compiler/bintoken.h>
21			#include <compiler/instance_enum.h>
22			#include <compiler/instquery.h>
23			#include <compiler/check.h>
24			#include <compiler/name.h>
25			#include <compiler/pending.h>
26
27			#include <utilities/readln.h>
28			#include <solver/mtx.h>
29			#include <solver/slv_types.h>
30			#include <solver/var.h>
31			#include <solver/rel.h>
32			#include <solver/discrete.h>
33			#include <solver/conditional.h>
34			#include <solver/logrel.h>
35			#include <solver/bnd.h>
36			#include <solver/calc.h>
37			#include <solver/relman.h>
38			#include <solver/slv_common.h>
39			#include <solver/linsol.h>
40			#include <solver/linsolqr.h>
41			#include <solver/slv_client.h>
42			#include <solver/system.h>
43			#include <solver/slv_interface.h>
44			#include <solver/slvDOF.h>
45			#include <solver/slv3.h>
46			#include <solver/slv_stdcalls.h>
47			}
48
49			#include "simulation.h"
50			#include "solver.h"
51	johnpye	208	#include "solverparameters.h"
52	johnpye	132	#include "name.h"
53	johnpye	233	#include "incidencematrix.h"
54	johnpye	132
55			/**
56			Create an instance of a type (call compiler etc)
57
58			@TODO fix mutex on compile command filenames
59			*/
60			Simulation::Simulation(Instance *i, const SymChar &name) : Instanc(i, name), simroot(GetSimulationRoot(i),SymChar("simroot")){
61	johnpye	153	is_built = false;
62	johnpye	132	// Create an Instance object for the 'simulation root' (we'll call
63	johnpye	190	// it the 'simulation model') and it can be fetched using 'getModel()'
64	johnpye	132	// any time later.
65			//simroot = Instanc(GetSimulationRoot(i),name);
66			}
67
68	johnpye	164	Simulation::Simulation(const Simulation &old) : Instanc(old), simroot(old.simroot){
69			is_built = old.is_built;
70			sys = old.sys;
71			bin_srcname = old.bin_srcname;
72			bin_objname = old.bin_objname;
73			bin_libname = old.bin_libname;
74			bin_cmd = old.bin_cmd;
75	johnpye	190	bin_rm = old.bin_rm;
76	johnpye	164	}
77
78			Simulation::~Simulation(){
79			//CONSOLE_DEBUG("Deleting simulation %s", getName().toString());
80			}
81
82	johnpye	132	Instanc &
83			Simulation::getModel(){
84			if(!simroot.getInternalType()){
85			throw runtime_error("Simulation::getModel: simroot.getInternalType()is NULL");
86			}
87			return simroot;
88			}
89
90			void
91			Simulation::checkDoF() const{
92			cerr << "CHECKING DOF..." << endl;
93			int dof, status;
94	johnpye	164	if(!sys){
95	johnpye	132	throw runtime_error("System not yet built");
96			}
97			slvDOF_status(sys, &status, &dof);
98			switch(status){
99	johnpye	190	case 1: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Underspecified; %d degrees of freedom",dof); break;
100			case 2: ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"Square"); break;
101			case 3: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Structurally singular"); break;
102			case 4: ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Overspecified"); break;
103	johnpye	132	case 5:
104			throw runtime_error("Unable to resolve degrees of freedom"); break;
105			default:
106			throw runtime_error("Invalid return status from slvDOF_status");
107			}
108			}
109
110			void
111			Simulation::run(const Method &method){
112			cerr << "RUNNING PROCEDURE " << method.getName() << endl;
113			Nam name = Nam(method.getSym());
114	johnpye	154	//cerr << "CREATED NAME '" << name.getName() << "'" << endl;
115	johnpye	132	Proc_enum pe;
116			pe = Initialize(
117	johnpye	164	&*(getModel().getInternalType()) ,name.getInternalType(), "__not_named__"
118	johnpye	132	,ASCERR
119			,0, NULL, NULL
120			);
121
122			if(pe == Proc_all_ok){
123	johnpye	190	ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Method '%s' was run (check above for errors)\n",method.getName());
124	johnpye	132	//cerr << "METHOD " << method.getName() << " COMPLETED OK" << endl;
125			}else{
126			stringstream ss;
127			ss << "Simulation::run: Method '" << method.getName() << "' returned error: ";
128			switch(pe){
129			case Proc_CallOK: ss << "Call OK"; break;
130			case Proc_CallError: ss << "Error occurred in call"; break;
131			case Proc_CallReturn: ss << "Request that caller return (OK)"; break;
132			case Proc_CallBreak: ss << "Break out of enclosing loop"; break;
133			case Proc_CallContinue: ss << "Skip to next iteration"; break;
134
135			case Proc_break: ss << "Break"; break;
136			case Proc_continue: ss << "Continue"; break;
137			case Proc_fallthru: ss << "Fall-through"; break;
138			case Proc_return: ss << "Return"; break;
139			case Proc_stop: ss << "Stop"; break;
140			case Proc_stack_exceeded: ss << "Stack exceeded"; break;
141			case Proc_stack_exceeded_this_frame: ss << "Stack exceeded this frame"; break;
142			case Proc_case_matched: ss << "Case matched"; break;
143			case Proc_case_unmatched: ss << "Case unmatched"; break;
144
145			case Proc_case_undefined_value: ss << "Undefined value in case"; break;
146			case Proc_case_boolean_mismatch: ss << "Boolean mismatch in case"; break;
147			case Proc_case_integer_mismatch: ss << "Integer mismatch in case"; break;
148			case Proc_case_symbol_mismatch: ss << "Symbol mismatch in case"; break;
149			case Proc_case_wrong_index: ss << "Wrong index in case"; break;
150			case Proc_case_wrong_value: ss << "Wrong value in case"; break;
151			case Proc_case_extra_values: ss << "Extra values in case"; break;
152			case Proc_bad_statement: ss << "Bad statement"; break;
153			case Proc_bad_name: ss << "Bad name"; break;
154			case Proc_for_duplicate_index: ss << "Duplicate index"; break;
155			case Proc_for_set_err: ss << "For set error"; break;
156			case Proc_for_not_set: ss << "For not set"; break;
157			case Proc_illegal_name_use: ss << "Illegal name use"; break;
158			case Proc_name_not_found: ss << "Name not found"; break;
159			case Proc_instance_not_found: ss << "Instance not found"; break;
160			case Proc_type_not_found: ss << "Type not found"; break;
161			case Proc_illegal_type_use: ss << "Illegal use"; break;
162			case Proc_proc_not_found: ss << "Method not found"; break;
163			case Proc_if_expr_error_typeconflict: ss << "Type conflict in 'if' expression"; break;
164			case Proc_if_expr_error_nameunfound: ss << "Name not found in 'if' expression"; break;
165			case Proc_if_expr_error_incorrectname: ss << "Incorrect name in 'if' expression"; break;
166			case Proc_if_expr_error_undefinedvalue: ss << "Undefined value in 'if' expression"; break;
167			case Proc_if_expr_error_dimensionconflict: ss << "Dimension conflict in 'if' expression"; break;
168			case Proc_if_expr_error_emptychoice: ss << "Empty choice in 'if' expression"; break;
169			case Proc_if_expr_error_emptyintersection: ss << "Empty intersection in 'if' expression"; break;
170			case Proc_if_expr_error_confused: ss << "Confused in 'if' expression"; break;
171			case Proc_if_real_expr: ss << "Real-valued result in 'if' expression"; break;
172			case Proc_if_integer_expr: ss << "Integeter-valued result in 'if' expression"; break;
173			case Proc_if_symbol_expr: ss << "Symbol-valued result in 'if' expression"; break;
174			case Proc_if_set_expr: ss << "Set-valued result in 'if' expression"; break;
175			case Proc_if_not_logical: ss << "If expression is not logical"; break;
176			case Proc_user_interrupt: ss << "User interrupt"; break;
177			case Proc_infinite_loop: ss << "Infinite loop"; break;
178			case Proc_declarative_constant_assignment: ss << "Declarative constant assignment"; break;
179			case Proc_nonsense_assignment: ss << "Nonsense assginment (bogus)"; break;
180			case Proc_nonconsistent_assignment: ss << "Inconsistent assignment"; break;
181			case Proc_nonatom_assignment: ss << "Non-atom assignment"; break;
182			case Proc_nonboolean_assignment: ss << "Non-boolean assignment"; break;
183			case Proc_noninteger_assignment: ss << "Non-integer assignment"; break;
184			case Proc_nonreal_assignment: ss << "Non-real assignment"; break;
185			case Proc_nonsymbol_assignment: ss << "Non-symbol assignment"; break;
186			case Proc_lhs_error: ss << "Left-hand-side error"; break;
187			case Proc_rhs_error: ss << "Right-hand-side error"; break;
188			case Proc_unknown_error: ss << "Unknown error"; break;
189			default:
190			ss << "Invalid error code";
191			}
192
193
194			ss << " (" << int(pe) << ")";
195			throw runtime_error(ss.str());
196			}
197			}
198
199			const bool
200			Simulation::check(){
201			cerr << "CHECKING SIMULATION" << endl;
202			Instance *i1 = getModel().getInternalType();
203	johnpye	164	CheckInstance(stderr, &*i1);
204	johnpye	132	cerr << "...DONE CHECKING" << endl;
205	johnpye	190	}
206	johnpye	132
207			void
208			Simulation::build(){
209			cerr << "BUILDING SIMULATION..." << endl;
210	johnpye	164	Instance *i1 = getModel().getInternalType();
211			sys = system_build(&*i1);
212			if(!sys){
213	johnpye	132	throw runtime_error("Unable to build system");
214			}
215	johnpye	153	is_built = true;
216	johnpye	132	cerr << "...DONE BUILDING" << endl;
217			}
218
219			vector<Variable>
220			Simulation::getFixableVariables(){
221			cerr << "GETTING FIXABLE VARIABLES..." << endl;
222			vector<Variable> vars;
223			vars.reserve(100);
224
225	johnpye	164	if(!sys){
226	johnpye	132	throw runtime_error("Simulation system not yet built");
227			}
228
229			int32 vip; / TODO ensure 32 bit integers are used */
230
231			// Get IDs of elegible variables in array at vip...
232			if(!slvDOF_eligible(sys,vip)){
233	johnpye	190	ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"No fixable variables found.");
234	johnpye	132	}else{
235			//cerr << "FIXABLE VARS FOUND" << endl;
236			struct var_variable **vp = slv_get_solvers_var_list(sys);
237
238			/struct var_variable first_var = vp[0];
239			char *first_var_name = var_make_name(sys,first_var);
240			cerr << "FIRST SYS VAR IS NAMED " << var_make_name(s,first_var) << endl;
241			ascfree(first_var_name);*/
242
243			if(vp==NULL){
244			throw runtime_error("Simulation variable list is null");
245			}
246	johnpye	190
247	johnpye	132	// iterate through this list until we find a -1:
248			int i=0;
249			int var_index = (*vip)[i];
250			while(var_index >= 0){
251			//cerr << "FOUND VARIABLE var_index = " << var_index << endl;
252			struct var_variable *var = vp[var_index];
253			//cerr << "VARIABLE " << var_index << " IS ELIGIBLE" << endl;
254			char *var_name = var_make_name(sys,var);
255			//cerr << "ELIGIBLE VAR: " << var_name << endl;
256			ascfree(var_name);
257			vars.push_back(Variable( sys,var ));
258			++i;
259			var_index = (*vip)[i];
260			}
261	johnpye	190	ERROR_REPORTER_NOLINE(ASC_USER_NOTE,"Found %d fixable variables.",i);
262	johnpye	132	//cerr << "END ELEGIBLE VARS LIST" << endl;
263			ascfree(*vip);
264			//cerr << "FREED VIP LIST" << endl;
265			}
266
267			//cerr << "FINISHED WITH FINDING ELEGIBLE VARIABLES" << endl;
268			return vars;
269			}
270
271	johnpye	196
272	johnpye	132	void
273			Simulation::solve(Solver solver){
274	johnpye	153	if(!is_built){
275			throw runtime_error("Simulation::solver: simulation is not yet built, can't start solving.");
276			}
277
278	johnpye	132	cerr << "SIMULATION::SOLVE STARTING..." << endl;
279			enum inst_t k = getModel().getKind();
280			if(k!=MODEL_INST)throw runtime_error("Can't solve: not an instance of type MODEL_INST");
281	johnpye	190
282	johnpye	164	Instance *i1 = getInternalType();
283			int npend = NumberPendingInstances(&*i1);
284	johnpye	190	if(npend)throw runtime_error("Can't solve: There are still %d pending instances");
285	johnpye	132
286			if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");
287	johnpye	190
288	johnpye	132	cerr << "SIMULATION::SOLVE: SET SOLVER..." << endl;
289			setSolver(solver);
290	johnpye	190
291
292	johnpye	207	cerr << "PRESOLVING SYSTEM...";
293	johnpye	132	slv_presolve(sys);
294	johnpye	207	cerr << "DONE" << endl;
295	johnpye	190
296	johnpye	132	cerr << "SOLVING SYSTEM..." << endl;
297	johnpye	196	// Add some stuff here for cleverer iteration....
298	johnpye	207	unsigned niter = 1000;
299	johnpye	196	double updateinterval = 0.02;
300	johnpye	132
301	johnpye	196	double starttime = tm_cpu_time();
302			double lastupdate = starttime;
303			slv_status_t status;
304			int solved_vars=0;
305			bool stop=false;
306
307	johnpye	207	for(int iter = 1; iter <= niter && !stop; ++iter){
308	johnpye	196	slv_get_status(sys,&status);
309			if(status.ready_to_solve){
310			slv_iterate(sys);
311			}
312			if(tm_cpu_time() - lastupdate > updateinterval && iter > 0){
313			if(solved_vars < status.block.previous_total_size){
314			solved_vars = status.block.previous_total_size;
315			CONSOLE_DEBUG("Solved %5d vars. Now block %3d/%-3d...", solved_vars, status.block.current_block, status.block.number_of);
316			}else{
317			CONSOLE_DEBUG("Total %5d iterations (%5d in current block of %d vars)", iter, status.block.iteration, status.block.current_size);
318			}
319			lastupdate = tm_cpu_time();
320			}
321			}
322	johnpye	230	double elapsed = tm_cpu_time() - starttime;
323	johnpye	196
324			if(status.ok){
325	johnpye	132	cerr << "... DONE SOLVING SYSTEM" << endl;
326			}else{
327	johnpye	190	ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Solver failed");
328	johnpye	132	}
329
330	johnpye	230	cerr << "SOLVER PERFORMED " << status.iteration << " ITERATIONS IN " << elapsed << "s" << endl;
331	johnpye	132
332	johnpye	196	if(status.iteration_limit_exceeded){
333	johnpye	190	ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Exceeded interation limit");
334	johnpye	132	}
335
336	johnpye	196	if(status.converged){
337	johnpye	230	ERROR_REPORTER_NOLINE(ASC_USER_SUCCESS,"Solver converged: %d iterations (%.2f s)"
338			,status.iteration,elapsed);
339	johnpye	132	}else{
340	johnpye	230	ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Solver not converged after %d iterations (%.2f s).",status.iteration,elapsed);
341	johnpye	132	}
342
343			}
344
345			void
346			Simulation::write(){
347			simroot.write();
348			}
349
350			//------------------------------------------
351			// ASSIGNING SOLVER TO SIMULATION
352	johnpye	190
353	johnpye	132	void
354			Simulation::setSolver(Solver &solver){
355	johnpye	151	cerr << "SETTING SOLVER ON SIMULATION TO " << solver.getName() << endl;
356
357	johnpye	132	if(!sys)throw runtime_error("Can't solve: Simulation system has not been built yet.");
358			// Update the solver object because sometimes an alternative solver can be returned, apparently.
359
360			int selected = slv_select_solver(sys, solver.getIndex());
361	johnpye	154	//cerr << "Simulation::setSolver: slv_select_solver returned " << selected << endl;
362	johnpye	132
363			if(selected<0){
364	johnpye	190	ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Failed to select solver");
365	johnpye	132	throw runtime_error("Failed to select solver");
366			}
367
368			if(selected!=solver.getIndex()){
369	johnpye	150	solver = Solver(slv_solver_name(selected));
370	johnpye	190	ERROR_REPORTER_NOLINE(ASC_PROG_NOTE,"Substitute solver '%s' (index %d) selected.\n", solver.getName().c_str(), selected);
371	johnpye	132	}
372	johnpye	190
373	johnpye	132	if( slv_eligible_solver(sys) <= 0){
374	johnpye	190	ERROR_REPORTER_NOLINE(ASC_PROG_ERROR,"Inelegible solver '%s'", solver.getName().c_str() );
375	johnpye	132	throw runtime_error("Inelegible solver");
376			}
377			}
378
379			const Solver
380			Simulation::getSolver() const{
381			int index = slv_get_selected_solver(sys);
382	johnpye	154	//cerr << "Simulation::getSolver: index = " << index << endl;
383	johnpye	132	if(index<0)throw runtime_error("No solver selected");
384
385	johnpye	150	return Solver(slv_solver_name(index));
386	johnpye	132	}
387
388
389	johnpye	208	/**
390			Get solver parameters struct wrapped up as a SolverParameters class.
391			*/
392			SolverParameters
393			Simulation::getSolverParameters() const{
394	johnpye	225	if(!sys)throw runtime_error("Can't getSolverParameters: Simulation system has not been built yet.");
395	johnpye	132
396	johnpye	208	slv_parameters_t p;
397			slv_get_parameters(sys,&p);
398			return SolverParameters(p);
399			}
400
401	johnpye	225	/**
402			Update the solver parameters by passing a new set back
403			*/
404			void
405			Simulation::setSolverParameters(SolverParameters &P){
406			if(!sys)throw runtime_error("Can't set solver parameters: simulation has not been built yet.");
407			slv_set_parameters(sys, &(P.getInternalType()));
408			}
409
410	johnpye	233	slv_system_structure *
411			Simulation::getSystem(){
412			if(!sys)throw runtime_error("Can't getSystem: simulation not yet built");
413			return sys;
414			}
415
416			IncidenceMatrix
417			Simulation::getIncidenceMatrix(){
418			return IncidenceMatrix(*this);
419			}