Skip to content

Commit

Permalink
BUG: r: Fix bug in the demo code.
Browse files Browse the repository at this point in the history 
In the demo script, ensure that parameter names are available
to be used in initial condition expressions.

[test r]
  • Loading branch information
@WarrenWeckesser
WarrenWeckesser committed Sep 30, 2024
1 parent 14ef357 commit cda1617
Showing 1 changed file with 27 additions and 18 deletions.
45 changes: 27 additions & 18 deletions src/vf_r.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,8 @@ using namespace std;
using namespace GiNaC;


static void vf_header_comment(ofstream& fout, string name, string indvar)
static void
vf_header_comment(ofstream& fout, string name, string indvar)
{
fout << "#" << endl;
fout << "# " << name << "(" << indvar << ", state, parameters)" << endl;
Expand All @@ -45,7 +46,8 @@ static void vf_header_comment(ofstream& fout, string name, string indvar)
return;
}

static void jac_header_comment(ofstream& fout, string name, string indvar)
static void
jac_header_comment(ofstream& fout, string name, string indvar)
{
fout << "#" << endl;
fout << "# " << name << "_jac(" << indvar << ", state, parameters)" << endl;
Expand Down Expand Up @@ -222,21 +224,24 @@ void VectorField::PrintRode(map<string, string> options)
if (HasPi) {
fout << "Pi = pi\n";
}
AssignNameValueLists(fout, " ", conname_list, "<-", convalue_list, "");
AssignNameValueLists(fout, "", conname_list, "<-", convalue_list, "");
fout << endl;
fout << "# Default parameter values" << endl;
AssignNameValueLists(fout, "", parname_list, "<-", pardefval_list, "");
fout << endl;
if (np > 0) {
fout << "# --- Parameters ---\n";
fout << "# Parameters list\n";
fout << "parameters = c(\n";
for (int i = 0; i < np; ++i) {
fout << " " << parname_list[i] << " = " << pardefval_list[i];
fout << " " << parname_list[i] << " = " << parname_list[i];
if (i < np - 1) {
fout << ",";
}
fout << "\n";
}
fout << ")\n\n";
}
fout << "# --- Initial conditions ---\n";
fout << "# Initial conditions\n";
fout << "state = c(\n";
for (int i = 0; i < nv; ++i) {
fout << " " << varname_list[i] << " = " << vardefic_list[i];
Expand All @@ -245,17 +250,17 @@ void VectorField::PrintRode(map<string, string> options)
}
fout << "\n";
}
fout << ")\n\n";
fout << ")\n";
fout << endl;
fout << "# --- Time values ---\n";
fout << "# Time values\n";
fout << "times = seq(0, 10, by = 0.02)\n";
fout << endl;
fout << "# --- Call the ODE solver ---" << endl;
fout << "# Call the ODE solver" << endl;
fout << "sol = ode(y = state, times = times, func = " << Name() << ", parms = parameters,\n";
fout << " jactype = \"fullusr\", jacfunc = " << Name() << "_jac,\n";
fout << " atol = 1e-8, rtol = 1e-6)\n";
fout << endl;
fout << "# --- Plot the solution ---" << endl;
fout << "# Plot the solution" << endl;
fout << "par(mfcol = c(" << nv << ", 1))" << endl;
fout << "t <- sol[, \"time\"]" << endl;
for (int i = 0; i < nv; ++i) {
Expand Down Expand Up @@ -509,21 +514,24 @@ void VectorField::PrintRdede(map<string, string> options)
if (HasPi) {
fout << "Pi = pi\n";
}
AssignNameValueLists(fout, " ", conname_list, "<-", convalue_list, "");
AssignNameValueLists(fout, "", conname_list, "<-", convalue_list, "");
fout << endl;
fout << "# Default parameter values" << endl;
AssignNameValueLists(fout, "", parname_list, "<-", pardefval_list, "");
fout << endl;
if (np > 0) {
fout << "# --- Parameters ---\n";
fout << "# Parameters list\n";
fout << "parameters = c(\n";
for (int i = 0; i < np; ++i) {
fout << " " << parname_list[i] << " = " << pardefval_list[i];
fout << " " << parname_list[i] << " = " << parname_list[i];
if (i < np - 1) {
fout << ",";
}
fout << "\n";
}
fout << ")\n\n";
}
fout << "# --- Initial conditions ---\n";
fout << "# Initial conditions\n";
fout << "state = c(\n";
for (int i = 0; i < nv; ++i) {
fout << " " << varname_list[i] << " = " << vardefic_list[i];
Expand All @@ -534,13 +542,14 @@ void VectorField::PrintRdede(map<string, string> options)
}
fout << ")\n";
fout << endl;
fout << "# --- Time values ---\n";
fout << "# Time values\n";
fout << "times = seq(0, 10, by = 0.02)\n";
fout << endl;
fout << "# --- Call the DDE solver ---" << endl;
fout << "sol = dede(y = state, times = times, func = " << Name() << ", parms = parameters)\n";
fout << "# Call the DDE solver" << endl;
fout << "sol = dede(y = state, times = times, func = " << Name() << ", parms = parameters,\n";
fout << " atol = 1e-12, rtol = 1e-8)\n";
fout << endl;
fout << "# --- Plot the solution ---" << endl;
fout << "# Plot the solution" << endl;
fout << "par(mfcol = c(" << nv << ", 1))" << endl;
fout << "t <- sol[, \"time\"]" << endl;
for (int i = 0; i < nv; ++i) {
Expand Down

0 comments on commit cda1617

Please sign in to comment.