minsky::RungeKutta Class Reference

#include <rungeKutta.h>

Inheritance diagram for minsky::RungeKutta:

Collaboration diagram for minsky::RungeKutta:

Public Member Functions
virtual	~RungeKutta ()=default
virtual bool	resetIfFlagged ()
	checks whether a reset is required, and resets the simulation if so More...
void	rkreset ()
	reset the simulation More...
void	rkstep ()
	step the equations (by n steps, default 1) evaluate the flow equations without stepping. More...
void	evalEquations ()

Public Attributes
double	t {0}
	time More...
bool	running =false
	controls whether simulation is running More...
bool	reverse =false
	reverse direction of simulation More...
EvalGodley	evalGodley
Public Attributes inherited from minsky::Simulation
double	stepMin {0}
double	stepMax {0.01}
int	nSteps {1}
double	epsRel {1e-2}
double	epsAbs {1e-3}
int	order {4}
bool	implicit {false}
int	simulationDelay {0}
std::string	timeUnit
double	tmax {INFINITY}
double	t0 {0}

Protected Member Functions
void	evalEquations (double result[], double, const double vars[])
void	evalJacobian (Matrix &, double, const double vars[])

Static Protected Member Functions
static int	RKfunction (double, const double y[], double f[], void *)
	function to be integrated (internal use) More...
static int	jacobian (double, const double y[], double *, double dfdt[], void *)
	compute jacobian (internal use) More...

Private Member Functions
	CLASSDESC_ACCESS (RungeKutta)

Friends
struct	RKdata

Additional Inherited Members
Static Public Attributes inherited from minsky::ValueVector
static std::vector< double, CIVITA_ALLOCATOR< double > >	stockVars
	vector of variables that are integrated via Runge-Kutta. These variables label the columns of the Godley table More...
static std::vector< double, CIVITA_ALLOCATOR< double > >	flowVars
	variables defined as a simple function of the stock variables, also known as lhs variables. These variables appear in the body of the Godley table More...

Detailed Description

Definition at line 45 of file rungeKutta.h.

Constructor & Destructor Documentation

~RungeKutta()

virtual minsky::RungeKutta::~RungeKutta ( )

virtualdefault

Member Function Documentation

CLASSDESC_ACCESS()

minsky::RungeKutta::CLASSDESC_ACCESS ( RungeKutta  )

private

evalEquations() [1/2]

void minsky::RungeKutta::evalEquations ( double  result[], double  t, const double  vars[]  )

protected

Definition at line 248 of file rungeKutta.cc.

References minsky::Minsky::displayErrorItem(), and pyminsky::minsky.

Referenced by minsky::VariableBase::onMouseMotion().

  {
    EvalOpBase::t=reverse? -t: t;
    const double reverseFactor=reverse? -1: 1;
    // firstly evaluate the flow variables. Initialise to flowVars so
    // that no input vars are correctly initialised
    auto flow(flowVars);
    for (size_t i=0; i<equations.size(); ++i)
      equations[i]->eval(flow.data(), flow.size(), vars);

    // then create the result using the Godley table
    for (size_t i=0; i<stockVars.size(); ++i) result[i]=0;
    evalGodley.eval(result, flow.data());

    // integrations are kind of a copy
    for (vector<Integral>::iterator i=integrals.begin(); i<integrals.end(); ++i)
      {
        if (i->input().idx()<0)
          {
            if (i->operation)
              minsky().displayErrorItem(*i->operation);
            throw error("integral not wired");
          }
        // enable element-wise integration of tensor variables. for feature 147
        assert(i->input().size()==i->stock->size());
        for (size_t j=0; j<i->input().size(); ++j)
          result[i->stock->idx()+j] = reverseFactor *
            (i->input().isFlowVar()? flow[i->input().idx()+j] : vars[i->input().idx()+j]);
      } 
  }

Here is the call graph for this function:

Here is the caller graph for this function:

evalEquations() [2/2]

void minsky::RungeKutta::evalEquations ( )

inline

Definition at line 70 of file rungeKutta.h.

                         {
      for (auto& eq: equations)
        eq->eval(ValueVector::flowVars.data(), ValueVector::flowVars.size(), ValueVector::stockVars.data());
    }

evalJacobian()

void minsky::RungeKutta::evalJacobian ( Matrix &  jac, double  t, const double  vars[]  )

protected

Definition at line 216 of file rungeKutta.cc.

  {
    EvalOpBase::t=reverse? -t: t;
    const double reverseFactor=reverse? -1: 1;
    // firstly evaluate the flow variables. Initialise to flowVars so
    // that no input vars are correctly initialised
    auto flow=flowVars;
    for (size_t i=0; i<equations.size(); ++i)
      equations[i]->eval(flow.data(), flow.size(), sv);

    // then determine the derivatives with respect to variable j
    for (size_t j=0; j<stockVars.size(); ++j)
      {
        vector<double> ds(stockVars.size()), df(flowVars.size());
        ds[j]=1;
        for (size_t i=0; i<equations.size(); ++i)
          equations[i]->deriv(df.data(), df.size(), ds.data(), sv, flow.data());
        vector<double> d(stockVars.size());
        evalGodley.eval(d.data(), df.data());
        for (vector<Integral>::iterator i=integrals.begin(); 
             i!=integrals.end(); ++i)
          {
            assert(i->stock->idx()>=0 && i->input().idx()>=0);
            d[i->stock->idx()] = 
              i->input().isFlowVar()? df[i->input().idx()]: ds[i->input().idx()];
          }
        for (size_t i=0; i<stockVars.size(); i++)
          jac(i,j)=reverseFactor*d[i];
      }
  
  }

