diff --git a/doc/dynare.texi b/doc/dynare.texi
index cb1629f213cdac0cc9ba3ae872c145094d8b728f..351ed5acd460649908733f5c33d88796c0d6887c 100644
--- a/doc/dynare.texi
+++ b/doc/dynare.texi
@@ -183,6 +183,7 @@ Steady state
* Finding the steady state with Dynare nonlinear solver::
* Using a steady state file::
+* Replace some equations during steady state computations::
Stochastic solution and simulation
@@ -2329,6 +2330,7 @@ providing it with a ``steady state file''.
@menu
* Finding the steady state with Dynare nonlinear solver::
* Using a steady state file::
+* Replace some equations during steady state computations::
@end menu
@node Finding the steady state with Dynare nonlinear solver
@@ -2677,6 +2679,50 @@ steady;
@end deffn
+@node Replace some equations during steady state computations
+@subsection Replace some equations during steady state computations
+
+When there is no steady state file, Dynare computes the steady state
+by solving the static model, @i{i.e.} the model from the @file{.mod}
+file from which leads and lags have been removed.
+
+In some specific cases, one may want to have more control over the way
+this static model is created. Dynare therefore offers the possibility
+to explicitly give the form of equations that should be in the static
+model.
+
+More precisely, if an equation is prepended by a @code{[static]} tag,
+then it will appear in the static model used for steady state
+computation, but that equation will not be used for other
+computations. For every equation tagged in this way, you must tag
+another equation with @code{[dynamic]}: that equation will not be used
+for steady state computation, but will be used for other computations.
+
+This functionality can be useful on models with a unit root, where
+there is an infinity of steady states. An equation (tagged
+@code{[dynamic]}) would give the law of motion of the nonstationary
+variable (like a random walk). To pin down one specific steady state,
+an equation tagged @code{[static]} would affect a constant value to
+the nonstationary variable.
+
+@examplehead
+
+This is a trivial example with two endogenous variables. The second equation
+takes a different form in the static model.
+
+@example
+var c k;
+varexo x;
+
+@dots{}
+
+model;
+c + k - aa*x*k(-1)^alph - (1-delt)*k(-1);
+[dynamic] c^(-gam) - (1+bet)^(-1)*(aa*alph*x(+1)*k^(alph-1) + 1 - delt)*c(+1)^(-gam);
+[static] k = ((delt+bet)/(x*aa*alph))^(1/(alph-1));
+end;
+@end example
+
@node Getting information about the model
@section Getting information about the model
diff --git a/matlab/evaluate_steady_state.m b/matlab/evaluate_steady_state.m
index 23936aaca5058f5b0942b8a01b3ea95c5ac9ba62..77d0c5897f1a8eeaa6d0fbfafb0e9b6cfb8e5286 100644
--- a/matlab/evaluate_steady_state.m
+++ b/matlab/evaluate_steady_state.m
@@ -22,7 +22,7 @@ function [ys,params,info] = evaluate_steady_state(ys_init,M,options,oo,steadysta
% SPECIAL REQUIREMENTS
% none
-% Copyright (C) 2001-2012 Dynare Team
+% Copyright (C) 2001-2013 Dynare Team
%
% This file is part of Dynare.
%
@@ -97,6 +97,31 @@ function [ys,params,info] = evaluate_steady_state(ys_init,M,options,oo,steadysta
return
end
+ % If some equations are tagged [static] or [dynamic], verify consistency
+ if M.static_and_dynamic_models_differ
+ % Evaluate residual of *dynamic* model using the steady state
+ % computed on the *static* one
+ z = repmat(ys,1,M.maximum_lead + M.maximum_lag + 1);
+ zx = repmat([oo.exo_simul oo.exo_det_simul],M.maximum_lead + M.maximum_lag + 1, 1);
+ if options.bytecode
+ [chck, r, junk]= bytecode('dynamic','evaluate', z, zx, M.params, ys, 1);
+ mexErrCheck('bytecode', chck);
+ elseif options.block
+ [r, data] = feval([M.fname '_dynamic'], z', zx, M.params, ys, M.maximum_lag+1, data);
+ else
+ iyv = M.lead_lag_incidence';
+ iyr0 = find(iyv(:));
+ xys = z(iyr0);
+ r = feval([M.fname '_dynamic'], z(iyr0), zx, M.params, ys, M.maximum_lag + 1);
+ end
+
+ % Fail if residual greater than tolerance
+ if max(abs(r)) > options.solve_tolf
+ info(1) = 25;
+ return
+ end
+ end
+
if ~isreal(ys)
info(1) = 21;
info(2) = sum(imag(ys).^2);
diff --git a/matlab/print_info.m b/matlab/print_info.m
index abf4f8faece7c54b2ff0b7c697b258e209391e16..83ac77c5cca35164142693393a5fc75e7d136079 100644
--- a/matlab/print_info.m
+++ b/matlab/print_info.m
@@ -10,7 +10,7 @@ function print_info(info,noprint)
% SPECIAL REQUIREMENTS
% none
-% Copyright (C) 2005-2012 Dynare Team
+% Copyright (C) 2005-2013 Dynare Team
%
% This file is part of Dynare.
%
@@ -78,6 +78,8 @@ if ~noprint
error('Some updated params are complex')
case 24
error('Some updated params contain NaN or Inf')
+ case 25
+ error('The solution to the static equations is not a steady state of the dynamic model: verify that the equations tagged by [static] and [dynamic] are consistent')
case 30
error('Variance can''t be computed')
case 41
diff --git a/preprocessor/DynamicModel.cc b/preprocessor/DynamicModel.cc
index b2fdc770569f86cba0ba6deb767ecdabf9d32217..142a2afd2d3579b177510ed69115ff1323d03bee 100644
--- a/preprocessor/DynamicModel.cc
+++ b/preprocessor/DynamicModel.cc
@@ -2354,6 +2354,12 @@ DynamicModel::writeOutput(ostream &output, const string &basename, bool block_de
<< equation_tags[i].second.second << "' ;" << endl;
output << "};" << endl;
+ /* Say if static and dynamic models differ (because of [static] and [dynamic]
+ equation tags) */
+ output << "M_.static_and_dynamic_models_differ = "
+ << (static_only_equations.size() > 0 ? "1" : "0")
+ << ";" << endl;
+
//In case of sparse model, writes the block_decomposition structure of the model
if (block_decomposition)
{
@@ -3403,6 +3409,11 @@ DynamicModel::cloneDynamic(DynamicModel &dynamic_model) const
for (deque<BinaryOpNode *>::const_iterator it = aux_equations.begin();
it != aux_equations.end(); it++)
dynamic_model.addAuxEquation((*it)->cloneDynamic(dynamic_model));
+
+ // Convert static_only equations
+ for (vector<BinaryOpNode *>::const_iterator it = static_only_equations.begin();
+ it != static_only_equations.end(); it++)
+ dynamic_model.addStaticOnlyEquation((*it)->cloneDynamic(dynamic_model));
}
void
@@ -3506,9 +3517,28 @@ DynamicModel::toStatic(StaticModel &static_model) const
static_model.AddLocalVariable(it->first, it->second->toStatic(static_model));
// Convert equations
- for (vector<BinaryOpNode *>::const_iterator it = equations.begin();
- it != equations.end(); it++)
- static_model.addEquation((*it)->toStatic(static_model));
+ int static_only_index = 0;
+ for (int i = 0; i < (int) equations.size(); i++)
+ {
+ // Detect if equation is marked [dynamic]
+ bool is_dynamic_only = false;
+ for (vector<pair<int, pair<string, string> > >::const_iterator it = equation_tags.begin();
+ it != equation_tags.end(); ++it)
+ if (it->first == i && it->second.first == "dynamic")
+ {
+ is_dynamic_only = true;
+ break;
+ }
+
+ // If yes, replace it by an equation marked [static]
+ if (is_dynamic_only)
+ {
+ static_model.addEquation(static_only_equations[static_only_index]->toStatic(static_model));
+ static_only_index++;
+ }
+ else
+ static_model.addEquation(equations[i]->toStatic(static_model));
+ }
// Convert auxiliary equations
for (deque<BinaryOpNode *>::const_iterator it = aux_equations.begin();
@@ -4150,3 +4180,32 @@ DynamicModel::isModelLocalVariableUsed() const
}
return used_local_vars.size() > 0;
}
+
+void
+DynamicModel::addStaticOnlyEquation(expr_t eq)
+{
+ BinaryOpNode *beq = dynamic_cast<BinaryOpNode *>(eq);
+ assert(beq != NULL && beq->get_op_code() == oEqual);
+
+ static_only_equations.push_back(beq);
+}
+
+size_t
+DynamicModel::staticOnlyEquationsNbr() const
+{
+ return static_only_equations.size();
+}
+
+size_t
+DynamicModel::dynamicOnlyEquationsNbr() const
+{
+ set<int> eqs;
+
+ for (vector<pair<int, pair<string, string> > >::const_iterator it = equation_tags.begin();
+ it != equation_tags.end(); ++it)
+ if (it->second.first == "dynamic")
+ eqs.insert(it->first);
+
+ return eqs.size();
+}
+
diff --git a/preprocessor/DynamicModel.hh b/preprocessor/DynamicModel.hh
index 0b73ff0b4c5d32fd8129eac5641c694a42fbfcd3..57bd30a79374fed0b3d5e8885b53da4f3d6e432b 100644
--- a/preprocessor/DynamicModel.hh
+++ b/preprocessor/DynamicModel.hh
@@ -31,6 +31,10 @@ using namespace std;
class DynamicModel : public ModelTree
{
private:
+ //! Stores equations declared as [static]
+ /*! They will be used in toStatic() to replace equations marked as [dynamic] */
+ vector<BinaryOpNode *> static_only_equations;
+
typedef map<pair<int, int>, int> deriv_id_table_t;
//! Maps a pair (symbol_id, lag) to a deriv ID
deriv_id_table_t deriv_id_table;
@@ -224,6 +228,15 @@ public:
//! Replaces the model equations in dynamic_model with those in this model
void replaceMyEquations(DynamicModel &dynamic_model) const;
+ //! Adds an equation marked as [static]
+ void addStaticOnlyEquation(expr_t eq);
+
+ //! Returns number of static only equations
+ size_t staticOnlyEquationsNbr() const;
+
+ //! Returns number of dynamic only equations
+ size_t dynamicOnlyEquationsNbr() const;
+
//! Writes LaTeX file with the equations of the dynamic model
void writeLatexFile(const string &basename) const;
diff --git a/preprocessor/DynareBison.yy b/preprocessor/DynareBison.yy
index fb62d7272c64ae80a025fb1926a18ea2afa00798..05a766efb9628672942b9e68976195b997abc826 100644
--- a/preprocessor/DynareBison.yy
+++ b/preprocessor/DynareBison.yy
@@ -603,6 +603,8 @@ tags_list : tags_list COMMA tag_pair
tag_pair : NAME EQUAL QUOTED_STRING
{ driver.add_equation_tags($1, $3); }
+ | NAME
+ { driver.add_equation_tags($1, new string()); }
;
hand_side : '(' hand_side ')'
diff --git a/preprocessor/ExprNode.cc b/preprocessor/ExprNode.cc
index 09368ecb91b0d1ee419577d4a2b760445f56ccb1..abe47a2ad7be7011f41dcb415080e9e5c5908c3c 100644
--- a/preprocessor/ExprNode.cc
+++ b/preprocessor/ExprNode.cc
@@ -469,6 +469,12 @@ NumConstNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
return const_cast<NumConstNode *>(this);
}
+bool
+NumConstNode::isInStaticForm() const
+{
+ return true;
+}
+
VariableNode::VariableNode(DataTree &datatree_arg, int symb_id_arg, int lag_arg) :
ExprNode(datatree_arg),
symb_id(symb_id_arg),
@@ -1314,6 +1320,12 @@ VariableNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
}
}
+bool
+VariableNode::isInStaticForm() const
+{
+ return lag == 0;
+}
+
UnaryOpNode::UnaryOpNode(DataTree &datatree_arg, UnaryOpcode op_code_arg, const expr_t arg_arg, int expectation_information_set_arg, int param1_symb_id_arg, int param2_symb_id_arg) :
ExprNode(datatree_arg),
arg(arg_arg),
@@ -2322,6 +2334,17 @@ UnaryOpNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
return buildSimilarUnaryOpNode(argsubst, datatree);
}
+bool
+UnaryOpNode::isInStaticForm() const
+{
+ if (op_code == oSteadyState || op_code == oSteadyStateParamDeriv
+ || op_code == oSteadyStateParam2ndDeriv
+ || op_code == oExpectation)
+ return false;
+ else
+ return arg->isInStaticForm();
+}
+
BinaryOpNode::BinaryOpNode(DataTree &datatree_arg, const expr_t arg1_arg,
BinaryOpcode op_code_arg, const expr_t arg2_arg) :
ExprNode(datatree_arg),
@@ -3558,6 +3581,12 @@ BinaryOpNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
}
+bool
+BinaryOpNode::isInStaticForm() const
+{
+ return arg1->isInStaticForm() && arg2->isInStaticForm();
+}
+
TrinaryOpNode::TrinaryOpNode(DataTree &datatree_arg, const expr_t arg1_arg,
TrinaryOpcode op_code_arg, const expr_t arg2_arg, const expr_t arg3_arg) :
ExprNode(datatree_arg),
@@ -4155,6 +4184,12 @@ TrinaryOpNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const
return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
}
+bool
+TrinaryOpNode::isInStaticForm() const
+{
+ return arg1->isInStaticForm() && arg2->isInStaticForm() && arg3->isInStaticForm();
+}
+
ExternalFunctionNode::ExternalFunctionNode(DataTree &datatree_arg,
int symb_id_arg,
const vector<expr_t> &arguments_arg) :
@@ -4732,6 +4767,16 @@ ExternalFunctionNode::removeTrendLeadLag(map<int, expr_t> trend_symbols_map) con
return buildSimilarExternalFunctionNode(arguments_subst, datatree);
}
+bool
+ExternalFunctionNode::isInStaticForm() const
+{
+ for (vector<expr_t>::const_iterator it = arguments.begin(); it != arguments.end(); ++it)
+ if (!(*it)->isInStaticForm())
+ return false;
+
+ return true;
+}
+
FirstDerivExternalFunctionNode::FirstDerivExternalFunctionNode(DataTree &datatree_arg,
int top_level_symb_id_arg,
const vector<expr_t> &arguments_arg,
diff --git a/preprocessor/ExprNode.hh b/preprocessor/ExprNode.hh
index 992ed9e49445b2b7e75f734230caae5b6fc2923c..526a26b168e289a5f9a8dfc12ce758ea7ed4b6cc 100644
--- a/preprocessor/ExprNode.hh
+++ b/preprocessor/ExprNode.hh
@@ -1,5 +1,5 @@
/*
- * Copyright (C) 2007-2012 Dynare Team
+ * Copyright (C) 2007-2013 Dynare Team
*
* This file is part of Dynare.
*
@@ -393,6 +393,9 @@ public:
//! Move a trend variable with lag/lead to time t by dividing/multiplying by its growth factor
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const = 0;
+
+ //! Returns true if the expression is in static form (no lead, no lag, no expectation, no STEADY_STATE)
+ virtual bool isInStaticForm() const = 0;
};
//! Object used to compare two nodes (using their indexes)
@@ -448,6 +451,7 @@ public:
virtual expr_t detrend(int symb_id, bool log_trend, expr_t trend) const;
virtual expr_t cloneDynamic(DataTree &dynamic_datatree) const;
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const;
+ virtual bool isInStaticForm() const;
};
//! Symbol or variable node
@@ -508,6 +512,7 @@ public:
virtual expr_t detrend(int symb_id, bool log_trend, expr_t trend) const;
virtual expr_t cloneDynamic(DataTree &dynamic_datatree) const;
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const;
+ virtual bool isInStaticForm() const;
};
//! Unary operator node
@@ -583,6 +588,7 @@ public:
virtual expr_t detrend(int symb_id, bool log_trend, expr_t trend) const;
virtual expr_t cloneDynamic(DataTree &dynamic_datatree) const;
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const;
+ virtual bool isInStaticForm() const;
};
//! Binary operator node
@@ -677,6 +683,7 @@ public:
expr_t addMultipliersToConstraints(int i);
//! Returns the non-zero hand-side of an equation (that must have a hand side equal to zero)
expr_t getNonZeroPartofEquation() const;
+ virtual bool isInStaticForm() const;
};
//! Trinary operator node
@@ -739,6 +746,7 @@ public:
virtual expr_t detrend(int symb_id, bool log_trend, expr_t trend) const;
virtual expr_t cloneDynamic(DataTree &dynamic_datatree) const;
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const;
+ virtual bool isInStaticForm() const;
};
//! External function node
@@ -812,6 +820,7 @@ public:
virtual expr_t detrend(int symb_id, bool log_trend, expr_t trend) const;
virtual expr_t cloneDynamic(DataTree &dynamic_datatree) const;
virtual expr_t removeTrendLeadLag(map<int, expr_t> trend_symbols_map) const;
+ virtual bool isInStaticForm() const;
};
class FirstDerivExternalFunctionNode : public ExternalFunctionNode
diff --git a/preprocessor/ModFile.cc b/preprocessor/ModFile.cc
index d6b7e2786ccc70fb292151feb2345ab0e7de119f..3c4933e40d952c4853225e572c6e76945bdc1eca 100644
--- a/preprocessor/ModFile.cc
+++ b/preprocessor/ModFile.cc
@@ -236,6 +236,19 @@ ModFile::checkPass()
<< "(deprecated) or the dsge_var option must be passed to the estimation statement (preferred)." << endl;
exit(EXIT_FAILURE);
}
+
+ if (dynamic_model.staticOnlyEquationsNbr() != dynamic_model.dynamicOnlyEquationsNbr())
+ {
+ cerr << "ERROR: the number of equations marked [static] must be equal to the number of equations marked [dynamic]" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ if (dynamic_model.staticOnlyEquationsNbr() > 0 &&
+ (mod_file_struct.ramsey_policy_present || mod_file_struct.discretionary_policy_present))
+ {
+ cerr << "ERROR: marking equations as [static] or [dynamic] is not possible with ramsey_policy or discretionary_policy" << endl;
+ exit(EXIT_FAILURE);
+ }
}
void
diff --git a/preprocessor/ModelTree.cc b/preprocessor/ModelTree.cc
index 33584f86892046a1b1ed3d105ebc1550cce99bb7..b6f3fd3acd918c1133bacca0518ed1528c7a3c9d 100644
--- a/preprocessor/ModelTree.cc
+++ b/preprocessor/ModelTree.cc
@@ -1385,9 +1385,12 @@ ModelTree::addEquation(expr_t eq)
}
void
-ModelTree::addEquationTags(int i, const string &key, const string &value)
+ModelTree::addEquation(expr_t eq, vector<pair<string, string> > &eq_tags)
{
- equation_tags.push_back(make_pair(i, make_pair(key, value)));
+ int n = equation_number();
+ for (size_t i = 0; i < eq_tags.size(); i++)
+ equation_tags.push_back(make_pair(n, eq_tags[i]));
+ addEquation(eq);
}
void
diff --git a/preprocessor/ModelTree.hh b/preprocessor/ModelTree.hh
index e8a7727fa4182fb2db89221fd629b55459332b82..4800ccc92500e1f924df662d4e3d78aacebf7917 100644
--- a/preprocessor/ModelTree.hh
+++ b/preprocessor/ModelTree.hh
@@ -297,8 +297,8 @@ public:
int mfs;
//! Declare a node as an equation of the model
void addEquation(expr_t eq);
- //! Adds tags to equation number i
- void addEquationTags(int i, const string &key, const string &value);
+ //! Declare a node as an equation of the model, also giving its tags
+ void addEquation(expr_t eq, vector<pair<string, string> > &eq_tags);
//! Declare a node as an auxiliary equation of the model, adding it at the end of the list of auxiliary equations
void addAuxEquation(expr_t eq);
//! Returns the number of equations in the model
diff --git a/preprocessor/ParsingDriver.cc b/preprocessor/ParsingDriver.cc
index 5d9ed99bb8de5c23a6b59af0040da8aa37c5af0c..ca7c515b4da7536bbfc176e345c7ac8627679700 100644
--- a/preprocessor/ParsingDriver.cc
+++ b/preprocessor/ParsingDriver.cc
@@ -249,8 +249,7 @@ ParsingDriver::add_predetermined_variable(string *name)
void
ParsingDriver::add_equation_tags(string *key, string *value)
{
- int n = model_tree->equation_number();
- model_tree->addEquationTags(n, *key, *value);
+ eq_tags.push_back(make_pair(*key, *value));
delete key;
delete value;
}
@@ -1933,7 +1932,28 @@ expr_t
ParsingDriver::add_model_equal(expr_t arg1, expr_t arg2)
{
expr_t id = model_tree->AddEqual(arg1, arg2);
- model_tree->addEquation(id);
+
+ // Detect if the equation is tagged [static]
+ bool is_static_only = false;
+ for (vector<pair<string, string> >::const_iterator it = eq_tags.begin();
+ it != eq_tags.end(); ++it)
+ if (it->first == "static")
+ {
+ is_static_only = true;
+ break;
+ }
+
+ if (is_static_only)
+ {
+ if (!id->isInStaticForm())
+ error("An equation tagged [static] cannot contain leads, lags, expectations or STEADY_STATE operators");
+
+ dynamic_model->addStaticOnlyEquation(id);
+ }
+ else
+ model_tree->addEquation(id, eq_tags);
+
+ eq_tags.clear();
return id;
}
diff --git a/preprocessor/ParsingDriver.hh b/preprocessor/ParsingDriver.hh
index 3ad6439ebda782f5e809ef5dee0b026995868198..279ba13d95d8601f8d8b893b3fd7758e19c01347 100644
--- a/preprocessor/ParsingDriver.hh
+++ b/preprocessor/ParsingDriver.hh
@@ -201,6 +201,8 @@ private:
expr_t add_model_variable(int symb_id, int lag);
//! For parsing the graph_format option
SymbolList graph_formats;
+ //! Temporary storage for equation tags
+ vector<pair<string, string> > eq_tags;
//! The mod file representation constructed by this ParsingDriver
ModFile *mod_file;
diff --git a/tests/Makefile.am b/tests/Makefile.am
index d9ea86cc35750eb3ba3b58cbad80f35b9a9e2df4..2cb8bd263d9ab60c405cd0b25dd1482c322b6a42 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -6,6 +6,7 @@ MODFILES = \
gsa/ls2003.mod \
ramst.mod \
ramst_a.mod \
+ ramst_static_tag.mod \
example1.mod \
example2.mod \
example1_use_dll.mod \
diff --git a/tests/ramst_static_tag.mod b/tests/ramst_static_tag.mod
new file mode 100644
index 0000000000000000000000000000000000000000..6b69f8a2e1ec9451e54df5ba8d7c556cdb4b035e
--- /dev/null
+++ b/tests/ramst_static_tag.mod
@@ -0,0 +1,37 @@
+var c k;
+varexo x;
+
+parameters alph gam delt bet aa;
+alph=0.5;
+gam=0.5;
+delt=0.02;
+bet=0.05;
+aa=0.5;
+
+
+model;
+c + k - aa*x*k(-1)^alph - (1-delt)*k(-1);
+[dynamic] c^(-gam) - (1+bet)^(-1)*(aa*alph*x(+1)*k^(alph-1) + 1 - delt)*c(+1)^(-gam);
+[static] k = ((delt+bet)/(x*aa*alph))^(1/(alph-1));
+end;
+
+initval;
+x = 1;
+k = ((delt+bet)/(1.0*aa*alph))^(1/(alph-1));
+c = aa*k^alph-delt*k;
+end;
+
+steady;
+
+check;
+
+shocks;
+var x;
+periods 1;
+values 1.2;
+end;
+
+simul(periods=200);
+
+rplot c;
+rplot k;