From b1e4884237376fdce06265b7f96fb63e0bd8b224 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Fri, 7 Apr 2023 15:45:44 +0200
Subject: [PATCH] =?UTF-8?q?=F0=9F=90=9B=20sign(0)=20was=20returning=201=20?=
=?UTF-8?q?instead=20of=200=20with=20use=5Fdll?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
The C99 copysign() function was used in the generated C output, but that
function does not correctly handle zero. Replace it by a custom sign()
function.
---
src/DataTree.cc | 12 ++++++++++--
src/DataTree.hh | 8 ++++----
src/ExprNode.cc | 10 ++++------
src/ModelTree.hh | 38 +++++++++++++++++---------------------
src/StaticModel.cc | 2 +-
5 files changed, 36 insertions(+), 34 deletions(-)
diff --git a/src/DataTree.cc b/src/DataTree.cc
index 129d1b85..96f16214 100644
--- a/src/DataTree.cc
+++ b/src/DataTree.cc
@@ -890,7 +890,7 @@ DataTree::minLagForSymbol(int symb_id) const
}
void
-DataTree::writePowerDeriv(ostream &output) const
+DataTree::writeCHelpersDefinition(ostream &output) const
{
if (isBinaryOpUsed(BinaryOpcode::powerDeriv))
output << "/*" << endl
@@ -909,13 +909,21 @@ DataTree::writePowerDeriv(ostream &output) const
<< " return dxp;" << endl
<< " }" << endl
<< "}" << endl;
+
+ if (isUnaryOpUsed(UnaryOpcode::sign))
+ output << "double sign(double x)" << endl
+ << "{" << endl
+ << " return (x > 0) ? 1 : ((x < 0) ? -1 : 0);" << endl
+ << "}" << endl;
}
void
-DataTree::writePowerDerivHeader(ostream &output) const
+DataTree::writeCHelpersDeclaration(ostream &output) const
{
if (isBinaryOpUsed(BinaryOpcode::powerDeriv))
output << "double getPowerDeriv(double x, double p, int k);" << endl;
+ if (isUnaryOpUsed(UnaryOpcode::sign))
+ output << "double sign(double x);" << endl;
}
vector<string>
diff --git a/src/DataTree.hh b/src/DataTree.hh
index bfc86938..108fb253 100644
--- a/src/DataTree.hh
+++ b/src/DataTree.hh
@@ -277,10 +277,10 @@ public:
//! Returns the minimum lag (as a negative number) of the given symbol in the whole data tree (and not only in the equations !!)
/*! Returns 0 if the symbol is not used */
int minLagForSymbol(int symb_id) const;
- //! Write getPowerDeriv in C (function body)
- void writePowerDeriv(ostream &output) const;
- //! Write getPowerDeriv in C (prototype)
- void writePowerDerivHeader(ostream &output) const;
+ //! Writes definitions of C function helpers (getPowerDeriv(), sign())
+ void writeCHelpersDefinition(ostream &output) const;
+ //! Writes declarations of C function helpers (getPowerDeriv(), sign())
+ void writeCHelpersDeclaration(ostream &output) const;
//! Thrown when trying to access an unknown variable by deriv_id
class UnknownDerivIDException
{
diff --git a/src/ExprNode.cc b/src/ExprNode.cc
index 12527dab..e69a177b 100644
--- a/src/ExprNode.cc
+++ b/src/ExprNode.cc
@@ -2950,10 +2950,10 @@ UnaryOpNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
output << "abs";
break;
case UnaryOpcode::sign:
- if (isCOutput(output_type))
- output << "copysign";
- else
- output << "sign";
+ /* C does not have a sign() function, and copysign() is not suitable
+ because it does not handle zero correctly, so we define our own sign()
+ helper function, see DataTree::writeCHelpersDefinition() */
+ output << "sign";
break;
case UnaryOpcode::steadyState:
ExprNodeOutputType new_output_type;
@@ -3054,8 +3054,6 @@ UnaryOpNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
&& arg->precedence(output_type, temporary_terms) < precedence(output_type, temporary_terms)))
{
output << LEFT_PAR(output_type);
- if (op_code == UnaryOpcode::sign && isCOutput(output_type))
- output << "1.0,";
close_parenthesis = true;
}
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index ddf72dc1..f1ab64ad 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -973,7 +973,7 @@ ModelTree::writeModelCFile(const string &basename, const string &mexext,
output << "#include <math.h>" << endl
<< R"(#include "mex.h")" << endl // Needed for calls to external functions
<< endl;
- writePowerDerivHeader(output);
+ writeCHelpersDeclaration(output);
output << endl
<< prototype_tt << endl
<< "{" << endl
@@ -1012,7 +1012,7 @@ ModelTree::writeModelCFile(const string &basename, const string &mexext,
output << "#include <math.h>" << endl
<< R"(#include "mex.h")" << endl // Needed for calls to external functions
<< endl;
- writePowerDerivHeader(output);
+ writeCHelpersDeclaration(output);
output << endl
<< prototype_main << endl
<< "{" << endl
@@ -1043,8 +1043,7 @@ ModelTree::writeModelCFile(const string &basename, const string &mexext,
output << "#include " << it.filename() << endl;
output << endl;
- // Write function definition if BinaryOpcode::powerDeriv is used
- writePowerDeriv(output);
+ writeCHelpersDefinition(output);
output << endl
<< "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
@@ -2645,18 +2644,13 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
NB: The prefix (static/dynamic) is added to the filename (even though it’s
the same source between static and dynamic) to avoid a race condition when
static and dynamic are compiled in parallel. */
- open_file(model_src_dir / (prefix + "power_deriv.h"));
- writePowerDerivHeader(output);
- output.close();
-
- filesystem::path power_deriv_src {model_src_dir / (prefix + "power_deriv.c")};
- open_file(power_deriv_src);
+ filesystem::path helpers_src {model_src_dir / (prefix + "helpers.c")};
+ open_file(helpers_src);
output << "#include <math.h>" << endl << endl;
- writePowerDeriv(output);
+ writeCHelpersDefinition(output);
output.close();
- auto power_deriv_object {compileMEX(model_src_dir, (prefix + "power_deriv"),
- mexext, { power_deriv_src },
- matlabroot, false)};
+ auto helpers_object {compileMEX(model_src_dir, (prefix + "helpers"),
+ mexext, { helpers_src }, matlabroot, false)};
size_t ttlen {0};
@@ -2724,8 +2718,9 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
open_file(source_tt);
output << "#include <math.h>" << endl
<< R"(#include "mex.h")" << endl // Needed for calls to external functions
- << R"(#include ")" << prefix << R"(power_deriv.h")" << endl
- << endl
+ << endl;
+ writeCHelpersDeclaration(output);
+ output << endl
<< prototype_tt << endl
<< "{" << endl
<< tt_sparse_output[i].str()
@@ -2747,7 +2742,7 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
output << "#include <math.h>" << endl
<< R"(#include "mex.h")" << endl // Needed for calls to external functions
<< endl;
- writePowerDerivHeader(output);
+ writeCHelpersDeclaration(output);
output << endl
<< prototype_main << endl
<< "{" << endl
@@ -2836,7 +2831,7 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
<< "}" << endl;
output.close();
- vector<filesystem::path> mex_input_files { power_deriv_object, main_object_file, source_mex };
+ vector<filesystem::path> mex_input_files { helpers_object, main_object_file, source_mex };
for (int j {0}; j <= i; j++)
mex_input_files.push_back(tt_object_files[j]);
compileMEX(mex_dir, funcname, mexext, mex_input_files, matlabroot);
@@ -2864,8 +2859,9 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
open_file(source_mex);
output << "#include <math.h>" << endl
<< R"(#include "mex.h")" << endl
- << R"(#include "../)" << prefix << R"(power_deriv.h")" << endl
- << endl
+ << endl;
+ writeCHelpersDeclaration(output);
+ output << endl
<< "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])" << endl
<< "{" << endl
<< " if (nrhs != " << nargin << ")" << endl
@@ -2928,7 +2924,7 @@ ModelTree::writeSparseModelCFiles(const string &basename, const string &mexext,
}
output << "}" << endl;
output.close();
- compileMEX(block_dir, funcname, mexext, { source_mex, power_deriv_object }, matlabroot);
+ compileMEX(block_dir, funcname, mexext, { source_mex, helpers_object }, matlabroot);
}
}
}
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index 7ff90b85..a923d322 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -920,7 +920,7 @@ StaticModel::writeRamseyMultipliersDerivativesCFile(const string &basename, cons
output << "#include <math.h>" << endl << endl
<< R"(#include "mex.h")" << endl // Needed for calls to external functions
<< endl;
- writePowerDeriv(output);
+ writeCHelpersDefinition(output);
output << endl
<< "void ramsey_multipliers_static_g1(const double *restrict y, const double *restrict x, const double *restrict params, double *restrict T, double *restrict g1m_v)" << endl
<< "{" << endl;
--
GitLab