jacobian()

int minsky::RungeKutta::jacobian ( double  t, const double  y[], double *  dfdy, double  dfdt[], void *  params  )

staticprotected

compute jacobian (internal use)

Definition at line 65 of file rungeKutta.cc.

  {
    if (params==NULL) return GSL_EBADFUNC;
    Matrix jac(stockVars.size(), dfdy);
    try
      {
        ((RungeKutta*)params)->evalJacobian(jac,t,y);
      }
    catch (std::exception& e)
      {
        ((RungeKutta*)params)->threadErrMsg=e.what();
        return GSL_EBADFUNC;
      }   
    return GSL_SUCCESS;
  }

resetIfFlagged()

virtual bool minsky::RungeKutta::resetIfFlagged ( )

inlinevirtual

checks whether a reset is required, and resets the simulation if so

Returns

whether simulation was reset

Reimplemented in minsky::Minsky.

Definition at line 65 of file rungeKutta.h.

{return false;}

RKfunction()

int minsky::RungeKutta::RKfunction ( double  t, const double  y[], double  f[], void *  params  )

staticprotected

function to be integrated (internal use)

checks if any GUI events are waiting, and proces an event if so TODO - make this a virtual method of RungeKutta?

Definition at line 50 of file rungeKutta.cc.

References f.

  {
    if (params==NULL) return GSL_EBADFUNC;
    try
      {
        ((RungeKutta*)params)->evalEquations(f,t,y);
      }
    catch (std::exception& e)
      {
        ((RungeKutta*)params)->threadErrMsg=e.what();
        return GSL_EBADFUNC;
      }
    return GSL_SUCCESS;
  }

rkreset()

void minsky::RungeKutta::rkreset

(

)

reset the simulation

Definition at line 122 of file rungeKutta.cc.

  {
    if (order==1 && !implicit)
      ode.reset(); // do explicit Euler
    else
      ode.reset(new RKdata(this)); // set up GSL ODE routines
  }

rkstep()

void minsky::RungeKutta::rkstep

(

)

step the equations (by n steps, default 1) evaluate the flow equations without stepping.

Exceptions

ecolab::error

if equations are illdefined

Definition at line 130 of file rungeKutta.cc.

References minsky::doOneEvent().

  {
    if (nSteps<1) return;
    running=true; // this ensure reset flag is cleared
    resetIfFlagged();
    running=true;
    
    // create a private copy for worker thread use
    auto stockVarsCopy(stockVars);
    RKThreadRunning=true;
    int err=GSL_SUCCESS;
    // run RK algorithm on a separate worker thread so as to not block UI. See ticket #6
    boost::thread rkThread([&]() {
      try
        { 
          double tp=reverse? -t: t;
          if (ode)
            {
              gsl_odeiv2_driver_set_nmax(ode->driver, nSteps);
              // we need to update Minsky's t synchronously to support the t operator
              // potentially means t and stockVars out of sync on GUI, but should still be thread safe
              err=gsl_odeiv2_driver_apply(ode->driver, &tp, numeric_limits<double>::max(), 
                                          stockVarsCopy.data());
            }
          else // do explicit Euler method
            {
              vector<double> d(stockVarsCopy.size());
              for (int i=0; i<nSteps; ++i, tp+=stepMax)
                {
                  evalEquations(d.data(), tp, stockVarsCopy.data());
                  for (size_t j=0; j<d.size(); ++j)
                    stockVarsCopy[j]+=d[j];
                }
            }
          t=reverse? -tp:tp;
        }
      catch (const std::exception& ex)
        {
          // catch any thrown exception, and report back to GUI thread
          threadErrMsg=ex.what();
        }
      catch (...)
        {
          threadErrMsg="Unknown exception thrown on ODE solver thread";
        }
      RKThreadRunning=false;
    });

    while (RKThreadRunning)
      {
        // while waiting for thread to finish, check and process any UI events
        usleep(1000);
        doOneEvent(false);
      }
    rkThread.join();

    if (!threadErrMsg.empty())
      {
        auto msg=std::move(threadErrMsg);
        threadErrMsg.clear(); // to be sure, to be sure
        // rethrow exception so message gets displayed to user
        throw runtime_error(msg);
      }

    if (resetIfFlagged())
      return; // in case reset() was called during the step evaluation

    switch (err)
      {
      case GSL_SUCCESS: case GSL_EMAXITER: break;
      case GSL_FAILURE:
        throw error("unspecified error GSL_FAILURE returned");
      case GSL_EBADFUNC: 
        gsl_odeiv2_driver_reset(ode->driver);
        throw error("Invalid arithmetic operation detected");
      default:
        throw error("gsl error: %s",gsl_strerror(err));
      }

    stockVars.swap(stockVarsCopy);

    // update flow variables
    evalEquations();

  }

Here is the call graph for this function:

Friends And Related Function Documentation

RKdata

friend struct RKdata

friend

Definition at line 55 of file rungeKutta.h.

Member Data Documentation

evalGodley

EvalGodley minsky::RungeKutta::evalGodley

Definition at line 60 of file rungeKutta.h.

reverse

bool minsky::RungeKutta::reverse =false

reverse direction of simulation

Definition at line 59 of file rungeKutta.h.

running

bool minsky::RungeKutta::running =false

controls whether simulation is running

Definition at line 58 of file rungeKutta.h.

t

double minsky::RungeKutta::t {0}

time

Definition at line 57 of file rungeKutta.h.

Referenced by minsky::VariableBase::draw().

---

The documentation for this class was generated from the following files:

• engine/rungeKutta.h
• engine/rungeKutta.cc