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

Diff of /trunk/pygtk/integrator.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 941 by johnpye, Fri Nov 24 10:46:32 2006 UTC revision 942 by johnpye, Sat Nov 25 05:26:47 2006 UTC
# Line 1  Line 1 
1  #include "integrator.h"  #include "integrator.h"
2  #include "integratorreporter.h"  #include "integratorreporter.h"
3    #include "solverparameters.h"
4  #include <stdexcept>  #include <stdexcept>
5  #include <sstream>  #include <sstream>
6    #include <cmath>
7  using namespace std;  using namespace std;
8    
9  /**  /**
# Line 32  Integrator::~Integrator(){ Line 34  Integrator::~Integrator(){
34      samplelist_free(samplelist);      samplelist_free(samplelist);
35  }  }
36    
37    SolverParameters
38    Integrator::getParameters() const{
39        SolverParameters params;
40        int res = integrator_params_get(blsys,&(params.getInternalType() ) );
41        if(res)throw runtime_error("Failed to get integrator parameters");
42        return params;
43    }
44    
45    void
46    Integrator::setParameters(const SolverParameters &params){
47        int res = integrator_params_set(blsys,&(params.getInternalTypeConst() ) );
48        if(res)throw runtime_error("Failed to set integrator parameters");
49    }
50    
51  void  void
52  Integrator::setReporter(IntegratorReporterCxx *reporter){  Integrator::setReporter(IntegratorReporterCxx *reporter){
# Line 78  Integrator::analyse(){ Line 93  Integrator::analyse(){
93    
94  /**  /**
95      @TODO what about root detection?      @TODO what about root detection?
96    
97        Integrate the function for the timesteps specified.
98    
99        This method will throw a runtime_error if integration fails.
100  */  */
101  int  void
102  Integrator::solve(){  Integrator::solve(){
103    
104      // check the integration limits      // check the integration limits
# Line 98  Integrator::solve(){ Line 117  Integrator::solve(){
117      simulation.processVarStatus();      simulation.processVarStatus();
118    
119      if(!res){      if(!res){
120          ERROR_REPORTER_NOLINE(ASC_USER_ERROR,"Failed integration");          throw runtime_error("Failed integration");
         return 0;  
121      }      }
   
     return 1;  
122  }  }
123    
124  void  void
# Line 161  Integrator::getEngineName() const{ Line 177  Integrator::getEngineName() const{
177      return f->second;      return f->second;
178  }        }      
179    
180    /**
181        @TODO what about conversion factors? Is an allowance being made?
182    */
183  void  void
184  Integrator::setLinearTimesteps(UnitsM units, double start, double end, unsigned long num){  Integrator::setLinearTimesteps(UnitsM units, double start, double end, unsigned long num){
185      if(samplelist!=NULL){      if(samplelist!=NULL){
# Line 176  Integrator::setLinearTimesteps(UnitsM un Line 195  Integrator::setLinearTimesteps(UnitsM un
195      }      }
196      integrator_set_samples(blsys,samplelist);      integrator_set_samples(blsys,samplelist);
197  }  }
198    
199    /**
200        @TODO what about conversion factors? Is an allowance being made?
201    */
202    void
203    Integrator::setLogTimesteps(UnitsM units, double start, double end, unsigned long num){
204        if(samplelist!=NULL){
205            ASC_FREE(samplelist);
206        }
207        const dim_type *d = units.getDimensions().getInternalType();
208    
209        if(start<=0)throw runtime_error("starting timestep needs to be > 0");
210        if(end<=0)throw runtime_error("end timestep needs to be > 0");
211        if(end <= start)throw runtime_error("end timestep needs to be > starting timestep");
212    
213        samplelist = samplelist_new(num+1, d);
214        double val = start;
215        double inc = exp((log(end)-log(start))/num);
216        for(unsigned long i=0;i<=num;++i){
217            samplelist_set(samplelist,i,val);
218            CONSOLE_DEBUG("samplelist[%lu] = %f",i,val);
219            val *= inc;
220        }
221        integrator_set_samples(blsys,samplelist);
222    }
223    
224  vector<double>  vector<double>
225  Integrator::getCurrentObservations(){  Integrator::getCurrentObservations(){

Legend:
Removed from v.941  
changed lines
  Added in v.942

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